#First install the glmnet package
library("glmnet")
set.seed(1)
n=100
X = matrix(rnorm(n*p), ncol = p)
y = X[,1] + rnorm(n) # only the first variable should matter

p=5
out = glmnet(X, y, int = F)

# Look at the output:
names(out)
par(mfrow = c(2,2))
plot(out, xvar = c("norm"))
plot(out, xvar = c("lambda"))
plot(out, xvar = c("dev"))
out$beta


# Now do it with p = 150
p=150
X = matrix(rnorm(n*p), ncol = p)
y = X[,1] + rnorm(n) # only the first variable should matter

out = glmnet(X, y, int = F)

# Look at the output:
out
par(mfrow = c(2,2))
plot(out, xvar = c("norm"))
plot(out, xvar = c("lambda"))
plot(out, xvar = c("dev"))

out$beta
theta = out$beta[,ncol(out$beta)]
sum(abs(theta))
win.graph()
par(mfrow = c(2,1))
plot(out)
plot(theta)