バックテストの過剰最適化を数値で検出する
Walk-Forward・Purged K-Fold・
Deflated Sharpe RatioのPython実装

"""overfit_demo.py — 期待値ゼロの戦略から高Sharpeを「発見」する実演
本コードは本人利用に限り自由に使用可。再配布・商用利用・販売は禁止。

目的: 「試行回数を増やせば、無意味な戦略でも良い成績が出る」
      ことを自分の手で確認する。
"""
from __future__ import annotations

import numpy as np


def annualised_sharpe(returns: np.ndarray, periods_per_year: int = 365) -> float:
    if returns.size < 2:
        return 0.0
    std = returns.std(ddof=1)
    if std == 0:
        return 0.0
    return float(returns.mean() / std * np.sqrt(periods_per_year))


def best_of_random(n_trials: int, n_days: int = 500,
                   seed: int = 42) -> tuple[float, np.ndarray]:
    """完全にランダムなリターン列をn_trials本作り、最良のSharpeを返す。

    真の期待値はゼロ。それでも試行回数を増やすとSharpeは上がる。
    """
    rng = np.random.default_rng(seed)
    # 日次1%のボラティリティ、期待値ちょうど0
    paths = rng.normal(loc=0.0, scale=0.01, size=(n_trials, n_days))
    sharpes = np.array([annualised_sharpe(p) for p in paths])
    best_idx = int(np.argmax(sharpes))
    return float(sharpes[best_idx]), sharpes


if __name__ == "__main__":
    print("真の期待値が0（=エッジなし）の戦略を、N回試して最良を選ぶ\n")
    print(f"{'試行回数':>8}{'最良Sharpe':>13}{'上位5%':>11}"
          f"{'中央値':>10}{'最悪':>10}")
    print("-" * 52)

    for n in (1, 10, 100, 1_000, 10_000):
        best, all_sharpes = best_of_random(n, n_days=500)
        print(f"{n:>8}{best:>13.2f}"
              f"{np.percentile(all_sharpes, 95):>11.2f}"
              f"{np.median(all_sharpes):>10.2f}"
              f"{all_sharpes.min():>10.2f}")

    print("\n--- 結論 ---")
    print("  エッジが一切なくても、1000回試せばSharpe 2.0前後が出る。")
    print("  つまり『Sharpe 2.0を見つけた』という報告には、")
    print("  『何回試したか』の情報が無いと意味がない。")
    print("\n  さらに悪いことに、この最良パスは見た目も良い:")

    best, sharpes = best_of_random(1_000, n_days=500)
    rng = np.random.default_rng(42)
    paths = rng.normal(0.0, 0.01, size=(1_000, 500))
    best_path = paths[int(np.argmax([annualised_sharpe(p) for p in paths]))]
    equity = np.cumprod(1 + best_path)

    print(f"    最終リターン : {(equity[-1] - 1) * 100:+.1f}%")
    print(f"    Sharpe       : {best:.2f}")
    peak = np.maximum.accumulate(equity)
    print(f"    最大DD       : {((equity / peak - 1).min()) * 100:.1f}%")
    print(f"    勝率(日次)   : {(best_path > 0).mean():.1%}")
    print("\n  これを見せられたら『良い戦略だ』と思ってしまう。")
    print("  しかし中身は乱数であり、将来の期待値はゼロ。")

時系列:  |---------------------- 全期間 ----------------------|

【通常のK-Fold（誤り）】
  訓練:  ####    ####################        ####
  テスト:     ####                    ####
         ^^^^ 未来のデータで学習して過去をテストしている

【Purged K-Fold（正しい）】
  訓練:  ######## [purge] ....テスト.... [embargo] ##########
                    ^^^^^^^              ^^^^^^^^^
                    ラベル期間が          系列相関の影響を
                    重なる部分を除去      避けるための空白

purge  : テスト期間のラベルが参照する未来と重なる訓練データを除く
embargo: テスト期間の直後、一定期間の訓練データを除く

"""purged_cv.py — 時系列リークを防ぐPurged K-Fold（embargo付き）
本コードは本人利用に限り自由に使用可。再配布・商用利用・販売は禁止。

参考: López de Prado (2018) "Advances in Financial Machine Learning" Ch.7

scikit-learnのKFold互換のインターフェースを持つので、
GridSearchCVなどにそのまま渡せる。
"""
from __future__ import annotations

from collections.abc import Iterator

import numpy as np
import pandas as pd


class PurgedKFold:
    """時系列用のK-Fold。テスト期間と重なる訓練データを除去する。

    n_splits    : 分割数
    t1          : 各サンプルのラベル終了時刻（Series: index=開始, value=終了）
                  Noneなら各サンプルは単一時点とみなす
    embargo_pct : テスト期間の直後に置く空白の割合（全期間に対する）
    """

    def __init__(self, n_splits: int = 5, t1: pd.Series | None = None,
                 embargo_pct: float = 0.01):
        if n_splits < 2:
            raise ValueError("n_splits must be >= 2")
        self.n_splits = n_splits
        self.t1 = t1
        self.embargo_pct = embargo_pct

    def get_n_splits(self, X=None, y=None, groups=None) -> int:
        return self.n_splits

    def split(self, X: pd.DataFrame, y=None, groups=None
              ) -> Iterator[tuple[np.ndarray, np.ndarray]]:
        if not isinstance(X.index, pd.DatetimeIndex):
            raise TypeError("XのindexはDatetimeIndexである必要があります")

        indices = np.arange(X.shape[0])
        embargo_size = int(X.shape[0] * self.embargo_pct)

        # ラベル終了時刻。指定が無ければサンプル自身の時刻
        t1 = self.t1 if self.t1 is not None else pd.Series(X.index,
                                                           index=X.index)
        t1 = t1.reindex(X.index).ffill()

        # 連続したブロックに分割（シャッフルしない）
        test_ranges = [(i[0], i[-1] + 1)
                       for i in np.array_split(indices, self.n_splits)]

        for start, end in test_ranges:
            test_idx = indices[start:end]
            test_start_time = X.index[start]
            test_end_time = t1.iloc[end - 1]

            # --- purge: テスト期間とラベルが重なる訓練データを除去 -------
            # 訓練サンプルのラベルが test_start_time 以降に及ぶものを除く
            train_mask = np.ones(X.shape[0], dtype=bool)
            train_mask[start:end] = False

            label_end = t1.to_numpy()
            overlap = label_end >= np.datetime64(test_start_time)
            before_test = indices < start
            train_mask[before_test & overlap] = False

            # --- embargo: テスト直後の一定期間を除去 ---------------------
            if embargo_size > 0:
                emb_end = min(end + embargo_size, X.shape[0])
                train_mask[end:emb_end] = False

            train_idx = indices[train_mask]
            if train_idx.size == 0:
                continue
            yield train_idx, test_idx


def make_t1(index: pd.DatetimeIndex, horizon_bars: int) -> pd.Series:
    """「N本先のリターンをラベルにする」場合のラベル終了時刻を作る。

    horizon_bars=5 なら、各サンプルのラベルは5本先まで参照する。
    このt1をPurgedKFoldに渡すと、正しくpurgeされる。
    """
    end_positions = np.minimum(np.arange(len(index)) + horizon_bars,
                               len(index) - 1)
    return pd.Series(index[end_positions], index=index)


def visualise_splits(X: pd.DataFrame, cv, max_show: int = 5) -> None:
    """分割の様子を文字で可視化する（実装の確認用）。"""
    n = X.shape[0]
    width = 80
    scale = width / n

    print(f"サンプル数 {n} / 分割 {cv.get_n_splits()}"
          f" / embargo {cv.embargo_pct:.1%}")
    print("  凡例: '#'=訓練  '.'=テスト  ' '=除外(purge/embargo)\n")

    for k, (tr, te) in enumerate(cv.split(X)):
        if k >= max_show:
            break
        line = [" "] * width
        for i in tr:
            pos = int(i * scale)
            if pos < width:
                line[pos] = "#"
        for i in te:
            pos = int(i * scale)
            if pos < width:
                line[pos] = "."
        excluded = n - len(tr) - len(te)
        print(f"  fold{k}: {''.join(line)}  "
              f"訓練{len(tr):>5} テスト{len(te):>4} 除外{excluded:>4}")


if __name__ == "__main__":
    # 1000本の1時間足を模したデータ
    idx = pd.date_range("2025-01-01", periods=1000, freq="1h", tz="UTC")
    X = pd.DataFrame({"feature": np.random.default_rng(0).normal(size=1000)},
                     index=idx)

    print("=== embargoなし・ラベル期間なし（通常の時系列分割） ===")
    visualise_splits(X, PurgedKFold(n_splits=5, embargo_pct=0.0))

    print("\n=== embargo 2% ===")
    visualise_splits(X, PurgedKFold(n_splits=5, embargo_pct=0.02))

    print("\n=== ラベルが24本先まで及ぶ場合（purgeが効く） ===")
    t1 = make_t1(idx, horizon_bars=24)
    visualise_splits(X, PurgedKFold(n_splits=5, t1=t1, embargo_pct=0.02))

    print("\n※ 除外数が増えているのがpurge/embargoの効果。")
    print("  これを省くと、テスト成績が楽観的になる（リークするため）。")

"""cpcv.py — Combinatorial Purged Cross-Validation
本コードは本人利用に限り自由に使用可。再配布・商用利用・販売は禁止。

参考: López de Prado (2018) Ch.12

N個のブロックからk個をテストに選ぶ全組み合わせを試すことで、
複数のOOS経路（バックテストパス）を生成する。
成績の「分布」が得られるため、1本の経路より判断が頑健になる。
"""
from __future__ import annotations

from collections.abc import Iterator
from dataclasses import dataclass
from itertools import combinations
from math import comb

import numpy as np
import pandas as pd


@dataclass
class CpcvSplit:
    train_idx: np.ndarray
    test_idx: np.ndarray
    test_groups: tuple[int, ...]


class CPCV:
    """組み合わせ的なPurged CV。

    n_groups : 全期間を分割するブロック数
    k_test   : テストに使うブロック数
    -> 分割数は C(n_groups, k_test)、OOS経路数は
       C(n_groups, k_test) * k_test / n_groups
    """

    def __init__(self, n_groups: int = 6, k_test: int = 2,
                 embargo_pct: float = 0.01, t1: pd.Series | None = None):
        if k_test >= n_groups:
            raise ValueError("k_test < n_groups である必要があります")
        self.n_groups = n_groups
        self.k_test = k_test
        self.embargo_pct = embargo_pct
        self.t1 = t1

    @property
    def n_splits(self) -> int:
        return comb(self.n_groups, self.k_test)

    @property
    def n_paths(self) -> int:
        """生成されるOOS経路の本数。"""
        return self.n_splits * self.k_test // self.n_groups

    def split(self, X: pd.DataFrame) -> Iterator[CpcvSplit]:
        n = X.shape[0]
        indices = np.arange(n)
        groups = np.array_split(indices, self.n_groups)
        embargo = int(n * self.embargo_pct)

        t1 = self.t1 if self.t1 is not None else pd.Series(X.index,
                                                           index=X.index)
        t1 = t1.reindex(X.index).ffill()
        label_end = t1.to_numpy()

        for combo in combinations(range(self.n_groups), self.k_test):
            test_idx = np.concatenate([groups[g] for g in combo])
            test_idx.sort()

            train_mask = np.ones(n, dtype=bool)
            train_mask[test_idx] = False

            # 各テストブロックについて purge と embargo を適用
            for g in combo:
                blk = groups[g]
                blk_start, blk_end = blk[0], blk[-1] + 1
                blk_start_time = np.datetime64(X.index[blk_start])

                # purge: ラベルがテスト期間に及ぶ訓練サンプルを除去
                overlap = (label_end >= blk_start_time) & (indices < blk_start)
                train_mask[overlap] = False

                # embargo: テストブロック直後を除去
                if embargo > 0:
                    train_mask[blk_end:min(blk_end + embargo, n)] = False

            train_idx = indices[train_mask]
            if train_idx.size == 0:
                continue
            yield CpcvSplit(train_idx=train_idx, test_idx=test_idx,
                            test_groups=combo)


def cpcv_paths(cv: CPCV, X: pd.DataFrame,
               oos_returns_by_split: list[pd.Series]) -> list[pd.Series]:
    """各分割のOOSリターンを、連続した経路に組み替える。

    同じブロックが複数の分割でテストになるため、
    それらを組み合わせて「全期間をカバーするOOS経路」を複数作れる。
    """
    n = X.shape[0]
    groups = np.array_split(np.arange(n), cv.n_groups)

    # group -> [そのgroupをテストにした分割のリターン, ...]
    by_group: dict[int, list[pd.Series]] = {g: [] for g in range(cv.n_groups)}
    for split, rets in zip(cv.split(X), oos_returns_by_split):
        for g in split.test_groups:
            idx = X.index[groups[g]]
            by_group[g].append(rets.reindex(idx).dropna())

    n_paths = min(len(v) for v in by_group.values()) if by_group else 0
    paths: list[pd.Series] = []
    for p in range(n_paths):
        pieces = [by_group[g][p] for g in range(cv.n_groups)
                  if len(by_group[g]) > p]
        if pieces:
            paths.append(pd.concat(pieces).sort_index())
    return paths


if __name__ == "__main__":
    idx = pd.date_range("2024-01-01", periods=2000, freq="1h", tz="UTC")
    X = pd.DataFrame({"f": np.random.default_rng(1).normal(size=2000)},
                     index=idx)

    print(f"{'n_groups':>9}{'k_test':>8}{'分割数':>8}{'OOS経路数':>11}")
    print("-" * 36)
    for ng, kt in ((5, 1), (6, 2), (8, 2), (10, 2), (10, 3)):
        cv = CPCV(n_groups=ng, k_test=kt)
        print(f"{ng:>9}{kt:>8}{cv.n_splits:>8}{cv.n_paths:>11}")

    print("\n=== CPCV(6, 2) の分割 ===")
    cv = CPCV(n_groups=6, k_test=2, embargo_pct=0.01)
    for i, s in enumerate(cv.split(X)):
        if i >= 6:
            print(f"  ... 他 {cv.n_splits - 6} 分割")
            break
        print(f"  分割{i}: テストブロック {s.test_groups}  "
              f"訓練{len(s.train_idx):>5} テスト{len(s.test_idx):>5}")

    print(f"\n※ 通常のWalk-ForwardではOOS経路が1本だが、")
    print(f"  CPCV(6,2)なら{cv.n_paths}本得られる。")
    print("  成績の分布が見えるため、「1本の経路が良かった」に騙されにくい。")

"""walk_forward.py — Walk-Forward最適化と効率比（WFE）
本コードは本人利用に限り自由に使用可。再配布・商用利用・販売は禁止。

考え方:
  1) IS（インサンプル）でパラメータを最適化
  2) 次の期間（OOS）でそのパラメータを検証
  3) 窓をずらして繰り返す
  4) OOS成績 / IS成績 = 効率比（WFE）。低いほど過剰最適化
"""
from __future__ import annotations

from collections.abc import Callable, Sequence
from dataclasses import dataclass, field

import numpy as np
import pandas as pd


@dataclass
class WfaWindow:
    fold: int
    is_start: pd.Timestamp
    is_end: pd.Timestamp
    oos_start: pd.Timestamp
    oos_end: pd.Timestamp
    best_params: dict
    is_score: float
    oos_score: float
    oos_returns: pd.Series = field(default_factory=lambda: pd.Series(dtype=float))

    @property
    def efficiency(self) -> float:
        """OOS / IS。1.0なら完全に再現、0以下なら機能していない。"""
        if self.is_score == 0:
            return 0.0
        return self.oos_score / self.is_score

    def line(self) -> str:
        return (f"  fold{self.fold:<3}"
                f"IS {self.is_start:%Y-%m-%d}〜{self.is_end:%m-%d}  "
                f"OOS {self.oos_start:%m-%d}〜{self.oos_end:%m-%d}  "
                f"IS {self.is_score:>+6.2f}  OOS {self.oos_score:>+6.2f}  "
                f"効率 {self.efficiency:>+6.2f}  {self.best_params}")


@dataclass
class WfaResult:
    windows: list[WfaWindow]

    @property
    def mean_efficiency(self) -> float:
        effs = [w.efficiency for w in self.windows]
        return float(np.mean(effs)) if effs else 0.0

    @property
    def oos_positive_ratio(self) -> float:
        if not self.windows:
            return 0.0
        return float(np.mean([w.oos_score > 0 for w in self.windows]))

    @property
    def param_stability(self) -> float:
        """選ばれたパラメータの安定性。同じ値が選ばれる割合。"""
        if len(self.windows) < 2:
            return 1.0
        keys = [tuple(sorted(w.best_params.items())) for w in self.windows]
        most_common = max(set(keys), key=keys.count)
        return keys.count(most_common) / len(keys)

    def stitched_oos(self) -> pd.Series:
        """全OOS期間をつないだリターン列。これが「本当の成績」。"""
        parts = [w.oos_returns for w in self.windows
                 if not w.oos_returns.empty]
        if not parts:
            return pd.Series(dtype=float)
        return pd.concat(parts).sort_index()

    def report(self) -> str:
        lines = ["=== Walk-Forward 結果 ==="]
        lines += [w.line() for w in self.windows]

        is_mean = np.mean([w.is_score for w in self.windows])
        oos_mean = np.mean([w.oos_score for w in self.windows])

        verdict = ("過剰最適化の疑いが強い" if self.mean_efficiency < 0.3 else
                   "やや過剰最適化" if self.mean_efficiency < 0.5 else
                   "許容範囲" if self.mean_efficiency < 0.8 else
                   "良好（IS/OOSの差が小さい）")

        lines += [
            "",
            f"  IS平均スコア      : {is_mean:>+7.3f}",
            f"  OOS平均スコア     : {oos_mean:>+7.3f}",
            f"  平均効率比(WFE)   : {self.mean_efficiency:>+7.3f}  <- {verdict}",
            f"  OOSがプラスの割合 : {self.oos_positive_ratio:>7.1%}",
            f"  パラメータ安定性  : {self.param_stability:>7.1%}",
        ]
        return "\n".join(lines)


def walk_forward(
    data: pd.DataFrame,
    param_grid: Sequence[dict],
    evaluate: Callable[[pd.DataFrame, dict], tuple[float, pd.Series]],
    is_days: int = 180,
    oos_days: int = 30,
    anchored: bool = False,
) -> WfaResult:
    """Walk-Forward最適化を実行する。

    evaluate(data, params) -> (スコア, リターン列) を返す関数を渡す。
    anchored=True なら IS の開始を固定（アンカード）、
    False ならローリング窓。
    """
    windows: list[WfaWindow] = []
    start = data.index[0]
    end = data.index[-1]

    is_delta = pd.Timedelta(days=is_days)
    oos_delta = pd.Timedelta(days=oos_days)

    fold = 0
    cursor = start
    while cursor + is_delta + oos_delta <= end:
        is_start = start if anchored else cursor
        is_end = cursor + is_delta
        oos_start = is_end
        oos_end = min(oos_start + oos_delta, end)

        is_data = data.loc[is_start:is_end]
        oos_data = data.loc[oos_start:oos_end]

        if len(is_data) < 30 or len(oos_data) < 5:
            cursor += oos_delta
            continue

        # ISで最良パラメータを探す
        best_score, best_params = -np.inf, param_grid[0]
        for params in param_grid:
            score, _ = evaluate(is_data, params)
            if score > best_score:
                best_score, best_params = score, params

        # OOSで検証（パラメータは固定）
        oos_score, oos_rets = evaluate(oos_data, best_params)

        windows.append(WfaWindow(
            fold=fold, is_start=is_start, is_end=is_end,
            oos_start=oos_start, oos_end=oos_end,
            best_params=best_params, is_score=best_score,
            oos_score=oos_score, oos_returns=oos_rets,
        ))
        fold += 1
        cursor += oos_delta

    return WfaResult(windows=windows)


# ---- 使用例: SMAクロスのWalk-Forward ------------------------------------
def sma_evaluate(bars: pd.DataFrame, params: dict) -> tuple[float, pd.Series]:
    """SMAクロス戦略のSharpeとリターン列を返す。

    ルックアヘッドを避けるため、シグナルを1本シフトしている
    （足Nの終値で判断し、足N+1で建てる）。
    """
    fast, slow = params["fast"], params["slow"]
    if fast >= slow or len(bars) < slow + 2:
        return -np.inf, pd.Series(dtype=float)

    close = bars["close"]
    sig = (close.rolling(fast).mean() > close.rolling(slow).mean()).astype(int)
    position = sig.shift(1).fillna(0)          # 1本シフト = 未来を見ない

    rets = close.pct_change().fillna(0) * position
    # 売買コスト（往復10bps相当）をポジション変化時に課す
    turnover = position.diff().abs().fillna(0)
    rets = rets - turnover * 0.0010

    std = rets.std(ddof=1)
    if std == 0 or np.isnan(std):
        return -np.inf, rets
    sharpe = float(rets.mean() / std * np.sqrt(24 * 365))
    return sharpe, rets


if __name__ == "__main__":
    from kline_store import KlineStore

    bars = KlineStore().load("BTCUSDT", "60")
    if bars.empty:
        raise SystemExit("先に記事07でデータを取得してください")
    bars = bars.loc["2024-01-01":]

    grid = [{"fast": f, "slow": s}
            for f in (5, 10, 20, 30, 50)
            for s in (40, 60, 100, 150, 200) if f < s]
    print(f"パラメータ候補: {len(grid)}通り\n")

    for label, anchored in (("ローリング窓", False), ("アンカード", True)):
        print(f"########## {label} ##########")
        res = walk_forward(bars, grid, sma_evaluate,
                           is_days=180, oos_days=30, anchored=anchored)
        print(res.report())
        print()

【1】期待される最大Sharpe（帰無仮説のもとで）
     N回試したとき、真の期待値がゼロでも最大Sharpeはこれくらいになる:

     SR* = sqrt(Var[SR_n]) * [ (1-γ)*Z^-1(1 - 1/N) + γ*Z^-1(1 - 1/(N*e)) ]

       γ = 0.5772...（オイラー・マスケローニ定数）
       Var[SR_n] = 試行間のSharpeの分散
       Z^-1 = 標準正規分布の逆累積分布関数

【2】Probabilistic Sharpe Ratio（PSR）でSR*を超える確率を出す

     PSR(SR*) = Z[ (SR_hat - SR*) * sqrt(T-1)
                   / sqrt(1 - γ3*SR_hat + (γ4-1)/4 * SR_hat^2) ]

       SR_hat = 観測されたSharpe（非年率化）
       T  = 観測数
       γ3 = 歪度（skewness）
       γ4 = 尖度（kurtosis, 非超過）

【3】DSR = PSR(SR*)

     DSR > 0.95 なら「偶然とは言いにくい」
     DSR < 0.5  なら「偶然の範囲」

"""deflated_sharpe.py — Probabilistic / Deflated Sharpe Ratio
本コードは本人利用に限り自由に使用可。再配布・商用利用・販売は禁止。

参考: Bailey, D. & Lopez de Prado, M. (2014)
      "The Deflated Sharpe Ratio: Correcting for Selection Bias,
       Backtest Overfitting and Non-Normality"
"""
from __future__ import annotations

from dataclasses import dataclass

import numpy as np
from scipy import stats

EULER_GAMMA = 0.5772156649015329


@dataclass
class SharpeStats:
    sharpe: float            # 非年率化（1期間あたり）
    sharpe_annual: float
    n_obs: int
    skew: float
    kurtosis: float          # 非超過（正規分布なら3）
    periods_per_year: int

    def line(self) -> str:
        return (f"Sharpe(期間) {self.sharpe:.4f} / "
                f"年率 {self.sharpe_annual:.3f} / "
                f"n={self.n_obs} / 歪度 {self.skew:+.3f} / "
                f"尖度 {self.kurtosis:.3f}")


def sharpe_stats(returns: np.ndarray | list[float],
                 periods_per_year: int = 365) -> SharpeStats:
    r = np.asarray(returns, dtype=float)
    r = r[~np.isnan(r)]
    if r.size < 3:
        raise ValueError("観測数が不足しています")

    std = r.std(ddof=1)
    sr = float(r.mean() / std) if std > 0 else 0.0

    return SharpeStats(
        sharpe=sr,
        sharpe_annual=sr * float(np.sqrt(periods_per_year)),
        n_obs=r.size,
        skew=float(stats.skew(r, bias=False)),
        kurtosis=float(stats.kurtosis(r, fisher=False, bias=False)),
        periods_per_year=periods_per_year,
    )


def probabilistic_sharpe_ratio(st: SharpeStats,
                               benchmark_sr: float = 0.0) -> float:
    """観測されたSharpeが benchmark_sr を上回る確率。

    非正規性（歪度・尖度）を織り込んだ標準誤差を使う。
    benchmark_sr は「1期間あたり」のSharpeで渡す。
    """
    sr = st.sharpe
    t = st.n_obs

    # Sharpeの標準誤差（Lo 2002 を非正規性に拡張した形）
    denom = 1.0 - st.skew * sr + (st.kurtosis - 1.0) / 4.0 * sr ** 2
    if denom <= 0:
        return float("nan")

    z = (sr - benchmark_sr) * np.sqrt(t - 1) / np.sqrt(denom)
    return float(stats.norm.cdf(z))


def expected_max_sharpe(n_trials: int, sharpe_variance: float) -> float:
    """N回試したときに期待される最大Sharpe（帰無仮説: 真の期待値ゼロ）。

    sharpe_variance は「試行間のSharpeの分散」（1期間あたりの単位）。
    """
    if n_trials < 2:
        return 0.0
    sd = np.sqrt(sharpe_variance)
    a = stats.norm.ppf(1.0 - 1.0 / n_trials)
    b = stats.norm.ppf(1.0 - 1.0 / (n_trials * np.e))
    return float(sd * ((1.0 - EULER_GAMMA) * a + EULER_GAMMA * b))


@dataclass
class DsrResult:
    dsr: float
    psr_vs_zero: float
    observed_sharpe_annual: float
    threshold_sharpe_annual: float   # SR*（年率換算）
    n_trials: int
    n_obs: int

    @property
    def verdict(self) -> str:
        if np.isnan(self.dsr):
            return "計算不能（分布が異常）"
        if self.dsr >= 0.95:
            return "偶然とは言いにくい（採用を検討できる）"
        if self.dsr >= 0.80:
            return "やや有望だが確信は持てない"
        if self.dsr >= 0.50:
            return "偶然の範囲を出ていない"
        return "ほぼ偶然（採用すべきでない）"

    def report(self) -> str:
        return "\n".join([
            f"  観測Sharpe(年率)     : {self.observed_sharpe_annual:>8.3f}",
            f"  必要水準SR*(年率)    : {self.threshold_sharpe_annual:>8.3f}"
            f"  <- {self.n_trials}回試行の帰無仮説",
            f"  PSR(vs 0)            : {self.psr_vs_zero:>8.3f}",
            f"  Deflated Sharpe(DSR) : {self.dsr:>8.3f}",
            f"  判定                 : {self.verdict}",
        ])


def deflated_sharpe_ratio(returns: np.ndarray | list[float],
                          n_trials: int,
                          trial_sharpes: np.ndarray | list[float] | None = None,
                          periods_per_year: int = 365) -> DsrResult:
    """Deflated Sharpe Ratio を計算する。

    n_trials      : 試したパラメータ・戦略の総数
    trial_sharpes : 各試行のSharpe（年率）。渡せば分散を実測できる。
                    渡さない場合は保守的な仮定を使う。
    """
    st = sharpe_stats(returns, periods_per_year)

    if trial_sharpes is not None and len(trial_sharpes) > 1:
        annual = np.asarray(trial_sharpes, dtype=float)
        # 年率Sharpeを1期間あたりに戻して分散を取る
        per_period = annual / np.sqrt(periods_per_year)
        var_sr = float(np.var(per_period, ddof=1))
    else:
        # 実測が無い場合の保守的な仮定:
        # 帰無仮説のもとでSharpeの標準誤差は約 1/sqrt(T)
        var_sr = 1.0 / st.n_obs

    sr_star = expected_max_sharpe(n_trials, var_sr)

    return DsrResult(
        dsr=probabilistic_sharpe_ratio(st, benchmark_sr=sr_star),
        psr_vs_zero=probabilistic_sharpe_ratio(st, benchmark_sr=0.0),
        observed_sharpe_annual=st.sharpe_annual,
        threshold_sharpe_annual=sr_star * float(np.sqrt(periods_per_year)),
        n_trials=n_trials,
        n_obs=st.n_obs,
    )


if __name__ == "__main__":
    rng = np.random.default_rng(7)

    print("=== 同じSharpeでも、試行回数で判定が変わる ===\n")

    # 年率Sharpe 1.5相当のリターン列を作る
    n_days = 500
    target_annual = 1.5
    daily_sr = target_annual / np.sqrt(365)
    rets = rng.normal(loc=daily_sr * 0.01, scale=0.01, size=n_days)

    st = sharpe_stats(rets)
    print(st.line())
    print()

    for n_trials in (1, 10, 100, 1_000, 10_000):
        res = deflated_sharpe_ratio(rets, n_trials=n_trials)
        print(f"--- 試行回数 {n_trials} ---")
        print(res.report())
        print()

    print("=== 実際の試行Sharpeの分布を渡す場合（より正確） ===")
    # グリッドサーチで得た全パラメータのSharpeを渡す想定
    trial_sharpes = rng.normal(0.2, 0.7, size=200)
    trial_sharpes[0] = 1.5          # 最良の1つ
    res = deflated_sharpe_ratio(rets, n_trials=200,
                                trial_sharpes=trial_sharpes)
    print(res.report())

    print("\n--- 実務での使い方 ---")
    print("  グリッドサーチの全結果（各パラメータのSharpe）を保存し、")
    print("  最良のリターン列と試行回数・Sharpe分布をDSRに渡す。")
    print("  DSR < 0.95 なら、その最良パラメータは採用すべきでない。")

"""pbo.py — Probability of Backtest Overfitting (CSCV法)
本コードは本人利用に限り自由に使用可。再配布・商用利用・販売は禁止。

参考: Bailey, Borwein, Lopez de Prado & Zhu (2017)
      "The Probability of Backtest Overfitting"

手続き:
  1) 各パラメータのリターン列を並べた行列 M (T x N) を作る
  2) 時間軸をS個のブロックに分割
  3) S/2個をIS、残りをOOSとする全組み合わせを試す
  4) 各組み合わせで「ISの最良パラメータ」を選び、そのOOS順位を見る
  5) OOS順位が中央値より下だった割合 = PBO
"""
from __future__ import annotations

from dataclasses import dataclass, field
from itertools import combinations
from math import comb

import numpy as np
import pandas as pd


@dataclass
class PboResult:
    pbo: float
    n_combinations: int
    n_params: int
    n_blocks: int
    logits: np.ndarray = field(default_factory=lambda: np.array([]))
    is_ranks: list[int] = field(default_factory=list)
    oos_ranks: list[int] = field(default_factory=list)
    is_scores: list[float] = field(default_factory=list)
    oos_scores: list[float] = field(default_factory=list)

    @property
    def verdict(self) -> str:
        if self.pbo <= 0.2:
            return "過剰最適化の兆候は小さい"
        if self.pbo <= 0.4:
            return "やや過剰最適化の疑い"
        if self.pbo <= 0.6:
            return "過剰最適化の疑いが強い"
        return "ほぼ過剰最適化（選択が機能していない）"

    @property
    def performance_degradation(self) -> float:
        """IS -> OOS のスコア劣化（回帰の傾き）。1に近ければ再現している。"""
        if len(self.is_scores) < 3:
            return float("nan")
        x = np.asarray(self.is_scores)
        y = np.asarray(self.oos_scores)
        if np.std(x) == 0:
            return float("nan")
        return float(np.polyfit(x, y, 1)[0])

    def report(self) -> str:
        return "\n".join([
            f"  パラメータ数       : {self.n_params}",
            f"  ブロック数         : {self.n_blocks}",
            f"  組み合わせ数       : {self.n_combinations}",
            f"  PBO                : {self.pbo:>8.3f}  <- {self.verdict}",
            f"  IS->OOS の傾き     : {self.performance_degradation:>8.3f}"
            f"  （1.0で完全再現、0以下で無意味）",
            f"  OOS順位の中央値    : "
            f"{np.median(self.oos_ranks) if self.oos_ranks else float('nan'):>8.1f}"
            f" / {self.n_params}",
        ])


def _sharpe(x: np.ndarray) -> float:
    if x.size < 2:
        return 0.0
    s = x.std(ddof=1)
    return float(x.mean() / s) if s > 0 else 0.0


def compute_pbo(returns_matrix: pd.DataFrame, n_blocks: int = 10,
                score_fn=None) -> PboResult:
    """PBOを計算する。

    returns_matrix : index=時刻, columns=パラメータID, values=リターン
                     （各列が1つのパラメータ設定のリターン列）
    n_blocks       : 時間軸の分割数（偶数推奨）。大きいほど組合せが増える
    """
    if score_fn is None:
        score_fn = _sharpe

    M = returns_matrix.dropna(how="all").fillna(0.0)
    n_params = M.shape[1]
    if n_params < 2:
        raise ValueError("パラメータが2つ以上必要です")
    if n_blocks % 2 != 0:
        n_blocks += 1

    blocks = np.array_split(np.arange(M.shape[0]), n_blocks)
    half = n_blocks // 2

    logits: list[float] = []
    is_ranks: list[int] = []
    oos_ranks: list[int] = []
    is_scores: list[float] = []
    oos_scores: list[float] = []

    for is_blocks in combinations(range(n_blocks), half):
        oos_blocks = [b for b in range(n_blocks) if b not in is_blocks]

        is_idx = np.concatenate([blocks[b] for b in is_blocks])
        oos_idx = np.concatenate([blocks[b] for b in oos_blocks])

        is_perf = np.array([score_fn(M.iloc[is_idx, j].to_numpy())
                            for j in range(n_params)])
        oos_perf = np.array([score_fn(M.iloc[oos_idx, j].to_numpy())
                             for j in range(n_params)])

        best = int(np.argmax(is_perf))

        # OOSでの順位（1が最下位、n_paramsが最上位）
        order = np.argsort(np.argsort(oos_perf)) + 1
        rank = int(order[best])

        # 相対順位を [0,1] に正規化してロジット変換
        omega = rank / (n_params + 1)
        omega = min(max(omega, 1e-6), 1 - 1e-6)
        logits.append(float(np.log(omega / (1 - omega))))

        is_ranks.append(int(np.argsort(np.argsort(is_perf))[best]) + 1)
        oos_ranks.append(rank)
        is_scores.append(float(is_perf[best]))
        oos_scores.append(float(oos_perf[best]))

    logit_arr = np.asarray(logits)
    # ロジットが負 = OOS順位が中央値未満 = 過剰最適化
    pbo = float((logit_arr < 0).mean()) if logit_arr.size else float("nan")

    return PboResult(
        pbo=pbo, n_combinations=comb(n_blocks, half), n_params=n_params,
        n_blocks=n_blocks, logits=logit_arr, is_ranks=is_ranks,
        oos_ranks=oos_ranks, is_scores=is_scores, oos_scores=oos_scores,
    )


def build_returns_matrix(bars: pd.DataFrame, param_grid: list[dict],
                         evaluate) -> pd.DataFrame:
    """全パラメータのリターン列を並べた行列を作る（PBOの入力）。"""
    cols: dict[str, pd.Series] = {}
    for i, params in enumerate(param_grid):
        _, rets = evaluate(bars, params)
        if rets.empty:
            continue
        key = "_".join(f"{k}{v}" for k, v in sorted(params.items()))
        cols[key] = rets
    return pd.DataFrame(cols)


if __name__ == "__main__":
    print("=== 1) 意味のないパラメータ群（全部ランダム）===")
    rng = np.random.default_rng(11)
    n_t, n_p = 1000, 50
    random_matrix = pd.DataFrame(
        rng.normal(0, 0.01, size=(n_t, n_p)),
        index=pd.date_range("2024-01-01", periods=n_t, freq="1h", tz="UTC"),
        columns=[f"p{i}" for i in range(n_p)],
    )
    res = compute_pbo(random_matrix, n_blocks=10)
    print(res.report())
    print("  -> エッジが無いのでPBOは高くなる（0.5前後）\n")

    print("=== 2) 本物のエッジがあるパラメータ群 ===")
    # 一部のパラメータに本当にプラスの期待値を持たせる
    edge_matrix = random_matrix.copy()
    for j in range(5):
        edge_matrix.iloc[:, j] += 0.0008        # 明確なエッジ
    res2 = compute_pbo(edge_matrix, n_blocks=10)
    print(res2.report())
    print("  -> エッジがあればPBOは下がる（選択が機能している）\n")

    print("=== 3) 実データでの使用例 ===")
    try:
        from kline_store import KlineStore
        from walk_forward import sma_evaluate

        bars = KlineStore().load("BTCUSDT", "60")
        if bars.empty:
            raise FileNotFoundError
        bars = bars.loc["2024-01-01":]
        grid = [{"fast": f, "slow": s}
                for f in (5, 10, 15, 20, 30, 40, 50)
                for s in (60, 80, 100, 120, 150, 200) if f < s]

        M = build_returns_matrix(bars, grid, sma_evaluate)
        print(f"  {M.shape[1]}パラメータ × {M.shape[0]}本")
        print(compute_pbo(M, n_blocks=10).report())
    except (ImportError, FileNotFoundError):
        print("  （記事07のデータが必要です）")

"""plateau.py — パラメータ空間のプラトー（頑健域）を探す
本コードは本人利用に限り自由に使用可。再配布・商用利用・販売は禁止。

思想: 「最良の1点」は運の可能性が高い。
      「周辺も良い領域の中心」を選べば、将来も機能しやすい。
"""
from __future__ import annotations

from dataclasses import dataclass

import numpy as np
import pandas as pd


@dataclass
class PlateauResult:
    grid: pd.DataFrame            # index/columns=パラメータ, values=スコア
    smoothed: pd.DataFrame        # 近傍平均を取ったスコア
    best_raw: tuple
    best_raw_score: float
    best_plateau: tuple
    best_plateau_score: float
    best_plateau_raw_score: float

    def report(self) -> str:
        return "\n".join([
            f"  素の最良      : {self.best_raw} "
            f"スコア {self.best_raw_score:+.3f}",
            f"  プラトー最良  : {self.best_plateau} "
            f"近傍平均 {self.best_plateau_score:+.3f} "
            f"(素のスコア {self.best_plateau_raw_score:+.3f})",
            "",
            "  -> プラトー最良を採用するほうが、将来の劣化が小さい傾向",
        ])


def score_grid(bars: pd.DataFrame, evaluate,
               param_a: str, values_a: list,
               param_b: str, values_b: list,
               extra: dict | None = None) -> pd.DataFrame:
    """2つのパラメータでグリッドを作り、スコア表を返す。"""
    extra = extra or {}
    rows = {}
    for a in values_a:
        row = {}
        for b in values_b:
            params = {param_a: a, param_b: b, **extra}
            try:
                score, _ = evaluate(bars, params)
            except Exception:
                score = np.nan
            row[b] = score if np.isfinite(score) else np.nan
        rows[a] = row
    return pd.DataFrame(rows).T.sort_index()


def find_plateau(grid: pd.DataFrame, window: int = 1) -> PlateauResult:
    """近傍平均でスコアを平滑化し、頑健な領域の中心を探す。"""
    arr = grid.to_numpy(dtype=float)
    n_r, n_c = arr.shape
    smoothed = np.full_like(arr, np.nan)

    for i in range(n_r):
        for j in range(n_c):
            r0, r1 = max(0, i - window), min(n_r, i + window + 1)
            c0, c1 = max(0, j - window), min(n_c, j + window + 1)
            block = arr[r0:r1, c0:c1]
            valid = block[np.isfinite(block)]
            if valid.size >= 2:
                # 近傍の平均から標準偏差を引く（ばらつきを罰する）
                smoothed[i, j] = valid.mean() - 0.5 * valid.std()

    sm = pd.DataFrame(smoothed, index=grid.index, columns=grid.columns)

    raw_flat = grid.stack(dropna=True)
    sm_flat = sm.stack(dropna=True)
    if raw_flat.empty or sm_flat.empty:
        raise ValueError("有効なスコアがありません")

    best_raw = raw_flat.idxmax()
    best_plat = sm_flat.idxmax()

    return PlateauResult(
        grid=grid, smoothed=sm,
        best_raw=best_raw, best_raw_score=float(raw_flat.max()),
        best_plateau=best_plat, best_plateau_score=float(sm_flat.max()),
        best_plateau_raw_score=float(grid.loc[best_plat[0], best_plat[1]]),
    )


def print_heatmap(grid: pd.DataFrame, title: str = "",
                  chars: str = " .:-=+*#%@") -> None:
    """スコアを文字でヒートマップ表示する（外部依存なし）。"""
    arr = grid.to_numpy(dtype=float)
    valid = arr[np.isfinite(arr)]
    if valid.size == 0:
        print("表示できるデータがありません")
        return
    lo, hi = valid.min(), valid.max()
    span = (hi - lo) or 1.0

    print(f"\n{title}  (低 '{chars[0]}' 〜 高 '{chars[-1]}',"
          f" 範囲 {lo:+.2f}〜{hi:+.2f})")
    header = "        " + "".join(f"{c:>6}" for c in grid.columns)
    print(header)
    for idx, row in grid.iterrows():
        cells = []
        for v in row:
            if not np.isfinite(v):
                cells.append("     -")
            else:
                level = int((v - lo) / span * (len(chars) - 1))
                cells.append(f"  {chars[level]}{v:>4.1f}"[-6:])
        print(f"{idx:>7} " + "".join(cells))


if __name__ == "__main__":
    from kline_store import KlineStore
    from walk_forward import sma_evaluate

    bars = KlineStore().load("BTCUSDT", "60")
    if bars.empty:
        raise SystemExit("先に記事07でデータを取得してください")
    bars = bars.loc["2024-01-01":]

    fasts = [5, 10, 15, 20, 25, 30, 40, 50]
    slows = [60, 80, 100, 120, 150, 180, 200]

    grid = score_grid(bars, sma_evaluate, "fast", fasts, "slow", slows)
    print_heatmap(grid, "素のSharpe（fast × slow）")

    res = find_plateau(grid, window=1)
    print_heatmap(res.smoothed, "平滑化スコア（近傍平均 − 0.5×標準偏差）")
    print()
    print(res.report())

    print("\n--- 判断の指針 ---")
    print("  ヒートマップが『まだら』なら過剰最適化の兆候。")
    print("  広い領域が連続して高いなら、そのパラメータには意味がある。")
    print("  最良の1点だけが飛び抜けている場合は、その点を採用しない。")

"""adoption_gate.py — 戦略を採用してよいかを機械的に判定する
本コードは本人利用に限り自由に使用可。再配布・商用利用・販売は禁止。

すべての検証を通し、基準を満たさない戦略は採用しない。
「なんとなく良さそう」で実資金を投入しないための仕組み。
"""
from __future__ import annotations

from dataclasses import dataclass, field

import numpy as np
import pandas as pd


@dataclass
class Criteria:
    """採用の最低基準。厳しすぎると何も通らないが、それが正常。"""

    min_trades: int = 100
    min_oos_sharpe: float = 0.5
    min_wfe: float = 0.5                 # Walk-Forward効率比
    min_dsr: float = 0.95                # Deflated Sharpe Ratio
    max_pbo: float = 0.4                 # 過剰最適化確率
    min_param_stability: float = 0.4
    max_dd_ratio: float = 0.5            # 最大DD / 年間リターン
    require_pessimistic_positive: bool = True
    min_oos_positive_ratio: float = 0.5  # OOS窓のうちプラスの割合


@dataclass
class Evidence:
    """各検証の結果をまとめて持つ。"""

    n_trades: int = 0
    oos_sharpe: float = float("nan")
    wfe: float = float("nan")
    dsr: float = float("nan")
    pbo: float = float("nan")
    param_stability: float = float("nan")
    max_dd: float = float("nan")
    annual_return: float = float("nan")
    pessimistic_net: float = float("nan")
    oos_positive_ratio: float = float("nan")
    lookahead_test_passed: bool = False
    n_trials_declared: int = 0


@dataclass
class Decision:
    approved: bool
    blockers: list[str] = field(default_factory=list)
    warnings: list[str] = field(default_factory=list)
    passed: list[str] = field(default_factory=list)

    def report(self) -> str:
        lines = [
            "=" * 60,
            f"  判定: {'採用可' if self.approved else '採用不可'}",
            "=" * 60,
        ]
        if self.passed:
            lines.append("通過した基準:")
            lines += [f"  [OK] {p}" for p in self.passed]
        if self.warnings:
            lines.append("警告:")
            lines += [f"  [!]  {w}" for w in self.warnings]
        if self.blockers:
            lines.append("不合格:")
            lines += [f"  [NG] {b}" for b in self.blockers]
        if not self.approved:
            lines += [
                "",
                "  実資金を投入してはいけません。",
                "  パラメータを再調整して再検証するのは"
                "『試行回数を増やす』ことなので、",
                "  DSRの前提が悪化する点に注意してください。",
            ]
        return "\n".join(lines)


def evaluate_adoption(ev: Evidence, cr: Criteria | None = None) -> Decision:
    cr = cr or Criteria()
    blockers: list[str] = []
    warnings: list[str] = []
    passed: list[str] = []

    def check(ok: bool | None, label: str, blocker: bool = True) -> None:
        if ok is None:
            warnings.append(f"{label}: 未検証")
        elif ok:
            passed.append(label)
        elif blocker:
            blockers.append(label)
        else:
            warnings.append(label)

    # 必須: ルックアヘッドテスト（記事31）
    check(ev.lookahead_test_passed or None if not ev.lookahead_test_passed
          else True,
          "ルックアヘッドテストに合格")

    check(ev.n_trades >= cr.min_trades if ev.n_trades else None,
          f"約定回数 {ev.n_trades} >= {cr.min_trades}")

    check(ev.oos_sharpe >= cr.min_oos_sharpe if np.isfinite(ev.oos_sharpe)
          else None,
          f"OOS Sharpe {ev.oos_sharpe:.2f} >= {cr.min_oos_sharpe}")

    check(ev.wfe >= cr.min_wfe if np.isfinite(ev.wfe) else None,
          f"Walk-Forward効率比 {ev.wfe:.2f} >= {cr.min_wfe}")

    check(ev.dsr >= cr.min_dsr if np.isfinite(ev.dsr) else None,
          f"Deflated Sharpe {ev.dsr:.3f} >= {cr.min_dsr}")

    check(ev.pbo <= cr.max_pbo if np.isfinite(ev.pbo) else None,
          f"PBO {ev.pbo:.3f} <= {cr.max_pbo}")

    check(ev.param_stability >= cr.min_param_stability
          if np.isfinite(ev.param_stability) else None,
          f"パラメータ安定性 {ev.param_stability:.1%} >= "
          f"{cr.min_param_stability:.0%}", blocker=False)

    check(ev.oos_positive_ratio >= cr.min_oos_positive_ratio
          if np.isfinite(ev.oos_positive_ratio) else None,
          f"OOSプラス率 {ev.oos_positive_ratio:.1%} >= "
          f"{cr.min_oos_positive_ratio:.0%}")

    if np.isfinite(ev.max_dd) and np.isfinite(ev.annual_return) \
            and ev.annual_return > 0:
        ratio = abs(ev.max_dd) / ev.annual_return
        check(ratio <= cr.max_dd_ratio,
              f"最大DD/年間リターン {ratio:.2f} <= {cr.max_dd_ratio}",
              blocker=False)

    if cr.require_pessimistic_positive:
        check(ev.pessimistic_net > 0 if np.isfinite(ev.pessimistic_net)
              else None,
              "悲観コストモデルでも収益がプラス")

    if ev.n_trials_declared < 10:
        warnings.append(
            f"申告した試行回数が {ev.n_trials_declared} と少なすぎます。"
            f"捨てたアイデアも数えていますか？"
        )

    return Decision(approved=not blockers, blockers=blockers,
                    warnings=warnings, passed=passed)


if __name__ == "__main__":
    print("### 例1: 過剰最適化された戦略 ###")
    bad = Evidence(
        n_trades=45, oos_sharpe=0.3, wfe=0.18, dsr=0.42, pbo=0.61,
        param_stability=0.15, max_dd=-0.35, annual_return=0.20,
        pessimistic_net=-120.0, oos_positive_ratio=0.38,
        lookahead_test_passed=True, n_trials_declared=2400,
    )
    print(evaluate_adoption(bad).report())

    print("\n\n### 例2: 検証を通過した戦略 ###")
    good = Evidence(
        n_trades=320, oos_sharpe=0.92, wfe=0.71, dsr=0.968, pbo=0.22,
        param_stability=0.55, max_dd=-0.14, annual_return=0.31,
        pessimistic_net=1850.0, oos_positive_ratio=0.68,
        lookahead_test_passed=True, n_trials_declared=180,
    )
    print(evaluate_adoption(good).report())