本內(nèi)容為【科研私家菜】R語言機器學(xué)習(xí)與臨床預(yù)測模型系列課程
你想要的R語言學(xué)習(xí)資料都在這里， 快來收藏關(guān)注【科研私家菜】
我們要建立的模型形式是Y = B0 + B1x1 + … + Bnxn + e，預(yù)測變量（特征）可以有1 ~ n個。
01 多變量線性回歸
檢查相關(guān)性方面的統(tǒng)計量，并繪出散點圖矩陣。
相關(guān)系數(shù)又稱Pearson’s r，可以用來測量兩個變量之間線性相關(guān)性的強度和方向肃拜。相關(guān)系數(shù)是一個數(shù)值，取值在-1和1之間雌团，-1表示完全負相關(guān)爆班，+1表示完全正相關(guān)。要計算相關(guān)系數(shù)辱姨，
需要使用兩個變量的協(xié)方差除以它們的標準差的乘積柿菩。如果取相關(guān)系數(shù)的平方，就可以得到R方雨涛。
data(water)
dim(water)
str(water)
head(water)
socal.water <- water[ ,-1]
head(socal.water)
library(corrplot)
water.cor <- cor(socal.water)
water.cor
corrplot(water.cor, method = "ellipse")
pairs(~ ., data = socal.water)
效果如下：
響應(yīng)變量與那些OP開頭的特征高度正相關(guān)枢舶，與OPBPC的相關(guān)系數(shù)是0.8857，與OPRC的相關(guān)系數(shù)是0.9196替久，與OPSLAKE的相關(guān)系數(shù)是0.9384凉泄。還可以看出，AP開頭的特征彼此之間高度相關(guān)蚯根，OP開頭的也是如此后众。我們會遇到多重共線性的問題。
02 特征選擇
前向逐步選擇
從一個零特征模型開始颅拦，然后每次添加一個特征蒂誉，直到所有特征添加完畢。在這個過程中距帅，被添加的選定特征建立的模型具有最小的RSS右锨。所以理論上，第一個選定的特征應(yīng)
該能最好地解釋響應(yīng)變量碌秸，依此類推绍移。
添加一個特征一定會使RSS減少，使R方增加讥电，但不一定能提高模型的擬合度和可解釋性蹂窖。
后向逐步回歸
從一個包含所有特征的模型開始，每次刪除一個起最小作用的特征《鞯校現(xiàn)在也有一種混合方法瞬测，這種算法先通過前向逐步回歸添加特征，然后檢查是否有特征不再對提高模型擬合度起作用，如果有則刪除涣楷。每次建模之后分唾，分析者都可以檢查模型輸出抗碰，并使用各種統(tǒng)計量選擇能提供最佳擬合的特征狮斗。
逐步回歸技術(shù)會遇到非常嚴重的問題。對于一個數(shù)據(jù)集弧蝇，你先用前向逐步回歸碳褒，然后再用后向逐步回歸，可能會得到兩個完全矛盾的模型看疗。最重要的一點是沙峻，逐步回歸會使回歸系數(shù)發(fā)生偏離.
最優(yōu)子集回歸
最優(yōu)子集回歸是逐步回歸的一個可接受的替代方案。在最優(yōu)子集回歸中两芳，算法使用所有可能的特征組合來擬合模型摔寨，所以，如果有3個特征怖辆，將生成7個模型是复。然后和逐步回
歸一樣，分析者需要應(yīng)用自己的判斷和統(tǒng)計分析來選擇最優(yōu)模型竖螃，模型選擇就是后面工作的關(guān)鍵淑廊。
如果數(shù)據(jù)集有多個特征，工作量就會非常大特咆。當(dāng)特征數(shù)多于觀測數(shù)時（p大于n）季惩，這個方法的效果就不會好。
可以使用lm()函數(shù)腻格。模型形式為fit = lm (y~x1 + x2 +x3 + … + xn)画拾。
library(leaps)
fit <- lm(BSAAM ~ ., data = socal.water)
summary(fit)
與單變量回歸一樣，我們要檢查F統(tǒng)計量的p值菜职，以檢驗是否至少有一個非零系數(shù)碾阁。確實，p
值是高度顯著的些楣。還可以看到脂凶，OPRC和OPSLAKE這兩個參數(shù)具有顯著的p值。有趣的是愁茁，OPBPC并不顯著蚕钦，盡管它與響應(yīng)變量高度相關(guān)。簡言之鹅很，當(dāng)我們控制其他OP開頭的特征時嘶居，OPBPC無法
對預(yù)測方差提供任何有意義的解釋。這就是說，模型中存在OPRC和OPSLAKE時邮屁，特征OPBPC從統(tǒng)計學(xué)角度來看沒有任何作用整袁。
03 多變量線性回歸模型構(gòu)建與模型評價
使用最優(yōu)子集法。我們使用leaps包中的regsubsets()函數(shù)建立一個sub.fit對象佑吝。summary()就生成了best.summary對象坐昙，幫助我們更深入地研究模型。對于R中的所有對象芋忿，都可以使用names()函數(shù)列出輸出結(jié)果炸客。
sub.fit <- regsubsets(BSAAM ~ ., data = socal.water)
best.summary <- summary(sub.fit)
names(best.summary)
which.min(best.summary$rss)
# 有6個特征的模型具有最小的RSS
which.min()和which.max()分別給出具有某
個輸出的最小值和最大值的模型.
增加特征必然會減少RSS！而且必然會增加R方戈钢。
4種用于特征選擇的統(tǒng)計方法：赤池信息量準則痹仙、馬洛斯的Cp、貝葉斯信息量準則和修正R方殉了。前三種方法的目標是追求統(tǒng)計量的值最小化开仰，修正R方的目標是追求統(tǒng)計量的值最大化。這些統(tǒng)計方法的目的是建立一個盡可能簡約的模型薪铜。
在線性模型中众弓，AIC和Cp成正比，所以我們只需關(guān)注Cp痕囱。BIC與Cp相比田轧，更傾向于選擇變量較少的模型，所以我們要對二者進行比較鞍恢。
par(mfrow = c(1,2))
plot(best.summary$cp, xlab = "number of features", ylab = "cp")
plot(sub.fit, scale = "Cp")
效果如下：
在左側(cè)的圖中傻粘，可以看出有3個特征的模型具有最小的Cp值。在右側(cè)的圖中帮掉，顯示了能給出
最小Cp的特征組合弦悉。這張圖應(yīng)該這么看：先在Y軸的最高點找到最小的Cp值，此處是1.2蟆炊；然后
向右在X軸上找到與之對應(yīng)的色塊稽莉。通過這張圖，我們可以看到這個具有最小Cp值的模型中的3個特征是APSLAKE涩搓、OPRC和OPSLAKE污秆。通過which.min()和which.max()函數(shù)，我們可以進行Cp與BIC和修正R方的比較昧甘。
which.min(best.summary$bic)
which.max(best.summary$adjr2)
可以看出良拼，在本例中，BIC充边、修正R方和Cp選擇的最優(yōu)模型是一致的∮雇疲現(xiàn)在，與單變量線性
回歸一樣，我們需要檢查模型并進行假設(shè)檢驗贬媒。正像之前做的那樣聋亡，我們需要生成一個線性模型對象并檢查統(tǒng)計圖。
用這三個變量構(gòu)建模型：
best.fit <- lm(BSAAM ~ APSLAKE + OPRC + OPSLAKE, data = socal.water)
summary(best.fit)
par(mfrow = c(2,2))
plot(best.fit)
關(guān)注R小鹽际乘，關(guān)注科研私家菜(VX_GZH: SciPrivate)坡倔，有問題請聯(lián)系R小鹽。讓我們一起來學(xué)習(xí) R語言機器學(xué)習(xí)與臨床預(yù)測模型
