如何为glm.nb预测建模设置随机种子？种子选择标准咨询

问题描述 #

我不确定该如何为glm.nb模型训练设置随机种子，不同种子会导致模型结果出现差异。想请教：有没有种子选择的黄金标准？还是应该尝试多个种子后选择最优结果？

目前我发现使用种子271时能得到最佳结果——最高的R-squared值，同时RMSE和MAE最低。

我的代码 #

library(MASS)        # For glm.nb function
library(dplyr)       # For data manipulation
library(tidyverse)   # For ggplot2 and other utilities
library(caret)       # For data splitting and cross-validation
library(Metrics)     # For RMSE, MAE, and R-squared metrics

# Step 1: Define file paths
microbial_abundance_file <- "glm_nb/asv_abundance_all.csv"
scfa_concentration_file <- "glm_nb/acetate_w2.csv"

# Step 2: Load your data
microbial_abundance <- read.csv(microbial_abundance_file)
scfa_concentration <- read.csv(scfa_concentration_file)

# Step 3: Combine the data by SampleID
data_combined <- inner_join(scfa_concentration, microbial_abundance, by = "SampleID")

# Step 4: Convert SampleID to a character vector to avoid factor issues
data_combined$SampleID <- as.character(data_combined$SampleID)

# Remove the first column for modeling
data_combined <- data_combined[, -1]

# Step 5: Split data into training and test sets
set.seed(271)  # For reproducibility
train_index <- createDataPartition(data_combined$Acetic.acid, p = 0.75, list = FALSE)
train_data <- data_combined[train_index, ]
test_data <- data_combined[-train_index, ]

# Define control for cross-validation within the training set
cv_control <- trainControl(method = "cv", number = 10)  # 10-fold cross-validation

# Define the columns to include in the model
columns <- c(
  "d__Bacteria.p__Firmicutes.c__Clostridia.o__Oscillospirales.f__Ruminococcaceae.g__Faecalibacterium.__",
  "d__Bacteria.p__Firmicutes.c__Negativicutes.o__Veillonellales.Selenomonadales.f__Selenomonadaceae.g__Megamonas.__"
)

# Create the formula dynamically
formula <- as.formula(paste("Acetic.acid ~", paste(columns, collapse = " + ")))

# Print the formula to verify
print(formula)

# Train the model with the updated formula
cv_model <- train(
  formula, 
  data = train_data, 
  method = "glm.nb", 
  trControl = cv_control
)

# Print cross-validation results
print(cv_model$results) 

# Display the final model details
print(cv_model$finalModel) 

# Step 7: Evaluate the final model on the test set
predicted_values <- predict(cv_model, newdata = test_data, type = "raw")
actual_values <- test_data$Acetic.acid
residuals <- actual_values - predicted_values

# Calculate performance metrics on the test data
test_rmse <- sqrt(mean(residuals^2))
test_mae <- mean(abs(residuals))
test_r_squared <- 1 - (sum(residuals^2) / sum((actual_values - mean(actual_values))^2))

cat("Test RMSE:", test_rmse, "\n")
cat("Test MAE:", test_mae, "\n")
cat("Test R-squared:", test_r_squared, "\n")

# Step 8: Plot predicted vs actual SCFA concentration for the test set
ggplot(data.frame(actual_values, predicted_values), aes(x = actual_values, y = predicted_values)) +
  geom_point() +
  geom_smooth(method = "lm", col = "red") +
  labs(title = "Predicted vs Actual SCFA Concentration (Test Set)",
       x = "Actual SCFA Concentration",
       y = "Predicted SCFA Concentration") +
  annotate("text", x = min(actual_values), y = max(predicted_values), 
           label = paste("RMSE:", round(test_rmse, 2)), hjust = 0, vjust = 1.5) +
  annotate("text", x = min(actual_values), y = max(predicted_values) - 0.5, 
           label = paste("MAE:", round(test_mae, 2)), hjust = 0, vjust = 1.5) +
  annotate("text", x = min(actual_values), y = max(predicted_values) - 1, 
           label = paste("R²:", round(test_r_squared, 2)), hjust = 0, vjust = 1.5)

解答 #

• 不存在种子选择的黄金标准：随机种子的核心作用是保证结果可复现，没有哪个种子天生"更优"。你用种子271得到的好结果，只是当前数据划分下的偶然情况。
• 不要仅选择最优种子的结果：如果换个种子模型性能就大幅波动，说明你的模型稳定性不足，根源可能是数据量偏小、特征噪声高，或者模型拟合状态有问题（过拟合/欠拟合）。
• 正确的实践方案：
  1. 固定种子用于复现：如果需要让他人复现你的实验结果，任选一个固定种子（比如你用的271）即可，把种子值明确写在代码里。
  2. 多种子验证稳定性：用10-20个不同的随机种子，分别完成数据划分、模型训练和测试评估，然后计算性能指标的均值和标准差。如果标准差较大，说明结果不可靠，需要优化数据或模型（比如增加样本量、筛选特征、换用更鲁棒的模型）。
  3. 严格避免数据泄露：绝对不能用测试集的性能来选择种子——这相当于把测试集的信息提前用到了模型选择环节，会导致最终的模型评估结果被高估。
• 针对你代码的优化建议：
  • 给交叉验证也固定种子：在trainControl中添加seeds参数，确保10折交叉验证的每一轮划分也可复现，进一步提升实验的严谨性。
  • 封装多种子测试函数：写一个循环函数，批量运行不同种子的实验，自动统计性能指标的分布情况，节省手动测试的时间。

内容的提问来源于stack exchange，提问作者Peter