r-计算线性回归的stderr、t值、p值和预测值

Question

我正在用MatrixModels:::lm.fit.sparse和MatrixModels::glm4(也是稀疏的(拟合线性模型
但是，这些函数仅返回coeff、residuals和fitted.values。

获取/计算其他值(如stderr、t-value、p-value、predict值(最快、最简单的方法是什么？

Answer 1

我使用MatrixModels:::lm.fit.sparse示例中的数据。我构建了一个自定义函数summary_sparse来为这个模型执行摘要。所有的矩阵运算都是用Matrix包来执行的。结果与稠密型模型进行了比较。

注意lm.fit.sparse必须用method = "chol"进行评估才能得到正确的结果。

功能：

summary_sparse <- function(l, X) {

XXinv <- Matrix::chol2inv(Matrix::chol(Matrix::crossprod(X)))
se <- sqrt(Matrix::diag(XXinv*sum(l$residuals**2)/(nrow(X)-ncol(X))))
ts <- l$coef/se
pvals <- 2*c(1 - pnorm(abs(ts)))

list(coef = l$coef, se = se, t = ts, p = pvals)

}
predict_sparse <- function(X, coef) {
X %*% coef
}

应用：

dd <- expand.grid(a = as.factor(1:3),
b = as.factor(1:4),
c = as.factor(1:2),
d= as.factor(1:8))
n <- nrow(dd <- dd[rep(seq_len(nrow(dd)), each = 10), ])
set.seed(17)
dM <- cbind(dd, x = round(rnorm(n), 1))
## randomly drop some
n <- nrow(dM <- dM[- sample(n, 50),])
dM <- within(dM, { A <- c(2,5,10)[a]
B <- c(-10,-1, 3:4)[b]
C <- c(-8,8)[c]
D <- c(10*(-5:-2), 20*c(0, 3:5))[d]
Y <- A + B + A*B + C + D + A*D + C*x + rnorm(n)/10
wts <- sample(1:10, n, replace=TRUE)
rm(A,B,C,D)
})
X <- Matrix::sparse.model.matrix( ~ (a+b+c+d)^2 + c*x, data = dM)
Xd <- as(X,"matrix")
fmDense <- lm(dM[,"Y"]~Xd-1)
ss <- summary(fmDense)
r1 <- MatrixModels:::lm.fit.sparse(X, y = dM[,"Y"], method = "chol")
f <- summary_sparse(r1, X)
all.equal(do.call(cbind, f), ss$coefficients, check.attributes = F)
#TRUE
all.equal(predict_sparse(X, r1$coef)@x, predict(fmDense), check.attributes = F, check.names=F)
#TRUE