对数正态蒙特卡洛模拟样本标准差偏低问题咨询

你遇到的这个问题其实挺常见的——对数正态分布蒙特卡洛模拟里样本标准差持续低估总体标准差，核心是两个因素叠加的结果：样本标准差本身的固有偏差，加上对数正态分布的强右偏特性放大了这个偏差，小样本下这种偏差会显得格外明显。

为什么会出现这个偏差？ #

1. 样本标准差的固有有偏性 #

不管是什么分布，样本标准差（除以n-1的那个）都是总体标准差的有偏估计。它的期望和总体标准差之间存在一个缩放系数：
E[s] = σ * sqrt(2/(n-1)) * Γ(n/2) / Γ((n-1)/2)
这个系数在n>1时始终小于1，意味着样本标准差的平均水平会天然低于总体标准差。比如：

• 当样本量n=5时，系数约为0.94，平均偏差约6%
• 当n=20时，系数约0.98，偏差约2%
• 只有当样本量极大时，这个系数才会趋近于1，偏差可以忽略。

2. 对数正态的右偏特性放大偏差 #

对数正态分布是典型的强右偏分布，总体中存在一些极端大值，但小样本中很难抽到这些值。这就导致样本的离散程度被“低估”，样本标准差比总体标准差低的情况更普遍，偏差也会比正态分布下更显著——这就是你感觉“偏差远超预期”的原因。

而样本均值没有这个问题，因为样本均值是总体均值的无偏估计，无论分布是否正态，只要样本独立同分布，大数定律就会保证样本均值收敛到总体均值。

用R代码验证这个偏差 #

我们可以通过蒙特卡洛模拟直观看到这个现象：

library(dplyr)

# 设置对数正态的参数：对数均值μ=0，对数标准差σ=1
mu_log <- 0
sigma_log <- 1
# 计算真实总体均值和标准差
pop_mean <- exp(mu_log + sigma_log^2 / 2)
pop_sd <- sqrt((exp(sigma_log^2) - 1) * exp(2*mu_log + sigma_log^2))

# 模拟设置：样本量n=20，模拟10000次
n <- 20
sim_times <- 10000

# 执行模拟
sim_data <- replicate(sim_times, {
  sample <- rlnorm(n, meanlog = mu_log, sdlog = sigma_log)
  tibble(
    sample_mean = mean(sample),
    sample_sd = sd(sample)
  )
}) %>% bind_rows()

# 查看模拟结果的平均水平
sim_summary <- sim_data %>%
  summarize(
    平均样本均值 = mean(sample_mean),
    真实总体均值 = pop_mean,
    平均样本标准差 = mean(sample_sd),
    真实总体标准差 = pop_sd
  )

print(sim_summary)

运行后你会发现：

• 平均样本均值几乎和真实总体均值重合
• 平均样本标准差明显低于真实总体标准差（比如n=20时，平均样本标准差约2.1，真实值约2.197，偏差约4%）

如何修正这个偏差？ #

如果你需要得到总体标准差的无偏估计，可以试试这几种方法：

• 增大样本量：这是最简单的方法，随着样本量增加，样本标准差的期望会快速趋近于总体标准差，偏差会大幅缩小。
• 用修正系数调整样本标准差：利用前面提到的缩放系数，对样本标准差进行修正：

  # 计算修正系数
  correction_factor <- gamma(n/2) / (gamma((n-1)/2) * sqrt(2/(n-1)))
  # 修正后的样本标准差
  sim_data <- sim_data %>%
    mutate(修正后样本标准差 = sample_sd * correction_factor)
  

  修正后的结果平均值会非常接近真实总体标准差。
• 转换到正态尺度计算：既然对数正态的对数是正态分布，我们可以先对样本取对数，计算对数样本的标准差（它是对数尺度标准差的无偏估计），再转换回原始尺度：

  # 以单次样本为例
  sample <- rlnorm(n, meanlog = mu_log, sdlog = sigma_log)
  log_sd <- sd(log(sample))
  log_mean <- mean(log(sample))
  # 转换为原始尺度的总体标准差估计
  estimated_pop_sd <- sqrt((exp(log_sd^2) - 1) * exp(2*log_mean + log_sd^2))
  

  这种方法利用了正态分布的良好统计性质，可靠性更高。

内容的提问来源于stack exchange，提问作者apax