【量化交易】风险模型：Barra 多因子、风险归因、压力测试

把组合的「风险」当作一个数字（年化波动、最大回撤、夏普）报给老板，是基金经理在路演 PPT 上做的事；把组合的「风险」当作一个工程系统拆给自己看，才是风控工程师每天上班要做的事。这两件事在多数团队里被混为一谈，结果就是基金经理在牛市里看夏普 3.0 自我感动，回撤来的那一周才发现，组合里 60% 的风险都集中在同一个风格因子上，行业暴露与市值暴露和指数没有任何区别，所谓的「α」不过是一个被加了杠杆的小盘 + 低估值组合。

风险模型这件事，从 1970 年代 Barr Rosenberg 在 Barra（旧称 BARRA Inc.）做出第一版多因子结构开始，到今天 MSCI Barra USE4、Axioma、Northfield、Bloomberg PORT、Wolfe Research、各家券商自研模型，几十年下来核心结构没有大变：把每只股票的收益拆成「国家因子 + 行业因子 + 风格因子 + 特异收益」四块，对每只股票做横截面回归（cross-sectional regression）估计因子收益，用因子协方差 + 特异波动合成股票协方差矩阵，再喂给优化器。这套结构在 A 股、港股、美股、日股都能跑得起来，区别只在于行业分类口径、风格因子定义、估计方法的细节。

本文不复述「Barra 因子是什么」这种百科式描述，而是把风险模型当作一个流水线工程系统讲清楚：从协方差矩阵估计的数值病根开始，到 Barra 多因子的横截面回归骨架，再到组合风险归因、VaR/ES 的三种估计法、压力测试情景的构造方法，最后落到工程实现上的每日刷新与回测嵌入。每一节都给出可以直接抄走的 Python 代码片段，关键结论附带数学推导与文献来源。

风险提示与适用范围：本文不构成任何投资建议。所有代码、数据与示例数值均用于说明算法，未经过完整生产审计，直接套用产生的资金损失由使用者自负。文中风险模型骨架与 MSCI Barra 商业产品的实际公式有差异，仅用于教学与原理复现；商用应使用厂商授权的官方因子库或自建经过严格回测的内部模型。VaR 与压力测试的监管报送（如银保监会《商业银行资本管理办法》、巴塞尔 FRTB、SEC Rule 18f-4、欧盟 CRR/CRD）有专门的口径要求，本文不替代任何监管文档。

一、风险模型的目的

写代码之前先要回答一个问题：风险模型这套东西到底为什么存在，做出来给谁用。多数人会脱口而出「估计组合波动」，这个回答是错的——估计组合波动只是风险模型的副产品，不是它的目的。如果只为了估计一个数字，市场上有 GARCH、EWMA、Realized Volatility 一堆方法，根本不需要 Barra 这套行业 + 风格 + 国家的三层结构。

风险模型真正存在的理由有三条，按工程上重要性从高到低排：

事前预测（ex ante forecasting）

第一条用途是事前给出组合的风险预算（risk budget）：在交易日开盘之前告诉基金经理，如果你按当前持仓一直拿到下个调仓日，未来 N 天可能损失多大、损失主要来自哪些因子、有没有超过任何一条硬约束。这件事如果靠「直接对组合历史净值算波动率」做不到——刚刚换过仓的组合，历史净值反映的是旧组合的风险，而不是新组合的风险。

事前预测必须建立在当前持仓 \(w_t\) 上，对每个成分股的未来协方差矩阵 \(\Sigma_{t+1}\) 做估计：

\[ \sigma_p^2 = w_t^\top \Sigma_{t+1} w_t \]

样本协方差 \(\hat\Sigma\) 在 \(N \gg T\)（股票数远多于历史长度）时退化为奇异矩阵，无法直接用。Barra 用因子分解把 \(\Sigma\) 的自由度从 \(\frac{N(N+1)}{2}\) 压到 \(K(K+1)/2 + N\)（\(K\) 为因子数，A 股大约 40 个），让协方差矩阵在 \(N=4000\) 这种规模下仍然可用。

事后归因（ex post attribution）

第二条用途是事后回答一个被无数 LP 在月度报告上质问过的问题：你这个月赚的钱（亏的钱）到底是怎么来的？是市场涨了你 beta 高跟着涨，还是你押对了某个行业，还是你某个风格因子（比如低市值）在风口上？还是说你确实选股选出了 alpha？

把组合实际收益 \(r_p\) 拆成因子贡献和特异贡献：

\[ r_p = \sum_{k=1}^{K} \beta_{p,k} f_k + \alpha_p \]

其中 \(\beta_{p,k} = w^\top X_{\cdot, k}\) 是组合在第 \(k\) 个因子上的暴露（就是把成分股的因子暴露按权重加权求和），\(f_k\) 是当期估计出的因子收益，\(\alpha_p\) 是无法被因子解释的残差，即所谓「真 alpha」。LP 看到「这个月组合赚了 2.3%，其中 1.8% 来自市值因子（小票上涨），0.4% 来自电子行业，0.1% 才是 alpha」之后，对基金经理的判断会立刻发生变化。

组合优化输入（optimizer input）

第三条用途，也是工程上最关键的一条：风险模型是组合优化器的输入。Markowitz 的均值-方差优化、风险平价、Black-Litterman、最大分散化（maximum diversification），所有这些方法的核心都是「把协方差矩阵塞进二次规划求解」。如果协方差矩阵不可信，优化器解出来的权重就是噪声放大器——理论上最优、实际上每天换手 80%、回测上夏普 5.0、实盘上扣完手续费亏钱。

风险模型给优化器提供的不是一个矩阵，而是一组结构化约束工具：

• 因子暴露上下界（如「市值因子暴露不超过基准 ±0.3 个标准差」）；
• 行业暴露上下界（如「电子行业相对基准超配不超过 5%」）；
• 跟踪误差（tracking error）目标（如「年化跟踪误差控制在 4%-6%」）；
• 特异风险占比下界（如「至少 40% 的预测方差来自特异收益」，防止组合变成纯因子押注）。

这些约束直接写在二次规划里，必须靠风险模型把每只股票的因子暴露 \(X\) 与特异波动 \(\delta\) 算出来才能用。

总结一句话：风险模型不是给基金经理报数字看的，是给组合优化器、风控系统、归因报告同时供数的统一基础设施。这条结论决定了后面所有工程决策——估计频率（每日 vs 每周）、刷新延迟（T+0 vs T+1）、覆盖范围（沪深 + 港股 + 中概 vs 仅 A 股）、版本切换（如何避免今天用新模型回算昨天报告）都按这条原则取舍。

下面这张图是后面整套讨论的几何骨架，展示 Barra 因子的四层分解结构。

二、协方差矩阵估计

风险模型的第一道关，是协方差矩阵 \(\Sigma\) 怎么算。这件事看起来是一行 np.cov(returns) 就能解决的小问题，实际上是几十年学术文献和工程实践堆起来的大坑。下面把主要病根一条一条讲清楚，再讲为什么 Barra 那套因子分解是目前为止最好的工程妥协。

样本协方差为什么不稳定

设有 \(N\) 只股票、\(T\) 天历史收益矩阵 \(R \in \mathbb{R}^{T \times N}\)，样本协方差为：

\[ \hat\Sigma = \frac{1}{T-1}(R - \bar R)^\top (R - \bar R) \]

它有几个工程上致命的问题：

第一，自由度不够。\(\hat\Sigma\) 有 \(\frac{N(N+1)}{2}\) 个独立参数，A 股 4000 只股票就是约 800 万个参数；而样本量是 \(N \times T\)，假设取 \(T=252\)（一年），样本量 100 万，参数比样本多 8 倍，矩阵秩为 \(\min(N, T-1) = 251\)，剩下 \(N - 251 \approx 3749\) 个特征值是数值零，矩阵奇异，求逆直接爆炸。

第二，特征值分布严重失真。即便 \(T \gg N\)，根据 Marčenko-Pastur 定理，当 \(N/T\) 不可忽略时，样本协方差的特征值会偏离真实分布，最大特征值被高估、最小特征值被低估到接近零。把这个矩阵直接喂给优化器（即便能求逆），优化器会过度配置到「最低估特征值方向」上，因为这个方向看起来是「最低风险方向」，但其实只是数值噪声。这个现象叫作 “errors maximization”（Michaud 1989），是均值-方差优化在实盘上失效的核心原因之一。

第三，时变性。股票之间的真实协方差结构本身就是时变的——牛市相关性低、危机时相关性接近 1。任何把过去 252 天平均的方法，本质上是用低相关性环境的历史去预测高相关性环境，结果在最需要风险模型的时候最不准。

shrinkage 估计：把样本协方差和先验混合

Ledoit & Wolf（2003、2004、2017）提出的 shrinkage 估计法，思路是把样本协方差 \(\hat\Sigma\) 与一个结构化先验 \(F\) 加权平均：

\[ \hat\Sigma_{\text{shrink}} = \delta F + (1 - \delta) \hat\Sigma, \quad \delta \in [0, 1] \]

先验 \(F\) 通常取三种之一：

1. 常数相关阵：\(F_{ii} = \hat\sigma_i^2\)，\(F_{ij} = \bar\rho \cdot \hat\sigma_i \hat\sigma_j\)，其中 \(\bar\rho\) 是所有股票对的平均相关系数。这是最简单也最常用的先验，A 股、美股都被验证有效。
2. 单因子模型：\(F = \sigma_m^2 \beta \beta^\top + D\)，其中 \(\beta\) 是各股票对市场指数的回归 beta，\(D\) 是对角残差方差矩阵。这是 Sharpe 单指数模型。
3. 常数方差对角阵：\(F = \bar\sigma^2 I\)。这是最强的收缩，几乎丢光所有相关结构信息，只在样本极少时才用。

收缩系数 \(\delta\) 由 Ledoit-Wolf 的最优公式估计，目标是最小化 \(\|\hat\Sigma_{\text{shrink}} - \Sigma\|_F^2\) 的期望（Frobenius 范数）。这套方法在 sklearn 里直接有 sklearn.covariance.LedoitWolf 实现，工程上几行就能跑：

from sklearn.covariance import LedoitWolf

lw = LedoitWolf().fit(returns)
sigma_shrunk = lw.covariance_
print(f"shrinkage intensity = {lw.shrinkage_:.4f}")

shrinkage 的优点是无须先验知识、数值稳定、对优化器友好；缺点是它把整个协方差结构当成「黑箱」收缩，不能告诉你「为什么两只股票相关性高」——也就是没有归因能力。这是为什么大型机构必须在 shrinkage 之外另上 Barra 这种因子模型：shrinkage 服务于优化器的数值稳定性，因子模型服务于归因与约束。

因子模型为什么是工程最优解

把每只股票的收益拆成因子部分和特异部分：

\[ r_i = \sum_{k=1}^{K} x_{ik} f_k + u_i \]

如果假设特异收益 \(u_i\) 在股票之间互相独立、且独立于因子，那么协方差矩阵可以写成：

\[ \Sigma = X \Sigma_f X^\top + \Delta \]

其中 \(X \in \mathbb{R}^{N \times K}\) 是暴露矩阵，\(\Sigma_f \in \mathbb{R}^{K \times K}\) 是因子收益协方差，\(\Delta = \mathrm{diag}(\delta_1^2, \dots, \delta_N^2)\) 是特异方差对角阵。这一拆，自由度从 \(\frac{N(N+1)}{2}\) 降到 \(\frac{K(K+1)}{2} + N\)。\(K\) 取 40，\(N\) 取 4000，参数从 800 万压到约 4820，压缩比超过 1600 倍。

这个压缩不是无代价的——它依赖「特异收益独立」这个假设，而这个假设在现实里经常被违反：同行业、同概念的股票特异收益强相关，比如 2020 年 7 月某天医美概念股集体涨停，特异收益矩阵里医美板块那个区块根本不对角。Barra 的工程对策是把行业因子尽可能细化（USE4 用 60 个细分行业，CNE6 用约 30 个），把绝大部分跨股票相关性吸进行业因子里，剩下的特异部分才当作独立。但这个对策在小盘股、概念股、ST 股上仍然不充分，导致风险模型在这些边缘资产上系统性低估风险。

估计的几个工程细节

实操中 \(\Sigma_f\) 与 \(\Delta\) 的估计有几条不写出来就做不对的细节：

第一，因子收益协方差的半衰期。\(\Sigma_f\) 不是简单对历史因子收益序列算协方差，而是用指数加权（EWMA），半衰期通常取 90 至 180 个交易日。Barra USE4 文档建议方差半衰期取 90 天、相关性半衰期取 480 天——方差变化快，相关性变化慢，分开估计能减小估计噪声。

第二，Newey-West 自相关修正。日频因子收益存在短期自相关（多数因子的 1 阶自相关在 0.05 至 0.15 之间），如果只用同期协方差会低估方差。Barra 用 Newey-West（1987）滞后加权修正：

\[ \hat\Sigma_f^{NW} = \hat\Sigma_f(0) + \sum_{l=1}^{L} \left(1 - \frac{l}{L+1}\right) \left(\hat\Sigma_f(l) + \hat\Sigma_f(l)^\top\right) \]

\(L\) 通常取 5 至 10。这一修正能让月频预测方差与实际波动的偏差缩小一个数量级。

第三，特异波动的截尾与平滑。残差 \(u_i\) 的样本方差在 \(T\) 短或股票最近停牌时不稳定，工程上要做：(1) 用 EWMA 平滑（半衰期 60 天）；(2) 用同行业、相近市值股票的特异波动做贝叶斯收缩；(3) 对极端值（前 0.5%、后 0.5%）做截尾。Barra 把这一步叫作 “structural specific risk model”，用回归 \(\log \delta_i = \gamma^\top z_i + \epsilon_i\) 把特异波动建模为公司基本面属性 \(z_i\)（市值、杠杆、流动性等）的函数，避免对单只股票样本依赖太重。

第四，Bayesian shrinkage of factor returns。Barra USE4 在 \(f_k\) 估计层面也做了一次贝叶斯收缩——把横截面回归得到的因子收益向「同类因子的均值」收缩，对小行业、小风格尤其重要。这一步一般叫作 “factor return shrinkage” 或 “eigen-adjusted covariance”，在 MSCI Barra 文档里专门有一节讲。

下面这段代码是协方差估计的工程骨架，把样本、shrinkage、因子模型三种估计放在一起对比：

import numpy as np
import pandas as pd
from sklearn.covariance import LedoitWolf

def sample_cov(returns: pd.DataFrame) -> np.ndarray:
    """样本协方差，作为基准"""
    return returns.cov().values

def ledoit_wolf_cov(returns: pd.DataFrame) -> tuple[np.ndarray, float]:
    """Ledoit-Wolf 收缩协方差，先验为常数相关阵"""
    lw = LedoitWolf().fit(returns.values)
    return lw.covariance_, lw.shrinkage_

def factor_model_cov(
    exposures: np.ndarray,
    factor_returns: pd.DataFrame,
    specific_var: np.ndarray,
    halflife_var: int = 90,
    halflife_cor: int = 480,
) -> np.ndarray:
    """
    Barra 风格因子模型协方差：Σ = X Σ_f X^T + Δ
    - exposures: (N, K) 因子暴露矩阵
    - factor_returns: (T, K) 因子日收益
    - specific_var: (N,) 特异方差
    - 用两种半衰期分别估方差与相关性
    """
    fr = factor_returns.values
    T, K = fr.shape

    weights_var = 0.5 ** (np.arange(T)[::-1] / halflife_var)
    weights_cor = 0.5 ** (np.arange(T)[::-1] / halflife_cor)
    weights_var /= weights_var.sum()
    weights_cor /= weights_cor.sum()

    mean = (fr * weights_var[:, None]).sum(axis=0)
    centered = fr - mean

    var_diag = (centered ** 2 * weights_var[:, None]).sum(axis=0)
    cor = np.zeros((K, K))
    std = np.sqrt(var_diag)
    for i in range(K):
        for j in range(K):
            cov_ij = (centered[:, i] * centered[:, j] * weights_cor).sum()
            cor[i, j] = cov_ij / (std[i] * std[j] + 1e-12)
    np.fill_diagonal(cor, 1.0)

    sigma_f = cor * np.outer(std, std)

    sigma = exposures @ sigma_f @ exposures.T + np.diag(specific_var)
    sigma = (sigma + sigma.T) / 2
    return sigma

这个骨架已经覆盖了 90% 的细节。剩下的 10%（Newey-West 修正、特异波动结构化建模、eigenfactor 调整）属于研究型团队的工作，工程上一开始可以忽略，等模型上线后再迭代。

三、Barra 多因子结构

Barra 的多因子模型把所有股票收益结构化拆成四块：国家因子（country）、行业因子（industry）、风格因子（style）、特异收益（idiosyncratic）。这一节把每一块的语义、数学、估计方法讲透。

国家因子

对单一市场（如 CNE6 只覆盖 A 股）模型，国家因子退化为常数：每只股票的国家暴露恒为 1，对应一个常数项。这相当于截距项，吸收所有股票共享的市场水平变动。对全球模型（如 Barra GEM3、MSCI ACWI），国家因子是一组哑变量，每只股票的国家暴露按主上市地或主营业务所在地确定。

国家因子的实际收益 \(f_{\text{country}}\) 解释力极强——A 股市场任何一天，这个因子常常解释 30% 至 60% 的横截面方差。它本质上是「今天市场涨跌」的统计回声。

行业因子

行业因子用哑变量编码：每只股票在所属行业上的暴露为 1，其它行业为 0。GICS 分到子行业一共约 158 个，CSI / 中信一级行业约 30 个，二级行业约 110 个。Barra USE4 用 60 个细分行业，CNE6 用约 32 个。

行业暴露在多元化集团股（如平安、中信、复星）上有歧义——一个公司同时有保险、银行、证券业务，暴露怎么分？两种主流做法：

1. 主行业法：按主营业务收入占比最高的行业打 1，其余为 0。优点简单，缺点是丢掉了多元化信息。
2. 加权暴露法：按主营业务收入占比拆分到多个行业，所有暴露之和仍为 1。Barra USE4 用这种方法。

行业因子的工程坑在于行业切换：一家公司从「化工」改名重组成「新能源」之后，历史时序上的行业暴露会跳变，导致因子收益序列出现伪信号。处理办法是 (1) 在 t-1 期用最新行业分类向前广播，避免数据泄漏；(2) 对行业切换日做 “industry rebalance” 标记，并在因子收益估计时排除该单只股票当天。

风格因子

风格因子是 Barra 模型的精髓，也是不同厂商差异最大的地方。MSCI Barra USE4 用 10 个风格因子（USE4-S）：

每个风格因子由若干个 descriptor 加权合成。每个 descriptor 在原始单位（市值是亿元、PB 是倍数）上没法直接当暴露用，必须做两步标准化：

1. 横截面标准化：在 t-1 期对所有股票，把 descriptor 减去截面均值（以市值平方根加权）、除以截面标准差，得到 z-score。
2. 截尾：把超过 ±3 倍标准差的截到 ±3，避免极端值（如刚 IPO 的微盘股）扭曲后面回归。

合成时用 descriptor 的加权和，权重通常取等权或按 IC（information coefficient，因子收益相关性）分配。

风格因子的工程坑在于多重共线性：Size、Beta、Volatility、Liquidity 之间相关性都比较高（市值小 → 流动性低 → 波动大 → beta 高）。Barra 的处理是 (1) Non-linear Size 与 Size 正交化；(2) Residual Volatility 与 Beta 正交化；(3) Liquidity 与 Size 正交化。剩下的相关性靠 WLS 回归本身的最小二乘性质处理。

横截面回归

每个交易日 \(t\)，在所有可交易股票（剔除停牌、新股、ST、涨跌停）上跑一次回归：

\[ r_{i,t} = \sum_{c} X^{\text{country}}_{i,c} f^{\text{country}}_{c,t} + \sum_{j} X^{\text{industry}}_{i,j} f^{\text{industry}}_{j,t} + \sum_{k} X^{\text{style}}_{i,k} f^{\text{style}}_{k,t} + u_{i,t} \]

注意暴露矩阵 \(X\) 用的是 t-1 期（前一日收盘）的值——绝对不能用当日的市值或当日的动量去解释当日收益，否则就是用未来信息预测过去，结果天上掉馅饼。

回归是加权最小二乘（WLS），权重取市值平方根 \(\sqrt{w_i}\)。为什么不是市值本身？因为大盘股的特异波动比小盘股小，按方差倒数加权应当用 \(1/\delta_i^2\)；但 \(\delta_i^2\) 与市值大致成 \(1/\sqrt{\text{cap}}\) 的关系（Heston-Rouwenhorst 1994），所以权重退化为 \(\sqrt{\text{cap}}\)。这一约定 Barra USE4 与 BARRA Cosmos 都遵循。

回归还需要额外约束以解决多重共线性——所有行业暴露之和恒为 1（每只股票主行业为 1，其余为 0，按加权暴露法时和也为 1），和国家因子的常数项 1 完全共线，矩阵不满秩。两种解决方案：

1. 舍弃一个行业作为基准：选第一个行业（按字母序或市值序）作为参考组，其余行业相对它估计。这是经济计量学传统做法，但解释起来不便（每个行业收益都是「相对基准的差」，不是「行业本身收益」）。
2. 加约束：要求行业因子收益的市值加权和为零（\(\sum_j w^{\text{ind}}_j f^{\text{industry}}_{j,t} = 0\)）。这是 Barra USE4 的做法，把约束加进 WLS 的拉格朗日方程。优点是每个行业收益都是「相对市场基准的超额」，归因清晰。

下面这段代码是横截面回归估计 Barra 风格因子收益的可运行骨架：

import numpy as np
import pandas as pd
from numpy.linalg import lstsq, pinv

STYLE_FACTORS = [
    "size", "beta", "momentum", "residual_volatility", "non_linear_size",
    "book_to_price", "liquidity", "earnings_yield", "growth", "leverage",
]

def winsorize(x: pd.Series, n_std: float = 3.0) -> pd.Series:
    mu = x.mean()
    sd = x.std()
    return x.clip(lower=mu - n_std * sd, upper=mu + n_std * sd)

def standardize(x: pd.Series, weights: pd.Series) -> pd.Series:
    """市值平方根加权标准化"""
    w = weights / weights.sum()
    mu = (x * w).sum()
    sd = np.sqrt(((x - mu) ** 2 * w).sum())
    return (x - mu) / (sd + 1e-12)

def build_exposure_matrix(
    descriptors: pd.DataFrame,
    industry: pd.Series,
    market_cap: pd.Series,
) -> tuple[pd.DataFrame, list[str], list[str]]:
    """
    构造一日的暴露矩阵 X：
    - 1 列国家因子（常数 1）
    - K_ind 列行业哑变量
    - K_style 列风格因子（标准化后的 descriptors）
    """
    weights = np.sqrt(market_cap)
    style_block = pd.DataFrame(index=descriptors.index)
    for col in STYLE_FACTORS:
        x = winsorize(descriptors[col].astype(float))
        style_block[col] = standardize(x, weights)

    industry_block = pd.get_dummies(industry, prefix="ind").astype(float)

    country_block = pd.DataFrame(
        {"country": np.ones(len(descriptors))}, index=descriptors.index
    )

    X = pd.concat([country_block, industry_block, style_block], axis=1)
    return X, list(industry_block.columns), STYLE_FACTORS

def cross_section_regression(
    returns: pd.Series,
    X: pd.DataFrame,
    market_cap: pd.Series,
    industry_cols: list[str],
) -> tuple[pd.Series, pd.Series]:
    """
    带行业市值加权和约束的 WLS：
        min ||W^{1/2}(r - X f)||^2
        s.t. sum_j w_j^{ind} f_j^{ind} = 0
    解析解通过拉格朗日法。
    """
    w = np.sqrt(market_cap.values)
    W_half = np.diag(w)

    Xv = X.values
    rv = returns.values

    A = W_half @ Xv
    b = W_half @ rv

    industry_idx = [X.columns.get_loc(c) for c in industry_cols]
    industry_cap = market_cap.groupby(X[industry_cols].idxmax(axis=1)).sum()
    industry_cap = industry_cap.reindex(industry_cols).fillna(0).values

    K = Xv.shape[1]
    R = np.zeros(K)
    R[industry_idx] = industry_cap / industry_cap.sum()

    AtA = A.T @ A
    Atb = A.T @ b
    KKT = np.block([
        [AtA, R[:, None]],
        [R[None, :], np.zeros((1, 1))],
    ])
    rhs = np.concatenate([Atb, [0.0]])
    sol = np.linalg.solve(KKT, rhs)
    f = sol[:K]

    factor_returns = pd.Series(f, index=X.columns)
    residuals = pd.Series(rv - Xv @ f, index=returns.index)
    return factor_returns, residuals

这段代码每天跑一次，得到当日的 1 + K_ind + K_style 个因子收益。把所有交易日的因子收益堆起来，就得到了 \(f_{k,t}\) 的时序矩阵；用前文的 EWMA 估计协方差，就完成了 \(\Sigma_f\) 的估计。

因子收益的解释

跑完一段时间，把每个因子的累积收益画出来，能直接看出市场结构。下面是 A 股 2016-2024 年的几个典型经验事实（来自笔者在某券商研究所内部的复现，量级与 MSCI Barra CNE6 报告一致）：

• Size 因子：2017 年大幅负收益（-25%，大盘股年）、2021 年 2 月「茅指数」崩盘后大幅正收益（小盘股回归），整体 9 年累积约 +15%，但波动巨大。
• Momentum 因子：A 股长期累积负收益（约 -20%），与美股长期为正完全相反。这反映 A 股散户主导的市场结构里，反转效应（reversal）大于动量效应。
• Volatility 因子：长期累积负收益（约 -40%），即「低波异象」（low volatility anomaly）在 A 股同样存在，且比美股更显著。
• Book-to-Price 因子：2014-2017 年正收益，2018-2021 年负收益（成长牛市），2022-2024 年再次正收益，呈现 5-7 年周期。

这些因子收益本身就是研究信号——一个长期负收益的因子，意味着「在该因子上空配能赚钱」；但要警惕这只是过去的样本，未来可能均值回归甚至反转。Barra 因子是风险模型的因子，不是直接拿来选股的 alpha 因子，请不要混淆。

四、特异风险与残差波动

横截面回归剩下的残差 \(u_{i,t}\) 是特异收益，对应的方差 \(\delta_i^2\) 是特异风险。这一节讲怎么把残差时序变成可用的特异波动估计。

朴素估计的问题

最朴素的做法：对每只股票，用过去 \(T\) 天残差的样本方差 \(\hat\delta_i^2 = \frac{1}{T-1}\sum_t u_{i,t}^2\)。这个估计有几个工程问题：

第一，停牌、新股、ST 摘帽这些事件让 \(T\) 实际可用样本变短，估计噪声大。

第二，单只股票的方差估计标准误是 \(\delta_i^2 \sqrt{2/(T-1)}\)，\(T=252\) 时相对误差约 9%，对个股来说能接受；但当组合优化器对 \(\delta_i^2\) 求倒数时（risk-parity 等），9% 的相对误差被放大，加权时小盘股的估计噪声更被放大。

第三，特异波动有时变性。市场进入高波动状态时，所有股票的特异波动同步上升，单只股票看不出这个共同冲击。

结构化特异风险模型

Barra 的对策是把特异波动建模为公司基本面属性的函数，让信息「跨股票借用」：

\[ \log \hat\delta_i^t = \gamma_0^t + \sum_l \gamma_l^t z_{i,l} + \epsilon_i^t \]

其中 \(z_{i,l}\) 是公司属性（log market cap、leverage、liquidity、book-to-price 等），\(\hat\delta_i^t\) 是 t 时刻样本估计的特异波动。把所有股票的 \(\log \hat\delta_i^t\) 与 \(z_{i,l}\) 跑一次横截面回归，得到系数 \(\gamma_l^t\)；然后用拟合值 \(\tilde\delta_i^t = \exp(\hat\gamma_0^t + \sum_l \hat\gamma_l^t z_{i,l})\) 作为「结构化估计」。

最终的特异波动是结构估计与样本估计的贝叶斯收缩：

\[ \delta_i^t = \theta_i \tilde\delta_i^t + (1 - \theta_i) \hat\delta_i^t \]

收缩权重 \(\theta_i\) 在样本量短（新股、停牌时间长）时偏向 1（更信任结构估计），样本充分时偏向 0（更信任样本估计）。Barra USE4 的具体公式给出了 \(\theta_i\) 是 \(h_i = T_i / 252\) 的函数。

时变波动：EWMA 与 GARCH

特异波动的时变性也要建模。最简单的方法是 EWMA：

\[ \hat\delta_i^{2,t} = \lambda \hat\delta_i^{2,t-1} + (1-\lambda) u_{i,t}^2 \]

\(\lambda\) 取 0.94（相当于半衰期约 11 天）适合短期，0.97（半衰期约 22 天）适合月频。RiskMetrics 1996 报告把 \(\lambda = 0.94\) 写进了行业标准。

更精细的是 GARCH(1,1)：

\[ \hat\delta_i^{2,t} = \omega_i + \alpha_i u_{i,t-1}^2 + \beta_i \hat\delta_i^{2,t-1} \]

GARCH 在单只股票上比 EWMA 略好（年化预测相对误差减少 5% 至 10%），但参数估计更费时间，且对 \(T\) 短的股票不稳定。工程上一般只在指数、ETF、大盘蓝筹这些样本充分的标的用 GARCH，一般个股用 EWMA。

特异波动的工程坑

第一，回归残差与特异收益不完全等价。回归残差包含了模型设定误差——如果某个未被模型覆盖的因子（比如「ESG 评分」、「龙虎榜资金流」）在某天有显著收益，所有暴露在这个因子上的股票残差会同向偏离，特异收益独立性假设被违反。这一点没有完美解法，工程上要定期检查残差矩阵是否有显著的低秩结构（用 PCA 看前几个主成分），如果有，要么把对应的隐藏因子加进模型，要么承认这部分跨股票相关性会被低估。

第二，停牌天数的处理。如果某只股票连续停牌 N 天，停牌期内 \(u_{i,t} = 0\)，把这个零计入样本会严重低估方差。正确做法：把停牌日从样本里剔除，重新计算有效 \(T_i\)；并把停牌期前后的「补涨补跌」也剔除（停牌期内市场跑了 -10%，复牌当天跌停，这个跌停不属于该日的特异收益）。

第三，新股、次新股的处理。上市不满 60 天的次新股没有足够样本，且经常一字涨停，特异波动失真。Barra 的做法是给次新股一个「同行业、相近市值的特异波动均值」作为先验，在样本积累足够之前完全用先验。

五、风险归因

风险模型最大的工程价值在归因——把组合的总风险拆成因子贡献和特异贡献，告诉基金经理「你的钱押在哪里」。这一节讲两类归因：事前风险归因（把预测方差拆开）和事后收益归因（把已实现收益拆开）。

事前风险归因：边际贡献与组合贡献

组合预测方差：

\[ \sigma_p^2 = w^\top \Sigma w = w^\top X \Sigma_f X^\top w + w^\top \Delta w = b^\top \Sigma_f b + \sum_i w_i^2 \delta_i^2 \]

其中 \(b = X^\top w\) 是组合的因子暴露向量。这一公式自然把方差拆成两块：因子方差 \(b^\top \Sigma_f b\) 和特异方差 \(\sum_i w_i^2 \delta_i^2\)。

但只拆到「因子 vs 特异」远远不够，基金经理还想知道哪个具体因子贡献最大。这里要用「边际风险贡献」（marginal contribution to risk，MCR）和「组合风险贡献」（component contribution to risk，CCR）：

\[ \text{MCR}_k = \frac{\partial \sigma_p}{\partial b_k} = \frac{(\Sigma_f b)_k}{\sigma_p} \]

\[ \text{CCR}_k = b_k \cdot \text{MCR}_k = \frac{b_k (\Sigma_f b)_k}{\sigma_p} \]

CCR 满足 \(\sum_k \text{CCR}_k + \sum_i w_i^2 \delta_i^2 / \sigma_p = \sigma_p\)，即所有 CCR 加起来恰好等于组合波动。这是欧拉齐次定理在风险度量上的应用——\(\sigma_p(w)\) 是 \(w\) 的一次齐次函数。

CCR 报告告诉基金经理：组合年化波动 12%，其中市值因子贡献 4%、电子行业贡献 3%、动量因子贡献 1.5%、其他因子合计 1.5%、特异部分 2%。如果市值因子贡献过高，说明组合实质上是一个「小盘股押注」，与号称的「选股模型」不符。

事后收益归因：Brinson-Fachler 与多因子

Brinson-Fachler（1985、1986）框架把收益归因拆成行业配置（allocation）和选股（selection）：

\[ r_p - r_b = \underbrace{\sum_j (w^p_j - w^b_j) (r^b_j - r^b)}_{\text{配置}} + \underbrace{\sum_j w^p_j (r^p_j - r^b_j)}_{\text{选股}} \]

其中 \(w^p_j, w^b_j\) 是组合与基准在行业 \(j\) 的权重，\(r^p_j, r^b_j\) 是组合与基准在行业 \(j\) 的实现收益。

Brinson 框架直观，但只拆到行业一层，做不到风格因子归因。多因子归因把组合实现收益拆成：

\[ r_p = \underbrace{\sum_k b_{p,k} f_k}_{\text{因子收益}} + \underbrace{\sum_i w_i u_i}_{\text{选股收益（特异）}} \]

每个因子的贡献等于「组合在该因子上的暴露 × 当期因子收益」。这种归因是 Barra 模型的标准输出，比 Brinson 多了风格归因这一层，业内已经成为标准。

下面这段代码做组合风险归因：

import numpy as np
import pandas as pd

def portfolio_risk_attribution(
    weights: pd.Series,
    exposures: pd.DataFrame,
    sigma_f: pd.DataFrame,
    specific_var: pd.Series,
) -> pd.DataFrame:
    """
    风险归因 component contribution to risk
    返回每个因子与特异部分的方差贡献、波动率贡献、占比
    """
    w = weights.reindex(exposures.index).fillna(0).values
    X = exposures.values
    Sf = sigma_f.values
    delta = specific_var.reindex(exposures.index).fillna(0).values

    b = X.T @ w
    factor_var = b @ Sf @ b
    specific_total = float((w ** 2 * delta).sum())
    total_var = factor_var + specific_total
    sigma_p = np.sqrt(total_var)

    Sf_b = Sf @ b
    factor_ccr_var = b * Sf_b
    factor_ccr_vol = factor_ccr_var / sigma_p

    rows = []
    for i, name in enumerate(exposures.columns):
        rows.append({
            "name": name,
            "exposure": float(b[i]),
            "var_contrib": float(factor_ccr_var[i]),
            "vol_contrib": float(factor_ccr_vol[i]),
            "share": float(factor_ccr_var[i] / total_var),
        })
    rows.append({
        "name": "specific",
        "exposure": float("nan"),
        "var_contrib": specific_total,
        "vol_contrib": specific_total / sigma_p,
        "share": specific_total / total_var,
    })
    rows.append({
        "name": "TOTAL",
        "exposure": float("nan"),
        "var_contrib": total_var,
        "vol_contrib": sigma_p,
        "share": 1.0,
    })
    return pd.DataFrame(rows)


def portfolio_return_attribution(
    weights: pd.Series,
    exposures: pd.DataFrame,
    factor_returns: pd.Series,
    asset_returns: pd.Series,
) -> pd.DataFrame:
    """
    收益归因：分解组合实现收益为因子贡献和特异贡献
    """
    w = weights.reindex(exposures.index).fillna(0).values
    X = exposures.values
    f = factor_returns.reindex(exposures.columns).fillna(0).values
    r = asset_returns.reindex(exposures.index).fillna(0).values

    b = X.T @ w
    factor_contrib = b * f
    factor_total = float(b @ f)

    r_explained_per_asset = X @ f
    specific_per_asset = r - r_explained_per_asset
    specific_total = float((w * specific_per_asset).sum())

    rows = []
    for i, name in enumerate(exposures.columns):
        rows.append({
            "name": name,
            "exposure": float(b[i]),
            "factor_return": float(f[i]),
            "contribution": float(factor_contrib[i]),
        })
    rows.append({
        "name": "specific",
        "exposure": float("nan"),
        "factor_return": float("nan"),
        "contribution": specific_total,
    })
    rows.append({
        "name": "TOTAL",
        "exposure": float("nan"),
        "factor_return": float("nan"),
        "contribution": factor_total + specific_total,
    })
    return pd.DataFrame(rows)

这两段代码是风险报告的核心引擎。生产环境里，每天收盘后跑一次，把当天的因子贡献、累计的因子贡献、年初至今的因子贡献堆成一张表，发给基金经理与风控。

归因的几个常见误用

第一，把因子归因当 alpha 来源。如果一个组合的 95% 收益来自因子（特别是 Size、Momentum），它本质上是一个「智能 beta」组合，而不是一个 alpha 组合。这种组合该收的是 0.5% 管理费而不是 2% 管理费。归因报告是 LP 鉴别「真假 alpha」的核心工具。

第二，归因结果不闭合。因子收益 + 特异收益应当与组合实际收益严格相等（不考虑交易成本和现金部分）。如果不等，说明 (1) 暴露矩阵 X 用错了日期（用了当日而不是 t-1）；(2) 因子收益估计与暴露不匹配；(3) 组合权重没用日初权重。一旦不闭合，整张报告作废，必须查到原因再报。

第三，跨期累计的非线性。日度归因可以加总到月度吗？答案是「对数收益可以、简单收益不行」。因为简单收益累乘不是累加，因子贡献的累乘需要做几何平均处理（Carino 1999 的 logarithmic linking）。生产里常见的偷懒做法是直接加总日度贡献，这在月度上误差能控制在 5% 以内；但在年度归因上误差累积，必须用 Carino 或 Menchero 链接公式。

六、VaR 与 ES

风险价值（value at risk，VaR）和预期短缺（expected shortfall，ES，又称 conditional VaR）是把组合风险压缩成一个数字的两种方式。它们在监管报送、对外披露、内控限额上有刚性需求，但在工程实现上有几个深坑。

VaR 的定义与三种估计法

VaR 的定义：在置信水平 \(\alpha\)（通常取 95% 或 99%）下，未来 \(h\) 天的损失不超过 VaR 的概率为 \(\alpha\)：

\[ \Pr(L_{t,t+h} \leq \text{VaR}_\alpha(h)) = \alpha \]

工业上有三类估计法：参数法、历史模拟法、Monte Carlo 法。

参数法（variance-covariance method）：假设组合收益服从均值 \(\mu_p\)、波动 \(\sigma_p\) 的正态分布（或 t 分布），VaR 直接由分位数公式给出：

\[ \text{VaR}_\alpha(h) = -\mu_p \cdot h - \sigma_p \cdot \sqrt{h} \cdot \Phi^{-1}(1-\alpha) \]

\(\Phi^{-1}\) 是标准正态分布反函数。99% 置信水平下，\(\Phi^{-1}(0.01) \approx -2.326\)。这种方法简单快、可以用 Barra 的 \(\Sigma\) 直接代入，但它假设收益正态——对单日股票收益严重不成立，对组合（中心极限定理稍稍生效）也只能近似。1998 年 LTCM 倒闭的核心原因之一就是参数 VaR 严重低估了尾部风险。

历史模拟法（historical simulation）：拿过去 \(T\) 天的收益序列直接当作未来 1 天收益的经验分布，VaR 取样本第 \((1-\alpha)\) 分位数：

\[ \text{VaR}_\alpha = -\text{Quantile}_{1-\alpha}(\{r_{p, t-T+1}, \dots, r_{p, t}\}) \]

这种方法不假设分布形式，能自然刻画偏度和厚尾；但它是「过去多少天就只能看多少天的最坏情况」，对结构性变化的反应慢，且 \(T\) 短时尾部估计不准。监管要求一般是 \(T \geq 250\) 天。

Monte Carlo 法：从拟合好的分布（如多元 t 分布、Cornish-Fisher 展开、copula 模型）里采样 \(M\) 条路径（\(M\) 通常取 10,000 至 100,000），按组合定价规则把每条路径的损失算出来，取分位数。对衍生品组合（含期权、可转债、CDS）是唯一可行的方法，因为这些工具的损益是非线性的。代价是计算成本极高——重定价 10,000 次需要几分钟到几小时。

三种方法的工程对比

工程上的标准做法是三种方法并行：参数法每分钟刷新一次（实时风控），历史模拟法每收盘后刷新（监管报送），Monte Carlo 每日收盘后跑一次（含期权头寸的全口径风险）。三者结果差异如果超过阈值（一般取 30%），自动触发告警，由风控人员人工复核——这是发现模型失效的第一道防线。

历史模拟法的实现

历史模拟法对工程友好，下面给出可直接用的实现：

import numpy as np
import pandas as pd

def historical_var(
    returns: pd.Series,
    alpha: float = 0.99,
    horizon: int = 1,
) -> dict:
    """
    历史模拟法 VaR / ES
    - returns: 组合日收益时序
    - alpha: 置信水平
    - horizon: 时间窗口（天）
    返回 VaR 与 ES（均为正数，表示损失）
    """
    if horizon > 1:
        rolling = returns.rolling(horizon).sum().dropna()
    else:
        rolling = returns.dropna()

    losses = -rolling.values
    var_quantile = np.quantile(losses, alpha)
    tail = losses[losses >= var_quantile]
    es = tail.mean() if len(tail) > 0 else var_quantile

    return {
        "VaR": float(var_quantile),
        "ES": float(es),
        "n_obs": int(len(losses)),
        "n_tail": int(len(tail)),
        "alpha": alpha,
        "horizon": horizon,
    }


def historical_var_with_volatility_scaling(
    returns: pd.Series,
    alpha: float = 0.99,
    halflife: int = 30,
) -> dict:
    """
    Hull-White (1998) 波动率调整历史模拟
    把每个历史观测按 当前波动率 / 历史波动率 重新缩放
    避免历史模拟在低波动期低估、高波动期高估
    """
    r = returns.dropna().values
    weights = 0.5 ** (np.arange(len(r))[::-1] / halflife)
    weights /= weights.sum()
    mean = (r * weights).sum()
    var_t = ((r - mean) ** 2 * weights).sum()
    sigma_now = np.sqrt(var_t)

    sigma_hist = np.zeros(len(r))
    sigma_hist[0] = sigma_now
    for t in range(1, len(r)):
        prev = r[:t]
        w_prev = 0.5 ** (np.arange(t)[::-1] / halflife)
        w_prev /= w_prev.sum()
        m = (prev * w_prev).sum()
        sigma_hist[t] = np.sqrt(((prev - m) ** 2 * w_prev).sum() + 1e-12)

    scaled = (r - mean) * (sigma_now / sigma_hist) + mean
    losses = -scaled
    var_quantile = np.quantile(losses, alpha)
    tail = losses[losses >= var_quantile]
    es = tail.mean() if len(tail) > 0 else var_quantile

    return {
        "VaR": float(var_quantile),
        "ES": float(es),
        "sigma_now": float(sigma_now),
        "alpha": alpha,
    }


def parametric_var(
    weights: np.ndarray,
    sigma: np.ndarray,
    alpha: float = 0.99,
    horizon: int = 1,
) -> dict:
    """
    参数法 VaR：基于因子模型协方差矩阵
    """
    from scipy.stats import norm
    var_p = float(weights @ sigma @ weights)
    sigma_p = np.sqrt(var_p)
    z = norm.ppf(1 - alpha)
    var_amount = -z * sigma_p * np.sqrt(horizon)
    es_amount = sigma_p * np.sqrt(horizon) * norm.pdf(z) / (1 - alpha)
    return {
        "VaR": float(var_amount),
        "ES": float(es_amount),
        "sigma_p": float(sigma_p),
        "alpha": alpha,
    }

ES 比 VaR 更稳健的几个原因

巴塞尔银行监管委员会（BCBS）从 2016 年开始（巴塞尔 III 修订版、FRTB 框架）就把市场风险资本计量从 99% VaR 切换到 97.5% ES。原因有三：

第一，VaR 不是相干风险度量（coherent risk measure）。Artzner et al.（1999）证明 VaR 不满足次可加性（subadditivity）—— 把两个组合合并后的 VaR 可能大于两个组合各自 VaR 之和，这违反了「分散化降低风险」的金融常识。ES 是相干的。

第二，VaR 对尾部分布不敏感。VaR 只问「分位数处的损失是多少」，不问「超过分位数时平均损失多大」。两个分布在 99% 分位数处 VaR 相同，但一个尾部短一个尾部极长，VaR 报同一个数；ES 能区分。LTCM、AIG、Bear Stearns 出事时，账面 VaR 没超阈值，但实际损失超过 VaR 数倍——这正是 ES 想解决的问题。

第三，ES 更适合监管资本计量。监管的目标是确保银行在最坏情况下还能履约，而不是最坏 1% 不履约时只损失 VaR 那么多。

工程上 ES 的估计同样用历史模拟即可：取尾部 \(1-\alpha\) 比例的损失均值。代码已经包含在上面 historical_var 的实现里。

Backtesting 与 Kupiec 检验

VaR 模型必须做回测验证：在过去 \(T\) 天里，实际损失超过 VaR 的次数应当近似等于 \(T(1-\alpha)\)。Kupiec（1995）的似然比检验给出统计判定：

\[ \text{LR}_{POF} = -2\log\left[\frac{(1-\alpha)^x \alpha^{T-x}}{(x/T)^x ((T-x)/T)^{T-x}}\right] \]

其中 \(x\) 是实际超阈值次数，\(T\) 是回测样本数。\(\text{LR}_{POF}\) 渐进服从 \(\chi^2(1)\)，95% 置信水平下，超过 3.84 拒绝「VaR 模型正确」的原假设。

巴塞尔委员会还要求做条件覆盖检验（Christoffersen 1998），即超阈值事件不仅频率正确、还要在时间上独立（不能扎堆）。模型失败时，监管会调高风险资本系数（traffic light approach：绿灯 ≤4 次、黄灯 5-9 次、红灯 ≥10 次，红灯一年内不能用内部模型）。

七、压力测试与情景分析

VaR 与 ES 都是从历史分布或参数分布里抽出来的「正常时期统计量」，对真正的极端事件——黑天鹅、流动性危机、政策转向——有系统性低估。压力测试（stress testing）是直接构造极端情景，逐情景重新定价组合，看损失会有多大。它和 VaR 是互补关系，不是替代关系。

下面这张图是压力测试的整体情景树，展示三类情景与共同冲击层的关系。

历史情景：复现真实事件

把历史上发生过的极端事件（雷曼、欧债、Covid、A 股 2015 股灾、2024/02 量化踩踏）当作模板，把当时的因子收益、行业收益、波动率、相关性矩阵当作冲击向量，加到当前组合上重定价。

历史情景的优点是真实可信——这件事真的发生过；缺点是「过去的极端不代表未来」，且每次危机的细节都不同。最重要的几个 A 股历史情景：

每个情景对应一个因子冲击向量 \(\Delta f\)，组合损失为：

\[ L_{\text{scen}} = -w^\top X \Delta f - w^\top \Delta u \]

其中 \(\Delta u\) 是该情景下个股的特异冲击。工程上把所有历史情景做成一张「stress library」，每周日例行跑一遍，看当前组合在哪个情景下损失最大。

假设情景：构造极端冲击

历史情景有局限——某些情景在历史上没发生过但理论可能发生。假设情景由风控人员手工构造：

• 利率冲击：10 年国债收益率瞬间上行 200 bp。对债券、银行、地产板块影响最大。
• 指数冲击：沪深 300 单日 -10%，全市场 beta 1.0 资产同步下跌。
• 风格反转：动量因子 -3σ、价值因子 +3σ，专门冲击「成长 + 动量」组合。
• 流动性枯竭：所有股票市场冲击成本 ×3，可成交规模 ×0.3。
• 行业政策冲击：某个行业（教培、地产、医药集采）单日 -25%，其余行业不变。
• 汇率冲击：人民币对美元 ±5%，影响出口、进口、外汇敞口。

每个情景对应 \(\Delta f\) 向量，组合风险按上述公式重算。假设情景的关键挑战是「哪些因子同时动、动多少」，这要求风控人员对市场结构有深入理解，不能拍脑袋。

反向情景：从损失反推风险源

反向情景（reverse stress test）思路相反：给定一个亏损目标（如 -10%），求解最有可能导致这个亏损的因子冲击组合。这相当于解一个约束优化：

\[ \min_{\Delta f} \frac{1}{2} \Delta f^\top \Sigma_f^{-1} \Delta f \quad \text{s.t.} \quad b^\top \Delta f \geq L^* \]

最小化「冲击的 Mahalanobis 距离」（即在因子收益分布下的不可能程度）使得损失达到 \(L^*\)。求解结果告诉风控人员：「要让组合亏 10%，最容易发生的方式是市值因子 +2σ、价值因子 -1.5σ、电子行业 -3σ」。这种情景特别适合监管报送（CCAR、SPRA 都要求做反向压力测试）。

尾部相关性

所有压力测试都要面对一个核心问题：极端时期资产相关性会上升。Longin-Solnik（2001）、Ang-Chen（2002）的实证表明，市场下跌 3σ 时，全球股票相关性从平均 0.3 抬升到 0.6 以上；下跌 5σ 时接近 0.85。这意味着「正常时期分散化」的组合在危机里完全不分散——你以为你买了 50 个不相关的资产，危机来时 50 个一起跌。

工程对策：在压力情景里，对相关矩阵做尾部上调：

\[ \rho^{\text{stress}} = \rho + \kappa (1 - \rho) \cdot \mathbb{1}[\text{tail}] \]

\(\kappa\) 取 0.3 至 0.5。或者更精细地，用 t-copula、Clayton copula 显式建模下尾相关性。这部分的细节属于学术前沿，工程上一般取经验上调系数即可。

监管框架：FRTB、CCAR、SPRA

巴塞尔 III 修订版（FRTB，Fundamental Review of the Trading Book，2019 年最终版）要求银行用 ES（97.5% 置信水平）替代 VaR，并对每个风险因子分类做「敏感度法」（standardized approach，SA）或「内部模型法」（internal models approach，IMA）。IMA 比 SA 资本占用低 30% 至 50%，但要求模型通过 PnL Attribution 测试和回测，且任何因子不通过都要回退到 SA。

美联储的 CCAR（Comprehensive Capital Analysis and Review）、欧央行的 SPRA、中国银保监会的资本管理办法都各自有压力测试要求。对资管机构（共同基金、私募）来说，监管压力相对小，但内部风控仍应建立类似框架。

八、工程实现

把上面所有数学和代码连成一个能上线的系统，是真正的工程挑战。这一节讲生产环境的几个关键决策。

每日刷新流水线

风险模型的每日流水线（pipeline）大致如下：

1. T 日 15:30 收盘后：等待官方收盘价、复权因子、行业分类、财务数据 ETL 完成（一般 17:00-19:00 完成）。
2. 19:30：构造 t-1 期暴露矩阵 X（用 t-1 期收盘后的市值、t-1 期最新财务、t-1 期累计动量等）。
3. 20:00：跑横截面回归，估计当日因子收益 \(f_t\) 与残差 \(u_t\)。
4. 20:15：把当日因子收益拼到历史时序上，重新估计因子协方差 \(\Sigma_f\) 与特异波动 \(\delta_i\)。
5. 20:30：合成股票协方差矩阵 \(V = X \Sigma_f X^\top + \Delta\)，写入风险数据库。
6. 20:45：对当前持仓跑风险归因、VaR、ES、压力测试，生成报告。
7. 21:00：报告通过邮件 / Slack / 内部 Web 推送给风控、基金经理。

整个流水线必须有依赖管理（哪一步在哪一步之后跑）、失败重试（数据没到怎么办）、版本控制（今天的模型版本是 V4.2.1）、回滚机制（V4.2.1 跑炸了切回 V4.2.0）。生产里通常用 Airflow、Prefect、Dagster 这类调度框架。

模型版本与回测一致性

风险模型迭代是常态——每 6 至 12 个月会有一次大版本升级（加新因子、调整权重、改估计方法）。但回测必须用「当时的模型」而不是「最新的模型」——这是不可妥协的工程纪律。

具体来说：

• 2019 年的回测用 V3.x，2020-2022 用 V4.x，2023 至今用 V5.x。
• 每个版本的因子定义、行业分类、估计参数、协方差快照都要存档，永不删除。
• 跨版本对比时，必须把同一区段用两个版本各跑一遍，对差异给出解释。

这套机制叫作 “point-in-time” 风险模型（Barra 内部叫 PIT model）。没有 PIT，回测结果不可信，模型迭代也不可控。工程实现上，每天的因子暴露、因子收益、协方差矩阵、特异波动都按日存盘（HDF5、Parquet 或专用风险数据库），并打上模型版本 tag。

与组合优化器对接

组合优化器的输入是 \(\Sigma\) 矩阵和因子暴露 \(X\)，输出是权重 \(w\)。三者间的接口要严格：

class RiskModelSnapshot:
    def __init__(
        self,
        date: pd.Timestamp,
        version: str,
        universe: list[str],
        exposures: pd.DataFrame,  # (N, K)
        factor_cov: pd.DataFrame,  # (K, K)
        specific_var: pd.Series,   # (N,)
    ):
        self.date = date
        self.version = version
        self.universe = universe
        self.exposures = exposures
        self.factor_cov = factor_cov
        self.specific_var = specific_var

    def asset_cov(self) -> np.ndarray:
        X = self.exposures.values
        Sf = self.factor_cov.values
        delta = self.specific_var.values
        Sigma = X @ Sf @ X.T + np.diag(delta)
        return (Sigma + Sigma.T) / 2

    def hash(self) -> str:
        import hashlib
        h = hashlib.sha256()
        h.update(self.date.isoformat().encode())
        h.update(self.version.encode())
        h.update(self.exposures.values.tobytes())
        h.update(self.factor_cov.values.tobytes())
        h.update(self.specific_var.values.tobytes())
        return h.hexdigest()[:16]

每个优化任务的输入快照都存档（带 hash），任何「同样输入得到不同输出」的现象都说明优化器或数据流有非确定性 bug。

在回测中嵌入风险模型

回测里嵌入风险模型有两种模式：

模式 A：在 alpha 信号之后做风控筛选。先用 alpha 信号选出 100 只候选股票，再用风险模型检查每只股票的因子暴露，剔除「过度暴露在不想要的因子上」的股票。这种模式简单，但在二次规划层面没有最优化。

模式 B：在二次规划中显式优化。把 alpha 信号当作期望收益 \(\mu\)，把风险模型协方差当作 \(\Sigma\)，求解：

\[ \max_w \quad \mu^\top w - \frac{\lambda}{2} w^\top \Sigma w \]

加上各种约束（行业偏离 ≤ 5%、风格偏离 ≤ 0.3、跟踪误差 ≤ 4%、单股 ≤ 2%）。这是 Markowitz 框架的标准形式，cvxpy / mosek 都能秒解。回测时每个调仓日跑一次，把当日的风险模型快照塞进去。

模式 B 的工程坑：

• 协方差矩阵规模：4000 只股票的 \(\Sigma\) 是 4000×4000，cvxpy 直接处理需要几十秒；用 \(\Sigma = X \Sigma_f X^\top + \Delta\) 的因子分解形式（cvxpy 1.4+ 支持），可以加速 10 倍以上。
• 优化数值稳定性：\(\Sigma\) 必须半正定（PSD）。因子分解形式自然 PSD（因为 \(\Sigma_f\) 半正定 + \(\Delta\) 正对角），直接用样本协方差时要做 PSD projection。
• 求解失败处理：对于约束矛盾（如 alpha 全负但要求满仓多头）的情况，优化器会返回 infeasible。回测必须对 infeasible 有降级策略——回退到上一日权重、或解放松版本（去掉一些约束）。

性能与基础设施

一个覆盖 A 股 4000 只股票、过去 10 年日频数据、含 40 个因子的风险模型，每天的计算量大约是：

• 横截面回归：4000 × 50 矩阵 WLS，单次约 50 ms。
• 协方差矩阵估计：50 × 50 EWMA，约 10 ms。
• 4000×4000 资产协方差合成：约 200 ms。
• 风险归因：单组合 20 ms，100 个组合并行 2 秒。
• VaR 历史模拟：单组合 5 ms。
• 压力测试 100 情景：单组合 100 ms。

整个每日流水线能在单台 16 核机器上 30 分钟内跑完。存储上，10 年历史快照（暴露 + 因子收益 + 协方差 + 特异波动）大约 50 GB。基础设施要求不高，但对数据完整性的要求极高——任何一天数据缺失或错误都会污染当天的因子估计、并通过协方差矩阵传播到后续几个月。

监控与告警

生产风险模型必须有持续监控：

• 因子收益的稳定性：每天的因子收益与过去 60 天均值偏差超过 ±3σ 触发告警（市场结构变化预警）。
• R²的稳定性：横截面回归的 \(R^2\) 突然下降（比如从 0.4 跌到 0.15）说明模型解释力崩塌，可能是数据错误或市场结构性突变。
• 特异波动分布：所有股票特异波动的分位数（25%、50%、75%）每天画一张图，分布漂移就告警。
• VaR 实际超出次数：每周看一次过去 20 天里实际损失超过 VaR 的次数，超过 1 次（99% VaR 期望 0.2 次）就重点关注。
• 协方差矩阵条件数：condition number 突然升高（>1e6）说明矩阵接近奇异，要么数据有问题、要么模型设定有问题。

这些监控组合起来构成风险模型的「自检系统」，是生产风险模型与玩具代码的核心区别。

九、写在最后

风险模型不是一个数学问题，是一个工程系统问题。Barra 多因子结构、shrinkage 协方差、historical-simulation VaR、stress testing 这套东西，每一块单独看都很简单，难的是把它们连成一个每天稳定跑、每年迭代一次、对 LP 报告闭合、对优化器友好、对回测一致的生产系统。

实战中真正决定胜负的几条：

第一，数据 hygiene 比模型先进。点对点数据时间戳错乱、复权因子算错、行业分类切换没标记、停牌股没剔除，任何一个环节出问题都会让最先进的模型变成噪声。先把数据质量搞到 99.9%，再考虑改模型。

第二，Point-in-time 是底线。任何「用今天的因子定义回算 5 年前因子收益」的操作都是数据泄漏，得到的回测结果完全不可信。所有版本、参数、定义、分类都必须按时间存档。

第三，VaR 不能替代压力测试。任何只看 VaR 不看压力测试的风控团队，都在为下一次黑天鹅交学费。VaR 是「正常时期的风险尺子」，压力测试是「极端情景的损失模拟器」，缺一不可。

第四，风险模型的输出是给人看的。归因报告、压力测试结果、VaR 监控图，最终都要让基金经理与风控人员看懂、信服、采取行动。一个跑得飞快但报告写成天书的系统，价值远低于一个慢一点但报告清晰的系统。

第五，不要追求 99% 完美。Barra 的官方文档每个细节展开都能讲一本书，但实际生产里 80% 的功能能解决 95% 的问题。先把基础流水线跑起来、跑稳、跑对，再迭代细节。

最后再次重申合规边界：本文所有代码与示例均为教学用途，未经过完整生产审计。商用风险模型应使用 MSCI Barra、Axioma、Bloomberg PORT 等厂商授权产品，或自建经过严格回测、版本管理、监管报备的内部模型。VaR 与压力测试涉及监管报送（FRTB、CCAR、SPRA、中国银保监会资本管理办法）的部分有专门的口径要求，本文不替代任何正式监管文档。任何把示例代码部署到生产环境的尝试都需要经过完整审计与法律审核，作者不对此承担任何责任。

参考文献

1. Rosenberg, B., & Marathe, V. (1976). “Common Factors in Security Returns: Microeconomic Determinants and Macroeconomic Correlates.” Proceedings of the Seminar on the Analysis of Security Prices, University of Chicago.
2. Grinold, R. C., & Kahn, R. N. (1999). Active Portfolio Management: A Quantitative Approach for Producing Superior Returns and Controlling Risk. McGraw-Hill, 2nd edition.
3. MSCI Barra (2011). “The Barra US Equity Model (USE4) Methodology Notes.” MSCI Research Notes.
4. MSCI Barra (2018). “The Barra China Equity Model (CNE6) Empirical Notes.” MSCI Research Notes.
5. Menchero, J., Orr, D. J., & Wang, J. (2011). “The Barra US Equity Model (USE4) Methodology Notes.” MSCI.
6. Ledoit, O., & Wolf, M. (2003). “Improved Estimation of the Covariance Matrix of Stock Returns with an Application to Portfolio Selection.” Journal of Empirical Finance, 10(5).
7. Ledoit, O., & Wolf, M. (2004). “Honey, I Shrunk the Sample Covariance Matrix.” Journal of Portfolio Management, 30(4).
8. Ledoit, O., & Wolf, M. (2017). “Nonlinear Shrinkage of the Covariance Matrix for Portfolio Selection: Markowitz Meets Goldilocks.” Review of Financial Studies, 30(12).
9. Newey, W. K., & West, K. D. (1987). “A Simple, Positive Semi-Definite, Heteroskedasticity and Autocorrelation Consistent Covariance Matrix.” Econometrica, 55(3).
10. Marčenko, V. A., & Pastur, L. A. (1967). “Distribution of Eigenvalues for Some Sets of Random Matrices.” Mathematics of the USSR-Sbornik, 1(4).
11. Michaud, R. O. (1989). “The Markowitz Optimization Enigma: Is ‘Optimized’ Optimal?” Financial Analysts Journal, 45(1).
12. Brinson, G. P., & Fachler, N. (1985). “Measuring Non-US Equity Portfolio Performance.” Journal of Portfolio Management, 11(3).
13. Brinson, G. P., Hood, L. R., & Beebower, G. L. (1986). “Determinants of Portfolio Performance.” Financial Analysts Journal, 42(4).
14. Cariño, D. (1999). “Combining Attribution Effects Over Time.” Journal of Performance Measurement, 3(4).
15. Menchero, J. (2000). “An Optimized Approach to Linking Attribution Effects Over Time.” Journal of Performance Measurement, 5(1).
16. Artzner, P., Delbaen, F., Eber, J. M., & Heath, D. (1999). “Coherent Measures of Risk.” Mathematical Finance, 9(3).
17. Acerbi, C., & Tasche, D. (2002). “On the Coherence of Expected Shortfall.” Journal of Banking & Finance, 26(7).
18. Hull, J., & White, A. (1998). “Incorporating Volatility Updating into the Historical Simulation Method for Value-at-Risk.” Journal of Risk, 1(1).
19. Kupiec, P. H. (1995). “Techniques for Verifying the Accuracy of Risk Measurement Models.” Journal of Derivatives, 3(2).
20. Christoffersen, P. F. (1998). “Evaluating Interval Forecasts.” International Economic Review, 39(4).
21. Longin, F., & Solnik, B. (2001). “Extreme Correlation of International Equity Markets.” Journal of Finance, 56(2).
22. Ang, A., & Chen, J. (2002). “Asymmetric Correlations of Equity Portfolios.” Journal of Financial Economics, 63(3).
23. Engle, R. F. (1982). “Autoregressive Conditional Heteroscedasticity with Estimates of the Variance of United Kingdom Inflation.” Econometrica, 50(4).
24. Bollerslev, T. (1986). “Generalized Autoregressive Conditional Heteroskedasticity.” Journal of Econometrics, 31(3).
25. RiskMetrics Group (1996). “RiskMetrics — Technical Document.” J.P. Morgan/Reuters, 4th edition.
26. Basel Committee on Banking Supervision (2019). “Minimum Capital Requirements for Market Risk.” BIS Standards（FRTB 最终版）.
27. Federal Reserve Board (2024). “Comprehensive Capital Analysis and Review 2024 Summary Instructions.”
28. Heston, S. L., & Rouwenhorst, K. G. (1994). “Does Industrial Structure Explain the Benefits of International Diversification?” Journal of Financial Economics, 36(1).
29. Connor, G. (1995). “The Three Types of Factor Models: A Comparison of Their Explanatory Power.” Financial Analysts Journal, 51(3).
30. Fama, E. F., & French, K. R. (2015). “A Five-Factor Asset Pricing Model.” Journal of Financial Economics, 116(1).

系列导航

• 上一篇：【量化交易】组合构建与风险预算：BL 模型、风险平价、Kelly
• 下一篇：【量化交易】仓位管理：Kelly 公式、固定分数、波动率目标