通过 Cholesky 前代法计算马氏距离

某风险口在做 copula 标定的合理性检查，需要量化每个历史样本相对中心点的"协方差校正后"距离。该原语就是马氏距离——一维 z-score 的多元推广。对样本 x、中心点 mu 与正定协方差 Sigma：

$$d_M(x, \mu, \Sigma) = \sqrt{(x - \mu)^T \, \Sigma^{-1} \, (x - \mu)}.$$

实现 solution(sample: list[float], mean: list[float], covariance: list[list[float]]) -> float，返回 d_M(sample, mean, covariance)。不要显式求 Sigma 的逆；做 Sigma = L * L^T 的下三角 Cholesky 分解后，把"逆二次型"化为"前代解 y 后取范数"——这是数值稳定的标准配方。

算法（不用 NumPy；纯 Python 嵌套循环；D <= 6 上限处仅 O(D^3)）：

1. 若 D == 0，返回 0.0（零维向量距离为 0）。
2. Cholesky 分解。构造下三角 L，使 L * L^T == covariance。标准递推：

• 对 i in 0..D-1：
• 对 j in 0..i-1：L[i][j] = (covariance[i][j] - sum_{k<j} L[i][k] * L[j][k]) / L[j][j]。
• diag = covariance[i][i] - sum_{k<i} L[i][k]^2。
• 若 diag <= 0，返回 float('nan')——Sigma 不严格正定。
• L[i][i] = sqrt(diag)。

1. 中心化残差 delta = [sample[i] - mean[i] for i in range(D)]。
2. 在下三角 L 上前代解 L * y = delta（top-down）：

• 对 i in 0..D-1：y[i] = (delta[i] - sum_{k=0..i-1} L[i][k] * y[k]) / L[i][i]。

1. 平方距离 m_sq = sum_i y[i]^2。把 m_sq 截断到 0（吸收任何微小负的 cancellation roundoff）。
2. 返回 sqrt(m_sq)。

化简之所以成立：Sigma = L * L^T 蕴含 Sigma^{-1} = L^{-T} * L^{-1}，从而

$$m_{sq} = (x - \mu)^T \, \Sigma^{-1} \, (x - \mu) = \|L^{-1} (x - \mu)\|^2 = \|y\|^2,$$

其中 y 是 L * y = delta 的唯一解。前代是 O(D^2)；不要显式求逆。

例

solution([3.0, 2.0], [1.0, 1.0], [[4.0, 1.0], [1.0, 2.0]]) 返回 sqrt(8/7) ≈ 1.0690449676496976。逐步：delta = [2.0, 1.0]。Cholesky：L[0][0] = sqrt(4) = 2；L[1][0] = covariance[1][0] / L[0][0] = 1 / 2 = 0.5；L[1][1] = sqrt(2 - 0.25) = sqrt(1.75)。前代解 L * y = delta：y[0] = 2 / 2 = 1；y[1] = (1 - 0.5 * 1) / sqrt(1.75) = 0.5 / sqrt(1.75)。平方：m_sq = 1 + 0.25 / 1.75 = 1 + 1/7 = 8/7。距离：sqrt(8/7)。交叉验证：Sigma^{-1} = (1 / 7) * [[2, -1], [-1, 4]]，delta^T * Sigma^{-1} * delta = (1 / 7) * (2 * 4 - 2 * 1 - 1 * 2 + 1 * 4) = (1 / 7) * 8 = 8 / 7。

三个常见坑

Sigma vs Sigma 逆。 马氏距离用的是协方差的逆。把 Sigma 直接代入二次型算 delta^T * Sigma * delta 是另一个标量（且不是真正意义上的"距离"）。对角算例 delta = [2, 3]、Sigma = [[4, 0], [0, 1]]，正解 sqrt(4/4 + 9/1) = sqrt(10)；bug 给 sqrt(16 + 9) = 5。两个都是正的有限浮点，单看一例容易混淆，因此测试集包含一个区分对，使两种输出至少差 1.5 倍以上。

平方 vs 开方。 契约返回的是距离，不是平方距离。平方距离是中间量，也是多元正态 log-PDF 期望的值（-0.5 * m_sq 是其二次型项）；距离带量纲（协方差校正后的标准差），是离群点 dashboard 阈值的入参。漏开方在 m_sq in {0, 1} 时巧合相等，其余都不一样。隐藏测试集包含 m_sq = 10 的算例，让该 bug 一击即穿。

前代 vs 回代。 L 是下三角。解 L * y = delta 必须 top-down：先 y[0]，再用 y[0] 算 y[1]，再用 y[0]、y[1] 算 y[2]，以此类推。若构造上三角因子 U = L^T 然后做 U * y = delta 的后代（先 y[D-1]，自下而上），得到的 y 的平方范数等于 delta^T * (U * U^T)^{-1} * delta，不等于 delta^T * Sigma^{-1} * delta。内层求和的 k 取 0..i-1（不含 i）；写成含 i 会把 L[i][i] 这次除丢，得偏差解。

Sentinel

D == 0：返回 0.0（零维，无距离可言）。Cholesky 非 PD 检测：当标准递推中 diag = covariance[i][i] - sum_{k<i} L[i][k]^2 <= 0（对某个 i），说明 Sigma 不严格正定（或在数值精度内 rank-deficient）；返回 float('nan')。契约把 PD 校验交给调用方；参考实现的 NaN 只是一个防御性 sentinel，恰好捕住"主元为 0 或为负"这种最显眼的失败模式，不抛异常。Roundoff 截断：近奇异但 PD 的 Sigma 可能因 cancellation 让 m_sq 变成微小负数；参考实现在最终 sqrt 前对 m_sq 做 0 截断。rho = 0.99 的算例就在测这个；参考实现对所给 fixture 都精确。

边界

• D == 1、Sigma = [[v]]：L = [[sqrt(v)]]，y = [delta[0] / sqrt(v)]，距离 |delta[0]| / sqrt(v)——一维 z-score 的绝对值。
• Identity Sigma：L = I，y = delta，距离 ||delta||_2——精确还原欧氏距离。
• sample == mean：delta = 0，y = 0，距离 0——与 Sigma 无关。
• 对角 Sigma：L = diag(sqrt(sigma_i^2))，距离 sqrt(sum_i delta[i]^2 / sigma_i^2)——对角马氏即 sigma 加权欧氏。

实践背景

马氏距离是一维 z-score 的多元推广：问"样本相对中心点偏离了多少个协方差校正后的标准差"。在 quant desk 上有三个具体用法：（1）离群点检测——把相对集群中心点马氏距离过大的历史样本标为异常（自由度 D 的卡方分布给出多元正态假设下"预期"马氏距离的分位阈值）。（2）copula 标定诊断——若样本来自相关矩阵给定的高斯 copula，平方马氏距离服从卡方分布；拟合优度检验就归约为"经验 m_sq 分布 vs 卡方分布"的比较。（3）多元正态 log-PDF——其二次型项就是 -0.5 * m_sq，且 Cholesky 前代配方"一次分解给两样东西"：m_sq 用作指数项，det(Sigma) = prod_i L[i][i]^2 用作归一化常数。本题是"诊断 / 距离"原语，与采样姊妹题 coding-gaussian-copula-correlated-sample-via-cholesky 区分明显——后者用 L 把相关性"加上去"，本题用 L 把相关性"扣掉"。

实现细节由 stubs/stub.py 提供。