if(!require(glmnet)){install.packages("glmnet"); library(glmnet)}
ifname<-"DataAgonist.csv"
dataAll <- read.delim(ifname, row.names=1,header=T, sep= ",", as.is=TRUE, 
                      strip.white=FALSE)
nf<-dim(dataAll)[[2]]
xmatrix<-as.matrix(dataAll[,2:nf])
yvector<-dataAll[,1]
#[1]- cross validation
alphaInput<-1
results.cvlasso<-cv.glmnet(x = xmatrix, y = yvector,alpha = alphaInput,
                           nlambda=50,nfolds=5)
#alpha=1 is the lasso penalty, and is alpha=0 the ridge penalty.
names(results.cvlasso)
plot(results.cvlasso$lambda,results.cvlasso$cvm,type="l",
     main="Lambda vs. mean cross-validated error")
bestlambda<-results.cvlasso$lambda.min; bestlambda
plot(results.cvlasso)
#[2]- beta validation
result.lasso<- glmnet(x = xmatrix, y = yvector, alpha = alphaInput,nlambda=50)
names(result.lasso)
ramda.value<-log(result.lasso$lambda); matplot(ramda.value,t(result.lasso$beta),pch=1:nf,col=1,main="Lambda vs. beta")
par(ps = 10,mgp=c(1, 0.1, 0.1)); par(mai=c(.8,.8,.8,1.5),xpd=TRUE)
featureName<-rownames(result.lasso$beta)
nf<-length(featureName)
legend(legend=featureName,pch=1:nf,col=1,"bottomright",ncol=2)
#[3]---------------------------------------------------------------------------------------
Best.lasso<-glmnet(x=xmatrix,y=yvector,alpha=alphaInput,lambda=bestlambda)
Best.lasso$beta
est.Y <- predict(Best.lasso, newx = xmatrix); est.Y
s <- paste("Lasso: r = ",sprintf("%5.3f",cor(yvector,est.Y)))
plot(yvector, est.Y,main=s)
cor(yvector,est.Y)
mse <- sum((yvector - est.Y)^2) / length(yvector); mse