【二分类因变量机器学习】R语言08-lightGBM机器学习模型（标准化代码）

# 设置工作目录和清理环境
rm(list = ls())
if (!is.null(dev.list())) dev.off()
setwd("C:/Users/hyy/Desktop/")

# 创建结果文件夹
if (!dir.exists("Results-gbm")) dir.create("Results-gbm")

# 加载必要的包
if (!require(pacman)) install.packages("pacman")
pacman::p_load(readxl, writexl, lightgbm, pROC, ggplot2, DALEX, DALEXtra, 
               officer, flextable, shapviz, DiagrammeR, data.table, caret)

# 设置全局中文字体支持
if (.Platform$OS.type == "windows") {
  windowsFonts(SimHei = windowsFont("SimHei"))
  windowsFonts(SimSun = windowsFont("SimSun"))
  chinese_font <- "SimHei"
} else {
  chinese_font <- "sans"
}

# 1. 数据准备
data <- read_excel("示例数据.xlsx")

# 检查数据结构
cat("数据结构:\n")
str(data)
cat("\n结局变量分布:\n")
table(data$结局)

# 将结局转换为因子
data$结局 <- as.factor(data$结局)

# 分层抽样划分训练集和测试集
set.seed(123)
train_index <- createDataPartition(data$结局, p = 0.7, list = FALSE)
Train <- data[train_index, ]
Test <- data[-train_index, ]

# 2. 数据预处理 - 确保所有特征都是数值型
# 检查并转换非数值型列
feature_cols <- setdiff(names(data), c("序号", "结局"))
non_numeric_cols <- sapply(data[, feature_cols], function(x) !is.numeric(x))

if (any(non_numeric_cols)) {
  cat("发现非数值型列，正在转换为数值型...\n")

# 使用独热编码或标签编码转换非数值型列
for (col in names(non_numeric_cols)[non_numeric_cols]) {
    if (is.character(data[[col]]) || is.factor(data[[col]])) {
      # 使用标签编码
      Train[[col]] <- as.numeric(factor(Train[[col]]))
      Test[[col]] <- as.numeric(factor(Test[[col]]))
      cat(paste("已将列", col, "转换为数值型\n"))
    }
  }
}

# 准备lightGBM数据
Train.matrix <- as.matrix(Train[, feature_cols])
Test.matrix <- as.matrix(Test[, feature_cols])
Train.labels <- as.numeric(Train$结局) - 1  # 转换为0/1

# 确保矩阵不包含NA或无限值
Train.matrix[is.na(Train.matrix)] <- 0
Train.matrix[is.infinite(Train.matrix)] <- 0
Test.matrix[is.na(Test.matrix)] <- 0
Test.matrix[is.infinite(Test.matrix)] <- 0

# 创建lightGBM数据集
lgb.Train <- lgb.Dataset(data = Train.matrix, label = Train.labels)

# 3. 设置LightGBM参数

# 5. 模型评价
# 预测测试集
pred_probs <- predict(lgb.model, Test.matrix)
pred_class <- ifelse(pred_probs > 0.5, 1, 0)

# 计算准确率
accuracy <- mean(pred_class == (as.numeric(Test$结局) - 1))
cat("模型准确率:", accuracy, "\n")

# 混淆矩阵
conf_matrix <- table(实际 = as.numeric(Test$结局) - 1, 预测 = pred_class)
conf_matrix_df <- as.data.frame.matrix(conf_matrix)

# ROC曲线
roc_obj <- roc(response = as.numeric(Test$结局) - 1, predictor = pred_probs)
auc_value <- auc(roc_obj)

# 绘制ROC曲线
jpeg("Results-gbm/gbm_roc_curve.jpg", width = 800, height = 600, quality = 100)
plot(roc_obj, main = "LightGBM ROC曲线", col = "blue", print.auc = TRUE)
dev.off()

pdf("Results-gbm/gbm_roc_curve.pdf", width = 8, height = 6)
plot(roc_obj, main = "LightGBM ROC曲线", col = "blue", print.auc = TRUE)
dev.off()

# 6. 特征重要性
importance <- lgb.importance(lgb.model, percentage = TRUE)

# 绘制特征重要性图
jpeg("Results-gbm/feature_importance.jpg", width = 800, height = 600, quality = 100)
lgb.plot.importance(importance, top_n = 10, measure = "Gain")
dev.off()

pdf("Results-gbm/feature_importance.pdf", width = 8, height = 6)
lgb.plot.importance(importance, top_n = 10, measure = "Gain")
dev.off()

# 7. SHAP分析 - 修正后的代码

pdf("Results-gbm/shap_waterfall.pdf", width = 8, height = 6)
sv_waterfall(shap, row_id = 1)
dev.off()

# 8. 保存结果表格
# 模型性能指标
performance_metrics <- data.frame(
  指标 = c("准确率", "AUC", "灵敏度", "特异度"),
  值 = c(round(accuracy, 4),
        round(auc_value, 4),
        round(conf_matrix[2,2]/sum(conf_matrix[2,]), 4),  # 灵敏度
        round(conf_matrix[1,1]/sum(conf_matrix[1,]), 4))  # 特异度
)

# 特征重要性表格
feature_importance_df <- as.data.frame(importance)

write_xlsx(list(
  混淆矩阵 = conf_matrix_df,
  性能指标 = performance_metrics,
  特征重要性 = feature_importance_df
), "Results-gbm/gbm_results.xlsx")

# 9. 生成Word报告
doc <- read_docx()

# 添加标题
doc <- doc %>% 
  body_add_par("LightGBM模型分析报告", style = "heading 1") %>%
  body_add_par(paste0("生成日期: ", Sys.Date()), style = "Normal") %>%
  body_add_par(" ", style = "Normal")

# 添加模型摘要
doc <- doc %>% 
  body_add_par("模型摘要", style = "heading 2") %>%
  body_add_par(paste("准确率:", round(accuracy, 4)), style = "Normal") %>%
  body_add_par(paste("AUC值:", round(auc_value, 4)), style = "Normal") %>%
  body_add_par(paste("最佳迭代次数:", lgb.model$best_iter), style = "Normal") %>%
  body_add_par(" ", style = "Normal")

# 添加图片
doc <- doc %>% 
  body_add_par("ROC曲线", style = "heading 3") %>%
  body_add_img("Results-gbm/gbm_roc_curve.jpg", width = 6, height = 4.5) %>%
  body_add_par(" ", style = "Normal") %>%
  body_add_par("特征重要性", style = "heading 3") %>%
  body_add_img("Results-gbm/feature_importance.jpg", width = 6, height = 4.5) %>%
  body_add_par(" ", style = "Normal") %>%
  body_add_par("SHAP摘要图", style = "heading 3") %>%
  body_add_img("Results-gbm/shap_summary.jpg", width = 6, height = 4.5) %>%
  body_add_par(" ", style = "Normal") %>%
  body_add_par("SHAP Waterfall图", style = "heading 3") %>%
  body_add_img("Results-gbm/shap_waterfall.jpg", width = 6, height = 4.5) %>%
  body_add_par(" ", style = "Normal")

# 添加性能表格
doc <- doc %>% 
  body_add_par("性能指标", style = "heading 3") %>%
  body_add_flextable(performance_metrics %>% 
                       flextable() %>% 
                       theme_box() %>%
                       set_caption("模型性能指标")) %>%
  body_add_par(" ", style = "Normal")

# 添加特征重要性表格
doc <- doc %>% 
  body_add_par("特征重要性", style = "heading 3") %>%
  body_add_flextable(feature_importance_df %>% 
                       head(10) %>%  # 只显示前10个最重要的特征
                       flextable() %>% 
                       theme_box() %>%
                       set_caption("特征重要性排名")) %>%
  body_add_par(" ", style = "Normal")

# 保存Word文档
print(doc, target = "Results-gbm/gbm_analysis_report.docx")

# 10. 保存工作空间
save.image("Results-gbm/gbm_analysis_workspace.RData")

cat("LightGBM分析完成!所有结果已保存到Results-gbm文件夹中。\n")