Ref

16 Dimensionality Reduction | Tidy Modeling with R

解读

1. 数据准备与划分

# 加载必要的包
library(tidymodels)
library(stringr)
library(beans)  # 包含豆类数据集
library(baguette)  # 装袋模型
library(discrim)   # 判别分析模型
library(embed)     # 降维方法
library(corrplot)  # 相关性可视化
library(patchwork) # 图形排版
library(ggforce)   # 高级图形功能
tidymodels_prefer() # 设置tidymodels优先使用函数
 
# 设置随机种子确保结果可重现
set.seed(1601)
 
# 创建训练集/验证集/测试集的三重分割
bean_split <- initial_validation_split(
  beans, 
  strata = class,  # 按类别分层抽样,保持各类别比例
  prop = c(0.75, 0.125)  # 75%训练集,12.5%验证集,12.5%测试集
)
 
bean_split
 
 
# 提取为数据框格式(用于手动分析)
bean_train <- training(bean_split)    # 训练集:10206行
bean_test <- testing(bean_split)      # 测试集:1703行  
bean_validation <- validation(bean_split) # 验证集:1702行
 
# 创建用于tune函数的重抽样对象
set.seed(1602)
bean_val <- validation_set(bean_split)
bean_val$splits[[1]]
#> <Training/Validation/Total>
#> <10206/1702/11908>

为什么用验证集而不是交叉验证?详细技术解释

计算效率考虑:数据集较大(13,611行),10折交叉验证需要训练10个模型,而单验证集只需1个

模型选择策略:

  • 验证集作用:快速筛选有希望的模型和参数组合(粗调)
  • 测试集作用:最终评估选出的最佳模型(精调)
  • 这种”验证集筛选+测试集确认”的策略在大型数据集上更实用

统计可靠性:

  • 虽然交叉验证通常更稳定,但当数据集足够大时,单次验证集也能提供可靠的性能估计
  • 验证集大小(1702行)足以获得稳定的性能评估

2. 数据探索与相关性分析

# 检查数据基本信息
glimpse(bean_train)
summary(bean_train)
 
# 可视化特征间的相关性
tmwr_cols <- colorRampPalette(c("#91CBD765", "#CA225E"))
bean_train %>% 
  select(-class) %>%  # 排除类别变量
  cor() %>%  # 计算相关系数矩阵
  corrplot(col = tmwr_cols(200), tl.col = "black")

3. 特征工程配方设计

# 基础预处理配方
bean_rec <-
  recipe(class ~ ., data = bean_train) %>%
  step_zv(all_numeric_predictors()) %>%        # 移除零方差特征
  step_orderNorm(all_numeric_predictors()) %>% # 有序正态化变换
  step_normalize(all_numeric_predictors())     # 标准化(均值为0,标准差为1)
 
# 在训练集上拟合预处理参数
bean_rec_trained <- prep(bean_rec, training = bean_train)
bean_rec_trained
 
# 查看预处理步骤的详细信息
tidy(bean_rec_trained)  # 显示所有步骤
tidy(bean_rec_trained, number = 2)  # 查看第2步的详细参数
 
# 将预处理应用于验证集
bean_val_processed <- bake(bean_rec_trained, new_data = bean_validation)
 
# 比较原始数据和处理后数据
original_summary <- bean_validation %>% 
  select(where(is.numeric)) %>% 
  summary()
 
processed_summary <- bean_val_processed %>% 
  select(where(is.numeric)) %>% 
  summary()
 
print("原始数据统计摘要:")
print(original_summary)
 
print("处理后数据统计摘要:")  
print(processed_summary)

为什么bake要用训练集的参数?深入解释数据泄露问题

技术原理:

  • prep() 在训练集上计算:均值、标准差、变换参数等
  • bake() 将这些参数应用于新数据

如果错误做法):

# 错误做法:验证集自己计算参数
validation_wrong <- bean_validation %>%
  recipe(class ~ .) %>%
  step_normalize(all_numeric_predictors()) %>%
  prep() %>%  # 在验证集上计算参数!
  bake(new_data = NULL)
 
# 正确做法:使用训练集的参数
validation_correct <- bake(bean_rec_trained, new_data = bean_validation)

为什么这是数据泄露?

  1. 测试集信息污染:验证集/测试集的分布信息”泄露”到预处理中
  2. 过于乐观的评估:模型在训练时”看到”了测试集的信息
  3. 泛化能力高估:在实际部署时,新数据没有对应的统计参数

真实场景类比:

  • 考试前知道了考题(数据泄露) vs. 真正遇到新问题(正确做法)
  • 预处理参数应该从”学习经验”(训练集)得来,而不是”考试内容”(测试集)

4. 数据分布可视化比较

# 比较处理前后的数据分布
p1 <- 
  bean_validation %>% 
  ggplot(aes(x = area)) + 
  geom_histogram(bins = 30, color = "white", fill = "blue", alpha = 1/3) + 
  ggtitle("Original Validation Set Data") +
  labs(x = "Area", y = "Count")
 
p2 <- 
  bean_val_processed %>% 
  ggplot(aes(x = area)) + 
  geom_histogram(bins = 30, color = "white", fill = "red", alpha = 1/3) + 
  ggtitle("Processed Validation Set Data") +
  labs(x = "Standardized Area", y = "Count")
 
# 并排显示
p1 + p2
 
# 多变量分布比较函数
plot_validation_results <- function(recipe, dat = bean_validation, title = "") {
  processed_data <- recipe %>%
    prep() %>%
    bake(new_data = dat)
  
  processed_data %>%
    ggplot(aes(x = .panel_x, y = .panel_y, color = class, fill = class)) +
    geom_point(alpha = 0.4, size = 0.5) +
    geom_autodensity(alpha = .3) +
    facet_matrix(vars(-class), layer.diag = 2) + 
    scale_color_brewer(palette = "Dark2") + 
    scale_fill_brewer(palette = "Dark2") +
    ggtitle(title)
}

5. 降维方法比较

# PCA降维(线性降维)
pca_plot <- bean_rec_trained %>%
  step_pca(all_numeric_predictors(), num_comp = 4) %>%
  plot_validation_results(title = "Principal Component Analysis (PCA)")
 
# UMAP降维(非线性降维)  
umap_plot <- bean_rec_trained %>%
  step_umap(all_numeric_predictors(), num_comp = 4) %>%
  plot_validation_results(title = "UMAP Visualization")
 
pca_plot
umap_plot

**为什么UMAP需要调参?

UMAP关键超参数及其影响:

  1. num_comp(嵌入维度)

    • 太低:信息损失严重,无法有效分离类别
    • 太高:计算成本增加,可能引入噪声
    • 调优范围:通常2-10维
  2. neighbors(近邻数)

    • 控制局部vs全局结构平衡
    • 太小:过度关注局部结构,可能分裂连续簇
    • 太大:忽略局部细节,簇边界模糊
    • 公式:建议值 = sqrt(样本数) ≈ 100-200
  3. min_dist(最小距离)

    • 控制簇内点的紧密程度
    • 太小:点过于聚集,难以区分
    • 太大:簇过于分散,结构不清晰
    • 典型范围:0.001-0.5

调参必要性:

# UMAP参数敏感度示例
umap_varied <- bean_rec_trained %>%
  step_umap(
    all_numeric_predictors(), 
    num_comp = tune(),      # 必须调优
    neighbors = tune(),     # 必须调优  
    min_dist = tune()       # 必须调优
  )
 
# 不同参数会产生完全不同的可视化效果
# 进而影响下游分类器性能

与PCA对比:

  • PCA:无需调参,确定性算法
  • UMAP:基于流形学习,参数敏感性强

6. 完整模型工作流定义

# 定义5种不同的分类模型
mlp_spec <-
  mlp(hidden_units = tune(), penalty = tune(), epochs = tune()) %>%
  set_engine('nnet') %>%
  set_mode('classification')
 
bagging_spec <-
  bag_tree() %>%
  set_engine('rpart') %>%
  set_mode('classification')
 
fda_spec <-
  discrim_flexible(prod_degree = tune()) %>%
  set_engine('earth')
 
rda_spec <-
  discrim_regularized(frac_common_cov = tune(), frac_identity = tune()) %>%
  set_engine('klaR')
 
bayes_spec <-
  naive_Bayes() %>%
  set_engine('klaR')
 
# 定义3种预处理流程
bean_rec <-
  recipe(class ~ ., data = bean_train) %>%
  step_zv(all_numeric_predictors()) %>%
  step_orderNorm(all_numeric_predictors()) %>%
  step_normalize(all_numeric_predictors())
 
pls_rec <- 
  bean_rec %>% 
  step_pls(all_numeric_predictors(), outcome = "class", num_comp = tune())
 
umap_rec <-
  bean_rec %>%
  step_umap(
    all_numeric_predictors(),
    outcome = "class",  # 监督式UMAP,利用类别信息
    num_comp = tune(),
    neighbors = tune(), 
    min_dist = tune()
  )

7. 大规模模型比较

# 设置并行计算控制参数
ctrl <- control_grid(
  parallel_over = "everything",  # 最大化并行
  save_pred = TRUE,             # 保存预测结果
  save_workflow = TRUE          # 保存工作流对象
)
 
# 创建工作流集(2种预处理 × 2种模型 = 4种组合)
bean_workflows <-
	workflow_set(
		preproc = list(
			basic = bean_rec,
			pls = pls_rec,
 
		),
		models = list(
			bayes = bayes_spec,
			fda = fda_spec,
 
		)
	)
 
# 执行网格搜索调优
bean_res <- 
  bean_workflows %>% 
  workflow_map(
    verbose = TRUE,
    seed = 1603,
    resamples = bean_val,        # 使用验证集
    grid = 10,                   # 每个工作流尝试10组参数
    metrics = metric_set(roc_auc), # 使用AUC作为评估指标
    control = ctrl
  )
 
# 结果排名和分析
rankings <- 
  rank_results(bean_res, select_best = TRUE) %>% 
  mutate(
    method = map_chr(wflow_id, ~ str_split(.x, "_", simplify = TRUE)[1]),
    preproc = map_chr(wflow_id, ~ str_split(.x, "_", simplify = TRUE)[2])
  )
 
# 显示前5名模型
filter(rankings, rank <= 5) %>% 
  dplyr::select(rank, mean, model, preproc, method)

8. 最终模型评估

# 选择最佳模型(基于验证集性能)
best_model_id <- rankings$wflow_id[which.min(rankings$rank)]
 
# 提取并最终化最佳工作流
final_workflow <- 
  bean_res %>% 
  extract_workflow(best_model_id) %>% 
  finalize_workflow(
    bean_res %>% 
      extract_workflow_set_result(best_model_id) %>% 
      select_best(metric = "roc_auc")
  )
 
# 在测试集上进行最终评估(无偏估计)
final_results <- 
  final_workflow %>% 
  last_fit(split = bean_split, metrics = metric_set(roc_auc, accuracy))
 
# 提取训练好的最终模型
final_fit <- extract_workflow(final_results)
 
# 查看测试集性能
test_performance <- collect_metrics(final_results)
cat("最终模型在测试集上的表现:\n")
print(test_performance)
 
# 可视化预测结果
test_predictions <- collect_predictions(final_results)
 
# 混淆矩阵
conf_mat(test_predictions, truth = class, estimate = .pred_class) %>%
  autoplot(type = "heatmap")