基于预定义标准的多元线性模型逐步选择
Multivariate linear model stepwise selection based on predefined criteria
我正在尝试编写自己的 R 函数,类似于向前逐步选择 step
,但我没有使用 AIC 作为选择标准,而是有一些标准需要在每次预测时进行评估添加了变量。模型的构建原理解释如下。该模型应从与因变量具有最高相关性的预测变量开始。然后每次根据新模型是否满足以下条件添加另一个预测变量。
- 调整后的r2值必须增加1%以上;
- 新增变量和现有变量的系数必须为正;
- 添加的变量必须显着,即 p 值 < 0.05。
重复此过程,直到没有剩余变量满足所有三个条件。我需要的输出可能只是最终模型中所有预测变量的名称、相应的系数和最终模型的 r2 值。
我的示例数据(y 是因变量,x1 - x6 是预测变量)
data = structure(list(y = c(23.6, 19.9, 40.7, 40.7, 40.7, 40.7, 40.2,
41.7, 41.7, 28.8), x1 = c(0.1, 0, 0.3, 0.3, 0.3,
0.3, 0.3, 0.3, 0.3, 0.1), x2 = c(0, 0.1, 0, 0, 0,
0, 0, 0.1, 0.1, 0), x3 = c(2277.6, 3038.1, 7797.9, 7797.9,
7797.9, 7797.9, 8392.2, 10127.2, 10127.2, 1799), x4 = c(34228.7,
49815, 76917.1, 76917.1, 76917.1, 76917.1, 75981.4, 74881.1,
74881.1, 56798.2), x5 = c(108786.5, 150465.5, 230397.1, 230397.1,
230397.1, 230397.1, 239300.9, 238493.8, 238493.8, 188799.5),
x6 = c(362.2, 198.2, 656.6, 656.6, 656.6, 656.6, 681,
655.3, 655.3, 222.3)), .Names = c("y", "x1",
"x2", "x3", "x4", "x5", "x6"), row.names = c(NA,
10L), class = "data.frame")
第一次尝试我的模型选择功能
modSel = function(data, var){
cor.result = cor(data[,var], df["y"]) #calculate correlation coeff for each variable against y
max.cor = rownames(cor.result)[which.max(cor.result)] #identify the variable with max cor
start.model = lm(as.formula(paste("y", max.cor, sep = "~")), data)
if #my criteria??
else #??]
没有太多的编程背景,我真的不知道如何将我的标准的评估重复未知次数。我意识到要实现这一点可能需要相当多的编码,但对于初学者来说,我将不胜感激任何关于整个框架应该是什么样子的指导。
干杯
希望这可以帮助您入门
对运行算法的函数
modSel <- function(data) {
# initial
cor.result <- cor(data$y, data[, -which(colnames(data) == "y")])
vars.model <- colnames(cor.result)[which.max(cor.result)]
vars.remaining <- colnames(data)[!colnames(data) %in% c("y", max.cor)]
start.model <- lm(as.formula(paste("y", vars.model, sep = "~")), data)
adj.rsq <- summary(start.model)$adj.r.squared
# algorithm
for (var in vars.remaining) {
# model
vars.test <- paste(vars.model, var, sep="+")
fit <- lm(as.formula(paste("y", vars.test, sep="~")), data)
new.rsq <- summary(fit)$adj.r.squared
# check adj rsq
cond1 <- new.rsq > adj.rsq + .01
# check coefficients
cond2 <- coefficients(fit)[var] > 0
# new var significant
cf <- summary(fit)$coefficients[, 4]
cond3 <- cf[var] < .05
if (cond1 & cond2 & cond3) {
vars.model <- vars.test
adj.rsq <- new.rsq
}
}
lm(as.formula(paste("y", vars.model, sep="~")), data)
}
调用 modSel
returns 算法中的最佳模型
bestfit <- modSel(data)
总结模型
summary(bestfit)
Call:
lm(formula = as.formula(paste("y", vars.model, sep = "~")), data = data)
Residuals:
Min 1Q Median 3Q Max
-0.8731 -0.3838 -0.3838 0.6273 1.5640
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 1.331e+01 1.941e+00 6.856 0.000241 ***
x1 5.417e+01 5.834e+00 9.285 3.48e-05 ***
x4 1.498e-04 4.460e-05 3.359 0.012099 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.9345 on 7 degrees of freedom
Multiple R-squared: 0.9904, Adjusted R-squared: 0.9876
F-statistic: 360.5 on 2 and 7 DF, p-value: 8.721e-08
我正在尝试编写自己的 R 函数,类似于向前逐步选择 step
,但我没有使用 AIC 作为选择标准,而是有一些标准需要在每次预测时进行评估添加了变量。模型的构建原理解释如下。该模型应从与因变量具有最高相关性的预测变量开始。然后每次根据新模型是否满足以下条件添加另一个预测变量。
- 调整后的r2值必须增加1%以上;
- 新增变量和现有变量的系数必须为正;
- 添加的变量必须显着,即 p 值 < 0.05。
重复此过程,直到没有剩余变量满足所有三个条件。我需要的输出可能只是最终模型中所有预测变量的名称、相应的系数和最终模型的 r2 值。
我的示例数据(y 是因变量,x1 - x6 是预测变量)
data = structure(list(y = c(23.6, 19.9, 40.7, 40.7, 40.7, 40.7, 40.2,
41.7, 41.7, 28.8), x1 = c(0.1, 0, 0.3, 0.3, 0.3,
0.3, 0.3, 0.3, 0.3, 0.1), x2 = c(0, 0.1, 0, 0, 0,
0, 0, 0.1, 0.1, 0), x3 = c(2277.6, 3038.1, 7797.9, 7797.9,
7797.9, 7797.9, 8392.2, 10127.2, 10127.2, 1799), x4 = c(34228.7,
49815, 76917.1, 76917.1, 76917.1, 76917.1, 75981.4, 74881.1,
74881.1, 56798.2), x5 = c(108786.5, 150465.5, 230397.1, 230397.1,
230397.1, 230397.1, 239300.9, 238493.8, 238493.8, 188799.5),
x6 = c(362.2, 198.2, 656.6, 656.6, 656.6, 656.6, 681,
655.3, 655.3, 222.3)), .Names = c("y", "x1",
"x2", "x3", "x4", "x5", "x6"), row.names = c(NA,
10L), class = "data.frame")
第一次尝试我的模型选择功能
modSel = function(data, var){
cor.result = cor(data[,var], df["y"]) #calculate correlation coeff for each variable against y
max.cor = rownames(cor.result)[which.max(cor.result)] #identify the variable with max cor
start.model = lm(as.formula(paste("y", max.cor, sep = "~")), data)
if #my criteria??
else #??]
没有太多的编程背景,我真的不知道如何将我的标准的评估重复未知次数。我意识到要实现这一点可能需要相当多的编码,但对于初学者来说,我将不胜感激任何关于整个框架应该是什么样子的指导。
干杯
希望这可以帮助您入门
对运行算法的函数
modSel <- function(data) {
# initial
cor.result <- cor(data$y, data[, -which(colnames(data) == "y")])
vars.model <- colnames(cor.result)[which.max(cor.result)]
vars.remaining <- colnames(data)[!colnames(data) %in% c("y", max.cor)]
start.model <- lm(as.formula(paste("y", vars.model, sep = "~")), data)
adj.rsq <- summary(start.model)$adj.r.squared
# algorithm
for (var in vars.remaining) {
# model
vars.test <- paste(vars.model, var, sep="+")
fit <- lm(as.formula(paste("y", vars.test, sep="~")), data)
new.rsq <- summary(fit)$adj.r.squared
# check adj rsq
cond1 <- new.rsq > adj.rsq + .01
# check coefficients
cond2 <- coefficients(fit)[var] > 0
# new var significant
cf <- summary(fit)$coefficients[, 4]
cond3 <- cf[var] < .05
if (cond1 & cond2 & cond3) {
vars.model <- vars.test
adj.rsq <- new.rsq
}
}
lm(as.formula(paste("y", vars.model, sep="~")), data)
}
调用 modSel
returns 算法中的最佳模型
bestfit <- modSel(data)
总结模型
summary(bestfit)
Call:
lm(formula = as.formula(paste("y", vars.model, sep = "~")), data = data)
Residuals:
Min 1Q Median 3Q Max
-0.8731 -0.3838 -0.3838 0.6273 1.5640
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 1.331e+01 1.941e+00 6.856 0.000241 ***
x1 5.417e+01 5.834e+00 9.285 3.48e-05 ***
x4 1.498e-04 4.460e-05 3.359 0.012099 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.9345 on 7 degrees of freedom
Multiple R-squared: 0.9904, Adjusted R-squared: 0.9876
F-statistic: 360.5 on 2 and 7 DF, p-value: 8.721e-08