使用dplyr按组计算加权统计量的R代码

Question

多年来，我一直使用Hmisc包和基本R语言计算加权统计摘要。通常情况下，我使用双重权重，其中一个是空间影响权重，另一个是数据支撑值，例如长度、体积、物理密度等。以“mtcars”数据集为例，其中mpg是所关注的变量，双重权重由汽车“wt”和“hp”构成，Hmisc + base R的工作流程通常如下。

require(Hmisc)

mtcars$Wt2 <- mtcars$wt * mtcars$hp               # double weight
mtcars$Acc <- mtcars$Wt2 * mtcars$mpg             # accumulation

min(mtcars$mpg)                                   # min
sqrt(wtd.var(mtcars$mpg, mtcars$mpg))             # wtd sd
wtd.quantile(mtcars$mpg, mtcars$mpg,0.05)         # wtd 5th 
wtd.quantile(mtcars$mpg, mtcars$mpg,0.50)         # wtd median
wtd.quantile(mtcars$mpg, mtcars$mpg,0.95)         # wtd 95th
max(mtcar$mpg)                                    # max

使用循环，可以对数据框中每个感兴趣的区域进行加权统计和筛选。 然而，由于我想学习如何使用dplyr，我在思考如何计算这些加权统计量。虽然有加权平均选项可用，但其他选项需要更多的工作。以下是我已经从头开始计算加权平均值并通过与dplyr内置函数的比较进行了检查的代码。然而，我被卡住了，不知道如何使用dplyr继续计算加权标准差和分位数，因为我需要将每个组的平方加权均值差（for each group）放入管道链中。

mtcars %>% 
  mutate(Car = row.names(mtcars)) %>%  # variable for car name
  mutate(Wt2 = wt * hp) %>%            # double weight
  mutate(Acc = Wt2 * mpg) %>%          # weighted consumption
  group_by(Car) %>%                    # group by car type
  summarise(n = n(),
            SmWt2 = sum(Wt2),                    # Sum of double weight
            SmAcc = sum(Acc),                    # Sum of accumulations
            WtMn = SmAcc/SmWt2,                  # Weighted mean
            WtMnChk = weighted.mean(mpg, Wt2)    # Check weighted mean
            )

- Markm0705

---

这个问题是否指引你朝着正确的方向？链接 - Kent Orr

在第一部分中，Acc 是 Wt2 * Wt2，在第二部分中是 Wt2 * mpg，哪一个是想要的？ - Jon Spring

我不太清楚您在单个汽车组时所说的加权平均值是什么意思。我认为您可能想要在一个组中获取加权平均值，比如所有四缸汽车的加权平均值... - Jon Spring

Jon - 感谢您发现了那个错误 - 我已经将示例更正为 mpg。您关于单个汽车加权平均数的评论也很明智 - 我应该像您建议的那样提供一个示例，因为那样更有意义。然而，问题的核心是关于使用 plyr 进行加权计算，目前仅限于加权平均数和加权残差。 - Markm0705

Answer 1

我不确定我完全理解你正在进行的方法，但这里是一个示例，通过 gear 找到加权平均值和加权标准差，使用wt作为加权：

library(dplyr)
datasets::mtcars %>% 
  group_by(gear) %>%
  summarize(n = n(),
            mpg_weighted_by_weight = sum(mpg*wt) / sum(wt),
            mpg_weighted_by_weight_check = weighted.mean(mpg, wt),
            
            mpg_sd = sqrt(sum(wt * ((mpg - mpg_weighted_by_weight)^2))/(sum(wt)-1)),
            mpg_sd_check = sqrt(Hmisc::wtd.var(mpg, wt)))


# A tibble: 3 x 6
   gear     n mpg_weighted_by_weight mpg_weighted_by_weight_check mpg_sd mpg_sd_check
* <dbl> <int>                  <dbl>                        <dbl>  <dbl>        <dbl>
1     3    15                   15.6                         15.6   3.32         3.32
2     4    12                   23.6                         23.6   4.81         4.81
3     5     5                   19.7                         19.7   5.63         5.63

我不熟悉加权标准差的公式，相反地使用了Hmisc::wtd.var中的公式。如果你在RStudio中控制点击公式名称，它会显示函数的底层代码。大部分是错误处理，一直到底部:

#Hmisc::wtd.var
function (x, weights = NULL, normwt = FALSE, na.rm = TRUE, method = c("unbiased", 
  "ML")) 
{
  # ...  skipping error handling
  sw <- sum(weights)
  # ...
  xbar <- sum(weights * x)/sw
  sum(weights * ((x - xbar)^2))/(sw - 1)
}

- Jon Spring

---

不错的解决方案 - 我没有意识到我可以在dplyr中调用Hsmic函数 - 也许这是避免从头开始工作的最佳方法。 - Markm0705