13  多重检验(Multiple Testing)

在前面的章节中，我们主要关注估计和预测。本章我们将转向假设检验（hypothesis testing），这是统计推断的核心内容。

13.1 学习闭环

先修与可观察目标

完成本章后，学生应能：

  1. 给定真值与拒绝集合计算 \(V,R\)、FWER、FDR 与功效，并准确区分概率、期望与单次实现。
  2. 手工执行 Bonferroni、Holm、BH 与 BY，并说明各自控制目标和依赖条件。
  3. 预注册检验族、效应估计与信息时点，记录计划/实际检验数及所有跳过理由。
  4. 为面板或重叠时序选择公司/年份双聚类、HAC 或时间块重抽样，使原始 p 值先对依赖稳健。
  5. 完成 M13 的冻结检验清单、稳健标准误、多重校正和独立复现合同，不把显著性等同于可交易性。

入口检查（5 分钟）

20 个独立真零假设各按 \(\alpha=0.05\) 检验，期望假阳性数是多少？“至少一次假阳性”的概率是否也等于 1？

展开入口评分、补救与异形复测

答案为 1 与否；后者为 \(1-0.95^{20}\approx0.642\)。若混淆期望计数与事件概率，先复习 式 13.1，再用 \(m=2\) 重测。

无提示检索、渐隐链与迁移

闭书写出 BH 的排序、阈值、最大满足索引和回填拒绝四步；随后仅给 p 值 \((0.006,0.021,0.049,0.20)\) 与 \(q=0.05\) 手算。渐隐第三轮改用 BY 并自行计算 \(c(4)\)。陌生迁移任务是为 30 个城市、12 个窗口的营销实验制定依赖稳健检验族，解释城市共同冲击和时间重叠如何改变原始 p 值与校正方法。

目标—评价映射

目标 评价证据 达标标准
1—2 练习 1—5、10、13 手算、拒绝集合与条件陈述无误
3 练习 6—8 计划数、实际数、跳过原因齐全
4 练习 7—8 每项均记录 estimand、\(n\)、依赖与 SE
5 M13、练习 6—8 冻结清单与复现合同完整

核心非项目练习共 66 分钟；M13 项目练习 55 分钟单列，不重复计入核心练习分钟。

多重检验问题

假设你要测试 1000 个量化因子是否能预测股票收益率。如果每个检验使用显著性水平\(\alpha = 0.05\)，即使所有因子实际上都没有预测能力，我们仍然期望约\(1000 \times 0.05 = 50\) 个因子被错误地认定为有效。

这就是多重检验问题（multiple testing problem）：当我们同时进行大量假设检验时，犯第一类错误的概率会显著增加。

13.2 假设检验回顾

13.2.1 基本概念

假设检验通常包含四个步骤：

  1. 设定假设：

    • 零假设\(H_0\)（默认状态）
    • 备择假设 \(H_a\)（与零假设对照的事前研究方向）

  2. 构造检验统计量：总结反对 \(H_0\) 的证据强度

  3. 计算 p 值：在 \(H_0\) 成立的条件下，获得至少与当前观察一样极端的结果的概率

  4. 做出决策：基于p 值决定是否拒绝\(H_0\)

统计检验只作“拒绝 \(H_0\)”或“不拒绝 \(H_0\)”的决定；前者表示当前设计下存在反对 \(H_0\) 的证据，后者表示证据不足，二者都不给任何命题判定真值。

13.2.2 两类错误

\(H_0\) 为真 \(H_0\) 为假
拒绝 \(H_0\) 第一类错误（Type I Error） 正确决策
不拒绝\(H_0\) 正确决策 第二类错误（Type II Error）
  • 第一类错误率：\(\alpha = \Pr(\text{拒绝 } H_0 | H_0 \text{ 为真})\)
  • 第二类错误率：\(\beta = \Pr(\text{不拒绝} H_0 | H_0 \text{ 为假})\)
  • 检验的功效（power）：\(1 - \beta\)

13.3 多重检验的问题

13.3.1 问题阐述

假设我们要同时检验\(m\) 个零假设 \(H_{01}, H_{02}, \ldots, H_{0m}\)。

令\(V\) 为错误拒绝的零假设的数量（假阳性，false positives），\(R\) 为总的拒绝数量。

若 \(m_0\) 个真零假设的 p 值相互独立，且每个检验都有精确的第一类错误率 \(\alpha\)，则

\[ \Pr(\text{至少犯一次第一类错误}) = 1 - (1 - \alpha)^{m_0}. \tag{13.1}\]

在“全部零假设都为真”的例子中 \(m_0=m\)。不独立时上式不一定成立；不需要独立性的通用上界是并集界 \(\mathrm{FWER}\le m_0\alpha\le m\alpha\)。

直观例子

假设你测试 100 个相互独立且全部为真的零假设，每个使用 \(\alpha = 0.05\)。

  • 单个检验的犯第一类错误概率：5%
  • 至少一次犯错的概率：\(1 - (1 - 0.05)^{100} \approx 99.4\%\)

这意味着几乎肯定会有假阳性发现！

13.3.2 中国案例：股票因子检验

在量化金融中，研究者经常测试大量潜在的交易因子。让我们模拟这个问题。

在量化投资的因子挖掘（Alpha Research）过程中，反复筛选海量股票基本面和量价特征会提高偶然出现亮眼历史回测的机会。一次研究可能找到若干候选，也可能一个都没有；关键问题是这些结果能否在预注册的独立样本中复现。 下面这段简短的 Python 破局代码，向你生动地展示了这种令人绝望的幻觉是如何产生的。我们设定了一个极端悲观且完全由噪音主导的虚拟市场环境：所有的 1000 个所谓“交易因子”（无论是市盈率、波动率还是某些极其复杂的非线性指标），在真实的未来都不具有任何一丁点的预测能力（即所有的零假设完全为真）。模型内部产生的检验统计量全都是毫无意义的标准正态分布随机数。 当 1000 个真零假设各按 \(\alpha=0.05\) 检验时，假阳性数的期望是 50，但单次模拟的实现值会波动。把这些偶然拒绝当作候选会提高实盘失效风险；它本身既不决定具体策略的损益，也不替代交易成本和独立样本评价。

表 13.1: 多重检验问题演示：股票因子测试 
import numpy as np  # 数值计算库
import pandas as pd  # 数据分析库
import matplotlib.pyplot as plt  # 绘图库
plt.rcParams['font.sans-serif'] = ['Source Han Serif SC', 'Arial Unicode MS']  # 设置中文字体
plt.rcParams['axes.unicode_minus'] = False  # 解决负号显示问题
from scipy import stats  # 科学计算统计模块

# 设置随机种子
np.random.seed(42)  # 固定随机种子以保证结果可复现

# 模拟测试 1000 个股票因子
test_count = 1000  # 检验数量
significance_level = 0.05  # 显著性水平

# 情况1：所有零假设都为真（没有真实因子）
# 生成检验统计量（标准正态分布）
null_test_statistics = np.random.normal(0, 1, test_count)  # 从标准正态分布抽取1000个检验统计量
null_p_values = 2 * (1 - stats.norm.cdf(np.abs(null_test_statistics)))  # 计算双侧检验的p值

# 统计显著的结果
significant_null_count = np.sum(null_p_values < significance_level)  # 统计p值低于阈值的假阳性数

print('多重检验问题演示')  # 输出标题
print('='*60)  # 输出分隔线
print(f'总检验数量: {test_count}')  # 输出总检验数
print(f'显著性水平: {significance_level}')  # 输出显著性水平
print(f'预期的假阳性数量: {test_count * significance_level:.1f}')  # 输出理论预期假阳性数
print(f'实际观察到的假阳性数量: {significant_null_count}')  # 输出实际观察到的假阳性数
print(f'假阳性率: {significant_null_count / test_count:.2%}')  # 输出假阳性占总检验数的比率
print(f'至少一次犯错的概率: {1 - (1-significance_level)**test_count:.2%}')  # 输出FWER：至少犯一次第一类错误的概率
多重检验问题演示
============================================================
总检验数量: 1000
显著性水平: 0.05
预期的假阳性数量: 50.0
实际观察到的假阳性数量: 43
假阳性率: 4.30%
至少一次犯错的概率: 100.00%

重复模拟在下一块估计长期错误率；它复用上块已经定义的检验数与显著性水平。

# 用重复实验估计长期错误率；单次实验中的 0/1 指示量不是“错误率”本身
simulation_repetitions = 2000
simulation_rng = np.random.default_rng(2025)
minimum_p_values = np.empty(simulation_repetitions)
false_discovery_proportions = np.empty(simulation_repetitions)
for repetition in range(simulation_repetitions):
    simulated_p_values = simulation_rng.uniform(size=test_count)  # 全部零假设为真
    rejected = simulated_p_values < significance_level
    minimum_p_values[repetition] = simulated_p_values.min()
    false_discovery_proportions[repetition] = 1.0 if rejected.any() else 0.0

empirical_fwer = np.mean(minimum_p_values < significance_level)
empirical_fdr = np.mean(false_discovery_proportions)
print(f'{simulation_repetitions} 次重复模拟的经验 FWER: {empirical_fwer:.3f}')
print(f'全零假设下的经验 FDR: {empirical_fdr:.3f}')
2000 次重复模拟的经验 FWER: 1.000
全零假设下的经验 FDR: 1.000

单次运行中的假阳性个数会随随机种子变化；其期望为 \(m\alpha=50\)。重复模拟估计的是长期 FWER，结果应接近理论值 \(1-(1-\alpha)^m\)。在全零假设下，只要发生拒绝就有 \(V/R=1\)，因此经验 FDR 与经验 FWER 相等；这并不表示任意混合真值情形下二者都相等。

13.3.3 错误度量

在多重检验中，我们需要定义新的错误度量：

符号 描述 频期
\(m\) 总检验数量 固定
\(m_0\) 真实零假设的数量 未知
\(V\) 假阳性数量（Type I error） 随机
\(S\) 真阳性数量（正确拒绝） 随机
\(R\) 总拒绝数量= \(V + S\) 随机
\(T\) 假阴性数量（Type II error） 随机
\(U\) 真阴性数量（正确不拒绝） 随机

13.3.4 FWER 与FDR

13.3.4.1 族错误率（Family-Wise Error Rate, FWER）

\[ \text{FWER} = \Pr(V \geq 1) \tag{13.2}\]

这是至少犯一次第一类错误的概率。FWER 是非常保守的标准。

13.3.4.2 错误发现率（False Discovery Rate, FDR）

\[ \text{FDR} = E\left[ \frac{V}{R} \right] = E\left[ \frac{V}{V + S} \right] \tag{13.3}\]

当\(R = 0\) 时，约定 \(V/R = 0\)。

FWER vs FDR

FWER（族错误率）：

  • 控制任何假阳性的概率
  • 非常保守，适合高风险应用（如重大投资决策、高频策略上线）
  • Bonferroni 校正控制 FWER

FDR（错误发现率）：

  • 控制假阳性在所有拒绝中的比例
  • 较宽松，允许一些假阳性
  • 适合量化策略回测、大规模筛选
  • Benjamini-Hochberg 方法控制 FDR

选择建议：

  • 需要严格控制错误：使用 FWER
  • 允许一些错误以获得更多发现：使用FDR

13.4 Bonferroni 校正

最简单且最保守的多重检验校正方法是 Bonferroni 校正。

13.4.1 方法原理

为了控制 FWER 在水平\(\alpha\)，我们使用调整后的显著性水平：

\[ \alpha_{\text{adjusted}} = \frac{\alpha}{m} \tag{13.4}\]

或者等价地，将 p 值乘以\(m\)。

\[ p_{\text{adjusted}} = \min(m \times p_{\text{original}}, 1) \tag{13.5}\]

13.4.2 中国案例：基于多重检验的强势股筛选(Bonferroni 校正)

下面用 Bonferroni 演示“检验数量增加会抬高至少一次误拒风险”。estimand 明确定义为 2020 年个股绝对日均收益率；它没有减市场、行业或无风险基准，因此不称为相对或超额收益。每只股票先用 HAC 标准误处理日序列相关和异方差，再对这些单项 p 值实施 Bonferroni。 右图把每个检验的阈值改为 \(0.05/m\)，所以拒绝数通常下降。真实市场样本中不知道哪些零假设为真，点的颜色不能称作“真 Alpha”或“假阳性”；Bonferroni 的功效代价必须在预先注入真值的重复模拟中评估。

import pandas as pd  # 数据分析库
import numpy as np  # 数值计算库
import matplotlib.pyplot as plt  # 绘图库
from scipy import stats  # 科学计算统计模块
import statsmodels.api as sm  # 计算个股均值的HAC标准误
import os  # 文件系统操作

np.random.seed(123)  # 设置随机种子以保证抽样结果可复现

from pathlib import Path  # 使用跨平台路径对象解析显式数据根
DATA_DIR = Path(os.environ['BOOK_DATA_DIR']).expanduser().resolve()  # 从必需环境变量取得数据根
if not DATA_DIR.is_dir():  # 在读取前验证数据根
    raise FileNotFoundError(f'BOOK_DATA_DIR 不存在或不是目录: {DATA_DIR}')  # 失败即停止，不回退固定路径
path_price = DATA_DIR / 'stock' / 'stock_price_post_adjusted.h5'  # 构建后复权股价文件路径
stock_price_history = pd.read_hdf(  # 在 HDF 存储层只读取 2020 年真实行情
    path_price,  # 指定后复权股价文件
    where="date>='2020-01-01' & date<'2021-01-01'",  # 固定案例估计年份并下推日期过滤
    columns=['date', 'order_book_id', 'close'],  # 只读取 HAC 均值检验所需字段
).reset_index()  # 恢复日期与证券索引列

# HDF 查询已固定为 2020 年；保留独立副本供收益构造
daily_stock_data = stock_price_history.copy()  # 使用存储层筛选后的2020年全年数据

# 计算每只股票的日收益率
daily_stock_data = daily_stock_data.sort_values(['order_book_id', 'date'])  # 按股票代码和日期排序
daily_stock_data['ret'] = daily_stock_data.groupby('order_book_id')['close'].pct_change()  # 按股票分组计算日收益率

# 获取至少 100天交易数据的有效股票
valid_stock_codes = daily_stock_data.groupby('order_book_id')['ret'].count()  # 统计每只股票的有效交易日数
valid_stock_codes = valid_stock_codes[valid_stock_codes >= 100].index.tolist()  # 筛选交易日不少于100天的股票

target_stock_count = 500  # 设定抽样目标数量为500只
if len(valid_stock_codes) < target_stock_count:
    raise ValueError({'status': 'stopped', 'reason': 'insufficient_planned_stocks', 'available': len(valid_stock_codes), 'required': target_stock_count})
sampled_stock_codes = np.random.choice(valid_stock_codes, target_stock_count, replace=False)  # 无放回随机抽样500只股票

对每只抽样股票进行单样本 t 检验，检验其2020年日均收益率是否显著大于零。

# 进行HAC稳健均值检验(H0: mean <= 0, H1: mean > 0)
calculated_p_values = []  # 初始化p值存储列表
strong_stock_records = []  # 保留效应、SE、带宽和样本量

for stock_code in sampled_stock_codes:  # 遍历每只抽样股票
    stock_returns = daily_stock_data.loc[daily_stock_data['order_book_id'] == stock_code, 'ret'].dropna()  # 获取该股票的日收益率序列
    hac_lags = max(1, int(np.sqrt(len(stock_returns))))  # 事前经验带宽规则
    mean_model = sm.OLS(stock_returns.to_numpy(), np.ones((len(stock_returns), 1))).fit(cov_type='HAC', cov_kwds={'maxlags': hac_lags})  # 稳健均值模型
    t_statistic = mean_model.params[0] / mean_model.bse[0]  # HAC t统计量
    p_value = stats.norm.sf(t_statistic)  # 正向单尾渐近p值
    calculated_p_values.append(p_value)  # 存储该股票的p值
    strong_stock_records.append({'order_book_id': stock_code, 'estimand': 'absolute mean daily return', 'estimate': mean_model.params[0], 'se': mean_model.bse[0], 'hac_lags': hac_lags, 'n': len(stock_returns), 'p_value': p_value})  # 保存单项推断证据

calculated_p_values = np.array(calculated_p_values)  # 将p值列表转换为numpy数组
# 颜色只标出原始 p 值较小的观测，不把同一 p 值反过来当作“真实效应”标签
small_p_highlight = calculated_p_values < 0.01

# 未校正的结果
significance_level = 0.05  # 设定显著性水平为5%
uncorrected_significant_flags = calculated_p_values < significance_level  # 未校正下的显著标记

# Bonferroni 校正
bonferroni_significance_level = significance_level / target_stock_count  # Bonferroni校正后的阈值 = α/m
bonferroni_significant_flags = calculated_p_values < bonferroni_significance_level  # Bonferroni校正后的显著标记
strong_stock_table = pd.DataFrame(strong_stock_records)  # 形成planned family完整单项证据
strong_stock_table['bonferroni_adjusted_p'] = np.minimum(strong_stock_table['p_value'] * target_stock_count, 1.0)  # 任意依赖FWER校正
print(strong_stock_table.to_string(index=False))  # 输出基准、estimand、HAC SE、带宽与校正p值
order_book_id                   estimand      estimate       se  hac_lags   n  p_value  bonferroni_adjusted_p
  000401.XSHE absolute mean daily return -3.591155e-04 0.001488        15 242 0.595372               1.000000
  300128.XSHE absolute mean daily return  1.203349e-04 0.002125        15 242 0.477425               1.000000
  600713.XSHG absolute mean daily return  4.033885e-05 0.000835        15 242 0.480741               1.000000
  600963.XSHG absolute mean daily return  1.080124e-03 0.001352        15 242 0.212102               1.000000
  600681.XSHG absolute mean daily return -7.974885e-04 0.000868        15 242 0.820875               1.000000
  002542.XSHE absolute mean daily return -6.461254e-04 0.001071        15 242 0.726800               1.000000
  600114.XSHG absolute mean daily return  1.584375e-04 0.002408        15 242 0.473768               1.000000
  300648.XSHE absolute mean daily return  5.024693e-03 0.002737        15 242 0.033166               1.000000
  300010.XSHE absolute mean daily return -5.262363e-04 0.002651        15 242 0.578669               1.000000
  601229.XSHG absolute mean daily return -5.628369e-04 0.000695        15 242 0.790869               1.000000
  600210.XSHG absolute mean daily return  1.029225e-03 0.001633        15 242 0.264203               1.000000
  300409.XSHE absolute mean daily return  1.006353e-03 0.001595        15 242 0.264052               1.000000
  300627.XSHE absolute mean daily return  2.141917e-03 0.001790        15 242 0.115781               1.000000
  300328.XSHE absolute mean daily return  3.823571e-06 0.001949        15 242 0.499218               1.000000
  600738.XSHG absolute mean daily return  1.445969e-03 0.001945        15 242 0.228625               1.000000
  300191.XSHE absolute mean daily return -6.784942e-04 0.001490        15 242 0.675523               1.000000
  601588.XSHG absolute mean daily return -1.210050e-03 0.001061        15 242 0.872923               1.000000
  300399.XSHE absolute mean daily return  1.483913e-03 0.002416        15 242 0.269572               1.000000
  603607.XSHG absolute mean daily return  2.729615e-05 0.001288        15 242 0.491549               1.000000
  600744.XSHG absolute mean daily return  1.082832e-03 0.001339        15 242 0.209427               1.000000
  000761.XSHE absolute mean daily return -7.916615e-04 0.000889        15 242 0.813354               1.000000
  600704.XSHG absolute mean daily return -4.014900e-04 0.000829        15 242 0.686010               1.000000
  600097.XSHG absolute mean daily return  4.161159e-04 0.001357        15 242 0.379534               1.000000
  002090.XSHE absolute mean daily return -7.688312e-04 0.001392        15 242 0.709676               1.000000
  002778.XSHE absolute mean daily return  8.003395e-04 0.001690        15 242 0.317865               1.000000
  002986.XSHE absolute mean daily return -2.830167e-03 0.001532        12 144 0.967640               1.000000
  600463.XSHG absolute mean daily return -3.331151e-04 0.001420        15 242 0.592719               1.000000
  603768.XSHG absolute mean daily return  1.841817e-03 0.001839        15 242 0.158305               1.000000
  603335.XSHG absolute mean daily return  2.826230e-05 0.001193        15 242 0.490554               1.000000
  300379.XSHE absolute mean daily return  5.721544e-04 0.001822        15 242 0.376733               1.000000
  000981.XSHE absolute mean daily return  6.711423e-04 0.002015        15 242 0.369515               1.000000
  300757.XSHE absolute mean daily return  1.171916e-03 0.002036        15 242 0.282475               1.000000
  603110.XSHG absolute mean daily return  2.743017e-03 0.001959        15 242 0.080724               1.000000
  600956.XSHG absolute mean daily return  7.383563e-03 0.007163        11 127 0.151318               1.000000
  000156.XSHE absolute mean daily return -3.110583e-04 0.001419        15 242 0.586758               1.000000
  601615.XSHG absolute mean daily return  2.171258e-03 0.001632        15 242 0.091649               1.000000
  300713.XSHE absolute mean daily return  2.597997e-03 0.003293        15 242 0.215077               1.000000
  000417.XSHE absolute mean daily return  2.770345e-04 0.001500        15 242 0.426718               1.000000
  601177.XSHG absolute mean daily return -4.272292e-04 0.001040        15 242 0.659348               1.000000
  600580.XSHG absolute mean daily return  1.597313e-03 0.001579        15 242 0.155933               1.000000
  300106.XSHE absolute mean daily return  3.976395e-03 0.003282        15 242 0.112829               1.000000
  002342.XSHE absolute mean daily return  1.600270e-03 0.001457        15 242 0.135977               1.000000
  600566.XSHG absolute mean daily return -3.210113e-04 0.001156        15 242 0.609416               1.000000
  002017.XSHE absolute mean daily return -1.805874e-03 0.001425        15 242 0.897440               1.000000
  300342.XSHE absolute mean daily return  2.063899e-03 0.003438        15 242 0.274132               1.000000
  002405.XSHE absolute mean daily return -1.369082e-04 0.001876        15 242 0.529087               1.000000
  300330.XSHE absolute mean daily return  4.542328e-04 0.001573        15 242 0.386412               1.000000
  600615.XSHG absolute mean daily return -1.840772e-03 0.001535        15 242 0.884805               1.000000
  601598.XSHG absolute mean daily return  5.293066e-04 0.001382        15 242 0.350891               1.000000
  002980.XSHE absolute mean daily return  5.320838e-03 0.005142        13 175 0.150380               1.000000
  603656.XSHG absolute mean daily return  1.793979e-04 0.000983        15 242 0.427600               1.000000
  002338.XSHE absolute mean daily return  2.234847e-03 0.002847        15 242 0.216227               1.000000
  603199.XSHG absolute mean daily return -5.484987e-04 0.000948        15 242 0.718569               1.000000
  300413.XSHE absolute mean daily return  3.167188e-03 0.001774        15 242 0.037131               1.000000
  002133.XSHE absolute mean daily return -8.227104e-05 0.000852        15 242 0.538463               1.000000
  300278.XSHE absolute mean daily return  1.415982e-03 0.003097        15 242 0.323770               1.000000
  600790.XSHG absolute mean daily return  2.768308e-04 0.000845        15 242 0.371640               1.000000
  000060.XSHE absolute mean daily return  7.630235e-04 0.001397        15 242 0.292421               1.000000
  603606.XSHG absolute mean daily return  3.822702e-03 0.001980        15 242 0.026770               1.000000
  600635.XSHG absolute mean daily return -4.715427e-04 0.001201        15 242 0.652664               1.000000
  002927.XSHE absolute mean daily return -6.157136e-04 0.001352        15 242 0.675549               1.000000
  300358.XSHE absolute mean daily return  2.873515e-03 0.002437        15 242 0.119201               1.000000
  300160.XSHE absolute mean daily return  4.709579e-03 0.004593        15 242 0.152603               1.000000
  002071.XSHE absolute mean daily return -4.673999e-03 0.003277        15 242 0.923112               1.000000
  603131.XSHG absolute mean daily return  3.902909e-03 0.002744        15 242 0.077492               1.000000
  603566.XSHG absolute mean daily return  9.275016e-04 0.001750        15 242 0.298029               1.000000
  600613.XSHG absolute mean daily return -1.275584e-03 0.001561        15 242 0.793021               1.000000
  603027.XSHG absolute mean daily return  4.236308e-03 0.001839        15 242 0.010612               1.000000
  603890.XSHG absolute mean daily return  2.337644e-03 0.002230        15 242 0.147234               1.000000
  002137.XSHE absolute mean daily return -2.644342e-04 0.001541        15 242 0.568126               1.000000
  300588.XSHE absolute mean daily return  6.382411e-04 0.002155        15 242 0.383533               1.000000
  600309.XSHG absolute mean daily return  2.469454e-03 0.001588        15 242 0.059936               1.000000
  002685.XSHE absolute mean daily return -1.718488e-03 0.002072        15 242 0.796555               1.000000
  688157.XSHG absolute mean daily return  7.581701e-05 0.002714        11 139 0.488857               1.000000
  601866.XSHG absolute mean daily return  9.532978e-04 0.001801        15 242 0.298250               1.000000
  000563.XSHE absolute mean daily return -3.809174e-04 0.001386        15 242 0.608253               1.000000
  300707.XSHE absolute mean daily return  2.390066e-04 0.002193        15 242 0.456604               1.000000
  603533.XSHG absolute mean daily return  3.717542e-03 0.003524        15 242 0.145711               1.000000
  300019.XSHE absolute mean daily return  3.172761e-03 0.002499        15 242 0.102150               1.000000
  300559.XSHE absolute mean daily return -1.078230e-03 0.001498        15 242 0.764139               1.000000
  300442.XSHE absolute mean daily return  4.055414e-03 0.002154        15 242 0.029897               1.000000
  601127.XSHG absolute mean daily return  2.513832e-03 0.003286        15 242 0.222104               1.000000
  002004.XSHE absolute mean daily return  8.409403e-04 0.001007        15 242 0.201855               1.000000
  002621.XSHE absolute mean daily return -1.160093e-03 0.001477        15 242 0.783942               1.000000
  603390.XSHG absolute mean daily return -1.962459e-03 0.001587        15 242 0.891909               1.000000
  002531.XSHE absolute mean daily return  1.534227e-03 0.001683        15 242 0.180942               1.000000
  603725.XSHG absolute mean daily return  4.354255e-05 0.001347        15 242 0.487105               1.000000
  002318.XSHE absolute mean daily return  1.063753e-03 0.001511        15 242 0.240694               1.000000
  600649.XSHG absolute mean daily return -7.140233e-05 0.001128        15 242 0.525246               1.000000
  600600.XSHG absolute mean daily return  3.145115e-03 0.001680        15 242 0.030635               1.000000
  300375.XSHE absolute mean daily return  2.183475e-03 0.002207        15 242 0.161290               1.000000
  300319.XSHE absolute mean daily return -8.764320e-04 0.001985        15 242 0.670576               1.000000
  000780.XSHE absolute mean daily return  9.644144e-04 0.001538        15 242 0.265340               1.000000
  600348.XSHG absolute mean daily return  5.573836e-04 0.001140        15 242 0.312503               1.000000
  600356.XSHG absolute mean daily return  4.756505e-04 0.001167        15 242 0.341804               1.000000
  000816.XSHE absolute mean daily return  4.627934e-03 0.004112        15 242 0.130179               1.000000
  300301.XSHE absolute mean daily return  2.022099e-03 0.004740        15 242 0.334838               1.000000
  002959.XSHE absolute mean daily return  3.711148e-03 0.002336        15 242 0.056104               1.000000
  300653.XSHE absolute mean daily return  1.158832e-03 0.001658        15 242 0.242353               1.000000
  002016.XSHE absolute mean daily return -8.524508e-04 0.001114        15 242 0.777978               1.000000
  002102.XSHE absolute mean daily return -9.901874e-04 0.000882        15 242 0.869092               1.000000
  688086.XSHG absolute mean daily return -2.638415e-03 0.002197        14 209 0.885089               1.000000
  002975.XSHE absolute mean daily return  7.809489e-03 0.004338        14 224 0.035898               1.000000
  601997.XSHG absolute mean daily return -4.704582e-04 0.000899        15 242 0.699704               1.000000
  600784.XSHG absolute mean daily return  5.825333e-04 0.001482        15 242 0.347147               1.000000
  002190.XSHE absolute mean daily return  1.692540e-03 0.001810        15 242 0.174813               1.000000
  600316.XSHG absolute mean daily return  6.909304e-03 0.003092        15 242 0.012722               1.000000
  002114.XSHE absolute mean daily return  6.798133e-04 0.001417        15 242 0.315751               1.000000
  603356.XSHG absolute mean daily return -3.918206e-04 0.001575        15 242 0.598207               1.000000
  605318.XSHG absolute mean daily return -1.963640e-03 0.003699        10 102 0.702245               1.000000
  002983.XSHE absolute mean daily return  3.848220e-03 0.005566        12 166 0.244680               1.000000
  002320.XSHE absolute mean daily return  7.966614e-04 0.001598        15 242 0.309055               1.000000
  603098.XSHG absolute mean daily return -9.487794e-04 0.001195        15 242 0.786314               1.000000
  002454.XSHE absolute mean daily return  1.127999e-03 0.001293        15 242 0.191422               1.000000
  002614.XSHE absolute mean daily return  1.526526e-03 0.002407        15 242 0.262991               1.000000
  002145.XSHE absolute mean daily return  1.281302e-03 0.001673        15 242 0.221819               1.000000
  603095.XSHG absolute mean daily return -2.463074e-03 0.001681        13 175 0.928601               1.000000
  002293.XSHE absolute mean daily return  1.671599e-03 0.001277        15 242 0.095228               1.000000
  600021.XSHG absolute mean daily return -3.606905e-04 0.000729        15 242 0.689604               1.000000
  002641.XSHE absolute mean daily return  1.848600e-03 0.002174        15 242 0.197613               1.000000
  300583.XSHE absolute mean daily return -1.288367e-03 0.001028        15 242 0.895033               1.000000
  300395.XSHE absolute mean daily return  4.527228e-03 0.002072        15 242 0.014450               1.000000
  603500.XSHG absolute mean daily return  1.929157e-04 0.003044        15 242 0.474733               1.000000
  002129.XSHE absolute mean daily return  3.799591e-03 0.002018        15 242 0.029865               1.000000
  002579.XSHE absolute mean daily return -1.609182e-04 0.001336        15 242 0.547924               1.000000
  605158.XSHG absolute mean daily return -6.329150e-04 0.002692        10 100 0.592939               1.000000
  002442.XSHE absolute mean daily return -6.788720e-04 0.001641        15 242 0.660427               1.000000
  000909.XSHE absolute mean daily return  5.259317e-04 0.001294        15 242 0.342223               1.000000
  600903.XSHG absolute mean daily return -3.827172e-04 0.001277        15 242 0.617782               1.000000
  603895.XSHG absolute mean daily return -1.052985e-03 0.001208        15 242 0.808381               1.000000
  603897.XSHG absolute mean daily return  2.807087e-05 0.001260        15 242 0.491111               1.000000
  002932.XSHE absolute mean daily return  3.108222e-03 0.003041        15 242 0.153338               1.000000
  300545.XSHE absolute mean daily return  6.137387e-04 0.001940        15 242 0.375892               1.000000
  000629.XSHE absolute mean daily return -9.900697e-04 0.001475        15 242 0.749005               1.000000
  002155.XSHE absolute mean daily return  6.518047e-04 0.001312        15 242 0.309672               1.000000
  601319.XSHG absolute mean daily return -3.692296e-04 0.001219        15 242 0.619032               1.000000
  300637.XSHE absolute mean daily return -1.065281e-03 0.001815        15 242 0.721328               1.000000
  300179.XSHE absolute mean daily return  1.124476e-03 0.002292        15 242 0.311844               1.000000
  600766.XSHG absolute mean daily return -9.392741e-04 0.001414        15 242 0.746797               1.000000
  300097.XSHE absolute mean daily return  1.128808e-03 0.002956        15 242 0.351300               1.000000
  600486.XSHG absolute mean daily return  3.039957e-03 0.001438        15 242 0.017282               1.000000
  603218.XSHG absolute mean daily return  3.277318e-03 0.001438        15 242 0.011326               1.000000
  603214.XSHG absolute mean daily return -8.099453e-04 0.001604        15 242 0.693169               1.000000
  002850.XSHE absolute mean daily return  3.801103e-03 0.002526        15 242 0.066190               1.000000
  600163.XSHG absolute mean daily return  8.344094e-04 0.001206        15 242 0.244497               1.000000
  600698.XSHG absolute mean daily return -5.263785e-04 0.001321        15 242 0.654831               1.000000
  600380.XSHG absolute mean daily return  1.631498e-03 0.001753        15 242 0.175993               1.000000
  300609.XSHE absolute mean daily return -2.552256e-03 0.002146        15 242 0.882825               1.000000
  688033.XSHG absolute mean daily return -2.197210e-03 0.001416        15 242 0.939594               1.000000
  300065.XSHE absolute mean daily return  1.173725e-03 0.002109        15 242 0.288916               1.000000
  600500.XSHG absolute mean daily return  2.797834e-04 0.001192        15 242 0.407190               1.000000
  002543.XSHE absolute mean daily return -2.984059e-04 0.001022        15 242 0.614804               1.000000
  300650.XSHE absolute mean daily return  3.334738e-04 0.002081        15 242 0.436349               1.000000
  600760.XSHG absolute mean daily return  4.375583e-03 0.003004        15 242 0.072629               1.000000
  603716.XSHG absolute mean daily return  1.080429e-04 0.001331        15 242 0.467657               1.000000
  002838.XSHE absolute mean daily return  4.720755e-03 0.004518        15 242 0.148019               1.000000
  002582.XSHE absolute mean daily return  2.198230e-03 0.001795        15 242 0.110392               1.000000
  002176.XSHE absolute mean daily return  3.713117e-04 0.003160        15 242 0.453226               1.000000
  603308.XSHG absolute mean daily return  3.537728e-03 0.002284        15 242 0.060707               1.000000
  002265.XSHE absolute mean daily return -5.699602e-04 0.001348        15 242 0.663826               1.000000
  300464.XSHE absolute mean daily return  1.697838e-03 0.004025        15 242 0.336585               1.000000
  000617.XSHE absolute mean daily return -7.023418e-04 0.001099        15 242 0.738681               1.000000
  300150.XSHE absolute mean daily return  8.341800e-04 0.001416        15 242 0.277830               1.000000
  600829.XSHG absolute mean daily return  7.443260e-04 0.001676        15 242 0.328488               1.000000
  600416.XSHG absolute mean daily return  5.276542e-03 0.003480        15 242 0.064748               1.000000
  002117.XSHE absolute mean daily return -9.996945e-04 0.001778        15 242 0.713040               1.000000
  002401.XSHE absolute mean daily return  9.736708e-04 0.001590        15 242 0.270128               1.000000
  600699.XSHG absolute mean daily return  2.051635e-03 0.002226        15 242 0.178322               1.000000
  600610.XSHG absolute mean daily return  2.363076e-03 0.002740        15 242 0.194231               1.000000
  002367.XSHE absolute mean daily return  1.516416e-03 0.001453        15 242 0.148363               1.000000
  600847.XSHG absolute mean daily return -2.797209e-04 0.001283        15 242 0.586312               1.000000
  002879.XSHE absolute mean daily return -1.410450e-04 0.001314        15 242 0.542755               1.000000
  603006.XSHG absolute mean daily return  8.176930e-04 0.001610        15 242 0.305774               1.000000
  002306.XSHE absolute mean daily return  4.326268e-04 0.001661        15 242 0.397266               1.000000
  300205.XSHE absolute mean daily return  2.837935e-04 0.001930        15 242 0.441546               1.000000
  600329.XSHG absolute mean daily return  1.272636e-03 0.001154        15 242 0.134982               1.000000
  600978.XSHG absolute mean daily return -4.771053e-03 0.002808        15 242 0.955331               1.000000
  300448.XSHE absolute mean daily return  1.242002e-03 0.002581        15 242 0.315204               1.000000
  603063.XSHG absolute mean daily return  3.436078e-03 0.001964        15 242 0.040138               1.000000
  002611.XSHE absolute mean daily return  6.162600e-04 0.001256        15 242 0.311876               1.000000
  688202.XSHG absolute mean daily return  4.788265e-03 0.001957        15 242 0.007217               1.000000
  600898.XSHG absolute mean daily return -6.061365e-04 0.002406        15 242 0.599459               1.000000
  300701.XSHE absolute mean daily return -9.785880e-05 0.001621        15 242 0.524066               1.000000
  600238.XSHG absolute mean daily return  1.891487e-03 0.002345        15 242 0.209901               1.000000
  600885.XSHG absolute mean daily return  2.384574e-03 0.001643        15 242 0.073358               1.000000
  600536.XSHG absolute mean daily return  9.419405e-04 0.002278        15 242 0.339652               1.000000
  603600.XSHG absolute mean daily return  1.022326e-03 0.002422        15 242 0.336463               1.000000
  300163.XSHE absolute mean daily return  3.315771e-04 0.001971        15 242 0.433204               1.000000
  002177.XSHE absolute mean daily return  4.022795e-04 0.001986        15 242 0.419728               1.000000
  688258.XSHG absolute mean daily return -4.657447e-04 0.002797        15 242 0.566129               1.000000
  601566.XSHG absolute mean daily return  5.081171e-04 0.001324        15 242 0.350557               1.000000
  300431.XSHE absolute mean daily return -1.116508e-02 0.006448        14 204 0.958317               1.000000
  000750.XSHE absolute mean daily return  1.169396e-03 0.001634        15 242 0.237150               1.000000
  300624.XSHE absolute mean daily return  1.998287e-03 0.002448        15 242 0.207149               1.000000
  002411.XSHE absolute mean daily return -3.953076e-03 0.002142        15 242 0.967502               1.000000
  002214.XSHE absolute mean daily return  4.526604e-03 0.002949        15 242 0.062370               1.000000
  002856.XSHE absolute mean daily return -8.364983e-04 0.001649        15 242 0.694034               1.000000
  002106.XSHE absolute mean daily return  1.162778e-03 0.001892        15 242 0.269467               1.000000
  300132.XSHE absolute mean daily return  2.667470e-03 0.002064        15 242 0.098121               1.000000
  688266.XSHG absolute mean daily return  3.281905e-04 0.002975        15 227 0.456073               1.000000
  002053.XSHE absolute mean daily return  3.033558e-04 0.001014        15 242 0.382416               1.000000
  601100.XSHG absolute mean daily return  5.417944e-03 0.001247        15 242 0.000007               0.003489
  300391.XSHE absolute mean daily return  2.978248e-03 0.003001        15 242 0.160515               1.000000
  300599.XSHE absolute mean daily return  1.200916e-03 0.001459        15 242 0.205179               1.000000
  603086.XSHG absolute mean daily return -3.111989e-04 0.001256        15 242 0.597814               1.000000
  600968.XSHG absolute mean daily return -6.926707e-04 0.001023        15 242 0.750755               1.000000
  000882.XSHE absolute mean daily return -6.265352e-04 0.000892        15 242 0.758865               1.000000
  000807.XSHE absolute mean daily return  2.048342e-03 0.002665        15 242 0.221034               1.000000
  002504.XSHE absolute mean daily return -7.999163e-04 0.001434        15 242 0.711514               1.000000
  600800.XSHG absolute mean daily return -4.801354e-04 0.001153        15 242 0.661479               1.000000
  002674.XSHE absolute mean daily return  1.429584e-03 0.002913        15 242 0.311800               1.000000
  603136.XSHG absolute mean daily return  3.014957e-04 0.001642        15 242 0.427142               1.000000
  603195.XSHG absolute mean daily return  4.520932e-03 0.003438        14 223 0.094240               1.000000
  002792.XSHE absolute mean daily return -1.297372e-03 0.001618        15 242 0.788654               1.000000
  002241.XSHE absolute mean daily return  3.206129e-03 0.002350        15 242 0.086196               1.000000
  000591.XSHE absolute mean daily return  3.491999e-03 0.002356        15 242 0.069178               1.000000
  603668.XSHG absolute mean daily return  2.191824e-04 0.001874        15 242 0.453455               1.000000
  300152.XSHE absolute mean daily return -2.458992e-05 0.003093        15 242 0.503172               1.000000
  600889.XSHG absolute mean daily return  1.951629e-04 0.002384        15 242 0.467380               1.000000
  300336.XSHE absolute mean daily return -3.646142e-04 0.002574        15 242 0.556316               1.000000
  002215.XSHE absolute mean daily return -1.928209e-04 0.001049        15 242 0.572935               1.000000
  002780.XSHE absolute mean daily return  7.709633e-04 0.001669        15 242 0.322102               1.000000
  600336.XSHG absolute mean daily return  4.194451e-03 0.002314        15 242 0.034936               1.000000
  300252.XSHE absolute mean daily return  6.591179e-04 0.001649        15 242 0.344706               1.000000
  000430.XSHE absolute mean daily return -7.656229e-05 0.001287        15 242 0.523712               1.000000
  603993.XSHG absolute mean daily return  2.026975e-03 0.002015        15 242 0.157189               1.000000
  601188.XSHG absolute mean daily return -8.912903e-05 0.000893        15 242 0.539730               1.000000
  300233.XSHE absolute mean daily return  4.798053e-04 0.002406        15 242 0.420978               1.000000
  300555.XSHE absolute mean daily return  8.968848e-04 0.002642        15 242 0.367130               1.000000
  600346.XSHG absolute mean daily return  2.735940e-03 0.001691        15 242 0.052887               1.000000
  600063.XSHG absolute mean daily return  9.422309e-05 0.001564        15 242 0.475983               1.000000
  300256.XSHE absolute mean daily return  5.382615e-04 0.001593        15 242 0.367752               1.000000
  688278.XSHG absolute mean daily return  8.662091e-04 0.003434        15 231 0.400421               1.000000
  002253.XSHE absolute mean daily return  8.186059e-05 0.001252        15 242 0.473934               1.000000
  002084.XSHE absolute mean daily return  1.277348e-03 0.002932        15 242 0.331567               1.000000
  300386.XSHE absolute mean daily return  1.790812e-03 0.002135        15 242 0.200829               1.000000
  603927.XSHG absolute mean daily return -1.817749e-03 0.001725        15 242 0.853971               1.000000
  002066.XSHE absolute mean daily return  7.546310e-04 0.001626        15 242 0.321280               1.000000
  002097.XSHE absolute mean daily return  1.644241e-03 0.001565        15 242 0.146700               1.000000
  600859.XSHG absolute mean daily return  4.350482e-03 0.004335        15 242 0.157797               1.000000
  601886.XSHG absolute mean daily return -5.671430e-04 0.000925        15 242 0.730095               1.000000
  600696.XSHG absolute mean daily return -2.616000e-05 0.002064        15 242 0.505056               1.000000
  688981.XSHG absolute mean daily return -2.660059e-03 0.002739        10 114 0.834309               1.000000
  002817.XSHE absolute mean daily return  2.062087e-04 0.001176        15 242 0.430378               1.000000
  603288.XSHG absolute mean daily return  3.619411e-03 0.001428        15 242 0.005642               1.000000
  300474.XSHE absolute mean daily return  1.150029e-03 0.001699        15 242 0.249264               1.000000
  002833.XSHE absolute mean daily return  2.706683e-03 0.001954        15 242 0.083034               1.000000
  002360.XSHE absolute mean daily return  1.955691e-03 0.001805        15 242 0.139363               1.000000
  300217.XSHE absolute mean daily return  2.876101e-03 0.003123        15 242 0.178535               1.000000
  002654.XSHE absolute mean daily return -1.022751e-03 0.001549        15 242 0.745454               1.000000
  600975.XSHG absolute mean daily return  3.516934e-04 0.001856        15 242 0.424848               1.000000
  300345.XSHE absolute mean daily return -2.299184e-03 0.002535        15 242 0.817748               1.000000
  600148.XSHG absolute mean daily return  6.982777e-04 0.001534        15 242 0.324478               1.000000
  300287.XSHE absolute mean daily return  1.082480e-03 0.002195        15 242 0.310930               1.000000
  688108.XSHG absolute mean daily return  3.499419e-05 0.003153        15 242 0.495572               1.000000
  000958.XSHE absolute mean daily return -4.199389e-04 0.001362        15 242 0.621059               1.000000
  000529.XSHE absolute mean daily return  2.454984e-04 0.001355        15 242 0.428128               1.000000
  300480.XSHE absolute mean daily return  6.859365e-04 0.001891        15 242 0.358435               1.000000
  688008.XSHG absolute mean daily return  1.362864e-03 0.002585        15 242 0.299001               1.000000
  603023.XSHG absolute mean daily return  8.754522e-04 0.001842        15 242 0.317317               1.000000
  600802.XSHG absolute mean daily return  8.483297e-04 0.001764        15 242 0.315302               1.000000
  600162.XSHG absolute mean daily return -3.457797e-04 0.001127        15 242 0.620465               1.000000
  600552.XSHG absolute mean daily return  1.034222e-03 0.001692        15 242 0.270511               1.000000
  601965.XSHG absolute mean daily return  2.956346e-03 0.001534        15 242 0.026999               1.000000
  002368.XSHE absolute mean daily return  5.112336e-05 0.001655        15 242 0.487675               1.000000
  300392.XSHE absolute mean daily return  1.892644e-03 0.002014        15 242 0.173725               1.000000
  002481.XSHE absolute mean daily return  2.922433e-03 0.002681        15 242 0.137804               1.000000
  600143.XSHG absolute mean daily return  4.130399e-03 0.002089        15 242 0.023983               1.000000
  002938.XSHE absolute mean daily return  9.378225e-04 0.001967        15 242 0.316729               1.000000
  603885.XSHG absolute mean daily return -8.091253e-04 0.001484        15 242 0.707146               1.000000
  300051.XSHE absolute mean daily return  1.134752e-04 0.003257        15 242 0.486104               1.000000
  600644.XSHG absolute mean daily return  4.830381e-04 0.001174        15 242 0.340345               1.000000
  002992.XSHE absolute mean daily return -1.112005e-03 0.004174        10 102 0.605050               1.000000
  601636.XSHG absolute mean daily return  4.294304e-03 0.002442        15 242 0.039355               1.000000
  300436.XSHE absolute mean daily return  9.607528e-04 0.001466        15 242 0.256127               1.000000
  603778.XSHG absolute mean daily return  2.155736e-04 0.001265        15 242 0.432331               1.000000
  600400.XSHG absolute mean daily return  6.829155e-05 0.001300        15 242 0.479055               1.000000
  002228.XSHE absolute mean daily return  5.991343e-04 0.001711        15 242 0.363131               1.000000
  002096.XSHE absolute mean daily return  1.038736e-03 0.003303        15 242 0.376584               1.000000
  002738.XSHE absolute mean daily return  2.393034e-03 0.001739        15 242 0.084372               1.000000
  300521.XSHE absolute mean daily return  3.495568e-03 0.003782        15 242 0.177684               1.000000
  002637.XSHE absolute mean daily return  1.251101e-03 0.001598        15 242 0.216826               1.000000
  002613.XSHE absolute mean daily return  2.770715e-03 0.003769        15 242 0.231159               1.000000
  000153.XSHE absolute mean daily return  1.718992e-03 0.002173        15 242 0.214413               1.000000
  300625.XSHE absolute mean daily return  3.462851e-04 0.001392        15 242 0.401792               1.000000
  688098.XSHG absolute mean daily return  6.150750e-04 0.002285        15 242 0.393912               1.000000
  688558.XSHG absolute mean daily return -3.722078e-03 0.001987        11 126 0.969513               1.000000
  603798.XSHG absolute mean daily return -1.289860e-03 0.001185        15 242 0.861751               1.000000
  300760.XSHE absolute mean daily return  3.907244e-03 0.001439        15 242 0.003306               1.000000
  002101.XSHE absolute mean daily return  1.841757e-04 0.001450        15 242 0.449471               1.000000
  600039.XSHG absolute mean daily return  1.397592e-03 0.001151        15 242 0.112379               1.000000
  600861.XSHG absolute mean daily return  3.654742e-03 0.002187        15 242 0.047353               1.000000
  688198.XSHG absolute mean daily return  3.308386e-03 0.003096        15 242 0.142599               1.000000
  002424.XSHE absolute mean daily return -1.340510e-04 0.001188        15 242 0.544914               1.000000
  688580.XSHG absolute mean daily return -4.364507e-03 0.003884        10 111 0.869444               1.000000
  605288.XSHG absolute mean daily return -1.145467e-04 0.001623        12 145 0.528126               1.000000
  600077.XSHG absolute mean daily return  4.188779e-04 0.001668        15 242 0.400865               1.000000
  002846.XSHE absolute mean daily return -7.341443e-04 0.002514        15 242 0.614886               1.000000
  601021.XSHG absolute mean daily return  1.238778e-03 0.001418        15 242 0.191160               1.000000
  600578.XSHG absolute mean daily return  1.270439e-04 0.000761        15 242 0.433745               1.000000
  601399.XSHG absolute mean daily return -2.786342e-03 0.002595        11 140 0.858503               1.000000
  002146.XSHE absolute mean daily return -1.277518e-03 0.001220        15 242 0.852495               1.000000
  688566.XSHG absolute mean daily return -3.183040e-03 0.001503        12 155 0.982931               1.000000
  600022.XSHG absolute mean daily return  1.405656e-04 0.000860        15 242 0.435053               1.000000
  600909.XSHG absolute mean daily return  9.733751e-04 0.001751        15 242 0.289140               1.000000
  002692.XSHE absolute mean daily return -1.092079e-03 0.000990        15 242 0.864902               1.000000
  002076.XSHE absolute mean daily return -5.775094e-04 0.003025        15 242 0.575696               1.000000
  600690.XSHG absolute mean daily return  2.010671e-03 0.001796        15 242 0.131517               1.000000
  600108.XSHG absolute mean daily return  1.357518e-03 0.001371        15 242 0.161111               1.000000
  002987.XSHE absolute mean daily return  2.926879e-03 0.003958        12 162 0.229805               1.000000
  600710.XSHG absolute mean daily return  3.687691e-04 0.001442        15 242 0.399079               1.000000
  600777.XSHG absolute mean daily return -9.885851e-04 0.001299        15 242 0.776733               1.000000
  000877.XSHE absolute mean daily return  1.595384e-03 0.002130        15 242 0.226884               1.000000
  000980.XSHE absolute mean daily return -2.843303e-03 0.002083        15 242 0.913925               1.000000
  000715.XSHE absolute mean daily return -3.753050e-04 0.001616        15 242 0.591824               1.000000
  600363.XSHG absolute mean daily return  2.520930e-03 0.002167        15 242 0.122317               1.000000
  600029.XSHG absolute mean daily return -6.231804e-04 0.001234        15 242 0.693250               1.000000
  002596.XSHE absolute mean daily return -4.797195e-04 0.001664        15 242 0.613456               1.000000
  002383.XSHE absolute mean daily return -8.911586e-04 0.001623        15 242 0.708554               1.000000
  601226.XSHG absolute mean daily return  1.413238e-04 0.000934        15 242 0.439842               1.000000
  002312.XSHE absolute mean daily return  1.649870e-03 0.001502        15 242 0.135988               1.000000
  300591.XSHE absolute mean daily return  1.687122e-03 0.002470        15 242 0.247255               1.000000
  300360.XSHE absolute mean daily return  7.651259e-05 0.001589        15 242 0.480800               1.000000
  600624.XSHG absolute mean daily return -6.517459e-04 0.001374        15 242 0.682388               1.000000
  300745.XSHE absolute mean daily return  6.499401e-05 0.001969        15 242 0.486832               1.000000
  300244.XSHE absolute mean daily return  2.179489e-03 0.001863        15 242 0.121017               1.000000
  300036.XSHE absolute mean daily return -3.531074e-05 0.001686        15 242 0.508357               1.000000
  300727.XSHE absolute mean daily return  4.140136e-03 0.003063        15 242 0.088260               1.000000
  300153.XSHE absolute mean daily return  6.545177e-04 0.002181        15 242 0.382026               1.000000
  601158.XSHG absolute mean daily return -1.255491e-04 0.000496        15 242 0.599977               1.000000
  000419.XSHE absolute mean daily return  1.747794e-04 0.001096        15 242 0.436641               1.000000
  002937.XSHE absolute mean daily return -4.585046e-04 0.001247        15 242 0.643422               1.000000
  002677.XSHE absolute mean daily return  1.149360e-03 0.001802        15 242 0.261823               1.000000
  000028.XSHE absolute mean daily return  3.762358e-04 0.001361        15 242 0.391087               1.000000
  300145.XSHE absolute mean daily return -2.630298e-05 0.001618        15 242 0.506486               1.000000
  605222.XSHG absolute mean daily return -1.753750e-03 0.001877        10 103 0.824971               1.000000
  601099.XSHG absolute mean daily return  7.241874e-04 0.001883        15 242 0.350268               1.000000
  600809.XSHG absolute mean daily return  6.413168e-03 0.001513        15 242 0.000011               0.005592
  002563.XSHE absolute mean daily return  5.883318e-04 0.001501        15 242 0.347587               1.000000
  002703.XSHE absolute mean daily return  1.333343e-03 0.001503        15 242 0.187470               1.000000
  600081.XSHG absolute mean daily return  2.166609e-03 0.001769        15 242 0.110331               1.000000
  002574.XSHE absolute mean daily return -2.607159e-04 0.000868        15 242 0.618040               1.000000
  002086.XSHE absolute mean daily return -2.614492e-03 0.001980        15 242 0.906678               1.000000
  300002.XSHE absolute mean daily return  2.656855e-03 0.002332        15 242 0.127262               1.000000
  002019.XSHE absolute mean daily return  1.230905e-03 0.001910        15 242 0.259634               1.000000
  603959.XSHG absolute mean daily return  1.444660e-03 0.002415        15 242 0.274828               1.000000
  600748.XSHG absolute mean daily return -3.753614e-04 0.001498        15 242 0.598924               1.000000
  600461.XSHG absolute mean daily return  7.466729e-04 0.000847        15 242 0.188911               1.000000
  300268.XSHE absolute mean daily return -2.581828e-04 0.001683        15 242 0.560955               1.000000
  000700.XSHE absolute mean daily return  2.244458e-03 0.004765        15 242 0.318799               1.000000
  002824.XSHE absolute mean daily return  2.546527e-03 0.004135        15 242 0.268989               1.000000
  300159.XSHE absolute mean daily return  1.313342e-03 0.002582        15 242 0.305521               1.000000
  300853.XSHE absolute mean daily return  4.523682e-03 0.006077        10 108 0.228313               1.000000
  002108.XSHE absolute mean daily return  1.917566e-03 0.002528        15 242 0.224022               1.000000
  002826.XSHE absolute mean daily return  1.957396e-04 0.001388        15 242 0.443913               1.000000
  600716.XSHG absolute mean daily return  1.131047e-04 0.000997        15 242 0.454834               1.000000
  002565.XSHE absolute mean daily return -1.236032e-03 0.001930        15 242 0.739104               1.000000
  603808.XSHG absolute mean daily return -9.653874e-05 0.001708        15 242 0.522540               1.000000
  603117.XSHG absolute mean daily return -3.821088e-05 0.000950        15 242 0.516046               1.000000
  000902.XSHE absolute mean daily return  3.244692e-03 0.001503        15 242 0.015460               1.000000
  600120.XSHG absolute mean daily return  1.272816e-04 0.001535        15 242 0.466947               1.000000
  600169.XSHG absolute mean daily return -1.562383e-04 0.001116        15 242 0.555680               1.000000
  002012.XSHE absolute mean daily return  1.633520e-04 0.001381        15 242 0.452918               1.000000
  002813.XSHE absolute mean daily return  1.168645e-04 0.001966        15 242 0.476305               1.000000
  601012.XSHG absolute mean daily return  5.735573e-03 0.002303        15 242 0.006380               1.000000
  002721.XSHE absolute mean daily return -2.186018e-03 0.001497        15 242 0.927836               1.000000
  002751.XSHE absolute mean daily return -1.047648e-03 0.001462        15 242 0.763116               1.000000
  002918.XSHE absolute mean daily return  2.559955e-03 0.002242        15 242 0.126760               1.000000
  600592.XSHG absolute mean daily return -1.168587e-03 0.001715        15 242 0.752172               1.000000
  688277.XSHG absolute mean daily return -4.948929e-03 0.004499        11 121 0.864326               1.000000
  300511.XSHE absolute mean daily return  3.346462e-03 0.002094        15 242 0.054989               1.000000
  002668.XSHE absolute mean daily return -1.184056e-03 0.001520        15 242 0.781979               1.000000
  300304.XSHE absolute mean daily return  1.288262e-03 0.002164        15 242 0.275808               1.000000
  002604.XSHE absolute mean daily return -1.555826e-02 0.009208        11 126 0.954447               1.000000
  000032.XSHE absolute mean daily return  1.420467e-03 0.002396        15 242 0.276665               1.000000
  300167.XSHE absolute mean daily return  1.220856e-04 0.001915        15 242 0.474589               1.000000
  300074.XSHE absolute mean daily return  2.960883e-04 0.001890        15 242 0.437760               1.000000
  002470.XSHE absolute mean daily return -2.702591e-03 0.001960        15 242 0.916037               1.000000
  600676.XSHG absolute mean daily return  7.707466e-04 0.001410        15 242 0.292258               1.000000
  002010.XSHE absolute mean daily return -1.435271e-03 0.000957        15 242 0.933220               1.000000
  603181.XSHG absolute mean daily return  2.360237e-03 0.001548        15 242 0.063635               1.000000
  000669.XSHE absolute mean daily return -1.907312e-03 0.002601        15 242 0.768272               1.000000
  603986.XSHG absolute mean daily return  1.699516e-03 0.002510        15 242 0.249154               1.000000
  002812.XSHE absolute mean daily return  4.761013e-03 0.002004        15 242 0.008769               1.000000
  002438.XSHE absolute mean daily return  2.432229e-03 0.001844        15 242 0.093526               1.000000
  600559.XSHG absolute mean daily return  4.936624e-03 0.002941        15 242 0.046601               1.000000
  002757.XSHE absolute mean daily return  1.913696e-03 0.002141        15 242 0.185680               1.000000
  603538.XSHG absolute mean daily return  2.398694e-03 0.002282        15 242 0.146595               1.000000
  002080.XSHE absolute mean daily return  3.215197e-03 0.001896        15 242 0.044946               1.000000
  002051.XSHE absolute mean daily return -1.076722e-03 0.001192        15 242 0.816868               1.000000
  600863.XSHG absolute mean daily return  5.994050e-06 0.000679        15 242 0.496478               1.000000
  300127.XSHE absolute mean daily return  1.845569e-04 0.001692        15 242 0.456576               1.000000
  300037.XSHE absolute mean daily return  4.944572e-03 0.002038        15 242 0.007636               1.000000
  000736.XSHE absolute mean daily return  1.233006e-03 0.002004        15 242 0.269187               1.000000
  603337.XSHG absolute mean daily return  2.047977e-03 0.001947        15 242 0.146491               1.000000
  002013.XSHE absolute mean daily return  2.542604e-03 0.001760        15 242 0.074319               1.000000
  603681.XSHG absolute mean daily return -7.955166e-04 0.001453        15 242 0.707958               1.000000
  002446.XSHE absolute mean daily return -8.878102e-04 0.001490        15 242 0.724308               1.000000
  600136.XSHG absolute mean daily return -2.226861e-03 0.001726        15 242 0.901561               1.000000
  002781.XSHE absolute mean daily return -6.456168e-08 0.001804        15 242 0.500014               1.000000
  002853.XSHE absolute mean daily return  1.725273e-03 0.002070        15 242 0.202271               1.000000
  002379.XSHE absolute mean daily return -5.929430e-04 0.001791        15 242 0.629726               1.000000
  603256.XSHG absolute mean daily return -1.693241e-03 0.001786        15 242 0.828472               1.000000
  300248.XSHE absolute mean daily return  8.682907e-04 0.002306        15 242 0.353255               1.000000
  000862.XSHE absolute mean daily return -1.213867e-04 0.001345        15 242 0.535952               1.000000
  002288.XSHE absolute mean daily return  2.561325e-03 0.002694        15 242 0.170843               1.000000
  688025.XSHG absolute mean daily return  8.391518e-04 0.002394        15 242 0.362962               1.000000
  002982.XSHE absolute mean daily return  3.441871e-03 0.005459        12 168 0.264170               1.000000
  002515.XSHE absolute mean daily return  6.822271e-04 0.001573        15 242 0.332250               1.000000
  603896.XSHG absolute mean daily return  9.693228e-04 0.001759        15 242 0.290787               1.000000
  002940.XSHE absolute mean daily return  1.948995e-03 0.001893        15 242 0.151615               1.000000
  002319.XSHE absolute mean daily return -8.272816e-04 0.002175        15 242 0.648170               1.000000
  601827.XSHG absolute mean daily return -1.040256e-03 0.001251        11 141 0.797196               1.000000
  688321.XSHG absolute mean daily return -1.237082e-03 0.002047        15 242 0.727170               1.000000
  002916.XSHE absolute mean daily return  6.651180e-04 0.002075        15 242 0.374298               1.000000
  002159.XSHE absolute mean daily return  2.077142e-04 0.001352        15 242 0.438959               1.000000
  000042.XSHE absolute mean daily return -5.364627e-04 0.000939        15 242 0.716089               1.000000
  300487.XSHE absolute mean daily return  1.322090e-03 0.002145        15 242 0.268819               1.000000
  000755.XSHE absolute mean daily return -3.696558e-04 0.000872        15 242 0.664207               1.000000
  002865.XSHE absolute mean daily return  1.018432e-03 0.001701        15 242 0.274658               1.000000
  002463.XSHE absolute mean daily return -5.274744e-04 0.001612        15 242 0.628262               1.000000
  300210.XSHE absolute mean daily return  1.771006e-03 0.002268        15 242 0.217461               1.000000
  002380.XSHE absolute mean daily return  3.331700e-04 0.001315        15 242 0.399982               1.000000
  300428.XSHE absolute mean daily return  7.272308e-04 0.001779        15 242 0.341328               1.000000
  002459.XSHE absolute mean daily return  6.136223e-03 0.002979        15 242 0.019708               1.000000
  002724.XSHE absolute mean daily return  5.083182e-04 0.001683        15 242 0.381323               1.000000
  600768.XSHG absolute mean daily return -1.226456e-03 0.001407        15 242 0.808299               1.000000
  300295.XSHE absolute mean daily return -2.260584e-05 0.001470        15 242 0.506137               1.000000
  603901.XSHG absolute mean daily return -3.159973e-04 0.002134        15 242 0.558853               1.000000
  601168.XSHG absolute mean daily return  2.980554e-03 0.002069        15 242 0.074899               1.000000
  002287.XSHE absolute mean daily return  1.260402e-03 0.002597        15 242 0.313731               1.000000
  603289.XSHG absolute mean daily return -1.480722e-04 0.001417        15 242 0.541605               1.000000
  300291.XSHE absolute mean daily return -5.860228e-04 0.001597        15 242 0.643174               1.000000
  300140.XSHE absolute mean daily return -9.744840e-04 0.001962        15 242 0.690303               1.000000
  688180.XSHG absolute mean daily return -4.518497e-03 0.003536        10 115 0.899371               1.000000
  603053.XSHG absolute mean daily return -1.472223e-03 0.001380        15 242 0.856916               1.000000
  603888.XSHG absolute mean daily return -1.243101e-05 0.001990        15 242 0.502492               1.000000
  002519.XSHE absolute mean daily return  4.934145e-04 0.001436        15 242 0.365590               1.000000
  002119.XSHE absolute mean daily return -5.966241e-04 0.001404        15 242 0.664603               1.000000
  600215.XSHG absolute mean daily return -3.401959e-04 0.002397        15 242 0.556423               1.000000
  002900.XSHE absolute mean daily return -7.883672e-04 0.001037        15 242 0.776440               1.000000
  603363.XSHG absolute mean daily return  4.722478e-04 0.002012        15 242 0.407231               1.000000
  600435.XSHG absolute mean daily return  6.752880e-04 0.001648        15 242 0.340969               1.000000
  600823.XSHG absolute mean daily return  6.640846e-04 0.001991        15 242 0.369383               1.000000
  601139.XSHG absolute mean daily return  3.672625e-06 0.001169        15 242 0.498747               1.000000
  601919.XSHG absolute mean daily return  3.864673e-03 0.002424        15 242 0.055411               1.000000
  002093.XSHE absolute mean daily return -6.751356e-05 0.001614        15 242 0.516686               1.000000
  300241.XSHE absolute mean daily return  4.039454e-04 0.002259        15 242 0.429048               1.000000
  000016.XSHE absolute mean daily return  2.468418e-03 0.003441        15 242 0.236548               1.000000
  002140.XSHE absolute mean daily return  2.501291e-04 0.001327        15 242 0.425250               1.000000
  002069.XSHE absolute mean daily return  2.382047e-03 0.002695        15 242 0.188390               1.000000
  600080.XSHG absolute mean daily return  7.766123e-05 0.001786        15 242 0.482658               1.000000
  603189.XSHG absolute mean daily return  4.540692e-04 0.002063        15 242 0.412880               1.000000
  300502.XSHE absolute mean daily return  3.440611e-03 0.002511        15 242 0.085308               1.000000
  600926.XSHG absolute mean daily return  2.325804e-03 0.001404        15 242 0.048792               1.000000
  600128.XSHG absolute mean daily return -4.549530e-04 0.000933        15 242 0.687102               1.000000
  002522.XSHE absolute mean daily return  5.056635e-04 0.001443        15 242 0.363024               1.000000
  603823.XSHG absolute mean daily return -5.559555e-04 0.001623        15 242 0.634070               1.000000
  002503.XSHE absolute mean daily return  1.623158e-04 0.003324        15 242 0.480525               1.000000
  603328.XSHG absolute mean daily return -1.139741e-03 0.001226        15 242 0.823777               1.000000
  000636.XSHE absolute mean daily return  4.063149e-03 0.002620        15 242 0.060492               1.000000
  000048.XSHE absolute mean daily return  7.799396e-04 0.002491        15 242 0.377084               1.000000
  300300.XSHE absolute mean daily return -2.661041e-03 0.002185        15 242 0.888414               1.000000
  002240.XSHE absolute mean daily return  5.544119e-03 0.002765        15 242 0.022483               1.000000
  603421.XSHG absolute mean daily return -1.004214e-03 0.001164        15 242 0.805950               1.000000
  600693.XSHG absolute mean daily return -6.552006e-04 0.001038        15 242 0.736015               1.000000
  603713.XSHG absolute mean daily return  5.395673e-03 0.002310        15 242 0.009739               1.000000
  603633.XSHG absolute mean daily return  7.736425e-04 0.002351        15 242 0.371058               1.000000
  300641.XSHE absolute mean daily return  8.340848e-04 0.001769        15 242 0.318603               1.000000
  601186.XSHG absolute mean daily return -8.805819e-04 0.001089        15 242 0.790572               1.000000
  002636.XSHE absolute mean daily return  6.263295e-04 0.001684        15 242 0.354990               1.000000
  603505.XSHG absolute mean daily return  1.729234e-03 0.001501        15 242 0.124726               1.000000
  000411.XSHE absolute mean daily return  2.672643e-03 0.003243        15 242 0.204901               1.000000
  600137.XSHG absolute mean daily return -3.347286e-04 0.001402        15 242 0.594350               1.000000
  002612.XSHE absolute mean daily return  4.727296e-03 0.003249        15 242 0.072840               1.000000
  300008.XSHE absolute mean daily return  3.939690e-03 0.003798        15 242 0.149772               1.000000
  002278.XSHE absolute mean daily return  6.994164e-04 0.001597        15 242 0.330672               1.000000
  002157.XSHE absolute mean daily return  7.093926e-04 0.002272        15 242 0.377408               1.000000
  600551.XSHG absolute mean daily return -2.050787e-04 0.001198        15 242 0.567971               1.000000
  603429.XSHG absolute mean daily return  1.342475e-03 0.002507        15 242 0.296123               1.000000
  300070.XSHE absolute mean daily return  2.663688e-04 0.001740        15 242 0.439160               1.000000
  300802.XSHE absolute mean daily return  6.047715e-04 0.002083        15 242 0.385799               1.000000
  300669.XSHE absolute mean daily return  2.385733e-03 0.001248        15 242 0.027991               1.000000
  688363.XSHG absolute mean daily return  2.822731e-03 0.001779        15 242 0.056250               1.000000
  600119.XSHG absolute mean daily return  1.704425e-03 0.002176        15 242 0.216743               1.000000
  600490.XSHG absolute mean daily return  3.388913e-04 0.001756        15 242 0.423489               1.000000
  002112.XSHE absolute mean daily return  1.192861e-03 0.001972        15 242 0.272633               1.000000
  601456.XSHG absolute mean daily return  1.336250e-02 0.008365        10 103 0.055075               1.000000
  002199.XSHE absolute mean daily return -9.734032e-04 0.001645        15 242 0.723019               1.000000
  300218.XSHE absolute mean daily return  1.593639e-03 0.002113        15 242 0.225338               1.000000
  002259.XSHE absolute mean daily return  1.857022e-03 0.001868        15 242 0.160139               1.000000
  002602.XSHE absolute mean daily return -1.044351e-03 0.001597        15 242 0.743435               1.000000
  300144.XSHE absolute mean daily return  3.498766e-04 0.001524        15 242 0.409215               1.000000
  000612.XSHE absolute mean daily return  2.656302e-03 0.002793        15 242 0.170752               1.000000
  002315.XSHE absolute mean daily return  9.302396e-04 0.002470        15 242 0.353257               1.000000
  600376.XSHG absolute mean daily return -9.308241e-04 0.001088        15 242 0.803783               1.000000
  600779.XSHG absolute mean daily return  2.482378e-03 0.001808        15 242 0.084890               1.000000
  002631.XSHE absolute mean daily return -2.983862e-05 0.001256        15 242 0.509478               1.000000
  601019.XSHG absolute mean daily return -6.043165e-04 0.000773        15 242 0.782936               1.000000
  603506.XSHG absolute mean daily return  1.568948e-04 0.001727        15 242 0.463813               1.000000

下面对未校正与 Bonferroni 校正的结果进行可视化对比。左图展示原始p值分布，右图展示校正后的显著性阈值变化。

图 13.1 左图展示未校正的 p 值分布，右图进一步叠加 Bonferroni 校正阈值线，直观对比校正前后显著性判定标准的差异。

fig, axes = plt.subplots(1, 2, figsize=(14, 5))  # 创建1行2列子图布局

sort_indices = np.argsort(calculated_p_values)  # 获取p值排序索引
sorted_p_values = calculated_p_values[sort_indices]  # 按p值升序排列
sorted_small_p_highlight = small_p_highlight[sort_indices]

# 左图：原始 p 值
ax1 = axes[0]  # 选取左侧子图
colors = ['red' if flag else 'gray' for flag in sorted_small_p_highlight]  # 红色仅表示原始 p<0.01
ax1.scatter(range(1, target_stock_count + 1), -np.log10(sorted_p_values),  # 在子图中绑制散点图
           c=colors, alpha=0.6, s=30)  # 绘制p值的负对数散点图
ax1.axhline(y=-np.log10(significance_level), color='blue', linestyle='--',  # 添加水平参考线
            linewidth=2, label=rf'$\alpha = {significance_level}$')  # 添加未校正的显著性阈值线
ax1.set_xlabel('排序后的标的', fontsize=11, fontproperties='Source Han Serif SC')  # 设置x轴标签
ax1.set_ylabel(r'$-\log_{10}(\text{p 值})$', fontsize=11)  # 设置y轴标签
ax1.set_title('未校正的 p 值(筛选收益显著标的)', fontsize=13, fontproperties='Source Han Serif SC')  # 设置标题
ax1.legend(prop={'family': 'Source Han Serif SC'})  # 添加图例
ax1.grid(True, alpha=0.3)  # 添加网格线

# 右图：Bonferroni 校正
ax2 = axes[1]  # 选取右侧子图
ax2.scatter(range(1, target_stock_count + 1), -np.log10(sorted_p_values),  # 在子图中绑制散点图
           c=colors, alpha=0.6, s=30)  # 绘制同样的p值散点图
ax2.axhline(y=-np.log10(significance_level), color='blue', linestyle='--',  # 添加水平参考线
            linewidth=2, label=rf'$\alpha = {significance_level}$')  # 添加原始阈值线作为参考
ax2.axhline(y=-np.log10(bonferroni_significance_level), color='green',  # 添加水平参考线
            linestyle='--', linewidth=2,  # 定义linestyle变量
            label=f'Bonferroni: $\\alpha = {bonferroni_significance_level:.6f}$')  # 添加Bonferroni校正后的阈值线
ax2.set_xlabel('排序后的标的', fontsize=11, fontproperties='Source Han Serif SC')  # 设置x轴标签
ax2.set_ylabel(r'$-\log_{10}(\text{p 值})$', fontsize=11)  # 设置y轴标签
ax2.set_title('Bonferroni 校正', fontsize=13, fontproperties='Source Han Serif SC')  # 设置标题
ax2.legend(prop={'family': 'Source Han Serif SC'})  # 添加图例
ax2.grid(True, alpha=0.3)  # 添加网格线

plt.tight_layout()  # 自动调整子图间距
plt.show()  # 显示图表
Font 'rm' does not have a glyph for '\u503c' [U+503c], substituting with a dummy symbol.
Font 'rm' does not have a glyph for '\u503c' [U+503c], substituting with a dummy symbol.
Font 'rm' does not have a glyph for '\u503c' [U+503c], substituting with a dummy symbol.
左右面板按 p 值排序显示负对数 p 值；左侧为原始阈值，右侧为更严格的 Bonferroni 阈值。
图 13.1: Bonferroni 校正效果演示：强势股筛选

真实市场样本里并不知道每只股票的零假设是否真的成立，因此不能从同一批 p 值构造“真值”，也不能据此估计实现的 FWER、FDR 或功效。下面先定义一个仅供随后“已知注入真值”的受控实验使用的度量函数。

# 计算错误度量的通用函数
def calculate_metrics(p_values_input, true_status_input, alpha_threshold):  # 定义函数calculate_metrics
    significant_flags = p_values_input < alpha_threshold  # 根据阈值判定是否显著
    actual_positives = np.sum((significant_flags == 1) & (true_status_input == 1))  # 真阳性数(TP)
    false_positives = np.sum((significant_flags == 1) & (true_status_input == 0))  # 假阳性数(FP)
    true_negatives = np.sum((significant_flags == 0) & (true_status_input == 0))  # 真阴性数(TN)
    false_negatives = np.sum((significant_flags == 0) & (true_status_input == 1))  # 假阴性数(FN)
    any_false_rejection = int(false_positives > 0)  # 单次族错误指示变量，不是FWER
    false_discovery_proportion = false_positives / max(actual_positives + false_positives, 1)  # 单次FDP，不是FDR
    true_positive_rate = actual_positives / (actual_positives + false_negatives) if (actual_positives + false_negatives) > 0 else 0  # TPR：真阳性检出率（统计功效）
    return {'TP': actual_positives, 'FP': false_positives, 'TN': true_negatives, 'FN': false_negatives,
            'AnyFalseRejection': any_false_rejection, 'FDP': false_discovery_proportion, 'TPR': true_positive_rate}

对真实样本，我们只报告两种规则的拒绝数量；“拒绝”不等于已证实的真阳性。

discovery_counts = pd.DataFrame({
    '规则': [f'未校正 p<{significance_level}', 'Bonferroni'],
    '拒绝数量': [uncorrected_significant_flags.sum(), bonferroni_significant_flags.sum()],
    '阈值': [significance_level, bonferroni_significance_level]
})
print(discovery_counts)
           规则  拒绝数量      阈值
0  未校正 p<0.05    34  0.0500
1  Bonferroni     2  0.0001

Bonferroni 通常减少拒绝数量，但这张真实样本表本身不能告诉我们哪些拒绝是错误的。FWER 与功效的经验比较必须转到零假设状态由设计者预先知道的重复模拟中；下一节的注入实验正是这种受控演示，而不是投资绩效证据。

13.5 Benjamini-Hochberg 方法

Benjamini-Hochberg (BH) 方法是控制 FDR 的经典 step-up 方法 (Benjamini 和 Hochberg 1995年)。在独立真零 p 值下有 \(\mathrm{FDR}=q m_0/m\le q\)；其保证可扩展到满足相应正依赖条件的情形，但任意依赖结构不能直接沿用这一结论。

13.5.1 算法步骤

给定 p 值\(p_1, p_2, \ldots, p_m\) 和目标FDR 水平 \(q\)。

  1. 将 p 值按从小到大排序：\(p_{(1)} \leq p_{(2)} \leq \cdots \leq p_{(m)}\)
  2. 找到最大的 \(k\) 使得：\(p_{(k)} \leq \frac{k}{m} q\)
  3. 拒绝所有对应\(p_{(1)}, p_{(2)}, \ldots, p_{(k)}\) 的假设

BH 方法的直观解释

BH 方法使用了一个自适应的阈值：不是固定的 \(\alpha\)，而是随着排序位置 \(k\) 线性增加。

这意味着：

  • 排名靠前的检验可以使用更宽松的阈值

这种自适应策略在控制FDR 的同时，比Bonferroni 有更高的功效。

数学推导：独立情形下的 BH 控制

设 \(H_0\) 是 \(m_0\) 个真零假设的集合，\(R\) 是 BH 的拒绝数，并约定 \(V/\max(R,1)=0\) 当 \(R=0\)。对每个 \(i\in H_0\)，令 \(I_i\) 表示该真零假设被拒绝。则 \[ \mathrm{FDR}=E\!\left[\frac{V}{\max(R,1)}\right] =\sum_{i\in H_0}\sum_{r=1}^{m}\frac1r P(I_i=1,R=r). \] 若 \(I_i=1\) 且 \(R=r\)，BH 的自洽阈值给出 \(P_i\le rq/m\)。把 \(P_i\) 暂时置为 0，并记其余 p 值决定的拒绝数为 \(R^{(-i)}\)；step-up 过程的单调性给出相应的 leave-one-out 事件 \(R^{(-i)}=r\)。在真零 p 值与其余 p 值独立且服从连续 \(U(0,1)\) 时， \[ P(P_i\le rq/m, R^{(-i)}=r) =\frac{rq}{m}P(R^{(-i)}=r). \] 因而 \[ E\!\left[\frac{I_i}{\max(R,1)}\right] =\sum_{r=1}^{m}\frac1r\frac{rq}{m}P(R^{(-i)}=r) =\frac qm. \] 对真零假设求和得到 \(\mathrm{FDR}=qm_0/m\le q\)；若真零 p 值只满足 super-uniform，上式相应变为不等式。任意依赖下不能使用这一步分解。全零假设时 \(V=R\)，所以 FDR=FWER；“恰等于 \(q\)”还需要独立、连续均匀等条件，不能无条件陈述。

13.5.2 中国案例：真实股票的描述性筛查与受控模拟

这里明确分开两种证据。真实 A 股样本没有“已知真 Alpha”标签：先从每日个股收益中减去当日等权市场收益，再对每只股票的平均市场调整收益做 HAC 均值检验。输出只能描述当前样本中的拒绝集合；它不能报告 TP、FP、FDP 或功效，也不能把拒绝解释成持续可交易的 Alpha。BH 的独立或正依赖保证未必适用于共同市场冲击，因此本例同时报告 BY 作为任意依赖下更保守的敏感性分析。

列表 13.1: 真实股票描述性筛查：数据与 HAC 推断设置 
import os
from pathlib import Path
import numpy as np
import pandas as pd
import statsmodels.api as sm
from scipy.stats import norm
from statsmodels.stats.multitest import multipletests
DATA_DIR = Path(os.environ['BOOK_DATA_DIR']).expanduser().resolve()
if not DATA_DIR.is_dir():
    raise FileNotFoundError(f'BOOK_DATA_DIR 不存在或不是目录: {DATA_DIR}')
stock_price_history = pd.read_hdf(  # 在存储层只读取 2021 年真实行情
    DATA_DIR / 'stock' / 'stock_price_post_adjusted.h5',  # 指定后复权股价文件
    where="date>='2021-01-01' & date<'2022-01-01'",  # 固定市场调整筛查年份
    columns=['date', 'order_book_id', 'close'],  # 只读取 HAC 调整所需字段
).reset_index()  # 恢复日期与证券索引列
列表 13.2: 真实股票描述性筛查：市场调整收益与资格规则 
stock_price_history['date'] = pd.to_datetime(stock_price_history['date'])
year_stock_data = stock_price_history.loc[stock_price_history['date'].dt.year.eq(2021), ['date', 'order_book_id', 'close']].copy()
year_stock_data['ret'] = year_stock_data.groupby('order_book_id')['close'].pct_change()
market_return = year_stock_data.groupby('date')['ret'].mean().rename('market_ret')
year_stock_data = year_stock_data.join(market_return, on='date')
year_stock_data['market_adjusted_ret'] = year_stock_data['ret'] - year_stock_data['market_ret']
valid_stock_distribution = year_stock_data.groupby('order_book_id')['market_adjusted_ret'].count()
valid_stock_codes = valid_stock_distribution.loc[valid_stock_distribution.ge(100)].index.sort_values().tolist()
列表 13.3: 真实股票描述性筛查：逐股票 HAC 均值检验 
finance_records = []
for stock_code in valid_stock_codes:
    stock_returns = year_stock_data.loc[year_stock_data['order_book_id'].eq(stock_code), 'market_adjusted_ret'].dropna()
    hac_fit = sm.OLS(stock_returns.to_numpy(), np.ones((len(stock_returns), 1))).fit(cov_type='HAC', cov_kwds={'maxlags': 5})
    finance_records.append({'order_book_id': stock_code, 'n': len(stock_returns), 'mean_market_adjusted_return': hac_fit.params[0], 'hac_se': hac_fit.bse[0], 'z_statistic': hac_fit.tvalues[0], 'p_value': 2 * norm.sf(abs(hac_fit.tvalues[0]))})
finance_results = pd.DataFrame(finance_records)
finance_p_values = finance_results['p_value'].to_numpy()
finance_t_statistics = finance_results['z_statistic'].to_numpy()
total_test_count = len(finance_results)
finance_results['bh_adjusted_p'] = multipletests(finance_p_values, method='fdr_bh')[1]
finance_results['by_adjusted_p'] = multipletests(finance_p_values, method='fdr_by')[1]
finance_results['bh_reject'] = finance_results['bh_adjusted_p'].le(0.05)
finance_results['by_reject'] = finance_results['by_adjusted_p'].le(0.05)

图 13.2 展示排序 p 值与 BH 边界，并按拒绝方向标出真实样本中的描述性发现。颜色编码的是统计决定，不是真值标签。

import matplotlib.pyplot as plt
sort_indices = np.argsort(finance_p_values)
sorted_p_values = finance_p_values[sort_indices]
bh_thresholds = np.arange(1, total_test_count + 1) / total_test_count * 0.05
fig, axes = plt.subplots(1, 2, figsize=(15, 6))
axes[0].scatter(np.arange(1, total_test_count + 1), -np.log10(sorted_p_values), color='gray', alpha=0.5, s=10)
axes[0].plot(np.arange(1, total_test_count + 1), -np.log10(bh_thresholds), color='green', label='BH 临界线')
axes[0].set(xlabel='排序后的标的', ylabel=r'$-\log_{10}(p)$', title='排序 p 值与 BH 边界')
axes[0].legend(prop={'family': 'Source Han Serif SC'})
plot_colors = np.where(finance_results['bh_reject'] & finance_results['z_statistic'].gt(0), 'red', np.where(finance_results['bh_reject'], 'blue', 'gray'))
axes[1].scatter(finance_t_statistics, -np.log10(finance_p_values), c=plot_colors, alpha=0.5, s=10)
axes[1].set(xlabel='HAC 统计量', ylabel=r'$-\log_{10}(p)$', title='BH 拒绝方向（非真值标签）')
plt.tight_layout()
plt.show()
左图比较排序 p 值与 BH 临界线；右图按 BH 拒绝的正负方向着色，颜色不表示真实 Alpha。
图 13.2: 真实 A 股市场调整收益的描述性多重检验
表 13.2: 真实股票样本中各规则的描述性拒绝数 
real_stock_rejection_summary = pd.DataFrame([
    {'method': '未校正', 'rejections': int((finance_p_values <= 0.05).sum())},
    {'method': 'Bonferroni', 'rejections': int((finance_p_values <= 0.05 / total_test_count).sum())},
    {'method': 'BH', 'rejections': int(finance_results['bh_reject'].sum())},
    {'method': 'BY', 'rejections': int(finance_results['by_reject'].sum())},
])
print(real_stock_rejection_summary.to_string(index=False))
    method  rejections
       未校正         173
Bonferroni           0
        BH           0
        BY           0

表 13.2 不含 TP、FP、FDP 或功效，因为真实样本没有可观测真值。要评价这些量，必须转入下面的独立受控实验。

列表 13.4: 已知真值模拟：共同冲击、AR(1) 误差与 HAC p 值 
def simulate_dependent_family(random_generator, test_count=100, time_count=160, nonnull_count=15):
    common_factor = random_generator.normal(size=time_count)
    innovations = random_generator.normal(size=(time_count, test_count))
    errors = np.zeros_like(innovations)
    for time_index in range(1, time_count):
        errors[time_index] = 0.45 * errors[time_index - 1] + innovations[time_index]
    means = np.zeros(test_count)
    means[:nonnull_count] = 0.32
    observations = means + 0.35 * common_factor[:, None] + errors
    p_values = np.empty(test_count)
    for test_index in range(test_count):
        fitted = sm.OLS(observations[:, test_index], np.ones((time_count, 1))).fit(cov_type='HAC', cov_kwds={'maxlags': 5})
        p_values[test_index] = norm.sf(fitted.tvalues[0])
    return p_values, np.arange(test_count) < nonnull_count
表 13.3: 依赖数据下重复模拟的错误发现与功效汇总 
simulation_rng = np.random.default_rng(2025)
simulation_records = []
for repetition in range(200):
    simulation_p, simulation_truth = simulate_dependent_family(simulation_rng)
    for method_name, method_code in [('BH', 'fdr_bh'), ('BY', 'fdr_by')]:
        simulation_reject = multipletests(simulation_p, alpha=0.05, method=method_code)[0]
        rejection_count = int(simulation_reject.sum())
        false_count = int((simulation_reject & ~simulation_truth).sum())
        true_count = int((simulation_reject & simulation_truth).sum())
        simulation_records.append({'repetition': repetition, 'method': method_name, 'fdp': false_count / max(rejection_count, 1), 'power': true_count / simulation_truth.sum(), 'null_rejection_rate': false_count / (~simulation_truth).sum()})
simulation_results = pd.DataFrame(simulation_records)
simulation_summary = simulation_results.groupby('method')[['fdp', 'power', 'null_rejection_rate']].agg(['mean', 'std'])
print(simulation_summary)
             fdp               power           null_rejection_rate          
            mean       std      mean       std                mean       std
method                                                                      
BH      0.134089  0.131555  0.444000  0.182720            0.015647  0.018465
BY      0.046225  0.101806  0.242667  0.153724            0.003118  0.007316

表 13.3 才能用已知真值计算每次 FDP、功效和零假设拒绝率，并以 200 次重复的均值与标准差描述长期表现。模拟含共同冲击与 AR(1) 误差，单项 p 值采用 HAC；它检验本数据生成机制下的经验行为，不是对所有依赖结构的普遍证明。

13.6 其他多重检验校正方法

除了 Bonferroni 和BH 方法，还有其他重要的校正方法。

13.6.1 Holm 方法

Holm 方法（也称 Holm-Bonferroni）是逐步下降的 FWER 程序 (Holm 1979年)。只要各真零 p 值边际有效，它在任意依赖下控制 FWER，并且不弱于单步 Bonferroni。

算法：

  1. 将p 值排序：\(p_{(1)} \leq p_{(2)} \leq \cdots \leq p_{(m)}\)
  2. 对于 \(k = 1, 2, \ldots, m\)：
    • 如果 \(p_{(k)} > \frac{\alpha}{m - k + 1}\)，停止检验
    • 否则，拒绝\(H_{(k)}\)，继续到下一个

13.6.2 Storey’s q-value

q-value 是p 值在 FDR 框架下的类比。

  • p 值：\(\Pr(\text{观察到的数据或更极端 } | H_0 \text{ 为真})\)
  • q-value：包含该检验的拒绝区域可达到的最小估计 FDR

Storey q-value 把每个有序检验与包含该检验的拒绝区域联系起来，可解释为在给定模型和估计规则下所能达到的最小估计 FDR (Storey 2002年)。它衡量错误发现风险，不把任何拒绝自动认证为有效信号。 下面的极简估算器使用高 p 值区域估计真零假设比例 \(\pi_0\)。这一估计依赖真零 p 值校准、阈值选择和依赖结构；有限样本中未必比 BH 更可靠，也不能确认被拒绝标的具有可复现的超额收益。

列表 13.5: q-value 计算示例 
def calculate_q_value(input_p_values):  # 定义函数calculate_q_value
    '''计算 q-value（Storey 方法），参数: input_p_values -- p值数组，返回: (q_values数组, pi0估计值)'''
    test_count = len(input_p_values)  # 获取检验总数

    # 估计 π0（真实零假设的比例），利用高p值区域的均匀分布特性
    lambda_threshold = 0.5  # 设置λ阈值为0.5
    # 高于λ的p值个数除以理论预期值，即为π0的估计
    # 计算总和
    null_proportion_estimate = np.sum(input_p_values > lambda_threshold) / (test_count * (1 - lambda_threshold))
    null_proportion_estimate = min(null_proportion_estimate, 1.0)  # π0不超过1

    sort_indices = np.argsort(input_p_values)  # 对p值从小到大排序
    sorted_p_values = input_p_values[sort_indices]  # 排序后的p值

    calculated_q_values = np.zeros(test_count)  # 初始化q值数组
    # 最大排名位置的q值 = π0 * 最大p值
    calculated_q_values[-1] = null_proportion_estimate * sorted_p_values[-1]  # 执行数据处理操作

    # 从倒数第二位开始逆序计算，q(i) = min(π0*m*p(i)/(i+1), q(i+1))
    for i in range(test_count-2, -1, -1):  # 遍历循环
        calculated_q_values[i] = min(  # 执行数据处理操作
            null_proportion_estimate * test_count * sorted_p_values[i] / (i + 1),  # 执行数据处理操作
            calculated_q_values[i + 1]  # 执行数据处理操作
        )  # 确保q值单调递增

    original_order_q_values = np.empty(test_count)  # 创建原始顺序q值数组
    original_order_q_values[sort_indices] = calculated_q_values  # 恢复原始排列顺序

    return original_order_q_values, null_proportion_estimate  # 返回q值和π0

使用上述 calculate_q_value 函数对全市场股票的 p 值进行 q-value 计算，展示结果。

# 对全市场股票p值计算q-value
example_q_values, example_pi0_estimate = calculate_q_value(finance_p_values)  # 调用Storey q-value函数

print('\nq-value 计算结果：')  # 打印结果标题
print('='*60)  # 打印分隔线
print(f'估计的 π0（真实零假设比例）: {example_pi0_estimate:.2%}')  # 打印π0估计值

print(f'\n前10个标的的 p-value 和 q-value：')  # 打印前10个标的标题
for i in range(10):  # 遍历前10个标的
    print(f'标的 {i+1:3d}: p-value={finance_p_values[i]:.4f}, q-value={example_q_values[i]:.4f}')  # 逐行打印p值和q值

q-value 计算结果：
============================================================
估计的 π0（真实零假设比例）: 100.00%

前10个标的的 p-value 和 q-value：
标的   1: p-value=0.5501, q-value=0.9996
标的   2: p-value=0.2065, q-value=0.9996
标的   3: p-value=0.7594, q-value=0.9996
标的   4: p-value=0.5313, q-value=0.9996
标的   5: p-value=0.0261, q-value=0.9996
标的   6: p-value=0.9507, q-value=0.9996
标的   7: p-value=0.9695, q-value=0.9996
标的   8: p-value=0.3698, q-value=0.9996
标的   9: p-value=0.4136, q-value=0.9996
标的  10: p-value=0.7319, q-value=0.9996

上方输出给出当前样本的 \(\hat\pi_0\) 与 q-value。较大的 q-value 只表示在本检验族、基准和样本期内没有达到预设错误率阈值；它不能证明真实 Alpha 不存在，也不能单独验证有效市场假说。

13.7 实际应用案例

13.7.1 P-hacking 与量化交易中的多重检验灾难

在金融量化研究中，经常存在一种被称为 P-hacking (P值操纵） 与数据窥探 (Data Snooping) 的现象。当研究人员测试成百上千个交易信号（例如不同均线组合、数百个财务因子、各类技术指标），但仅报告那些在历史回测中获得显著高收益（\(p < 0.05\)）的个别策略时，多重检验灾难就发生了：

如果使用了未经校正的经典假设检验：

  • 即使这1000 个策略本质上都只是随机数生成的无意义策略（零假设全部为真），按照 \(\alpha=0.05\) 的阈值，依然会有大约 50 个策略在统计上表现出“显著的高收益”。
  • 若把这些偶然拒绝未经独立验证就上线，样本外失效风险会升高；在全零模拟中其期望超额收益为零，但单次实盘损益仍受市场路径、成本和持仓规则影响，不能预先确定损益方向。

本章把 \(t>3.0\) 仅作为课堂敏感性阈值，用来观察提高单项证据门槛后结论如何变化；它不是 Bonferroni 校正，也不自动控制特定 planned family 的 FWER 或 FDR。正式研究应先冻结检验族与依赖结构，再选用 Holm、BH、BY 等有明确错误度量的方法；若要把某一经验阈值归于特定文献，还须在参考文献表提供可核验正式来源。

13.7.2 案例 1：A股市场技术指标有效性检验

均线规则包含大量相邻且高度重叠的窗口；若只报告历史样本中最小的 p 值，就会产生数据窥探风险。下面以海康威视（002415.XSHE）为教学样本，预先冻结 5 日到 54 日共 50 个窗口，逐项估计“次日持有或空仓”相对同期买入持有的平均收益差。这个案例只比较错误控制规则，不预设任何均线有效或无效。

表 13.4: A股技术指标（均线策略）有效性检验 
import pandas as pd  # 导入pandas数据分析库
import numpy as np  # 导入numpy数值计算库
from scipy import stats  # 导入scipy统计模块
import statsmodels.api as sm  # HAC协方差用于时间依赖收益
from statsmodels.stats.multitest import multipletests  # 比较依赖稳健主校正与条件性对照
import os  # 导入os模块用于跨平台路径处理

# 根据操作系统自动选择本地数据根目录
# 根据操作系统设置数据根目录路径
from pathlib import Path  # 使用跨平台路径对象解析显式数据根
DATA_DIR = Path(os.environ['BOOK_DATA_DIR']).expanduser().resolve()  # 从必需环境变量取得数据根
if not DATA_DIR.is_dir():  # 在读取前验证数据根
    raise FileNotFoundError(f'BOOK_DATA_DIR 不存在或不是目录: {DATA_DIR}')  # 失败即停止，不回退固定路径
path_price = DATA_DIR / 'stock' / 'stock_price_post_adjusted.h5'  # 后复权股价路径

target_stock_data = pd.read_hdf(  # 在存储层直接读取海康威视真实行情
    path_price,  # 指定后复权股价文件
    where="order_book_id='002415.XSHE'",  # 固定唯一教学标的
    columns=['date', 'order_book_id', 'close'],  # 只读取均线检验所需字段
).reset_index()  # 恢复日期与证券索引列
target_stock_data = target_stock_data.sort_values('date').set_index('date')  # 按日期排序并设为索引
stock_close_prices = target_stock_data['close']  # 提取收盘价序列
daily_returns = stock_close_prices.pct_change().shift(-1)  # 计算次日收益率（T日信号预测T+1日收益）

遍历 50 种均线周期（5 日线到 54 日线），构建“收盘价 > MA(n) 则次日持有”策略。每条规则只从均线首次可用日开始，并对相对同期买入持有的平均超额收益使用 HAC 标准误的单尾检验。

indicator_p_values = []  # 存储各策略的p值
indicator_names = []  # 存储各策略名称

# 遍历50种均线周期：MA(5)到MA(54)
for moving_average_window in range(5, 55):  # 遍历循环
    # 计算n日移动平均线
    # 计算均值
    moving_average_series = stock_close_prices.rolling(window=moving_average_window).mean()
    valid_observation = moving_average_series.notna() & daily_returns.notna()
    trading_signal = stock_close_prices.loc[valid_observation] > moving_average_series.loc[valid_observation]
    benchmark_returns = daily_returns.loc[valid_observation]
    strategy_returns = benchmark_returns.where(trading_signal, 0.0)
    strategy_excess_returns = strategy_returns - benchmark_returns

    if len(strategy_excess_returns) > 10:
        hac_lag_count = max(1, int(np.sqrt(len(strategy_excess_returns))))
        hac_mean_model = sm.OLS(
            strategy_excess_returns.to_numpy(), np.ones((len(strategy_excess_returns), 1))
        ).fit(cov_type='HAC', cov_kwds={'maxlags': hac_lag_count})
        t_statistic = float(hac_mean_model.tvalues[0])
        two_sided_p_value = float(hac_mean_model.pvalues[0])
        p_value = two_sided_p_value / 2 if t_statistic > 0 else 1 - two_sided_p_value / 2
    else:  # 默认分支
        p_value = 1.0  # 样本量不足则p值设为1（不显著）

    indicator_p_values.append(p_value)  # 记录p值
    indicator_names.append(f'MA({moving_average_window})')  # 记录策略名称

对 50 个高度依赖的均线策略，以任意依赖下有效的 BY 为 FDR 主分析、Bonferroni 为 FWER 稳健性分析。BH 只作为“尚未证明独立或 PRDS 条件，因此不具保证”的教学对照。

indicators_p_values_array = np.array(indicator_p_values)  # 将p值列表转为numpy数组
indicator_count = len(indicators_p_values_array)  # 获取策略总数（50）

# Bonferroni校正：将每个p值乘以检验总数，上限截断为1
# 定义bonferroni_adjusted_p_values变量
bonferroni_adjusted_p_values = np.minimum(indicators_p_values_array * indicator_count, 1)

by_adjusted_p_values = multipletests(indicators_p_values_array, method='fdr_by')[1]  # 任意依赖下的 FDR 主分析
bh_unverified_p_values = multipletests(indicators_p_values_array, method='fdr_bh')[1]  # 保留无条件保证的教学对照

汇总展示 50 个均线策略的原始 p 值、Bonferroni、BY 主分析和无保证的 BH 教学对照。

# 构建包含三种p值的结果数据框
indicator_test_results = pd.DataFrame({  # 构建DataFrame数据表
    'Strategy': indicator_names,  # 策略名称
    'Original p': indicators_p_values_array,  # 原始p值
    'Bonferroni p': bonferroni_adjusted_p_values,  # Bonferroni校正p值
    'BY p (primary)': by_adjusted_p_values,  # 任意依赖下的主分析
    'BH p (conditions unverified)': bh_unverified_p_values  # 未证明条件的教学对照
})  # 完成构建

print('均线策略有效性检验结果（前10个）：')  # 打印标题
print('='*60)  # 打印分隔线
print(indicator_test_results.head(10).round(4).to_string(index=False))  # 打印前10个策略结果

significance_level = 0.05  # 设定显著性水平
significant_original_count = np.sum(indicators_p_values_array < significance_level)  # 未校正显著策略数
significant_bonferroni_count = np.sum(bonferroni_adjusted_p_values < significance_level)  # Bonferroni显著数
significant_by_count = np.sum(by_adjusted_p_values < significance_level)  # BY 主分析拒绝数

print(f'\n在 {indicator_count} 个策略中发现显著策略数量 (α={significance_level})：')  # 打印汇总标题
print(f'未校正: {significant_original_count}')  # 打印未校正显著数
print(f'Bonferroni: {significant_bonferroni_count}')  # 打印Bonferroni显著数
print(f'BY（主分析）: {significant_by_count}')  # 打印依赖稳健主分析拒绝数

if significant_original_count > significant_by_count:  # 比较未校正与依赖稳健主分析
    print('\n结论：未校正拒绝数更高，说明选择性报告风险需要单独控制。')  # 仅陈述现场差异
均线策略有效性检验结果（前10个）：
============================================================
Strategy  Original p  Bonferroni p  BY p (primary)  BH p (conditions unverified)
   MA(5)      0.9856           1.0             1.0                        0.9963
   MA(6)      0.9861           1.0             1.0                        0.9963
   MA(7)      0.9844           1.0             1.0                        0.9963
   MA(8)      0.9864           1.0             1.0                        0.9963
   MA(9)      0.9939           1.0             1.0                        0.9963
  MA(10)      0.9950           1.0             1.0                        0.9963
  MA(11)      0.9963           1.0             1.0                        0.9963
  MA(12)      0.9923           1.0             1.0                        0.9963
  MA(13)      0.9769           1.0             1.0                        0.9963
  MA(14)      0.9852           1.0             1.0                        0.9963

在 50 个策略中发现显著策略数量 (α=0.05)：
未校正: 0
Bonferroni: 0
BY（主分析）: 0

本次运行会现场报告未校正、Bonferroni 和 BY 主分析的拒绝数，不预写固定结果；BH 调整值只展示条件未验证时的差异，不能承载保证。每个检验对比“收盘后形成信号、次日持有或空仓”与同期买入持有的超额收益，并用 HAC 协方差允许收益时间依赖。多重校正只是研究内部防线；交易成本、预注册和独立测试期仍不可缺少。

13.7.3 案例 2：行业平均收益的多重检验

这种对“伪阿尔法”的讨伐同样适用于更为宏观的行业配置研究。 下面以数据中最后一个交易日为锚点，取往前 3 年的可用样本，计算各行业等权日收益。这里的零基准是“平均原始收益为 0”，并非市场或因子模型调整后的 alpha；每个行业使用 HAC 标准误的均值检验。 三幅图比较同一检验族在未校正、Bonferroni 与 Benjamini–Yekutieli（BY）规则下的拒绝结果。HAC 处理每个行业序列的时间依赖，BY 则允许行业 p 值因共同市场冲击而具有任意依赖。校正后的拒绝数由现场输出决定；未拒绝不证明效应为零，拒绝也不证明可交易。

列表 13.6: 行业等权收益检验的数据加载与合同检查 
# 目标：检验 30+ 个一级行业的日均收益率是否显著异于 0
# 这是一个典型的多重检验问题：检验越多，偶然出现小 p 值的机会越高

import pandas as pd  # 数据分析库
import numpy as np  # 数值计算库
import matplotlib.pyplot as plt  # 绘图库
from scipy import stats  # 科学计算统计模块
import statsmodels.api as sm  # 每个行业的 HAC 均值检验

# 1. 加载数据
import os  # 文件系统操作
# 根据操作系统设置数据根目录路径
from pathlib import Path  # 使用跨平台路径对象解析显式数据根
industry_data_root_value = os.environ.get('BOOK_DATA_DIR', '').strip()
if not industry_data_root_value:
    raise RuntimeError({'status': 'stopped', 'reason': 'BOOK_DATA_DIR_missing'})
DATA_DIR = Path(industry_data_root_value).expanduser().resolve()  # 从必需环境变量取得数据根
if not DATA_DIR.is_dir():  # 在读取前验证数据根
    raise FileNotFoundError(f'BOOK_DATA_DIR 不存在或不是目录: {DATA_DIR}')  # 失败即停止，不回退固定路径
path_price = DATA_DIR / 'stock' / 'stock_price_post_adjusted.h5'  # 后复权股价数据路径
path_basic = DATA_DIR / 'stock' / 'stock_basic_data.h5'  # 上市公司基本信息路径
industry_missing_files = [str(path) for path in (path_price, path_basic) if not path.is_file()]
if industry_missing_files:
    raise FileNotFoundError({'status': 'stopped', 'reason': 'missing_input_files', 'paths': industry_missing_files})

stock_basic_data = pd.read_hdf(path_basic)  # 读取上市公司基本信息
with pd.HDFStore(path_price, mode='r') as industry_price_store:  # 打开只读 HDF 元数据入口
    industry_price_row_count = industry_price_store.get_storer('data').nrows  # 读取总行数而不载入全表
    industry_tail = industry_price_store.select('data', start=max(0, industry_price_row_count - 10000), columns=['close'])  # 只读末尾一万行及索引确定日期上界
available_end_date = pd.to_datetime(industry_tail.index.get_level_values('date')).max()  # 从真实尾部索引取得样本最后日期
available_start_date = available_end_date - pd.DateOffset(years=3)  # 冻结往前三年支持域
stock_price_history = pd.read_hdf(  # 在 HDF 存储层读取三年真实行业窗口
    path_price,  # 指定后复权股价文件
    where=f"date>='{available_start_date:%Y-%m-%d}' & date<='{available_end_date:%Y-%m-%d}'",  # 下推三年日期过滤
    columns=['date', 'order_book_id', 'close'],  # 只读取行业收益所需字段
).reset_index()  # 恢复日期与证券索引列

接下来，我们将股价数据与行业分类信息合并，计算各行业的等权日均收益率，为后续的多重检验做数据准备。

列表 13.7 
# 合并行业信息到股价数据
merged_stock_data = pd.merge(  # 合并数据表
    stock_price_history[['order_book_id', 'date', 'close']],  # 选取股票代码、日期和收盘价
    stock_basic_data[['order_book_id', 'industry_name']],  # 选取股票代码和行业分类
    on='order_book_id', how='inner'  # 按股票代码内连接
)  # 完成构建
merged_stock_data['date'] = pd.to_datetime(merged_stock_data['date'])
merged_stock_data = merged_stock_data.sort_values(['order_book_id', 'date'])
assert merged_stock_data['date'].between(available_start_date, available_end_date).all()  # 核验 HDF 日期下推没有越界
merged_stock_data['ret'] = merged_stock_data.groupby('order_book_id', sort=False)['close'].pct_change()  # 公司内按日期计算收益
industry_daily_returns = merged_stock_data.groupby(['date', 'industry_name'])['ret'].mean().unstack()  # 聚合到行业层面（等权平均）
industry_daily_returns = industry_daily_returns.dropna(axis=1, thresh=len(industry_daily_returns)*0.8).dropna()  # 去除空值过多的行业
print(f'分析 {len(industry_daily_returns.columns)} 个行业，共 {len(industry_daily_returns)} 个交易日')  # 输出数据概况

代码会现场报告行业数与交易日数，不预写固定的“82 个行业”或“1210 日”。行业收益共享市场冲击且在时间上相关，因此不能称为相互独立的检验；下面先对每个行业使用 HAC 标准误，再对行业检验族使用允许任意依赖的 BY 校正。

有了行业日收益率数据后，我们用 HAC 协方差检验每个行业的平均日收益是否异于零，然后分别应用 Bonferroni 和 BY 方法校正同一检验族。

列表 13.8 
# 对每个行业检验 H0: 平均日收益率 = 0，用HAC允许时间相关
industry_test_results = []  # 存储各行业检验结果的列表
for industry_name in industry_daily_returns.columns:  # 遍历每个行业
    industry_return_series = industry_daily_returns[industry_name]  # 提取该行业日收益率序列
    hac_lag_count = max(1, int(np.sqrt(len(industry_return_series))))
    industry_mean_model = sm.OLS(
        industry_return_series.to_numpy(), np.ones((len(industry_return_series), 1))
    ).fit(cov_type='HAC', cov_kwds={'maxlags': hac_lag_count})
    t_statistic = float(industry_mean_model.tvalues[0])
    p_value = float(industry_mean_model.pvalues[0])
    industry_test_results.append({  # 记录行业名、p值和平均收益
        # 字典条目定义
        'Industry': industry_name, 'p_value': p_value, 'Mean_Ret': industry_return_series.mean()
    })  # 完成构建
test_results_df = pd.DataFrame(industry_test_results)  # 转换为DataFrame
industry_p_values = test_results_df['p_value'].values  # 提取所有p值
industry_count = len(industry_p_values)  # 行业总数（即检验次数）
test_results_df['p_Bonf'] = np.minimum(industry_p_values * industry_count, 1)  # Bonferroni校正：p值乘以检验次数
列表 13.9: 共享冲击下行业检验族的 BY 调整值 
def benjamini_yekutieli(p_values):
    '''BY 调整值：在任意检验依赖下控制 FDR。'''
    test_count = len(p_values)  # 检验总数
    if test_count == 0 or not np.isfinite(p_values).all():
        raise ValueError({'status': 'stopped', 'reason': 'invalid_industry_p_values'})
    sort_indices = np.argsort(p_values)  # 按p值升序排列的索引
    sorted_p_values = p_values[sort_indices]  # 排序后的p值
    harmonic_factor = np.sum(1 / np.arange(1, test_count + 1))
    by_adjusted_p_values = np.ones(test_count)
    by_adjusted_p_values[-1] = min(sorted_p_values[-1] * harmonic_factor, 1.0)
    for i in range(test_count-2, -1, -1):  # 从后向前逐步调整
        by_adjusted_p_values[i] = min(
            sorted_p_values[i] * test_count * harmonic_factor / (i + 1),
            by_adjusted_p_values[i + 1],
            1.0,
        )
    original_order_p_values = np.empty(test_count)  # 恢复原始顺序的数组
    original_order_p_values[sort_indices] = by_adjusted_p_values  # 将校正后的p值映射回原始位置
    return original_order_p_values  # 返回校正后的p值

test_results_df['p_BY'] = benjamini_yekutieli(industry_p_values)

最后，我们将原始 p 值、Bonferroni 校正 p 值和 BY 校正 p 值并排可视化，展示三种规则在当前检验族中的拒绝差异。

fig, axes = plt.subplots(1, 3, figsize=(16, 5))  # 创建1行3列子图
sorted_test_results_df = test_results_df.sort_values('p_value')  # 按原始p值排序
industry_rank_range = range(1, industry_count + 1)  # 行业排名序列
significance_threshold = 0.05  # 显著性水平
correction_methods = [  # 三种校正方法及其对应列名
    ('Original p-value', 'p_value'),  # 执行数据处理操作
    ('Bonferroni p-value', 'p_Bonf'),  # 执行数据处理操作
    ('BY adjusted p-value', 'p_BY')
]  # 执行数据处理操作
for ax, (title, col) in zip(axes, correction_methods):  # 遍历三个子图
    colors = ['red' if p < significance_threshold else 'gray' for p in sorted_test_results_df[col]]  # 显著为红色，不显著为灰色
    ax.scatter(industry_rank_range, -np.log10(sorted_test_results_df[col]), c=colors, s=30)  # 绘制散点图（y轴为-log10(p)）
    ax.axhline(-np.log10(significance_threshold), color='blue', linestyle='--', label='alpha=0.05')  # 添加显著性阈值线
    ax.set_title(title, fontsize=12, fontproperties='Source Han Serif SC')  # 设置子图标题
    ax.set_xlabel('Rank', fontsize=10)  # x轴标签
    ax.set_ylabel('-log10(p)', fontsize=10)  # y轴标签
    ax.grid(True, alpha=0.3)  # 添加网格线
plt.tight_layout()  # 自动调整子图间距
plt.show()  # 显示图形
significant_industries_df = test_results_df[test_results_df['p_BY'] < 0.05]
print(f'\nBY 校正后拒绝零假设的行业({len(significant_industries_df)}个):')
if len(significant_industries_df) > 0:  # 如果存在显著行业
    print(significant_industries_df[['Industry', 'Mean_Ret', 'p_BY']].head(10))
else:  # 如果没有显著行业
    print('无显著行业')  # 输出提示信息
图 13.3

图 13.3 展示当前有效行业数下三种规则的拒绝集合，阈值由代码按实际检验数计算。图只支持对本次检验族的描述，不能由拒绝数推出市场有效性或持续可交易收益。

13.8 实践指南与建议

13.8.1 选择合适的校正方法

在见识了多重检验能给量化回测带来多么恐怖的幸存者偏差之后，作为一个理性的宽客，我们需要在不同的防守策略之间做出抉择。

在决策树可视化之后，我们用一张汇总表来系统地对比 Bonferroni、Holm 和 BH 三种方法的核心差异，为实际研究中的方法选择提供简明参考。

图 13.4 把 Bonferroni、Holm、BH、BY 与 Storey q-value 放在同一选择框架中，但它是检查清单而不是自动决策器。 当一次假阳性的代价很高时，可预先选择控制 FWER 的 Bonferroni 或 Holm；探索性大规模检验可选择控制 FDR 的方法，但 BH 需要独立或适当正依赖，任意依赖时可用 BY。任何校正都不修复无效的单项 p 值，也不保证收益最大化。

在绘制了 FWER 控制方法的相关文本标签后，在决策树右侧绘制 FDR 控制方法的具体方法列表和典型应用场景。完成决策树所有文本标签和连接线的绘制后，设置坐标轴范围并显示最终的方法选择指南图，同时以文本格式输出方法对比表供读者参考。

# 创建决策树可视化
fig, ax = plt.subplots(figsize=(12, 8))  # 创建子图布局

# 绘制决策树
ax.text(0.5, 0.95, '如何选择多重检验校正方法？',  # 添加文本标签
        # 定义fontsize变量
        fontsize=16, fontproperties='Source Han Serif SC', ha='center', weight='bold')

# 第一层：目标
ax.annotate('你的首要目标是什么？', xy=(0.5, 0.85),  # 添加注释标注
            xytext=(0.5, 0.85), fontsize=12, fontproperties='Source Han Serif SC',  # 设置标注文字的偏移位置
            ha='center', va='center',  # 定义ha变量
            bbox=dict(boxstyle='round,pad=0.5', fc='lightblue', alpha=0.5))  # 定义bbox变量

# 第二层：严格性
# 添加文本标签
ax.text(0.25, 0.70, '控制至少一次误拒的概率\n(FWER)', fontsize=11, fontproperties='Source Han Serif SC',
        ha='center', va='center',  # 定义ha变量
        bbox=dict(boxstyle='round,pad=0.5', fc='lightyellow', alpha=0.5))  # 定义bbox变量

# 添加文本标签
ax.text(0.75, 0.70, '允许一些错误\n(最大化发现)', fontsize=11, fontproperties='Source Han Serif SC',
        ha='center', va='center',  # 定义ha变量
        bbox=dict(boxstyle='round,pad=0.5', fc='lightgreen', alpha=0.5))  # 定义bbox变量

# 第三层：具体方法
# 添加文本标签
ax.text(0.25, 0.50, '使用 FWER 控制方法：', fontsize=11, fontproperties='Source Han Serif SC',
        ha='center', weight='bold')  # 定义ha变量
ax.text(0.25, 0.42, '• 任意依赖：Bonferroni / Holm\n• 独立或特定正依赖：Hochberg',  # 添加文本标签
        fontsize=10, fontproperties='Source Han Serif SC', ha='center',  # 定义fontsize变量
        bbox=dict(boxstyle='round,pad=0.5', fc='white', alpha=0.8))  # 定义bbox变量

# 添加文本标签
ax.text(0.75, 0.50, '使用 FDR 控制方法：', fontsize=11, fontproperties='Source Han Serif SC',
        ha='center', weight='bold')  # 定义ha变量
# 添加文本标签
ax.text(0.75, 0.42, '• 独立/PRDS：BH\n• 任意依赖：BY\n• Storey：需额外条件',
        fontsize=10, fontproperties='Source Han Serif SC', ha='center',  # 定义fontsize变量
        bbox=dict(boxstyle='round,pad=0.5', fc='white', alpha=0.8))  # 定义bbox变量

# 第四层：应用场景
ax.text(0.25, 0.25, '典型应用：', fontsize=11, fontproperties='Source Han Serif SC',  # 添加文本标签
        ha='center', weight='bold')  # 定义ha变量
ax.text(0.25, 0.15, '• 高频交易策略验证\n• 重大投资决策\n• 高风险决策',  # 添加文本标签
        fontsize=10, fontproperties='Source Han Serif SC', ha='center',  # 定义fontsize变量
        bbox=dict(boxstyle='round,pad=0.5', fc='white', alpha=0.8))  # 定义bbox变量

ax.text(0.75, 0.25, '典型应用：', fontsize=11, fontproperties='Source Han Serif SC',  # 添加文本标签
        ha='center', weight='bold')  # 定义ha变量
ax.text(0.75, 0.15, '• 多因子选股\n• 量化策略回测\n• 大规模筛选',  # 添加文本标签
        fontsize=10, fontproperties='Source Han Serif SC', ha='center',  # 定义fontsize变量
        bbox=dict(boxstyle='round,pad=0.5', fc='white', alpha=0.8))  # 定义bbox变量

# 连接线
ax.annotate('', xy=(0.35, 0.70), xytext=(0.45, 0.82),  # 添加注释标注
            arrowprops=dict(arrowstyle='->', lw=2))  # 定义arrowprops变量
ax.annotate('', xy=(0.65, 0.70), xytext=(0.55, 0.82),  # 添加注释标注
            arrowprops=dict(arrowstyle='->', lw=2))  # 定义arrowprops变量

ax.set_xlim(0, 1)  # 设置子图X轴范围
ax.set_ylim(0, 1)  # 设置子图Y轴范围
ax.axis('off')  # 执行数据处理操作
plt.tight_layout()  # 自动调整子图间距
plt.show()  # 显示图形
自上而下的决策图先区分 FWER 与 FDR 目标，再连接到 Bonferroni、Holm、BH 与 q-value 及典型场景。
图 13.4: 多重检验方法选择指南
表 13.5: 多重检验校正方法的控制目标与适用条件 
method_comparison_table = pd.DataFrame([
    {'方法': 'Bonferroni', '控制目标': 'FWER', '依赖条件': '任意依赖', '相对功效': '低'},
    {'方法': 'Holm', '控制目标': 'FWER', '依赖条件': '任意依赖', '相对功效': '中等'},
    {'方法': 'BH', '控制目标': 'FDR', '依赖条件': '独立或 PRDS', '相对功效': '较高'},
    {'方法': 'BY', '控制目标': 'FDR', '依赖条件': '任意依赖', '相对功效': '较低'},
    {'方法': 'Storey q', '控制目标': 'FDR', '依赖条件': '需验证额外条件', '相对功效': '较高'},
])  # 将非绘图输出与图形块分离
print(method_comparison_table.to_string(index=False))
        方法 控制目标     依赖条件 相对功效
Bonferroni FWER     任意依赖    低
      Holm FWER     任意依赖   中等
        BH  FDR 独立或 PRDS   较高
        BY  FDR     任意依赖   较低
  Storey q  FDR  需验证额外条件   较高

方法选择必须同时看错误控制目标和依赖条件。FWER 方法控制的是“至少一次误拒”的概率，不保证某次研究零假阳性；检验个数少也不自动允许取消校正。Bonferroni 与 Holm 不要求独立，Hochberg、BH 与 Storey 类方法需要各自的独立或正依赖等条件；任意依赖下的 FDR 可使用 BY，但前提仍是每个输入 p 值本身有效。少簇、面板共同冲击或重叠时序应先采用相应稳健推断；无法得到有效单项 p 值时必须停止，而不是换一种多重校正补救。

13.9 本章小结

13.9.1 核心要点

  1. 多重检验问题：当同时进行大量假设检验时，假阳性会累积
  2. 两类主要错误度量：
    • FWER（族错误率）：至少一次假阳性的概率
    • FDP：一次实现中假阳性占所有拒绝的比例；FDR 是 FDP 的重复抽样期望
  3. 三种主要校正方法：
    • Bonferroni：最简单但最保守，控制FWER
    • Holm：逐步下降，比 Bonferroni 更有功效
    • Benjamini-Hochberg：控制FDR，适合大规模研究

13.9.2 应用建议

多重检验最佳实践

  1. 预注册研究计划：在看到数据前确定主要假设
  2. 分层检验：将检验分为主要检验和探索性检验
  3. 报告完整结果：不要只报告显著结果
  4. 复制验证：使用独立数据集验证发现
  5. 效应量考虑：不要只看p 值，也要考虑实际意义

13.10 理论来源与前沿

多重检验问题的根源在于：当同时进行 \(m\) 个检验时，即使每个检验都控制在显著性水平\(\alpha\)，至少一次假阳性发生的概率会随\(m\) 增大而迅速上升。经典的 Bonferroni 思想可以用并集上界直接推导出 FWER 控制，但它往往过于保守。随着Holm 的逐步法在保证 FWER 控制的同时提高功效；Benjamini-Hochberg (BH) 则把目标从‘避免任何假阳性’转向‘控制假阳性比例’，成为大规模量化策略回测的工作马。

近年来的前沿主要集中在：

  1. 依赖性检验的校正：当检验统计量相关且无法验证 BH 所需条件时，Benjamini–Yekutieli（BY）用调和级数因子修正阈值，在任意依赖下控制 FDR (Benjamini 和 Yekutieli 2001年)；代价通常是功效下降。更精细的依赖建模只有在其假设可审计时才可替代 BY。
  2. 自适应 FDR 与局部FDR：利用\(\pi_0\)（原假设比例）估计来提升功效，并用经验贝叶斯框架解释‘显著性’。
  3. 与可重复性预注册结合：把多重检验控制嵌入研究流程，强调独立样本复现与透明报告，减少p-hacking。

13.11 M13：冻结检验族与独立复现

本里程碑与 小节 6.3 一一对应，且只在本章实施。

输入与依赖：承接 M12 的冻结公司与日期边界，但不使用 M12 的聚类标签挑选假设。M13 必须逐字节复核 M12 preregistration、三层 roster、排除账本、六份固定输出和稳定 handoff；M12 成功 telemetry 由运行器另行提供，并用独立声明的 SHA-256 核验，不属于稳定 handoff 绑定的内容。缺一项即停止。六份输出的唯一名称是 loadings、scores、silhouettes、cluster-profile、stability-seeds、lineage。数据根来自 BOOK_DATA_DIR；检验清单逐项写明 estimand、单位、方向、形成日、标签窗、最小样本、依赖处理与原始 SE 方法。快照编号、冻结日和校验值必须来自课程冻结清单；缺一项即停止，不得补写虚构值。

冻结规则：在查看 p 值前冻结完整检验族、主分析与稳健性分析、\(q=0.05\)、跳过规则、BY/Holm 选择以及独立复现区间。探索期严格截止于复现起点之前；只有探索期 BY 拒绝的原始逐成员合同进入冻结发现集。复现期不得重新选择行业、方向或窗口，失败候选也须保留在 Holm 家族中并以 \(p=1\) 校正。

运行器必须同时提供月初且严格有序的 BOOK_M13_EXPLORATION_START、BOOK_M13_REPLICATION_ORIGIN 与 BOOK_M13_REPLICATION_END_EXCLUSIVE。探索起点月仅建立基准收盘价，不产生可评价收益；价格、完整月历、复现配对和结果均服从最后一个边界的开区间上界。设计、血缘、handoff 与独立 telemetry 都逐字段携带同一边界和暖启动合同。

收益构造先建立“冻结总体中每只证券 × 完整日历月”的面板，再为每行同时保留相邻的 formation_month 与 label_month；只有两个紧邻日历月的月末收盘价都有限且形成月价格非零时才计算收益。内部缺月会使相邻收益保持缺失，不能用下一次观测跨月补算。探索与复现的含下界、开上界筛选以及结果表的起止日期一律作用于 label_month，所以等于复现起点的标签月属于复现期，其形成月可位于探索期末但不得被误报为复现结果日期。

输出：提交计划/实际/跳过检验清单、显式最大交易日选出的月末价审计表、效应与稳健 SE 表、原始 p 值及 BY 调整值、冻结拒绝集合、带 Holm 调整值和失败原因的非重叠复现表、数据血缘、逐次运行 telemetry 与最终 handoff。期末价表、完整家族表、BY 表、发现集、复现表和血缘稳定序列化，写入同一 preregistration hash，并由稳定 handoff 绑定实际 SHA-256；含时刻、耗时和峰值内存的 telemetry 独立交付及核验，不进入 handoff 的哈希投影。每个跳过项保留在清单中并记录原因，不能从分母静默删除。

教师启动器：许可数据路线先由教师创建并导出一个专属 BOOK_M13_ARTIFACT_DIR，再运行 conda run -n peter python scripts/m13_preflight.py --self-test 与 conda run -n peter python scripts/m13_instructor.py --manifest tests/fixtures/m13/classroom-manifest-v2.json。第二条命令是许可数据命令；未提供制品目录时按 CLI 合同退出 2 并报告 artifact_directory_missing，不是无挂载替代活动。启动器从清单读取必需环境、精确输出名及本章规范标签链，在全新 Python 进程中依次执行；它先验证唯一制品目录，再注册退出守卫并签发 run_id，签发后缺参则以该 ID 发布唯一 stopped。成功时要求十一个最终名称逐项相等；墙钟 wall_seconds 与代码活动 active_seconds 分开记录。

无挂载学生—教师路线：两条命令都从仓库根复制；先删除或改名任何同名旧目录/旧成绩文件，以免把冲突重跑混为新提交。学生运行：

conda run --no-capture-output -n peter python scripts/m13_instructor.py --no-mount-student --output-dir m13-no-mount-submission

预期退出 0，且 m13-no-mount-submission/ 只有 m13-process-evidence-v1.json 与 m13-process-terminal-v1.json。教师收到这两份学生交接文件后，从独立干净进程运行：

conda run --no-capture-output -n peter python scripts/m13_instructor.py --verify-no-mount-submission m13-no-mount-submission --gradebook-output m13-no-mount-grade-v1.json

预期退出 0 并只新增一份确定性成绩记录 m13-no-mount-grade-v1.json。教师端会重跑仓库内 CC0 preflight、复核 manifest 与两份提交的规范字节和 hash，任何增删、篡改或实证完成冒称都在生成成绩前失败。流程终态固定为 status=completed、activity=contract_fixture_only、grade_status=eligible_process_substitute，同时固定 empirical_outputs_created=false 与 empirical_milestone_status=blocked_pending_licensed_snapshot；它绝不是 m13-run-<run_id>-v5.json，也不生成设计、期末价、检验族、发现或复现制品。

当学校未向整班提供许可挂载、教师事先指定本替代路线时，四项各 5 分：边界/hash/schema preflight、行序与推断守卫、fresh-process 生命周期、提交完整性与流程终态，共 20 分。四项全部由教师核验才写入 completed_process_substitute、20/20、课程权重 4%；失败则不写成绩记录，修正后重交。该 20 分是本队列 M13 的公平应急成绩，不因机构数据缺失扣分；以后取得许可数据的实跑只作形成性证据，不加分、不减分、不替换该成绩。成绩完成与实证里程碑状态分离，后者持续 blocked_pending_licensed_snapshot，不得据此声称有效应、BY 发现或独立复现。30 秒字段只约束无数据 preflight；真实数据耗时必须从许可运行 telemetry 读取。

验收与失败条件：面板用公司与年度双聚类或预注册块重抽样；重叠时序用 HAC 或时间块重抽样；任意依赖下默认 BY，只有写清并验证适用条件才可用 BH。每个单项检验须同时满足效应、完整协方差与统计量有限，SE 有限且严格大于 0，原始 p 值有限且位于 \([0,1]\)；否则登记具体失败原因、标记 skipped，并以 \(p=1\) 留在完整校正家族。若有效检验为零、协方差不可估、复现期不足或数据边界失败，结论为“不可评价”，不得改用朴素 p 值。

量规（20 分）：预注册与时点 5 分，依赖稳健单项推断 5 分，多重校正 4 分，独立复现 3 分，结论边界与复现审计 3 分。任何看过结果后删改检验族，整项最高 8 分。评分只认下表的“代码输出—制品文件”配对。

可评分动词 必须由代码打印或断言的输出 必交制品文件 分值
复核 M12 并预注册 M12 全制品 hash、探索/复现边界、planned \(m\) m13-design-v3.json、m13-preregistration-v3.json 5
冻结总体与失败规则 日度键唯一、显式期末日、多行业数、逐证券排除原因、两期键 hash m13-period-end-closes-v1.csv、m13-universe-exclusion-ledger-v3.csv、m13-lineage-v3.json 2
估计稳健单项检验 estimand、有效 \(n\)、HAC lag、效应、SE、原始 p 值 m13-family-v3.csv 3
执行 BY 并冻结发现 planned/estimated/skipped 数、BY 调整值、发现数 m13-by-results-v3.csv、m13-frozen-discoveries-v1.json 4
独立复现并执行 Holm 不重叠断言、冻结候选数、状态/失败原因、Holm p 值 m13-replication-v3.csv 3
记录与限定结论 命令、版本、起止时刻、耗时、峰值内存、加载行数、输入/输出 hash；非因果、非可交易边界 m13-run-<run_id>-v5.json、m13-handoff-v5.json、决策说明 3

上表 20 分只适用于学校已提供许可快照的实证路线。若教师在布置前因机构挂载不可用而统一指定无挂载路线，则改用上一段四项各 5 分的流程量规，并以 m13-no-mount-grade-v1.json 作为唯一成绩簿导入记录；两套量规不可混计。无挂载记录的 20/20 只结清本课程 M13 的 4% 公平替代成绩，不改变实证状态，也不能充当上表任何实证制品。

13.12 练习

13.12.1 概念题

  1. [核心｜难度:1｜时间:8分钟｜分值:4｜项目:无] 解释 FWER 与FDR 的区别。在什么情况下应该控制 FWER 而不是FDR。

  2. [核心｜难度:2｜时间:8分钟｜分值:4｜项目:无] Bonferroni 校正为什么非常保守？在什么情况下这会导致问题？

  3. [核心｜难度:2｜时间:10分钟｜分值:5｜项目:无] Benjamini-Hochberg 方法的核心思想是什么？为什么它比Bonferroni 更有功效。

  4. [拓展｜难度:2｜时间:10分钟｜分值:5｜项目:无] 解释 q-value 的含义。它和p 值有什么关系？

  5. [核心｜难度:2｜时间:10分钟｜分值:5｜项目:无] 什么是“p-hacking”？多重检验校正如何帮助减少 p-hacking 的问题？

13.12.2 应用题

  1. [拓展｜难度:3｜时间:40分钟｜分值:15｜项目:无] 基于本地数据的技术指标有效性检验 使用本地 stock_price_post_adjusted.h5 数据：

    • 选择一只代表性股票（如海康威视）：
    • 构建 50 个不同参数的均线策略（MA5, MA6, …, MA54）。
    • 以“策略减同股买入持有”的日均超额收益为 estimand，用至少覆盖均线窗口的 HAC 标准误形成双侧原始 p 值；不可估项以 \(p=1\) 留在 50 项家族。
    • 以 BY 作为共享日期与嵌套窗口任意依赖下的 FDR 主分析，并以 Bonferroni 作为 FWER 敏感性分析；BH 只有在另交可核验设计证据证明其依赖条件后才可运行，否则标记为 not_run。
    • 量规：HAC 单项推断 5 分，完整家族与失败日志 3 分，BY 主分析 3 分，Bonferroni 敏感性 2 分，BH 条件证据或 not_run 审计 2 分。

  2. [拓展｜难度:3｜时间:50分钟｜分值:18｜项目:无] 财务指标选股的多重检验 使用 financial_statement.h5 和stock_price_post_adjusted.h5：

    • 预注册 5 个可由本地字段稳定构造的比率：ROE、净利率、负债率、资产周转率、流动比率。
    • 采用年报期末后至少一年的保守可得期，并在公司内构造下一年度收益率。
    • 对 5 个因子的秩相关 estimand 使用公司与年度双向聚类标准误；若有效公司或年度簇不足则跳过并登记原因。
    • 因因子与收益检验共享公司—年度冲击，默认使用 BY 控制任意依赖下的 FDR；只有另行以设计证据验证适用依赖条件才可改用 BH。
    • 逐项报告 estimand、有效 \(n\)、公司数、年度数、依赖处理、标准误、原始/调整 p 值及信息时点限制，不把显著性称为“有效选股能力”。

  3. [核心｜难度:3｜时间:55分钟｜分值:20｜项目:M13] 行业动量策略的多重检验 使用 stock_basic_data.h5 和stock_price_post_adjusted.h5：

    • 使用数据中实际可得的行业分类，计算等权月度收益率。
    • 测试“过去N 个月收益率”预测“下月收益率”的能力（动量效应）。
    • 在看任何 p 值前冻结探索期、非重叠复现期、方向、\(N=1,2,\ldots,12\)、逐成员 skip 规则和完整“有效行业 × 12”家族。
    • 每项以标准化回归斜率表示 Pearson 相关 estimand，并使用至少覆盖回看窗口的 HAC 标准误；跨行业、跨窗口默认用 BY 校正。
    • 探索期对完整家族执行 BY，随后把 BY 拒绝项逐字节冻结；复现期只检验该发现集，以 Holm 控制 FWER，失败候选以 \(p=1\) 留在复现家族。
    • 由代码报告计划数、实际检验数、缺失原因、冻结发现数、复现成功/失败数，并按上方证据矩阵提交家族、BY、复现、血缘、运行日志与 handoff 文件；不得只交屏幕输出。

13.12.3 理论题

  1. [拓展｜难度:3｜时间:25分钟｜分值:12｜项目:无] 证明 Benjamini-Hochberg 方法在独立检验的情况下控制FDR。

  2. [核心｜难度:2｜时间:15分钟｜分值:8｜项目:无] 推导 Bonferroni 校正的 FWER 控制性质。

  3. [拓展｜难度:3｜时间:25分钟｜分值:12｜项目:无] 比较 Holm 方法和Hochberg 方法的性质，证明它们都控制 FWER。

  4. [拓展｜难度:2｜时间:25分钟｜分值:10｜项目:无] 研究并总结以下主题：

    • 自适应 FDR 控制方法
    • 依赖性检验的多重校正
    • 贝叶斯多重检验方法

  5. [核心｜难度:2｜时间:15分钟｜分值:12｜项目:无] 四种校正的同族手算：给定已排序 p 值 \(p_{(1:4)}=(0.004,0.011,0.021,0.200)\)，令 FWER 水平与 FDR 目标均为 \(0.05\)。分别写出 Bonferroni、Holm、BH 与 BY 的逐项阈值、最大可拒绝序号及拒绝集合；计算 \(c(4)=\sum_{j=1}^4 1/j\)，并逐项说明四种方法保证所需的依赖条件。答案必须使用“拒绝/不拒绝”与证据强弱措辞，不得给任何假设命题判定真值。

13.13 练习参考解答

展开完整解答、评分点与常见失败模式

评分以题面分值为准；手算调整值允许 \(10^{-6}\) 绝对误差。应用题若缺 estimand、有效 \(n\)、依赖处理或 SE 方法，每项扣 20%；若把朴素 iid p 值直接送入校正、静默删除失败检验或把拒绝解释为可交易，推断与结论项不得分。

13.13.1 概念题解答

  1. FWER vs FDR

    • FWER：\(P(V\ge 1)\)，强调‘一次假阳性都不想要’。适合监管、关键决策、confirmatory 研究（例如临床主要终点）。
    • FDR：\(E\big[\frac{V}{\max(R,1)}\big]\)，允许一定比例的假阳性以换取发现能力，适合探索性大规模研究（量化因子挖掘）。

  2. Bonferroni 为什么保守

    Bonferroni 用并集上界： \[ P\Big(\cup_{j=1}^m \{p_j \le \alpha/m\}\Big) \le \sum_{j=1}^m P(p_j \le \alpha/m) = \alpha. \] 上界一般不紧，尤其当检验正相关且\(m\) 很大时，阈值\(\alpha/m\) 太小导致功效显著下降。

  3. BH 方法的核心思想

    将p 值从小到大排序\(p_{(1)}\le\dots\le p_{(m)}\)，寻找最大的 \(k\) 使得 \[ p_{(k)} \le \frac{k}{m}q, \] 然后拒绝前\(k\) 个假设。阈值随排名放宽，既‘保护尾部’又‘奖励强信号’，因此通常比Bonferroni 更有功效。

  4. q-value 的含义

    q-value 可理解为：在包含某个检验的拒绝区域中，给定模型与估计规则所能达到的最小估计 FDR。它不是该假设为假的后验概率，也不等同于局部 fdr；解释依赖 p 值校准、\(\pi_0\) 估计和依赖条件。

  5. p-hacking 与校正作用

    p-hacking 指研究者在多个模型/指标/子样本上反复试验并只报告显著结果。多重检验校正把这种“搜索成本”显式计入阈值，从而降低‘偶然显著’的概率，但更根本的解决是预注册与复现。

13.13.2 应用题解答

  1. 技术指标有效性检验(完整代码)

    下例把每个参数视为预注册假设，信号在 \(t\) 日收盘后形成并对应 \(t+1\) 日收益。超额收益定义为策略相对同一股票买入持有基准的差；HAC p 值进入 BY 主分析与 Bonferroni 敏感性分析。由于未提供可核验的 BH 依赖条件证据，BH 明确记为 not_run；结果仍未计交易成本、重叠持仓或模型选择后的独立复现。

当前交付 HDF 已现场核验为 table 格式，且 order_book_id、date 均为 data_column；下面将公司过滤直接下推到 read_hdf(where=...)，并只读取日期与收盘价。

列表 13.10: 习题6：环境数据根与单股票输入 
import os  # 读取必需的数据根环境变量
from pathlib import Path  # 以跨平台路径对象拼接文件
import pandas as pd  # 读取并整理行情数据
import numpy as np  # 构造 HAC 回归所需数组
import statsmodels.api as sm  # 估计 HAC 均值检验
from statsmodels.stats.multitest import multipletests  # 执行家族校正
book_data_dir = Path(os.environ['BOOK_DATA_DIR']).expanduser().resolve()  # 必须由运行环境显式指定数据根
if not book_data_dir.is_dir():  # 在读取前验证数据根存在
    raise FileNotFoundError(f'BOOK_DATA_DIR 不存在或不是目录: {book_data_dir}')  # 失败即停止而非回退固定路径
price_path = book_data_dir / 'stock' / 'stock_price_post_adjusted.h5'  # 定位后复权行情文件
one_stock = (pd.read_hdf(  # 在存储层直接读取海康威视真实行情
    price_path,  # 指定后复权股价文件
    where="order_book_id='002415.XSHE'",  # 下推唯一证券过滤
    columns=['date', 'close'],  # 只读取均线家族所需字段
).reset_index().sort_values('date')[['date', 'close']].copy())  # 恢复日期索引并固定顺序
one_stock['future_return'] = one_stock['close'].pct_change().shift(-1)  # 对齐信号形成后的下一日收益
列表 13.11: 习题6：50个均线规则的 HAC 检验与跳过日志 
ma_records, ma_skip_log = [], []  # 保留成功检验和失败原因
for window in range(5, 55):  # 执行完整预注册窗口族
    moving_average = one_stock['close'].rolling(window, min_periods=window).mean()  # 排除 warm-up 期
    valid_observation = moving_average.notna() & one_stock['future_return'].notna()  # 只留完整信号与标签
    valid_count = int(valid_observation.sum())  # 记录该规则的有效样本量
    minimum_count = max(60, 4 * window)  # 使有效样本相对窗口和 HAC 带宽足够
    if valid_count < minimum_count:  # 在拟合前执行最小样本保护
        ma_skip_log.append({'window': window, 'n': valid_count, 'reason': f'需要 n>={minimum_count}'})  # 明确登记跳过原因
        continue  # 不以朴素 t 检验替代
    benchmark_return = one_stock.loc[valid_observation, 'future_return']  # 提取同日期买入持有基准
    strategy_return = benchmark_return.where(
        one_stock.loc[valid_observation, 'close'] > moving_average.loc[valid_observation], 0.0
    )  # 构造长仓或现金策略收益
    strategy_excess = strategy_return - benchmark_return  # 定义策略相对同股基准的 estimand
    if strategy_excess.nunique() < 2:  # 防止零方差导致协方差不可估
        ma_skip_log.append({'window': window, 'n': valid_count, 'reason': '超额收益为常数'})  # 登记数值失败
        continue  # 跳过不可估规则
    hac_lag_count = max(window, int(np.sqrt(valid_count)))  # HAC 同时覆盖窗口和经验带宽
    fitted_mean = sm.OLS(
        strategy_excess.to_numpy(), np.ones((len(strategy_excess), 1))
    ).fit(cov_type='HAC', cov_kwds={'maxlags': hac_lag_count})  # 使用 HAC 推断均值
    estimate, standard_error = strategy_excess.mean(), fitted_mean.bse[0]  # 提取效应与HAC标准误
    covariance_values = np.asarray(fitted_mean.cov_params(), dtype=float)  # 读取完整HAC协方差
    test_statistic, raw_p_value = fitted_mean.tvalues[0], fitted_mean.pvalues[0]  # 提取双侧统计量和p值
    inference_valid = np.isfinite(estimate) and np.isfinite(covariance_values).all() and np.isfinite(standard_error) and standard_error > 0 and np.isfinite(test_statistic) and np.isfinite(raw_p_value) and 0 <= raw_p_value <= 1  # 校正前验证全部推断量
    if not inference_valid:
        ma_skip_log.append({'window': window, 'n': valid_count, 'reason': '效应/协方差/统计量/p值须有限，SE须有限且>0，p须在[0,1]'})  # 数值失败不得冒充估计成功
        continue  # 保留计划成员并在完整家族中赋p=1
    ma_records.append({'window': window, 'status': 'estimated', 'failure_reason': None, 'estimand': '日均策略减买入持有收益', 'estimate': estimate, 'n': valid_count, 'dependence': '重叠均线与日序列相关', 'se_method': f'HAC(maxlags={hac_lag_count})', 'se': standard_error, 'statistic': test_statistic, 'p_value': raw_p_value, 'start': one_stock.loc[valid_observation, 'date'].min(), 'end': one_stock.loc[valid_observation, 'date'].max()})  # 留存完整推断字段
列表 13.12: 习题6：家族校正与最小样本审计 
ma_results = pd.DataFrame(ma_records, columns=['window', 'status', 'failure_reason', 'estimand', 'estimate', 'n', 'dependence', 'se_method', 'se', 'statistic', 'p_value', 'start', 'end'])  # 汇总成功估计的规则并固定空表schema
ma_planned_family = pd.DataFrame({'window': list(range(5, 55))})  # 在校正前保留全部50个计划假设
ma_skip_table = pd.DataFrame(ma_skip_log).rename(columns={'reason': 'failure_reason'})  # 统一失败字段
ma_family_table = ma_planned_family.merge(ma_results, on='window', how='left', validate='one_to_one')  # 以计划家族为权威左表
if not ma_skip_table.empty:  # 将跳过原因写回完整计划表
    ma_family_table = ma_family_table.merge(ma_skip_table[['window', 'failure_reason']], on='window', how='left', suffixes=('', '_skip'), validate='one_to_one')  # 一对一合并失败日志
    ma_family_table['failure_reason'] = ma_family_table['failure_reason'].fillna(ma_family_table.pop('failure_reason_skip'))  # 保留唯一失败字段
ma_family_table['status'] = ma_family_table['status'].fillna('skipped')  # 明确每个计划假设的状态
ma_family_table['p_value'] = ma_family_table['p_value'].where(ma_family_table['status'].eq('estimated'), 1.0)  # 跳过项在完整家族中显式保存p=1
ma_correction_p = ma_family_table['p_value']  # 全部50项都有合法校正输入
ma_family_table['by_adjusted_p'] = multipletests(ma_correction_p, method='fdr_by')[1]  # 任意依赖下的FDR主分析
ma_family_table['bonferroni_adjusted_p'] = multipletests(ma_correction_p, method='bonferroni')[1]  # 任意依赖下的FWER稳健性分析
ma_bh_condition_verified = False  # 嵌套窗口尚无BH正依赖或独立性设计证据
ma_bh_condition_evidence = None  # 未提交证据时不得计算BH调整值
if ma_bh_condition_verified and ma_bh_condition_evidence:  # 只有双重条件成立才运行BH
    ma_family_table['bh_adjusted_p'] = multipletests(ma_correction_p, method='fdr_bh')[1]  # 条件性BH补充分析
    ma_bh_status = 'run_with_verified_condition'  # 记录BH已经授权运行
else:
    ma_family_table['bh_adjusted_p'] = np.nan  # 不用未验证BH数值暗示有效推断
    ma_bh_status = 'not_run'  # 明示依赖条件证据缺失
print({'planned': len(ma_family_table), 'actual': len(ma_results), 'skipped': len(ma_skip_log), 'primary': 'BY', 'sensitivity': 'Bonferroni', 'bh_status': ma_bh_status})  # 审计固定家族和方法层级
assert len(ma_family_table) == len(ma_results) + len(ma_skip_log) == 50  # planned=estimated+skipped且固定为50
print(ma_family_table.to_string(index=False))  # 显示全部计划项、状态、依赖与调整值
{'planned': 50, 'actual': 50, 'skipped': 0, 'primary': 'BY', 'sensitivity': 'Bonferroni', 'bh_status': 'not_run'}
 window    status failure_reason    estimand  estimate    n dependence       se_method       se  statistic  p_value      start        end  by_adjusted_p  bonferroni_adjusted_p  bh_adjusted_p
      5 estimated           None 日均策略减买入持有收益 -0.000500 3784 重叠均线与日序列相关 HAC(maxlags=61) 0.000229  -2.187065 0.028738 2010-06-03 2025-12-30       0.583796               1.000000            NaN
      6 estimated           None 日均策略减买入持有收益 -0.000515 3783 重叠均线与日序列相关 HAC(maxlags=61) 0.000234  -2.200103 0.027800 2010-06-04 2025-12-30       0.583796               1.000000            NaN
      7 estimated           None 日均策略减买入持有收益 -0.000499 3782 重叠均线与日序列相关 HAC(maxlags=61) 0.000232  -2.155263 0.031141 2010-06-07 2025-12-30       0.583796               1.000000            NaN
      8 estimated           None 日均策略减买入持有收益 -0.000502 3781 重叠均线与日序列相关 HAC(maxlags=61) 0.000227  -2.207743 0.027262 2010-06-08 2025-12-30       0.583796               1.000000            NaN
      9 estimated           None 日均策略减买入持有收益 -0.000560 3780 重叠均线与日序列相关 HAC(maxlags=61) 0.000224  -2.505455 0.012229 2010-06-09 2025-12-30       0.583796               0.611471            NaN
     10 estimated           None 日均策略减买入持有收益 -0.000572 3779 重叠均线与日序列相关 HAC(maxlags=61) 0.000222  -2.574323 0.010044 2010-06-10 2025-12-30       0.583796               0.502182            NaN
     11 estimated           None 日均策略减买入持有收益 -0.000621 3778 重叠均线与日序列相关 HAC(maxlags=61) 0.000232  -2.674733 0.007479 2010-06-11 2025-12-30       0.583796               0.373944            NaN
     12 estimated           None 日均策略减买入持有收益 -0.000570 3777 重叠均线与日序列相关 HAC(maxlags=61) 0.000235  -2.424173 0.015343 2010-06-17 2025-12-30       0.583796               0.767164            NaN
     13 estimated           None 日均策略减买入持有收益 -0.000467 3776 重叠均线与日序列相关 HAC(maxlags=61) 0.000234  -1.993651 0.046190 2010-06-18 2025-12-30       0.711912               1.000000            NaN
     14 estimated           None 日均策略减买入持有收益 -0.000515 3775 重叠均线与日序列相关 HAC(maxlags=61) 0.000237  -2.174691 0.029653 2010-06-21 2025-12-30       0.583796               1.000000            NaN
     15 estimated           None 日均策略减买入持有收益 -0.000407 3774 重叠均线与日序列相关 HAC(maxlags=61) 0.000230  -1.765796 0.077430 2010-06-22 2025-12-30       0.805895               1.000000            NaN
     16 estimated           None 日均策略减买入持有收益 -0.000399 3773 重叠均线与日序列相关 HAC(maxlags=61) 0.000231  -1.730323 0.083573 2010-06-23 2025-12-30       0.805895               1.000000            NaN
     17 estimated           None 日均策略减买入持有收益 -0.000359 3772 重叠均线与日序列相关 HAC(maxlags=61) 0.000228  -1.575394 0.115165 2010-06-24 2025-12-30       0.901427               1.000000            NaN
     18 estimated           None 日均策略减买入持有收益 -0.000395 3771 重叠均线与日序列相关 HAC(maxlags=61) 0.000230  -1.717010 0.085977 2010-06-25 2025-12-30       0.805895               1.000000            NaN
     19 estimated           None 日均策略减买入持有收益 -0.000394 3770 重叠均线与日序列相关 HAC(maxlags=61) 0.000224  -1.757624 0.078811 2010-06-28 2025-12-30       0.805895               1.000000            NaN
     20 estimated           None 日均策略减买入持有收益 -0.000489 3769 重叠均线与日序列相关 HAC(maxlags=61) 0.000224  -2.183268 0.029016 2010-06-29 2025-12-30       0.583796               1.000000            NaN
     21 estimated           None 日均策略减买入持有收益 -0.000519 3768 重叠均线与日序列相关 HAC(maxlags=61) 0.000226  -2.290467 0.021994 2010-06-30 2025-12-30       0.583796               1.000000            NaN
     22 estimated           None 日均策略减买入持有收益 -0.000442 3767 重叠均线与日序列相关 HAC(maxlags=61) 0.000225  -1.960388 0.049951 2010-07-01 2025-12-30       0.711912               1.000000            NaN
     23 estimated           None 日均策略减买入持有收益 -0.000363 3766 重叠均线与日序列相关 HAC(maxlags=61) 0.000225  -1.612085 0.106944 2010-07-02 2025-12-30       0.891038               1.000000            NaN
     24 estimated           None 日均策略减买入持有收益 -0.000313 3765 重叠均线与日序列相关 HAC(maxlags=61) 0.000220  -1.421240 0.155247 2010-07-05 2025-12-30       0.957727               1.000000            NaN
     25 estimated           None 日均策略减买入持有收益 -0.000288 3764 重叠均线与日序列相关 HAC(maxlags=61) 0.000217  -1.323716 0.185597 2010-07-06 2025-12-30       1.000000               1.000000            NaN
     26 estimated           None 日均策略减买入持有收益 -0.000315 3763 重叠均线与日序列相关 HAC(maxlags=61) 0.000213  -1.478490 0.139277 2010-07-07 2025-12-30       0.956411               1.000000            NaN
     27 estimated           None 日均策略减买入持有收益 -0.000251 3762 重叠均线与日序列相关 HAC(maxlags=61) 0.000211  -1.187929 0.234861 2010-07-08 2025-12-30       1.000000               1.000000            NaN
     28 estimated           None 日均策略减买入持有收益 -0.000190 3761 重叠均线与日序列相关 HAC(maxlags=61) 0.000209  -0.911470 0.362048 2010-07-09 2025-12-30       1.000000               1.000000            NaN
     29 estimated           None 日均策略减买入持有收益 -0.000261 3760 重叠均线与日序列相关 HAC(maxlags=61) 0.000202  -1.291889 0.196396 2010-07-12 2025-12-30       1.000000               1.000000            NaN
     30 estimated           None 日均策略减买入持有收益 -0.000220 3759 重叠均线与日序列相关 HAC(maxlags=61) 0.000200  -1.101420 0.270714 2010-07-13 2025-12-30       1.000000               1.000000            NaN
     31 estimated           None 日均策略减买入持有收益 -0.000263 3758 重叠均线与日序列相关 HAC(maxlags=61) 0.000194  -1.352457 0.176229 2010-07-14 2025-12-30       0.995202               1.000000            NaN
     32 estimated           None 日均策略减买入持有收益 -0.000274 3757 重叠均线与日序列相关 HAC(maxlags=61) 0.000193  -1.416858 0.156524 2010-07-15 2025-12-30       0.957727               1.000000            NaN
     33 estimated           None 日均策略减买入持有收益 -0.000345 3756 重叠均线与日序列相关 HAC(maxlags=61) 0.000198  -1.744059 0.081149 2010-07-16 2025-12-30       0.805895               1.000000            NaN
     34 estimated           None 日均策略减买入持有收益 -0.000267 3755 重叠均线与日序列相关 HAC(maxlags=61) 0.000197  -1.350187 0.176956 2010-07-19 2025-12-30       0.995202               1.000000            NaN
     35 estimated           None 日均策略减买入持有收益 -0.000299 3754 重叠均线与日序列相关 HAC(maxlags=61) 0.000194  -1.542210 0.123022 2010-07-20 2025-12-30       0.922506               1.000000            NaN
     36 estimated           None 日均策略减买入持有收益 -0.000272 3753 重叠均线与日序列相关 HAC(maxlags=61) 0.000192  -1.412089 0.157924 2010-07-21 2025-12-30       0.957727               1.000000            NaN
     37 estimated           None 日均策略减买入持有收益 -0.000254 3752 重叠均线与日序列相关 HAC(maxlags=61) 0.000192  -1.317505 0.187669 2010-07-22 2025-12-30       1.000000               1.000000            NaN
     38 estimated           None 日均策略减买入持有收益 -0.000235 3751 重叠均线与日序列相关 HAC(maxlags=61) 0.000192  -1.226480 0.220018 2010-07-23 2025-12-30       1.000000               1.000000            NaN
     39 estimated           None 日均策略减买入持有收益 -0.000289 3750 重叠均线与日序列相关 HAC(maxlags=61) 0.000196  -1.474681 0.140298 2010-07-26 2025-12-30       0.956411               1.000000            NaN
     40 estimated           None 日均策略减买入持有收益 -0.000248 3749 重叠均线与日序列相关 HAC(maxlags=61) 0.000196  -1.268641 0.204569 2010-07-27 2025-12-30       1.000000               1.000000            NaN
     41 estimated           None 日均策略减买入持有收益 -0.000297 3748 重叠均线与日序列相关 HAC(maxlags=61) 0.000197  -1.510380 0.130946 2010-07-28 2025-12-30       0.950250               1.000000            NaN
     42 estimated           None 日均策略减买入持有收益 -0.000346 3747 重叠均线与日序列相关 HAC(maxlags=61) 0.000201  -1.722035 0.085063 2010-07-29 2025-12-30       0.805895               1.000000            NaN
     43 estimated           None 日均策略减买入持有收益 -0.000361 3746 重叠均线与日序列相关 HAC(maxlags=61) 0.000200  -1.802977 0.071392 2010-07-30 2025-12-30       0.805895               1.000000            NaN
     44 estimated           None 日均策略减买入持有收益 -0.000388 3745 重叠均线与日序列相关 HAC(maxlags=61) 0.000199  -1.954570 0.050634 2010-08-02 2025-12-30       0.711912               1.000000            NaN
     45 estimated           None 日均策略减买入持有收益 -0.000431 3744 重叠均线与日序列相关 HAC(maxlags=61) 0.000198  -2.174442 0.029672 2010-08-03 2025-12-30       0.583796               1.000000            NaN
     46 estimated           None 日均策略减买入持有收益 -0.000395 3743 重叠均线与日序列相关 HAC(maxlags=61) 0.000197  -2.004734 0.044991 2010-08-04 2025-12-30       0.711912               1.000000            NaN
     47 estimated           None 日均策略减买入持有收益 -0.000360 3742 重叠均线与日序列相关 HAC(maxlags=61) 0.000200  -1.798996 0.072019 2010-08-05 2025-12-30       0.805895               1.000000            NaN
     48 estimated           None 日均策略减买入持有收益 -0.000340 3741 重叠均线与日序列相关 HAC(maxlags=61) 0.000200  -1.697411 0.089619 2010-08-06 2025-12-30       0.806429               1.000000            NaN
     49 estimated           None 日均策略减买入持有收益 -0.000327 3740 重叠均线与日序列相关 HAC(maxlags=61) 0.000199  -1.646557 0.099649 2010-08-09 2025-12-30       0.862195               1.000000            NaN
     50 estimated           None 日均策略减买入持有收益 -0.000317 3739 重叠均线与日序列相关 HAC(maxlags=61) 0.000202  -1.570906 0.116204 2010-08-10 2025-12-30       0.901427               1.000000            NaN
     51 estimated           None 日均策略减买入持有收益 -0.000289 3738 重叠均线与日序列相关 HAC(maxlags=61) 0.000204  -1.411706 0.158036 2010-08-11 2025-12-30       0.957727               1.000000            NaN
     52 estimated           None 日均策略减买入持有收益 -0.000284 3737 重叠均线与日序列相关 HAC(maxlags=61) 0.000203  -1.399117 0.161778 2010-08-12 2025-12-30       0.957727               1.000000            NaN
     53 estimated           None 日均策略减买入持有收益 -0.000256 3736 重叠均线与日序列相关 HAC(maxlags=61) 0.000203  -1.262257 0.206856 2010-08-13 2025-12-30       1.000000               1.000000            NaN
     54 estimated           None 日均策略减买入持有收益 -0.000221 3735 重叠均线与日序列相关 HAC(maxlags=61) 0.000206  -1.076033 0.281912 2010-08-16 2025-12-30       1.000000               1.000000            NaN
  1. 财务指标选股的多重检验*

本题使用真实上市公司财务数据估计多个财务因子与未来收益的秩相关，并对双向聚类 p 值实施 BY 校正。

列表 13.13: 习题7：加载财务报表与股价数据 
import pandas as pd  # 数据分析库
import numpy as np  # 数值计算库
import os  # 文件系统操作
import hashlib  # 为事前计划生成稳定内容指纹
import json  # 稳定序列化预注册计划
import statsmodels.api as sm  # 在 fresh kernel 中显式提供 HAC 回归
from statsmodels.stats.multitest import multipletests  # 在 fresh kernel 中显式提供家族校正

# 构建跨平台数据路径
# 根据操作系统设置数据根目录路径
from pathlib import Path  # 使用跨平台路径对象解析显式数据根
DATA_DIR = Path(os.environ['BOOK_DATA_DIR']).expanduser().resolve()  # 从必需环境变量取得数据根
if not DATA_DIR.is_dir():  # 在读取前验证数据根
    raise FileNotFoundError(f'BOOK_DATA_DIR 不存在或不是目录: {DATA_DIR}')  # 失败即停止，不回退固定路径
m13_factor_financial_slice = Path(os.environ['BOOK_M13_FACTOR_FINANCIAL_SLICE']).expanduser().resolve()  # 运行器冻结的年报列/年份切片
m13_factor_price_slice = Path(os.environ['BOOK_M13_FACTOR_PRICE_SLICE']).expanduser().resolve()  # 冻结年末价切片
factor_slice_contract = [(m13_factor_financial_slice, 'BOOK_M13_FACTOR_FINANCIAL_SLICE_SHA256'), (m13_factor_price_slice, 'BOOK_M13_FACTOR_PRICE_SLICE_SHA256')]  # 两个切片合同
for slice_path, hash_name in factor_slice_contract:  # 逐切片流式核验字节
    factor_digest = hashlib.sha256()  # 初始化当前切片摘要
    with slice_path.open('rb') as source_file:  # 打开字节流
        for block in iter(lambda: source_file.read(1024 * 1024), b''):  # 固定块读取到EOF
            factor_digest.update(block)  # 累加摘要
    if factor_digest.hexdigest() != os.environ[hash_name].lower():
        raise RuntimeError({'status': 'stopped', 'reason': 'factor_slice_sha256_mismatch', 'path': str(slice_path)})
financial_statement_raw = pd.read_hdf(m13_factor_financial_slice)  # 哈希通过后只读取冻结年报切片
annual_financial_data = financial_statement_raw[financial_statement_raw['quarter'].str.endswith('q4')].copy()  # 只保留年报
annual_financial_data['report_year'] = annual_financial_data['quarter'].str[:4].astype(int)  # 提取报告年
annual_financial_data['available_year'] = annual_financial_data['report_year'] + 1  # 保守形成年
annual_financial_data['info_date'] = pd.to_datetime(annual_financial_data['info_date'])  # 解析披露/修订日
annual_financial_data['formation_date'] = pd.to_datetime(annual_financial_data['available_year'].astype(str) + '-12-31')  # 年末形成日
annual_financial_data = annual_financial_data[annual_financial_data['info_date'].le(annual_financial_data['formation_date'])].copy()  # 排除未来重述
annual_financial_data = annual_financial_data.sort_values(['order_book_id', 'available_year', 'info_date']).groupby(['order_book_id', 'available_year']).tail(1)  # 逐形成日取最后可得版本

完成年报数据筛选和年份提取后，接下来基于财务报表原始字段计算五个核心财务因子（ROE、净利润率、资产负债率、总资产周转率和流动比率），作为后续相关性检验的解释变量。

列表 13.14: 习题7：计算财务因子 
# 计算五个核心财务因子
# ROE = 净利润 / 归属母公司股东权益
annual_financial_data['ROE'] = annual_financial_data['net_profit'] / annual_financial_data['equity_parent_company']
# 净利润率 = 净利润 / 营业收入
annual_financial_data['Net_Margin'] = annual_financial_data['net_profit'] / annual_financial_data['operating_revenue']
# 资产负债率 = (流动负债+非流动负债) / 总资产
annual_financial_data['Debt_Ratio'] = (
    annual_financial_data['current_liabilities'].fillna(0) + annual_financial_data['non_current_liabilities'].fillna(0)
) / annual_financial_data['total_assets']  # 填充缺失值后计算负债率
# 总资产周转率 = 营业收入 / 总资产
annual_financial_data['Asset_Turnover'] = annual_financial_data['operating_revenue'] / annual_financial_data['total_assets']
# 流动比率 = 流动资产 / 流动负债
annual_financial_data['Current_Ratio'] = annual_financial_data['current_assets'] / annual_financial_data['current_liabilities']

财务因子计算完毕后，下一步从后复权股价数据中提取每年末收盘价、计算下一年度收益率，并与财务因子表进行合并，构建因子-收益率配对数据集。

列表 13.15: 习题7：计算年度收益率并合并 
# 读取后复权股价数据
stock_price_annual = pd.read_hdf(m13_factor_price_slice)  # 只读取冻结年末价切片
# 提取每年最后一个交易日收盘价（用于计算年度收益率）
stock_price_annual = stock_price_annual.reset_index()  # 重置DataFrame索引
stock_price_annual['date'] = pd.to_datetime(stock_price_annual['date'], errors='coerce')  # 严格解析日期后才能定义期末
stock_price_annual['close'] = pd.to_numeric(stock_price_annual['close'], errors='coerce')  # 价格必须可验证为有限数
invalid_price_rows = stock_price_annual['date'].isna() | ~np.isfinite(stock_price_annual['close']) | stock_price_annual['close'].le(0)
if invalid_price_rows.any() or stock_price_annual.duplicated(['order_book_id', 'date']).any():
    raise RuntimeError({'status': 'stopped', 'reason': 'factor_price_date_key_or_close_invalid', 'invalid_rows': int(invalid_price_rows.sum()), 'duplicate_keys': int(stock_price_annual.duplicated(['order_book_id', 'date']).sum())})
stock_price_annual = stock_price_annual.sort_values(['order_book_id', 'date'], kind='mergesort').reset_index(drop=True)  # 物理HDF行序不是统计合同
stock_price_annual['year'] = stock_price_annual['date'].dt.year  # 提取交易年份
year_end_dates = stock_price_annual.groupby(['order_book_id', 'year'])['date'].transform('max')  # 显式求每年最大交易日
year_end_price = stock_price_annual.loc[stock_price_annual['date'].eq(year_end_dates), ['order_book_id', 'year', 'date', 'close']].rename(columns={'date': 'year_end_date'}).copy()
if year_end_price.duplicated(['order_book_id', 'year']).any():
    raise RuntimeError({'status': 'stopped', 'reason': 'factor_year_end_not_unique'})
year_end_price = year_end_price.sort_values(['order_book_id', 'year']).reset_index(drop=True)  # 显式固定shift的公司内年序
# 严格在公司组内计算从当年末到下一年末的收益，禁止跨公司 shift
year_end_price['next_year'] = year_end_price.groupby('order_book_id')['year'].shift(-1)
year_end_price['next_year_return'] = (year_end_price.groupby('order_book_id')['close'].shift(-1) / year_end_price['close'] - 1).where(year_end_price['next_year'].eq(year_end_price['year'] + 1))  # 缺年不得跨期冒充一年收益
year_end_price['next_year_end_date'] = year_end_price.groupby('order_book_id')['year_end_date'].shift(-1)  # 保留收益窗口的真实期末日
year_end_selection_bytes = year_end_price[['order_book_id', 'year', 'year_end_date', 'close']].to_csv(index=False, date_format='%Y-%m-%d', lineterminator='\n').encode('utf-8')
print({'year_end_selection_sha256': hashlib.sha256(year_end_selection_bytes).hexdigest(), 'selected_rows': len(year_end_price), 'selection_rule': 'unique security-date; stable sort; explicit max date'})
# 将财务因子与下一年收益率合并（基于股票代码和年份）
factor_return_merged = annual_financial_data.merge(
    year_end_price[['order_book_id', 'year', 'year_end_date', 'next_year_end_date', 'next_year_return']],  # 收益率与真实起止日一起交付
    left_on=['order_book_id', 'available_year'],  # 至少滞后到下一年后再启动收益窗
    right_on=['order_book_id', 'year'],  # 右表合并键
    how='inner'  # 内连接仅保留匹配记录
)  # 合并财务因子与未来收益率数据

得到因子—收益配对后，把两个变量转换为秩并标准化；标准化秩回归的斜率就是 Spearman 相关 estimand。推断使用公司与年度双向聚类协方差，以同时允许公司内持续性和同年度共同冲击。五个因子共享样本且相关，故默认用 BY，而不是把朴素 p 值送入 BH。

列表 13.16: 习题7：稳定两维簇键与数值验收 
def stable_cluster_codes(frame):
    company_values = frame['order_book_id'].astype(str)  # 不把字符串直接交给statsmodels协方差实现
    company_levels = sorted(company_values.unique())
    year_levels = sorted(frame['year'].unique())
    return np.column_stack([pd.Categorical(company_values, categories=company_levels).codes, pd.Categorical(frame['year'], categories=year_levels).codes])
def validated_cluster_result(model):
    estimate, standard_error = model.params.iloc[1], model.bse.iloc[1]
    covariance_values = np.asarray(model.cov_params(), dtype=float)
    statistic, raw_p_value = model.tvalues.iloc[1], model.pvalues.iloc[1]
    valid = np.isfinite(estimate) and np.isfinite(covariance_values).all() and np.isfinite(standard_error) and standard_error > 0 and np.isfinite(statistic) and np.isfinite(raw_p_value) and 0 <= raw_p_value <= 1
    return estimate, standard_error, statistic, raw_p_value, valid
factor_names_to_test = ['ROE', 'Net_Margin', 'Debt_Ratio', 'Asset_Turnover', 'Current_Ratio']  # 冻结五项检验族
factor_test_results, factor_skip_log = [], []  # 分别保留有效检验与跳过原因
列表 13.17: 习题7：秩相关 estimand、公司/年度双向聚类推断与 BY 校正 
for factor_name in factor_names_to_test:  # 逐项执行同一预注册程序
    required_columns = ['order_book_id', 'year', factor_name, 'next_year_return']  # 保留 estimand 与聚类键
    factor_pair = factor_return_merged[required_columns].replace([np.inf, -np.inf], np.nan).dropna().copy()  # 删除无效配对
    company_count, year_count = factor_pair['order_book_id'].nunique(), factor_pair['year'].nunique()  # 审计两维簇数
    if len(factor_pair) < 100 or company_count < 20 or year_count < 20:  # 普通渐近双聚类要求两维均有足够簇
        factor_skip_log.append({'factor': factor_name, 'n': len(factor_pair), 'companies': company_count, 'years': year_count, 'reason': '普通渐近双向聚类要求 n>=100、公司>=20、年度>=20；少簇需另行预注册 wild cluster bootstrap，本例停止'})
        continue  # 不以朴素检验替代失败的稳健检验
    ranked_factor = factor_pair[factor_name].rank(method='average')  # 把因子转换为平均秩
    ranked_return = factor_pair['next_year_return'].rank(method='average')  # 把未来收益转换为平均秩
    if ranked_factor.nunique() < 2 or ranked_return.nunique() < 2 or not np.isfinite(ranked_factor).all() or not np.isfinite(ranked_return).all():
        factor_skip_log.append({'factor': factor_name, 'n': len(factor_pair), 'companies': company_count, 'years': year_count, 'reason': '因子秩与收益秩必须有限且非常数'})
        continue
    standardized_factor = (ranked_factor - ranked_factor.mean()) / ranked_factor.std(ddof=0)  # 标准化使斜率等于秩相关
    standardized_return = (ranked_return - ranked_return.mean()) / ranked_return.std(ddof=0)  # 标准化被解释变量
    cluster_groups = stable_cluster_codes(factor_pair)  # 稳定整数编码避免object簇键失败
    try:
        rank_model = sm.OLS(standardized_return, sm.add_constant(standardized_factor)).fit(cov_type='cluster', cov_kwds={'groups': cluster_groups})  # 双向聚类拟合
        estimate, standard_error, statistic, raw_p_value, inference_valid = validated_cluster_result(rank_model)
    except Exception as error:
        factor_skip_log.append({'factor': factor_name, 'n': len(factor_pair), 'companies': company_count, 'years': year_count, 'reason': f'双向聚类拟合失败:{type(error).__name__}'})
        continue
    if not inference_valid:
        factor_skip_log.append({'factor': factor_name, 'n': len(factor_pair), 'companies': company_count, 'years': year_count, 'reason': '效应/协方差/统计量/p值须有限，SE须有限且>0，p须在[0,1]'})
        continue
    factor_test_results.append({'factor': factor_name, 'estimand': 'Spearman rho（标准化秩回归斜率）', 'estimate': estimate, 'n': len(factor_pair), 'companies': company_count, 'years': year_count, 'dependence': '公司内与年度共同冲击', 'se_method': 'company/year two-way cluster', 'se': standard_error, 'statistic': statistic, 'p_value': raw_p_value})  # 保存完整审计字段
factor_result_columns = ['factor', 'estimand', 'estimate', 'n', 'companies', 'years', 'dependence', 'se_method', 'se', 'statistic', 'p_value']  # 固定空表schema
factor_results_df = pd.DataFrame(factor_test_results, columns=factor_result_columns)  # 汇总所有可估检验
factor_family_df = pd.DataFrame({'factor': factor_names_to_test})  # 权威五项计划家族
factor_family_df = factor_family_df.merge(factor_results_df, on='factor', how='left', validate='one_to_one')  # 空结果也必须得到完整五行schema
factor_skip_df = pd.DataFrame(factor_skip_log, columns=['factor', 'n', 'companies', 'years', 'reason']).rename(columns={'reason': 'failure_reason'})  # 定型失败表
if not factor_skip_df.empty:
    factor_family_df = factor_family_df.merge(factor_skip_df[['factor', 'failure_reason']], on='factor', how='left', validate='one_to_one')  # 保留唯一失败原因
else:
    factor_family_df['failure_reason'] = None
factor_family_df['status'] = np.where(factor_family_df['p_value'].notna(), 'estimated', 'skipped')  # 标记状态
factor_family_df['p_value'] = factor_family_df['p_value'].where(factor_family_df['status'].eq('estimated'), 1.0)  # 跳过项以1保留m=5
factor_correction_p = factor_family_df['p_value']
factor_family_df['by_adjusted_p'] = multipletests(factor_correction_p, alpha=0.05, method='fdr_by')[1]  # 固定五项BY
factor_family_df['by_reject'] = factor_family_df['status'].eq('estimated') & factor_family_df['by_adjusted_p'].le(0.05)  # 跳过不可拒绝
print({'planned': 5, 'estimated': int(factor_family_df['status'].eq('estimated').sum()), 'skipped': int(factor_family_df['status'].eq('skipped').sum())})  # 审计分母
assert len(factor_family_df) == int(factor_family_df['status'].eq('estimated').sum()) + int(factor_family_df['status'].eq('skipped').sum()) == 5
print(factor_family_df.to_string(index=False))  # 输出完整五行
{'planned': 5, 'estimated': 0, 'skipped': 5}
        factor estimand estimate   n companies years dependence se_method  se statistic p_value failure_reason  status by_adjusted_p  by_reject
           ROE      NaN      NaN NaN       NaN   NaN        NaN       NaN NaN       NaN     1.0           None skipped             1      False
    Net_Margin      NaN      NaN NaN       NaN   NaN        NaN       NaN NaN       NaN     1.0           None skipped             1      False
    Debt_Ratio      NaN      NaN NaN       NaN   NaN        NaN       NaN NaN       NaN     1.0           None skipped             1      False
Asset_Turnover      NaN      NaN NaN       NaN   NaN        NaN       NaN NaN       NaN     1.0           None skipped             1      False
 Current_Ratio      NaN      NaN NaN       NaN   NaN        NaN       NaN NaN       NaN     1.0           None skipped             1      False

上方表逐项记录 estimand、有效 \(n\)、公司/年度簇数、依赖结构和 SE 方法。相关方向、效应量和拒绝数必须从当前输出读取；即使 BY 后拒绝零假设，也不等于具备样本外选股能力，还需独立测试期、基准收益与交易成本审计。

  1. 行业动量策略的多重检验*

这是一个典型的 “Data Mining” 场景：对多个行业 × 多个回看窗口的动量策略进行检验，演示多重检验校正如何过滤虚假发现。

# 步骤1：加载数据并构建行业月度收益率
import pandas as pd  # 数据分析库
import numpy as np  # 数值计算库
import os  # 文件系统操作
import hashlib  # 核验M12血缘、数据切片与计划字节
import json  # 规范化计划和发现集
import statsmodels.api as sm  # 本题独立执行HAC回归
from statsmodels.stats.multitest import multipletests  # BY与Holm校正
from scipy.stats import norm  # 单尾稳健p值
from pathlib import Path  # 使用跨平台路径对象解析显式数据根
列表 13.18: M13 在签发运行身份前验证唯一可写目录与发布依赖 
import sys  # 记录解释器版本以复核运行环境
import sklearn  # 记录统计学习库版本
import time  # 记录端到端墙钟耗时
import tracemalloc  # 记录本进程 Python 内存峰值
import uuid  # 为追加式 telemetry 生成唯一运行标识
import atexit  # 进程异常结束时补写唯一stopped终态
from datetime import datetime  # 解析运行器冻结的运行时刻
from zoneinfo import ZoneInfo  # 将运行日志统一到中国标准时间
from scripts.m13_preflight import atomic_write_same, publish_terminal, terminal_record, validate_artifact_directory  # 复用可独立测试的原子发布实现
m13_started_perf = time.perf_counter()  # 必须早于入口哈希和数据读取
tracemalloc.start()  # 从 fresh-kernel 设置块开始计量峰值
m13_failure_artifact_value = os.environ.get('BOOK_M13_ARTIFACT_DIR', '').strip()  # 先读取唯一授权目录
if not m13_failure_artifact_value:  # 没有可发布位置时绝不签发run_id
    raise RuntimeError({'status': 'stopped', 'reason': 'required_environment_missing', 'name': 'BOOK_M13_ARTIFACT_DIR'})
m13_artifact_dir = validate_artifact_directory(m13_failure_artifact_value)  # 真实创建、刷新并清理权限探针
m13_runtime_state = {'run_id': None, 'terminal': False}  # 守卫注册前运行身份保持未签发
列表 13.19: M13 终态 schema、运行标识、状态与字节哈希核验 
def m13_validate_terminal(terminal_bytes, expected_sha256=None):  # 对候选和回读终态执行同一验证
    expected_run_id = m13_runtime_state['run_id']  # 只消费守卫持有的唯一身份
    if expected_run_id is None:
        raise RuntimeError({'status': 'terminal_unchanged', 'reason': 'run_id_not_issued'})  # 未签发时禁止形成终态
    try:
        validated_record = terminal_record(terminal_bytes, expected_run_id)  # 核验规范JSON、schema、身份与终态
    except ValueError as error:
        raise RuntimeError({'status': 'terminal_unchanged', 'reason': str(error), 'run_id': expected_run_id}) from error
    actual_sha256 = hashlib.sha256(terminal_bytes).hexdigest()  # 核验持久化字节而非内存对象
    if expected_sha256 and actual_sha256 != expected_sha256:  # 发布后回读必须完全相同
        raise RuntimeError({'status': 'terminal_unchanged', 'reason': 'persisted_terminal_hash_mismatch', 'run_id': expected_run_id})
    return validated_record, actual_sha256  # 向调用者返回已验证状态与摘要
列表 13.20: M13 稳定与终态候选共享同目录临时文件政策 
def m13_stable_payload(record):  # 为所有JSON稳定制品统一字节边界
    return (json.dumps(record, ensure_ascii=False, sort_keys=True, separators=(',', ':')) + '\n').encode('utf-8')  # 固定键序与单换行
m13_temporary_output_policy = '.<authorized-final-name>.<uuid>.tmp; flush+fsync; non-overwriting link; always remove; never submit'  # 与manifest逐项对齐
列表 13.21: M13 完整落盘后原子且排他发布单一终态 
def m13_write_terminal(record, fail_after_flush=False):  # completed、stopped与退出守卫共用入口
    expected_run_id = m13_runtime_state['run_id']  # 读取守卫持有身份
    if expected_run_id is None:
        raise RuntimeError({'status': 'terminal_unchanged', 'reason': 'run_id_not_issued'})  # 未签发不得发布
    if record.get('run_id') != expected_run_id:
        raise RuntimeError({'status': 'terminal_unchanged', 'reason': 'terminal_run_id_mismatch', 'run_id': expected_run_id})  # 候选不得串用身份
    terminal_bytes = (json.dumps(record, ensure_ascii=False, sort_keys=True, separators=(',', ':')) + '\n').encode('utf-8')  # 固定候选字节
    _, terminal_sha256 = m13_validate_terminal(terminal_bytes)  # 写前拒绝非法候选
    if fail_after_flush:
        raise OSError('m13_fault_after_terminal_flush')  # 故障探针由可再分发preflight覆盖
    try:
        terminal_path = publish_terminal(m13_artifact_dir, record)  # 完整刷新后原子且排他发布
    except ValueError as error:
        raise RuntimeError({'status': 'terminal_unchanged', 'reason': str(error), 'run_id': expected_run_id}) from error
    persisted_bytes = terminal_path.read_bytes()  # 回读最终路径而非信任写入缓冲
    m13_validate_terminal(persisted_bytes, terminal_sha256)  # 核验持久化schema、ID、状态与哈希
    m13_runtime_state['terminal'] = True  # 只有回读成功才关闭退出守卫
    return terminal_path, terminal_bytes
列表 13.22: M13 异常退出时的唯一 stopped 终态守卫 
def m13_finalize_incomplete_run():  # fresh-kernel退出前补全未落盘终态
    expected_run_id = m13_runtime_state['run_id']  # 守卫只依赖注册前已有闭包状态
    if expected_run_id is None or m13_runtime_state['terminal']:
        return  # 未签发或已严格核验终态时无需动作
    terminal_path = m13_artifact_dir / f'm13-run-{expected_run_id}-v5.json'  # 定位当前ID唯一终态
    try:
        if terminal_path.exists():  # 正常完成或显式停止后只做严格回读验证
            m13_validate_terminal(terminal_path.read_bytes())  # 截断与非法schema不得静默通过
            m13_runtime_state['terminal'] = True  # 合法先到终态关闭守卫
            return  # 已验证终态不得被退出守卫重写
        record = {'schema_version': '5', 'run_id': expected_run_id, 'status': 'stopped', 'reason': 'process_exit_without_terminal', 'finished_at': datetime.now(ZoneInfo('Asia/Shanghai')).isoformat()}
        m13_write_terminal(record)  # 与正常路径共用完整落盘和排他发布
    except Exception as terminal_error:
        print({'m13_exit_guard': str(terminal_error)}, file=sys.stderr)  # 明示截断或竞争但不覆盖既有字节
atexit.register(m13_finalize_incomplete_run)
列表 13.23: M13 在发布器和退出守卫就绪后签发身份并定义结构化停止 
m13_run_id = os.environ.get('BOOK_RUN_ID', '').strip() or uuid.uuid4().hex  # 运行器可注入ID，否则签发唯一ID
m13_runtime_state['run_id'] = m13_run_id  # 守卫已注册且目录已验证后才公开身份
def m13_stop(reason, **details):  # 每条失败路径只尝试写一份stopped终态
    failure_record = {'schema_version': '5', 'run_id': m13_run_id, 'status': 'stopped', 'reason': reason, 'finished_at': datetime.now(ZoneInfo('Asia/Shanghai')).isoformat(), 'elapsed_seconds': round(time.perf_counter() - m13_started_perf, 6), 'peak_memory_bytes': int(tracemalloc.get_traced_memory()[1]), 'loaded_rows': details.pop('loaded_rows', {}), 'details': details}
    m13_write_terminal(failure_record)  # 已有终态时保留原文件且报告ID复用
    raise RuntimeError(failure_record)
def m13_require_env(name):  # 在任何直接索引前统一形成结构化缺参证据
    value = os.environ.get(name, '').strip()
    if not value:
        m13_stop('required_environment_missing', name=name)
    return value
列表 13.24: M13 输入路径与流式哈希准备 
def m13_stream_sha256(path, block_size=1024 * 1024):  # 流式计算大文件指纹
    digest = hashlib.sha256()  # 初始化摘要
    with path.open('rb') as source_file:  # 打开真实字节流
        for block in iter(lambda: source_file.read(block_size), b''):  # 固定块读取
            digest.update(block)  # 更新摘要
    return digest.hexdigest()  # 返回完整指纹
DATA_DIR = Path(m13_require_env('BOOK_DATA_DIR')).expanduser().resolve()  # 从必需环境变量取得数据根
if not DATA_DIR.is_dir():  # 在读取前验证数据根
    m13_stop('BOOK_DATA_DIR_not_directory', path=str(DATA_DIR))  # 失败即停止，不回退固定路径
path_price = DATA_DIR / 'stock' / 'stock_price_post_adjusted.h5'  # 后复权股价路径
path_basic = DATA_DIR / 'stock' / 'stock_basic_data.h5'  # 上市公司基本信息路径
m13_price_slice = Path(m13_require_env('BOOK_M13_PRICE_SLICE')).expanduser().resolve()  # 冻结日期/列/公司切片
m13_basic_slice = Path(m13_require_env('BOOK_M13_BASIC_SLICE')).expanduser().resolve()  # 冻结行业映射切片
m12_contract_path = Path(m13_require_env('BOOK_M12_HANDOFF')).expanduser().resolve()  # 上章稳定机器合同
m12_run_log_path = Path(m13_require_env('BOOK_M12_RUN_LOG')).expanduser().resolve()  # 上章某次成功 telemetry
列表 13.25: M13 稳定序列化、幂等写入与确定性入口 hash 核验 
def m13_stable_csv(frame):  # 固定账本与映射表字节
    return frame.to_csv(index=False, date_format='%Y-%m-%d', lineterminator='\n').encode('utf-8')  # 固定编码和换行
def m13_write_same(path, payload):  # 稳定制品也必须完整刷新后原子发布
    try:
        persisted = atomic_write_same(path, payload)  # 同目录临时名、flush、fsync与非覆盖链接
    except ValueError as error:
        m13_stop('artifact_changed', path=str(path), error=str(error))  # 同名不同字节统一停止
    if hashlib.sha256(persisted).digest() != hashlib.sha256(payload).digest():
        m13_stop('persisted_artifact_hash_mismatch', path=str(path))  # 只信最终路径回读字节
    return persisted  # 后续hash可绑定已落盘权威字节
m13_expected_hashes = {'price_sha256': m13_require_env('BOOK_STOCK_PRICE_SHA256').lower(), 'basic_sha256': m13_require_env('BOOK_STOCK_BASIC_SHA256').lower(), 'price_slice_sha256': m13_require_env('BOOK_M13_PRICE_SLICE_SHA256').lower(), 'basic_slice_sha256': m13_require_env('BOOK_M13_BASIC_SLICE_SHA256').lower(), 'm12_handoff_sha256': m13_require_env('BOOK_M12_HANDOFF_SHA256').lower()}  # 运行前冻结确定性入口字节
m13_paths = {'price': path_price, 'basic': path_basic, 'price_slice': m13_price_slice, 'basic_slice': m13_basic_slice, 'm12_handoff': m12_contract_path}  # 建立入口血缘路径
if not all(path.is_file() for path in m13_paths.values()):  # 所有入口必须存在
    m13_stop('m13_input_missing', paths={name: str(path) for name, path in m13_paths.items() if not path.is_file()})
m13_actual_hashes = {name + '_sha256': m13_stream_sha256(path) for name, path in m13_paths.items()}  # 流式复算实际字节
if m13_actual_hashes != m13_expected_hashes:  # 比较完整expected/actual映射
    m13_stop('m13_input_sha256_mismatch', expected=m13_expected_hashes, actual=m13_actual_hashes)
列表 13.26: M12 独立交付 telemetry 的 SHA-256 与合同核验 
m12_run_log_expected_sha256 = m13_require_env('BOOK_M12_RUN_LOG_SHA256').lower()  # telemetry 独立核验而不进入稳定计划
if not m12_run_log_path.is_file() or m13_stream_sha256(m12_run_log_path) != m12_run_log_expected_sha256:
    m13_stop('m12_run_log_missing_or_hash_mismatch', path=str(m12_run_log_path))
m12_handoff = json.loads(m12_contract_path.read_text(encoding='utf-8'))  # 读取已核验的M12边界
if m12_handoff.get('schema_version') != '4' or m12_handoff.get('freeze_date') != m13_require_env('BOOK_DATA_FREEZE_DATE'):  # 只接受稳定多层roster合同
    m13_stop('m12_handoff_schema_or_freeze_mismatch')
m12_run_log = json.loads(m12_run_log_path.read_text(encoding='utf-8'))
if m12_run_log.get('schema_version') != '4' or m12_run_log.get('status') != 'completed' or m12_run_log.get('preregistration_sha256') != m12_handoff.get('preregistration_sha256'):
    m13_stop('m12_run_log_schema_or_contract_mismatch')
if m12_run_log.get('handoff_sha256') != m13_actual_hashes['m12_handoff_sha256']:
    m13_stop('m12_run_log_handoff_sha256_mismatch', expected=m13_actual_hashes['m12_handoff_sha256'], actual=m12_run_log.get('handoff_sha256'))
列表 13.27: M13 可再分发三边界、哈希与唯一终态探针 
from scripts.m13_preflight import run_self_test  # 载入纯标准库课堂preflight
m13_preflight_result = run_self_test()  # 执行四类边界、hash、schema和终态并发探针
assert m13_preflight_result['matching_handoff'] == 'completed'
assert m13_preflight_result['mismatched_telemetry'] == 'stopped'
assert m13_preflight_result['tampered_handoff'] == 'stopped'
assert m13_preflight_result['schema_mismatch'] == 'stopped'
assert m13_preflight_result['pre_start_invariant'] and m13_preflight_result['warm_up_invariant']
assert m13_preflight_result['first_eligible_sensitive'] and m13_preflight_result['post_end_invariant']
assert m13_preflight_result['fault_safe'] and m13_preflight_result['truncated_detected']
assert m13_preflight_result['exit_guard_atomic'] and m13_preflight_result['idempotent']
assert m13_preflight_result['concurrent_single_terminal'] and m13_preflight_result['consecutive_setup_exit']
assert m13_preflight_result['stable_artifact_atomic'] and m13_preflight_result['temporary_candidates_clean']
assert m13_preflight_result['calendar_adjacent'] and m13_preflight_result['internal_missing_month_skipped']
assert m13_preflight_result['replication_origin_inclusive'] and m13_preflight_result['inference_reconciled']
assert m13_preflight_result['period_end_selection_correct'] and m13_preflight_result['period_end_permutation_invariant']
assert m13_preflight_result['annual_period_end_permutation_invariant']
assert m13_preflight_result['invalid_price_rows_rejected'] and m13_preflight_result['stable_two_way_cluster_codes']
assert m13_preflight_result['output_contract_exact']
assert m13_preflight_result['no_mount_fixture_status']['milestone_status'] == 'blocked_pending_licensed_snapshot'
print(m13_preflight_result)  # 输出不含外部数据的课堂验收证据
{'matching_handoff': 'completed', 'mismatched_telemetry': 'stopped', 'tampered_handoff': 'stopped', 'schema_mismatch': 'stopped', 'pre_start_invariant': True, 'warm_up_invariant': True, 'first_eligible_sensitive': True, 'post_end_invariant': True, 'stable_artifact_atomic': True, 'fault_safe': True, 'idempotent': True, 'truncated_detected': True, 'exit_guard_atomic': True, 'consecutive_setup_exit': True, 'concurrent_single_terminal': True, 'temporary_candidates_clean': True, 'calendar_adjacent': True, 'internal_missing_month_skipped': True, 'replication_origin_inclusive': True, 'period_end_selection_correct': True, 'period_end_permutation_invariant': True, 'annual_period_end_permutation_invariant': True, 'invalid_price_rows_rejected': True, 'stable_two_way_cluster_codes': True, 'inference_reconciled': True, 'reconciliation': {'planned': 4, 'estimated': 1, 'skipped': 3, 'correction_p_values': [0.0227501319, 1.0, 1.0, 1.0]}, 'discoveries': ['A'], 'environment_inputs': 21, 'required_outputs': 11, 'output_contract_exact': True, 'no_mount_fixture_status': {'launcher_terminal': 'completed', 'activity': 'contract_fixture_only', 'empirical_outputs_created': False, 'milestone_status': 'blocked_pending_licensed_snapshot', 'grading_scope': 'process_contract_evidence_only'}, 'no_mount_student_outputs': ['m13-process-evidence-v1.json', 'm13-process-terminal-v1.json'], 'no_mount_gradebook_contract': {'schema_version': '1', 'record_name': 'm13-no-mount-grade-v1.json', 'status': 'completed_process_substitute', 'points_earned': 20, 'points_possible': 20, 'course_weight_percent': 4, 'empirical_milestone_status': 'blocked_pending_licensed_snapshot', 'replacement_policy': 'For a cohort assigned the institutional no-mount route, this verified 20/20 process substitute is the final M13 4% grade; a later licensed run is formative evidence only and neither adds to nor reduces this score.'}}
列表 13.28: M13 逐字节复核 M12 完整证据包 
m12_artifact_root = m12_contract_path.parent  # 从已核验handoff解析同目录制品
m12_expected_output_names = {'loadings', 'scores', 'silhouettes', 'cluster-profile', 'stability-seeds', 'lineage'}  # 与M12生产者完全一致
if set(m12_handoff.get('output_names', [])) != m12_expected_output_names or set(m12_handoff.get('outputs', {})) != m12_expected_output_names:
    m13_stop('m12_output_schema_mismatch', expected=sorted(m12_expected_output_names), actual=sorted(m12_handoff.get('outputs', {})))
m12_artifact_paths = {'preregistration': m12_artifact_root / 'm12-preregistration-v2.json', 'universe': m12_artifact_root / m12_handoff['rosters']['universe']['file'], 'analysis_subset': m12_artifact_root / m12_handoff['rosters']['analysis_subset']['file'], 'clustering': m12_artifact_root / m12_handoff['rosters']['clustering']['file'], 'exclusions': m12_artifact_root / m12_handoff['exclusions']['file']}  # 只定位稳定制品
m12_artifact_expected = {'preregistration': m12_handoff['preregistration_sha256'], 'universe': m12_handoff['rosters']['universe']['sha256'], 'analysis_subset': m12_handoff['rosters']['analysis_subset']['sha256'], 'clustering': m12_handoff['rosters']['clustering']['sha256'], 'exclusions': m12_handoff['exclusions']['sha256']}  # 恢复稳定指纹
for output_name, output_hash in m12_handoff['outputs'].items():  # 纳入全部六份M12正式输出
    m12_artifact_paths[f'output_{output_name}'] = m12_artifact_root / f'm12-{output_name}-v2.csv'  # 从固定接口定位输出
    m12_artifact_expected[f'output_{output_name}'] = output_hash  # 恢复handoff登记指纹
if not all(path.is_file() for path in m12_artifact_paths.values()):  # 禁止缺制品时继续
    m13_stop('m12_bound_artifact_missing')  # 完整证据包缺失即停止
m12_artifact_actual = {name: m13_stream_sha256(path) for name, path in m12_artifact_paths.items()}  # 复算全部M12制品字节
if m12_artifact_actual != m12_artifact_expected:  # 逐一比较handoff绑定值
    m13_stop('m12_bound_artifact_mismatch', expected=m12_artifact_expected, actual=m12_artifact_actual)  # 任一篡改即停止
if m12_run_log.get('output_hashes') != m12_handoff['outputs']:
    m13_stop('m12_run_log_output_hash_mismatch')
列表 13.29: M13 重核 M12 三层名册与行业映射 
universe_roster = pd.read_csv(m12_artifact_paths['universe'], dtype={'order_book_id': str})  # M13权威总体只取课程级多行业roster
analysis_subset_roster = pd.read_csv(m12_artifact_paths['analysis_subset'], dtype={'order_book_id': str})  # 读取但不作为M13总体
clustering_roster = pd.read_csv(m12_artifact_paths['clustering'], dtype={'order_book_id': str})  # 读取但不作为M13总体
if not set(clustering_roster['order_book_id']).issubset(analysis_subset_roster['order_book_id']):
    m13_stop('m12_clustering_roster_outside_analysis_subset')
if not set(analysis_subset_roster['order_book_id']).issubset(universe_roster['order_book_id']):
    m13_stop('m12_analysis_subset_outside_universe')
industry_mapping = universe_roster[['order_book_id', 'industry_name']].copy()  # 直接消费M12课程级行业映射
mapping_columns = ['order_book_id', 'industry_name', 'industry_mapping_version']  # 恢复映射字节列序
if hashlib.sha256(m13_stable_csv(universe_roster[mapping_columns])).hexdigest() != m12_handoff['industry_mapping']['sha256']:  # 核对映射真实字节
    m13_stop('industry_mapping_hash_mismatch')  # 映射漂移即停止
列表 13.30: M13 价格与行业切片的模式及总体审计 
stock_price_raw = pd.read_hdf(m13_price_slice).reset_index()  # 只读运行器冻结的小切片
stock_basic_info = pd.read_hdf(m13_basic_slice).reset_index(drop=True)  # 只读冻结行业映射切片
required_price_columns = {'order_book_id', 'date', 'close'}  # 冻结价量schema
required_basic_columns = {'order_book_id', 'citics_2019_l1_name'}  # 冻结行业schema
if not required_price_columns.issubset(stock_price_raw) or not required_basic_columns.issubset(stock_basic_info):  # 核验字段完整性
    m13_stop('m13_slice_schema_missing')  # 缺字段即停止
stock_price_raw['order_book_id'] = stock_price_raw['order_book_id'].astype(str)  # 统一公司键类型
stock_basic_info['order_book_id'] = stock_basic_info['order_book_id'].astype(str)  # 统一映射键类型
stock_price_raw['date'] = pd.to_datetime(stock_price_raw['date'], errors='coerce')  # 严格日期是期末选择的前提
stock_price_raw['close'] = pd.to_numeric(stock_price_raw['close'], errors='coerce')  # 将不可解析价格归为数据边界失败
invalid_price_rows = stock_price_raw['date'].isna() | ~np.isfinite(stock_price_raw['close']) | stock_price_raw['close'].le(0)
duplicate_price_keys = stock_price_raw.duplicated(['order_book_id', 'date'], keep=False)
if invalid_price_rows.any() or duplicate_price_keys.any():
    m13_stop('price_date_key_or_close_invalid', invalid_rows=int(invalid_price_rows.sum()), duplicate_keys=int(duplicate_price_keys.sum()))
stock_price_raw = stock_price_raw.sort_values(['order_book_id', 'date'], kind='mergesort').reset_index(drop=True)  # 任意HDF行序不得影响期末价
universe_ids = set(universe_roster['order_book_id'])  # 冻结课程级公司集合
if not set(stock_price_raw['order_book_id']).issubset(universe_ids):  # 禁止总体外公司混入
    m13_stop('price_slice_outside_universe')  # 总体外键即停止
basic_mapping = stock_basic_info[['order_book_id', 'citics_2019_l1_name']].rename(columns={'citics_2019_l1_name': 'industry_name'}).sort_values('order_book_id').reset_index(drop=True)  # 定型切片行业映射
if not basic_mapping.equals(industry_mapping.sort_values('order_book_id').reset_index(drop=True)):  # 要求切片映射与M12总体逐行相同
    m13_stop('basic_mapping_not_equal_to_m12_universe')  # 增删改均失败
列表 13.31: 习题8：在查看 p 值前冻结探索与复现计划 
m13_replication_origin = pd.Timestamp(m13_require_env('BOOK_M13_REPLICATION_ORIGIN'))  # 读取事前登记的复现起点
m13_exploration_start = pd.Timestamp(m13_require_env('BOOK_M13_EXPLORATION_START'))  # 读取事前登记的探索起点
m13_replication_end_exclusive = pd.Timestamp(m13_require_env('BOOK_M13_REPLICATION_END_EXCLUSIVE'))  # 读取复现开区间上界
m13_artifact_dir = Path(m13_require_env('BOOK_M13_ARTIFACT_DIR')).expanduser().resolve()  # 定位课程授权制品目录
m13_artifact_dir.mkdir(parents=True, exist_ok=True)  # 仅创建授权的项目制品目录
if m13_require_env('BOOK_M13_INDUSTRY_VERSION') != m12_handoff['industry_mapping']['version']:  # M13不得另换行业版本
    m13_stop('industry_mapping_version_mismatch')  # 版本不一致即停止
if not m13_exploration_start < m13_replication_origin < m13_replication_end_exclusive:  # 三边界必须严格有序
    m13_stop('invalid_period_boundary')  # 禁止重叠、空探索期或空复现期
period_starts = [m13_exploration_start, m13_replication_origin, m13_replication_end_exclusive]  # 三边界必须落在月初
if any(value != value.to_period('M').start_time for value in period_starts):
    m13_stop('period_boundary_not_month_start')  # 避免月内日期被静默归并
m13_warm_up_policy = 'no pre-start prices or outcomes; first exploration month establishes baseline close and has no return'
m13_first_eligible_month = (m13_exploration_start.to_period('M') + 1).start_time  # 暖启动后首个可评价结果月
m13_period_boundaries = {'exploration_start_inclusive': m13_exploration_start.date().isoformat(), 'first_eligible_month': m13_first_eligible_month.date().isoformat(), 'exploration_end_exclusive': m13_replication_origin.date().isoformat(), 'replication_start_inclusive': m13_replication_origin.date().isoformat(), 'replication_end_exclusive': m13_replication_end_exclusive.date().isoformat()}  # 形成三边界与暖启动合同
m13_design = {'schema_version': '3', 'analysis_unit': 'industry-month equal-weight simple return in decimal units', 'formation_date_rule': 'signal at calendar month t final trading date sums returns from t-window+1 through t', 'label_window': 'calendar month t+1 industry return; all boundaries and reported dates use label_month', 'price_key_rule': 'parseable unique (order_book_id,date); close finite and strictly positive; stable sort by order_book_id,date', 'period_end_selection_rule': 'explicit row at maximum trading date per security and calendar month; persisted in m13-period-end-closes-v1.csv', 'security_month_rule': 'complete security-month calendar; return requires finite closes in adjacent formation_month and label_month', 'estimand': 'Pearson rho via standardized slope with intercept', 'alternative': 'positive', 'lookback_windows': list(range(1, 13)), **m13_period_boundaries, 'period_boundaries': m13_period_boundaries, 'warm_up_months': 1, 'warm_up_policy': m13_warm_up_policy, 'security_coverage_rule': 'observed exploration returns after warm-up >=80%', 'security_minimum_coverage': 0.8, 'industry_eligibility_rule': 'at least 2 preregistration-eligible securities', 'industry_minimum_securities': 2, 'replication_eligibility_rule': 'replication missingness never changes industry or planned family', 'exploration_minimum_rule': 'max(36,4*(window+1))', 'replication_minimum_rule': 'max(24,4*(window+1))', 'hac_rule': 'maxlags=max(1,window)', 'skip_rule': 'skip if n is below minimum, signal/label is constant, effect/covariance/statistic is nonfinite, SE is nonfinite or <=0, or p is outside [0,1]; retain member with p=1', 'primary_correction': 'BY', 'replication_correction': 'Holm over frozen discoveries', 'alpha': 0.05, 'industry_mapping_version': m12_handoff['industry_mapping']['version'], 'input_hashes': m13_actual_hashes, 'm12_handoff_sha256': m13_actual_hashes['m12_handoff_sha256'], 'm12_preregistration_sha256': m12_handoff['preregistration_sha256'], 'm12_universe_sha256': m12_handoff['rosters']['universe']['sha256']}  # 保留完整日历、期末选择、标签边界与数值失败合同
m13_design_bytes = (json.dumps(m13_design, ensure_ascii=False, sort_keys=True, separators=(',', ':')) + '\n').encode('utf-8')  # 规范化基础设计
m13_design_sha256 = hashlib.sha256(m13_design_bytes).hexdigest()  # 计算基础设计指纹
m13_write_same(m13_artifact_dir / 'm13-design-v3.json', m13_design_bytes)  # 在构造任何统计量前幂等冻结

加载股价和基本信息数据后，接下来将两者合并以获取每只股票的行业归属，并在此基础上计算各行业的等权月度收益率序列。

# 步骤2：合并行业信息，计算行业月度收益率
merged_data = pd.merge(  # 合并股价与行业分类数据
    stock_price_raw[['order_book_id', 'date', 'close']],  # 选取核心字段
    industry_mapping,
    on='order_book_id', how='inner'  # 按股票代码内连接
)
merged_data['date'] = pd.to_datetime(merged_data['date'])  # 确保日期为datetime类型
merged_data = merged_data[merged_data['date'].ge(m13_exploration_start) & merged_data['date'].lt(m13_replication_end_exclusive)].copy()  # 价格严格服从总区间
if merged_data['date'].lt(m13_exploration_start).any():  # 禁止任何起点前价格进入暖启动
    m13_stop('pre_start_price_leakage')
if merged_data['date'].ge(m13_replication_end_exclusive).any():  # 禁止开区间上界及以后价格进入分析
    m13_stop('post_replication_end_price_leakage')
merged_data['year_month'] = merged_data['date'].dt.to_period('M')  # 提取年月标识
列表 13.32: M13 与输入行序无关的月末价选择与审计表 
# 先显式选择月内最大交易日，再补齐每只证券的完整日历月
merged_data = merged_data.sort_values(['order_book_id', 'date'], kind='mergesort').reset_index(drop=True)  # 稳定行序便于审计
month_end_dates = merged_data.groupby(['order_book_id', 'year_month'])['date'].transform('max')  # 不使用当前行序的last
observed_monthly_close = merged_data.loc[merged_data['date'].eq(month_end_dates), ['order_book_id', 'year_month', 'date', 'close']].rename(columns={'date': 'month_end_date'}).copy()
if observed_monthly_close.duplicated(['order_book_id', 'year_month']).any():
    m13_stop('month_end_selection_not_unique')  # 日键已唯一时期末行也必须唯一
period_end_audit = observed_monthly_close.rename(columns={'year_month': 'label_month'}).sort_values(['order_book_id', 'label_month']).reset_index(drop=True)
period_end_bytes = m13_stable_csv(period_end_audit[['order_book_id', 'label_month', 'month_end_date', 'close']])
period_end_sha256 = hashlib.sha256(period_end_bytes).hexdigest()
m13_write_same(m13_artifact_dir / 'm13-period-end-closes-v1.csv', period_end_bytes)  # 持久化每只证券实际选中的月末日
列表 13.33: M13 完整证券—日历月收益面板 
calendar_months = pd.period_range(m13_exploration_start.to_period('M'), m13_replication_end_exclusive.to_period('M') - 1, freq='M')  # 固定含基线月的完整月历
security_month_index = pd.MultiIndex.from_product([sorted(universe_ids), calendar_months], names=['order_book_id', 'label_month'])  # 建立完整证券—标签月笛卡尔积
monthly_close = period_end_audit.set_index(['order_book_id', 'label_month']).reindex(security_month_index).reset_index()  # 内部缺月显式保留NA
monthly_close['industry_name'] = monthly_close['order_book_id'].map(industry_mapping.set_index('order_book_id')['industry_name'])  # 从冻结总体恢复行业
monthly_close['formation_month'] = monthly_close['label_month'] - 1  # 每项收益明确保留相邻形成月
previous_close = monthly_close.groupby('order_book_id', sort=False)['close'].shift(1)  # 只取日历紧邻前月，绝不跳过NA
monthly_close['monthly_ret'] = (monthly_close['close'] / previous_close - 1).replace([np.inf, -np.inf], np.nan)  # 两个相邻月收盘价均有限才有收益
if observed_monthly_close.empty:  # 空价格切片不能形成资格账本
    m13_stop('empty_monthly_close')  # 明确停止而非生成空家族
if not (monthly_close['label_month'].astype(int) - monthly_close['formation_month'].astype(int)).eq(1).all():
    m13_stop('nonadjacent_security_month_return')  # 形成月与标签月必须严格相邻

在获取个股月度收益率后，将数据按行业聚合为等权月度收益率矩阵，每一列代表一个行业的时间序列，为后续动量策略检验提供标准化输入。

# 步骤3：按冻结规则建立逐证券排除账本
exploration_origin_month = m13_exploration_start.to_period('M')  # 转换探索起始月份
replication_origin_month = m13_replication_origin.to_period('M')  # 转换复现起始月份
replication_end_month = m13_replication_end_exclusive.to_period('M')  # 转换复现开区间上界月份
expected_exploration_months = len(pd.period_range(exploration_origin_month, replication_origin_month - 1, freq='M'))  # 冻结覆盖率分母
expected_exploration_return_months = expected_exploration_months - m13_design['warm_up_months']  # 首月只形成基准价
if expected_exploration_return_months < 1:  # 暖启动后分母必须为正
    m13_stop('empty_expected_exploration_months')  # 空探索窗即停止
exploration_security_rows = monthly_close[monthly_close['label_month'].ge(m13_first_eligible_month.to_period('M')) & monthly_close['label_month'].lt(replication_origin_month)]  # 覆盖率边界只作用于收益标签月
security_observations = exploration_security_rows.groupby('order_book_id')['monthly_ret'].count()  # 统计探索期有效收益月数
m13_exclusion_ledger = universe_roster[['order_book_id', 'industry_name', 'industry_mapping_version']].copy()  # 从完整课程总体开始审计
m13_exclusion_ledger['design_sha256'] = m13_design_sha256  # 把资格账本绑定统计前设计合同
m13_exclusion_ledger['exploration_n'] = m13_exclusion_ledger['order_book_id'].map(security_observations).fillna(0).astype(int)  # 不让无价格公司静默消失
m13_exclusion_ledger['exploration_coverage'] = m13_exclusion_ledger['exploration_n'] / expected_exploration_return_months  # 排除首个暖启动月
coverage_eligible = m13_exclusion_ledger['exploration_coverage'].ge(m13_design['security_minimum_coverage'])  # 执行预注册证券覆盖门槛
eligible_industry_counts = m13_exclusion_ledger[coverage_eligible].groupby('industry_name')['order_book_id'].nunique()  # 统计每行业预合格证券
validated_industries = sorted(eligible_industry_counts[eligible_industry_counts.ge(m13_design['industry_minimum_securities'])].index)  # 冻结满足证券数门槛的多行业集合
m13_exclusion_ledger['reason'] = np.select([m13_exclusion_ledger['exploration_n'].eq(0), ~coverage_eligible, ~m13_exclusion_ledger['industry_name'].isin(validated_industries)], ['no_exploration_return', 'exploration_coverage_below_80pct', 'industry_has_fewer_than_2_eligible_securities'], default='included')  # 按优先级给唯一原因
m13_exclusion_ledger = m13_exclusion_ledger.sort_values('order_book_id').reset_index(drop=True)  # 固定账本行序
m13_exclusion_bytes = m13_stable_csv(m13_exclusion_ledger)  # 固定完整资格账本字节
m13_exclusion_sha256 = hashlib.sha256(m13_exclusion_bytes).hexdigest()  # 计算排除证据指纹
m13_write_same(m13_artifact_dir / 'm13-universe-exclusion-ledger-v3.csv', m13_exclusion_bytes)  # 在检验前冻结严格排除账本
列表 13.34: M13 从排除账本生成行业月度收益矩阵 
eligible_security_ids = set(m13_exclusion_ledger.loc[m13_exclusion_ledger['reason'].eq('included'), 'order_book_id'])  # 恢复严格纳入证券
eligible_monthly_close = monthly_close[monthly_close['order_book_id'].isin(eligible_security_ids)].copy()  # 仅在账本冻结后清洗总体
industry_monthly_returns = eligible_monthly_close.groupby(['label_month', 'industry_name'])['monthly_ret'].mean().unstack()  # 按明确标签月计算行业等权简单收益
complete_month_index = pd.period_range(exploration_origin_month, replication_end_month - 1, freq='M')  # 日历键严格止于复现开区间上界前
industry_monthly_returns = industry_monthly_returns.reindex(index=complete_month_index)  # 缺月保留NA而非跳到下一观察月
industry_monthly_returns.index = industry_monthly_returns.index.to_timestamp()  # 将月份转换为时间戳
industry_monthly_returns = industry_monthly_returns.sort_index().reindex(columns=validated_industries)  # 固定全时期顺序与行业全集
if len(validated_industries) < 2:  # 正式M13必须是多行业家族
    m13_stop('fewer_than_two_validated_industries')  # 禁止退回M12单行业子集
print({'validated_industries': validated_industries, 'universe_n': len(universe_roster), 'eligible_security_n': len(eligible_security_ids), 'exclusion_ledger_sha256': m13_exclusion_sha256})  # 输出资格审计
列表 13.35: 习题8：建立不重叠探索期与复现期 
exploration_mask = industry_monthly_returns.index.to_series().between(m13_exploration_start, m13_replication_origin, inclusive='left')
exploration_raw = industry_monthly_returns.loc[exploration_mask.to_numpy()].copy()  # 同时执行探索下界和上界
replication_mask = industry_monthly_returns.index.to_series().between(m13_replication_origin, m13_replication_end_exclusive, inclusive='left')
replication_raw = industry_monthly_returns.loc[replication_mask.to_numpy()].copy()  # 同时执行复现下界和开区间上界
exploration_returns = exploration_raw.reindex(columns=validated_industries).copy()  # 探索只用账本冻结行业
replication_returns = replication_raw.reindex(columns=validated_industries).copy()  # 复制缺失不得删除行业
if exploration_returns.loc[m13_exploration_start].notna().any():  # 探索起点月只能提供基准收盘价
    m13_stop('warm_up_outcome_present')  # 禁止基线月结果进入发现统计量
if exploration_returns.index.min() < m13_exploration_start:  # 起点前结果不得进入检验或发现集
    m13_stop('exploration_lower_bound_violated')
if exploration_returns.empty or replication_returns.empty:  # 两个时期都必须存在
    m13_stop('empty_exploration_or_replication_period')  # 禁止用同一时期替代
if exploration_returns.index.max() >= replication_returns.index.min():  # 明确验证时期不重叠
    m13_stop('period_overlap')  # 边界冲突即停止
if replication_returns.index.max() >= m13_replication_end_exclusive:  # 复现结果不得越过声明上界
    m13_stop('replication_upper_bound_violated')
列表 13.36: M13 冻结多行业乘十二窗口的完整计划族 
planned_rows = [{'industry': industry_name, 'window': lookback_months, 'strategy': f'{industry_name}_M{lookback_months}', 'analysis_unit': m13_design['analysis_unit'], 'formation_date_rule': m13_design['formation_date_rule'], 'label_window': m13_design['label_window'], 'exploration_minimum_n': max(36, 4 * (lookback_months + 1)), 'replication_minimum_n': max(24, 4 * (lookback_months + 1)), 'hac_maxlags': max(1, lookback_months), 'skip_rule': m13_design['skip_rule'], 'alternative': m13_design['alternative']} for industry_name in validated_industries for lookback_months in m13_design['lookback_windows']]  # 为每个家族成员冻结全部推断字段
for planned_member in planned_rows:  # 为逐项追溯生成独立指纹
    planned_member['member_sha256'] = hashlib.sha256(json.dumps(planned_member, ensure_ascii=False, sort_keys=True, separators=(',', ':')).encode('utf-8')).hexdigest()  # 哈希该成员全部字段
if len(planned_rows) != len(validated_industries) * 12:  # planned family必须是多行业乘12
    m13_stop('planned_family_not_industries_times_12')  # 禁止缺成员
m13_plan = {**m13_design, 'design_sha256': m13_design_sha256, 'industries': validated_industries, 'planned_family': planned_rows, 'period_end_closes_sha256': period_end_sha256, 'period_end_selected_rows': len(period_end_audit), 'exclusion_ledger_sha256': m13_exclusion_sha256, 'exploration_key_sha256': hashlib.sha256(exploration_returns.to_csv(lineterminator='\n').encode()).hexdigest(), 'replication_key_sha256': hashlib.sha256(replication_returns.index.to_series().to_csv(index=False, lineterminator='\n').encode()).hexdigest()}  # 绑定期末日、总体、账本、两期键与完整成员
m13_plan_bytes = (json.dumps(m13_plan, ensure_ascii=False, sort_keys=True, separators=(',', ':')) + '\n').encode()  # 规范化计划
m13_plan_sha256 = hashlib.sha256(m13_plan_bytes).hexdigest()  # 计算完整计划指纹
m13_plan_path = m13_artifact_dir / 'm13-preregistration-v3.json'  # 固定版本文件名
m13_write_same(m13_plan_path, m13_plan_bytes)  # 幂等冻结完整计划族
print({'preregistration_sha256': m13_plan_sha256, 'industries': validated_industries, 'planned_family_n': len(planned_rows), 'exploration': [exploration_returns.index.min(), exploration_returns.index.max()], 'replication': [replication_returns.index.min(), replication_returns.index.max()]})  # 输出时期审计

行业月度收益率矩阵准备就绪后，对每个“行业 × 回看窗口”组合估计标准化回归斜率；它在含截距的一元回归中等于 Pearson 相关系数。HAC 最大滞后至少覆盖该回看窗口，以处理重叠信号诱发的序列相关。

列表 13.37: 习题8：行业动量相关 estimand 与 HAC 推断 
momentum_test_results, skipped_combinations = [], []  # 保存完整结果与跳过日志
for planned_member in m13_plan['planned_family']:  # 只遍历已冻结的逐项合同
    industry_name, lookback_months = planned_member['industry'], int(planned_member['window'])  # 恢复行业与窗口
    momentum_signal = industry_monthly_returns[industry_name].rolling(lookback_months).sum()  # 在形成月t末汇总截至t的收益
    future_return = industry_monthly_returns[industry_name].shift(-1)  # 只对齐日历紧邻的t+1标签收益
    aligned_df = pd.concat({'signal': momentum_signal, 'future_return': future_return}, axis=1)  # 先保留完整形成月历
    aligned_df['formation_month'] = aligned_df.index  # 明确保留信号形成月
    aligned_df['label_month'] = aligned_df['formation_month'] + pd.offsets.MonthBegin(1)  # 标签必须恰为下一日历月
    aligned_df = aligned_df[aligned_df['label_month'].ge(m13_first_eligible_month) & aligned_df['label_month'].lt(m13_replication_origin)].dropna(subset=['signal', 'future_return'])  # 探索边界只作用于标签月
    minimum_months = int(planned_member['exploration_minimum_n'])  # 从逐成员合同读取探索门槛
    if len(aligned_df) < minimum_months or aligned_df.nunique().min() < 2:  # 执行冻结skip规则
        skipped_combinations.append({'industry': industry_name, 'window': lookback_months, 'member_sha256': planned_member['member_sha256'], 'n': len(aligned_df), 'reason': f'需要 n>={minimum_months} 且两变量非恒定'})  # 保留失败成员
        continue  # 禁止退回其他窗口或朴素检验
    standardized_signal = (aligned_df['signal'] - aligned_df['signal'].mean()) / aligned_df['signal'].std(ddof=0)  # 按冻结口径标准化信号
    standardized_return = (aligned_df['future_return'] - aligned_df['future_return'].mean()) / aligned_df['future_return'].std(ddof=0)  # 按冻结口径标准化标签
    hac_lags = int(planned_member['hac_maxlags'])  # 从逐成员合同读取HAC带宽
    momentum_model = sm.OLS(standardized_return, sm.add_constant(standardized_signal)).fit(cov_type='HAC', cov_kwds={'maxlags': hac_lags})  # 估计冻结estimand
    estimate, standard_error = momentum_model.params.iloc[1], momentum_model.bse.iloc[1]  # 提取效应与HAC标准误
    covariance_values = np.asarray(momentum_model.cov_params(), dtype=float)  # 取得完整HAC协方差矩阵
    test_statistic = estimate / standard_error  # 显式形成冻结单尾统计量
    one_sided_p = norm.sf(test_statistic)  # 执行冻结正向单尾备择
    inference_valid = np.isfinite(estimate) and np.isfinite(covariance_values).all() and np.isfinite(standard_error) and standard_error > 0 and np.isfinite(test_statistic) and np.isfinite(one_sided_p) and 0 <= one_sided_p <= 1  # 所有推断量必须有限且在定义域内
    if not inference_valid:
        skipped_combinations.append({'industry': industry_name, 'window': lookback_months, 'member_sha256': planned_member['member_sha256'], 'n': len(aligned_df), 'reason': '效应/协方差/统计量/p值须有限，SE须有限且>0，p须在[0,1]'})  # 数值失败保留在完整家族
        continue  # 失败项不得进入estimated计数
    momentum_test_results.append({'strategy': planned_member['strategy'], 'member_sha256': planned_member['member_sha256'], 'estimand': m13_plan['estimand'], 'estimate': estimate, 'n': len(aligned_df), 'dependence': f'{lookback_months}月重叠信号与月度序列相关', 'se_method': f'HAC(maxlags={hac_lags})', 'se': standard_error, 'statistic': test_statistic, 'p_value': one_sided_p, 'formation_start': aligned_df['formation_month'].min(), 'formation_end': aligned_df['formation_month'].max(), 'label_start': aligned_df['label_month'].min(), 'label_end': aligned_df['label_month'].max()})  # 报告边界只引用标签月并保留形成月

在收集稳健 p 值后，行业之间共享市场冲击、窗口之间高度嵌套，因此采用任意依赖下有效的 BY。失败检验仍保留在计划清单和跳过日志中。

列表 13.38: 习题8：依赖稳健 p 值的 BY 校正与跳过审计 
momentum_results_df = pd.DataFrame(momentum_test_results, columns=['strategy', 'member_sha256', 'estimand', 'estimate', 'n', 'dependence', 'se_method', 'se', 'statistic', 'p_value', 'formation_start', 'formation_end', 'label_start', 'label_end'])  # 汇总成功估计的策略并固定空表schema
momentum_family_df = pd.DataFrame(planned_rows).merge(momentum_results_df, on=['strategy', 'member_sha256'], how='left', validate='one_to_one')  # 以逐成员计划合同为权威左表
skip_df = pd.DataFrame(skipped_combinations).rename(columns={'reason': 'failure_reason'})  # 统一跳过原因字段
if not skip_df.empty:  # 将失败状态写回完整家族
    momentum_family_df = momentum_family_df.merge(skip_df[['industry', 'window', 'member_sha256', 'failure_reason']], on=['industry', 'window', 'member_sha256'], how='left', validate='one_to_one')  # 保留计划成员与失败追溯
else:
    momentum_family_df['failure_reason'] = None  # 即使无人跳过也固定完整家族schema
momentum_family_df['status'] = np.where(momentum_family_df['p_value'].notna(), 'estimated', 'skipped')  # 标记每项是否可估
momentum_family_df['p_value'] = momentum_family_df['p_value'].where(momentum_family_df['status'].eq('estimated'), 1.0)  # 失败项在完整家族中显式保存p=1
correction_p_values = momentum_family_df['p_value']  # 所有计划项已有合法校正输入
momentum_family_df['by_adjusted_p'] = multipletests(correction_p_values, alpha=m13_plan['alpha'], method='fdr_by')[1]  # 任意依赖BY主分析
momentum_family_df['by_reject'] = momentum_family_df['status'].eq('estimated') & momentum_family_df['by_adjusted_p'].le(m13_plan['alpha'])  # 只允许有效检验被拒绝
frozen_discoveries = momentum_family_df.loc[momentum_family_df['by_reject'], ['industry', 'window', 'strategy', 'member_sha256', 'replication_minimum_n', 'hac_maxlags', 'alternative']].copy()  # 原样冻结复制所需逐成员字段
print({'planned_tests': len(momentum_family_df), 'actual_tests': len(momentum_results_df), 'skipped_tests': len(skipped_combinations), 'by_rejections': len(frozen_discoveries)})  # 输出家族审计
assert len(momentum_family_df) == len(momentum_results_df) + len(skipped_combinations)  # planned=estimated+skipped必须逐次对账
print(momentum_family_df.to_string(index=False))  # 逐项显示状态、依赖和调整值
列表 13.39: M13 完整家族与 BY 结果的稳定序列化 
family_columns = ['industry', 'window', 'strategy', 'member_sha256', 'analysis_unit', 'formation_date_rule', 'label_window', 'exploration_minimum_n', 'replication_minimum_n', 'hac_maxlags', 'skip_rule', 'alternative', 'status', 'failure_reason', 'estimand', 'estimate', 'n', 'dependence', 'se_method', 'se', 'statistic', 'p_value', 'formation_start', 'formation_end', 'label_start', 'label_end']  # 固定形成月、标签月与推断证据列序
m13_family_artifact = momentum_family_df.reindex(columns=family_columns).copy()  # 保留计划、估计与失败但不混入校正列
m13_family_artifact['preregistration_sha256'] = m13_plan_sha256  # 把完整家族绑定预注册
m13_family_bytes = m13_stable_csv(m13_family_artifact)  # 固定完整家族字节
m13_write_same(m13_artifact_dir / 'm13-family-v3.csv', m13_family_bytes)  # 幂等写入完整家族
m13_by_artifact = momentum_family_df.copy()  # 保留BY调整值与拒绝状态
m13_by_artifact['preregistration_sha256'] = m13_plan_sha256  # 绑定同一完整计划
m13_by_bytes = m13_stable_csv(m13_by_artifact)  # 固定BY结果字节
m13_write_same(m13_artifact_dir / 'm13-by-results-v3.csv', m13_by_bytes)  # 幂等写入BY结果
列表 13.40: 习题8：在进入复现期前冻结探索发现集 
discovery_payload = {'preregistration_sha256': m13_plan_sha256, 'candidates': frozen_discoveries.to_dict(orient='records')}  # 绑定计划与探索发现
discovery_bytes = (json.dumps(discovery_payload, ensure_ascii=False, sort_keys=True, separators=(',', ':')) + '\n').encode('utf-8')  # 稳定序列化发现集
discovery_sha256 = hashlib.sha256(discovery_bytes).hexdigest()  # 计算不可变发现集指纹
discovery_path = m13_artifact_dir / 'm13-frozen-discoveries-v1.json'  # 使用固定路径防止选择性改名
m13_write_same(discovery_path, discovery_bytes)  # 冲突统一进入m13_stop且不改写原字节
print({'discovery_sha256': discovery_sha256, 'candidate_count': len(frozen_discoveries)})  # 输出复现入口证据
列表 13.41: 习题8：冻结候选在不重叠复现期的 HAC 检验 
replication_records = []  # 保存每个冻结候选的独立时期证据
for candidate in frozen_discoveries.to_dict(orient='records'):  # 不允许在复现期重新选择行业或窗口
    industry_name, lookback_months = candidate['industry'], int(candidate['window'])  # 读取冻结方向和窗口
    full_signal = industry_monthly_returns[industry_name].rolling(lookback_months).sum()  # 允许用起点前历史形成首个复现信号
    full_future = industry_monthly_returns[industry_name].shift(-1)  # 对齐下一月复现结果
    replication_pair = pd.concat({'signal': full_signal, 'future_return': full_future}, axis=1)  # 先保留完整形成月历
    replication_pair['formation_month'] = replication_pair.index  # 保留冻结信号形成月
    replication_pair['label_month'] = replication_pair['formation_month'] + pd.offsets.MonthBegin(1)  # 标签恰为下一日历月
    replication_pair = replication_pair[replication_pair['label_month'].ge(m13_replication_origin) & replication_pair['label_month'].lt(m13_replication_end_exclusive)].dropna(subset=['signal', 'future_return'])  # 复现边界只作用于标签月且含起点
    minimum_months = int(candidate['replication_minimum_n'])  # 从冻结成员读取复现最小样本
    if len(replication_pair) < minimum_months or replication_pair.nunique().min() < 2:  # 在拟合前执行失败守卫
        replication_records.append({'strategy': candidate['strategy'], 'member_sha256': candidate['member_sha256'], 'status': 'skipped', 'failure_reason': f'需要 n>={minimum_months} 且两变量非恒定', 'n': len(replication_pair), 'formation_start': replication_pair['formation_month'].min(), 'formation_end': replication_pair['formation_month'].max(), 'label_start': replication_pair['label_month'].min(), 'label_end': replication_pair['label_month'].max(), 'estimate': np.nan, 'se': np.nan, 'statistic': np.nan, 'p_value': 1.0})  # 失败候选仍以p=1留在完整复现家族
        continue  # 不用其他窗口或朴素检验替代
    replication_signal = (replication_pair['signal'] - replication_pair['signal'].mean()) / replication_pair['signal'].std(ddof=0)  # 使用冻结标准化 estimand
    replication_return = (replication_pair['future_return'] - replication_pair['future_return'].mean()) / replication_pair['future_return'].std(ddof=0)  # 标准化复现结果
    replication_model = sm.OLS(replication_return, sm.add_constant(replication_signal)).fit(cov_type='HAC', cov_kwds={'maxlags': int(candidate['hac_maxlags'])})  # 使用冻结成员HAC规则
    replication_estimate, replication_se = replication_model.params.iloc[1], replication_model.bse.iloc[1]  # 提取复现效应与标准误
    replication_covariance = np.asarray(replication_model.cov_params(), dtype=float)  # 读取完整HAC协方差
    replication_statistic = replication_estimate / replication_se  # 形成冻结单尾统计量
    replication_p = norm.sf(replication_statistic)  # 计算冻结正向单尾p值
    replication_valid = np.isfinite(replication_estimate) and np.isfinite(replication_covariance).all() and np.isfinite(replication_se) and replication_se > 0 and np.isfinite(replication_statistic) and np.isfinite(replication_p) and 0 <= replication_p <= 1  # 执行统一数值有效性合同
    if not replication_valid:
        replication_records.append({'strategy': candidate['strategy'], 'member_sha256': candidate['member_sha256'], 'status': 'skipped', 'failure_reason': '效应/协方差/统计量/p值须有限，SE须有限且>0，p须在[0,1]', 'n': len(replication_pair), 'formation_start': replication_pair['formation_month'].min(), 'formation_end': replication_pair['formation_month'].max(), 'label_start': replication_pair['label_month'].min(), 'label_end': replication_pair['label_month'].max(), 'estimate': np.nan, 'se': np.nan, 'statistic': np.nan, 'p_value': 1.0})  # 数值失败仍占Holm家族
        continue  # 不以无效统计量冒充复现结果
    replication_records.append({'strategy': candidate['strategy'], 'member_sha256': candidate['member_sha256'], 'status': 'estimated', 'failure_reason': None, 'n': len(replication_pair), 'formation_start': replication_pair['formation_month'].min(), 'formation_end': replication_pair['formation_month'].max(), 'label_start': replication_pair['label_month'].min(), 'label_end': replication_pair['label_month'].max(), 'estimate': replication_estimate, 'se': replication_se, 'statistic': replication_statistic, 'p_value': replication_p})  # 保存形成月、标签月与合法推断
replication_columns = ['strategy', 'member_sha256', 'status', 'failure_reason', 'n', 'formation_start', 'formation_end', 'label_start', 'label_end', 'estimate', 'se', 'statistic', 'p_value']  # 固定空表与非空表列序
replication_results_df = pd.DataFrame(replication_records, columns=replication_columns)  # 形成独立复现表
表 13.6: 习题8：冻结探索发现集的独立时期复现表 
if not replication_results_df.empty:  # 仅在探索期冻结了候选时评价复现
    replication_correction_p = replication_results_df['p_value']  # 跳过项已经显式以1保留在冻结复现家族
    replication_results_df['holm_adjusted_p'] = multipletests(replication_correction_p, alpha=0.05, method='holm')[1]  # 控制冻结发现集复现FWER
    replication_results_df['replicated'] = replication_results_df['status'].eq('estimated') & replication_results_df['estimate'].gt(0) & replication_results_df['holm_adjusted_p'].le(0.05)  # 同方向且校正后拒绝才记复现
else:  # 无探索发现也是可审计结果
    replication_results_df = pd.DataFrame(columns=replication_columns + ['holm_adjusted_p', 'replicated'])  # 输出空但定型的复现表
assert len(replication_results_df) == len(frozen_discoveries)  # 每个冻结候选必须恰有estimated或skipped记录
assert int(replication_results_df['status'].eq('estimated').sum()) + int(replication_results_df['status'].eq('skipped').sum()) == len(frozen_discoveries)  # planned=estimated+skipped
print({'preregistration_sha256': m13_plan_sha256, 'discovery_sha256': discovery_sha256, 'period_boundaries': m13_period_boundaries, 'replication_observed_end': replication_returns.index.max()})  # 绑定三边界和两份冻结制品
print(replication_results_df.to_string(index=False))  # 输出效应、SE、p值、状态与失败原因
列表 13.42: M13 复现与血缘稳定制品的哈希冻结 
replication_artifact = replication_results_df.copy()  # 不改写现场展示表
replication_artifact['preregistration_sha256'] = m13_plan_sha256  # 绑定复现所依据的完整计划
replication_artifact['discovery_sha256'] = discovery_sha256  # 绑定复现入口发现集
m13_replication_bytes = m13_stable_csv(replication_artifact)  # 固定复现与失败表字节
m13_write_same(m13_artifact_dir / 'm13-replication-v3.csv', m13_replication_bytes)  # 幂等冻结复现表
m13_lineage = {'schema_version': '3', 'm12_handoff_sha256': m13_actual_hashes['m12_handoff_sha256'], 'm12_bound_artifacts': m12_artifact_actual, 'input_hashes': m13_actual_hashes, 'price_key_rule': m13_design['price_key_rule'], 'period_end_selection_rule': m13_design['period_end_selection_rule'], 'period_end_closes_sha256': period_end_sha256, 'period_end_selected_rows': len(period_end_audit), 'period_boundaries': m13_period_boundaries, 'warm_up_policy': m13_warm_up_policy, 'design_sha256': m13_design_sha256, 'preregistration_sha256': m13_plan_sha256, 'exclusion_ledger_sha256': m13_exclusion_sha256, 'exploration_key_sha256': m13_plan['exploration_key_sha256'], 'replication_key_sha256': m13_plan['replication_key_sha256']}  # 串联上章、输入、期末日、三边界和两期键
m13_lineage_bytes = (json.dumps(m13_lineage, ensure_ascii=False, sort_keys=True, separators=(',', ':')) + '\n').encode('utf-8')  # 稳定序列化血缘
m13_write_same(m13_artifact_dir / 'm13-lineage-v3.json', m13_lineage_bytes)  # 幂等冻结血缘
列表 13.43: M13 稳定 handoff 的确定性制品哈希投影 
m13_artifact_hashes = {'design': m13_design_sha256, 'preregistration': m13_plan_sha256, 'period_end_closes': period_end_sha256, 'exclusion_ledger': m13_exclusion_sha256, 'family': hashlib.sha256(m13_family_bytes).hexdigest(), 'by_results': hashlib.sha256(m13_by_bytes).hexdigest(), 'discoveries': discovery_sha256, 'replication': hashlib.sha256(m13_replication_bytes).hexdigest(), 'lineage': hashlib.sha256(m13_lineage_bytes).hexdigest()}  # 汇总全部必需稳定制品指纹
m13_stable_final_names = {'m13-design-v3.json', 'm13-preregistration-v3.json', 'm13-period-end-closes-v1.csv', 'm13-universe-exclusion-ledger-v3.csv', 'm13-family-v3.csv', 'm13-by-results-v3.csv', 'm13-frozen-discoveries-v1.json', 'm13-replication-v3.csv', 'm13-lineage-v3.json', 'm13-handoff-v5.json'}  # 固定十份稳定最终名称
m13_handoff = {'schema_version': '5', 'preregistration_sha256': m13_plan_sha256, 'period_boundaries': m13_period_boundaries, 'warm_up_policy': m13_warm_up_policy, 'planned_family_n': len(momentum_family_df), 'estimated_family_n': int(momentum_family_df['status'].eq('estimated').sum()), 'skipped_family_n': int(momentum_family_df['status'].eq('skipped').sum()), 'frozen_discovery_n': len(frozen_discoveries), 'stable_final_names': sorted(m13_stable_final_names), 'artifacts': m13_artifact_hashes, 'telemetry_contract': {'schema_version': '5', 'filename_pattern': 'm13-run-<run_id>-v5.json'}}  # 声明计数、精确名称和独立telemetry接口
m13_handoff_bytes = (json.dumps(m13_handoff, ensure_ascii=False, sort_keys=True, separators=(',', ':')) + '\n').encode('utf-8')  # 稳定序列化最终交接
m13_write_same(m13_artifact_dir / 'm13-handoff-v5.json', m13_handoff_bytes)  # 完成态telemetry前先核验或写入handoff
visible_before_terminal = {path.name for path in m13_artifact_dir.iterdir() if path.is_file()}  # 枚举授权目录实际文件
if visible_before_terminal != m13_stable_final_names:
    m13_stop('artifact_final_name_set_mismatch', expected=sorted(m13_stable_final_names), actual=sorted(visible_before_terminal))  # 缺失、额外与临时残留均停止
列表 13.44: M13 每个 run_id 的唯一 completed 或 stopped telemetry 
required_run_fields = ['BOOK_RUN_STARTED_AT', 'BOOK_RUN_COMMAND', 'BOOK_GIT_COMMIT']  # 冻结运行身份必需字段
missing_run_fields = [field for field in required_run_fields if not os.environ.get(field, '').strip()]  # 检查环境提供情况
if missing_run_fields:  # 缺字段时不能声称可复现
    m13_stop('run_identity_incomplete', missing=missing_run_fields, loaded_rows={'price_slice': len(stock_price_raw), 'basic_slice': len(stock_basic_info), 'planned_family': len(momentum_family_df)})
run_started_at = datetime.fromisoformat(m13_require_env('BOOK_RUN_STARTED_AT'))  # 解析运行器固定的起始时刻
if run_started_at.tzinfo is None:  # 禁止无时区运行日志
    m13_stop('run_started_at_timezone_missing')  # 要求明确时区
m13_finished_at = datetime.now(ZoneInfo('Asia/Shanghai'))  # 记录完成时刻
m13_active_seconds = time.perf_counter() - m13_started_perf  # Python链实际活动计时
m13_wall_seconds = (m13_finished_at - run_started_at.astimezone(ZoneInfo('Asia/Shanghai'))).total_seconds()  # 从运行器起点计算含启动开销墙钟时间
if not np.isfinite(m13_wall_seconds) or m13_wall_seconds < 0 or m13_wall_seconds + 1e-6 < m13_active_seconds:
    m13_stop('run_timing_invalid', wall_seconds=m13_wall_seconds, active_seconds=m13_active_seconds)  # 两类计时必须独立且物理一致
m13_peak_memory_bytes = tracemalloc.get_traced_memory()[1]  # 读取 fresh-kernel Python 峰值
m13_loaded_rows = {'price_slice': int(len(stock_price_raw)), 'basic_slice': int(len(stock_basic_info)), 'm12_universe_roster': int(len(universe_roster)), 'eligible_securities': int(len(eligible_security_ids)), 'industry_months': int(len(industry_monthly_returns)), 'planned_family': int(len(momentum_family_df))}  # 明示关键输入和分析表行数
m13_run_log = {'schema_version': '5', 'run_id': m13_run_id, 'status': 'completed', 'started_at': run_started_at.astimezone(ZoneInfo('Asia/Shanghai')).isoformat(), 'finished_at': m13_finished_at.isoformat(), 'wall_seconds': round(m13_wall_seconds, 6), 'active_seconds': round(m13_active_seconds, 6), 'elapsed_seconds': round(m13_active_seconds, 6), 'peak_memory_bytes': int(m13_peak_memory_bytes), 'loaded_rows': m13_loaded_rows, 'command': os.environ['BOOK_RUN_COMMAND'], 'git_commit': os.environ['BOOK_GIT_COMMIT'], 'python': sys.version.split()[0], 'packages': {'pandas': pd.__version__, 'numpy': np.__version__, 'scikit_learn': sklearn.__version__}, 'period_boundaries': m13_period_boundaries, 'warm_up_policy': m13_warm_up_policy, 'preregistration_sha256': m13_plan_sha256, 'handoff_sha256': hashlib.sha256(m13_handoff_bytes).hexdigest(), 'input_hashes': m13_actual_hashes, 'm12_run_log_sha256': m12_run_log_expected_sha256, 'artifact_hashes': m13_artifact_hashes}  # 分开记录墙钟与活动时间并绑定稳定handoff
m13_run_log_path, m13_run_log_bytes = m13_write_terminal(m13_run_log)  # 稳定制品完成后排他写入唯一完成态
expected_completed_names = m13_stable_final_names | {f'm13-run-{m13_run_id}-v5.json'}  # 展开本次运行唯一终态名称
actual_completed_names = {path.name for path in m13_artifact_dir.iterdir() if path.is_file()}  # 完成后重新枚举实际文件
if actual_completed_names != expected_completed_names:
    raise RuntimeError({'status': 'terminal_unchanged', 'reason': 'completed_output_set_mismatch', 'expected': sorted(expected_completed_names), 'actual': sorted(actual_completed_names)})  # 完成态后不得发现未授权文件
print({'m13_handoff_sha256': hashlib.sha256(m13_handoff_bytes).hexdigest(), 'm13_run_log': m13_run_log_path.name, 'm13_run_log_sha256': hashlib.sha256(m13_run_log_bytes).hexdigest(), 'artifact_hashes': m13_handoff['artifacts']})  # 分别输出稳定入口与本次遥测

代码动态报告完整计划族、实际可估数、跳过数与 BY 拒绝数，再把探索期冻结发现集原样带入不重叠复现期。HAC 处理单一序列的时间依赖，BY 处理探索检验间未建模依赖，Holm 控制冻结发现集的复现 FWER；未拒绝只表示当前证据不足，复现结果也不给因果或可交易命题判定真值，仍需基准、成本与执行审计。

13.13.3 理论题解答

  1. BH 在独立检验下控制 FDR

    令 \(I_i\) 表示真零假设 \(i\) 被拒绝，\(R\) 为总拒绝数。按 \(R=1,\ldots,m\) 分解 \[ E\!\left[\frac{I_i}{\max(R,1)}\right] =\sum_{r=1}^m\frac1rP\!\left(P_i\le \frac{rq}{m},R^{(-i)}=r\right), \] 其中 \(R^{(-i)}\) 是把 \(P_i\) 置零后的 step-up 拒绝数。独立且连续均匀的真零 p 值给出 \[ \sum_{r=1}^m\frac1r\frac{rq}{m}P(R^{(-i)}=r)=\frac qm. \] 对 \(m_0\) 个真零假设求和，\(\mathrm{FDR}=qm_0/m\le q\)。若只满足 super-uniform，则等号改为上界；任意依赖下该分解不成立。

  2. Bonferroni 控制 FWER（完整推导）

    令\(A_j=\{p_j\le\alpha/m\}\)，则 \[ \mathrm{FWER} = P\Big(\cup_{j\in H_0} A_j\Big) \le \sum_{j\in H_0} P(A_j) \le \sum_{j\in H_0} \alpha/m \le \alpha. \] 不需要独立性。

  3. Holm 与 Hochberg 的 FWER 控制

    Holm 中若发生至少一个真零误拒绝，设排序中第一个真零位于 \(k\)，则在到达它时至少还有 \(m_0\) 个真零未处理，必有 \(P_{(k)}\le\alpha/(m-k+1)\le\alpha/m_0\)。对真零 p 值用并集上界，概率至多 \(m_0(\alpha/m_0)=\alpha\)，无需独立性。

    Hochberg 是 step-up 过程。发生真零误拒绝时，某个真零有序 p 值必须满足相应 Simes 边界；在独立或满足所需正依赖条件时，Simes 不等式把该事件概率界在 \(\alpha\)。因此 Hochberg 控制 FWER，但不能像 Holm 一样在任意依赖下无条件保证；任意依赖时应改用 Holm 等方法。

  4. 延伸主题写作要点

    • 自适应 FDR：说明如何估计\(\pi_0\) 以及对阈值的影响。
    • 依赖性校正：解释 BY 校正相当于把阈值除以调和级数常数，从而更保守。
    • 贝叶斯方法：强调先验、后验与局部fdr 的关系，以及与经验贝叶斯的链接。

  5. 四种校正的同族手算

    令 \(H_{(i)}\) 对应 \(p_{(i)}\)。Bonferroni 的共同阈值为 \(0.05/4=0.0125\)，故拒绝 \(\{H_{(1)},H_{(2)}\}\)。Holm 阈值依次为 \((0.0125,0.0167,0.0250,0.0500)\)；前三步通过、第四步停止，拒绝 \(\{H_{(1)},H_{(2)},H_{(3)}\}\)。

    BH 阈值为 \((0.0125,0.0250,0.0375,0.0500)\)，最大满足序号为 \(k=3\)，拒绝前三项。BY 有 \(c(4)=1+1/2+1/3+1/4=25/12\)，阈值为 \((0.0060,0.0120,0.0180,0.0240)\)，最大满足序号为 \(k=2\)，拒绝前两项。

    Bonferroni 与 Holm 在边际 p 值有效时允许任意依赖并控制 FWER；BH 控制 FDR 需要独立或适当正依赖（PRDS）；BY 在边际 p 值有效时允许任意依赖并控制 FDR。其余假设均记为“不拒绝”；所有结论只描述当前设计下反对零假设的证据强弱，不给任何假设命题判定真值。

13.14 章末闭环

逐项目标自检：定义 planned family 与错误度量；为每个假设选择依赖稳健的单项标准误；说明 Holm、BH 与 BY 的适用条件；在看 p 值前冻结跳过规则和复现期；把探索发现集带到不重叠时期复现。另设三项不可抵消出口：完整的“验证多行业 × 12”计划族、明确的单项时间依赖与家族依赖控制、探索与复现严格不重叠且复制缺失不改家族；任一项失败都不得完成 M13。

禁用情境：单项 p 值无效、planned family 事后缩小或探索/复现时期重叠时，不应报告“经过校正的发现”；没有成本、基准和执行假设时，不应把拒绝零假设称为可交易收益。常见误区是把未拒绝解释为无效应，或把 FDR 控制误解为每个已拒绝假设为真的概率保证。

无提示检索：1）任意依赖下为何通常选择 BY 而不是直接套用 BH？2）跳过的预注册假设为何仍须留在 planned family 审计中？3）多重校正与独立时期复现分别解决什么问题？

展开检索反馈与全书出口决策

1）经典 BH 保证依赖独立或适当正依赖条件；BY 通过调和级数修正提供任意依赖下的控制。2）事后删除失败项会改变家族并带来选择性报告；应保留状态和原因，按预注册规则处理分母。3）校正控制同一家族内的错误累积，独立复现检查探索选择在新时期能否重现，两者不能互相替代。检索题每题 1 分；错第 1 题回到“依赖结构与 BY”并用相关 p 值矩阵异形题复测，错第 2 题回到 lst-ex8-by-correction 并用两项不可估的新家族复测，错第 3 题回到 lst-ex8-period-freeze 与 lst-ex8-independent-replication 并用改变边界月的新时间轴复测，答对错项才重入。三项出口另各 1 分且必须 3 分：计划族逐项可追溯、HAC+BY/Holm 分工正确、非重叠复现与冻结家族证据完整；任一缺项均从 lst-ex8-preregistration 重跑 fresh-kernel 链，概念题得分不得抵消。

全书综合出口：提交 M02—M09、A10、A11、M12、M13 的冻结合同、制品 hash、失败日志、共同测试比较、无监督结构审计和独立复现表；再由同伴从干净环境抽取一个里程碑复算，并在答辩中解释哪些结论只适用于当前数据、损失与验证设计。无法复算的环节进入补修清单，不得以展示性图表替代证据。

Benjamini, Yoav, 和 Yosef Hochberg. 1995年. 《Controlling the False Discovery Rate: A Practical and Powerful Approach to Multiple Testing》. Journal of the Royal Statistical Society: Series B (Methodological) 57 (1): 289～300. https://doi.org/10.1111/j.2517-6161.1995.tb02031.x.
Benjamini, Yoav, 和 Daniel Yekutieli. 2001年. 《The Control of the False Discovery Rate in Multiple Testing under Dependency》. The Annals of Statistics 29 (4): 1165～88. https://doi.org/10.1214/aos/1013699998.
Holm, Sture. 1979年. 《A Simple Sequentially Rejective Multiple Test Procedure》. Scandinavian Journal of Statistics 6 (2): 65～70. https://www.jstor.org/stable/4615733.
Storey, John D. 2002年. 《A Direct Approach to False Discovery Rates》. Journal of the Royal Statistical Society: Series B (Statistical Methodology) 64 (3): 479～98. https://doi.org/10.1111/1467-9868.00346.