ixxmu
commented
4 months ago 101 种机器学习算法组合筛选最优预后风险模型（Mime 包） by 生信碱移
生信碱移
▲ DOI: 10.1016/j.csbj.2024.06.035
github 链接如下：
https://github.com/l-magnificence/Mime
一、安装 R 包
if (!requireNamespace("BiocManager", quietly = TRUE)) install.packages("BiocManager")
depens<-c('GSEABase', 'GSVA', 'cancerclass', 'mixOmics', 'sparrow', 'sva' , 'ComplexHeatmap' )
for(i in 1:length(depens)){
depen<-depens[i]
if (!requireNamespace(depen, quietly = TRUE)) BiocManager::install(depen,update = FALSE)
}
if (!requireNamespace("CoxBoost", quietly = TRUE))
devtools::install_github("binderh/CoxBoost")
if (!requireNamespace("fastAdaboost", quietly = TRUE))
devtools::install_github("souravc83/fastAdaboost")
if (!requireNamespace("Mime", quietly = TRUE))
devtools::install_github("mumdark/MimeDL")
Mime包似乎与devtools的依赖包mime名称一致，导致安装失败，可能需要改一下R包名称（大家遇到了再说吧）。二、简单的示例
① 加载示例数据（这个治疗反应数据在本文用不上，在官网是用作二分类的的示例，感兴趣的铁子可以自行看看）。：
# 01 加载示例数据
load("./External data/Example.cohort.Rdata") # 生存数据与基因表达信息
list_train_vali_Data[["Dataset1"]][1:5,1:5]
# ID OS.time OS MT-CO1 MT-CO3
#60 TCGA.DH.A66B.01 1281.65322 0 13.77340 13.67931
#234 TCGA.HT.7607.01 96.19915 1 14.96535 14.31857
#42 TCGA.DB.A64Q.01 182.37755 0 13.90659 13.65321
#126 TCGA.DU.8167.01 471.97707 0 14.90695 14.59776
#237 TCGA.HT.7610.01 1709.53901 0 15.22784 14.62756
load("./External data/Example.ici.Rdata") # 治疗反应数据
list_train_vali_Data[["training"]][1:5,1:5]
# ID Var FTH1 EEF1A1 ACTB
#Mariathasan_UC_pre_aPDL1_combo_tpm.1 SAMf2ce197162ce N 10.114846 4.817746 11.230180
#Gide_SKCM_pre_combo.66 ERR2208915 Y 2.044180 5.038854 3.977902
#Bruan_RCC_pre_aPD1_combo_tpm.3 G138701_RCCBMS-00141-T_v1_RNA_OnPrem Y 5.406008 5.341635 5.366668
#Mariathasan_UC_pre_aPDL1_combo_tpm.203 SAMe41b1e773582 N 9.215794 4.707360 11.412721
#Mariathasan_UC_pre_aPDL1_combo_tpm.58 SAM5ffd7e4cd794 N 9.003710 3.908884 10.440559
load("./External data/genelist.Rdata") # 基因列表
genelist
# [1] "MYC" "CTNNB1" "JAG2" "NOTCH1" "DLL1" "AXIN2" "PSEN2" "FZD1" "NOTCH4" "LEF1" "AXIN1" "NKD1" "WNT5B" "CUL1"
#[15] "JAG1" "MAML1" "KAT2A" "GNAI1" "WNT6" "PTCH1" "NCOR2" "DKK4" "HDAC2" "DKK1" "TCF7" "WNT1" "NUMB" "ADAM17"
#[29] "DVL2" "PPARD" "NCSTN" "HDAC5" "CCND2" "FRAT1" "CSNK1E" "RBPJ" "FZD8" "TP53" "SKP2" "HEY2" "HEY1" "HDAC11"
list_train_vali_Data与genelist，介绍如下：list_train_vali_Data：行是样本，第一列是样本ID，第二列是生存时间（单位year），第三列是生存状态，随后的列都是相应基因的表达量； genelist：是用于建模的基因列表。
② 构建预后模型：
load("./External data/Example.cohort.Rdata")
load("./External data/genelist.Rdata")
res <- ML.Dev.Prog.Sig(train_data = list_train_vali_Data$Dataset1,
list_train_vali_Data = list_train_vali_Data,
unicox.filter.for.candi = T,
unicox_p_cutoff = 0.05,
candidate_genes = genelist,
mode = 'all',
nodesize =5,
seed = 123)
③ cindex 评估预后模型：
cindex_dis_all(res,
validate_set = "Dataset2", # 指定测试集的名称
order =names(list_train_vali_Data),
width = 0.35)
cindex_dis_select(res,
model="StepCox[forward] + plsRcox", # 指定选择模型
order= names(list_train_vali_Data))
④ ROC下面积AUC 评估预后模型（下面展示的是1年的 AUC，通过AUC_time参数指定）：
all.auc.1y <- cal_AUC_ml_res(res.by.ML.Dev.Prog.Sig = res,
train_data = list_train_vali_Data[["Dataset1"]],
inputmatrix.list = list_train_vali_Data,
mode = 'all',
AUC_time = 1, # 一年生存
auc_cal_method="KM")
auc_dis_all(all.auc.1y,
dataset = names(list_train_vali_Data),
validate_set="Dataset1",
order= names(list_train_vali_Data),
width = 0.35,
year=1)
model_name参数指定模型。下面给大家介绍几个评估方式。首先，可以看看指定模型的ROC下面积AUC（下面展示的是1年的 AUC，通过AUC_time参数指定）：all.auc.1y <- cal_AUC_ml_res(res.by.ML.Dev.Prog.Sig = res,
train_data = list_train_vali_Data[["Dataset1"]],
inputmatrix.list = list_train_vali_Data,
mode = 'all',
AUC_time = 1, # 一年生存
auc_cal_method="KM")
auc_dis_all(all.auc.1y,
dataset = names(list_train_vali_Data),
validate_set="Dataset1",
order= names(list_train_vali_Data),
width = 0.35,
year=1)
roc_vis(all.auc.1y,
model_name = "StepCox[forward] + plsRcox",
dataset = names(list_train_vali_Data),
order= names(list_train_vali_Data),
anno_position=c(0.65,0.55),
year=1)
all.auc.3y <- cal_AUC_ml_res(res.by.ML.Dev.Prog.Sig = res,train_data = list_train_vali_Data[["Dataset1"]],
inputmatrix.list = list_train_vali_Data,mode = 'all',AUC_time = 3,
auc_cal_method="KM")
all.auc.5y <- cal_AUC_ml_res(res.by.ML.Dev.Prog.Sig = res,train_data = list_train_vali_Data[["Dataset1"]],
inputmatrix.list = list_train_vali_Data,mode = 'all',AUC_time = 5,
auc_cal_method="KM")
auc_dis_select(list(all.auc.1y,all.auc.3y,all.auc.5y),
model_name="StepCox[forward] + plsRcox",
dataset = names(list_train_vali_Data),
order= names(list_train_vali_Data),
year=c(1,3,5))
rs_sur(res,
model_name = "StepCox[forward] + plsRcox", # 选择你的模型
dataset = "Dataset1",
#color=c("blue","green"),
median.line = "hv",
cutoff = 0.5,
conf.int = T,
xlab="Day",pval.coord=c(1000,0.9))
unicox.rs.res <- cal_unicox_ml_res(res.by.ML.Dev.Prog.Sig = res,
optimal.model = "StepCox[forward] + plsRcox",
type ='categorical')
metamodel <- cal_unicox_meta_ml_res(input = unicox.rs.res)
meta_unicox_vis(metamodel,
dataset = names(list_train_vali_Data))
⑥ 与其他已发表文章的模型比较（作者仅仅提供了胶质瘤发表文章的模型）:
首先是Cindex的比较：
cc.glioma.lgg.gbm <- cal_cindex_pre.prog.sig(use_your_own_collected_sig = F,
type.sig = c('Glioma','LGG','GBM'),
list_input_data = list_train_vali_Data)
cindex_comp(cc.glioma.lgg.gbm,
res,
model_name="StepCox[forward] + plsRcox",
dataset=names(list_train_vali_Data))
然后是 AUC的比较：
auc.glioma.lgg.gbm.1 <- cal_auc_pre.prog.sig(use_your_own_collected_sig = F,
type.sig = c('Glioma','LGG','GBM'),
list_input_data = list_train_vali_Data,
AUC_time = 1,
auc_cal_method = 'KM')
auc_comp(auc.glioma.lgg.gbm.1,
all.auc.1y,
model_name="StepCox[forward] + plsRcox",
dataset=names(list_train_vali_Data))
这个r包功能还有其他功能，包括其他潜在的参数大家都可以自行看看文档。
欢迎点击关注
仅供粉丝老铁们参考
如有侵权或错误，请联系删除改正～
ZimLea
https://mp.weixin.qq.com/s/B_Yvb1jqAI1H2Up411X2Sw