"""M07 教学图：French 五因子与动量在三段窗口里的算术溢价。

五因子月频从 1963-07 起。动量来自 M06 已下载的 Mom。只读月频块。
"""

from __future__ import annotations

import math
import urllib.request
import zipfile
from io import BytesIO
from pathlib import Path

import matplotlib.pyplot as plt
import numpy as np
import pandas as pd

plt.rcParams["font.sans-serif"] = ["Microsoft YaHei", "SimHei", "DejaVu Sans"]
plt.rcParams["axes.unicode_minus"] = False

BASE = "https://mba.tuck.dartmouth.edu/pages/faculty/ken.french/ftp/"
CACHE = Path(__file__).resolve().parent / "_cache"
MOM_ZIP = (
    Path(__file__).resolve().parents[1]
    / "M06"
    / "_cache"
    / "F-F_Momentum_Factor_CSV.zip"
)
OUT = Path(__file__).resolve().parent / "ff5_windows.png"
NAMES = ["Mkt-RF", "SMB", "HML", "RMW", "CMA", "Mom"]


def fetch_ff5() -> bytes:
    path = CACHE / "F-F_Research_Data_5_Factors_2x3_CSV.zip"
    if not path.exists():
        CACHE.mkdir(parents=True, exist_ok=True)
        url = BASE + "F-F_Research_Data_5_Factors_2x3_CSV.zip"
        path.write_bytes(urllib.request.urlopen(url, timeout=60).read())
    return path.read_bytes()


def monthly_block(text: str, marker: str) -> pd.DataFrame:
    lines = text.splitlines()
    start = next(i for i, line in enumerate(lines) if marker in line)
    if lines[start].startswith(","):
        header_line = lines[start]
        data_start = start + 1
    else:
        header_line = lines[start + 1]
        data_start = start + 2
    header = [part.strip() for part in header_line.split(",")]
    rows: list[list[str]] = []
    for line in lines[data_start:]:
        if not line.strip():
            break
        parts = [part.strip() for part in line.split(",")]
        if len(parts[0]) != 6:
            break
        rows.append(parts)
    frame = pd.DataFrame(rows, columns=header)
    frame = frame.rename(columns={header[0]: "date"})
    for column in frame.columns[1:]:
        frame[column] = pd.to_numeric(frame[column]) / 100.0
    frame["date"] = pd.to_datetime(frame["date"], format="%Y%m")
    return frame.set_index("date")


def stats(series: pd.Series) -> dict[str, float]:
    values = series.dropna()
    n = len(values)
    mean = float(values.mean())
    std = float(values.std(ddof=1))
    nav = (1.0 + values).cumprod()
    mdd = float((nav / nav.cummax() - 1.0).min())
    t_stat = mean / (std / math.sqrt(n))
    return {
        "n": float(n),
        "ann": mean * 12.0,
        "sharpe": mean / std * math.sqrt(12.0),
        "mdd": mdd,
        "t": t_stat,
    }


def ols_table(y: np.ndarray, x: np.ndarray) -> tuple[np.ndarray, np.ndarray, float]:
    """经典 OLS。标准误不带 Newey–West，截距的 t 值只供和斜率对照。"""
    design = np.column_stack([np.ones(len(y)), x])
    coef = np.linalg.lstsq(design, y, rcond=None)[0]
    resid = y - design @ coef
    df = len(y) - design.shape[1]
    sigma2 = float(resid @ resid) / df
    cov = sigma2 * np.linalg.inv(design.T @ design)
    se = np.sqrt(np.diag(cov))
    ss_tot = float(((y - y.mean()) ** 2).sum())
    r2 = 1.0 - float(resid @ resid) / ss_tot
    return coef, se, r2


def main() -> None:
    raw = zipfile.ZipFile(BytesIO(fetch_ff5())).read(
        "F-F_Research_Data_5_Factors_2x3.csv"
    )
    ff5 = monthly_block(raw.decode("latin1"), ",Mkt-RF")
    mom_raw = zipfile.ZipFile(MOM_ZIP).read("F-F_Momentum_Factor.csv")
    mom = monthly_block(mom_raw.decode("latin1"), ",Mom")
    data = ff5.join(mom.rename(columns={mom.columns[0]: "Mom"}), how="inner")
    print("sample", data.index.min().date(), data.index.max().date(), len(data))
    windows = {
        "full": data,
        "y1963_2006": data.loc["1963-07":"2006-12"],
        "y2007_2020": data.loc["2007-01":"2020-12"],
        "y2021_end": data.loc["2021-01":],
        "y2016_2025": data.loc["2016-01":"2025-12"],
    }
    for name, frame in windows.items():
        print("WINDOW", name, frame.index.min().date(), frame.index.max().date(), len(frame))
        for column in NAMES:
            row = stats(frame[column])
            print(
                column,
                "ann", f"{row['ann']:.4f}",
                "sharpe", f"{row['sharpe']:.2f}",
                "mdd", f"{row['mdd']:.4f}",
                "t", f"{row['t']:.2f}",
            )
    corr = data[NAMES].corr()
    print("corr_hml_mom", f"{float(corr.loc['HML', 'Mom']):.3f}")
    print("corr_hml_rmw", f"{float(corr.loc['HML', 'RMW']):.3f}")
    print("corr_hml_cma", f"{float(corr.loc['HML', 'CMA']):.3f}")
    y = data["HML"].to_numpy(dtype=float)
    x = data[["Mkt-RF", "SMB", "RMW", "CMA"]].to_numpy(dtype=float)
    coef, se, r2 = ols_table(y, x)
    names = ["alpha", "Mkt-RF", "SMB", "RMW", "CMA"]
    for name, value, err in zip(names, coef, se):
        print("hml_on", name, f"{float(value):.6f}", "t", f"{float(value / err):.2f}")
    print("hml_r2", f"{r2:.3f}", "df", len(y) - len(names))
    hml = (1.0 + data["HML"]).cumprod()
    # 切片含 2006-12。高点若不是这一月，跌幅仍相对该峰值，不改用更早的定义。
    peak = hml.loc[: "2006-12-01"].idxmax()
    trough = (hml / hml.cummax() - 1.0).loc["2006-12-01":].idxmin()
    level_peak = float(hml.loc[peak])
    print(
        "hml_peak", peak.date(),
        "trough", trough.date(),
        "dd", f"{float(hml.loc[trough] / level_peak - 1.0):.4f}",
        "vs_peak_end", f"{float(hml.iloc[-1] / level_peak - 1.0):.4f}",
    )
    year = data.loc["2022"]
    for column in NAMES:
        compounded = float((1.0 + year[column]).prod() - 1.0)
        print("c2022", column, f"{compounded:.4f}")
    mix = 0.5 * data["HML"] + 0.5 * data["Mom"]
    print(
        "vol_hml", f"{float(data['HML'].std(ddof=1) * math.sqrt(12)):.4f}",
        "vol_mom", f"{float(data['Mom'].std(ddof=1) * math.sqrt(12)):.4f}",
        "vol_mix", f"{float(mix.std(ddof=1) * math.sqrt(12)):.4f}",
    )
    for label, series in (("HML", data["HML"]), ("Mom", data["Mom"]), ("mix", mix)):
        row = stats(series)
        print(label, "mixann", f"{row['ann']:.4f}", "sharpe", f"{row['sharpe']:.2f}")

    order = ["y1963_2006", "y2007_2020", "y2021_end"]
    labels = ["1963–2006", "2007–2020", "2021–2026.08"]
    fig, axes = plt.subplots(1, 2, figsize=(11.2, 4.8))
    ax = axes[0]
    xpos = np.arange(len(NAMES))
    width = 0.25
    colors = ["#1f4e79", "#c1121f", "#2a9d8f"]
    for offset, key, label, color in zip((-width, 0.0, width), order, labels, colors):
        ann = [stats(windows[key][column])["ann"] * 100 for column in NAMES]
        ax.bar(xpos + offset, ann, width=width, color=color, label=label)
    ax.axhline(0.0, color="#333333", lw=0.6)
    ax.set_xticks(xpos)
    ax.set_xticklabels(["市场", "规模", "价值", "盈利", "投资", "动量"], fontsize=8)
    ax.set_ylabel("算术年化（%）")
    ax.legend(frameon=False, fontsize=7)
    ax.grid(True, axis="y", alpha=0.3)
    ax = axes[1]
    ax.plot((1.0 + data["HML"]).cumprod(), color="#1f4e79", lw=1.3, label="HML")
    ax.plot((1.0 + data["Mom"]).cumprod(), color="#c1121f", lw=1.2, label="Mom")
    ax.plot((1.0 + mix).cumprod(), color="#2a9d8f", lw=1.4, label="各一半")
    ax.set_yscale("log")
    ax.set_ylabel("累计倍数（对数刻度）")
    ax.legend(frameon=False, fontsize=8)
    ax.grid(True, alpha=0.3)
    fig.suptitle("M07 French 五因子与 Mom（月频，未扣成本）", fontsize=12)
    fig.tight_layout()
    fig.savefig(OUT, dpi=140, bbox_inches="tight")
    print("wrote", OUT)


if __name__ == "__main__":
    main()