原文鏈接：http://tecdat.cn/?p=11387

盡管線性模型是最簡(jiǎn)單的機(jī)器學(xué)習(xí)技術(shù)之一，但它們?nèi)匀皇沁M(jìn)行預(yù)測(cè)的強(qiáng)大工具。這尤其是由于線性模型特別容易解釋這一事實(shí)。在這里，我將討論使用空氣質(zhì)量數(shù)據(jù)集的普通最小二乘回歸示例解釋線性模型時(shí)最重要的方面。

空氣質(zhì)量數(shù)據(jù)集

空氣質(zhì)量數(shù)據(jù)集包含對(duì)在紐約獲得的以下四個(gè)空氣質(zhì)量指標(biāo)的154次測(cè)量：

• 臭氧：平均臭氧水平，以十億分之一為單位

• Solar.R：太陽(yáng)輻射?

• 風(fēng)：平均風(fēng)速，每小時(shí)英里

• 溫度：每日最高溫度，以華氏度為單位

我們將通過(guò)刪除所有NA?并排除??Month?和Day?列來(lái)清理數(shù)據(jù)集? ，這些列和? 列不應(yīng)充當(dāng)預(yù)測(cè)變量。

1. data(airquality)

2. ozone <- subset(na.omit(airquality),

3. select = c("Ozone", "Solar.R", "Wind", "Temp"))

數(shù)據(jù)探索和準(zhǔn)備

預(yù)測(cè)任務(wù)如下：根據(jù)太陽(yáng)輻射，風(fēng)速和溫度，我們可以預(yù)測(cè)臭氧水平嗎？要查看線性模型的假設(shè)是否適合手頭的數(shù)據(jù)，我們將計(jì)算變量之間的相關(guān)性：

1. # 散點(diǎn)圖矩陣

2. plot(ozone)

?

?

1. # 成對(duì)變量的相關(guān)性

2. cors <- cor(ozone)

3. print(cors)

1. ## ? ? ? ? ? ? ?Ozone ? ?Solar.R ? ? ? Wind ? ? ? Temp

2. ## Ozone ? ?1.0000000 ?0.3483417 -0.6124966 ?0.6985414

3. ## Solar.R ?0.3483417 ?1.0000000 -0.1271835 ?0.2940876

4. ## Wind ? ?-0.6124966 -0.1271835 ?1.0000000 -0.4971897

5. ## Temp ? ? 0.6985414 ?0.2940876 -0.4971897 ?1.0000000

1. # 哪些變量高度相關(guān)，排除自相關(guān)

2. 
   

3. print(cor.names)

1. ## [1] "Wind+Ozone: -0.61" "Temp+Ozone: 0.7" ? "Ozone+Wind: -0.61"

2. ## [4] "Ozone+Temp: 0.7"

由于臭氧參與兩個(gè)線性相互作用，即：

• 臭氧與溫度呈正相關(guān)

• 臭氧與風(fēng)負(fù)相關(guān)

這表明應(yīng)該有可能使用其余特征來(lái)形成預(yù)測(cè)臭氧水平的線性模型。

分為訓(xùn)練和測(cè)試集

我們將抽取70％的樣本進(jìn)行訓(xùn)練，并抽取30％的樣本進(jìn)行測(cè)試：

1. set.seed(123)

2. N.train <- ceiling(0.7 * nrow(ozone))

3. N.test <- nrow(ozone) - N.train

4. trainset <- sample(seq_len(nrow(ozone)), N.train)

5. testset <- setdiff(seq_len(nrow(ozone)), trainset)

研究線性模型

為了說(shuō)明解釋線性模型的最重要方面，我們將通過(guò)以下方式訓(xùn)練訓(xùn)練數(shù)據(jù)的普通最小二乘模型：

為了解釋模型，我們使用以下??summary?函數(shù)：

1. model.summary <- summary(model)

2. print(model.summary)

1. ##

2. ## Call:

4. ## Residuals:

5. ## ? ? Min ? ? ?1Q ?Median ? ? ?3Q ? ? Max

6. ## -36.135 -12.670 ?-2.221 ? 9.420 ?65.914

7. ##

8. ## Coefficients:

9. ## ? ? ? ? ? ? ?Estimate Std. Error t value Pr(>|t|)

10. ## (Intercept) -65.76604 ? 22.52381 ?-2.920 0.004638 **

11. ## Solar.R ? ? ? 0.05309 ? ?0.02305 ? 2.303 0.024099 *

12. ## Temp ? ? ? ? ?1.56320 ? ?0.25530 ? 6.123 4.03e-08 ***

13. ## Wind ? ? ? ? -2.61904 ? ?0.68921 ?-3.800 0.000295 ***

14. ## ---

15. ## Signif. codes: ?0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

16. ##

17. ## Residual standard error: 18.17 on 74 degrees of freedom

18. ## Multiple R-squared: ?0.5924, Adjusted R-squared: ?0.5759

19. ## F-statistic: 35.85 on 3 and 74 DF, ?p-value: 2.039e-14

殘差

我們獲得的第一條信息是殘差。?

殘留中值表明，該模型通常預(yù)測(cè)的臭氧值略高于觀測(cè)值。但是，最大值很大，表明某些離群值預(yù)測(cè)也太低了。查看數(shù)字可能有點(diǎn)抽象，因此讓我們根據(jù)觀察值繪制模型的預(yù)測(cè)：

1. res <- residuals(model)

3. # 向圖中添加殘差

4. segments(obs, pred, obs, pred + res)

?系數(shù)

現(xiàn)在我們了解了殘差，讓我們看一下系數(shù)。我們可以使用該??coefficients?函數(shù)來(lái)獲取模型的擬合系數(shù)：

1. ## ?(Intercept) ? ? ?Solar.R ? ? ? ? Temp ? ? ? ? Wind

2. ## -65.76603538 ? 0.05308965 ? 1.56320267 ?-2.61904128

請(qǐng)注意，模型的截距值非常低。這是在所有獨(dú)立值均為零的情況下模型將預(yù)測(cè)的值。 低系數(shù)??Solar.R?表示太陽(yáng)輻射對(duì)預(yù)測(cè)臭氧水平?jīng)]有重要作用，這不足為奇，因?yàn)樵谖覀兊奶剿餍苑治鲋?，它與臭氧水平?jīng)]有很大的相關(guān)性。 系數(shù)??Temp?表示溫度高時(shí)臭氧水平高（因?yàn)槌粞鯐?huì)更快形成）。 系數(shù)??Wind?告訴我們快風(fēng)時(shí)臭氧水平會(huì)降低（因?yàn)槌粞鯐?huì)被吹走）。

與系數(shù)關(guān)聯(lián)的其他值提供有關(guān)估計(jì)的統(tǒng)計(jì)確定性的信息。

1. ## ? ? ? ? ? ? ? ? Estimate ?Std. Error ? t value ? ? Pr(>|t|)

2. ## (Intercept) -65.76603538 22.52380940 -2.919845 4.638426e-03

3. ## Solar.R ? ? ? 0.05308965 ?0.02305379 ?2.302860 2.409936e-02

4. ## Temp ? ? ? ? ?1.56320267 ?0.25530453 ?6.122894 4.034064e-08

5. ## Wind ? ? ? ? -2.61904128 ?0.68920661 -3.800081 2.946349e-04

• Std. Error?是系數(shù)估計(jì)的標(biāo)準(zhǔn)誤差

• t value?以標(biāo)準(zhǔn)誤差表示系數(shù)的值

• Pr(>|t|)?是t檢驗(yàn)的p值，表示檢驗(yàn)統(tǒng)計(jì)量的重要性

標(biāo)準(zhǔn)誤差

系數(shù)的標(biāo)準(zhǔn)誤差定義為特征方差的標(biāo)準(zhǔn)偏差：

?

?

在R中，可以通過(guò)以下方式計(jì)算模型估計(jì)的標(biāo)準(zhǔn)誤差：

1. ## ? ? ? ? ? ? ?(Intercept) ? ? ? Solar.R ? ? ? ? Temp ? ? ? ? ?Wind

2. ## (Intercept) 507.32198977 ?2.893612e-02 -5.345957524 -9.940961e+00

3. ## Solar.R ? ? ? 0.02893612 ?5.314773e-04 -0.001667748 ?7.495211e-05

4. ## Temp ? ? ? ? -5.34595752 -1.667748e-03 ?0.065180401 ?6.715467e-02

5. ## Wind ? ? ? ? -9.94096142 ?7.495211e-05 ?0.067154670 ?4.750058e-01

1. ## (Intercept) ? ? Solar.R ? ? ? ?Temp ? ? ? ?Wind

2. ## 22.52380940 ?0.02305379 ?0.25530453 ?0.68920661

現(xiàn)在， 你可能想知道這些值的??vcov?來(lái)源。它定義為設(shè)計(jì)矩陣的方差-協(xié)方差矩陣，該矩陣按誤差的方差標(biāo)準(zhǔn)化：

1. ## ? ? ? ? ? ? ?(Intercept) ? ? ? Solar.R ? ? ? ? Temp ? ? ? ? ?Wind

2. ## (Intercept) 507.32198977 ?2.893612e-02 -5.345957524 -9.940961e+00

3. ## Solar.R ? ? ? 0.02893612 ?5.314773e-04 -0.001667748 ?7.495211e-05

4. ## Temp ? ? ? ? -5.34595752 -1.667748e-03 ?0.065180401 ?6.715467e-02

5. ## Wind ? ? ? ? -9.94096142 ?7.495211e-05 ?0.067154670 ?4.750058e-01

用于標(biāo)準(zhǔn)化的方差-協(xié)方差矩陣的方差是誤差的估計(jì)方差，其定義為

?

? ?cov.unscaled?參數(shù)是簡(jiǎn)單地所有的方差-協(xié)方差矩陣? ：

1. # 通過(guò)'model.matrix'將截距作為特征

2. X <- model.matrix(model) # design matrix

3. # 求解X ^ T％*％X = I求X ^ T * X的逆

4. unscaled.var.matrix <- solve(crossprod(X), diag(4))

5. print(paste("Is this the same?", isTRUE(all.equal(unscaled.var.matrix, model.summary$cov.unscaled, check.attributes = FALSE))))

## [1] "Is this the same? TRUE"

t值

t值定義為

?

在R中?

1. ## (Intercept) ? ? Solar.R ? ? ? ?Temp ? ? ? ?Wind

2. ## ? -2.919845 ? ?2.302860 ? ?6.122894 ? -3.800081

p值

在所有系數(shù)βi=0 的假設(shè)下計(jì)算p值。t值遵循t分布

model.df <- df.residual(model)

自由程度。線性模型的自由度定義為

其中n 是樣本數(shù)，p 是特征數(shù)（包括inctercept）。p值表示獲得的系數(shù)估計(jì)純粹是偶然地與零不同的可能性。因此，低p值表明變量與結(jié)果之間存在顯著關(guān)聯(lián)。

進(jìn)一步統(tǒng)計(jì)

該summary?函數(shù)提供以下附加統(tǒng)計(jì)信息? ：多個(gè)R平方，調(diào)整后的R平方和F統(tǒng)計(jì)。?

殘留標(biāo)準(zhǔn)誤差

顧名思義，殘留標(biāo)準(zhǔn)誤差是模型的平均RSS（MSE）的平方根：

## [1] 18.16979

殘留標(biāo)準(zhǔn)誤差僅表示模型的平均精度。在這種情況下，該值非常低，表明該模型具有良好的擬合度。

多個(gè)R平方

R平方的倍數(shù)表示確定系數(shù)。它定義為估計(jì)值與觀察到的結(jié)果之間的相關(guān)性的平方：

## [1] 0.5924073

與[-1,1] [-1,1]中的相關(guān)性相反，R平方在[0,1] [0,1]中。

調(diào)整后的R平方

調(diào)整后的R平方值會(huì)根據(jù)模型的復(fù)雜性來(lái)調(diào)整R平方：

?

其中nn是觀察數(shù)，pp是特征數(shù)。因此，調(diào)整后的R平方可以像這樣計(jì)算：

1. n <- length(trainset) # 樣本數(shù)

2. 


3. print(r.squared.adj)

## [1] 0.5758832

如果R平方和調(diào)整后的R平方之間存在相當(dāng)大的差異，則表明可以考慮減少特征空間。

F統(tǒng)計(jì)

F統(tǒng)計(jì)量定義為已解釋方差與無(wú)法解釋方差的比率。為了進(jìn)行回歸，F(xiàn)統(tǒng)計(jì)量始終指示兩個(gè)模型之間的差異，其中模型1（p1p1）由模型2（p2p2）的特征子集定義：

?

F統(tǒng)計(jì)量描述模型2的預(yù)測(cè)性能（就RSS而言）優(yōu)于模型1的程度。報(bào)告的默認(rèn)F統(tǒng)計(jì)量是指訓(xùn)練后的模型與僅截距模型之間的差異：

1. ##

2. ## Call:

3. ##

4. ## Coefficients:

5. ## (Intercept)

6. ## ? ? ? 36.76

?

? ?因此，測(cè)試的零假設(shè)是唯一的截距-模型的擬合和指定的模型是相等的。如果可以拒絕原假設(shè)，則意味著指定模型比原模型具有更好的擬合度。

讓我們通過(guò)手工計(jì)算得出這個(gè)想法：

1. rss <- function(model) {

2. return(sum(model$residuals^2))

3. }

4. 


5. # 比較僅截距模型和具有3個(gè)模型

6. f.statistic(null.model, model)

## [1] 35.85126

在這種情況下，F(xiàn)統(tǒng)計(jì)量具有較大的值，這表明我們訓(xùn)練的模型明顯優(yōu)于僅攔截模型。

置信區(qū)間

置信區(qū)間是解釋線性模型的有用工具。默認(rèn)情況下，??confint?計(jì)算95％置信區(qū)間（±1.96σ^±1.96σ^）：

ci <- confint(model)

1. ## ? ? ? ? ? ? ? ?(Intercept) ? ? ? ? ? ? ? ? ? ?Solar.R

2. ## "95% CI: [-110.65,-20.89]" ? ? ? "95% CI: [0.01,0.1]"

3. ## ? ? ? ? ? ? ? ? ? ? ? Temp ? ? ? ? ? ? ? ? ? ? ? Wind

4. ## ? ? ?"95% CI: [1.05,2.07]" ? ?"95% CI: [-3.99,-1.25]"

這些值表明模型對(duì)截距的估計(jì)不確定。這可能表明需要更多數(shù)據(jù)才能獲得更好的擬合度。

檢索估計(jì)值的置信度和預(yù)測(cè)間隔

通過(guò)提供自interval?變量，可以將線性模型的預(yù)測(cè)轉(zhuǎn)換為間隔? 。這些間隔給出了對(duì)預(yù)測(cè)值的置信度。間隔有兩種類型：置信間隔和預(yù)測(cè)間隔。讓我們將模型應(yīng)用于測(cè)試集，使用不同的參數(shù)作為??interval?參數(shù)，以查看兩種間隔類型之間的差異：

1. # 計(jì)算預(yù)測(cè)的置信區(qū)間（CI）

2. preds.ci <- predict(model, newdata = ozone[testset,], interval = "confidence")

1. ## ? ? ?fit ? ? ? ? ? ? ? ? lwr ? ? ? ? ? ? ? ? upr ? ? ? ? ? ? ? ?Method

2. ## [1,] "-4.42397219873667" "-13.4448767773931" "4.59693237991976" "CI"

3. ## [2,] "-22.0446990408131" "-61.0004555440645" "16.9110574624383" "PI"

置信區(qū)間是窄區(qū)間，而預(yù)測(cè) 區(qū)間是寬區(qū)間。它們的值基于level?參數(shù)指定的提供的公差/重要性水平? （默認(rèn)值：0.95）。

它們的定義略有不同。給定新的觀測(cè)值xx，配置項(xiàng)和PI定義如下

?

其中tα/ 2，dftα/ 2，df是df = 2df = 2自由度且顯著性水平為αα的t值，σerrσerr是殘差的標(biāo)準(zhǔn)誤差，σ2xσx2是獨(dú)立特征的方差， x（x）表示特征的平均值。?

最受歡迎的見(jiàn)解

1.R語(yǔ)言多元Logistic邏輯回歸 應(yīng)用案例

2.面板平滑轉(zhuǎn)移回歸(PSTR)分析案例實(shí)現(xiàn)

3.matlab中的偏最小二乘回歸（PLSR）和主成分回歸（PCR）

4.R語(yǔ)言泊松Poisson回歸模型分析案例

5.R語(yǔ)言回歸中的Hosmer-Lemeshow擬合優(yōu)度檢驗(yàn)

6.r語(yǔ)言中對(duì)LASSO回歸，Ridge嶺回歸和Elastic Net模型實(shí)現(xiàn)

7.在R語(yǔ)言中實(shí)現(xiàn)Logistic邏輯回歸

8.python用線性回歸預(yù)測(cè)股票價(jià)格

9.R語(yǔ)言如何在生存分析與Cox回歸中計(jì)算IDI，NRI指標(biāo)