与R中列名称对应的指标矩阵

Question

我要解决的问题是“如何创建一系列自动化代码，这些代码将从数据集中提取所需的列标题名称以适合一般的线性化模型（glm）？”我有一个包含8个变量的数据集；但是，我只想使用其中的3个来交叉验证并找到“最佳”模型。这是我想出的：

library(boot)

if (! file.exists("salary_data.csv"))
  download.file("http://homepage.stat.uiowa.edu/~rdecook/stat6220/datasets/salary_data.csv", "salary_data.csv")
salary <- read.csv("salary_data.csv")
vars <- colnames(salary[c(2,3,7)])
nvars <- length(vars)
list.to.expand = vector(mode = "list", length = nvars)
for (i in 1:nvars){
  list.to.expand[[i]]=c(0,1)
}
model.spec.matrix <- expand.grid(list.to.expand)

vars
model.spec.matrix

names(model.spec.matrix) <- vars
vars.to.use <- model.spec.matrix[1,]
vars.to.use <- as.numeric(vars.to.use)
model.vars <- c()
count = 1
for (i in vars.to.use){
  if(i==1){model.vars <- append(model.vars, vars[count])
  } else {
    count + 1
  }
}
paste(model.vars, collapse = "+")
glm.out = glm(paste("LogACG~",paste(model.vars,collapse = "+"),sep = ""), family = gaussian, data = salary)
cv.err = cv.glm(salary, glm.out, K = 10)$delta

我的问题在于for循环。我试图构造一个循环，将“ vars”中的值附加到“ model.use”中，但是我似乎无法让它读取到矩阵的第二行之外。有什么建议么？谢谢

Answer 1

似乎这里发生了几件事。

您已将LogACG设置为倒数第二行中的解释变量（"LogACG~，但它也是vars之一，由于model.vars而最终以vars <- colnames(salary[c(2,3,7)])结尾，因此不正确。

接下来您的第二个for循环应遍历model.spec.matrix的各行，即

for(i in 1:nrow(model.spec.matrix)){

并以编程方式捕获该行所指示的列名（变量），您可以这样做>

cn <- colnames(model.spec.matrix[sapply(model.spec.matrix[i,], function(x) x > 0)])

在循环内。您还应该在循环内移动glm和cv.glm。

因此它应该看起来更接近

for(i in 1:nrow(model.spec.matrix)){
  cn <- colnames(model.spec.matrix[sapply(model.spec.matrix[i,], function(x) x > 0)])
  if(length(cn) == 0){cn <- "."}  
  glm.out = glm(as.formula(paste("LogACG~",paste(cn,collapse = "+"),sep = "")), family = gaussian, data = salary)
  cv.err = cv.glm(salary, glm.out, K = 10)$delta
}
但是，每次都会覆盖glm.out和cv.err，因此您需要将它们创建为空列表并在每次迭代中追加该列表。

因此最终产品将看起来像这样（但请注意，该示例需要对因素进行一些琐碎的功能设计，而这些因素可能只针对该示例）：

# Since you can't use LogACG to explain itself, 
#    suppose you meant to use Engineering as a candidate X
vars <- colnames(salary[c(2,3,8)])

# Make your grid
model.spec.matrix        <- expand.grid(list.to.expand)
names(model.spec.matrix) <- vars

glm.out <- rep(NA, nrow(model.spec.matrix))
cv.err  <- rep(NA, nrow(model.spec.matrix))

for(i in 1:nrow(model.spec.matrix)){
  cn <- colnames(model.spec.matrix[sapply(model.spec.matrix[i,], function(x) x > 0)])
  if(length(cn) == 0){cn <- "."}  
  tmp        <- glm(as.formula(paste("LogACG~",paste(cn,collapse = "+"),sep = "")), family = gaussian, data = salary)
  glm.out[i] <- tmp
  cv.err[i]  <-  cv.glm(salary, tmp, K = 10)$delta
}
此外，最好使用bestglm之类的包，或者仅使用默认包装step中的stats函数（请参见?step）。

[如果我误解了任何一个，请告诉我。