如何使用ggplot在直方图上叠加任意参数分布?

9

如何使用ggplot对直方图进行任意参数分布的叠加?

我根据Quick-R样例尝试了一下,但我不理解其中的缩放因子从何而来。这种方法合理吗?如何修改以便在ggplot中使用?

以下是使用此方法叠加正态分布和对数正态分布的示例:

## Get a log-normalish data set: the number of characters per word in "Alice in Wonderland"
alice.raw <- readLines(con = "http://www.gutenberg.org/cache/epub/11/pg11.txt", 
                       n = -1L, ok = TRUE, warn = TRUE,
                       encoding = "UTF-8")

alice.long <- paste(alice.raw, collapse=" ")
alice.long.noboilerplate <- strsplit(alice.long, split="\\*\\*\\*")[[1]][3]
alice.words <- strsplit(alice.long.noboilerplate, "[[:space:]]+")[[1]]
alice.nchar <- nchar(alice.words)
alice.nchar <- alice.nchar[alice.nchar > 0]

# Now we want to plot both the histogram and then log-normal probability dist
require(MASS)
h <- hist(alice.nchar, breaks=1:50, xlab="Characters in word", main="Count")
xfit <- seq(1, 50, 0.1)

# Plot a normal curve
yfit<-dnorm(xfit,mean=mean(alice.nchar),sd=sd(alice.nchar))
yfit <- yfit * diff(h$mids[1:2]) * length(alice.nchar) 
lines(xfit, yfit, col="blue", lwd=2) 

# Now plot a log-normal curve
params <- fitdistr(alice.nchar, densfun="lognormal")
yfit <- dlnorm(xfit, meanlog=params$estimate[1], sdlog=params$estimate[1])
yfit <- yfit * diff(h$mids[1:2]) * length(alice.nchar)
lines(xfit, yfit, col="red", lwd=2)

这将产生以下图表: 由上述代码生成的图表,显示单词长度的直方图,叠加了正态分布曲线和对数正态分布曲线 为了澄清,我希望y轴上有计数,而不是密度估计。

请注意,正态分布在单词长度都大于0且值为离散整数的情况下是没有意义的,因为正态分布是连续的。 - David LeBauer
同意 - 这只是一个玩具示例,带有方便的数据集。正态曲线可能不太合适。 - fmark
1个回答

12

请查看stat_function()函数

alice.raw <- readLines(con = "http://www.gutenberg.org/cache/epub/11/pg11.txt", 
                       n = -1L, ok = TRUE, warn = TRUE,
                       encoding = "UTF-8")

alice.long <- paste(alice.raw, collapse=" ")
alice.long.noboilerplate <- strsplit(alice.long, split="\\*\\*\\*")[[1]][3]
alice.words <- strsplit(alice.long.noboilerplate, "[[:space:]]+")[[1]]
alice.nchar <- nchar(alice.words)
alice.nchar <- alice.nchar[alice.nchar > 0]
dataset <- data.frame(alice.nchar = alice.nchar)
library(ggplot2)
ggplot(dataset, aes(x = alice.nchar)) + geom_histogram(aes(y = ..density..)) +
  stat_function(fun = dnorm, 
    args = c(
      mean = mean(dataset$alice.nchar), 
      sd = sd(dataset$alice.nchar)), 
    colour = "red")

enter image description here

如果您想在y轴上像示例一样显示计数,则需要编写一个将密度转换为计数的函数:
dnorm.count <- function(x, mean = 0, sd = 1, log = FALSE, n = 1, binwidth = 1){
  n * binwidth * dnorm(x = x, mean = mean, sd = sd, log = log) 
}

ggplot(dataset, aes(x = alice.nchar)) + geom_histogram(binwidth=1.6) + 
  stat_function(fun = dnorm.count, 
                args = c(
                  mean = mean(dataset$alice.nchar), 
                  sd = sd(dataset$alice.nchar), 
                  n = nrow(dataset), binwidth=1.6), 
                colour = "red")

enter image description here


不错。我认为stat_function一定是新的。这是对我的以前方法的巨大改进,以前我是先创建一个x、dnorm(x, , )数据框架。 - David LeBauer
1
@David stat_function 已经存在了很久,至少我还记得! :) - joran
这很棒 - 可以将y轴上的计数显示出来,而不是像上面的例子中一样显示密度吗? - fmark
@fmark:你可以这样做。你需要一个将密度转换为计数的函数。 - Thierry
@Thierry 谢谢你的回答。我稍微修改了你的答案,因为你需要在转换函数中包含 binwidth,并添加了结果图。 - fmark

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接