<p><div><div><div><div data-show-faces="false" data-height="21" data-width="100" data-layout="button_count" data-send="true" data-href="https://www.r-bloggers.com/a-deep-dive-into-glmnet-penalty-factor/"></div></div><div></div></div><div></div></div>
<div>
(This article was first published on <strong><a href="https://statisticaloddsandends.wordpress.com/2018/11/13/a-deep-dive-into-glmnet-penalty-factor/"> R – Statistical Odds &amp; Ends</a></strong>, and kindly contributed to <a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers)</a> </div>
<p>The <code>glmnet</code> function (from the package of the same name) is probably the most used function for fitting the elastic net model in R. (It also fits the lasso and ridge regression, since they are special cases of elastic net.) The <code>glmnet</code> function is very powerful and has several function options that users may not know about. In a series of posts, I hope to shed some light on what these options do.</p>
<p>Here is the full signature of the <code>glmnet</code> function:</p>
<pre>
glmnet(x, y, family=c(&quot;gaussian&quot;,&quot;binomial&quot;,&quot;poisson&quot;,&quot;multinomial&quot;,&quot;cox&quot;,&quot;mgaussian&quot;),
    weights, offset=NULL, alpha = 1, nlambda = 100,
    lambda.min.ratio = ifelse(nobs&lt;nvars,0.01,0.0001), lambda=NULL,
    standardize = TRUE, intercept=TRUE, thresh = 1e-07,  dfmax = nvars + 1,
    pmax = min(dfmax * 2+20, nvars), exclude, penalty.factor = rep(1, nvars),
    lower.limits=-Inf, upper.limits=Inf, maxit=100000,
    type.gaussian=ifelse(nvars&lt;500,&quot;covariance&quot;,&quot;naive&quot;),
    type.logistic=c(&quot;Newton&quot;,&quot;modified.Newton&quot;),
    standardize.response=FALSE, type.multinomial=c(&quot;ungrouped&quot;,&quot;grouped&quot;))
</pre>
<p>In this post, we will focus on the <strong>penalty.factor</strong> option.</p>
<p>Unless otherwise stated, <img alt="n" title="n" src="https://s0.wp.com/latex.php?latex=n&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> will denote the number of observations, <img alt="p" title="p" src="https://s0.wp.com/latex.php?latex=p&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> will denote the number of features, and <code>fit</code> will denote the output/result of the <code>glmnet</code> call.</p>
<p><span><strong>penalty.factor</strong></span></p>
<p>When this option is not set, for each value of <img alt="lambda" title="lambda" src="https://s0.wp.com/latex.php?latex=%5Clambda&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> in <code>lambda</code>, <code>glmnet</code> is minimizing the following objective function:</p>
<p><img alt="begin{aligned} underset{beta}{text{minimize}} quad frac{1}{2}frac{text{RSS}}{n} + lambda displaystylesum_{j=1}^p left(frac{1 - alpha}{2}|beta_j |_2^2 + alpha |beta_j |_1 right). end{aligned}" title="begin{aligned} underset{beta}{text{minimize}} quad frac{1}{2}frac{text{RSS}}{n} + lambda displaystylesum_{j=1}^p left(frac{1 - alpha}{2}|beta_j |_2^2 + alpha |beta_j |_1 right). end{aligned}" src="https://s0.wp.com/latex.php?latex=%5Cbegin%7Baligned%7D+%5Cunderset%7B%5Cbeta%7D%7B%5Ctext%7Bminimize%7D%7D+%5Cquad+%5Cfrac%7B1%7D%7B2%7D%5Cfrac%7B%5Ctext%7BRSS%7D%7D%7Bn%7D+%2B+%5Clambda+%5Cdisplaystyle%5Csum_%7Bj%3D1%7D%5Ep+%5Cleft%28%5Cfrac%7B1+-+%5Calpha%7D%7B2%7D%5C%7C%5Cbeta_j+%5C%7C_2%5E2+%2B+%5Calpha+%5C%7C%5Cbeta_j+%5C%7C_1+%5Cright%29.+%5Cend%7Baligned%7D&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"></p>
<p>When the option is set to a vector <code>c(c_1, ..., c_p)</code>, <code>glmnet</code> minimizes the following objective instead:</p>
<p><img alt="begin{aligned} underset{beta}{text{minimize}} quad frac{1}{2}frac{text{RSS}}{n} + lambda displaystylesum_{j=1}^p c_j left(frac{1 - alpha}{2}|beta_j |_2^2 + alpha |beta_j |_1 right). end{aligned}" title="begin{aligned} underset{beta}{text{minimize}} quad frac{1}{2}frac{text{RSS}}{n} + lambda displaystylesum_{j=1}^p c_j left(frac{1 - alpha}{2}|beta_j |_2^2 + alpha |beta_j |_1 right). end{aligned}" src="https://s0.wp.com/latex.php?latex=%5Cbegin%7Baligned%7D+%5Cunderset%7B%5Cbeta%7D%7B%5Ctext%7Bminimize%7D%7D+%5Cquad+%5Cfrac%7B1%7D%7B2%7D%5Cfrac%7B%5Ctext%7BRSS%7D%7D%7Bn%7D+%2B+%5Clambda+%5Cdisplaystyle%5Csum_%7Bj%3D1%7D%5Ep+c_j+%5Cleft%28%5Cfrac%7B1+-+%5Calpha%7D%7B2%7D%5C%7C%5Cbeta_j+%5C%7C_2%5E2+%2B+%5Calpha+%5C%7C%5Cbeta_j+%5C%7C_1+%5Cright%29.+%5Cend%7Baligned%7D&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"></p>
<p>In the documentation, it is stated that “the penalty factors are internally rescaled to sum to nvars and the lambda sequence will reflect this change.” However, from my own experiments, it seems that the penalty factors are internally rescaled to sum to nvars <strong>but</strong> the lambda sequence remains the same. Let’s generate some data:</p>
<pre>
n &lt;- 100; p &lt;- 5; true_p &lt;- 2
set.seed(11)
X &lt;- matrix(rnorm(n * p), nrow = n)
beta &lt;- matrix(c(rep(1, true_p), rep(0, p - true_p)), ncol = 1)
y &lt;- X %*% beta + 3 * rnorm(n)
</pre>
<p>We fit two models, <code>fit</code> which uses the default options for <code>glmnet</code>, and <code>fit2</code> which has <code>penalty.factor = rep(2, 5)</code>:</p>
<pre>
fit &lt;- glmnet(X, y)
fit2 &lt;- glmnet(X, y, penalty.factor = rep(2, 5))
</pre>
<p>What we find is that these two models have the exact same <code>lambda</code> sequence and produce the same <code>beta</code> coefficients.</p>
<pre>
sum(fit$lambda != fit2$lambda)
# [1] 0
sum(fit$beta != fit2$beta)
# [1] 0
</pre>
<p>The same thing happens when we supply our own <code>lambda</code> sequence:</p>
<pre>
fit3 &lt;- glmnet(X, y, lambda = c(1, 0.1, 0.01), penalty.factor = rep(10, 5))
fit4 &lt;- glmnet(X, y, lambda = c(1, 0.1, 0.01), penalty.factor = rep(1, 5))
sum(fit3$lambda != fit4$lambda)
# [1] 0
sum(fit3$beta != fit4$beta)
# [1] 0
</pre>
<p>Hence, my conclusion is that if <code>penalty.factor</code> is set to <code>c(c_1, ..., c_p)</code>, <code>glmnet</code> is really minimizing</p>
<p><img alt="begin{aligned} underset{beta}{text{minimize}} quad frac{1}{2}frac{text{RSS}}{n} + lambda displaystylesum_{j=1}^p frac{c_j}{bar{c}} left(frac{1 - alpha}{2}|beta_j |_2^2 + alpha |beta_j |_1 right), end{aligned}" title="begin{aligned} underset{beta}{text{minimize}} quad frac{1}{2}frac{text{RSS}}{n} + lambda displaystylesum_{j=1}^p frac{c_j}{bar{c}} left(frac{1 - alpha}{2}|beta_j |_2^2 + alpha |beta_j |_1 right), end{aligned}" src="https://s0.wp.com/latex.php?latex=%5Cbegin%7Baligned%7D+%5Cunderset%7B%5Cbeta%7D%7B%5Ctext%7Bminimize%7D%7D+%5Cquad+%5Cfrac%7B1%7D%7B2%7D%5Cfrac%7B%5Ctext%7BRSS%7D%7D%7Bn%7D+%2B+%5Clambda+%5Cdisplaystyle%5Csum_%7Bj%3D1%7D%5Ep+%5Cfrac%7Bc_j%7D%7B%5Cbar%7Bc%7D%7D+%5Cleft%28%5Cfrac%7B1+-+%5Calpha%7D%7B2%7D%5C%7C%5Cbeta_j+%5C%7C_2%5E2+%2B+%5Calpha+%5C%7C%5Cbeta_j+%5C%7C_1+%5Cright%29%2C+%5Cend%7Baligned%7D&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"></p>
<p>where <img alt="bar{c} = frac{1}{p}sum_{j=1}^p c_j" title="bar{c} = frac{1}{p}sum_{j=1}^p c_j" src="https://s0.wp.com/latex.php?latex=%5Cbar%7Bc%7D+%3D+%5Cfrac%7B1%7D%7Bp%7D%5Csum_%7Bj%3D1%7D%5Ep+c_j&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex">.</p>
	<p>
<div>
<div>To <strong>leave a comment</strong> for the author, please follow the link and comment on their blog: <strong><a href="https://statisticaloddsandends.wordpress.com/2018/11/13/a-deep-dive-into-glmnet-penalty-factor/"> R – Statistical Odds &amp; Ends</a></strong>.</div>
<hr>
<a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers.com</a> offers <strong><a rel="nofollow" href="https://feedburner.google.com/fb/a/mailverify?uri=RBloggers">daily e-mail updates</a></strong> about <a rel="nofollow" href="https://www.r-project.org/" title="The R Project for Statistical Computing">R</a> news and <a rel="nofollow" href="https://www.r-bloggers.com/search/tutorial" title="R tutorials">tutorials</a> on topics such as: <a rel="nofollow" href="https://www.r-bloggers.com/search/data%20science" title="Data science">Data science</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/Big%20Data" title="Big Data">Big Data, <a rel="nofollow" href="https://www.r-users.com/" title="R jobs">R jobs</a>, visualization (<a rel="nofollow" href="https://www.r-bloggers.com/search/ggplot2" title="ggplot and ggplot2 tutorials">ggplot2</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/boxplot" title="Boxplots using lattice and ggplot2 tutorials">Boxplots</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/map" title="Maps and gis">maps</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/animation" title="Animation in R">animation</a>), programming (<a rel="nofollow" href="https://www.r-bloggers.com/search/RStudio" title="RStudio IDE for R">RStudio</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/sweave" title="Sweave and literate programming">Sweave</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/LaTeX" title="LaTeX in R">LaTeX</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/SQL" title="SQL and databases">SQL</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/eclipse" title="Eclipse IDE for R">Eclipse</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/git" title="git and github, Version Control System">git</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/hadoop" title="Large data in R using Hadoop">hadoop</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/Web+Scraping" title="Web Scraping of google, facebook, yahoo, twitter and more using R">Web Scraping</a>) statistics (<a rel="nofollow" href="https://www.r-bloggers.com/search/regression" title="Regressions and ANOVA analysis tutorials">regression</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/PCA" title="principal component analysis tutorial">PCA</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/time+series" title="Time series">time series</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/trading" title="finance trading">trading</a>) and more...
<img width="1" alt="" src="http://feeds.feedburner.com/~r/RBloggers/~4/TVDpdvRl0nw" height="1"></a></div>