【R语言】二项分布，正态分布，极大似然估计实现

二项分布

生成二项分布概率

s <- 0:60
prob <- dbinom(s, size = 60, prob = 1/6)

• s <- 0:60：生成 0 到 60 之间的整数，表示可能的成功次数。

• dbinom(s, size = 60, prob = 1/6)
  

  • dbinom(x, size, prob) 计算二项分布的概率质量函数（PMF），即：
    $$P(X=s)=\left(\genfrac{}{}{0ex}{}{60}{s}\right)(1/6{)}^s(5/6{)}^{60-s}$$

  • size = 60：表示总试验次数（即 60 次）。

  • prob = 1/6：每次试验的成功概率（如掷骰子得到 1 的概率为 1/6）。

---

2. 用 Base R 绘图

plot(s, prob, type = "h")

• plot()：绘制散点图。
• type = "h"：表示绘制竖直线条，类似于直方图。

---

3. 用 ggplot2 绘图

dat <- data.frame(x = s, y = prob)
ggplot(dat, aes(x = s, y = y, xend = s, yend = 0)) +
  geom_segment() + 
  theme_minimal()

代码解析

• data.frame(x = s, y = prob)：创建数据框，包含 x 轴的成功次数 s 和 y 轴的概率 prob。

• ggplot(dat, aes(x = s, y = y, xend = s, yend = 0))
  

  • x = s, y = y：将 s 作为横坐标，prob 作为纵坐标。
  • xend = s, yend = 0：绘制从 (s, prob) 到 (s, 0) 的竖直线段。

• geom_segment()：绘制竖线（类似于 Base R 的 type = "h"）。

• theme_minimal()：应用简洁主题，去掉多余背景。

正态分布概率密度计算

dnorm(10, 10, 0.1)

• 用于计算

  正态分布的概率密度函数

  • x 是你想要计算的点。
  • mean 是正态分布的均值（平均值）。
  • sd 是正态分布的标准差。

所以，dnorm(10, 10, 0.1) 计算的是在均值为 10，标准差为 0.1 的正态分布下，点 x = 10 的概率密度值。

数学公式：

正态分布的概率密度函数（PDF）公式为：
$$f(x)=\frac{1}{\sqrt{2\pi {\sigma }^2}}\exp \left(-\frac{(x-\mu {)}^2}{2{\sigma }^2}\right)$$
其中：

• x 是你要求概率密度的值。
• μ 是正态分布的均值。
• σ 是标准差。

3.989423

pdf是可以大于1的，因为是概率密度而不是概率

生成正态分布密度范围

x <- seq(from = mu - 4 * sig, to = mu + 4 * sig,
         length.out = 128)

• mu 代表正态分布的均值（mean）。

• sig 代表标准差（standard deviation）。

• seq(from = ..., to = ..., length.out = ...)
  

  是 R 中用于生成数值序列的函数。

  • from 是序列的起始值。
  • to 是序列的结束值。
  • length.out 指定生成的序列包含的元素个数。

功能：

• 该代码生成了一个从 mu - 4 * sig 到 mu + 4 * sig 的序列。通常，正态分布的95%数据都落在均值 ± 2倍标准差的范围内，所以 mu ± 4 * sig 给出了一个比 95% 范围更广的范围（约 99.99%）。
• length.out = 128 指定生成的序列包含128个点，这些点均匀分布在这个范围内。

极大似然估计二项分布

n <- 50  
y <- 10  

neg_log_likelihood <- function(p) {
  -dbinom(y, n, p, log = TRUE)
}

initial_guess <- 0.5

result <- optim(initial_guess, neg_log_likelihood, method = "Brent", lower = 0, upper = 1)

mle_proportion_boys <- result$par
mle_proportion_boys

详细解释：

1. 给定的数据：
   • n <- 50：总共50个人。
   • y <- 10：其中10个是男孩。
2. 负对数似然函数的定义：
   • neg_log_likelihood：该函数定义了给定男孩比例 ppp 下的负对数似然函数，表示为 −log⁡(P(y boys))，其中 P(y boys) 是二项分布的概率质量函数 dbinom(y, n, p)。
   • dbinom(y, n, p, log = TRUE) 返回的是概率值的对数。
3. 优化：
   • 使用 optim() 函数进行最大化，初始猜测是 initial_guess <- 0.5，即假设男孩的比例是50%。
   • optim() 方法使用的是 “Brent” 方法，这是一个适合一维问题的优化算法，lower = 0 和 upper = 1 限制了男孩比例 ppp 的取值范围在 [0, 1] 之间。
4. 得到的结果：
   • result$par 提供了使负对数似然最小的比例，这就是我们想要的最大似然估计值（MLE）。

核密度估计

density.cbd <- density(x = distance.cbd, bw = 1.25, from = 0, to = 50)

代码解释：

1. density(x = distance.cbd, bw = 1.25, from = 0, to = 50)
   • x = distance.cbd：distance.cbd 是输入数据，可能是一个包含观测值的向量。
   • bw = 1.25：带宽（bandwidth）是KDE中的一个重要参数，决定了估计曲线的平滑度。较小的带宽会使曲线更接近数据点，而较大的带宽会使曲线更加平滑。bw = 1.25 指定了带宽为1.25。
   • from = 0 和 to = 50：这些参数定义了估计的范围，即从0到50的区间内进行核密度估计。