量化交易中文教材

第 09 章 参数推断

本章对应 Wasserman 原书第 9 章。前面几章(经验分布、统计泛函、Bootstrap)走的是非参数路线:几乎不对分布做假设。本章转向参数模型:假设数据来自一个由有限个参数刻画的分布族,然后估计这几个参数。参数推断的主角是最大似然估计(MLE),它是金融计量里 GARCH、ARMA、状态空间模型、跳跃扩散、Hawkes 过程等几乎所有模型的估计方法。本章把 MLE 的来龙去脉讲透:怎么求、为什么好、标准误怎么算、什么时候会坏、算不出闭式解时怎么数值求解。

学习目标

  1. 会对常见分布写出似然函数并求 MLE,能识别"支撑依赖参数"这类不能求导的情形(如均匀分布)。
  2. 理解 MLE 的五条性质:相合、同变、渐近正态、渐近有效、近似贝叶斯;能用 KL 距离解释相合性。
  3. 会计算得分函数与 Fisher 信息,用 \(\widehat{\text{se}}=1/\sqrt{I_n(\hat\theta)}\) 构造 Wald 型置信区间;多参数时会用 Fisher 信息矩阵的逆和多元 Delta 方法。
  4. 会用参数 Bootstrap 求任意参数函数的标准误,知道它与 Delta 方法各自的优劣。
  5. 理解充分统计量、因子分解定理、Rao–Blackwell 定理和指数族的基本结论。
  6. 会用 Newton–Raphson(拟牛顿)与 EM 算法数值求 MLE,并能把 EM 用于两状态收益混合模型。

读前导读

这一章在解决什么问题。 你在 CFA 二级学过 GARCH 的结论,也在软件里见过"Maximum Likelihood"字样,但多半没追问过:软件到底在最大化什么,报告出来的标准误从哪里来,为什么能信。这一章回答这三件事。最大似然的思想很朴素:在所有候选参数里,挑一个"让已经发生的这批数据看起来最不意外"的参数。例如观察到 400 笔交易赚了 220 笔,胜率 0.55 让这个结果最可能发生,所以 MLE 就是 0.55。

本章最重要的结论是:MLE 的标准误由对数似然函数在峰顶的弯曲程度决定。峰越尖,说明参数稍微偏离一点数据就"很不满意",参数就被定位得越准。这和你熟悉的凸性是同一个数学对象——二阶导数。债券的凸性衡量价格曲线有多弯;Fisher 信息衡量对数似然曲线有多弯。

需要先想起来的数学。

  • 对数与求导。 \(\log(ab)=\log a+\log b\),所以乘积形式的似然取对数后变成求和;\((\log u)'=u'/u\);最大值点处一阶导为 0、二阶导为负。例:\(\ell(p)=220\log p+180\log(1-p)\),\(\ell'(p)=220/p-180/(1-p)=0\) 解出 \(p=0.55\)。见 第 00 册第 02 章 导数与泰勒展开。
  • 一阶泰勒展开用于求根。 \(h(\hat\theta)\approx h(\theta)+h'(\theta)(\hat\theta-\theta)\);令左边为 0 就能解出 \(\hat\theta-\theta\)。渐近正态的证明和 Newton–Raphson 迭代都是这一步。你用过的"根据久期反推收益率变动"是同一个操作。
  • 偏导、Hessian 与矩阵求逆。 多个参数时,二阶导排成矩阵(Hessian,\(H_{jk}=\partial^2\ell/\partial\theta_j\partial\theta_k\)),"除以二阶导"变成"乘以 Hessian 的逆矩阵"。2×2 矩阵的逆:\(\begin{pmatrix}a&b\\b&d\end{pmatrix}^{-1}=\frac{1}{ad-b^2}\begin{pmatrix}d&-b\\-b&a\end{pmatrix}\)。见 第 00 册第 05 章 与 第 06 章 线性代数速成。
  • 积分号下求导。 \(\frac{d}{d\theta}\int f(x;\theta)\,dx=\int\frac{\partial f}{\partial\theta}dx\),在"正则条件"下成立。得分均值为零的引理只靠这一步。见 第 00 册第 03 章 积分。
  • \(\arg\max\) 记号。 \(\arg\max_\theta J(\theta)\) 是"使 \(J\) 最大的那个 \(\theta\)",不是最大值本身。见 第 00 册第 08 章 读懂数学证明与符号。

怎么读这一章。 核心必读:9.3(似然与 MLE 的求法)、9.7(得分、Fisher 信息、渐近正态,这是本章的心脏)、9.9–9.10(Delta 方法与多参数情形)、9.11(参数 bootstrap)。9.5 的相合性证明和 9.13–9.14 的充分统计量、指数族属于原书附录,第一次可以只看直觉和结论。9.2 矩估计读一遍即可,记住它主要用于提供初值。9.15 的 EM 对处理市场状态模型很实用,可以和 9.16.2 的代码对照读。


9.1 参数模型与"感兴趣参数"

参数模型是一族分布

\[\mathfrak F=\{f(x;\theta):\theta\in\Theta\},\qquad\Theta\subset\mathbb R^k,\ \theta=(\theta_1,\dots,\theta_k).\]
推断问题于是化为:根据 IID 样本 \(X_1,\dots,X_n\) 估计 \(\theta\)。

学生常问:我怎么知道数据真的来自某个参数模型?Wasserman 的回答很坦率——我们很少能知道,所以非参数方法往往更可取。但学参数方法仍有两个理由:一是背景知识有时表明某个参数模型是合理近似(例如交通事故计数近似服从 Poisson);二是参数模型里的推断概念(似然、Fisher 信息、渐近正态)是理解很多非参数方法的基础。对量化来说还有第三个理由:在数据有限的情况下,一个结构合理的参数模型(比如带厚尾的 t 分布、GARCH)能比纯非参数方法更稳定地外推到尾部。

感兴趣参数与冗余参数。 我们常常只关心 \(\theta\) 的某个函数 \(T(\theta)\)。例如 \(X\sim N(\mu,\sigma^2)\) 而目标是 \(\mu\),则 \(\mu\) 是感兴趣参数(parameter of interest),\(\sigma\) 是冗余参数(nuisance parameter)。

例 9.1(原书)。 血液检测值 \(X_i\sim N(\mu,\sigma^2)\),关心检测值大于 1 的人群比例

\[\tau=\mathbb P(X>1)=1-\Phi\!\left(\frac{1-\mu}{\sigma}\right).\]
把"1"换成"−3%",这就是"日收益跌破 −3% 的概率",和 VaR 是同一类感兴趣参数。

例 9.2(原书)。 Gamma\((\alpha,\beta)\) 分布 \(f(x;\alpha,\beta)=\frac{1}{\beta^\alpha\Gamma(\alpha)}x^{\alpha-1}e^{-x/\beta}\) 常用来描述寿命,平均寿命 \(T(\alpha,\beta)=\alpha\beta\) 是感兴趣参数。


9.2 矩估计法

直觉。 总体矩是参数的函数;样本矩是总体矩的相合估计。让二者相等,解出参数。

定义第 \(j\) 阶总体矩与样本矩:

\[\alpha_j(\theta)=\mathbb E_\theta(X^j)=\int x^j\,dF_\theta(x),\qquad\hat\alpha_j=\frac1n\sum_{i=1}^nX_i^j.\]

定义 9.3(矩估计量,method of moments estimator)。 \(\hat\theta_n\) 是使

\[\alpha_j(\hat\theta_n)=\hat\alpha_j,\qquad j=1,\dots,k\]
成立的 \(\theta\)——\(k\) 个方程解 \(k\) 个未知数。

例 9.4。 Bernoulli\((p)\):\(\alpha_1=p\),故 \(\hat p_n=\bar X_n\)。

例 9.5。 \(N(\mu,\sigma^2)\):\(\alpha_1=\mu\),\(\alpha_2=\sigma^2+\mu^2\)。解方程得

\[\hat\mu=\bar X_n,\qquad\hat\sigma^2=\frac1n\sum(X_i-\bar X_n)^2.\]

定理 9.6。 在适当条件下,矩估计量 (1) 以趋于 1 的概率存在;(2) 相合;(3) 渐近正态:\(\sqrt n(\hat\theta_n-\theta)\rightsquigarrow N(0,\Sigma)\),其中 \(\Sigma=g\,\mathbb E_\theta(YY^T)g^T\),\(Y=(X,X^2,\dots,X^k)^T\),\(g_j=\partial\alpha_j^{-1}(\theta)/\partial\theta\)。

第 (3) 条可以用来算标准误,但实践中用 Bootstrap 更省事。矩估计一般不是最优的(方差比 MLE 大),但它计算简单,最大的用处是给数值 MLE 提供初值——后面的量化实战就是这样做的:先用样本峰度反推 t 分布的自由度,再交给优化器。

量化旁注:金融计量里的广义矩估计(GMM)是矩估计的推广,方程个数可以多于参数个数,用加权最小二乘去"尽量满足"所有矩条件。资产定价模型检验(随机贴现因子的欧拉方程)大量使用 GMM,详见第 05 册与第 06 册相关章节。


9.3 最大似然

定义 9.7。 IID 样本、密度 \(f(x;\theta)\) 时,似然函数(likelihood function)与对数似然(log-likelihood)为

\[\mathcal L_n(\theta)=\prod_{i=1}^nf(X_i;\theta),\qquad\ell_n(\theta)=\log\mathcal L_n(\theta).\]

似然就是数据的联合密度,只不过把数据固定、把它看作参数的函数。务必记住:似然不是 \(\theta\) 的概率密度,对 \(\theta\) 积分一般不等于 1。把似然当成"参数的概率"是初学者最常见的误解,第 11 章的贝叶斯推断会说明,要得到参数的概率分布必须再乘上先验。

定义 9.8。 使 \(\mathcal L_n(\theta)\) 最大的 \(\hat\theta_n\) 称为最大似然估计量(maximum likelihood estimator, MLE)。

白话解释:似然回答的问题是"如果参数是 \(\theta\),看到手上这批数据的可能性有多大"。固定数据、让 \(\theta\) 变化,就得到一条曲线;MLE 取曲线的最高点。它和"参数的概率"不同:问"胜率为 0.55 时看到 220/400 的可能性"是似然;问"看到 220/400 后胜率为 0.55 的可能性"才是关于参数的概率,后者需要先验。 和你熟悉的隐含波动率对比:隐含波动率是"让模型价格等于市场价格"的那个 \(\sigma\);MLE 是"让模型对整批数据的拟合度(联合密度)最大"的那个参数。两者都是"反推参数",只是匹配的标准不同。

\(\ell_n\) 与 \(\mathcal L_n\) 的最大值点相同,实践中总是对数似然更好处理(乘积变求和,数值上也不会下溢)。另外,\(\mathcal L_n\) 乘以与 \(\theta\) 无关的正常数不改变 MLE,所以写似然时可以随意丢掉常数(注 9.9)。

例 9.10(Bernoulli)。 \(f(x;p)=p^x(1-p)^{1-x}\),记 \(S=\sum X_i\),

\[\mathcal L_n(p)=p^S(1-p)^{n-S},\qquad\ell_n(p)=S\log p+(n-S)\log(1-p).\]
令导数为 0 得 \(\hat p_n=S/n\)。原书图 9.1 画出 \(n=20\)、\(S=12\) 时的似然,最大值在 \(\hat p=0.6\)。

例 9.11(正态)。 \(\theta=(\mu,\sigma)\)。利用分解 \(\sum(X_i-\mu)^2=nS^2+n(\bar X-\mu)^2\)(把 \(X_i-\mu\) 写成 \((X_i-\bar X)+(\bar X-\mu)\) 展开,交叉项为 0),其中 \(S^2=n^{-1}\sum(X_i-\bar X)^2\),忽略常数:

\[\ell(\mu,\sigma)=-n\log\sigma-\frac{nS^2}{2\sigma^2}-\frac{n(\bar X-\mu)^2}{2\sigma^2}.\]

推导拆解:正态密度取对数是 \(\log f(x_i)=-\log\sigma-\frac12\log(2\pi)-\frac{(x_i-\mu)^2}{2\sigma^2}\),对 \(i\) 求和并丢掉常数 \(-\frac n2\log(2\pi)\),得 \(-n\log\sigma-\frac{1}{2\sigma^2}\sum(X_i-\mu)^2\),再代入原文的平方和分解即得上式。 对 \(\mu\) 求偏导(\(\sigma\) 当常数):只有最后一项含 \(\mu\),\(\frac{\partial\ell}{\partial\mu}=\frac{n(\bar X-\mu)}{\sigma^2}=0\Rightarrow\hat\mu=\bar X\)。 对 \(\sigma\) 求偏导:\(-\frac n\sigma+\frac{nS^2}{\sigma^3}+\frac{n(\bar X-\mu)^2}{\sigma^3}\),代入 \(\mu=\bar X\) 后第三项为 0,令其为 0 得 \(\sigma^2=S^2\)。 交叉项为 0 的原因:\(\sum(X_i-\bar X)(\bar X-\mu)=(\bar X-\mu)\sum(X_i-\bar X)\),而离差之和恒为 0。

解 \(\partial\ell/\partial\mu=0\)、\(\partial\ell/\partial\sigma=0\) 得 \(\hat\mu=\bar X\),\(\hat\sigma=S\)。注意 MLE 的方差估计除以 \(n\) 而不是 \(n-1\),因此有偏,但偏差是 \(O(1/n)\),大样本下无关紧要。

例 9.12(一个让很多人困惑的例子)。 \(X_i\sim\text{Unif}(0,\theta)\),\(f(x;\theta)=1/\theta\),\(0\le x\le\theta\)。如果某个 \(X_i>\theta\),则 \(f(X_i;\theta)=0\),整个似然为 0。所以

\[\mathcal L_n(\theta)=\begin{cases}\theta^{-n},&\theta\ge X_{(n)}\\0,&\theta<X_{(n)}\end{cases}\]
其中 \(X_{(n)}=\max_iX_i\)。似然在 \([X_{(n)},\infty)\) 上严格递减,故 \(\hat\theta_n=X_{(n)}\)。这里似然在最大值点处不连续,求导令为 0 的方法完全失效。教训:分布的支撑依赖参数时,不能机械地求导。这类模型也不满足后文的正则条件,渐近正态性不成立(\(X_{(n)}\) 的收敛速度是 \(1/n\) 而不是 \(1/\sqrt n\),极限分布是指数型的)。

多元正态与多项分布的 MLE 见原书第 14 章(定理 14.5、14.3)。


9.4 MLE 的五条性质(总览)

在模型满足正则条件(regularity conditions,本质上是 \(f(x;\theta)\) 关于 \(\theta\) 足够光滑、支撑不依赖参数等)时,MLE 有以下性质:

  1. 相合(consistent):\(\hat\theta_n\xrightarrow{P}\theta_\star\),\(\theta_\star\) 为真值;
  2. 同变(equivariant):若 \(\hat\theta_n\) 是 \(\theta\) 的 MLE,则 \(g(\hat\theta_n)\) 是 \(g(\theta)\) 的 MLE;
  3. 渐近正态(asymptotically normal):\((\hat\theta_n-\theta_\star)/\widehat{\text{se}}\rightsquigarrow N(0,1)\),且 \(\widehat{\text{se}}\) 常可解析计算;
  4. 渐近有效(asymptotically efficient):粗略地说,在所有表现良好的估计量中,MLE 的渐近方差最小;
  5. 近似为贝叶斯估计量(第 11 章定理 11.5)。

Wasserman 特别强调:在足够复杂的问题中,这些性质不再成立,MLE 也不再是好估计。高维问题就是典型(第 12 章的"多个正态均值"与 Stein 悖论)。下面逐条展开。


9.5 相合性:为什么最大化似然能找到真值

Kullback–Leibler 距离。 两个密度 \(f,g\) 之间

\[D(f,g)=\int f(x)\log\frac{f(x)}{g(x)}\,dx.\]
由 Jensen 不等式可证 \(D(f,g)\ge0\),且 \(D(f,f)=0\)。它不对称,所以不是严格意义上的距离。记 \(D(\theta,\psi)=D(f(\cdot;\theta),f(\cdot;\psi))\)。若 \(\theta\ne\psi\Rightarrow D(\theta,\psi)>0\),即不同参数对应不同分布,称模型可识别(identifiable)。以下总假设模型可识别。

可识别性在量化里不是摆设。两状态混合模型交换两个状态的标签,似然完全不变("标签切换");因子模型中因子载荷乘以可逆矩阵、因子乘以其逆,似然也不变。这些情形必须加约束(例如规定危机状态方差更大)才能谈估计的相合性。

直觉。 最大化 \(\ell_n(\theta)\) 等价于最大化

\[M_n(\theta)=\frac1n\sum_{i=1}^n\log\frac{f(X_i;\theta)}{f(X_i;\theta_\star)},\]
因为两者只差一个与 \(\theta\) 无关的项 \(\ell_n(\theta_\star)/n\)。由大数定律,
\[M_n(\theta)\to\mathbb E_{\theta_\star}\log\frac{f(X;\theta)}{f(X;\theta_\star)}=-D(\theta_\star,\theta)\equiv M(\theta).\]
\(M(\theta)\) 在 \(\theta_\star\) 处取最大值 0,其他地方为负。所以 \(M_n\) 的最大值点应当趋向 \(\theta_\star\)。换句话说:MLE 是在寻找与真实分布 KL 距离最近的模型成员。

白话解释:KL 距离 \(D(f,g)\) 可以理解为"真实分布是 \(f\),你却用 \(g\) 来描述它,平均每个观测要付出多少'意外'"。\(\log\frac{f(x)}{g(x)}\) 是在 \(x\) 处真实密度比模型密度高出的对数倍数,按真实分布加权平均就是 \(D\)。\(D\ge0\) 的证明用 Jensen:\(-D=\mathbb E_f\log\frac gf\le\log\mathbb E_f\frac gf=\log\int g=\log1=0\)(\(\log\) 是凹函数)。 上面那条链的逻辑是:MLE 最大化的样本平均 \(M_n(\theta)\),按大数定律收敛到 \(-D(\theta_\star,\theta)\);后者在真值处最大。于是"最大化样本版本"最终等于"最小化 KL 距离"。

这个视角在模型设定错误时尤其有用:如果真实分布不在模型族里,MLE 会收敛到族里与真分布 KL 距离最小的那个成员("伪真值")。用正态分布去拟合厚尾收益,得到的 \(\hat\sigma\) 仍收敛到真实标准差,但你据此算出的尾部概率是错的——这正是后面量化实战要演示的。

定理 9.13。 令 \(M(\theta)=-D(\theta_\star,\theta)\)。若

\[\sup_{\theta\in\Theta}|M_n(\theta)-M(\theta)|\xrightarrow{P}0\quad\text{(一致收敛)}\]
且对每个 \(\epsilon>0\),
\[\sup_{\theta:|\theta-\theta_\star|\ge\epsilon}M(\theta)<M(\theta_\star)\quad\text{(真值是"分离良好"的最大值点)},\]
则 \(\hat\theta_n\xrightarrow{P}\theta_\star\)。

证明思路。 因为 \(\hat\theta_n\) 最大化 \(M_n\),有 \(M_n(\hat\theta_n)\ge M_n(\theta_\star)\),于是

\[M(\theta_\star)-M(\hat\theta_n)\le[M_n(\hat\theta_n)-M(\hat\theta_n)]+[M(\theta_\star)-M_n(\theta_\star)]\le\sup_\theta|M_n-M|+[M(\theta_\star)-M_n(\theta_\star)]\xrightarrow{P}0.\]
所以 \(M(\hat\theta_n)\) 依概率逼近最大值 \(M(\theta_\star)\)。由分离条件,离 \(\theta_\star\) 超过 \(\epsilon\) 的点 \(M\) 值至少低 \(\delta\),因此 \(\mathbb P(|\hat\theta_n-\theta_\star|>\epsilon)\le\mathbb P(M(\hat\theta_n)<M(\theta_\star)-\delta)\to0\)。

只有逐点的大数定律不够,必须一致收敛——否则 \(M_n\) 可能在远离真值的某处偶然冒出一个尖峰。


9.6 同变性

定理 9.14。 设 \(\tau=g(\theta)\),\(\hat\theta_n\) 是 \(\theta\) 的 MLE,则 \(\hat\tau_n=g(\hat\theta_n)\) 是 \(\tau\) 的 MLE。

证明(\(g\) 可逆时)。 令 \(h=g^{-1}\)。对任意 \(\tau\),以 \(\tau\) 为参数的似然 \(\mathcal L(\tau)=\prod f(x_i;h(\tau))=\mathcal L(\theta)\),\(\theta=h(\tau)\)。于是 \(\mathcal L_n(\tau)=\mathcal L(\theta)\le\mathcal L(\hat\theta)=\mathcal L_n(\hat\tau)\)。

例 9.15。 \(X_i\sim N(\theta,1)\),\(\hat\theta=\bar X_n\),则 \(\tau=e^\theta\) 的 MLE 为 \(e^{\bar X_n}\)。

同变性的实用价值很大:VaR、尾部概率、半衰期、年化波动率,只要是参数的函数,MLE 就是把参数的 MLE 代进去。注意无偏性没有这个性质:\(\bar X\) 对 \(\theta\) 无偏,\(e^{\bar X}\) 对 \(e^\theta\) 却有偏(Jensen 不等式)。


9.7 得分函数、Fisher 信息与渐近正态性

这是本章的核心。

定义 9.16。 得分函数(score function)与 Fisher 信息(Fisher information):

\[s(X;\theta)=\frac{\partial\log f(X;\theta)}{\partial\theta},\qquad I_n(\theta)=\mathbb V_\theta\Big(\sum_{i=1}^ns(X_i;\theta)\Big)=\sum_{i=1}^n\mathbb V_\theta\big(s(X_i;\theta)\big).\]
\(n=1\) 时记 \(I(\theta)\)。

引理(得分均值为零)。 \(\mathbb E_\theta[s(X;\theta)]=0\)。证明:对 \(1=\int f(x;\theta)dx\) 两边关于 \(\theta\) 求导(假设可交换求导与积分),

\[0=\int\frac{\partial f}{\partial\theta}dx=\int\frac{\partial f/\partial\theta}{f}f\,dx=\int\frac{\partial\log f}{\partial\theta}f\,dx=\mathbb E_\theta s(X;\theta).\]
因此 \(\mathbb V_\theta(s)=\mathbb E_\theta(s^2)\)。

定理 9.17。 \(I_n(\theta)=nI(\theta)\),并且

\[I(\theta)=-\mathbb E_\theta\left(\frac{\partial^2\log f(X;\theta)}{\partial\theta^2}\right).\]

推导拆解:从引理 \(\int\frac{\partial\log f}{\partial\theta}f\,dx=0\) 出发,两边再对 \(\theta\) 求导。被积函数是乘积 \(s\cdot f\),用乘积法则: \(0=\int\frac{\partial^2\log f}{\partial\theta^2}f\,dx+\int\frac{\partial\log f}{\partial\theta}\cdot\frac{\partial f}{\partial\theta}dx\)。 第二项里 \(\frac{\partial f}{\partial\theta}=\frac{\partial\log f}{\partial\theta}\cdot f\)(这是 \((\log f)'=f'/f\) 反过来用),所以第二项 \(=\int s^2f\,dx=\mathbb E s^2\)。移项得 \(\mathbb E s^2=-\mathbb E\frac{\partial^2\log f}{\partial\theta^2}\)。 \(I_n=nI\) 则来自独立性:独立变量之和的方差等于方差之和。数据量翻倍,信息量翻倍,标准误缩小到 \(1/\sqrt2\)。

(对引理的等式再求一次导即得。)于是 Fisher 信息有两个等价的面孔:得分的方差,以及对数似然曲率的期望。曲率越大,似然峰越尖,数据对参数的"定位"越精确,估计方差越小。

定理 9.18(MLE 的渐近正态性)。 令 \(\text{se}=\sqrt{\mathbb V(\hat\theta_n)}\)。在正则条件下:

  1. \(\text{se}\approx\sqrt{1/I_n(\theta)}\),且 \(\dfrac{\hat\theta_n-\theta}{\text{se}}\rightsquigarrow N(0,1)\);
  2. 令 \(\widehat{\text{se}}=\sqrt{1/I_n(\hat\theta_n)}\),则 \(\dfrac{\hat\theta_n-\theta}{\widehat{\text{se}}}\rightsquigarrow N(0,1)\)。

证明(原书附录,值得掌握)。 把得分在真值 \(\theta\) 处做一阶展开:

\[0=\ell'(\hat\theta)\approx\ell'(\theta)+(\hat\theta-\theta)\ell''(\theta)\ \Rightarrow\ \sqrt n(\hat\theta-\theta)\approx\frac{\frac{1}{\sqrt n}\ell'(\theta)}{-\frac1n\ell''(\theta)}=\frac{\text{TOP}}{\text{BOTTOM}}.\]
令 \(Y_i=\partial\log f(X_i;\theta)/\partial\theta\),则 \(\mathbb EY_i=0\)、\(\mathbb VY_i=I(\theta)\),由 CLT,\(\text{TOP}=\sqrt n\bar Y\rightsquigarrow N(0,I(\theta))\)。令 \(A_i=-\partial^2\log f(X_i;\theta)/\partial\theta^2\),\(\mathbb EA_i=I(\theta)\),由大数定律 \(\text{BOTTOM}=\bar A\xrightarrow{P}I(\theta)\)。由 Slutsky 定理(本册第 01 章 1.5.4 节),
\[\sqrt n(\hat\theta-\theta)\rightsquigarrow\frac{N(0,I(\theta))}{I(\theta)}=N\Big(0,\frac{1}{I(\theta)}\Big).\]
第 2 条再用一次 Slutsky:\(I(\hat\theta_n)/I(\theta)\xrightarrow{P}1\)。

推导拆解:这个证明只有四个动作。 (1) MLE 满足一阶条件 \(\ell'(\hat\theta)=0\);把 \(\ell'\) 在真值 \(\theta\) 处做一阶泰勒展开,就得到一个关于 \(\hat\theta-\theta\) 的线性方程,解出 \(\hat\theta-\theta\approx-\ell'(\theta)/\ell''(\theta)\)。 (2) 分子分母同乘一个合适的 \(n\) 的幂:分子 \(\ell'(\theta)=\sum Y_i\) 是 \(n\) 个均值为 0 的独立项之和,量级 \(\sqrt n\),所以除以 \(\sqrt n\);分母 \(\ell''(\theta)\) 是 \(n\) 项之和,量级 \(n\),所以除以 \(n\)。两边合起来正好是 \(\sqrt n(\hat\theta-\theta)\)。 (3) 分子用 CLT("零均值项的和除以 \(\sqrt n\)"趋于正态),分母用大数定律("平均值"趋于期望 \(I(\theta)\))。 (4) Slutsky:分子正态、分母趋于常数,比值就是正态除以常数,方差 \(I/I^2=1/I\)。 金融直觉:\(1/I\) 的含义可以用收益率曲线类比:对数似然在峰顶的曲率越大(类似高凸性),参数稍有偏离"拟合度"就急剧下降,数据把参数钉得越牢,标准误越小;曲线平坦时,一大片参数都差不多好,标准误就大。

定理 9.19(渐近置信区间)。

\[C_n=\big(\hat\theta_n-z_{\alpha/2}\widehat{\text{se}},\ \hat\theta_n+z_{\alpha/2}\widehat{\text{se}}\big),\qquad\mathbb P_\theta(\theta\in C_n)\to1-\alpha.\]
\(\alpha=0.05\) 时就是熟悉的 \(\hat\theta_n\pm2\widehat{\text{se}}\)。报纸上"民调误差 ±1 个百分点,95% 可信"说的就是这种区间。

例 9.20(Bernoulli)。 \(\log f=x\log p+(1-x)\log(1-p)\),

\[s(X;p)=\frac Xp-\frac{1-X}{1-p},\qquad-s'(X;p)=\frac X{p^2}+\frac{1-X}{(1-p)^2},\qquad I(p)=\frac1p+\frac1{1-p}=\frac1{p(1-p)}.\]
故 \(\widehat{\text{se}}=\sqrt{\hat p(1-\hat p)/n}\),近似 95% 区间 \(\hat p\pm2\sqrt{\hat p(1-\hat p)/n}\)。量化读法:一个信号做了 400 笔交易、胜率 55%,标准误 \(\sqrt{0.55\times0.45/400}\approx2.5\%\),95% 区间约 \([50\%,60\%]\)——刚好擦着 50% 的边。

例 9.21(正态,\(\sigma\) 已知)。 \(s=(X-\theta)/\sigma^2\),\(I=1/\sigma^2\),\(\bar X\sim N(\theta,\sigma^2/n)\) 精确成立。

例 9.22(Poisson)。 \(\hat\lambda=\bar X\),\(I(\lambda)=1/\lambda\),\(\widehat{\text{se}}=\sqrt{\hat\lambda/n}\)。量化读法:每分钟到达的大单数、每天的跳跃次数都常用 Poisson 建模,其强度的区间就是 \(\hat\lambda\pm z_{\alpha/2}\sqrt{\hat\lambda/n}\)。

观测信息 vs 期望信息。 实践中常用观测 Fisher 信息(observed information)\(-\ell_n''(\hat\theta)\) 代替 \(I_n(\hat\theta)=nI(\hat\theta)\)——优化器输出的负对数似然 Hessian 就是它。二者都给出相合的标准误;复杂模型(GARCH 等)里期望信息算不出来,只能用观测信息。

模型错设时的标准误。 上面的推导假设模型正确。若模型错设(例如用正态似然估计 GARCH,但真实新息厚尾),"TOP 的方差"不再等于"BOTTOM 的极限",此时 \(\sqrt n(\hat\theta-\theta^*)\) 的渐近方差是三明治形式 \(A^{-1}BA^{-1}\),其中 \(A=-\mathbb E\ell_i''\)、\(B=\mathbb V(\ell_i')\)。这就是准极大似然(QMLE)的稳健标准误(Bollerslev–Wooldridge 标准误),arch 包中 cov_type="robust" 即此。原书未讲,但在金融计量中是默认做法。


9.8 最优性与渐近相对效率

正态模型中,样本均值 \(\bar X_n\)(MLE)满足 \(\sqrt n(\bar X_n-\theta)\rightsquigarrow N(0,\sigma^2)\);样本中位数 \(\tilde\theta_n\) 也收敛到真值,但 \(\sqrt n(\tilde\theta_n-\theta)\rightsquigarrow N(0,\sigma^2\pi/2)\),方差更大。

定义(渐近相对效率,ARE)。 若 \(\sqrt n(T_n-\theta)\rightsquigarrow N(0,t^2)\),\(\sqrt n(U_n-\theta)\rightsquigarrow N(0,u^2)\),则 \(\text{ARE}(U,T)=t^2/u^2\)。

正态下 \(\text{ARE}(\text{中位数},\text{均值})=2/\pi\approx0.63\):用中位数相当于只用了约 63% 的数据。

定理 9.23。 若 \(\hat\theta_n\) 是 MLE,\(\tilde\theta_n\) 是任一其他(表现良好的)估计量,则 \(\text{ARE}(\tilde\theta_n,\hat\theta_n)\le1\)。MLE 因此被称为有效的(efficient)或渐近最优的。这与经典的 Cramér–Rao 下界 \(\mathbb V(\tilde\theta)\ge1/I_n(\theta)\)(对无偏估计)是同一件事的两种说法。

关键前提是模型正确。 收益分布厚尾时,"正态 MLE = 样本均值"不再最优;在 t 分布、拉普拉斯分布下,中位数或截尾均值的效率可以超过均值。这就是为什么稳健统计量(中位数、截尾均值、Huber 估计)在估计期望收益、因子截面均值时常常更可靠。更一般的最优性讨论见第 12 章决策理论。


9.9 Delta 方法(MLE 版)

Delta 方法的一般形式与 Sharpe 比率标准误的推导已在本册第 01 章 1.5.7 节讲过,这里只给 MLE 语境下的表述。

定理 9.24。 设 \(\tau=g(\theta)\),\(g\) 可微且 \(g'(\theta)\ne0\),\(\hat\tau_n=g(\hat\theta_n)\),则

\[\frac{\hat\tau_n-\tau}{\widehat{\text{se}}(\hat\tau_n)}\rightsquigarrow N(0,1),\qquad\widehat{\text{se}}(\hat\tau_n)=|g'(\hat\theta_n)|\,\widehat{\text{se}}(\hat\theta_n).\]
(原书 (9.15) 分子多写了 \(\sqrt n\),因为 \(\widehat{\text{se}}\) 已含 \(1/\sqrt n\),这里按正确形式给出。)

例 9.25(对数几率)。 Bernoulli 中 \(\psi=\log\frac{p}{1-p}\),\(g'(p)=\frac1{p(1-p)}\),所以

\[\widehat{\text{se}}(\hat\psi)=\frac{1}{\hat p(1-\hat p)}\sqrt{\frac{\hat p(1-\hat p)}{n}}=\frac{1}{\sqrt{n\hat p(1-\hat p)}}.\]

例 9.26(对数标准差)。 \(N(\mu,\sigma^2)\),\(\mu\) 已知。\(\log f=-\log\sigma-\frac{(X-\mu)^2}{2\sigma^2}\),二阶导为 \(\frac{1}{\sigma^2}-\frac{3(X-\mu)^2}{\sigma^4}\),取期望并变号得 \(I(\sigma)=2/\sigma^2\),\(\widehat{\text{se}}(\hat\sigma)=\hat\sigma/\sqrt{2n}\)。对 \(\psi=\log\sigma\),\(g'=1/\sigma\),于是

\[\widehat{\text{se}}(\hat\psi)=\frac{1}{\sqrt{2n}},\]
与 \(\sigma\) 无关——对数变换是标准差的方差稳定化变换。量化读法:用 252 个日收益估计年化波动率,\(\log\hat\sigma\) 的标准误约 \(1/\sqrt{504}\approx4.5\%\),即波动率估计的相对误差约 ±9%(95%),无论波动率本身是 15% 还是 40%。(这是正态假设下的结果;厚尾时更大。)


9.10 多参数模型

设 \(\theta=(\theta_1,\dots,\theta_k)\),\(H_{jk}=\partial^2\ell_n/\partial\theta_j\partial\theta_k\)。Fisher 信息矩阵

\[I_n(\theta)=-\big[\mathbb E_\theta(H_{jk})\big]_{j,k=1}^k,\qquad J_n(\theta)=I_n^{-1}(\theta).\]

定理 9.27。 正则条件下 \(\hat\theta-\theta\approx N(0,J_n)\)。第 \(j\) 个分量 \(\frac{\hat\theta_j-\theta_j}{\widehat{\text{se}}_j}\rightsquigarrow N(0,1)\),\(\widehat{\text{se}}_j^2=J_n(j,j)\);且 \(\text{Cov}(\hat\theta_j,\hat\theta_k)\approx J_n(j,k)\)。

注意是先求逆再取对角元,而不是对角元取倒数。两者只在信息矩阵对角时相同;参数之间相关时,后者会低估标准误——这正是"冗余参数要付出代价"的数学表达。

推导拆解:用 2×2 的例子看清差别。设 \(I_n=\begin{pmatrix}a&b\\b&d\end{pmatrix}\),其逆的左上角是 \(\frac{d}{ad-b^2}=\frac{1}{a-b^2/d}\)。而"对角元取倒数"给的是 \(\frac1a\)。只要 \(b\ne0\),\(a-b^2/d<a\),所以真实方差 \(\frac{1}{a-b^2/d}\) 更大。 金融直觉:\(1/a\) 对应"假装另一个参数已知"时的方差;\(\frac{1}{a-b^2/d}\) 对应"另一个参数也要一起估"时的方差。CAPM 回归里估 alpha 时,市场收益均值不为零使 alpha 与 beta 的估计相关(\(b\ne0\)),beta 的不确定性就会"漏进"alpha 的标准误。这也是多重共线性让回归系数标准误膨胀的同一机制。

定理 9.28(多参数 Delta 方法)。 \(\tau=g(\theta)\),梯度 \(\nabla g\) 在 \(\hat\theta\) 处非零,则

\[\widehat{\text{se}}(\hat\tau)=\sqrt{(\hat\nabla g)^T\hat J_n(\hat\nabla g)}.\]

例 9.29(变异系数)。 \(X_i\sim N(\mu,\sigma^2)\),\(\tau=\sigma/\mu\)。可算得(原书习题 8)

\[I_n(\mu,\sigma)=\begin{pmatrix}n/\sigma^2&0\\0&2n/\sigma^2\end{pmatrix},\quad J_n=\frac1n\begin{pmatrix}\sigma^2&0\\0&\sigma^2/2\end{pmatrix},\quad\nabla g=\begin{pmatrix}-\sigma/\mu^2\\1/\mu\end{pmatrix},\]
\[\widehat{\text{se}}(\hat\tau)=\frac{1}{\sqrt n}\left(\frac{\hat\sigma^4}{\hat\mu^4}+\frac{\hat\sigma^2}{2\hat\mu^2}\right)^{1/2}.\]
(原书印刷第一项作 \(1/\hat\mu^4\),按推导应为 \(\hat\sigma^4/\hat\mu^4\)。)变异系数 \(\sigma/\mu\) 恰好是 Sharpe 比率的倒数;对 \(SR=\mu/\sigma\) 做同样计算,梯度为 \((1/\sigma,-\mu/\sigma^2)\),得 \(\text{se}(\widehat{SR})=\sqrt{(1+SR^2/2)/n}\),与第 01 章的结果一致。


9.11 参数 Bootstrap

Bootstrap(见本册第 08 章,原书第 8 章)在参数模型中也能用,唯一变化是抽样来源:

  • 非参数 Bootstrap:从经验分布 \(\hat F_n\) 有放回抽样;
  • 参数 Bootstrap:从拟合好的模型 \(f(x;\hat\theta_n)\) 抽样。

例 9.30。 续例 9.29:模拟 \(X_1^*,\dots,X_n^*\sim N(\hat\mu,\hat\sigma^2)\),计算 \(\hat\tau^*=\hat\sigma^*/\hat\mu^*\),重复 \(B\) 次,用 \(\hat\tau^*_b\) 的标准差作为 \(\widehat{\text{se}}_{\text{boot}}\)。

Bootstrap 比 Delta 方法容易得多——不用求导,不用求逆矩阵,任何复杂的参数函数(VaR、ES、期权价格、最优仓位)都照此办理。Delta 方法的优势是给出闭式标准误、计算快。参数 Bootstrap 的风险在于它继承了模型的全部假设:模型错了,Bootstrap 分布也跟着错。


9.12 检查模型假设

既然假设了参数模型,就应当检查它。非正式方法是画图:直方图明显双峰,正态假设就可疑;QQ 图尾部明显偏离直线,就说明厚尾。正式方法是拟合优度检验(第 10a 章 10.7 节)。但要记住:检验不拒绝不代表模型正确,可能只是功效不足。金融收益在大样本下几乎总会拒绝正态。


9.13 充分统计量(原书附录)

直觉。 数据中与参数有关的全部信息,能否压缩成少数几个数?能压缩到的那个"摘要"叫充分统计量。

定义 9.32。 若 \(f(x^n;\theta)=c\,f(y^n;\theta)\)(\(c\) 可以依赖 \(x^n,y^n\) 但不依赖 \(\theta\)),记 \(x^n\leftrightarrow y^n\)(两组数据似然形状相同)。若 \(T(x^n)=T(y^n)\) 蕴含 \(x^n\leftrightarrow y^n\),则称 \(T\) 充分(sufficient)。粗略地说:只知道 \(T\) 就能算出似然函数。

例 9.33–9.34。 Bernoulli 中 \(S=\sum X_i\) 充分;正态中 \((\bar X,S)\) 充分,因为

\[f(x^n;\mu,\sigma)=\Big(\frac{1}{\sigma\sqrt{2\pi}}\Big)^n\exp\Big\{-\frac{nS^2}{2\sigma^2}\Big\}\exp\Big\{-\frac{n(\bar X-\mu)^2}{2\sigma^2}\Big\}\]
只通过 \((\bar X,S)\) 依赖数据。充分统计量远非唯一:全部数据 \(T_1=(X_1,\dots,X_n)\) 充分;\(T_2=(\bar X,S)\) 充分;\(T_3=\bar X\) 不充分(算不出 \(\mathcal L(\mu,\sigma)\));\(T_4=(\bar X,S,X_3)\) 充分但冗余。

定义 9.35 与定理 9.36(最小充分)。 \(T\) 是最小充分统计量,若它充分且是任何其他充分统计量的函数。判别准则:"\(T(x^n)=T(y^n)\) 当且仅当 \(x^n\leftrightarrow y^n\)"。

例 9.37(划分视角)。 两次掷硬币,\(V=X_1\)、\(T=X_1+X_2\)、\(U=(T,X_1)\) 在结果集 \(\{(0,0),(0,1),(1,0),(1,1)\}\) 上诱导的划分分别是:\(V\):\(\{(0,0),(0,1)\},\{(1,0),(1,1)\}\);\(T\):\(\{(0,0)\},\{(0,1),(1,0)\},\{(1,1)\}\);\(U\):四个单点。\(V\) 不充分,\(T\) 和 \(U\) 充分,\(T\) 最小充分;\(U\) 不最小,因为 \((1,0)\leftrightarrow(0,1)\) 但 \(U\) 把它们分开了。

通常的定义。 若给定 \(T=t\) 时数据的条件分布不依赖 \(\theta\),则 \(T\) 充分(例 9.39:给定 \(T=1\),\((0,1)\) 与 \((1,0)\) 各 1/2,与 \(p\) 无关)。

白话解释:充分统计量就像一份合格的财务摘要。如果你只关心公司的盈利能力,一份正确的利润表就"充分"了,逐笔流水账不再提供额外信息;但只有"营业收入"一个数就不充分。对正态模型,\((\bar X,S)\) 是那份摘要:两组数据只要均值和标准差相同,不管具体数值怎么排布,它们给出的似然曲线形状完全一样,对 \((\mu,\sigma)\) 的推断也就完全一样。"最小充分"则是最精简、不能再压缩的摘要。

定理 9.40(因子分解定理,Factorization Theorem)。 \(T\) 充分当且仅当

\[f(x^n;\theta)=g\big(t(x^n),\theta\big)\,h(x^n).\]

定理 9.42(Rao–Blackwell)。 设 \(\hat\theta\) 为估计量,\(T\) 充分,\(\tilde\theta=\mathbb E(\hat\theta\mid T)\),则对所有 \(\theta\),\(R(\theta,\tilde\theta)\le R(\theta,\hat\theta)\)(\(R\) 为均方误差)。估计量应当只依赖充分统计量,否则对它取条件期望就能改进。 例 9.43:\(\hat\theta=X_1\) 无偏但浪费信息,\(\mathbb E(X_1\mid\sum X_i)=\bar X\)。

量化旁注:充分统计量是流式计算的理论基础。正态模型只需维护 \((n,\sum X_i,\sum X_i^2)\) 就能随时更新均值和方差;指数加权版本(EWMA 波动率)维护的也是加权的 \(\sum X_i^2\)。


9.14 指数族(原书附录)

单参数指数族:

\[f(x;\theta)=h(x)\exp\{\eta(\theta)T(x)-B(\theta)\}.\]
由因子分解定理 \(T(X)\) 充分,称为自然充分统计量。IID 样本的联合密度仍是指数族,\(\sum_iT(X_i)\) 充分。

  • 例 9.44(Poisson):\(\frac{\theta^xe^{-\theta}}{x!}=\frac1{x!}e^{x\log\theta-\theta}\),\(\eta=\log\theta\),\(T(x)=x\),\(B=\theta\)。
  • 例 9.45(Binomial):\(\binom nx\exp\{x\log\frac{\theta}{1-\theta}+n\log(1-\theta)\}\),\(\eta=\log\frac\theta{1-\theta}\),\(B(\theta)=-n\log(1-\theta)\)。(原书误印为 \(-n\log\theta\)。)
  • 例 9.46(Uniform\((0,\theta)\)):\(f(x^n;\theta)=\theta^{-n}I(x_{(n)}\le\theta)\),充分统计量是 \(\max X_i\),不能写成 \(\sum T(X_i)\)——支撑依赖参数的分布不是指数族。

自然参数形式:\(f(x;\eta)=h(x)e^{\eta T(x)-A(\eta)}\),\(A(\eta)=\log\int h(x)e^{\eta T(x)}dx\) 为对数配分函数。

定理 9.47。 \(\mathbb E(T(X))=A'(\eta)\),\(\mathbb V(T(X))=A''(\eta)\)。(以 Poisson 为例:\(A(\eta)=e^\eta\),\(A'=A''=e^\eta=\theta\),均值与方差都是 \(\theta\)。)

多参数指数族:\(f(x;\theta)=h(x)\exp\{\sum_j\eta_j(\theta)T_j(x)-B(\theta)\}\)。例 9.48:正态分布 \(\eta_1=\mu/\sigma^2\)、\(T_1=x\),\(\eta_2=-1/(2\sigma^2)\)、\(T_2=x^2\),故 \((\sum X_i,\sum X_i^2)\) 充分。

指数族的意义:对数似然关于自然参数是凹函数,MLE 唯一、数值上好求;它也是广义线性模型(logistic 回归、Poisson 回归)的基础。


9.15 数值计算 MLE:Newton–Raphson 与 EM

绝大多数实用模型的 MLE 没有闭式解。两类迭代法都生成 \(\theta^0,\theta^1,\dots\),在理想条件下收敛到 MLE;好的初值非常重要,矩估计常是好初值。

9.15.1 Newton–Raphson

在当前点 \(\theta^j\) 展开得分:\(0=\ell'(\hat\theta)\approx\ell'(\theta^j)+(\hat\theta-\theta^j)\ell''(\theta^j)\),得迭代

\[\theta^{j+1}=\theta^j-\frac{\ell'(\theta^j)}{\ell''(\theta^j)},\qquad\text{多参数:}\ \theta^{j+1}=\theta^j-H^{-1}\ell'(\theta^j).\]

金融直觉:Newton–Raphson 就是你在计算债券 YTM 时用过的迭代:猜一个收益率,用价格对收益率的导数(久期)修正猜测,重复直到价格吻合。这里要找的是"得分 \(\ell'\) 等于 0"的根,修正量用的是得分的导数 \(\ell''\)。一元数值例:\(\ell'(p)=220/p-180/(1-p)\),从 \(p^0=0.5\) 出发,\(\ell'(0.5)=80\),\(\ell''(0.5)=-220/0.25-180/0.25=-1600\),一步得 \(p^1=0.5+80/1600=0.55\),已到 MLE。

它与渐近正态性证明里的展开是同一个式子。一个副产品:收敛时的 Hessian \(H\) 就是观测信息,\(-H^{-1}\) 直接给出协方差矩阵。实践中常用拟牛顿法(BFGS)用梯度信息逐步逼近 Hessian,避免显式计算二阶导;详见第 04 册数值最优化。

9.15.2 EM 算法

适用场景。 观测数据 \(Y\) 的似然难以最大化,但若能"看到"某个隐变量 \(Z\),完整数据 \((Y,Z)\) 的似然就很好处理。即感兴趣模型是一个简单模型的边缘:\(f(y;\theta)=\int f(y,z;\theta)dz\)。\(Z\) 称为隐藏(hidden)、潜在(latent)或缺失(missing)数据。

例 9.49(正态混合)。 观测以概率 \(p\) 来自第二个正态、以概率 \(1-p\) 来自第一个,但不知道每个观测来自哪个:

\[f(y;\theta)=(1-p)\phi(y;\mu_0,\sigma_0)+p\,\phi(y;\mu_1,\sigma_1),\qquad\theta=(\mu_0,\sigma_0,\mu_1,\sigma_1,p).\]
直接最大化 \(\prod_i[\cdots]\) 的对数很难(对数里有和)。若知道标签 \(Z_i\in\{0,1\}\),就只是两组分别求均值方差。

EM 步骤。 选初值 \(\theta^0\),对 \(j=0,1,2,\dots\) 重复:

  1. E 步:计算 \(J(\theta\mid\theta^j)=\mathbb E_{\theta^j}\Big[\log\frac{f(Y^n,Z^n;\theta)}{f(Y^n,Z^n;\theta^j)}\,\Big|\,Y^n=y^n\Big]\),期望对缺失数据 \(Z^n\) 在当前参数下的条件分布求;
  2. M 步:令 \(\theta^{j+1}=\arg\max_\theta J(\theta\mid\theta^j)\)。

单调性定理:似然永不下降。 由 \(f(y^n,z^n;\theta)=f(z^n\mid y^n;\theta)f(y^n;\theta)\),

\[\log\frac{\mathcal L(\theta^{j+1})}{\mathcal L(\theta^j)}=J(\theta^{j+1}\mid\theta^j)+K(f_j,f_{j+1}),\]
其中 \(f_j=f(z^n\mid y^n;\theta^j)\),\(K\) 为 KL 距离。M 步保证 \(J(\theta^{j+1}\mid\theta^j)\ge J(\theta^j\mid\theta^j)=0\),又 \(K\ge0\),故 \(\mathcal L(\theta^{j+1})\ge\mathcal L(\theta^j)\)。注意单调不降只保证收敛到局部极大或鞍点,多个初值试算是标准做法。

例 9.50(简化版:\(p=1/2\)、\(\sigma_0=\sigma_1=1\))。 完整数据对数似然(略常数)

\[\tilde\ell=-\frac12\sum_i(1-z_i)(y_i-\mu_0)^2-\frac12\sum_iz_i(y_i-\mu_1)^2.\]
(原书此处漏了平方号。)它关于 \(z_i\) 是线性的,所以 E 步只需把 \(z_i\) 换成条件期望,由贝叶斯定理
\[\tau_i\equiv\mathbb E(Z_i\mid y^n,\theta^j)=\frac{\phi(y_i;\mu_1^j,1)}{\phi(y_i;\mu_1^j,1)+\phi(y_i;\mu_0^j,1)}.\]
M 步对 \(\mu_0,\mu_1\) 求导令为 0:
\[\mu_1^{j+1}=\frac{\sum_i\tau_iy_i}{\sum_i\tau_i},\qquad\mu_0^{j+1}=\frac{\sum_i(1-\tau_i)y_i}{\sum_i(1-\tau_i)}.\]

推导拆解:\(\tau_i\) 的公式就是贝叶斯定理,和第 01 章医学检测的例子结构相同。先验:两种状态各 1/2;"似然":在状态 1 下看到 \(y_i\) 的密度 \(\phi(y_i;\mu_1,1)\),状态 0 下是 \(\phi(y_i;\mu_0,1)\)。后验 \(\mathbb P(Z_i=1\mid y_i)=\frac{\frac12\phi_1}{\frac12\phi_1+\frac12\phi_0}\),1/2 约掉即得。 M 步:对 \(\mu_1\) 求导,\(\frac{\partial\tilde\ell}{\partial\mu_1}=\sum_i\tau_i(y_i-\mu_1)=0\)(E 步已把 \(z_i\) 换成 \(\tau_i\)),移项得加权平均。 金融直觉:可以把 \(\tau_i\) 想成"这一天有多大比例算作危机日"。一个 −5% 的日子几乎全部计入危机组(\(\tau\approx1\)),一个 +0.1% 的日子几乎全部计入平静组;然后各组按这些"部分计入"的权重重算均值和波动。重复下去,分组和参数互相校正,直到稳定。

一句话概括:E 步算每个观测属于各成分的"责任"\(\tau_i\),M 步做加权平均。一般情形(5 个参数全未知)同理,M 步还要更新加权方差 \(\sigma_1^2=\sum\tau_i(y_i-\mu_1)^2/\sum\tau_i\) 与 \(p=\bar\tau\)。


9.16 量化实战

9.16.1 用 t 分布 MLE 估计 VaR,并给出它的置信区间

日收益厚尾,用正态 MLE 估 VaR 会低估尾部风险。下面用 Student-t 分布做 MLE:矩估计(由峰度反推自由度)给初值,BFGS 最大化似然,用观测信息求标准误,再用多参数 Delta 方法和参数 Bootstrap 两种方式求 1% VaR 的标准误。参数做了变换 \(\theta=(\mu,\log s,\log(\nu-2))\) 以保证 \(s>0\)、\(\nu>2\)——这利用了同变性:在变换后的参数上求 MLE,再变回去即可。

import numpy as np
from scipy import stats, optimize
from statsmodels.tools.numdiff import approx_hess3, approx_fprime

rng = np.random.default_rng(42)
# 模拟 4 年日收益:t(4) 厚尾,位置 0.04%,尺度 1%
n, mu0, s0, nu0 = 1000, 0.0004, 0.01, 4.0
x = mu0 + s0 * rng.standard_t(nu0, size=n)

# 参数化 theta = (mu, log s, log(nu-2)),保证 s>0、nu>2
def unpack(th):
    return th[0], np.exp(th[1]), 2 + np.exp(th[2])

def negll(th, x=x):
    mu, s, nu = unpack(th)
    return -np.sum(stats.t.logpdf(x, df=nu, loc=mu, scale=s))

# 矩估计作初值:均值、标准差、由峰度反推 nu(超额峰度 = 6/(nu-4))
m, sd = x.mean(), x.std()
ek = max(stats.kurtosis(x), 0.5)
nu_mm = 4 + 6 / ek
s_mm = sd * np.sqrt((nu_mm - 2) / nu_mm)
th0 = np.array([m, np.log(s_mm), np.log(nu_mm - 2)])
res = optimize.minimize(negll, th0, method="BFGS")
th = res.x
mu_h, s_h, nu_h = unpack(th)

# 观测 Fisher 信息 = 负对数似然的 Hessian;其逆是渐近协方差(变换后参数)
H = approx_hess3(th, negll)
J = np.linalg.inv(H)
# 用 Delta 方法把协方差转换回 (mu, s, nu)
G = np.diag([1.0, s_h, nu_h - 2])           # d(mu,s,nu)/d(theta)
J_nat = G @ J @ G.T
se = np.sqrt(np.diag(J_nat))
print(f"矩估计初值: mu={m:.5f}, s={s_mm:.5f}, nu={nu_mm:.2f}")
print(f"t-MLE     : mu={mu_h:.5f}({se[0]:.5f}), s={s_h:.5f}({se[1]:.5f}), nu={nu_h:.2f}({se[2]:.2f})")

# 感兴趣参数:1% VaR = -(mu + s * t_nu^{-1}(0.01))
def var1(th):
    mu, s, nu = unpack(th)
    return -(mu + s * stats.t.ppf(0.01, nu))
v_h = var1(th)
grad = approx_fprime(th, var1)
se_delta = np.sqrt(grad @ J @ grad)
print(f"1% VaR (t-MLE) = {v_h:.4%}, Delta 方法 se = {se_delta:.4%}, "
      f"95% CI = [{v_h-1.96*se_delta:.4%}, {v_h+1.96*se_delta:.4%}]")

# 参数 bootstrap:从拟合的 t 分布重新抽样、重新估计
B, vb = 300, []
for b in range(B):
    xb = mu_h + s_h * rng.standard_t(nu_h, size=n)
    rb = optimize.minimize(negll, th, args=(xb,), method="BFGS")
    vb.append(var1(rb.x))
vb = np.array(vb)
print(f"参数 bootstrap se = {vb.std(ddof=1):.4%}, 百分位区间 = "
      f"[{np.quantile(vb,0.025):.4%}, {np.quantile(vb,0.975):.4%}]")

# 对比:错误地假设正态
v_norm = -(m + sd * stats.norm.ppf(0.01))
v_true = -(mu0 + s0 * stats.t.ppf(0.01, nu0))
print(f"正态 MLE 的 1% VaR = {v_norm:.4%};真实 1% VaR = {v_true:.4%};样本 1% 分位 = {-np.quantile(x,0.01):.4%}")

输出:

矩估计初值: mu=-0.00016, s=0.01157, nu=5.37
t-MLE     : mu=0.00003(0.00037), s=0.00968(0.00041), nu=3.25(0.39)
1% VaR (t-MLE) = 4.1315%, Delta 方法 se = 0.2669%, 95% CI = [3.6084%, 4.6545%]
参数 bootstrap se = 0.2720%, 百分位区间 = [3.6588%, 4.6910%]
正态 MLE 的 1% VaR = 3.4129%;真实 1% VaR = 3.7069%;样本 1% 分位 = 3.9864%

读法:(1) 矩估计给出的 \(\nu=5.37\) 偏离较大(样本峰度本身方差极大,厚尾时尤甚),但作为初值足够;MLE 更有效。(2) Delta 方法与参数 Bootstrap 的标准误几乎相同(0.267% vs 0.272%),说明样本量 1000 时渐近正态近似已经不错。(3) 正态 MLE 的 VaR 3.41% 明显低于真实值 3.71%,而 t-MLE 的区间 [3.61%, 4.65%] 覆盖了真值——模型设定错了,再"精确"的 MLE 也是对错误问题的精确回答。(4) 1% VaR 的 95% 区间宽度约 ±0.5 个百分点,相对误差约 ±13%:风险报告里的 VaR 数字远没有它的小数位看起来那么精确。

9.16.2 用 EM 估计两状态收益混合模型

市场常被描述为在"平静"与"危机"两种状态间切换。把每天的状态当作隐变量,就是例 9.49 的正态混合。下面手写 EM(5 个参数全部未知),在每步检查似然单调不降,并与 sklearn 的结果对照。

import numpy as np
from scipy import stats
from sklearn.mixture import GaussianMixture

rng = np.random.default_rng(7)
# 两状态收益:平静(85%)与危机(15%),日收益,单位 %
n = 2000
z = rng.random(n) < 0.15
y = np.where(z, rng.normal(-0.20, 2.5, n), rng.normal(0.05, 0.8, n))

def loglik(y, p, m0, s0, m1, s1):
    return np.sum(np.log((1-p)*stats.norm.pdf(y, m0, s0) + p*stats.norm.pdf(y, m1, s1)))

# 初值:两个成分均值相同、方差一小一大
p, m0, s0, m1, s1 = 0.5, y.mean(), 0.5*y.std(), y.mean(), 2*y.std()
ll_old = -np.inf
for it in range(500):
    # E 步:每个观测属于"危机"成分的后验概率(责任)tau_i
    a1 = p * stats.norm.pdf(y, m1, s1)
    a0 = (1-p) * stats.norm.pdf(y, m0, s0)
    tau = a1 / (a0 + a1)
    # M 步:加权的 MLE
    p = tau.mean()
    m1 = np.sum(tau*y)/tau.sum();      s1 = np.sqrt(np.sum(tau*(y-m1)**2)/tau.sum())
    m0 = np.sum((1-tau)*y)/(1-tau).sum(); s0 = np.sqrt(np.sum((1-tau)*(y-m0)**2)/(1-tau).sum())
    ll = loglik(y, p, m0, s0, m1, s1)
    assert ll >= ll_old - 1e-9        # EM 单调性:似然不降
    if ll - ll_old < 1e-8: break
    ll_old = ll
print(f"EM 迭代 {it} 次, 对数似然 {ll:.2f}")
print(f"平静: mu={m0:.3f}, sigma={s0:.3f};  危机: mu={m1:.3f}, sigma={s1:.3f}, p={p:.3f}")

gm = GaussianMixture(2, n_init=5, tol=1e-8, max_iter=2000, random_state=0).fit(y.reshape(-1,1))
k = np.argmax(gm.covariances_.ravel())
print(f"sklearn 对照(对数似然 {gm.score(y.reshape(-1,1))*n:.2f}): 危机 mu={gm.means_[k,0]:.3f}, sigma={np.sqrt(gm.covariances_.ravel()[k]):.3f}, p={gm.weights_[k]:.3f}")

# 责任 tau 作为"危机概率"信号:与真实状态的吻合度
pred = tau > 0.5
print(f"tau>0.5 判为危机的天数 {pred.sum()},真实危机天数 {z.sum()},判对率 {np.mean(pred==z):.3f}")

输出:

EM 迭代 99 次, 对数似然 -3044.73
平静: mu=0.052, sigma=0.828;  危机: mu=-0.186, sigma=2.581, p=0.143
sklearn 对照(对数似然 -3044.73): 危机 mu=-0.186, sigma=2.580, p=0.143
tau>0.5 判为危机的天数 138,真实危机天数 294,判对率 0.914

读法:EM 准确恢复了两个状态的参数(真值:平静 0.05/0.8,危机 −0.20/2.5,\(p=0.15\))。值得注意的是最后一行:参数估计得很准,但单日状态分类并不准——真实危机日有 294 天,只有 138 天被判为危机,因为危机状态的大部分日收益落在 ±1% 内,与平静日无法区分。混合模型适合刻画"收益分布是什么样",用来逐日判断"今天是不是危机"则要加上状态的持续性(隐马尔可夫模型,HMM,其 Baum–Welch 算法就是 EM 的推广;见第 06 册马尔可夫转换模型相关内容)。

另一个工程细节:编写本例时,第一次用 sklearn 默认设置(单个初值、较松的容差)时得到了不同的局部解,加上 n_init=5 与更严的容差后才与手写 EM 一致。EM 只保证似然不降,不保证全局最优,多初值是必需的。

9.16.3 其他应用速记

  • GARCH 族的估计就是本章方法的直接应用:写出条件正态(或条件 t)对数似然,BFGS 最大化,Hessian 求标准误;arch 包默认报告 QMLE 稳健标准误(见 9.7 节旁注)。详见第 06 册。
  • Fisher 信息与样本量规划:\(\text{se}\approx1/\sqrt{nI(\theta)}\) 告诉你要把某个参数估计到给定精度需要多少数据。例如估计日收益均值需要 \(I=1/\sigma^2\),于是 \(\text{se}(\hat\mu)=\sigma/\sqrt n\)——年化 Sharpe 0.5 的策略需要约 16 年数据才能让 \(t\approx2\)。
  • logit 变换与比例的区间:胜率、违约率、成交率接近 0 或 1 时,直接对 \(\hat p\) 用 Wald 区间会越界,先对 \(\log\frac{p}{1-p}\) 构造区间(例 9.25)再变换回去更可靠。
  • 模型错设:原书习题 6 的结论——若数据不是正态,依赖正态模型的 \(\Phi(\bar X)\) 估计 \(\mathbb P(X>0)\) 不相合,而非参数的 \(\bar Y=\frac1n\sum I(X_i>0)\) 仍相合——对估计"上涨天数占比"这类量同样成立。

本章小结

参数推断的主线是:写出似然 → 求 MLE → 用 Fisher 信息求标准误 → 用 Delta 方法或参数 Bootstrap 处理参数函数。MLE 之所以是默认方法,是因为它在正则条件下相合(它在寻找与真分布 KL 距离最小的模型)、同变(参数函数的 MLE 直接代入)、渐近正态且渐近方差达到下界 \(1/I_n(\theta)\)。这些好性质有两个前提:模型正确、参数维数相对样本量较小。模型错设时 MLE 收敛到"伪真值",标准误要用三明治形式;高维时 MLE 可能远非最优(第 12 章)。充分统计量和指数族解释了"哪些数据摘要就够了";Newton–Raphson 与 EM 解决"算不出闭式解"的问题,其中 EM 是处理隐状态(市场状态、缺失数据)的标准工具。

概念 公式 / 结论
矩估计 解 \(\alpha_j(\hat\theta)=\frac1n\sum X_i^j\),\(j=1..k\)
似然 / 对数似然 \(\mathcal L_n(\theta)=\prod f(X_i;\theta)\),\(\ell_n=\log\mathcal L_n\)
KL 距离 \(D(f,g)=\int f\log(f/g)\ge0\)
得分 \(s=\partial_\theta\log f(X;\theta)\),\(\mathbb Es=0\)
Fisher 信息 \(I(\theta)=\mathbb E s^2=-\mathbb E\,\partial^2_\theta\log f\),\(I_n=nI\)
渐近正态 \(\hat\theta\approx N(\theta,1/I_n(\theta))\),\(\widehat{\text{se}}=1/\sqrt{I_n(\hat\theta)}\)
Wald 区间 \(\hat\theta\pm z_{\alpha/2}\widehat{\text{se}}\)
ARE \(t^2/u^2\);正态下中位数对均值 \(2/\pi\)
Delta 方法 \(\widehat{\text{se}}(g(\hat\theta))=\lvert g'(\hat\theta)\rvert\widehat{\text{se}}(\hat\theta)\);多元 \(\sqrt{\nabla g^T J_n\nabla g}\)
多参数 \(J_n=I_n^{-1}\),\(\widehat{\text{se}}_j=\sqrt{J_n(j,j)}\)
常见 Fisher 信息 Bernoulli \(\frac1{p(1-p)}\);Poisson \(\frac1\lambda\);正态 \(\mu\):\(\frac1{\sigma^2}\),\(\sigma\):\(\frac2{\sigma^2}\)
因子分解 \(T\) 充分 ⇔ \(f(x^n;\theta)=g(t(x^n),\theta)h(x^n)\)
Rao–Blackwell \(\mathbb E(\hat\theta\mid T)\) 的 MSE 不大于 \(\hat\theta\)
指数族 \(f=h(x)e^{\eta T(x)-A(\eta)}\),\(\mathbb ET=A'\),\(\mathbb VT=A''\)
Newton–Raphson \(\theta^{j+1}=\theta^j-H^{-1}\ell'(\theta^j)\)
EM E 步求 \(\mathbb E[\log f(Y,Z;\theta)\mid y,\theta^j]\),M 步最大化;似然单调不降
QMLE 稳健协方差 \(A^{-1}BA^{-1}/n\),\(A=-\mathbb E\ell_i''\),\(B=\mathbb V\ell_i'\)

练习

基础

  1. 求 Gamma\((\alpha,\beta)\) 的矩估计。(提示:\(\mathbb EX=\alpha\beta\),\(\mathbb VX=\alpha\beta^2\),得 \(\hat\beta=\hat\sigma^2/\bar X\),\(\hat\alpha=\bar X^2/\hat\sigma^2\)。)
  2. Poisson\((\lambda)\) 样本:求矩估计、MLE 与 Fisher 信息,并写出 \(\lambda\) 的 95% Wald 区间。(答案:二者都是 \(\bar X\);\(I(\lambda)=1/\lambda\);\(\bar X\pm1.96\sqrt{\bar X/n}\)。)
  3. 指数分布 \(f(x;\beta)=\frac1\beta e^{-x/\beta}\) 描述订单间隔时间。求 \(\beta\) 的 MLE、Fisher 信息,以及"间隔超过 1 秒的概率"\(e^{-1/\beta}\) 的 MLE 与 Delta 方法标准误。(答案:\(\hat\beta=\bar X\),\(I(\beta)=1/\beta^2\),\(\text{se}(\hat\beta)=\hat\beta/\sqrt n\);\(g'(\beta)=e^{-1/\beta}/\beta^2\),\(\text{se}=e^{-1/\hat\beta}/(\hat\beta\sqrt n)\)。)
  4. Uniform\((0,\theta)\) 中证明 \(\hat\theta=X_{(n)}\) 相合。(提示:\(\mathbb P(X_{(n)}<\theta-\epsilon)=(1-\epsilon/\theta)^n\to0\)。)
  5. 某信号 400 笔交易中 220 笔盈利。求胜率的 95% Wald 区间与对数几率 \(\psi\) 的 95% 区间,再把后者变换回 \(p\) 的尺度,与前者比较。(提示:\(\hat\psi=\log(220/180)\approx0.20\),\(\text{se}=1/\sqrt{400\times0.55\times0.45}\approx0.10\)。)

进阶

  1. 设收益 \(X_i\sim N(\mu,\sigma^2)\),求 5% 分位数 \(\tau=\mu-1.645\sigma\)(即 VaR 的相反数)的 MLE 与标准误。(答案:\(\text{se}^2=\frac{\sigma^2}{n}(1+1.645^2/2)\)。这是原书习题 3 的量化版本。)
  2. 两种执行算法,各 200 笔订单,在限定时间内全部成交的分别为 160 与 148 笔。用多参数 Delta 方法求 \(p_1-p_2\) 的 90% 区间,再用参数 Bootstrap 计算并比较。(原书习题 7。)
  3. 推导 Student-t 位置-尺度族的得分函数,说明为什么它对极端观测的"权重"有上界,而正态模型没有。由此解释 t-MLE 比样本均值更稳健的原因。(提示:\(\partial_\mu\log f\propto\frac{(\nu+1)(x-\mu)}{\nu s^2+(x-\mu)^2}\),\(|x-\mu|\to\infty\) 时趋于 0。)
  4. 修改 9.16.2 节代码:把危机比例改为 5%、危机均值改为 −0.5,样本量改为 500,多次换随机种子运行,观察 EM 估计的稳定性与局部最优问题。
  5. 证明 EM 的单调性公式 \(\log\frac{\mathcal L(\theta^{j+1})}{\mathcal L(\theta^j)}=J(\theta^{j+1}\mid\theta^j)+K(f_j,f_{j+1})\)。

原书推荐习题:第 9 章 3(分位数 MLE 与 Delta 方法,必做)、6(模型错设与 ARE)、7(两比例之差)、9、10(Delta、参数与非参数 Bootstrap 的比较及失效情形)、2;另建议完成例 9.50 的 EM 编程。

原书对照

本章内容 原书章节 PDF 页码
参数模型、感兴趣参数 第 9 章引言、9.1 p.131–132
矩估计 9.2 p.132–134
最大似然、例 9.10–9.12 9.3 p.134–136
MLE 性质总览 9.4 p.136–138
相合性、KL 距离 9.5 p.138–139
同变性 9.6 p.139–140
得分、Fisher 信息、渐近正态 9.7 p.140–142
最优性、ARE 9.8 p.142–143
Delta 方法 9.9 p.143–145
多参数模型、多元 Delta 9.10 p.145–146
参数 Bootstrap、检查假设 9.11–9.12 p.146–147
附录:证明、充分性与 Rao–Blackwell、指数族、Newton–Raphson 与 EM 9.13 p.147–158
习题 9.14 p.158–160

(PDF 页码 = 原书正文页码 + 17。)