在R中自行实现稳定分布极大似然拟合：FFT与对数似然函数问题

核心问题说明 #

稳定分布没有闭合形式的概率密度函数（PDF），没法像t分布那样直接写出对数似然表达式。要构建极大似然估计的目标函数，必须通过特征函数的傅里叶逆变换数值计算PDF，而FFT是高效完成这一逆变换的工具。

稳定分布特征函数（对应你提供的形式） #

当$\alpha \neq 1$时：
$$
\phi(t) = \exp\left{ i\mu t - \sigma^\alpha |t|^\alpha \left(1 + i\beta \text{sign}(t) \tan\left(\frac{\pi\alpha}{2}\right)\right) \right}
$$
当$\alpha = 1$时（柯西分布特例），特征函数调整为：
$$
\phi(t) = \exp\left{ i\mu t - \sigma |t| \left(1 + i\beta \text{sign}(t) \frac{2}{\pi}\log|t|\right) \right}
$$
其中四个参数的取值范围：

• $\alpha$：特征指数，$(0,2]$，控制尾部厚度（α=2时为正态分布）
• $\beta$：偏度参数，$[-1,1]$，控制分布对称性
• $\sigma$：尺度参数，$(0,+\infty)$
• $\mu$：位置参数，$(-\infty,+\infty)$

R代码实现步骤 #

1. 定义稳定分布的特征函数 #

stable_char_func <- function(t, alpha, beta, sigma, mu) {
  if (alpha == 1) {
    # alpha=1时的特征函数分支
    term <- -sigma * abs(t) * (1 + 1i * beta * (2/pi) * sign(t) * log(abs(t)))
  } else {
    # alpha≠1时的特征函数分支
    term <- -sigma^alpha * abs(t)^alpha * (1 + 1i * beta * sign(t) * tan(pi * alpha / 2))
  }
  exp(1i * mu * t + term)
}

2. 用FFT反演得到PDF的数值近似 #

通过FFT快速计算傅里叶逆变换，得到PDF在采样点上的值，再通过插值匹配到你的数据点：

stable_pdf_fft <- function(x, alpha, beta, sigma, mu, N = 2^14, L = 10) {
  # N: FFT采样点数（选2的幂次计算更快），L: 积分区间半长，控制计算范围
  dt <- 2 * L / N
  t <- seq(-L, L - dt, by = dt)
  
  # 计算特征函数在所有t点的取值
  cf <- stable_char_func(t, alpha, beta, sigma, mu)
  
  # 执行逆FFT并归一化，得到PDF的数值近似
  fft_vals <- fft(Conj(cf), inverse = TRUE) / (2 * pi)
  
  # 生成FFT对应的x轴采样点
  dx <- 2 * pi / (N * dt)
  x_fft <- seq(-pi / dt, pi / dt - dx, by = dx)
  
  # 用线性插值得到输入数据x对应的PDF值
  approx(x_fft, Re(fft_vals), xout = x)$y
}

3. 构建负对数似然函数 #

和你实现t分布的逻辑一致，只是把PDF的计算换成FFT反演的结果：

nll_stable <- function(pars, xs) {
  alpha <- pars[1]
  beta <- pars[2]
  sigma <- pars[3]
  mu <- pars[4]
  
  # 参数合法性校验，非法参数返回极大值让优化器避开
  if (alpha <= 0 || alpha > 2 || abs(beta) > 1 || sigma <= 0) {
    return(Inf)
  }
  
  # 计算每个数据点的PDF值
  pdf_vals <- stable_pdf_fft(xs, alpha, beta, sigma, mu)
  
  # 处理极小值，避免log(0)产生NaN或无穷大
  pdf_vals[pdf_vals <= 1e-15] <- 1e-15
  
  # 返回负对数似然
  -sum(log(pdf_vals))
}

4. 调用优化器进行参数估计 #

使用L-BFGS-B方法，因为它支持参数上下界约束，适配稳定分布的参数范围：

# 替换为你的数据
xs <- dx$DAX1$Tagesrendite[2:441]

# 初始参数建议：alpha选1.5（介于正态和柯西之间），beta=0（先假设对称），sigma和mu参考数据的标准差和均值
init_pars <- c(alpha = 1.5, beta = 0, sigma = sd(xs), mu = mean(xs))

# 执行极大似然估计
mle_stable <- optim(init_pars, fn = nll_stable, method = "L-BFGS-B",
                    lower = c(1e-5, -1, 1e-5, -Inf),  # 参数下界
                    upper = c(2, 1, Inf, Inf),        # 参数上界
                    xs = xs)

# 查看估计结果
print(mle_stable)

关键注意事项 #

• FFT参数调整：N和L决定PDF的计算精度，N越大、L越大，精度越高但计算越慢，一般2^14（16384）的采样点数足够应对大多数金融数据场景。
• 参数约束：必须严格限制参数的取值范围，否则特征函数会出现无意义的计算，优化器也容易陷入局部最优。
• 数值稳定性：PDF值可能出现极小值，必须替换为一个很小的正数（如1e-15），避免对数计算出错。

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