使用regsubsets筛选最优模型后lm结果调整R平方不符的技术求助

• 使用regsubsets工具寻找使**调整R平方（adjusted r squared）**最大化的自变量组合，设置nvmax为自变量总数
• 提取到的最优模型调整R平方，与用选中变量构建lm回归模型后的调整R平方差距极大，后者表现很差
• 存在模型识别错误：可视化显示18个预测变量的模型调整R平方最高，但代码提取了19个预测变量的模型信息，且lm结果也不匹配该19变量模型的调整R平方值

用户原代码：

##############################################

Op = V[1:60,c(2:60)]  
attach(Op)  
dependent_var <- Op$Prequin.VC #select the VC Index for dependant variable  
independent_vars <- Op[, c(11:27,30:46,48:50)] #define all the columns of predictors  

#test all combinations for independent variables
result <- regsubsets(dependent_var ~ ., data = independent_vars, nvmax = dim(independent_vars)[2], really.big = TRUE)  
 ##Stock adjusted r squared for each models
adjr2_values <- summary(result)$adjr2  
adjr2_values  
#Stock the model with the best Adjust R squared
best_model <- which.max(adjr2_values)  

#Identify variables of best model
cat("Best variables combinations : ", paste(colnames(independent_vars) 
 [which(coef(result, id = best_model) != 0)], collapse = ", "), "
")  
cat("Number of variables : ", paste(num_selected_variables <- length(coef(result, id = best_model))-1))  
cat("Best adjusted Rsquare : ", adjr2_values[best_model], "
")  

#plot summary lm of the best model
selected_vars <- colnames(independent_vars)[which(coef(result, id = best_model) != 0)]  
formula <- as.formula(paste("dependent_var ~", paste(selected_vars, collapse = " + ")))  
lm_model <- lm(formula)  
summary(lm_model)  

1. attach()引发的变量环境冲突 #

attach(Op)会将数据框的列加载到全局环境，后续构建lm模型时，变量可能被误读为全局环境中的对象，而非independent_vars数据框内的对应变量，导致拟合数据错误。

• 解决：直接移除attach(Op)，全程使用数据框引用明确指定变量来源。

2. 模型索引与变量数量的匹配错误 #

summary(result)$adjr2的索引顺序是从1个变量的模型开始，索引值直接对应模型的变量数量（如索引18对应18个变量的模型）。原代码中which.max(adjr2_values)返回的索引值就是最优模型的变量数，但后续提取系数时未严格对应，导致选到错误模型。

• 解决：用which.max的结果直接作为模型的变量数量参数，传递给coef(result, id = ...)。

3. lm模型未指定数据来源 #

原代码中lm(formula)未指定data参数，lm会在全局环境中查找变量，但selected_vars是independent_vars的列名，全局环境中不存在这些变量，导致模型拟合的数据集错误。

• 解决：将因变量与自变量合并到同一数据框，构建lm时明确指定data参数。

修正后的代码 #

##############################################
# 移除attach，避免环境冲突
Op <- V[1:60, c(2:60)]  
dependent_var <- Op$Prequin.VC # 选择VC指数作为因变量
independent_vars <- Op[, c(11:27, 30:46, 48:50)] # 定义所有预测变量列

# 测试所有自变量组合
result <- regsubsets(dependent_var ~ ., data = independent_vars, 
                     nvmax = ncol(independent_vars), really.big = TRUE)  

# 提取各模型的调整R平方
adjr2_values <- summary(result)$adjr2
print(adjr2_values)

# 找到调整R平方最大的模型对应的变量数量
best_nvar <- which.max(adjr2_values)
cat("最优模型的变量数量：", best_nvar, "\n")

# 提取最优模型的变量
best_coef <- coef(result, id = best_nvar)
selected_vars <- names(best_coef)[-1] # 移除截距项
cat("最优变量组合：", paste(selected_vars, collapse = ", "), "\n")
cat("最优调整R平方：", adjr2_values[best_nvar], "\n")

# 构建并拟合正确的lm模型
model_data <- cbind(independent_vars, dependent_var) # 合并因变量与自变量
formula <- as.formula(paste("dependent_var ~", paste(selected_vars, collapse = " + ")))
lm_model <- lm(formula, data = model_data)
summary(lm_model)

内容的提问来源于stack exchange，提问作者Cédric Delaveau