4

对于 R 中的多个(约 100 万)响应变量,我必须used lm()拟合多个回归模型。

allModels <- lm(t(responseVariablesMatrix) ~ modelMatrix)

这将返回一个“mlm”类的对象,它就像一个包含所有模型的巨大对象。我想获得每个模型中第一个系数的t 统计量,我可以使用该summary(allModels)函数来完成,但是在这个大数据上它非常慢并且也返回了很多不需要的信息。

是否有更快的t-statistic手动计算方法,可能比使用summary()函数更快

谢谢!

4

1 回答 1

1

您可以破解 summary.lm() 函数以获取您需要的位并留下其余部分。

如果你有

nVariables <- 5
nObs <- 15

y <- rnorm(nObs)
x <- matrix(rnorm(nVariables*nObs),nrow=nObs)

allModels <-lm(y~x)

然后这是来自 lm.summary() 函数的代码,但删除了所有多余的行李(请注意,所有错误处理也已删除)。

p <- allModels$rank
rdf <- allModels$df.residual
Qr <- allModels$qr
n <- NROW(Qr$qr)
p1 <- 1L:p
r <- allModels$residuals
f <- allModels$fitted.values
w <- allModels$weights
mss <- if (attr(allModels$terms, "intercept")) 
sum((f - mean(f))^2) else sum(f^2)
rss <- sum(r^2)
resvar <- rss/rdf
R <- chol2inv(Qr$qr[p1, p1, drop = FALSE])
se <- sqrt(diag(R) * resvar)
est <- allModels$coefficients[Qr$pivot[p1]]
tval <- est/se

tval现在是 t 统计量的向量,也由

summary(allModels)$coefficients[,3]

如果您在大型模型上遇到问题,您可能需要重写代码,以便通过将多行/分配复合成更少的行来保留更少的对象。

我知道的哈克解决方案。但它会尽可能快。我想将所有代码行也放入一个函数中会更整洁。

于 2013-05-16T02:05:33.167 回答