Cross Validated is a question and answer site for people interested in statistics, machine learning, data analysis, data mining, and data visualization. It only takes a minute to sign up.
Sign up to join this community
I'm dealing with a large p small n problem (p>1000, n=150) and have decided to use glmnet and lars to select possible explanatory variables (as I understand is their functionality) rather than do this manually.
I have gotten out both a lars, and a glmnet object from the functions, but now can't discern how to extract from those objects, the names of the coefficients / the optimal model, that the functions have chosen.
Any help would be greatly appreciated!
In order to determine the optimal model you'll need to estimate an optimal penalization parameter $\lambda$. This can be done using the built in cross-validation function in glmnet called cv.glmnet. Once you've fit your cross validated model with cv.glmnet the model object will contain all sorts of useful values. These are described under the Values header in the cv.glmnet help file. They can be extracted with the typical extractors, for example coef(), or using the $ operator on the cv.glmnet object.
A word of warning though: be very careful with the default glmnet arguments. Make sure you understand what each one is doing and why it's appropriate for your specific problem.
So almost verbatim from the cv.glmnet help file:
> ## simulate data
> set.seed(1010)
> n=1000;p=100
> nzc=trunc(p/10)
> x=matrix(rnorm(n*p),n,p)
> beta=rnorm(nzc)
> fx= x[,seq(nzc)] %*% beta
> eps=rnorm(n)*5
> y=drop(fx+eps)
> ## run cross-validation with all the default argumnets
> cvob1=cv.glmnet(x,y)
> ## one way to print the coefficients of the optimal model
> coef(cvob1)
101 x 1 sparse Matrix of class "dgCMatrix"
1
(Intercept) -0.1162737
V1 -0.2171531
V2 0.3237422
V3 .
V4 -0.2190339
V5 -0.1856601
V6 0.2530652
V7 0.1874832
V8 -1.3574323
V9 1.0162046
V10 0.1558299
V11 .
V12 .
V13 .
V14 .
V15 .
V16 .
V17 .
V18 .
V19 .
V20 .
V21 .
V22 .
V23 .
V24 .
V25 .
V26 .
V27 .
V28 .
V29 .
V30 .
V31 .
V32 .
V33 .
V34 .
V35 .
V36 .
V37 .
V38 .
V39 .
V40 .
V41 .
V42 .
V43 .
V44 .
V45 .
V46 .
V47 .
V48 .
V49 .
V50 .
V51 .
V52 .
V53 .
V54 .
V55 .
V56 .
V57 .
V58 .
V59 .
V60 .
V61 .
V62 .
V63 .
V64 .
V65 .
V66 .
V67 .
V68 .
V69 .
V70 .
V71 .
V72 .
V73 .
V74 .
V75 -0.1420966
V76 .
V77 .
V78 .
V79 .
V80 .
V81 .
V82 .
V83 .
V84 .
V85 .
V86 .
V87 .
V88 .
V89 .
V90 .
V91 .
V92 .
V93 .
V94 .
V95 .
V96 .
V97 .
V98 .
V99 .
V100 .
> ## a view of the cvob1 values
> str(cvob1)
List of 10
$ lambda : num [1:74] 1.78 1.62 1.48 1.35 1.23 ...
$ cvm : num [1:74] 31.3 30.8 30.3 29.8 29.3 ...
$ cvsd : num [1:74] 0.994 1.008 1.022 1.045 1.054 ...
$ cvup : num [1:74] 32.2 31.8 31.4 30.9 30.3 ...
$ cvlo : num [1:74] 30.3 29.8 29.3 28.8 28.2 ...
$ nzero : Named int [1:74] 0 1 1 2 2 2 2 2 2 2 ...
..- attr(*, "names")= chr [1:74] "s0" "s1" "s2" "s3" ...
$ name : Named chr "Mean-Squared Error"
..- attr(*, "names")= chr "mse"
$ glmnet.fit:List of 12
..$ a0 : Named num [1:74] -0.157 -0.156 -0.155 -0.153 -0.149 ...
.. ..- attr(*, "names")= chr [1:74] "s0" "s1" "s2" "s3" ...
..$ beta :Formal class 'dgCMatrix' [package "Matrix"] with 6 slots
.. .. ..@ i : int [1:4266] 7 7 7 8 7 8 7 8 7 8 ...
.. .. ..@ p : int [1:75] 0 0 1 2 4 6 8 10 12 14 ...
.. .. ..@ Dim : int [1:2] 100 74
.. .. ..@ Dimnames:List of 2
.. .. .. ..$ : chr [1:100] "V1" "V2" "V3" "V4" ...
.. .. .. ..$ : chr [1:74] "s0" "s1" "s2" "s3" ...
.. .. ..@ x : num [1:4266] -0.1565 -0.2991 -0.4299 0.0459 -0.5505 ...
.. .. ..@ factors : list()
..$ df : int [1:74] 0 1 1 2 2 2 2 2 2 2 ...
..$ dim : int [1:2] 100 74
..$ lambda : num [1:74] 1.78 1.62 1.48 1.35 1.23 ...
..$ dev.ratio: num [1:74] 0 0.0172 0.0315 0.0475 0.0676 ...
..$ nulldev : num 31247
..$ npasses : int 298
..$ jerr : int 0
..$ offset : logi FALSE
..$ call : language glmnet(x = x, y = y)
..$ nobs : int 1000
..- attr(*, "class")= chr [1:2] "elnet" "glmnet"
$ lambda.min: num 0.23
$ lambda.1se: num 0.402
- attr(*, "class")= chr "cv.glmnet"
> cvob1$lambda.min
[1] 0.2300227
