R：向 ggplot2 中的分组直方图添加正态拟合

我正在寻找最优雅的方法来将正态分布拟合叠加到分组直方图中ggplot2。我知道这个问题之前已经被问过很多次了，但是没有一个建议的选项，比如this one https://stackoverflow.com/questions/7182556/how-to-add-gaussian-curve-to-histogram-created-with-qplot or this one https://stackoverflow.com/questions/1376967/using-stat-function-and-facet-wrap-together-in-ggplot2-in-r我觉得非常优雅，至少不是，除非stat_function可以对数据的每个特定部分进行处理。

我遇到的一种将正态分布拟合叠加到非分组直方图上的相对优雅的方法是使用geom_smooth and method="nls"（除了它不是自启动函数并且必须指定起始值这一事实之外）：

library(ggplot2)
myhist = data.frame(size = 10:27, counts = c(1L, 3L, 5L, 6L, 9L, 14L, 13L, 23L, 31L, 40L, 42L, 22L, 14L, 7L, 4L, 2L, 2L, 1L) )
ggplot(data=myhist, aes(x=size, y=counts)) + geom_point() + 
     geom_smooth(method="nls", formula = y ~ N * dnorm(x, m, s), se=F, 
                 start=list(m=20, s=5, N=300)) 

我想知道这种方法是否也可以用于将正态分布拟合添加到分组直方图，如

library(devtools)
install_github("tomwenseleers/easyGgplot2",type="source")
library("easyGgplot2") # load weight data
ggplot(weight,aes(x = weight)) + 
+     geom_histogram(aes(y = ..count.., colour=sex, fill=sex),alpha=0.5,position="identity")

我还想知道是否有任何包可以定义+ stat_distrfit() or + stat_normfit() for ggplot2有没有机会（有分组的可能性）？ （我真的找不到任何东西，但这似乎是一个很常见的任务，所以我只是想知道）

我希望代码尽可能短的原因是这是一个课程，并且我想让事情尽可能简单......

PS geom_density不适合我的目标，我还想绘制计数/频率而不是密度。我还想将它们放在同一个面板中，并避免使用facet_wrap

像这样？

## simulate your dataset - could not get easyGplot2 to load....
set.seed(1)     # for reproducible example
weight <- data.frame(sex=c("Female","Male"), weight=rnorm(1000,mean=c(65,67),sd=1))

library(ggplot2)
library(MASS)       # for fitdistr(...)
get.params <- function(z) with(fitdistr(z,"normal"),estimate[1:2])
df <- aggregate(weight~sex, weight, get.params)
df <- data.frame(sex=df[,1],df[,2])
x  <- with(weight, seq(min(weight),max(weight),len=100))
gg <- data.frame(weight=rep(x,nrow(df)),df)
gg$y <- with(gg,dnorm(x,mean,sd))
gg$y <- gg$y * aggregate(weight~sex, weight,length)$weight * diff(range(weight$weight))/30

ggplot(weight,aes(x = weight, colour=sex)) + 
  geom_histogram(aes(y = ..count.., fill=sex), alpha=0.5,position="identity") +
  geom_line(data=gg, aes(y=y))  

我想“优雅”是情人眼里出西施。使用时出现的问题stat_function(...)那是args=...无法使用列表进行映射aes(...)，正如评论中的帖子所解释的那样。所以你必须创建一个辅助data.frame（gg在此示例中）具有拟合分布的 x 和 y 值，并使用geom_line(...).

上面的代码使用fitdistr(...) in the MASS包来计算数据的平均值和标准差的最大似然估计，根据正态性假设按性别分组（如果有意义的话，您可以使用不同的分布）。然后，它通过将范围划分来创建 x 轴weight分成 100 个增量，并计算dnorm(x,...)以获得适当的均值和标准差。由于结果是密度，我们必须使用以下方法进行调整：

gg$y <- gg$y * aggregate(weight~sex, weight,length)$weight * diff(range(weight$weight))/30

因为您想将其映射到计数数据。请注意，这假设您使用 geom_histogram 中的默认分箱（它将 x 中的范围划分为 30 个相等的增量）。最后，我们添加调用geom_line(...) using gg作为特定于图层的数据集。