首页 > 解决方案 > 使用带有 lm 的样本权重时更正 dfs

问题描述

我试图弄清楚权重lm实际上是如何工作的,我看到了这个 7.5 年前的问题,它提供了一些关于权重如何工作的见解。该问题的数据部分复制并在下面扩展。

我在 Cross Validated 上发布了这个相关的问题。

library(plyr)
set.seed(100)
df <- data.frame(uid=1:200,
                      bp=sample(x=c(100:200),size=200,replace=TRUE),
                      age=sample(x=c(30:65),size=200,replace=TRUE),
                      weight=sample(c(1:10),size=200,replace=TRUE),
                      stringsAsFactors=FALSE)

set.seed(100)
df.double_weights <- data.frame(uid=1:200,
                      bp=sample(x=c(100:200),size=200,replace=TRUE),
                      age=sample(x=c(30:65),size=200,replace=TRUE),
                      weight=2*df$weight,
                      stringsAsFactors=FALSE)

df.expand <- ddply(df,
                        c("uid"),
                        function(df) {
                          data.frame(bp=rep(df[,"bp"],df[,"weight"]),
                                     age=rep(df[,"age"],df[,"weight"]),
                                     stringsAsFactors=FALSE)})

df.lm <- lm(bp~age,data=df,weights=weight)
df.double_weights.lm <- lm(bp~age,data=df.double_weights,weights=weight)
df.expand.lm <- lm(bp~age,data=df.expand)

summary(df.lm)
summary(df.double_weights.lm)
summary(df.expand.lm)

这三个 data.frames 包含完全相同的数据。然而;

其中df有 200 个观察值,加权加起来为 1178 sum(df.$weight) == 1178,.

df.double_weights中,权重只是翻了一番sum(df.double_weights$weight) == 2356

df.expand中,有 1178 个未加权观测值,而不是 200 个加权观测值。

两者的系数summary(df.lm)相同summary(df.double_weights.lm),重要性也相同(这意味着,如果权重工作正常,权重的绝对大小无关紧要)。编辑:但是,绝对大小似乎确实很重要,请参阅底部结果。

但是,对于summary(df.lm)summary(df.expand.lm),系数相同,但显着性不同。

summary(df.lm)

            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 165.6545    10.3850  15.951   <2e-16 ***
age          -0.2852     0.2132  -1.338    0.183    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 98.84 on 198 degrees of freedom
Multiple R-squared:  0.008956,  Adjusted R-squared:  0.003951 
F-statistic: 1.789 on 1 and 198 DF,  p-value: 0.1825

summary(df.expand.lm)

             Estimate Std. Error t value Pr(>|t|)    
(Intercept) 165.65446    4.26123   38.88  < 2e-16 ***
age          -0.28524    0.08749   -3.26  0.00115 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 28.68 on 1176 degrees of freedom
Multiple R-squared:  0.008956,  Adjusted R-squared:  0.008114 
F-statistic: 10.63 on 1 and 1176 DF,  p-value: 0.001146

根据@IRTFM 的说法,自由度没有被正确添加,提供了以下代码来修复它:

df.lm.aov <- anova(df.lm)
df.lm.aov$Df[length(df.lm.aov$Df)] <- 
        sum(df.lm$weights)-   
        sum(df.lm.aov$Df[-length(df.lm.aov$Df)]  ) -1
df.lm.aov$`Mean Sq` <- df.lm.aov$`Sum Sq`/df.lm.aov$Df
df.lm.aov$`F value`[1] <- df.lm.aov$`Mean Sq`[1]/
                                        df.lm.aov$`Mean Sq`[2]
df.lm.aov$`Pr(>F)`[1] <- pf(df.lm.aov$`F value`[1], 1, 
                                      df.lm.aov$Df, lower.tail=FALSE)[2]
df.lm.aov

Analysis of Variance Table

Response: bp
            Df Sum Sq Mean Sq F value   Pr(>F)   
age          1   8741  8740.5  10.628 0.001146 **
Residuals 1176 967146   822.4                    

现在,将近 8 年后,显然这个问题仍然存在(这是否意味着几乎所有使用加权变量与lmfrom结合的研究的R显着性值都太低了?)更实际地,我遇到的问题是我几乎不明白 IRTFM 是什么正在做什么,或者它与多元回归分析(甚至是在后台使用lm的其他函数)有何关系?

问题:有没有更通用的方法来解决这个问题,可以应用于多元回归?

编辑:

如果我们在 上运行 IRTFM 的解决方案df.double_weights.lm,我们会得到不同的结果,因此显然权重的绝对大小确实很重要。

Analysis of Variance Table

Response: bp
            Df  Sum Sq Mean Sq F value    Pr(>F)    
age          1   17481 17481.0  21.274 4.194e-06 ***
Residuals 2354 1934293   821.7                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

标签: rlmweighted

解决方案


如果我正确理解您的问题,那么您在权重列中的内容通常称为“频率权重”。它们用于通过指示您对每个协变量组合有多少观察值来节省数据集中的空间。

要使用“聚合”数据集估计模型并获得正确的标准误差,您需要做的就是更正模型中的自由度数lm

正确的自由度数是观察总数减去模型中的参数数。这可以通过获取weights变量的总和或通过查看“完整”数据中的观察总数并减去估计的参数数量(即系数)来计算。

这是一个更简单的例子,我认为这更清楚了这一点:

library(dplyr)
library(modelsummary)

set.seed(1024)

# individual (true) dataset
x <- round(rnorm(1e5))
y <- round(x + x^2 + rnorm(1e5))
ind <- data.frame(x, y)

# aggregated dataset
agg <- ind %>%
  group_by(x, y) %>%
  summarize(freq = n())

models <- list( 
  "True"                = lm(y ~ x, data = ind),
  "Aggregated"          = lm(y ~ x, data = agg),
  "Aggregated & W"      = lm(y ~ x, data = agg, weights=freq),
  "Aggregated & W & DF" = lm(y ~ x, data = agg, weights=freq)
)

现在我们要更正列表中最后一个模型的自由度数。我们通过对我们的freq列求和来做到这一点。我们也可以使用nrow(ind),因为它们是相同的:

# correct degrees of freedom
models[[4]]$df.residual <- sum(agg$freq) - length(coef(models[[4]]))

最后,我们总结了使用该modelsummary包的所有 5 个模型。请注意,第一个模型和最后一个模型完全相同,即使第一个模型是使用完整的单个数据集估计的,而最后一个模型是使用聚合数据估计的:

modelsummary(models, fmt=5)
真的 聚合的 聚合 & W 聚合 & W & DF
(截距) 1.08446 5.51391 1.08446 1.08446
(0.00580) (0.71710) (0.22402) (0.00580)
X 1.00898 0.91001 1.00898 1.00898
(0.00558) (0.30240) (0.21564) (0.00558)
编号Obs。 1e+05 69 69 69
R2 0.246 0.119 0.246 0.246
R2 调整。 0.246 0.106 0.235 0.999
AIC 405058.1 446.0 474.1 474.1
BIC 405086.7 452.7 480.8 480.8
日志.Lik。 -202526.074 -219.977 -234.046 -234.046
F 32676.664 9.056 21.894 32676.664

推荐阅读