Bagaimana cara menginterpretasikan glmnet?

36

Saya mencoba untuk menyesuaikan model regresi linier multivariat dengan sekitar 60 variabel prediktor dan 30 pengamatan, jadi saya menggunakan paket glmnet untuk regresi yang diatur karena p> n.

Saya telah melalui dokumentasi dan pertanyaan lain tetapi saya masih belum dapat menginterpretasikan hasilnya, berikut ini contoh kode (dengan 20 prediktor dan 10 pengamatan untuk disederhanakan):

Saya membuat matriks x dengan num rows = num observasi dan num cols = num prediktor dan vektor y yang mewakili variabel respons

> x=matrix(rnorm(10*20),10,20)
> y=rnorm(10)

Saya cocok dengan model glmnet yang menjadikan alpha sebagai default (= 1 untuk penalti laso)

> fit1=glmnet(x,y)
> print(fit1)

Saya mengerti saya mendapatkan prediksi yang berbeda dengan penurunan nilai lambda (yaitu penalti)

Call:  glmnet(x = x, y = y) 

        Df    %Dev   Lambda
  [1,]  0 0.00000 0.890700
  [2,]  1 0.06159 0.850200
  [3,]  1 0.11770 0.811500
  [4,]  1 0.16880 0.774600
   .
   .
   .
  [96,] 10 0.99740 0.010730
  [97,] 10 0.99760 0.010240
  [98,] 10 0.99780 0.009775
  [99,] 10 0.99800 0.009331
 [100,] 10 0.99820 0.008907

Sekarang saya memperkirakan nilai Beta saya memilih, misalnya, nilai lambda terkecil yang diberikan glmnet

> predict(fit1,type="coef", s = 0.008907)

21 x 1 sparse Matrix of class "dgCMatrix"
                  1
(Intercept) -0.08872364
V1           0.23734885
V2          -0.35472137
V3          -0.08088463
V4           .         
V5           .         
V6           .         
V7           0.31127123
V8           .         
V9           .         
V10          .         
V11          0.10636867
V12          .         
V13         -0.20328200
V14         -0.77717745
V15          .         
V16         -0.25924281
V17          .         
V18          .         
V19         -0.57989929
V20         -0.22522859

Jika sebaliknya saya memilih lambda dengan

cv <- cv.glmnet(x,y)
model=glmnet(x,y,lambda=cv$lambda.min)

Semua variabel akan menjadi (.).

Keraguan dan pertanyaan:

  1. Saya tidak yakin tentang cara memilih lambda.
  2. Haruskah saya menggunakan variabel non (.) Agar sesuai dengan model lain? Dalam kasus saya, saya ingin menyimpan variabel sebanyak mungkin.
  3. Bagaimana saya tahu nilai-p, yaitu variabel mana yang secara signifikan memprediksi respons?

Saya minta maaf atas pengetahuan statistik saya yang buruk! Dan terima kasih atas bantuannya.

Alice
sumber
Mungkin kita lihat paket CRAN hdi , yang memberikan inferensi untuk model dimensi tinggi ...
Tom Wenseleers
Untuk penjelasan lengkap tentang metode yang digunakan, saya merujuk Anda ke makalah ini: projecteuclid.org/euclid.ss/1449670857
Tom Wenseleers

Jawaban:

40

Ini fakta yang tidak intuitif - Anda sebenarnya tidak seharusnya memberi glmnet nilai tunggal lambda. Dari dokumentasi di sini :

Jangan berikan nilai tunggal untuk lambda (untuk prediksi setelah penggunaan CV memprediksi () sebagai gantinya). Alih-alih, berikan urutan penurunan nilai lambda. glmnet bergantung pada penghangatnya yang dimulai untuk kecepatan, dan seringkali lebih cepat untuk menyesuaikan seluruh lintasan daripada menghitung satu kecocokan.

cv.glmnetakan membantu Anda memilih lambda, seperti yang Anda singgung dalam contoh Anda. Penulis paket glmnet menyarankan cv$lambda.1sebukan cv$lambda.min, tetapi dalam praktiknya saya telah sukses dengan yang terakhir.

Setelah menjalankan cv.glmnet, Anda tidak perlu menjalankan kembali glmnet! Setiap lambda di kisi ( cv$lambda) telah dijalankan. Teknik ini disebut "Mulai Hangat" dan Anda dapat membaca lebih lanjut di sini . Mengutip dari pengantar, teknik Mulai Hangat mengurangi waktu berjalan dari metode berulang dengan menggunakan solusi dari masalah optimasi yang berbeda (misalnya, glmnet dengan lambda yang lebih besar) sebagai nilai awal untuk masalah optimasi nanti (misalnya, glmnet dengan lambda yang lebih kecil) ).

Untuk mengekstrak run yang diinginkan cv.glmnet.fit, coba ini:

small.lambda.index <- which(cv$lambda == cv$lambda.min)
small.lambda.betas <- cv$glmnet.fit$beta[, small.lambda.index]

Revisi (28/1/2017)

Tidak perlu kembali ke objek glmnet seperti yang saya lakukan di atas; ikuti saran @ alex23lemm di bawah ini dan berikan nomor s = "lambda.min", s = "lambda.1se"atau nomor lain (misalnya, s = .007) kepada keduanya coefdan predict. Perhatikan bahwa koefisien dan prediksi Anda bergantung pada nilai ini yang ditetapkan dengan validasi silang. Gunakan benih untuk reproduksibilitas! Dan jangan lupa bahwa jika Anda tidak menyediakan "s"in coefdan predict, Anda akan menggunakan default dari s = "lambda.1se". Saya telah melakukan pemanasan ke default itu setelah melihatnya berfungsi lebih baik dalam situasi data kecil.s = "lambda.1se"juga cenderung memberikan lebih banyak regularisasi, jadi jika Anda bekerja dengan alpha> 0, itu juga akan cenderung ke arah model yang lebih pelit. Anda juga dapat memilih nilai numerik s dengan bantuan plot.glmnet untuk sampai ke suatu tempat di antaranya (jangan lupa untuk mengekspansi nilai dari sumbu x!).

Ben Ogorek
sumber
1
Terima kasih! Ini membantu ... apakah Anda mungkin punya jawaban untuk pertanyaan 2 dan 3?
Alice
3
Ha jangan khawatir. Tanda (.) Mewakili nol. Karena Anda menggunakan Lasso, Anda telah menentukan bahwa Anda menginginkan solusi "jarang" (yaitu, banyak nol). Jika Anda ingin semuanya memiliki nilai, atur alpha = 0. Sekarang Anda telah beralih dari Lasso ke regresi Ridge. nilai-p untuk glmnet secara konsep sulit. Jika Anda mencari "nilai-p untuk laso" di Google, misalnya, Anda akan melihat banyak penelitian dan debat terbaru. Saya bahkan membaca satu akun (sumber amnesia) di mana penulis berpendapat bahwa nilai-p tidak masuk akal untuk regresi bias seperti laso dan regresi ridge.
Ben Ogorek
6
Cara alternatif untuk mengekstrak koefisien yang terkait dengan nilai lambda yang memberikan cvm minimum adalah sebagai berikut:small.lambda.betas <- coef(cv, s = "lambda.min")
alex23lemm
1
@ Benengorek, pembaruan luar biasa! Referensi lain yang bermanfaat adalah Friedman J, Hastie T, Hoefling H, optimasi koordinat Tibshirani R. Pathwise. Statistik Statistik Terapan. 2007; 2 (1): 302–332. ( arxiv.org/pdf/0708.1485.pdf )
dv_bn
1
@erosennin, periksa argumen lambda dari cv.glmnet: "Urutan lambda yang disediakan oleh pengguna; standarnya adalah NULL, dan glmnet memilih urutannya sendiri." Anda akan ingin menggunakan prinsip awal hangat dan memulai urutan dengan beberapa nilai lambda yang lebih besar sebelum menurun ke kisaran yang Anda minati.
Ben Ogorek
2

Q1) Saya tidak yakin tentang cara memilih lambda. T2) Haruskah saya menggunakan variabel non (.) Agar sesuai dengan model lain? Dalam kasus saya, saya ingin menyimpan variabel sebanyak mungkin.

Sebagai per @ BenOgorek jawaban yang bagus, biasanya Anda membiarkan pas menggunakan seluruh urutan lambda, maka ketika mengekstraksi koefisien optimal menggunakan nilai lambda.1se (tidak seperti apa yang Anda lakukan).

Selama Anda mengikuti tiga peringatan di bawah ini, maka jangan melawan regularisasi atau mengubah model: jika suatu variabel dihilangkan, maka itu karena memberikan penalti keseluruhan yang lebih rendah. Peringatan adalah:

  1. Agar koefisien yang diatur menjadi bermakna, pastikan Anda secara eksplisit menormalkan mean dan stdev variabel sebelumnya scale(); jangan mengandalkan glmnet(standardize=T). Untuk pembenaran lihat, apakah standardisasi sebelum Lasso benar-benar diperlukan? ; pada dasarnya variabel dengan nilai-nilai besar mungkin dihukum secara tidak adil dalam regularisasi.

  2. Agar dapat direproduksi, jalankan dengan set.seedbeberapa benih acak dan periksa koefisien yang diatur untuk stabilitas.

  3. Jika Anda ingin regularisasi yang kurang keras yaitu lebih banyak variabel yang disertakan, gunakan alpha <1 (ie elastic-net) daripada bubungan sederhana. Saya sarankan Anda menyapu alpha dari 0 ke 1. Jika Anda akan melakukannya, maka untuk menghindari overfitting alpha hyperparameter dan kesalahan regresi, Anda harus menggunakan crossvalidation, yaitu menggunakan cv.glmnet()daripada sederhana glmnet():

.

for (alpha in c(0,.1,.3,.5,.7,.9,1)) {
  fit <- cv.glmnet(..., alpha=alpha, nfolds=...)
  # Look at the CVE at lambda.1se to find the minimum for this alpha value...
}

Jika Anda ingin mengotomatiskan pencarian jaringan seperti itu dengan CV, Anda dapat mengkodekannya sendiri atau menggunakan paket caret di atas glmnet; caret melakukan ini dengan baik. Untuk cv.glmnet nfoldsnilai parameter, pilih 3 (minimum) jika dataset Anda kecil, atau 5 atau 10 jika itu besar.

T3) Bagaimana saya tahu nilai-p, yaitu variabel mana yang secara signifikan memprediksi respons?

Jangan, itu tidak berarti . Seperti yang dijelaskan secara terperinci dalam Mengapa tidak disarankan untuk mendapatkan informasi ringkasan statistik untuk koefisien regresi dari model glmnet?

Biarkan saja cv.glmnet()lakukan pemilihan variabel secara otomatis. Dengan peringatan di atas. Dan tentu saja distribusi variabel respons harus normal (dengan asumsi Anda menggunakan family='gaussian').

smci
sumber
Terima kasih atas komentarnya yang sangat membantu! Saya juga mengalami bahwa standarisasi variabel itu sendiri tampaknya berfungsi daripada menggunakan glmnet (standardize = T).
Michelle
Saya punya pertanyaan @smci, tentang nilai beta yang dikembalikan oleh cvglmnet. Saya mengerti bahwa mereka adalah nilai beta pada setiap titik kisi dari nilai lambda yang dicoba. Namun, adalah nilai beta yang dikembalikan untuk setiap nilai lambda (1) nilai koefisien rata-rata dari 10 lipatan (dengan asumsi saya menggunakan 10 kali lipatCV), (2) nilai beta dari lipatan yang memberikan akurasi terbaik, atau (3) koefisien dari menjalankan kembali model pada seluruh dataset?
Michelle