我使用小鼠估算了回归模型和估算数据。
model1 <- with(imp, lm(outcome~ predictor1+ predictor2+ predictor3+ predictor4))).
在输出中我得到一些信息
summary(pool(model1), conf.int = TRUE)
如估计值、标准误差和 p 值。现在我想知道整个模型的F值和R^2。
对于 R^2,我找到了以下代码:
pool.r.squared(model1)
。但我仍在寻找显示 F 值的代码。有人有这方面的经验吗
传统的 F 统计量,我们通过对
anova
的 F 值进行平均得到,比较:
mean(anova(aov(bmi ~ hyp + chl, nhanes))[, 4], na.rm=TRUE)
summary(lm(bmi ~ hyp + chl, nhanes))$fstatistic[1]
对于汇总分析,我们可以使用
miceadds::mi.anova
来获取 R^2 和 F 统计量。
library('miceadds')
nul <- capture.output(
aov_fit <- miceadds::mi.anova(mi.res=imp, formula="bmi ~ hyp + chl" )
)
(
capture.output
不一定需要,但可以防止控制台混乱。)
所需的信息现已存储在对象中
aov_fit
。
aov_fit$r.squared ## R-squared
# [1] 0.1158705
(fval <- mean(round(aov_fit$anova.table$`F value`, 2), na.rm=TRUE) ) ## F-statistic
# [1] 0.97
df_mod <- aov_fit$anova.table$df1[- nrow(aov_fit$anova.table)] ## DF model
df_res <- el(fit$analyses)$df.residual ## DF residual
c(df_mod, df_res)
# [1] 1 1 22
可以使用 F 分布的分布函数通过右尾检验来计算模型 p 值
pf()
。
pf(q=fval, df1=sum(df_mod), df_2=df_res, lower.tail=FALSE) ## p-value
# [1] 0.3947152
我们现在可以使用
sprintf
来类似于 lm()
的 GOF 指标:
sprintf('Pooled R-squared: %s', round(aov_fit$r.squared, 4))
# [1] "Pooled R-squared: 0.1159"
tmp <- aov_fit$anova.table
sprintf('Pooled F-statistic: %s on %s and %s DF, p-value: %s',
mean(round(tmp$`F value`, 2), na.rm=TRUE),
round(sum(tmp$df1[- nrow(aov_fit$anova.table)]), 2),
round(el(fit$analyses)$df.residual, 2),
format.pval(pf(fval, sum(df_mod), df_res, lower.tail=FALSE)))
# [1] "Pooled F-statistic: 0.97 on 2 and 22 DF, p-value: 0.39472"
要获得 r2adj.,我们可以使用通常的公式,
adjR2 <- \(r2, n, p) {
1 - (n - 1)/(n - p - 1)*(1 - r2)
}
adjR2(aov_fit$r.squared, nrow(nhanes), sum(aov_fit$anova.table$df1, na.rm=TRUE))
# [1] 0.03549512
其中 n= 观测值数量,p= 参数数量。
数据:
使用
nhanes
包的 mice
数据集。
library('mice')
set.seed(42)
imp <- mice(nhanes, m=100, printFlag=FALSE)
fit <- with(data=imp, exp=lm(bmi ~ hyp + chl))
我无法发表评论,但对于那些尝试计算 p 值的人,您需要将 df_2 更改为 df2。我挣扎了好久,终于成功了。 @jay.sf - 你能确认是这样吗?