How to export all coefficients of a penlized regression model from package `penalized`? Need them for reporting rolling regression estimate

Question

I have been able to run regression with some coefficients constrained to positive territory, but I'm doing alot of rolling regressions where I face the problem. Here is my sample code:

library(penalized)
set.seed(1)
x1=rnorm(100)*10
x2=rnorm(100)*10
x3=rnorm(100)*10
y=sin(x1)+cos(x2)-x3+rnorm(100)

data <- data.frame(y, x1, x2, x3)

win <- 10

coefs <- matrix(NA, ncol=4, nrow=length(y))

for(i in 1:(length(y)-win)) {
  d <- data[(1+i):(win+i),]
  p <- win+i
  # Linear Regression
  coefs[p,] <- as.vector(coef(penalized(y, ~ x1 + x2 + x3, ~1, 
                                      lambda1=0, lambda2=0, positive = c(F, F, T), data=data)))}

This is how I usually populate matrix with coefs from rolling regression and now I receive error:

Error in coefs[p, ] <- as.vector(coef(penalized(y, ~x1 + x2 + x3, ~1,  : 
  number of items to replace is not a multiple of replacement length

I assume that this error is produced because there is not always Intercept + 3 coefficients coming out of that penalized regression function. Is there away to get penalized function to show 0 coefs as well? or other way to populated matrix / data.frame?

Actually I would like to see those 0 coefs, as now I don't have a clue what coefs belong to what variable. I will be doing some timeseries analysis, so building y back from these will be pretty hard. — Hakki
The example I'm trying to reproduce ordinary linear model is used and trying to mimic that. So can't say purposely, but just mimicking example. Only diffenrence is that they constrain some coefs to be positive for fundamental reason. — Hakki

Zheyuan Li Zheyuan Li · Accepted Answer · 2016-09-28T17:51:12

Perhaps you are unaware of the which argument for coef for "penfit" object. Have a look at:

getMethod(coef, "penfit")

#function (object, ...) 
#{
#    .local <- function (object, which = c("nonzero", "all", "penalized", 
#        "unpenalized"), standardize = FALSE) 
#    {
#        coefficients(object, which, standardize)
#    }
#    .local(object, ...)
#}
#<environment: namespace:penalized>

We can set which = "all" to report all coefficients. The default is which = "nonzero", which is causing the "replacement length differs" issue.

The following works:

library(penalized)
set.seed(1)
x1 = rnorm(100)*10
x2 = rnorm(100)*10
x3 = rnorm(100)*10
y = sin(x1) + cos(x2) - x3 + rnorm(100)

data <- data.frame(y, x1, x2, x3)

win <- 10

coefs <- matrix(NA, ncol=4, nrow=length(y))

for(i in 1:(length(y)-win)) {
  d <- data[(1+i):(win+i),]
  p <- win + i
  pen <- penalized(y, ~ x1 + x2 + x3, ~1, lambda1 = 0, lambda2 = 0,
                   positive = c(F, F, T), data = data)
  beta <- coef(pen, which = "all")
  coefs[p,] <- unname(beta)
  }

How to export all coefficients of a penlized regression model from package `penalized`? Need them for reporting rolling regression estimate

1 Answers