生存随机、随机森林变量重要性R包大评测
随机森林（Random Forest，RF）是一种集成学习算法，可以应用于分类和回归问题，ML文章中屡见不鲜。RF通过构建多个决策树并将它们的预测结果进行集成，提高整体模型robust。随着ML领域的发展，目前R包RF众多，包括不限于以下9个包下面随着我一起全面总结一下RF重要性与预测模型可解释。本文相对重要性为全局可解释性，仅抛砖引玉了局部解释，感兴趣ML可解释性问题的同学可以深入学习A课程。
一、9大R包
比较重要的R包是ranger、randomForest、randomForestSRC、cforest。
二、变量重要性VI
随机森林中，特征重要性通常是通过计算每个特征在所有决策树中的平均不纯度减少（分类y）或平均误差减少（回归y）来度量的。不纯度可以是基尼不纯度（Gini impurity）或信息熵（information entropy），用于衡量数据集的混乱程度。回归问题使用平均误差减少，包括平均绝对误差减少量（mean decrease in accuracy）或称为平均绝对误差（mean decrease in MSE）。
三、y01随机森林R代码
y01二分类非生存数据，RF手到擒来。
setwd('D:/聂个人文件/../A款2024课件/17.机器学习/分享rar')
#######训练集和测试集划分
##data = read.csv("data.csv",check.names = F,encoding = "GBK")
# 将id列移动到数据的第一列
##sampled_data = sampled_data[, c("id", setdiff(names(sampled_data), "id"))]
# 保存多个对象
##save(sampled_data, file = "seerdata.RData")
load("seerdata.RData")
##remotes::install_github("zheng-gt/deR")
pacman::p_load(rms, DynNom, regplot,dplyr,deR,tidyverse,survival,survminer,SmartEDA)
#----- 1???? RF随机森林与VI重要性------
##--ranger包5种重要性----
# 计算基于不纯度减少（impurity）的变量重要性
r1 <- ranger::ranger(status ~ ., data = rt, importance = "impurity", num.trees = 100)
ranger.impurity <- ranger::importance(r1)
# 备注: 计算基于 Gini 指标（分类任务）或方差减少（回归任务）的变量重要性。计算效率高，但容易受变量尺度影响，可能对连续变量偏向。
# 计算基于修正后不纯度减少（impurity_corrected）的变量重要性
r1 <- ranger::ranger(status ~ ., data = rt, importance = "impurity_corrected", num.trees = 100)
ranger.impurity_corrected <- ranger::importance(r1)
# 备注: 在原始 impurity 重要性计算的基础上进行了 bias correction（偏差校正），减少对高维数据中特征数目的偏倚影响。
# 计算 Janitza 修正的变量重要性 p 值
ranger.janitza <- ranger::importance_pvalues(r1, method = "janitza")
# 备注: Janitza 方法基于 Bias-Corrected Permutation Importance，适用于分类问题（分类变量重要性计算）。该方法修正了随机森林默认方法对变量重要性高估的问题。
# 计算基于置换误差增加（permutation）的变量重要性
r1 <- ranger::ranger(status ~ ., data = rt, importance = "permutation", num.trees = 100)
ranger.permutation <- ranger::importance(r1)
# 备注: 通过随机打乱单个变量并观察模型性能下降程度来计算变量重要性，适用于回归和分类任务。相比于 impurity 方法，该方法更鲁棒但计算开销较大。
# 计算 Altman 修正的变量重要性 p 值
ranger.altman <- ranger::importance_pvalues(r1, method = "altman", formula = status ~ ., data = rt)
# 备注: Altman 方法适用于回归和分类问题，可调整变量重要性的统计显著性，特别适用于共线性变量评估。
##----?? grf包------
# 确保 status 变量是01因子型
rt <- rt %>% mutate(status = as.factor(status))
# 使用 model.matrix() 进行 One-Hot 编码（去掉第一列截距项）
gfr.matrix <- model.matrix(~ . -1, data = rt[, -c(1,2)])
p.forest <- grf::probability_forest(gfr.matrix, rt$status,num.trees = 100)
grf.V = grf::variable_importance(p.forest,decay.exponent = 1.5,max.depth = 8)
##----randomForest包------
#MeanDecreaseAccuracy MeanDecreaseGini;MDA MDI
set.seed(2025)
# 确保 status 变量是01因子型!
rt <- rt %>% mutate(status = as.factor(status))
rf_model <- randomForest::randomForest(status ~ ., data = rt, ntree = 100, importance = TRUE)
rF <- randomForest::importance(rf_model)
rF
#保存随机森林结果
## write.csv(rF ,"rF随机森林结果.csv")
#---?? y01 基于基尼指数的特征重要性
importance_gini <- importance(rf_model, type = 2)
#--y连续，基于准确度的特征重要性
importance_accuracy <- importance(rf_model, type = 1)
# 获取特征重要性并Gini排序
importance_scores <- rF
feature_importance <- data.frame(var = rownames(importance_scores), Importance = importance_scores[, "MeanDecreaseGini"])
#--VI降序
ordered_features <- feature_importance[order(-feature_importance$Importance), ]
# 选出重要性排名前15的基因
top_var15 <- head(ordered_features, 15)
##---- type1 蓝色/渐进色VI
#定义蓝色柱func
create_rfplot <- function(rf, type){
imp <- importance(rf, type = type, scale = F)
featureImportance <- data.frame(Feature = row.names(imp), Importance = imp[,1])
p <- ggplot(featureImportance, aes(x = reorder(Feature, Importance), y = Importance)) +
geom_bar(stat = "identity", fill = "#53cfff", width = 0.65) +
coord_flip() +
theme_light(base_size = 20) +
theme(axis.title.x = element_text(size = 15, color = "black"),
axis.title.y = element_blank(),
axis.text.x = element_text(size = 15, color = "black"),
axis.text.y = element_text(size = 15, color = "black"))
return(p)
}
#importance为IncNodePurity，gini系数
create_rfplot(rf_model, type = 2)
##----3 randomForest包------
library(randomForest)
#MeanDecreaseAccuracy MeanDecreaseGini;MDA MDI
set.seed(2025)
# 确保 status 变量是01因子型!
rt <- rt %>% mutate(status = as.factor(status))
rf_model <- randomForest::randomForest(status ~ ., data = rt, ntree = 100, importance = TRUE)
rF <- randomForest::importance(rf_model)
rF
#保存随机森林结果
## write.csv(rF ,"rF随机森林结果.csv")
#---?? y01 基于基尼指数的特征重要性
importance_gini <- importance(rf_model, type = 2)
#--y连续，基于准确度的特征重要性
importance_accuracy <- importance(rf_model, type = 1)
# 获取特征重要性并Gini排序
importance_scores <- rF
feature_importance <- data.frame(var = rownames(importance_scores), Importance = importance_scores[, "MeanDecreaseGini"])
#--VI降序
ordered_features <- feature_importance[order(-feature_importance$Importance), ]
# 选出重要性排名前15的基因
top_var15 <- head(ordered_features, 15)
##---- type1 蓝色/渐进色VI
#定义蓝色柱func
create_rfplot <- function(rf, type){
imp <- importance(rf, type = type, scale = F)
featureImportance <- data.frame(Feature = row.names(imp), Importance = imp[,1])
p <- ggplot(featureImportance, aes(x = reorder(Feature, Importance), y = Importance)) +
geom_bar(stat = "identity", fill = "#53cfff", width = 0.65) +
coord_flip() +
theme_light(base_size = 20) +
theme(axis.title.x = element_text(size = 15, color = "black"),
axis.title.y = element_blank(),
axis.text.x = element_text(size = 15, color = "black"),
axis.text.y = element_text(size = 15, color = "black"))
return(p)
}
#importance为IncNodePurity，gini系数
create_rfplot(rf_model, type = 2)
##---??4 randomForestSRC包-----
#3种重要性"anti", "permute", "random"
#https://www.randomforestsrc.org/reference/rfsrc.html
#---> 自动安装需要>R 4.3.0 ，?? 手动安装旧版
##install.packages("data.tree")
## install.packages("D:/聂个人文件/副业项目/广州智荟数据有限公司/A款2024课件/17.机器学习/手动安装R包/randomForestSRC_3.3.0.tar.gz", repos = NULL, type = "source")
## install.packages("randomForestSRC", repos = "https://cran.us.r-project.org")
pacman::p_load(randomForestSRC, ggplot2, dplyr)
rt$status <- as.numeric(as.character(rt$status))
#---?? 后续可以跟进survex shap！参考智荟可解释.R
rc <- randomForestSRC::rfsrc(status~., data = rt)
# 计算变量重要性（anti: 反对称重要性）
src.anti <- randomForestSRC::vimp(rc, importance = c("anti"))
# - 反对称重要性（Antisymmetric importance）
# - 计算变量重要性时考虑变量之间的交互效应，以减少偏差
# - 适用于特征之间可能具有交互作用的情况
src.random <- randomForestSRC::vimp(rc, importance = c("random"))
# - 随机变量重要性（Random importance）
# - 通过随机变量的影响来调整重要性度量
# - 适用于变量分布较复杂或高维数据场景
src.permute <- randomForestSRC::vimp(rc, importance = c("permute"))
# - 置换重要性（Permutation importance）
# - 通过随机置换每个变量，观察预测误差的变化
# - 是最常用、最稳健的变量重要性计算方法，适用于回归、分类和生存分析任务
##----5 party包cforest-----
library(party)
#条件推断森林（Conditional Inference Forest, cforest）
rc <- party::cforest(status ~., data = rt, control = cforest_unbiased(ntree = 100))
cf.vip = varimp(rc)
##----5 汇总rf包结果----
# 查看每个变量的行数
sapply(list(
ranger.impurity,
ranger.impurity_corrected,
ranger.permutation,
ranger.janitza[,1],
ranger.altman[,1],
grf.V,
rF[,3],
rF[,4],
src.anti[["importance"]][,1],
src.random[["importance"]][,1],
src.permute[["importance"]][,1],
cf.vip
), length)
#不考虑gfr包
imp_score = data.frame(ranger.impurity,
ranger.impurity_corrected,
ranger.permutation,
ranger.janitza = ranger.janitza[,1],
ranger.altman = ranger.altman[,1],
rf.MDA = rF[,3],
rf.MDI = rF[,4],
src.anti = src.anti[["importance"]][,1],
src.random = src.random[["importance"]][,1],
src.permute = src.permute[["importance"]][,1],
cf = cf.vip )
t(which(imp_score<0,arr.ind = TRUE))
pacman::p_load(dplyr,ggplot2,forcats)
# 计算绝对值的重要性得分
imp_score <- abs(imp_score)
# 归一化处理
norm_score <- as.matrix(imp_score) %*% diag(1 / colSums(imp_score))
# 构造数据框
ggdt <- data.frame(
method = rep(colnames(imp_score), each = nrow(imp_score)),
Xi = rep(rownames(imp_score), times = ncol(imp_score)),
score = as.vector(norm_score) )
# 数据排序与分组
ggdt <- ggdt %>%
arrange(method, score) %>%
mutate(group_order = forcats::fct_inorder(interaction(Xi, score)))
# 自定义颜色TheBestColors
TheBestColors::Best100(30)
#custom_colors <- c(
# "dodgerblue2", "#E31A1C", "green4", "#6A3D9A", "#FF7F00",
# "black", "gold1", "skyblue2", "#FB9A99", "palegreen2",
# "#CAB2D6", "#FDBF6F", "gray70", "khaki2", "maroon",
# "orchid1", "deeppink1", "blue1", "steelblue4",
# "darkorange4", "brown"
#)
# 更新颜色方案
custom_colors <- c(
"#795F35", "#FBC1AF", "#EF2C30", "#706D6E", "#F4D49D",
"#FFEE00", "#FBC41E", "#B45241", "#5A69B2", "#1AB3EF",
"#DD9785", "#938ABB", "#FFF4A4", "#86817E", "#D8D7D7",
"#F58971", "#E56D4C", "#F2833F", "#28B45D", "#FDE09F",
"#EE1A95" )
#画图
ggplot(ggdt, aes(x = method, y = score, fill = Xi, group = group_order)) +
geom_bar(position = "stack", stat = "identity") +
scale_fill_manual(values = custom_colors, breaks = rownames(imp_score)) +
coord_flip() +
theme_minimal(base_size = 14) + # 统一美化主题
theme(
axis.text.x = element_text(angle = 45, hjust = 1), # x轴文字倾斜
legend.position = "bottom", # 调整图例位置
legend.title = element_blank(), # 去掉图例标题
panel.grid.major = element_blank(), # 去掉大网格
panel.grid.minor = element_blank() # 去掉小网格
) +
labs(
x = "变量重要性计算方法",
y = "归一化变量重要性得分",
title = "不同方法计算的变量重要性对比"
)四、生存随机森林R代码
对于生存数据需要用RandomForestSRC包处理！RandomForestSRC包可以基于3种相对重要性评估，diy生存重要性图，后续可以与surviz包联合进行ML可解释性的系列分析！十分方便。
#----- 2???? 生存RF随机森林与VI重要性------
##----??randomForestSRC----
#---?? 变成 Surv(time, status==1) ~ .
rc <- randomForestSRC::rfsrc(Surv(time, status==1) ~ ., data = rt,
ntree = 100, block.size = 1)
pacman::p_load(survex)
rsf_exp <- survex::explain(rc)
mp_rsf <- model_performance(rsf_exp)
plot(mp_rsf)
#--三重要性指标
src.anti <- randomForestSRC::vimp(rc,importance = c("anti"))
src.random <- randomForestSRC::vimp(rc,importance = c("random"))
src.permute <- randomForestSRC::vimp(rc,importance = c( "permute"))
imp_score = data.frame(
src.anti = src.anti[["importance"]],
src.random = src.random[["importance"]],
src.permute = src.permute[["importance"]] )
#即使你的 status 只有 0/1，也会以 竞争风险格式处理它；event.1 event.2
#--- ?? 不存在012竞争，仅保留 event.1 变量！！
imp_score <- imp_score[, grepl("event.1", colnames(imp_score))]
# 获取变量数量
num_variables <- nrow(imp_score) # 变量个数
num_methods <- ncol(imp_score) # 方法个数
# 构造数据框，自动适应变量数量
ggdt <- data.frame(
method = rep(colnames(imp_score), each = num_variables), # 根据方法数重复
Xi = rep(rownames(imp_score), times = num_methods), # 根据变量数重复
score = as.vector(norm_score) )
# 加载 dplyr 进行数据排序
library(dplyr)
ggdt <- ggdt %>%
arrange(method, score) %>%
mutate(group_order = forcats::fct_inorder(interaction(Xi, score)))
# 颜色列表（自动扩展）
custom_colors <- c(
"#795F35", "#FBC1AF", "#EF2C30", "#706D6E", "#F4D49D",
"#FFEE00", "#FBC41E", "#B45241", "#5A69B2", "#1AB3EF",
"#DD9785", "#938ABB", "#FFF4A4", "#86817E", "#D8D7D7",
"#F58971", "#E56D4C", "#F2833F", "#28B45D", "#FDE09F",
"#EE1A95")
# 画图，不限制 event.2 变量，自动适应变量数量
library(ggplot2)
ggplot(ggdt, aes(fill = Xi, y = score, x = method, group = group_order)) +
geom_bar(position = "stack", stat = "identity") +
scale_fill_manual(
values = rep(custom_colors, length.out = num_variables), # 颜色自动扩展
breaks = rownames(imp_score)
) +
coord_flip() +
theme_minimal(base_size = 14) +
theme(
axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "bottom",
legend.title = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank()
) +
labs(
x = "变量重要性计算方法",
y = "归一化变量重要性得分",
title = "不同方法计算的变量重要性对比"
)
五、发散思维
1.RF变量显著性p值：rfPermute 随机置换检验（permutation test）来评估变量的重要性，并提供 p 值，以判断某个变量的重要性是否显著高于随机噪声。rfPermute与randomForest包gini重要性结果一致！
2.RF的包不同，算法不同，重要性VI不同。选什么是个问题啊！seed不同、分区不同、model不同、VI不同、shap不同结果大不相同！ML真是很玄。
六、小结
“随机森林R包众多，推荐randomForestSRC 支持生存数据，可后续可解释性。
本公众号建立了学术交流群，仅供SCI学术交流，人数有限需要实名制。入群请加笔者微信 popnie，请备注说明：姓名-学校(单位)-专业。