目录

中文EN

理论篇：因子模型与风格

1、单因子模型代码2、多因子模型理论基础3、经典风格因子4、因子版图扩展

操作篇：因子生产与管理

5、因子构建与中性化6、因子测试框架代码7、因子组合优化与归因8、因子择时与拥挤

拓展篇：宏观量化与资产配置

9、宏观数据与因子10、宏观状态识别11、大类资产配置代码12、宏观量化十年笔记代码

1、单因子模型

以下先用几小节交代因子投资的基本概念（原第1章导论），再进入 CAPM 与 A 股单因子复现。

这是一个我在学习量化投资研究过程中持续整理的个人笔记和代码，主要记录了从单因子到多因子模型的知识框架、实操方法和阶段性思考。内容不追求标准答案，更注重把学习过程中的理解、踩坑和复盘真实地沉淀下来。希望这些笔记能为同样在路上的你提供一点参考，也非常欢迎交流讨论、互相启发。联系邮箱：adingminjie@163.com

导论.1 什么是因子？

因子，通俗来讲，就是能解释股票收益的共同特征。为什么有些股票涨得好？可能是因为它们便宜，可能是因为它们盈利能力强，也可能是因为市场情绪好。这些选择股票的逻辑，就是因子。

举个例子：

• 价值因子 ：买便宜的股票
• 规模因子 ：买小市值的股票
• 盈利因子 ：买盈利能力强、盈利更稳健的股票
• 动量因子 ：买涨得好的股票
• 质量因子 ：买赚钱能力强的股票
• 低波因子 ：买波动小的股票

每个因子都代表一种风险溢价。投资者承担了某种风险，市场就给予相应的回报。这就是因子的核心逻辑。另外，因子不是凭空捏造的，它背后一定有经济学逻辑支撑。没有逻辑的因子，再漂亮也是垃圾。

导论.2 因子投资的起源与发展

1960年代 ，夏普提出了CAPM模型。那时候大家觉得，股票的收益只跟市场风险有关。你承担了市场风险，就能获得市场收益。

1970-80年代 ，Fama和French开始质疑这个观点。他们发现，小盘股和低市净率的股票，长期表现更好。于是1993年，他们提出了著名的三因子模型：市场、规模、价值。

1990年代至今 ，因子家族不断壮大：

• 动量因子（Jegadeesh & Titman, 1993）
• 质量因子（Novy-Marx, 2013）
• 低波因子（Ang et al., 2006）
• 红利因子、成长因子、情绪因子...

现在主流的因子库，少说有上百个。但真正有效的，其实也就那么十来个。建议先聚焦5-10个经典因子，把它们的逻辑吃透，再慢慢扩展。

因子投资已成为大学金融学教育中的核心内容之一。以北京大学光华管理学院为例，本科大二的专业课《证券投资学》就已开始系统介绍因子投资的完整脉络；在更高层次培养中，清华大学五道口金融学院博士阶段也开设了《实证资产定价》和《理论资产定价》等相关课程。

导论.3 因子投资的核心逻辑

因子投资的优势，可以归纳成三句话：

1. 收益可分解 ：任何股票的收益，都可以拆解成因子暴露和因子收益的乘积
2. 风险可管理 ：通过控制因子暴露，你可以精确管理投资组合的风险
3. 收益可预测 ：因子收益在长期是稳定的，短期虽有波动，但方向不变

但理论成立，不等于实操简单。市场环境会持续变化，因子有效性也会随时间演化，这正是因子投资最核心的挑战。

基于以上三点，落地到实操层面，通常遵循三条原则：

1. 分散化：避免单一因子暴露，多因子组合更稳健。
2. 风险控制：将因子暴露保持在可承受范围内，防止单一因子拖累整体组合。
3. 动态调整：跟踪因子生命周期，定期检视并更新配置。

此外，因子研究一定要警惕过拟合。策略在历史数据上表现再好，也不代表具备稳定的实盘能力。应通过样本外检验、滚动验证和简化模型结构来提升策略的稳健性。

在因子投资研究中，CAPM是无法回避的基础框架，也是金融学本科大二阶段较早系统接触的核心模型之一。

初学时，许多人会质疑其设定过于简化，尤其是仅以一个$\beta$系数刻画资产收益与市场风险的关系，似乎难以覆盖现实市场中的复杂异质性。然而，随着研究与实践的深入可以发现，CAPM 的价值并不在于解释一切，而在于提供了风险定价的最小理论起点。这种高度简化的结构，为后续多因子模型的扩展与检验奠定了方法论基础。

1.1 CAPM模型介绍

CAPM最早由Sharpe（1964）提出，并由Lintner（1965）与Mossin（1966）进一步完善。公式如下：

$$ E(R_i) = R_f + \beta_i \bigl(E(R_m) - R_f\bigr) $$

意思是： 一只股票的预期收益 = 无风险利率 + 它的$\beta$ × 市场风险溢价

$\beta$(beta)是股票对市场的敏感度。$\beta=1$，市场涨1%它涨1%；$\beta=1.5$，市场涨1%它涨1.5%。

CAPM认为，只有系统性风险（市场风险）才被定价。非系统性风险可以通过分散化消除，所以市场不会给投资者额外补偿。

高$\beta$股票在上行市场中通常表现出更高的收益弹性，但在下行阶段也往往伴随更大的回撤幅度。这一现象从经验层面揭示了CAPM的关键局限：模型仅以市场因子解释资产收益，难以充分刻画由风格、行业与基本面等非市场维度所驱动的收益差异。

1.2 单因子回归分析

如何计算$\beta$？最直接的方法就是用线性回归。对于个股，建议使用滚动窗口算$\beta$，窗口长度取252个交易日，也就是1年。为什么选取1年，因为太短了噪声大，太长了反应迟钝。

import numpy as np
import pandas as pd
import statsmodels.api as sm

def calculate_beta(stock_returns, market_returns, window=252):
    """
    滚动计算beta值
    stock_returns: 个股日收益率序列
    market_returns: 市场日收益率序列
    window: 滚动窗口长度
    """
    betas = pd.Series(index=stock_returns.index, dtype=float)

    for i in range(window, len(stock_returns)):
        y = stock_returns.iloc[i - window:i]
        x = market_returns.iloc[i - window:i]
        x = sm.add_constant(x)

        model = sm.OLS(y, x).fit
        betas.iloc[i] = model.params.iloc[1]  # 斜率就是beta

    return betas

def calculate_beta_newey_west(stock_returns, market_returns, window=36, nw_lags=3, return_stats=True):
    """
    滚动计算 beta，并使用 Newey-West(HAC) 稳健标准误。
    stock_returns: 个股收益率序列（建议与 market_returns 同频）
    market_returns: 市场收益率序列
    window: 滚动窗口长度（月度数据常用 24/36/60）
    nw_lags: Newey-West 的滞后阶数
    return_stats: True 返回 alpha/beta/t/p；False 仅返回 beta 序列
    """
    df = pd.concat([stock_returns, market_returns], axis=1)
    df.columns = ["ret_i", "ret_m"]
    df = df.dropna.copy

    if nw_lags is None:
        nw_lags = max(1, int(window ** 0.25))

    result = pd.DataFrame(
        index=df.index,
        columns=["alpha", "beta", "beta_se_nw", "beta_t_nw", "beta_p_nw"],
        dtype=float,
    )

    for i in range(window, len(df) + 1):
        win = df.iloc[i - window:i]
        y = win["ret_i"]
        X = sm.add_constant(win["ret_m"])

        model = sm.OLS(y, X).fit(cov_type="HAC", cov_kwds={"maxlags": nw_lags})
        t = win.index[-1]

        result.loc[t, "alpha"] = model.params.get("const", np.nan)
        result.loc[t, "beta"] = model.params.get("ret_m", np.nan)
        result.loc[t, "beta_se_nw"] = model.bse.get("ret_m", np.nan)
        result.loc[t, "beta_t_nw"] = model.tvalues.get("ret_m", np.nan)
        result.loc[t, "beta_p_nw"] = model.pvalues.get("ret_m", np.nan)

    if return_stats:
        return result
    return result["beta"]

# 月度数据示例：
# beta_df = calculate_beta_newey_west(stock_returns, market_returns, window=36, nw_lags=3, return_stats=True)
# beta_series = calculate_beta_newey_west(stock_returns, market_returns, window=36, nw_lags=3, return_stats=False)

提一下计量经济学中的基本知识，这里的回归需要加截距项。否则，回归会强制过原点，算出来的beta会偏。

细节上，注意选择收益率的形式，log(1+r)≈r。简单收益率和对数收益率都可以用于$\beta$估计，但个股、市场与无风险利率必须使用一致的收益定义。

1.3 因子收益与因子暴露

因子收益和因子暴露是两个容易混淆的概念。

因子暴露：是在这个因子上的敞口/载荷。比如beta=1.2，这就是对市场因子的暴露。

因子收益：是这个因子本身赚了多少钱。比如这个月市场涨了3%，那市场因子的收益就是3%。

比如，一只股票的收益 = 因子暴露 × 因子收益 + 残差。

1.4 单因子模型的局限

CAPM并不能完全解释股票涨跌。一个典型例子是2017年的A股白马行情。贵州茅台全年涨幅超过100%，而同期市场指数仅上涨约6%。若完全依据CAPM解释，这会对应一个异常高的$\alpha$。但这种超额收益并不一定是定价偏差，更可能反映了价值、质量等其他系统性因子的暴露。

这也是资产定价模型持续扩展的原因：从Fama-French三因子到五因子，再到更高维的多因子框架。本质上，这一演进是在CAPM基础上逐步纳入被遗漏的风险来源。

因此，CAPM并非无效，而是边界明确。它适合刻画最基础的市场风险补偿关系，却难以覆盖现实市场中的全部收益结构。单因子模型的真正价值在于提供了因子定价的分析起点；在此基础上，多因子体系才得以形成更强的解释力与应用价值。

1.5 A股数据与 CAPM 拟合

本节完成 A 股月频数据准备与 CAPM 回归（全样本 + 个股层面），为 2.6 节滚动 beta 估计提供数据口径对照。

读入数据

import pandas as pd
from urllib.parse import urlparse, parse_qs, urlencode, urlunparse
from IPython.display import display

# 1) Dropbox 分享链接（原始分享页 URL）
dropbox_urls = {
    "market_return": "https://www.dropbox.com/scl/fi/wd4zxa3safkrjyt4xjy3y/AShare_Market_Return_Monthly.xlsx?rlkey=glixs89u4unkwu3j0w478rdmh&st=38c5588m&dl=0",
    "riskfree_rate": "https://www.dropbox.com/scl/fi/hfl350d4t9fy3mw6dd9cp/AShare_Riskfree_Rate_Monthly.xlsx?rlkey=jrcfsu9vk856kdoigodfr03o7&st=ex1vxhiy&dl=0",
    "stock_return": "https://www.dropbox.com/scl/fi/roaossdn05bvxyr4b8n1i/AShare_Stocks_Return_Monthly.xlsx?rlkey=tw72id38dgyxd7qfso090hjgg&st=27ej4fqj&dl=0",
}

# 2) 将分享链接改为可直接下载链接（dl=1）
def to_direct_download(url: str) -> str:
    u = urlparse(url)
    q = parse_qs(u.query)
    q["dl"] = ["1"]  # 强制下载原文件而不是网页预览
    new_query = urlencode(q, doseq=True)
    return urlunparse((u.scheme, u.netloc, u.path, u.params, new_query, u.fragment))

direct_urls = {k: to_direct_download(v) for k, v in dropbox_urls.items}

# 3) 读入 Excel 到字典 dfs，后续统一从 dfs 取表
dfs = {}
for name, url in direct_urls.items:
    dfs[name] = pd.read_excel(url)

# 4) 快速检查：看每张表的大小和前几行
for name, df in dfs.items:
    print(f"\n{name} shape = {df.shape}")
    display(df.head)

d:\pythonprojects\venv\Lib\site-packages\openpyxl\styles\stylesheet.py:237: UserWarning: Workbook contains no default style, apply openpyxl's default
  warn("Workbook contains no default style, apply openpyxl's default")
d:\pythonprojects\venv\Lib\site-packages\openpyxl\styles\stylesheet.py:237: UserWarning: Workbook contains no default style, apply openpyxl's default
  warn("Workbook contains no default style, apply openpyxl's default")
d:\pythonprojects\venv\Lib\site-packages\openpyxl\styles\stylesheet.py:237: UserWarning: Workbook contains no default style, apply openpyxl's default
  warn("Workbook contains no default style, apply openpyxl's default")

market_return shape = (8119, 4)

	Markettype	Trdmnt	Cmretwdos	Cmretwdtl
0	市场类型	交易月份	考虑现金红利再投资的综合月市场回报率(流通市值加权平均法)	考虑现金红利再投资的综合月市场回报率(总市值加权平均法)
1	没有单位	没有单位	没有单位	没有单位
2	5	1990-12	NaN	NaN
3	15	1990-12	NaN	NaN
4	21	1990-12	NaN	NaN

riskfree_rate shape = (13016, 3)

	Nrr1	Clsdt	Nrrmtdt
0	无风险利率基准	统计日期	月度化无风险利率(%)
1	没有单位	没有单位	没有单位
2	NRI01	1990-04-15	0.8035
3	NRI01	1990-08-21	0.693
4	NRI01	1990-12-19	0.693

stock_return shape = (942618, 8)

	Stkcd	Trdmnt	Clsdt	Mclsprc	Msmvosd	Msmvttl	Mretwd	Markettype
0	证券代码	交易月份	月收盘日期	月收盘价	月个股流通市值	月个股总市值	考虑现金红利再投资的月个股回报率	市场类型
1	没有单位	没有单位	没有单位	元/股	千元	千元	没有单位	没有单位
2	000001	1991-04	30	43.68	1157520	2118487.47	NaN	4
3	000001	1991-05	31	38.34	1016010	1859496.56	-0.122253	4
4	000001	1991-06	28	33.99	900735	1648520.81	-0.113459	4

CAPM检验（中国市场，Markettype=53）

下面代码完成四件事：

1. 自动识别三张表的关键字段（股票/月度收益、市场收益、无风险利率）
2. 过滤市场类型为53（沪深A股+创业板+科创板）
3. 构造个股与市场超额收益
4. 给出CAPM回归结果（整体样本 + 个股截面汇总）

import numpy as np
import statsmodels.api as sm

# 三张原始表（来自上一个单元的 dfs）
stock_raw = dfs["stock_return"].copy
market_raw = dfs["market_return"].copy
rf_raw = dfs["riskfree_rate"].copy

# 手工映射列名（按你提供的 CSMAR 字段）
# 说明：固定列名比自动识别更稳，也更便于复现。
STOCK_CODE_COL = "Stkcd"
STOCK_MONTH_COL = "Trdmnt"
STOCK_RET_COL = "Mretwd"

MARKET_TYPE_COL = "Markettype"
MARKET_MONTH_COL = "Trdmnt"
MARKET_RET_COL = "Cmretwdos"  # 含红利再投资 + 流通市值加权

RF_MONTH_COL = "Clsdt"
RF_COL = "Nrrmtdt"  # 月度化无风险利率(%)

def to_month_str(x):
    """将日期列统一为 YYYY-MM 字符串，便于跨表按月合并。"""
    return pd.to_datetime(x, errors="coerce").dt.to_period("M").astype("string")

def to_return_decimal(s):
    """
    将收益率列尽量转换为小数口径：
    - 若中位数明显 > 0.5，通常表示百分数口径（如 3.2），自动 /100
    - 否则视作已是小数口径（如 0.032）
    """
    s = pd.to_numeric(s, errors="coerce")
    s_clean = s.dropna
    if len(s_clean) == 0:
        return s
    if s_clean.abs.median > 0.5:
        return s / 100.0
    return s

# 检查关键列是否存在，避免后续 silently fail
required_cols = {
    "stock": [STOCK_CODE_COL, STOCK_MONTH_COL, STOCK_RET_COL],
    "market": [MARKET_TYPE_COL, MARKET_MONTH_COL, MARKET_RET_COL],
    "rf": [RF_MONTH_COL, RF_COL],
}

missing_cols = {
    "stock": [c for c in required_cols["stock"] if c not in stock_raw.columns],
    "market": [c for c in required_cols["market"] if c not in market_raw.columns],
    "rf": [c for c in required_cols["rf"] if c not in rf_raw.columns],
}

if any(missing_cols[k] for k in missing_cols):
    raise ValueError(f"列名不匹配，请检查: {missing_cols}")

print("列名映射已固定，不再自动识别。")
print("stock:", required_cols["stock"])
print("market:", required_cols["market"])
print("rf:", required_cols["rf"])

列名映射已固定，不再自动识别。
stock: ['Stkcd', 'Trdmnt', 'Mretwd']
market: ['Markettype', 'Trdmnt', 'Cmretwdos']
rf: ['Clsdt', 'Nrrmtdt']

# 1) 个股数据：取所需列 -> 标准命名 -> 时间统一 -> 收益率转小数
START_MONTH = "2000-01" # 仅保留2000年及以后样本
stock_df = stock_raw[[STOCK_CODE_COL, STOCK_MONTH_COL, STOCK_RET_COL]].copy
stock_df.columns = ["stock_code", "month", "ret_i"]
stock_df["month"] = to_month_str(stock_df["month"])
stock_df["ret_i"] = to_return_decimal(stock_df["ret_i"])

stock_df = stock_df.dropna(subset=["stock_code", "month", "ret_i"])

stock_df["stock_code"] = stock_df["stock_code"].astype(str).str.strip
stock_df = stock_df[stock_df["stock_code"].str.match(r"^\d{6}$", na=False)]
stock_df = stock_df[stock_df["month"] >= START_MONTH].copy

market_df = market_raw[[MARKET_TYPE_COL, MARKET_MONTH_COL, MARKET_RET_COL]].copy
market_df.columns = ["market_type", "month", "ret_m"]
market_df["month"] = to_month_str(market_df["month"])
market_df["ret_m"] = to_return_decimal(market_df["ret_m"])

market_df["market_type_num"] = pd.to_numeric(market_df["market_type"], errors="coerce")

market_df = market_df[(market_df["market_type_num"] == 53) | (market_df["market_type"].astype(str).str.strip == "53")].copy
market_df = market_df[["month", "ret_m"]].dropna.drop_duplicates(subset=["month"])
market_df = market_df[market_df["month"] >= START_MONTH].copy

rf_df = rf_raw[[RF_MONTH_COL, RF_COL]].copy
rf_df.columns = ["month", "rf"]
rf_df["month"] = to_month_str(rf_df["month"])

rf_df["rf"] = pd.to_numeric(rf_df["rf"], errors="coerce") / 100.0
rf_df = rf_df.dropna.drop_duplicates(subset=["month"])
rf_df = rf_df[rf_df["month"] >= START_MONTH].copy

print(f"样本窗口: {START_MONTH} 及以后")
print("样本规模:")
print("stock_df:", stock_df.shape)
print("market_df (type=53):", market_df.shape)
print("rf_df:", rf_df.shape)
display(stock_df.head)
display(market_df.head)
display(rf_df.head)
print("\nrf描述统计（应为小数口径，例如0.006092表示0.6092%）：")
display(rf_df["rf"].describe)

<ipython-input-4-379ecc8467f1>:24: UserWarning: Could not infer format, so each element will be parsed individually, falling back to `dateutil`. To ensure parsing is consistent and as-expected, please specify a format.
  return pd.to_datetime(x, errors="coerce").dt.to_period("M").astype("string")

样本窗口: 2000-01 及以后
样本规模:
stock_df: (890537, 3)
market_df (type=53): (319, 2)
rf_df: (320, 2)

<ipython-input-4-379ecc8467f1>:24: UserWarning: Could not infer format, so each element will be parsed individually, falling back to `dateutil`. To ensure parsing is consistent and as-expected, please specify a format.
  return pd.to_datetime(x, errors="coerce").dt.to_period("M").astype("string")
<ipython-input-4-379ecc8467f1>:24: UserWarning: Could not infer format, so each element will be parsed individually, falling back to `dateutil`. To ensure parsing is consistent and as-expected, please specify a format.
  return pd.to_datetime(x, errors="coerce").dt.to_period("M").astype("string")

	stock_code	month	ret_i
107	000001	2000-01	0.061891
108	000001	2000-02	-0.011333
109	000001	2000-03	0.002729
110	000001	2000-04	0.037017
111	000001	2000-05	-0.055118

	month	ret_m
2066	2000-01	0.160838
2085	2000-02	0.122024
2104	2000-03	0.055926
2123	2000-04	0.013014
2142	2000-05	0.027691

	month	rf
3304	2000-01	0.001856
3335	2000-02	0.001856
3364	2000-03	0.001856
3395	2000-04	0.001856
3425	2000-05	0.001856

rf描述统计（应为小数口径，例如0.006092表示0.6092%）：

count    320.000000
mean       0.001752
std        0.000612
min        0.000788
25%        0.001241
50%        0.001635
75%        0.002060
max        0.003386
Name: rf, dtype: float64

# 3) 按 month 合并三张表，并构造超额收益
capm_df = (
    stock_df
    .merge(market_df, on="month", how="inner")
    .merge(rf_df, on="month", how="inner")
)
capm_df = capm_df.dropna(subset=["ret_i", "ret_m", "rf", "stock_code", "month"]).copy
capm_df["excess_i"] = capm_df["ret_i"] - capm_df["rf"]
capm_df["excess_m"] = capm_df["ret_m"] - capm_df["rf"]

print("CAPM合并后样本:", capm_df.shape)
print("股票数:", capm_df["stock_code"].nunique, "月份数:", capm_df["month"].nunique)
display(capm_df.head)

# 4) Pooled CAPM 回归：全样本合并后的一元回归
X_pool = sm.add_constant(capm_df[["excess_m"]])
y_pool = capm_df["excess_i"]
pool_model = sm.OLS(y_pool, X_pool, missing="drop").fit(cov_type="HAC", cov_kwds={"maxlags": 3})

print("\nPooled CAPM (HAC, maxlags=3):")
print(pool_model.summary)

# 5) 个股层面 CAPM 回归：每只股票单独估计 alpha/beta
rows = []
min_obs = 24  # 最少24个月，避免样本过短导致估计不稳

for code, g in capm_df.groupby("stock_code"):
    g = g.dropna(subset=["excess_i", "excess_m"])

    if len(g) < min_obs:
        continue

    X = sm.add_constant(g[["excess_m"]])
    y = g["excess_i"]
    m = sm.OLS(y, X).fit(cov_type="HAC", cov_kwds={"maxlags": 3})

    rows.append({
        "stock_code": code,
        "n_obs": len(g),
        "alpha": m.params.get("const", np.nan),
        "beta": m.params.get("excess_m", np.nan),
        "alpha_t": m.tvalues.get("const", np.nan),
        "beta_t": m.tvalues.get("excess_m", np.nan),
        "alpha_p": m.pvalues.get("const", np.nan),
        "beta_p": m.pvalues.get("excess_m", np.nan),
        "r2": m.rsquared,
    })

capm_by_stock = pd.DataFrame(rows)
print("\n个股回归结果:", capm_by_stock.shape)
display(capm_by_stock.head)

# 6) 汇总指标：用于快速判断 CAPM 在样本中的表现
if len(capm_by_stock) > 0:
    summary_df = pd.DataFrame({
        "metric": ["alpha_mean", "beta_mean", "r2_mean", "alpha_sig_5pct_ratio"],
        "value": [
            capm_by_stock["alpha"].mean,
            capm_by_stock["beta"].mean,
            capm_by_stock["r2"].mean,
            (capm_by_stock["alpha_p"] < 0.05).mean,
        ],
    })
    print("\n汇总统计:")
    display(summary_df)
else:
    print("未达到最小样本要求（min_obs=24）的股票，请检查时间区间或字段映射。")

CAPM合并后样本: (890537, 7)
股票数: 5977 月份数: 319

	stock_code	month	ret_i	ret_m	rf	excess_i	excess_m
0	000001	2000-01	0.061891	0.160838	0.001856	0.060035	0.158982
1	000001	2000-02	-0.011333	0.122024	0.001856	-0.013189	0.120168
2	000001	2000-03	0.002729	0.055926	0.001856	0.000873	0.054070
3	000001	2000-04	0.037017	0.013014	0.001856	0.035161	0.011158
4	000001	2000-05	-0.055118	0.027691	0.001856	-0.056974	0.025835

Pooled CAPM (HAC, maxlags=3):
                            OLS Regression Results                            
==============================================================================
Dep. Variable:               excess_i   R-squared:                       0.207
Model:                            OLS   Adj. R-squared:                  0.207
Method:                 Least Squares   F-statistic:                 1.622e+05
Date:                Sat, 08 Aug 2026   Prob (F-statistic):               0.00
Time:                        21:00:30   Log-Likelihood:             4.6529e+05
No. Observations:              890537   AIC:                        -9.306e+05
Df Residuals:                  890535   BIC:                        -9.306e+05
Df Model:                           1                                         
Covariance Type:                  HAC                                         
==============================================================================
                 coef    std err          z      P>|z|      [0.025      0.975]
------------------------------------------------------------------------------
const          0.0045      0.000     30.965      0.000       0.004       0.005
excess_m       1.1083      0.003    402.792      0.000       1.103       1.114
==============================================================================
Omnibus:                  1713749.296   Durbin-Watson:                   2.060
Prob(Omnibus):                  0.000   Jarque-Bera (JB):      47072717237.924
Skew:                          14.038   Prob(JB):                         0.00
Kurtosis:                    1128.976   Cond. No.                         15.1
==============================================================================

Notes:
[1] Standard Errors are heteroscedasticity and autocorrelation robust (HAC) using 3 lags and without small sample correction

个股回归结果: (5729, 9)

	stock_code	n_obs	alpha	beta	alpha_t	beta_t	alpha_p	beta_p	r2
0	000001	317	0.001696	1.000171	0.435371	17.537020	0.663293	7.474557e-69	0.509825
1	000002	313	0.005495	0.872107	0.983000	8.214996	0.325607	2.121721e-16	0.254732
2	000004	309	-0.001135	1.055256	-0.141812	10.194937	0.887228	2.088800e-24	0.224673
3	000005	275	-0.001002	1.338899	-0.108143	7.944917	0.913882	1.943202e-15	0.216694
4	000006	314	0.004957	1.240507	0.904831	13.150238	0.365555	1.696403e-39	0.407076

汇总统计:

	metric	value
0	alpha_mean	0.004637
1	beta_mean	1.163326
2	r2_mean	0.237669
3	alpha_sig_5pct_ratio	0.035259

CAPM 回归结果解读（Pooled + 个股截面）

下面对 Pooled CAPM (HAC, maxlags=3) 与个股回归结果做简要解读：

1. 市场因子显著有效

• excess_m 系数约为 1.1083，且 z=402.792、p<0.001。
• 说明个股超额收益与市场超额收益存在非常显著的正相关关系。
• 系数大于 1，意味着样本平均意义上，股票组合对市场波动的暴露偏高。

2. 截距（alpha）为正且显著

• const=0.0045（若是月度小数口径，约为每月 0.45%）。
• 在控制市场因子后仍存在显著正的平均超额收益。
• 这通常提示：单因子 CAPM 不能完全解释收益，可能存在遗漏因子（规模、价值、质量、动量等）。

3. 模型解释力有限但合理

• R-squared=0.207，即 CAPM 解释了约 20.7% 的个股超额收益波动。
• 其余约 79.3% 来自其他风险因子、个股特质和噪声。
• 这与CAPM 是基础框架但解释力有限的理论预期一致。

4. 统计显著性很强，但需区分经济显著性

• 样本量 N=890,537，F 统计量非常高。
• 在超大样本下，参数显著并不意外；更关键是系数大小、可交易性和稳健性。

5. 残差呈明显厚尾分布

• Skew=14.038、Kurtosis=1128.976、JB 极大且显著。
• 说明残差远非正态，存在极端值与厚尾现象（股票收益中常见）。
• 使用 HAC（Newey-West）标准误是合理的，可缓解异方差和自相关对标准误的影响。

6. Durbin-Watson 接近 2

• DW=2.060，表明残差一阶自相关不明显。

7. 个股回归覆盖面较广

• 个股回归结果: (5729, 9) 表示对 5729 只股票估计了 alpha/beta/t/p/R2 等指标。
• 后续可重点看：
• alpha 的均值/中位数与显著比例；
• beta 的分布稳定性；
• 是否存在行业或风格暴露导致的伪 alpha。

一句话总结： CAPM 在样本中方向上成立（beta 显著），但解释不充分（R²有限且alpha显著），这正是从单因子走向多因子建模的经验依据。

1.6 滚动窗口 Beta 估计（A 股完整复现）

本节用 A 股数据 完整复现实证资产定价中标准的市场 beta 估计流程：先做单只股票的 CAPM 时序回归，再对全市场做闭式解滚动窗口估计，并比较月频与日频结果。

对应关系（字段命名尽量与常见学术代码一致）：

本节字段	含义	2.5 节对应
permno	股票代码	stock_code
date	月末/月首日期	month
ret_excess	个股超额收益	excess_i
mkt_excess	市场超额收益	excess_m

CAPM 时序回归式为

$$ r_{i,t}-r_{f,t}=\alpha_i+\beta_i(r_{m,t}-r_{f,t})+\varepsilon_{i,t}. $$

数据优先读取本地 202608_02_传统量化/Data 下的 CSMAR 导出文件（字段口径与 2.5 节 Dropbox 数据一致）。若此前已运行 2.5 并生成 capm_df，月频部分也可直接复用其合并逻辑。

1.6.1 环境与依赖

以下使用 pandas / numpy / statsmodels / matplotlib 实现。滚动估计函数 roll_capm_estimation 的输入输出结构，与常见 Python 复现代码中的同名函数对齐（permno, date, coefficient, estimate, t_statistic）。

import warnings
from pathlib import Path

import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
from statsmodels.regression.linear_model import OLS
from statsmodels.tools.tools import add_constant
from IPython.display import display

warnings.filterwarnings("ignore", category=UserWarning)

# 中文字体（与本章前文回测图一致，按本机可用字体回退）
from matplotlib import font_manager

_cn_candidates = ["Microsoft YaHei", "SimHei", "PingFang SC", "Noto Sans CJK SC", "Arial Unicode MS"]
_available = {f.name for f in font_manager.fontManager.ttflist}
for _name in _cn_candidates:
    if _name in _available:
        plt.rcParams["font.sans-serif"] = [_name]
        break
plt.rcParams["axes.unicode_minus"] = False

DATA_DIR = Path(r"D:/A_Topics/202607_03_TidyFinanceAShare/202608_02_传统量化/Data")
CACHE_DIR = DATA_DIR / "_cache_beta"
CACHE_DIR.mkdir(parents=True, exist_ok=True)

print("DATA_DIR =", DATA_DIR)
print("exists =", DATA_DIR.exists)

DATA_DIR = D:\A_Topics\202607_03_TidyFinanceAShare\202608_02_传统量化\Data
exists = True

1.6.2 读取并整理月频面板

步骤与 2.5 节相同：

1. 读取个股月收益、市场月收益、无风险利率；
2. 市场类型取 53（沪深 A 股 + 创业板 + 科创板，流通市值加权含红利再投资）；
3. 统一为小数口径，并构造超额收益。

随后把列名映射为 permno / date / ret_excess / mkt_excess，便于后面直接套用滚动估计函数。

def read_csmar_excel(path: Path) -> pd.DataFrame:
    """读取 CSMAR 风格 Excel：第 1 行为字段名，随后两行常为中文说明/单位。"""
    df = pd.read_excel(path, header=0)
    # 固定跳过两行注释（中文名/单位）；勿按首列是否数值过滤（无风险利率首列为 NRI01）
    if len(df) >= 2:
        df = df.iloc[2:].copy
    return df.reset_index(drop=True)

def to_month_period(x) -> pd.Series:
    return pd.to_datetime(x, errors="coerce").dt.to_period("M")

def to_return_decimal(s: pd.Series) -> pd.Series:
    s = pd.to_numeric(s, errors="coerce")
    clean = s.dropna
    if len(clean) == 0:
        return s
    if clean.abs.median > 0.5:
        return s / 100.0
    return s

def load_ashare_monthly(data_dir: Path, start_month: str = "2000-01") -> pd.DataFrame:
    stock_raw = read_csmar_excel(data_dir / "AShare_Stocks_Return_Monthly.xlsx")
    market_raw = read_csmar_excel(data_dir / "AShare_Market_Return_Monthly.xlsx")
    rf_raw = read_csmar_excel(data_dir / "AShare_Riskfree_Rate_Monthly.xlsx")

    stock = stock_raw[["Stkcd", "Trdmnt", "Mretwd", "Msmvosd"]].copy
    stock.columns = ["permno", "month", "ret", "mktcap"]
    stock["permno"] = stock["permno"].astype(str).str.zfill(6)
    stock["month"] = to_month_period(stock["month"])
    stock["ret"] = to_return_decimal(stock["ret"])
    stock["mktcap"] = pd.to_numeric(stock["mktcap"], errors="coerce")

    market = market_raw[["Markettype", "Trdmnt", "Cmretwdos"]].copy
    market.columns = ["market_type", "month", "ret_m"]
    market["market_type"] = pd.to_numeric(market["market_type"], errors="coerce")
    market = market[market["market_type"] == 53].copy
    market["month"] = to_month_period(market["month"])
    market["ret_m"] = to_return_decimal(market["ret_m"])
    market = market[["month", "ret_m"]].drop_duplicates("month")

    rf = rf_raw.copy
    if "Nrr1" in rf.columns:
        rf = rf[rf["Nrr1"].astype(str).str.upper.str.contains("NRI01", na=False)]
    rf = rf[["Clsdt", "Nrrmtdt"]].copy
    rf.columns = ["month", "rf"]
    rf["month"] = to_month_period(rf["month"])
    rf["rf"] = pd.to_numeric(rf["rf"], errors="coerce") / 100.0  # Nrrmtdt 为百分数
    rf = rf.dropna.drop_duplicates("month")

    out = (
        stock.merge(market, on="month", how="inner")
        .merge(rf, on="month", how="inner")
        .dropna(subset=["permno", "month", "ret", "ret_m", "rf"])
    )
    out = out[out["month"] >= pd.Period(start_month, freq="M")].copy
    out["date"] = out["month"].dt.to_timestamp(how="start")
    out["ret_excess"] = out["ret"] - out["rf"]
    out["mkt_excess"] = out["ret_m"] - out["rf"]
    # 上期市值，供后文按规模分组描述
    out = out.sort_values(["permno", "date"])
    out["mktcap_lag"] = out.groupby("permno")["mktcap"].shift(1)
    return out.reset_index(drop=True)

cache_monthly = CACHE_DIR / "crsp_monthly_ashare.parquet"
if cache_monthly.exists:
    crsp_monthly = pd.read_parquet(cache_monthly)
    print("loaded cache:", cache_monthly)
else:
    crsp_monthly = load_ashare_monthly(DATA_DIR, start_month="2000-01")
    crsp_monthly.to_parquet(cache_monthly, index=False)
    print("saved cache:", cache_monthly)

print(crsp_monthly.shape)
print("stocks:", crsp_monthly["permno"].nunique, "months:", crsp_monthly["date"].nunique)
display(crsp_monthly.head)

loaded cache: D:\A_Topics\202607_03_TidyFinanceAShare\202608_02_传统量化\Data\_cache_beta\crsp_monthly_ashare.parquet
(890537, 10)
stocks: 5977 months: 319

	permno	month	ret	mktcap	ret_m	rf	date	ret_excess	mkt_excess	mktcap_lag
0	000001	2000-01	0.061891	19843822.88	0.160838	0.001856	2000-01-01	0.060035	0.158982	NaN
1	000001	2000-02	-0.011333	19618933.36	0.122024	0.001856	2000-02-01	-0.013189	0.120168	19843822.88
2	000001	2000-03	0.002729	19672478.48	0.055926	0.001856	2000-03-01	0.000873	0.054070	19618933.36
3	000001	2000-04	0.037017	20400692.17	0.013014	0.001856	2000-04-01	0.035161	0.011158	19672478.48
4	000001	2000-05	-0.055118	19276244.57	0.027691	0.001856	2000-05-01	-0.056974	0.025835	20400692.17

1.6.3 单只股票的 CAPM 回归示例

先对一只代表性股票做全样本时序回归，确认回归式与系数含义。这里用 600519（贵州茅台） 对应原文中的 Apple 示例。

example_permno = "600519"

model_fit = OLS(
    crsp_monthly.loc[crsp_monthly["permno"] == example_permno, "ret_excess"],
    add_constant(
        crsp_monthly.loc[crsp_monthly["permno"] == example_permno, "mkt_excess"]
    ),
    missing="drop",
).fit

coefficients = pd.DataFrame(
    {
        "Estimate": model_fit.params,
        "Std. Error": model_fit.bse,
        "t value": model_fit.tvalues,
        "Pr(>|t|)": model_fit.pvalues,
    }
)
coefficients.index = ["Intercept", "mkt_excess"]
display(coefficients)

print(
    "贵州茅台全样本 beta ≈ {:.3f}；alpha ≈ {:.4f}/月".format(
        model_fit.params["mkt_excess"], model_fit.params["const"]
    )
)

	Estimate	Std. Error	t value	Pr(>|t|)
Intercept	0.018621	0.005104	3.648361	3.118157e-04
mkt_excess	0.633779	0.067782	9.350247	2.229768e-18

贵州茅台全样本 beta ≈ 0.634；alpha ≈ 0.0186/月

1.6.4 闭式解滚动窗口估计

对全市场逐股、逐月拟合 CAPM 会很慢。对单回归元模型，滚动 OLS 有闭式解。记 $x$ 为市场超额收益、$y$ 为个股超额收益，在窗口内

$$ \hat\beta_i=\frac{S_{xy}}{S_{xx}},\qquad \hat\alpha_i=\bar y-\hat\beta_i\bar x, $$

其中

$$ S_{xy}=\sum_t x_t y_t-\frac{1}{n}\sum_t x_t\sum_t y_t,\qquad S_{xx}=\sum_t x_t^2-\frac{1}{n}\left(\sum_t x_t\right)^2. $$

因此只需滚动计算 $\sum x,\sum y,\sum x^2,\sum y^2,\sum xy$，即可向量化得到全部股票的滚动 beta 与 $t$ 统计量。默认窗口 60 个月，至少 48 个有效观测（与常见设定一致）。

def roll_capm_estimation(
    data: pd.DataFrame,
    look_back: int = 60,
    min_obs: int = 48,
) -> pd.DataFrame:
    """闭式解滚动 CAPM 估计（结构对齐常见 roll_capm_estimation 实现）。

    Parameters
    ----------
    data : 需含 permno, date, ret_excess, mkt_excess
    look_back : 滚动窗口长度（按 date 聚合后的期数，通常为月）。
    min_obs : 窗口内最少原始观测数。月频时约等于最少月数；
        日频截断到月后，n 为窗口内交易日总数（例如 200），
        不能把它直接当作 pandas rolling 的 min_periods（否则日频
        look_back=12、min_obs=200 会触发 min_periods > window）。
    """
    df = data[["permno", "date", "ret_excess", "mkt_excess"]].dropna.copy
    df["date"] = pd.to_datetime(df["date"], errors="coerce")
    df = df.dropna(subset=["date"])
    df["xx"] = df["mkt_excess"] ** 2
    df["yy"] = df["ret_excess"] ** 2
    df["xy"] = df["mkt_excess"] * df["ret_excess"]

    cumulants = (
        df.groupby(["permno", "date"], as_index=False)
        .agg(
            n=("mkt_excess", "size"),
            sum_x=("mkt_excess", "sum"),
            sum_y=("ret_excess", "sum"),
            sum_xx=("xx", "sum"),
            sum_yy=("yy", "sum"),
            sum_xy=("xy", "sum"),
        )
        .sort_values(["permno", "date"])
    )

    pieces = []
    for permno, g in cumulants.groupby("permno", sort=False):
        g = g.sort_values("date").copy
        # rolling 的 min_periods 是“至少多少个月/期”，与 min_obs（原始观测数）分离
        rolled = g[["n", "sum_x", "sum_y", "sum_xx", "sum_yy", "sum_xy"]].rolling(
            window=look_back, min_periods=1
        ).sum
        est = g[["permno", "date"]].copy
        est = est.join(rolled)
        est = est[est["n"] >= min_obs].copy
        if est.empty:
            continue

        est["s_xx"] = est["sum_xx"] - est["sum_x"] ** 2 / est["n"]
        est["s_xy"] = est["sum_xy"] - est["sum_x"] * est["sum_y"] / est["n"]
        est["s_yy"] = est["sum_yy"] - est["sum_y"] ** 2 / est["n"]
        # 避免窗口内市场收益无变异导致除零
        est = est[est["s_xx"].abs > 1e-18].copy
        if est.empty:
            continue
        est["beta"] = est["s_xy"] / est["s_xx"]
        est["alpha"] = (est["sum_y"] - est["beta"] * est["sum_x"]) / est["n"]
        est["sigma2"] = (est["s_yy"] - est["beta"] * est["s_xy"]) / (est["n"] - 2)
        est.loc[est["sigma2"] <= 0, "sigma2"] = np.nan
        est["t_alpha"] = est["alpha"] / np.sqrt(
            est["sigma2"]
            * (1.0 / est["n"] + (est["sum_x"] / est["n"]) ** 2 / est["s_xx"])
        )
        est["t_beta"] = est["beta"] / np.sqrt(est["sigma2"] / est["s_xx"])
        pieces.append(est)

    if not pieces:
        return pd.DataFrame(
            columns=["permno", "date", "coefficient", "estimate", "t_statistic"]
        )

    estimates = pd.concat(pieces, ignore_index=True)
    alpha = estimates.assign(
        coefficient="alpha", estimate=estimates["alpha"], t_statistic=estimates["t_alpha"]
    )[["permno", "date", "coefficient", "estimate", "t_statistic"]]
    beta = estimates.assign(
        coefficient="mkt_excess",
        estimate=estimates["beta"],
        t_statistic=estimates["t_beta"],
    )[["permno", "date", "coefficient", "estimate", "t_statistic"]]
    out = pd.concat([alpha, beta], ignore_index=True)
    return out.sort_values(["permno", "date", "coefficient"]).reset_index(drop=True)

print("roll_capm_estimation ready")

roll_capm_estimation ready

1.6.5 示例股票的滚动 Beta

先在少数知名股票上验证函数：贵州茅台、宁德时代、中国平安、五粮液。窗口为 60 个月。

examples = pd.DataFrame(
    {
        "permno": ["600519", "300750", "601318", "000858"],
        "company": ["贵州茅台", "宁德时代", "中国平安", "五粮液"],
    }
)

capm_examples = roll_capm_estimation(
    crsp_monthly[crsp_monthly["permno"].isin(examples["permno"])],
    look_back=60,
    min_obs=48,
)
display(capm_examples.head(10))
print("rows:", len(capm_examples))

	permno	date	coefficient	estimate	t_statistic
0	000858	2003-12-01	alpha	0.003234	0.328274
1	000858	2003-12-01	mkt_excess	0.609631	3.699961
2	000858	2004-01-01	alpha	0.003228	0.334763
3	000858	2004-01-01	mkt_excess	0.609517	3.798020
4	000858	2004-02-01	alpha	0.005334	0.550281
5	000858	2004-02-01	mkt_excess	0.650284	4.045008
6	000858	2004-03-01	alpha	0.005596	0.588700
7	000858	2004-03-01	mkt_excess	0.652414	4.108064
8	000858	2004-04-01	alpha	0.006662	0.710796
9	000858	2004-04-01	mkt_excess	0.623492	4.052087

rows: 1504

图 1：示例股票的月频滚动 Beta

只保留 coefficient == "mkt_excess" 的斜率估计，画时间序列。可以看到 beta 随时间变化，且多数时期为正。

beta_examples = (
    capm_examples[capm_examples["coefficient"] == "mkt_excess"]
    .merge(examples, on="permno", how="left")
)

fig, ax = plt.subplots(figsize=(10, 4.5))
for company, g in beta_examples.groupby("company"):
    ax.plot(g["date"], g["estimate"], label=company, linewidth=1.5)
ax.set_title("示例股票的月频滚动 Beta（5 年窗口）")
ax.set_xlabel("")
ax.set_ylabel("beta")
ax.legend(frameon=False, ncol=2)
ax.grid(True, alpha=0.3)
plt.tight_layout
plt.show

1.6.6 全样本月频滚动估计

对全部 A 股月频面板一次性估计。闭式解下通常数十秒内可完成（视机器而定），无需按股票显式循环回归。

cache_capm_m = CACHE_DIR / "capm_monthly_ashare.parquet"
if cache_capm_m.exists:
    capm_monthly = pd.read_parquet(cache_capm_m)
    print("loaded cache:", cache_capm_m)
else:
    capm_monthly = roll_capm_estimation(crsp_monthly, look_back=60, min_obs=48)
    capm_monthly.to_parquet(cache_capm_m, index=False)
    print("saved cache:", cache_capm_m)

print(capm_monthly.shape)
display(capm_monthly.head)

beta_monthly = (
    capm_monthly[capm_monthly["coefficient"] == "mkt_excess"]
    [["permno", "date", "estimate"]]
    .rename(columns={"estimate": "beta"})
    .assign(return_type="monthly")
)
print("beta_monthly:", beta_monthly.shape)

loaded cache: D:\A_Topics\202607_03_TidyFinanceAShare\202608_02_传统量化\Data\_cache_beta\capm_monthly_ashare.parquet
(1246512, 5)

	permno	date	coefficient	estimate	t_statistic
0	000001	2003-12-01	alpha	-0.008803	-0.999836
1	000001	2003-12-01	mkt_excess	0.945177	6.418647
2	000001	2004-01-01	alpha	-0.008221	-0.951793
3	000001	2004-01-01	mkt_excess	0.956954	6.657053
4	000001	2004-02-01	alpha	-0.007131	-0.835241

beta_monthly: (623256, 4)

1.6.7 按市值分组看 Beta 分布

原文按行业画箱线图。本地月频表未附带行业分类，这里改用 上期流通市值五分组 做同样的截面描述：先对每只股票在样本内取平均 beta，再按平均市值分组画箱线。

预期上，小市值组的 beta 离散程度或中枢可能与大市值组不同；同时几乎很少出现长期为负的市场 beta。

size_panel = (
    beta_monthly.merge(
        crsp_monthly[["permno", "date", "mktcap_lag"]],
        on=["permno", "date"],
        how="left",
    )
    .dropna(subset=["beta", "mktcap_lag"])
)

firm_beta = (
    size_panel.groupby("permno", as_index=False)
    .agg(beta=("beta", "mean"), mktcap_lag=("mktcap_lag", "mean"))
)
firm_beta["size_group"] = pd.qcut(
    firm_beta["mktcap_lag"],
    5,
    labels=["Q1 小市值", "Q2", "Q3", "Q4", "Q5 大市值"],
)

fig, ax = plt.subplots(figsize=(8, 4.5))
order = ["Q1 小市值", "Q2", "Q3", "Q4", "Q5 大市值"]
data_box = [firm_beta.loc[firm_beta["size_group"] == k, "beta"].values for k in order]
ax.boxplot(data_box, labels=order, showfliers=False)
ax.axhline(1.0, color="gray", ls="--", lw=1)
ax.set_title("按市值分组的个股平均 Beta 分布")
ax.set_ylabel("平均 beta")
ax.grid(True, axis="y", alpha=0.3)
plt.tight_layout
plt.show

print(firm_beta.groupby("size_group", observed=True)["beta"].describe[["count", "mean", "50%"]])

             count      mean       50%
size_group                            
Q1 小市值      1026.0  1.083252  1.065115
Q2          1025.0  1.136182  1.118821
Q3          1025.0  1.156031  1.135611
Q4          1025.0  1.192340  1.162915
Q5 大市值      1026.0  1.151092  1.129414

<ipython-input-14-8c79d19563c4>:23: MatplotlibDeprecationWarning: The 'labels' parameter of boxplot() has been renamed 'tick_labels' since Matplotlib 3.9; support for the old name will be dropped in 3.11.
  ax.boxplot(data_box, labels=order, showfliers=False)

1.6.8 Beta 截面分位数随时间的变化

每个月计算全市场 beta 的 10%–90% 分位，观察风险暴露的共同起伏与离散度变化。

quantiles = np.arange(0.1, 1.0, 0.1)
rows = []
for dt, g in beta_monthly.groupby("date"):
    qs = g["beta"].quantile(quantiles)
    for q, v in qs.items:
        rows.append({"date": dt, "quantile": int(round(q * 100)), "beta": v})
beta_quantiles = pd.DataFrame(rows)

fig, ax = plt.subplots(figsize=(10, 4.5))
for q, g in beta_quantiles.groupby("quantile"):
    ax.plot(g["date"], g["beta"], label=str(q), linewidth=1.2)
ax.set_title("月频估计 Beta 的截面分位数（10%–90%）")
ax.set_xlabel("")
ax.set_ylabel("beta")
ax.legend(title="分位", ncol=5, frameon=False, fontsize=8)
ax.grid(True, alpha=0.3)
plt.tight_layout
plt.show

1.6.9 使用日频收益估计 Beta

日频数据直接使用本地现成文件，不再自行合成市场收益或把月度无风险摊成日度：

• 个股日收益：AShare_Stocks_Return_Daily_1.xlsx … _7.xlsx（Dretwd）
• 市场日收益：AShare_Market_Return_Daily.xlsx（Markettype=53 的 Cdretwdos，流通市值加权含红利再投资）
• 无风险日利率：AShare_Riskfree_Rate_Daily.xlsx（NRI01 的 Nrrdaydt，百分数转小数）

因个股日频样本约从 2021-08 起偏短，滚动窗口取 12 个月、min_obs=200；若日后日频历史足够长，可将 look_back / min_obs 改回 60 / 1000。将交易日截断到月初后，调用同一个 roll_capm_estimation。

def load_ashare_daily_panel(data_dir: Path) -> pd.DataFrame:
    """读取现成日频个股、市场收益与无风险利率，构造超额收益面板。"""
    files = sorted(data_dir.glob("AShare_Stocks_Return_Daily_*.xlsx"))
    if not files:
        raise FileNotFoundError("未找到 AShare_Stocks_Return_Daily_*.xlsx")

    parts = []
    for f in files:
        raw = read_csmar_excel(f)
        keep = raw[["Stkcd", "Trddt", "Dretwd", "Trdsta"]].copy
        keep.columns = ["permno", "date", "ret", "trdsta"]
        parts.append(keep)
        print("read", f.name, "rows", len(keep))

    daily = pd.concat(parts, ignore_index=True)
    daily["permno"] = daily["permno"].astype(str).str.zfill(6)
    daily["date"] = pd.to_datetime(daily["date"], errors="coerce")
    daily["ret"] = pd.to_numeric(daily["ret"], errors="coerce")
    daily["trdsta"] = pd.to_numeric(daily["trdsta"], errors="coerce")
    daily = daily.dropna(subset=["permno", "date", "ret"])
    # 交易状态 1 = 正常交易（CSMAR 常见口径）
    daily = daily[(daily["trdsta"].isna) | (daily["trdsta"] == 1)].copy

    # 现成日度市场收益：类型 53，流通市值加权含红利再投资
    market_raw = read_csmar_excel(data_dir / "AShare_Market_Return_Daily.xlsx")
    market = market_raw[["Markettype", "Trddt", "Cdretwdos"]].copy
    market.columns = ["market_type", "date", "ret_m"]
    market["market_type"] = pd.to_numeric(market["market_type"], errors="coerce")
    market = market[market["market_type"] == 53].copy
    market["date"] = pd.to_datetime(market["date"], errors="coerce")
    market["ret_m"] = pd.to_numeric(market["ret_m"], errors="coerce")
    market = market[["date", "ret_m"]].dropna.drop_duplicates("date")

    # 现成日度无风险利率（百分数 -> 小数）
    rf_raw = read_csmar_excel(data_dir / "AShare_Riskfree_Rate_Daily.xlsx")
    rf = rf_raw.copy
    if "Nrr1" in rf.columns:
        rf = rf[rf["Nrr1"].astype(str).str.upper.str.contains("NRI01", na=False)]
    rf = rf[["Clsdt", "Nrrdaydt"]].copy
    rf.columns = ["date", "rf_daily"]
    rf["date"] = pd.to_datetime(rf["date"], errors="coerce")
    rf["rf_daily"] = pd.to_numeric(rf["rf_daily"], errors="coerce") / 100.0
    rf = rf.dropna.drop_duplicates("date")

    out = (
        daily.merge(market, on="date", how="inner")
        .merge(rf, on="date", how="inner")
        .dropna(subset=["ret", "ret_m", "rf_daily"])
    )
    out["ret_excess"] = out["ret"] - out["rf_daily"]
    out["mkt_excess"] = out["ret_m"] - out["rf_daily"]
    return out[["permno", "date", "ret_excess", "mkt_excess"]].reset_index(drop=True)

cache_daily = CACHE_DIR / "crsp_daily_ashare.parquet"
# 数据口径已改为现成日度市场/无风险，旧缓存作废
for _stale in [
    cache_daily,
    CACHE_DIR / "capm_daily_ashare.parquet",
]:
    if _stale.exists:
        _stale.unlink
        print("removed old cache:", _stale)

crsp_daily = load_ashare_daily_panel(DATA_DIR)
crsp_daily.to_parquet(cache_daily, index=False)
print("saved cache:", cache_daily)

print(crsp_daily.shape)
print(
    "daily range:",
    crsp_daily["date"].min.date,
    "->",
    crsp_daily["date"].max.date,
    "stocks:",
    crsp_daily["permno"].nunique,
)
display(crsp_daily.head)

removed old cache: D:\A_Topics\202607_03_TidyFinanceAShare\202608_02_传统量化\Data\_cache_beta\crsp_daily_ashare.parquet
removed old cache: D:\A_Topics\202607_03_TidyFinanceAShare\202608_02_传统量化\Data\_cache_beta\capm_daily_ashare.parquet
read AShare_Stocks_Return_Daily_1.xlsx rows 1000000
read AShare_Stocks_Return_Daily_2.xlsx rows 1000000
read AShare_Stocks_Return_Daily_3.xlsx rows 1000000
read AShare_Stocks_Return_Daily_4.xlsx rows 1000000
read AShare_Stocks_Return_Daily_5.xlsx rows 1000000
read AShare_Stocks_Return_Daily_6.xlsx rows 1000000
read AShare_Stocks_Return_Daily_7.xlsx rows 370309
saved cache: D:\A_Topics\202607_03_TidyFinanceAShare\202608_02_传统量化\Data\_cache_beta\crsp_daily_ashare.parquet
(6167307, 4)
daily range: 2021-08-05 -> 2026-08-04 stocks: 5782

	permno	date	ret_excess	mkt_excess
0	000001	2021-08-05	-0.007340	-0.005102
1	000001	2021-08-06	0.010140	-0.001528
2	000001	2021-08-09	0.067148	0.009096
3	000001	2021-08-10	0.035111	0.011574
4	000001	2021-08-11	0.004014	0.001276

# 将日频日期截断到月初，使滚动窗口仍按“月”对齐
# 注意：look_back=12 表示 12 个月；min_obs=200 表示窗口内至少 200 个交易日，
# 由 roll_capm_estimation 对汇总后的 n 过滤，而不是 pandas rolling 的 min_periods。
crsp_daily_month = crsp_daily.copy
crsp_daily_month["date"] = pd.to_datetime(
    crsp_daily_month["date"], errors="coerce"
).dt.to_period("M").dt.to_timestamp
crsp_daily_month = crsp_daily_month.dropna(subset=["date", "ret_excess", "mkt_excess"])

cache_capm_d = CACHE_DIR / "capm_daily_ashare.parquet"
# 若此前因参数错误生成了空缓存，自动重算
need_recompute = True
if cache_capm_d.exists:
    capm_daily = pd.read_parquet(cache_capm_d)
    if len(capm_daily) > 0:
        need_recompute = False
        print("loaded cache:", cache_capm_d)
    else:
        print("empty cache detected, recomputing...")

if need_recompute:
    capm_daily = roll_capm_estimation(
        crsp_daily_month,
        look_back=12,   # 本地日频偏短，使用 12 个月
        min_obs=200,    # 窗口内交易日总数下限
    )
    capm_daily.to_parquet(cache_capm_d, index=False)
    print("saved cache:", cache_capm_d)

print(capm_daily.shape)
display(capm_daily.head)

beta_daily = (
    capm_daily[capm_daily["coefficient"] == "mkt_excess"]
    [["permno", "date", "estimate"]]
    .rename(columns={"estimate": "beta"})
    .assign(return_type="daily")
)
beta = pd.concat([beta_monthly, beta_daily], ignore_index=True)
print("beta stacked:", beta.shape, beta["return_type"].value_counts.to_dict)

saved cache: D:\A_Topics\202607_03_TidyFinanceAShare\202608_02_传统量化\Data\_cache_beta\capm_daily_ashare.parquet
(511308, 5)

	permno	date	coefficient	estimate	t_statistic
0	000001	2022-06-01	alpha	-0.000465	-0.358986
1	000001	2022-06-01	mkt_excess	0.966083	8.728630
2	000001	2022-07-01	alpha	-0.000928	-0.765064
3	000001	2022-07-01	mkt_excess	0.953419	8.980394
4	000001	2022-08-01	alpha	-0.000816	-0.701891

beta stacked: (878910, 4) {'monthly': 623256, 'daily': 255654}

图 4：示例股票的日频 vs 月频 Beta

日频估计通常更平滑（有效样本更多），但趋势应与月频大体一致。因日频样本起点较晚，对比主要落在重叠区间。

# 为公平对比，月频也在同一重叠区间用 12 个月窗口重估一版（仅用于图 4）
capm_examples_m12 = roll_capm_estimation(
    crsp_monthly[crsp_monthly["permno"].isin(examples["permno"])],
    look_back=12,
    min_obs=10,
)
beta_m12 = (
    capm_examples_m12[capm_examples_m12["coefficient"] == "mkt_excess"]
    [["permno", "date", "estimate"]]
    .rename(columns={"estimate": "beta"})
    .assign(return_type="monthly")
)
beta_d12 = beta_daily[beta_daily["permno"].isin(examples["permno"])].copy
beta_comparison = pd.concat([beta_m12, beta_d12], ignore_index=True).merge(
    examples, on="permno", how="inner"
)

companies = examples["company"].tolist
fig, axes = plt.subplots(len(companies), 1, figsize=(10, 8), sharex=True)
if len(companies) == 1:
    axes = [axes]
for ax, company in zip(axes, companies):
    g = beta_comparison[beta_comparison["company"] == company]
    for rt, gg in g.groupby("return_type"):
        ax.plot(gg["date"], gg["beta"], label=rt, linewidth=1.4)
    ax.set_ylabel("beta")
    ax.set_title(company)
    ax.grid(True, alpha=0.3)
    ax.legend(frameon=False, loc="best")
fig.suptitle("日频 vs 月频 Beta（均为 12 个月窗口）", y=1.01)
plt.tight_layout
plt.show

1.6.10 覆盖率与描述统计

最后做两类体检：

1. 每月有多少股票成功估出 beta（覆盖率）；
2. 日频/月频 beta 的分布摘要。

若覆盖率出现大段空洞，通常意味着数据对齐或窗口参数有问题。

# 覆盖率：以月频股票-月份宇宙为分母
universe = crsp_monthly[["permno", "date"]].drop_duplicates
coverage_rows = []
for rt, b in [("monthly", beta_monthly), ("daily", beta_daily)]:
    tmp = universe.merge(
        b[["permno", "date", "beta"]],
        on=["permno", "date"],
        how="left",
    )
    cov = tmp.groupby("date")["beta"].apply(lambda s: s.notna.mean).reset_index
    cov.columns = ["date", "share"]
    cov["return_type"] = rt
    coverage_rows.append(cov)
beta_coverage = pd.concat(coverage_rows, ignore_index=True)

fig, ax = plt.subplots(figsize=(10, 4))
for rt, g in beta_coverage.groupby("return_type"):
    ax.plot(g["date"], g["share"], label=rt, linewidth=1.5)
ax.set_ylim(0, 1)
ax.set_title("月末：具备 Beta 估计的证券占比")
ax.set_ylabel("share")
ax.legend(frameon=False)
ax.grid(True, alpha=0.3)
plt.tight_layout
plt.show

summary = (
    beta.groupby("return_type")["beta"]
    .agg(
        count="count",
        mean="mean",
        std="std",
        min="min",
        q05=lambda s: s.quantile(0.05),
        q50="median",
        q95=lambda s: s.quantile(0.95),
        max="max",
    )
    .round(3)
)
display(summary)

# 持久化，供后续章节（组合排序等）直接读取
out_path = DATA_DIR / "beta_ashare.parquet"
beta.to_parquet(out_path, index=False)
print("wrote", out_path)

	count	mean	std	min	q05	q50	q95	max
return_type
daily	255654	1.228	0.451	-18.801	0.564	1.192	1.992	37.881
monthly	623256	1.118	0.377	-7.364	0.585	1.099	1.707	12.785

wrote D:\A_Topics\202607_03_TidyFinanceAShare\202608_02_传统量化\Data\beta_ashare.parquet

1.6.11 要点小结

1. 市场 beta 来自个股超额收益对市场超额收益的时序回归斜率；滚动窗口可捕捉时变风险暴露。
2. 闭式解滚动估计把“百万次回归”变成滚动求和，是全样本可行的关键。
3. 月频长窗口适合刻画中长期暴露；日频在有效样本更多时通常更平滑，但对数据长度与微观结构噪声更敏感。
4. A 股落地时注意：市场组合口径（本节月频用类型 53）、无风险利率频率转换、ST/停牌（日频用 Trdsta）都会影响估计。
5. 得到的 beta_ashare.parquet 可直接用于后续单维排序、低 beta 异象检验等（参见第 6 章分层回测与单维排序）。