Lasso回归

#install.packages("glmnet")
#set.seed(123)
library(glmnet)             
setwd("C:\\Users\\76325\\Desktop\\immune_gene")      

data=read.table("limma_expr_control+treat.csv" , header=T, sep=",", check.names=F, row.names=1)
geneRT=read.table("immune_gene.csv" , header=T, sep=",", check.names=F)
data=data[as.vector(geneRT[,1]),]#交集基因的表达量



#lasso回归
#参数family规定了回归模型的类型
# family = "binomial"适用于二元离散因变量，我们这里结局指标是二分类变量，所以使用bionomial
x=as.matrix(t(data))#交集基因的表达量（行名样本名，列名基因名）
y=c(rep("0",15),rep("1",23))#样本类型

#对数据集进行交叉验证
cvfit=cv.glmnet(x, y, family="binomial", alpha=1,type.measure='mse',nfolds = 10)
plot(cvfit)

#lasso回归
fit=glmnet(x, y, family = "binomial", nlambda = 1000,alpha = 1)#alpha=1为lasso回归，alpha=0是岭回归
coef(fit,s=c(cvfit$lambda.min,cvfit$lambda.1se))

#输出疾病的特征基因
coef=coef(fit, s = cvfit$lambda.min)
index=which(coef != 0)
lassoGene=row.names(coef)[index]
lassoGene=lassoGene[-1]
write.table(lassoGene, file="lasso.gene.txt", sep="\t", quote=F, row.names=F, col.names=F)

#输出疾病特征基因的表达数据
outTab=data[lassoGene,]
write.table(outTab, file="lasso.geneExp.txt", sep="\t", quote=F, row.names=F)