Glmnet是一個(gè)通過懲罰最大似然關(guān)系擬合廣義線性模型的軟件包。正則化路徑是針對正則化參數(shù)λ的值網(wǎng)格處的lasso或Elastic Net(彈性網(wǎng)絡(luò))懲罰值計(jì)算的。該算法非??欤⑶铱梢岳幂斎刖仃囍械南∈栊?nbsp;x
。它適合線性,邏輯和多項(xiàng)式,泊松和Cox回歸模型??梢詮臄M合模型中做出各種預(yù)測。它也可以擬合多元線性回歸。
glmnet
解決以下問題
在覆蓋整個(gè)范圍的λ值網(wǎng)格上。這里l(y,η)是觀察i的負(fù)對數(shù)似然貢獻(xiàn);例如對于高斯分布是
眾所周知,嶺懲罰使相關(guān)預(yù)測因子的系數(shù)彼此縮小,而套索傾向于選擇其中一個(gè)而丟棄其他預(yù)測因子。彈性網(wǎng)絡(luò)則將這兩者混合在一起。
glmnet
算法使用循環(huán)坐標(biāo)下降法,該方法在每個(gè)參數(shù)固定不變的情況下連續(xù)優(yōu)化目標(biāo)函數(shù),并反復(fù)循環(huán)直到收斂,我們的算法可以非??焖俚赜?jì)算求解路徑。
代碼可以處理稀疏的輸入矩陣格式,以及系數(shù)的范圍約束,還包括用于預(yù)測和繪圖的方法,以及執(zhí)行K折交叉驗(yàn)證的功能。
首先,我們加載 glmnet
包:
library(glmnet)
該命令 從此保存的R數(shù)據(jù)中加載輸入矩陣 x
和因向量 y
。
我們擬合模型 glmnet
。
fit = glmnet(x, y)
plot
函數(shù)來可視化系數(shù) :plot(fit)
每條曲線對應(yīng)一個(gè)變量。它顯示了當(dāng)λ變化時(shí),其系數(shù)相對于整個(gè)系數(shù)向量的?1范數(shù)的路徑。上方的軸表示當(dāng)前λ處非零系數(shù)的數(shù)量,這是套索的有效自由度(df)。用戶可能還希望對曲線進(jìn)行注釋。這可以通過label = TRUE
在plot命令中進(jìn)行設(shè)置來完成 。
glmnet
如果我們只是輸入對象名稱或使用print
函數(shù),則會顯示每個(gè)步驟的路徑 摘要 :
print(fit)
##
## Call: glmnet(x = x, y = y)
##
## Df %Dev Lambda
## [1,] 0 0.0000 1.63000
## [2,] 2 0.0553 1.49000
## [3,] 2 0.1460 1.35000
## [4,] 2 0.2210 1.23000
## [5,] 2 0.2840 1.12000
## [6,] 2 0.3350 1.02000
## [7,] 4 0.3900 0.93300
## [8,] 5 0.4560 0.85000
## [9,] 5 0.5150 0.77500
## [10,] 6 0.5740 0.70600
## [11,] 6 0.6260 0.64300
## [12,] 6 0.6690 0.58600
## [13,] 6 0.7050 0.53400
## [14,] 6 0.7340 0.48700
## [15,] 7 0.7620 0.44300
## [16,] 7 0.7860 0.40400
## [17,] 7 0.8050 0.36800
## [18,] 7 0.8220 0.33500
## [19,] 7 0.8350 0.30600
## [20,] 7 0.8460 0.27800
它從左到右顯示了非零系數(shù)的數(shù)量(Df
),解釋的(零)偏差百分比(%dev
)和λ(Lambda
)的值。
我們可以在序列范圍內(nèi)獲得一個(gè)或多個(gè)λ處的實(shí)際系數(shù):
coef(fit,s=0.1)
## 21 x 1 sparse Matrix of class "dgCMatrix"
## 1
## (Intercept) 0.150928
## V1 1.320597
## V2 .
## V3 0.675110
## V4 .
## V5 -0.817412
## V6 0.521437
## V7 0.004829
## V8 0.319416
## V9 .
## V10 .
## V11 0.142499
## V12 .
## V13 .
## V14 -1.059979
## V15 .
## V16 .
## V17 .
## V18 .
## V19 .
## V20 -1.021874
還可以使用新的輸入數(shù)據(jù)在特定的λ處進(jìn)行預(yù)測:
predict(fit,newx=nx,s=c(0.1,0.05))
## 1 2
## [1,] 4.4641 4.7001
## [2,] 1.7509 1.8513
## [3,] 4.5207 4.6512
## [4,] -0.6184 -0.6764
## [5,] 1.7302 1.8451
## [6,] 0.3565 0.3512
## [7,] 0.2881 0.2662
## [8,] 2.7776 2.8209
## [9,] -3.7016 -3.7773
## [10,] 1.1546 1.1067
glmnet
返回一系列模型供用戶選擇。交叉驗(yàn)證可能是該任務(wù)最簡單,使用最廣泛的方法。cv.glmnet
是交叉驗(yàn)證的主要函數(shù)。
cv.glmnet
返回一個(gè) cv.glmnet
對象,此處為“ cvfit”,其中包含交叉驗(yàn)證擬合的所有成分的列表。
我們可以繪制對象。
它包括交叉驗(yàn)證曲線(紅色虛線)和沿λ序列的上下標(biāo)準(zhǔn)偏差曲線(誤差線)。垂直虛線表示兩個(gè)選定的λ。
我們可以查看所選的λ和相應(yīng)的系數(shù)。例如,
cvfit$lambda.min
## [1] 0.08307
lambda.min
是給出最小平均交叉驗(yàn)證誤差的λ值。保存的另一個(gè)λ是 lambda.1se
,它給出了的模型,使得誤差在最小值的一個(gè)標(biāo)準(zhǔn)誤差以內(nèi)。我們只需要更換 lambda.min
到lambda.1se
以上。
coef(cvfit, s = "lambda.min")
## 21 x 1 sparse Matrix of class "dgCMatrix"
## 1
## (Intercept) 0.14936
## V1 1.32975
## V2 .
## V3 0.69096
## V4 .
## V5 -0.83123
## V6 0.53670
## V7 0.02005
## V8 0.33194
## V9 .
## V10 .
## V11 0.16239
## V12 .
## V13 .
## V14 -1.07081
## V15 .
## V16 .
## V17 .
## V18 .
## V19 .
## V20 -1.04341
可以根據(jù)擬合的cv.glmnet
對象進(jìn)行預(yù)測 。讓我們看一個(gè)示例。
# 1
## [1,] -1.3647
## [2,] 2.5686
## [3,] 0.5706
## [4,] 1.9682
## [5,] 1.4964
newx
與新的輸入矩陣 s
相同,如前所述,是預(yù)測的λ值。
這里的線性回歸是指兩個(gè)模型系列。一個(gè)是 gaussian
正態(tài)分布,另一個(gè)是 mgaussian
多元正態(tài)分布。
假設(shè)我們有觀測值xi∈Rp并且yi∈R,i = 1,...,N。目標(biāo)函數(shù)是
其中λ≥0是復(fù)雜度參數(shù),0≤α≤1在嶺回歸(α=0)和套索LASSO(α=1)之間。
應(yīng)用坐標(biāo)下降法解決該問題。具體地說,通過計(jì)算βj=β?j處的梯度和簡單的演算,更新為
其中
當(dāng)x
變量標(biāo)準(zhǔn)化為具有單位方差(默認(rèn)值)時(shí),以上公式適用 。
glmnet
提供各種選項(xiàng)供用戶自定義。我們在這里介紹一些常用的選項(xiàng),它們可以在glmnet
函數(shù)中指定 。
alpha
表示彈性網(wǎng)混合參數(shù)α,范圍α∈[0,1]。α=1是套索(默認(rèn)),α=0是Ridge。
weights
用于觀察權(quán)重。每個(gè)觀察值的默認(rèn)值為1。
nlambda
是序列中λ值的數(shù)量。默認(rèn)值為100。
lambda
可以提供,但通常不提供,程序會構(gòu)建一個(gè)序列。自動生成時(shí),λ序列由lambda.max
和 確定 lambda.min.ratio
。
standardize
是x
在擬合模型序列之前進(jìn)行變量標(biāo)準(zhǔn)化的邏輯標(biāo)志 。
例如,我們設(shè)置α=0.2,并對后半部分的觀測值賦予兩倍的權(quán)重。為了避免在此處顯示太長時(shí)間,我們將其設(shè)置 nlambda
為20。但是,實(shí)際上,建議將λ的數(shù)量設(shè)置為100(默認(rèn)值)或更多。
然后我們可以輸出glmnet
對象。
print(fit)
##
## Call: glmnet(x = x, y = y, weights = c(rep(1, 50), rep(2, 50)), alpha = 0.2, nlambda = 20)
##
## Df %Dev Lambda
## [1,] 0 0.000 7.94000
## [2,] 4 0.179 4.89000
## [3,] 7 0.444 3.01000
## [4,] 7 0.657 1.85000
## [5,] 8 0.785 1.14000
## [6,] 9 0.854 0.70300
## [7,] 10 0.887 0.43300
## [8,] 11 0.902 0.26700
## [9,] 14 0.910 0.16400
## [10,] 17 0.914 0.10100
## [11,] 17 0.915 0.06230
## [12,] 17 0.916 0.03840
## [13,] 19 0.916 0.02360
## [14,] 20 0.916 0.01460
## [15,] 20 0.916 0.00896
## [16,] 20 0.916 0.00552
## [17,] 20 0.916 0.00340
fit
以及帶有列Df
(非零系數(shù)的數(shù)量), %dev
(解釋的偏差百分比)和Lambda
(對應(yīng)的λ值) 的三列矩陣 。我們可以繪制擬合的對象。
讓我們針對log-lambda值標(biāo)記每個(gè)曲線來繪制“擬合”。
這是訓(xùn)練數(shù)據(jù)中的偏差百分比。我們在這里看到的是,在路徑末端時(shí),該值變化不大,但是系數(shù)有點(diǎn)“膨脹”。這使我們可以將注意力集中在重要的擬合部分上。
我們可以提取系數(shù)并在某些特定值的情況下進(jìn)行預(yù)測。兩種常用的選項(xiàng)是:
左列是,exact = TRUE
右列是 FALSE
。從上面我們可以看到,0.01不在序列中,因此盡管沒有太大差異,但還是有一些差異。如果沒有特殊要求,則線性插補(bǔ)就足夠了。
s 指定進(jìn)行提取的λ值。
exact 指示是否需要系數(shù)的精確值。
一個(gè)簡單的例子是:
## 21 x 2 sparse Matrix of class "dgCMatrix"
## 1 1
## (Intercept) 0.19657 0.199099
## V1 1.17496 1.174650
## V2 . .
## V3 0.52934 0.531935
## V4 . .
## V5 -0.76126 -0.760959
## V6 0.46627 0.468209
## V7 0.06148 0.061927
## V8 0.38049 0.380301
## V9 . .
## V10 . .
## V11 0.14214 0.143261
## V12 . .
## V13 . .
## V14 -0.91090 -0.911207
## V15 . .
## V16 . .
## V17 . .
## V18 . 0.009197
## V19 . .
## V20 -0.86099 -0.863117
用戶可以根據(jù)擬合的對象進(jìn)行預(yù)測。除中的選項(xiàng)外 coef
,主要參數(shù)是 newx
的新值矩陣 x
。type
選項(xiàng)允許用戶選擇預(yù)測類型:*“鏈接”給出擬合值
因變量與正態(tài)分布的“鏈接”相同。
“系數(shù)”計(jì)算值為的系數(shù) s
例如,
## 1
## [1,] -0.9803
## [2,] 2.2992
## [3,] 0.6011
## [4,] 2.3573
## [5,] 1.7520
給出在λ=0.05時(shí)前5個(gè)觀測值的擬合值。如果提供的多個(gè)值, s
則會生成預(yù)測矩陣。
用戶可以自定義K折交叉驗(yàn)證。除所有 glmnet
參數(shù)外, cv.glmnet
還有特殊的參數(shù),包括 nfolds
(次數(shù)), foldid
(用戶提供的次數(shù)), type.measure
(用于交叉驗(yàn)證的損失):*“ deviance”或“ mse”
“ mae”使用平均絕對誤差
舉個(gè)例子,
cvfit = cv.glmnet(x, y, type.measure = "mse", nfolds = 20)
根據(jù)均方誤差標(biāo)準(zhǔn)進(jìn)行20折交叉驗(yàn)證。
并行計(jì)算也受 cv.glmnet
。為我們在這里給出一個(gè)簡單的比較示例。
system.time(cv.glmnet(X, Y))
## user system elapsed
## 3.591 0.103 3.724
system.time(cv.glmnet(X, Y, parallel = TRUE))
## user system elapsed
## 4.318 0.391 2.700
從上面的建議可以看出,并行計(jì)算可以大大加快計(jì)算過程。
“ lambda.min”:達(dá)到最小MSE的λ。在這里,我們使用相同的k折,為α選擇一個(gè)值。
cvfit$lambda.min
## [1] 0.08307
## 21 x 1 sparse Matrix of class "dgCMatrix"
## 1
## (Intercept) 0.14936
## V1 1.32975
## V2 .
## V3 0.69096
## V4 .
## V5 -0.83123
## V6 0.53670
## V7 0.02005
## V8 0.33194
## V9 .
## V10 .
## V11 0.16239
## V12 .
## V13 .
## V14 -1.07081
## V15 .
## V16 .
## V17 .
## V18 .
## V19 .
## V20 -1.04341
將它們?nèi)糠胖迷谕焕L圖上:
我們看到lasso(alpha=1
)在這里表現(xiàn)最好。
假設(shè)我們要擬合我們的模型,但將系數(shù)限制為大于-0.7且小于0.5。這可以通過upper.limits
和 lower.limits
參數(shù)實(shí)現(xiàn) :
通常,我們希望系數(shù)為正,因此我們只能lower.limit
將其設(shè)置 為0。
此參數(shù)允許用戶將單獨(dú)的懲罰因子應(yīng)用于每個(gè)系數(shù)。每個(gè)參數(shù)的默認(rèn)值為1,但可以指定其他值。特別是,任何penalty.factor
等于零的變量 都不會受到懲罰
在許多情況下,某些變量可能是重要,我們希望一直保留它們,這可以通過將相應(yīng)的懲罰因子設(shè)置為0來實(shí)現(xiàn):
我們從標(biāo)簽中看到懲罰因子為0的三個(gè)變量始終保留在模型中,而其他變量遵循典型的正則化路徑并最終縮小為0。
有時(shí),尤其是在變量數(shù)量很少的情況下,我們想在圖上添加變量標(biāo)簽。
我們首先生成帶有10個(gè)變量的一些數(shù)據(jù),然后,我們擬合glmnet模型,并繪制標(biāo)準(zhǔn)圖。
我們希望用變量名標(biāo)記曲線。在路徑的末尾放置系數(shù)的位置。
使用family = "mgaussian"
option 獲得多元正態(tài)分布glmnet
。
顯然,顧名思義,y不是向量,而是矩陣。結(jié)果,每個(gè)λ值的系數(shù)也是一個(gè)矩陣。
在這里,我們解決以下問題:
這里,βj是p×K系數(shù)矩陣β的第j行,對于單個(gè)預(yù)測變量xj,我們用每個(gè)系數(shù)K向量βj的組套索罰分代替每個(gè)單一系數(shù)的絕對罰分。
我們使用預(yù)先生成的一組數(shù)據(jù)進(jìn)行說明。
我們擬合數(shù)據(jù),并返回對象“ mfit”。
mfit = glmnet(x, y, family = "mgau
ssian")
如果為 standardize.response = TRUE
,則將因變量標(biāo)準(zhǔn)化。
為了可視化系數(shù),我們使用 plot
函數(shù)。
注意我們設(shè)置了 type.coef = "2norm"
。在此設(shè)置下,每個(gè)變量繪制一條曲線,其值等于?2范數(shù)。默認(rèn)設(shè)置為 type.coef = "coef"
,其中為每個(gè)因變量創(chuàng)建一個(gè)系數(shù)圖。
通過使用該函數(shù)coef
,我們可以提取要求的λ值的系數(shù), 并通過進(jìn)行預(yù)測 。
## , , 1
##
## y1 y2 y3 y4
## [1,] -4.7106 -1.1635 0.6028 3.741
## [2,] 4.1302 -3.0508 -1.2123 4.970
## [3,] 3.1595 -0.5760 0.2608 2.054
## [4,] 0.6459 2.1206 -0.2252 3.146
## [5,] -1.1792 0.1056 -7.3353 3.248
##
## , , 2
##
## y1 y2 y3 y4
## [1,] -4.6415 -1.2290 0.6118 3.780
## [2,] 4.4713 -3.2530 -1.2573 5.266
## [3,] 3.4735 -0.6929 0.4684 2.056
## [4,] 0.7353 2.2965 -0.2190 2.989
## [5,] -1.2760 0.2893 -7.8259 3.205
預(yù)測結(jié)果保存在三維數(shù)組中,其中前兩個(gè)維是每個(gè)因變量的預(yù)測矩陣,第三個(gè)維表示因變量。
我們還可以進(jìn)行k折交叉驗(yàn)證。
我們繪制結(jié)果 cv.glmnet
對象“ cvmfit”。
顯示選定的λ最佳值
cvmfit$lambda.min
## [1] 0.04732
cvmfit$lambda.1se
## [1] 0.1317
邏輯回歸
當(dāng)因變量是分類的時(shí),邏輯回歸是另一個(gè)廣泛使用的模型。如果有兩個(gè)可能的結(jié)果,則使用二項(xiàng)式分布,否則使用多項(xiàng)式。
對于二項(xiàng)式模型,假設(shè)因變量的取值為G = {1,2} 。表示yi = I(gi = 1)。我們建模
可以用以下形式寫
懲罰邏輯回歸的目標(biāo)函數(shù)使用負(fù)二項(xiàng)式對數(shù)似然
我們的算法使用對數(shù)似然的二次逼近,然后對所得的懲罰加權(quán)最小二乘問題進(jìn)行下降。這些構(gòu)成了內(nèi)部和外部循環(huán)。
出于說明目的,我們 從數(shù)據(jù)文件加載預(yù)生成的輸入矩陣 x
和因變量 y
。
對于二項(xiàng)式邏輯回歸,因變量y可以是兩個(gè)級別的因子,也可以是計(jì)數(shù)或比例的兩列矩陣。
glmnet
二項(xiàng)式回歸的其他可選參數(shù)與正態(tài)分布的參數(shù) 幾乎相同。不要忘記將family
選項(xiàng)設(shè)置 為“ binomial”。
fit = glmnet(x, y, family = "binomial")
邏輯回歸略有不同,主要體現(xiàn)在選擇上 type
。“鏈接”和“因變量”不等價(jià),“類”僅可用于邏輯回歸??傊?“鏈接”給出了線性預(yù)測變量
“因變量”給出合適的概率
“類別”產(chǎn)生對應(yīng)于最大概率的類別標(biāo)簽。
“系數(shù)”計(jì)算值為的系數(shù) s
在下面的示例中,我們在λ=0.05,0.01的情況下對類別標(biāo)簽進(jìn)行了預(yù)測。
## 1 2
## [1,] "0" "0"
## [2,] "1" "1"
## [3,] "1" "1"
## [4,] "0" "0"
## [5,] "1" "1"
對于邏輯回歸,type.measure
:
“偏差”使用實(shí)際偏差。
“ mae”使用平均絕對誤差。
“class”給出錯(cuò)誤分類錯(cuò)誤。
“ auc”(僅適用于兩類邏輯回歸)給出了ROC曲線下的面積。
例如,
它使用分類誤差作為10倍交叉驗(yàn)證的標(biāo)準(zhǔn)。
我們繪制對象并顯示λ的最佳值。
cvfit$lambda.min
## [1] 0.01476
cvfit$lambda.1se
## [1] 0.02579
coef
并且 predict
類似于正態(tài)分布案例,因此我們省略了細(xì)節(jié)。我們通過一些例子進(jìn)行回顧。
## 31 x 1 sparse Matrix of class "dgCMatrix"
## 1
## (Intercept) 0.24371
## V1 0.06897
## V2 0.66252
## V3 -0.54275
## V4 -1.13693
## V5 -0.19143
## V6 -0.95852
## V7 .
## V8 -0.56529
## V9 0.77454
## V10 -1.45079
## V11 -0.04363
## V12 -0.06894
## V13 .
## V14 .
## V15 .
## V16 0.36685
## V17 .
## V18 -0.04014
## V19 .
## V20 .
## V21 .
## V22 0.20882
## V23 0.34014
## V24 .
## V25 0.66310
## V26 -0.33696
## V27 -0.10570
## V28 0.24318
## V29 -0.22445
## V30 0.11091
如前所述,此處返回的結(jié)果僅針對因子因變量的第二類。
## 1
## [1,] "0"
## [2,] "1"
## [3,] "1"
## [4,] "0"
## [5,] "1"
## [6,] "0"
## [7,] "0"
## [8,] "0"
## [9,] "1"
## [10,] "1"
多項(xiàng)式模型
對于多項(xiàng)式模型,假設(shè)因變量變量的K級別為G = {1,2,…,K}。在這里我們建模
設(shè)Y為N×K指標(biāo)因變量矩陣,元素yi?= I(gi =?)。然后彈性網(wǎng)懲罰的負(fù)對數(shù)似然函數(shù)變?yōu)?/p>
β是系數(shù)的p×K矩陣。βk指第k列(對于結(jié)果類別k),βj指第j行(變量j的K個(gè)系數(shù)的向量)。最后一個(gè)懲罰項(xiàng)是||βj|| q ,我們對q有兩個(gè)選擇:q∈{1,2}。當(dāng)q = 1時(shí),這是每個(gè)參數(shù)的套索懲罰。當(dāng)q = 2時(shí),這是對特定變量的所有K個(gè)系數(shù)的分組套索懲罰,這使它們在一起全為零或非零。
對于多項(xiàng)式情況,用法類似于邏輯回歸,我們加載一組生成的數(shù)據(jù)。
glmnet
除少數(shù)情況外,多項(xiàng)式邏輯回歸中的可選參數(shù) 與二項(xiàng)式回歸基本相似。
多項(xiàng)式回歸的一個(gè)特殊選項(xiàng)是 type.multinomial
,如果允許,則允許使用分組的套索罰分 type.multinomial = "grouped"
。這將確保變量的多項(xiàng)式系數(shù)全部一起輸入或輸出,就像多元因變量一樣。
我們繪制結(jié)果。
我們還可以進(jìn)行交叉驗(yàn)證并繪制返回的對象。
預(yù)測最佳選擇的λ:
## 1
## [1,] "3"
## [2,] "2"
## [3,] "2"
## [4,] "1"
## [5,] "1"
## [6,] "3"
## [7,] "3"
## [8,] "1"
## [9,] "1"
## [10,] "2"
Poisson回歸用于在假設(shè)Poisson誤差的情況下對計(jì)數(shù)數(shù)據(jù)進(jìn)行建模,或者在均值和方差成比例的情況下使用非負(fù)數(shù)據(jù)進(jìn)行建模。泊松也是指數(shù)分布族的成員。我們通常以對數(shù)建模:
和以前一樣,我們優(yōu)化了懲罰對數(shù):
Glmnet使用外部牛頓循環(huán)和內(nèi)部加權(quán)最小二乘循環(huán)(如邏輯回歸)來優(yōu)化此標(biāo)準(zhǔn)。
首先,我們加載一組泊松數(shù)據(jù)。
再次,繪制系數(shù)。
像以前一樣,我們可以 分別使用coef
和 提取系數(shù)并在特定的λ處進(jìn)行預(yù)測 predict
。
例如,我們可以
## 21 x 1 sparse Matrix of class "dgCMatrix"
## 1
## (Intercept) 0.61123
## V1 0.45820
## V2 -0.77061
## V3 1.34015
## V4 0.04350
## V5 -0.20326
## V6 .
## V7 .
## V8 .
## V9 .
## V10 .
## V11 .
## V12 0.01816
## V13 .
## V14 .
## V15 .
## V16 .
## V17 .
## V18 .
## V19 .
## V20 .
## 1 2
## [1,] 2.4944 4.4263
## [2,] 10.3513 11.0586
## [3,] 0.1180 0.1782
## [4,] 0.9713 1.6829
## [5,] 1.1133 1.9935
我們還可以使用交叉驗(yàn)證來找到最佳的λ,從而進(jìn)行推斷。
選項(xiàng)幾乎與正態(tài)族相同,不同之處在于 type.measure
*,“ mse”代表均方誤差*,“ mae”代表均值絕對誤差。
我們可以繪制 cv.glmnet
對象。
我們還可以顯示最佳的λ和相應(yīng)的系數(shù)。
## 21 x 2 sparse Matrix of class "dgCMatrix"
## 1 2
## (Intercept) 0.031263 0.18570
## V1 0.619053 0.57537
## V2 -0.984550 -0.93212
## V3 1.525234 1.47057
## V4 0.231591 0.19692
## V5 -0.336659 -0.30469
## V6 0.001026 .
## V7 -0.012830 .
## V8 . .
## V9 . .
## V10 0.015983 .
## V11 . .
## V12 0.030867 0.02585
## V13 -0.027971 .
## V14 0.032750 .
## V15 -0.005933 .
## V16 0.017506 .
## V17 . .
## V18 0.004026 .
## V19 -0.033579 .
## V20 0.012049 0.00993
Cox比例風(fēng)險(xiǎn)模型通常用于研究預(yù)測變量與生存時(shí)間之間的關(guān)系。
Cox比例風(fēng)險(xiǎn)回歸模型,它不是直接考察
式中,
由于Cox回歸模型對
公式可以轉(zhuǎn)化為:
我們使用一組預(yù)先生成的樣本數(shù)據(jù)。用戶可以加載自己的數(shù)據(jù)并遵循類似的過程。在這種情況下,x必須是協(xié)變量值的n×p矩陣-每行對應(yīng)一個(gè)患者,每列對應(yīng)一個(gè)協(xié)變量。y是一個(gè)n×2矩陣。
## time status
## [1,] 1.76878 1
## [2,] 0.54528 1
## [3,] 0.04486 0
## [4,] 0.85032 0
## [5,] 0.61488 1
Surv
包中的 函數(shù) survival
可以創(chuàng)建這樣的矩陣。
我們計(jì)算默認(rèn)設(shè)置下的求解路徑。
繪制系數(shù)。
提取特定值λ處的系數(shù)。
## 30 x 1 sparse Matrix of class "dgCMatrix"
## 1
## V1 0.37694
## V2 -0.09548
## V3 -0.13596
## V4 0.09814
## V5 -0.11438
## V6 -0.38899
## V7 0.24291
## V8 0.03648
## V9 0.34740
## V10 0.03865
## V11 .
## V12 .
## V13 .
## V14 .
## V15 .
## V16 .
## V17 .
## V18 .
## V19 .
## V20 .
## V21 .
## V22 .
## V23 .
## V24 .
## V25 .
## V26 .
## V27 .
## V28 .
## V29 .
## V30 .
函數(shù) cv.glmnet
可用于計(jì)算Cox模型的k折交叉驗(yàn)證。
擬合后,我們可以查看最佳λ值和交叉驗(yàn)證的誤差圖,幫助評估我們的模型。
如前所述,圖中的左垂直線向我們顯示了CV誤差曲線達(dá)到最小值的位置。右邊的垂直線向我們展示了正則化的模型,其CV誤差在最小值的1個(gè)標(biāo)準(zhǔn)偏差之內(nèi)。我們還提取了最優(yōu)λ。
cvfit$lambda.min
## [1] 0.01594
cvfit$lambda.1se
## [1] 0.04869
我們可以檢查模型中的協(xié)變量并查看其系數(shù)。
index.min
## [1] 0.491297 -0.174601 -0.218649 0.175112 -0.186673 -0.490250 0.335197
## [8] 0.091587 0.450169 0.115922 0.017595 -0.018365 -0.002806 -0.001423
## [15] -0.023429 0.001688 -0.008236
coef.min
## 30 x 1 sparse Matrix of class "dgCMatrix"
## 1
## V1 0.491297
## V2 -0.174601
## V3 -0.218649
## V4 0.175112
## V5 -0.186673
## V6 -0.490250
## V7 0.335197
## V8 0.091587
## V9 0.450169
## V10 0.115922
## V11 .
## V12 .
## V13 0.017595
## V14 .
## V15 .
## V16 .
## V17 -0.018365
## V18 .
## V19 .
## V20 .
## V21 -0.002806
## V22 -0.001423
## V23 .
## V24 .
## V25 -0.023429
## V26 .
## V27 0.001688
## V28 .
## V29 .
## V30 -0.008236
我們的程序包支持稀疏的輸入矩陣,該矩陣可以高效地存儲和操作大型矩陣,但只有少數(shù)幾個(gè)非零條目。
我們加載一組預(yù)先創(chuàng)建的樣本數(shù)據(jù)。
加載100 * 20的稀疏矩陣和 y
因向量。
## [1] "dgCMatrix"
## attr(,"package")
## [1] "Matrix"
我們可以像以前一樣擬合模型。
fit = glmnet(x, y)
預(yù)測新輸入矩陣 。例如,
## 1
## [1,] 0.3826
## [2,] -0.2172
## [3,] -1.6622
## [4,] -0.4175
## [5,] -1.3941
Jerome Friedman, Trevor Hastie and Rob Tibshirani. (2008).
Regularization Paths for Generalized Linear Models via Coordinate Descent