首页 > 解决方案 > 将正态分布拟合到分组数据,给出预期频率

问题描述

我有观察的频率分布,分组为类间隔内的计数。我想拟合正态(或其他连续)分布,并根据该分布在每个区间中找到预期频率。

例如,假设以下,我想计算另一列,expected给出胸围在 给出的区间内的士兵的预期数量chest,其中假设这些以名义值为中心。例如,35= 34.5 <= y < 35.5。我看到的一项分析给出了这个单元格中的预期频率为 72.5 与观察到的 81。

> data(ChestSizes, package="HistData")
> 
> ChestSizes
   chest count
1     33     3
2     34    18
3     35    81
4     36   185
5     37   420
6     38   749
7     39  1073
8     40  1079
9     41   934
10    42   658
11    43   370
12    44    92
13    45    50
14    46    21
15    47     4
16    48     1
> 

> # ungroup to a vector of values
> chests <- vcdExtra::expand.dft(ChestSizes, freq="count")

这个问题有很多变体,其中大部分与在直方图顶部绘制正常密度有关,按比例缩放以表示计数而不是密度。但没有一个明确显示预期频率的计算。一个很接近的问题是R: add normal fit to grouped histograms in ggplot2

我可以很好地完成标准图(下图),但对于其他事情,比如卡方检验或vcd::rootogram图,我需要相同类别间隔中的预期频率。

> bw <- 1
n_obs <- nrow(chests)
xbar <- mean(chests$chest)
std <- sd(chests$chest)

plt <-
ggplot(chests, aes(chest))  + 
  geom_histogram(color="black", fill="lightblue",  binwidth = bw) + 
  stat_function(fun = function(x) 
    dnorm(x, mean = xbar, sd = std) * bw * n_obs,
    color = "darkred", size = 1)

plt

在此处输入图像描述

标签: rggplot2histogramnormal-distribution

解决方案


以下是假设正态性的情况下如何计算每个组的预期频率。

xbar <- with(ChestSizes, weighted.mean(chest, count))
sdx <- with(ChestSizes, sd(rep(chest, count)))
transform(ChestSizes, Expected = diff(pnorm(c(32, chest) + .5, xbar, sdx)) * sum(count))

   chest count     Expected
1     33     3    4.7600583
2     34    18   20.8822328
3     35    81   72.5129162
4     36   185  199.3338028
5     37   420  433.8292832
6     38   749  747.5926687
7     39  1073 1020.1058521
8     40  1079 1102.2356155
9     41   934  943.0970605
10    42   658  638.9745241
11    43   370  342.7971793
12    44    92  145.6089948
13    45    50   48.9662992
14    46    21   13.0351612
15    47     4    2.7465640
16    48     1    0.4579888

推荐阅读