Performing ridge regression using `optim`. Not sure where this extra term comes from?

Question

Doing problems from ISLR and I've taken up the task of trying to do linear regression (and by extension lasso and ridge regression) using R's optim function. Recall equation 6.5 of page 215, which states that ridge regression is the minimization of:

$$ \sum_{i = 1}^n \left( y_i - \beta_0 - \sum_{j = 1}^p \beta_j x_{ij} \right)^2 + \lambda \sum_{j = 1}^p \beta_j^2 $$

I created a minimization function based on this. The coefficients it generates do not match up very well with those given by glmnet.

In frustration, I continued searching and I found an existing working implementation:

https://github.com/m-clark/Miscellaneous-R-Code/blob/master/ModelFitting/ridge.R

The part that confuses me is the minimization function:

ridge <- function(w, X, y, lambda = .1) {
  # X: model matrix
  # y: target; lambda: penalty parameter
  # w: the weights/coefficients

  crossprod(y - X%*%w) + lambda*length(y)*crossprod(w)
}

This implementation works as expected, but notice the length(y) term used in the penalty. How did that get there? It is not present in the original equation!

Am I missing something in the original equation? Is this extra term necessary when using optim? If so, why?

============== UPDATE ==============

Here is the simulated_data code:

simulated_data <- function(seed, sample_size, predictors) {
  set.seed(seed)

  X <- scale(matrix(rnorm(sample_size*predictors), ncol = predictors))
  b <- rep(1, predictors)
  y <- scale(X %*% b + rnorm(sample_size))

  return(list(X = X, y = y))
}

Here is the ridge regression code I've written which uses glmnet:

tuple <- simulated_data(1, 100, 3)

# Load the glmnet library
library(glmnet)

# Build the ridge regression model.
# Note that alpha = 0 is for regression while alpha = 1 is for lasso.
ridge.model <- glmnet(tuple$X, tuple$y, 
                      alpha = 0, 
                      lambda = c(10, 1, .1),
                      thresh = 1e-12,
                      intercept = F)
ridge.cf <- coef(ridge.model, s = 0.1)
print(ridge.cf)

This outputs:

## 4 x 1 sparse Matrix of class "dgCMatrix"
##                     1
## (Intercept) .        
## V1          0.4357297
## V2          0.4322341
## V3          0.5070639

And here is the optim code for ridge regression:

minimize.ridge <- function(par, X, y, lambda = 0.1) {
  rss <- crossprod((X %*% par) - y)
  #penalty <- lambda*length(y)*crossprod(par)
  penalty <- lambda * (par %*% par)
  return(rss + penalty)
}

op.result <- optim(rep(0, 3), 
                   fn = minimize.ridge,
                   method = 'BFGS',
                   X = tuple$X,
               y = tuple$y,
                   lambda = 0.1)
op.cf <- op.result$par
print(op.cf)

This outputs

## [1] 0.4781309 0.4778045 0.5589391

Answer 1

I continued to dig and I found the glmnet paper. As it turns out, this version of the minimization function is not for the vanilla ridge regression that is stated in ISLR but for the elastic-net version of ridge regression that is implemented.

Kudos to Kornel for catching that length(y) is actually $n$ not $p$.