从 R 中的 lm 脚本获取正确的输出

Question

我正在尝试计算 50 多个 SNP 的 Beta 误差和标准误差，以及它们对绝经年龄的个体影响。我的脚本不太正确。我想使用 country 作为协变量并对绝经年龄 (a_menopause) 和每个单独的 SNP 进行线性回归。

当我通过SNP计算SNP时，我得到了正确的答案。这是我使用的命令：

rs <- lm(a_menopause~rs17465637_metabo + country, data="x")
summary(rs)

输出看起来像这样（我只对第二行 (rs) 感兴趣，特别是估计值和标准误差：

Coefficients:
                      Estimate Std. Error t value Pr(>|t|)    
(Intercept)          49.515274   0.340216 145.541   <2e-16 ***
rs                   -0.004258   0.171557  -0.025    0.980    
RS$country            0.019684   0.047994   0.410    0.682

我现在有以下脚本来尝试每次用新变量替换 rs17465637_metabo。然而，当我查看我的输出时，所有的值都非常相似（并且与我一次只做一个 SNP 时非常不同）。我也只得到（截距）值而不是 RS 值。

这是我的脚本：

read.data <- function(file="RS.csv") {
  x <- read.csv(file, sep=";", dec=",") 
  x$country <- factor(x$country)        
  x
}

## Fit the LMs
fit.lms <- function(x) {
  rsid <- colnames(x)[-(1:3)]         # all column names except the first three
  res <- vector("list", length=length(rsid))
  i <- 1
  for (rs in rsid) {
    modeltxt <- paste("lm(a_menopause ~ 0 + ", rsid, " + country, data=x)") # build model expression text
    r <- eval(parse(text=modeltxt)) 
    s <- summary(r)$coefficients    
    res[[i]] <- c(s[2,], s[2,])     
    i <- i+1
  }
  data.frame(do.call(rbind, res), rsid=rsid)
}

run <- function() {
  x <- read.data()
  lms <- fit.lms(x)
  write.table(lms, file="fits.csv", sep=";", dec=",", row.names=F)
}

我哪里错了？如何更改输出以提供第二行数据？

提前致谢！

Answer 1

不use eval(parse())但是as.formula(paste())

my.form <- as.formula(paste("a_menopause ~ 0 +", rs, "+ country")))
r <- lm(my.form, data=x)

这是提取第二行的方法。

output <- lapply(rsid, function(rs){
    my.form <- as.formula(paste("a_menopause ~ 0 +", rs, "+ country")))
    r <- lm(my.form, data=x)
    summary(r)$coefficients[2, c("Estimate", "Std. Error")]
})
do.call(rbind, output)

请注意，由于您未能提供样本数据，因此未经测试

从 R 中的 lm 脚本获取正确的输出

Getting the right output from lm script in R

r

lm