小白勿入，小果教你做lasso回归分析

尔云间  一个专门做科研的团队

我们在通过Cox单因素分析，得到数量较多的基因与生存数据相关，将全部的基因纳入后续分析是不合适的，我们需要通过lasso回归进一步的筛选基因，用于后续的模型构建，所以今天小果开始着手进行lasso回归，开始敲代码

 

代码如下： 

 

remotes::install_github("rpkgs/Ipaper")#安装Ipaper包

library("Ipaper")

library('glmnet')

 

 

data=read.table("training_logFC0.585.CSV",sep = ",",header = T,row.names = 1)

 

x_all <- subset(data, select = -c(OS, OS.time)) #基因表达矩阵

 

y_all <- subset(data, select = c(OS, OS.time))#生存数据

cvfit = cv.glmnet(as.matrix(x_all), Surv(y_all$OS.time,y_all$OS),nfold=10,family = "cox")

 

png(filename = "lasso1.png", height = 450, width = 600)

plot(cvfit, las =1)

dev.off()

pdf(file = "lasso1.pdf", height = 5)

plot(cvfit, las =1)

dev.off()

plot(cvfit)

fit <- glmnet(as.matrix(x_all), Surv(y_all$OS.time,y_all$OS),

              family = "cox")

plot(fit, xvar = "norm", label = TRUE)

 

help("glmnet")

png(filename = "lasso2.png", height = 450, width = 600)

plot(fit, xvar = "lambda",label = TRUE, las=1)

 

dev.off()

pdf(file = "lasso2.pdf", height = 5)

plot(fit, xvar = "lambda",label = TRUE, las=1)

dev.off()

 

 

 

coef.min = coef(cvfit, s = "lambda.min")  ## lambda.min & lambda.1se 取一个

 

active.min = which(coef.min@i != 0) ## 找出那些回归系数没有被惩罚为0的

 

lasso_geneids <- coef.min@Dimnames[[1]][coef.min@i+1] ## 提取基因名称

 

write(lasso_geneids, "lasso_genes.csv")

 

 结果展示

lasso1.pdf如下，两条虚竖线间的基因就是我们筛选到的基因

Lasso2.pdf如下，我们关注的基因是随着 -log lambda 的增加，coefficients值，最晚归零的基因

lasso_genes.csv 如下

 

好了，以上就是今天的内容，欢迎大家留言讨论交流！