import numpy as np  # 用数组按同一叶均值口径复算所有候选分裂
import pandas as pd  # 将候选阈值与复算结果整理为检查表
tree_feature_values = np.array([2.2, 3.5, 6.0, 9.8, 15.0])  # 固定题面五个特征值
tree_target_values = np.array([10.0, 12.0, 25.0, 80.0, 90.0])  # 固定题面五个目标值
candidate_thresholds = np.array([2.85, 4.75, 7.90, 12.40])  # 由相邻特征中点生成候选阈值
def calculate_split_rss(split_threshold):  # 对一个候选按两叶各自均值计算总平方误差
    left_targets = tree_target_values[tree_feature_values < split_threshold]  # 取阈值左侧目标
    right_targets = tree_target_values[tree_feature_values >= split_threshold]  # 取阈值右侧目标
    return float(np.square(left_targets - left_targets.mean()).sum() + np.square(right_targets - right_targets.mean()).sum())  # 合并两叶 RSS
candidate_rss = np.array([calculate_split_rss(threshold) for threshold in candidate_thresholds])  # 一次复算四个候选
expected_rss = np.array([4556.75, 2452.00, 182.6666666667, 3236.75])  # 确定可手算的正确答案
assert np.allclose(candidate_rss, expected_rss, atol=1e-8) and candidate_thresholds[candidate_rss.argmin()] == 7.90  # 确定四值与最优阈值
display(pd.DataFrame({'分割点 s': candidate_thresholds, '总 RSS': candidate_rss}).round(2))  # 输出与后续手算和图形同源的比较表

# 为“基尼 vs 熵：有何区别？”导入 `numpy` 并绑定 `np`，用于执行当前任务的数组、数值或随机机制计算。
import numpy as np
# 为“基尼 vs 熵：有何区别？”导入 `matplotlib.pyplot` 并绑定 `plt`，用于构建当前任务的坐标轴并呈现比较结果。
import matplotlib.pyplot as plt

# 建立 `p` 的有序取值网格，用于展示“二分类问题中，基尼不纯度与信息熵的比较”随参数变化的比较结果。
p = np.linspace(0.001, 0.999, 200)
# 按二分类定义计算各概率下的基尼不纯度。
gini = 2 * p * (1-p)
# 对 `p * np.log2(p) + (1 - p) * np.log2(1 - p)` 应用一元运算，得到 `entropy`。
entropy = - (p * np.log2(p) + (1-p) * np.log2(1-p))


# 创建 `fig, ax` 画布，承载“二分类问题中，基尼不纯度与信息熵的比较”的并排视觉比较。
fig, ax = plt.subplots(figsize=(7, 3.5))
# 以 `p` 为横轴、`gini` 为纵轴绘制曲线，展示“二分类问题中，基尼不纯度与信息熵的比较”。
ax.plot(p, gini, label='Gini Impurity (基尼不纯度)', color='royalblue', linewidth=2.5)
# 以 `p` 为横轴、`entropy` 为纵轴绘制曲线，展示“二分类问题中，基尼不纯度与信息熵的比较”。
ax.plot(p, entropy, label='Entropy (信息熵)', color='crimson', linewidth=2.5)
# 将图题设为“基尼不纯度 vs. 信息熵 (二分类)”，直接说明当前图形的比较目的。
ax.set_title('基尼不纯度 vs. 信息熵 (二分类)', fontsize=16)
# 将横轴标为“类别 1 的比例 (p)”，明确横向编码的变量。
ax.set_xlabel('类别 1 的比例 (p)', fontsize=12)
# 将纵轴标为“不纯度”，明确纵向编码的变量。
ax.set_ylabel('不纯度', fontsize=12)
# 在 `0.5` 处添加垂直参考线，标出“二分类问题中，基尼不纯度与信息熵的比较”的基准或阈值。
ax.axvline(0.5, color='grey', linestyle='--', lw=1)
# 以箭头标注“不纯度最高”，指出“二分类问题中，基尼不纯度与信息熵的比较”中的候选拐点或阈值。
ax.annotate('不纯度最高', xy=(0.5, 1.0), xytext=(0.55, 0.9),
            # 指定 `arrowprops` 为注释箭头样式，细化“二分类问题中，基尼不纯度与信息熵的比较”的输出。
            arrowprops=dict(facecolor='black', shrink=0.05, width=1, headwidth=8),
            # 将“最大差异”注释设为 12pt，保证箭头说明在图中可读。
            fontsize=12)
# 显示“基尼 vs 熵：有何区别？”图例，使颜色或线型与比较对象一一对应。
ax.legend()
# 显示基尼不纯度与熵随正类比例变化的曲线，比较两者在纯节点处为零、在均衡节点处最大的共同性质。
plt.show()

# 为“步骤一：导入必要的库”导入 `pandas` 并绑定 `pd`，用于整理当前任务的表格、字段与时间索引。
import pandas as pd
# 为“步骤一：导入必要的库”导入 `numpy` 并绑定 `np`，用于执行当前任务的数组、数值或随机机制计算。
import numpy as np

# 为“步骤一：导入必要的库”，机器学习库
from sklearn.tree import DecisionTreeRegressor
# 为“步骤一：导入必要的库”，从 `sklearn.metrics` 导入`mean_squared_error` 用于计算真实值与预测值的均方误差。
from sklearn.metrics import mean_squared_error

# 为“步骤一：导入必要的库”导入 `matplotlib.pyplot` 并绑定 `plt`，用于构建当前任务的坐标轴并呈现比较结果。
import matplotlib.pyplot as plt
# 为“步骤一：导入必要的库”，导入决策树绘图函数，展示分裂特征与阈值。
from sklearn.tree import plot_tree

from pathlib import Path  # 管理 Linux 与 Windows 的公开数据路径

from urllib.request import urlretrieve  # 复用本章隐藏设置单元安装的浏览器标识下载器
# 按 Linux 共享数据、Windows 共享数据、项目缓存的顺序选择财务报表文件。
financial_path = next((candidate_path for candidate_path in [Path('/home/ubuntu/r2_data_mount/data/stock/financial_statement.h5'), Path('C:/qiufei/data/stock/financial_statement.h5'), Path('data/course/financial_statement.h5')] if candidate_path.exists()), Path('data/course/financial_statement.h5'))
if not financial_path.exists():
    financial_path.parent.mkdir(parents=True, exist_ok=True)
    urlretrieve('https://assets.qiufei.site/data/stock/financial_statement.h5', financial_path)
yangtze_company_ids = ['600276.XSHG', '600009.XSHG', '600585.XSHG', '603259.XSHG', '002415.XSHE', '600104.XSHG', '600019.XSHG', '600406.XSHG']  # 固定沪苏浙皖公司池
financial_columns = ['order_book_id', 'quarter', 'info_date', 'total_assets', 'total_liabilities', 'net_profit_parent_company', 'basic_earnings_per_share']  # 限定读取字段
company_frames = []  # 收集公司级年度财报
for company_id in yangtze_company_ids:  # 避免载入全市场全部列
    company_frame = pd.read_hdf(financial_path, key='financial_data', where=f'order_book_id == "{company_id}"', columns=financial_columns)  # 选择性读取单家公司
    company_frames.append(company_frame)  # 保存当前公司数据
annual_raw = pd.concat(company_frames, ignore_index=True)  # 合并但保留全部披露版本
annual_raw = annual_raw.query("quarter.str.endswith('q4')", engine='python').copy()  # 只保留年度报告版本
annual_raw['info_date'] = pd.to_datetime(annual_raw['info_date'])  # 统一信息可得日类型
annual_raw['feature_year'] = annual_raw['quarter'].str[:4].astype(int)  # 提取报告财年
annual_raw = annual_raw.query('2012 <= feature_year <= 2024').copy()  # 确定财年范围
display(annual_raw[['order_book_id', 'quarter', 'info_date', 'basic_earnings_per_share']].rename(columns={
    'order_book_id': '证券代码', 'quarter': '报告期', 'info_date': '披露日',
    'basic_earnings_per_share': '基本每股收益'
}).head(8))

def last_asof(raw, fiscal_year, cutoff):  # 每个截止日都从原始版本重建数据版本
    all_versions = raw.loc[raw['feature_year'].eq(fiscal_year)].copy()  # 保留该财年全部版本
    eligible = all_versions.loc[all_versions['info_date'].le(cutoff)].copy()  # 先按截止日过滤
    if eligible.empty:  # 空窗口不得真空通过
        raise RuntimeError(f'{fiscal_year}@{cutoff.date()} 空窗')  # 发现点时空窗即中止，避免静默使用未来版本
    snapshot = eligible.sort_values(['order_book_id', 'feature_year', 'info_date']).drop_duplicates(['order_book_id', 'feature_year'], keep='last')  # 再选点时版本
    assert snapshot['info_date'].le(cutoff).all()  # 确定可得日不变量
    flow = {'raw_versions': len(all_versions), 'eligible_versions': len(eligible), 'excluded_after_cutoff': len(all_versions) - len(eligible), 'unique_companies': snapshot['order_book_id'].nunique(), 'max_info_date': snapshot['info_date'].max()}  # 记录样本数量变化
    return snapshot, flow  # 返回数据版本与检查结果

def build_eps_year(target_year):  # 构造固定目标年的特征与事后评分标签
    decision = pd.Timestamp(f'{target_year + 2}-06-30')  # 按公开披露滞后在目标年后第二个六月末预测
    score = pd.Timestamp(f'{target_year + 3}-06-30')  # 再过一年评分
    features, feature_flow = last_asof(annual_raw, target_year - 1, decision)  # 取上一财年点时特征
    label_all = annual_raw.loc[annual_raw['feature_year'].eq(target_year)].copy()  # 标签版本先不去重
    assert len(label_all) > 0 and label_all['info_date'].min() > decision  # 任何目标版本均不得在决策日可见
    labels, label_flow = last_asof(annual_raw, target_year, score)  # 到评分日选择标签版本
    features = features.rename(columns={'info_date': 'feature_info_date'})  # 明示特征可得日
    labels = labels[['order_book_id', 'info_date', 'basic_earnings_per_share']].rename(columns={'info_date': 'target_info_date', 'basic_earnings_per_share': 'next_eps'})  # 明示标签可得日
    panel = features.merge(labels, on='order_book_id', how='inner').assign(target_year=target_year, decision_date=decision, score_cutoff=score)  # 合并公司年
    panel['log_assets'] = np.log(panel['total_assets'])  # 构造规模
    panel['debt_ratio'] = panel['total_liabilities'] / panel['total_assets']  # 构造负债率
    panel['roa'] = panel['net_profit_parent_company'] / panel['total_assets']  # 构造盈利率
    assert len(panel) > 0 and panel['feature_info_date'].le(decision).all() and panel['target_info_date'].gt(decision).all() and panel['target_info_date'].le(score).all()  # 阻止泄漏
    return panel, {'target_year': target_year, 'feature': feature_flow, 'label': label_flow, 'merged': len(panel)}  # 返回样本与样本数量变化

point_time_parts, version_flows = [], []  # 只收集开发期逐目标年样本与检查结果
for target_year in range(2016, 2023):  # 构造 2016—2022 开发目标年，不提前构造 2023 标签
    year_panel, year_flow = build_eps_year(target_year)  # 执行固定决策日与评分日说明
    point_time_parts.append(year_panel)  # 保存当年面板
    version_flows.append(year_flow)  # 保存版本选择证据
development_panel = pd.concat(point_time_parts, ignore_index=True).replace([np.inf, -np.inf], np.nan).dropna(subset=['log_assets', 'debt_ratio', 'roa', 'next_eps'])  # 合并开发期有效样本
assert not development_panel.duplicated(['order_book_id', 'target_year']).any()  # 每公司目标年唯一
assert development_panel['target_year'].max() == 2022  # 保证最终测试目标年尚未构造或查看
assert set(development_panel['target_year']) == set(range(2016, 2023))  # 核对七个开发目标年完整
feature_names = ['log_assets', 'debt_ratio', 'roa']  # 确定特征字典
X_development, y_development = development_panel[feature_names], development_panel['next_eps']  # 构造开发期特征和目标
flow_table = pd.json_normalize(version_flows)[['target_year', 'feature.raw_versions', 'feature.excluded_after_cutoff', 'label.raw_versions', 'label.excluded_after_cutoff', 'merged']].set_axis(['年', '特征版本', '特征排除', '标签版本', '标签排除', '合并'], axis=1)  # 用短列名展示完整样本数量变化
display(flow_table)  # 输出全部目标年的版本与排除数
print({'开发期X': X_development.shape, 'feature_max': development_panel.feature_info_date.max().date(), 'label_min': development_panel.target_info_date.min().date(), '最终测试状态': '2023 标签尚未构造'})  # 输出紧凑不变量摘要

demo_validation_year = 2022  # 确定开发期伪测试目标年
demo_training_mask = development_panel['target_year'] < demo_validation_year  # 较早开发目标年用于训练
demo_validation_mask = development_panel['target_year'] == demo_validation_year  # 2022 只作开发期评分
X_demo_train = X_development.loc[demo_training_mask]  # 提取开发训练特征
X_demo_valid = X_development.loc[demo_validation_mask]  # 提取开发期伪测试特征
y_demo_train = y_development.loc[demo_training_mask]  # 对齐开发训练目标
y_demo_valid = y_development.loc[demo_validation_mask]  # 对齐开发期伪测试目标
print({'开发训练': X_demo_train.shape, '2022伪测试': X_demo_valid.shape})  # 报告开发窗样本形状
print(f'开发训练披露期: {development_panel.loc[demo_training_mask, "target_info_date"].min().date()}—{development_panel.loc[demo_training_mask, "target_info_date"].max().date()}')  # 报告训练窗口
print(f'2022伪测试披露期: {development_panel.loc[demo_validation_mask, "target_info_date"].min().date()}—{development_panel.loc[demo_validation_mask, "target_info_date"].max().date()}')  # 报告开发评分窗口

# 建立最大深度为二的回归树并固定随机状态，限制教学模型复杂度。
cart_regressor = DecisionTreeRegressor(max_depth=2, random_state=42)

# 仅用训练期特征与 EPS 目标估计回归树的分裂阈值和叶节点均值。
# .fit() 方法是所有scikit-learn模型的训练入口
cart_regressor.fit(X_demo_train, y_demo_train)

# 报告示范树已只用较早开发窗拟合，2023 仍未构造。
print('2022 开发期伪测试示范树训练完成；2023 仍最终测试。')

# 对 2022 开发期伪测试特征生成预测。
demo_predictions = cart_regressor.predict(X_demo_valid)

# 将开发期真实值与预测值整理为可核对表格。
demo_results = pd.DataFrame({'2022开发真实值': y_demo_valid, '示范树预测值': demo_predictions})
# 输出前五条开发期结果，明确它不是 2023 最终测试证据。
display(demo_results.head())

# 计算 2022 开发期伪测试 MSE，只用于流程示范。
demo_tree_mse = mean_squared_error(y_demo_valid, demo_predictions)
# 只用较早开发训练期均值建立同窗朴素基准。
demo_mean_predictions = np.repeat(y_demo_train.mean(), len(y_demo_valid))
# 计算同一 2022 开发窗的训练均值基准误差。
demo_baseline_mse = mean_squared_error(y_demo_valid, demo_mean_predictions)
# 报告开发证据并再次声明 2023 尚未开启。
print({'2022示范树MSE': round(demo_tree_mse, 4), '2022训练均值MSE': round(demo_baseline_mse, 4), '最终测试状态': '2023 未开启'})

# 为“步骤七：可视化与解读决策树”建立投影画布，给坐标、图例与标注预留空间。
plt.figure(figsize=(22, 10))

# 绘出回归树的分裂条件、样本量、节点误差与预测值。
plot_tree(cart_regressor,
          # 固定 `feature_names` 的特征名称，使图表标签与模型输入列一一对应。
          feature_names=feature_names,
          filled=True, # 用颜色梯度编码各回归节点的预测值
          rounded=True, # 使用圆角矩形
          impurity=False, # 省略课堂解读不需要的节点误差字段
          proportion=True, # 用样本比例替代较难比较的原始计数
          # 指定 `fontsize` 为投影字号，细化“训练好的每股收益预测决策树 (最大深度=2)”的输出。
          fontsize=16,
          precision=3) # 显示3位小数

# 将图题设为“决策树模型：预测每股收益 (EPS)”，直接说明当前图形的比较目的。
plt.title('2022 开发期伪测试示范树：预测每股收益 (EPS)', fontsize=20)
# 显示已拟合决策树的分裂变量、阈值、样本数与叶节点预测，核对 EPS 预测路径。
plt.show()

expected_development_years = [2019, 2020, 2021, 2022]  # 确定四个开发验证窗
candidate_counts = candidate_table.groupby('validation_year').size()  # 统计每窗候选覆盖
assert candidate_counts.reindex(expected_development_years).eq(9).all()  # 每窗必须九规格齐全
assert len(candidate_table) == 36  # 四窗乘九规格必须形成三十六条比较结果
candidate_wide = candidate_table.pivot(index=['max_depth', 'min_samples_leaf'], columns='validation_year', values='mse').reindex(columns=expected_development_years).reset_index()  # 压缩为九规格乘四窗
assert candidate_wide.shape == (9, 6)  # 两个参数列加四个开发窗
development_baseline_by_year = candidate_table.groupby('validation_year')['mean_baseline_mse'].first().reindex(expected_development_years)  # 提取四窗训练均值基线
baseline_consistency = candidate_table.groupby('validation_year')['mean_baseline_mse'].nunique().reindex(expected_development_years)  # 核对同窗基线未随树规格变化
assert baseline_consistency.eq(1).all()  # 每个开发窗只能有一个训练均值基线
print({'每窗候选数': candidate_counts.to_dict(), '候选总行数': len(candidate_table), '训练均值基线MSE': development_baseline_by_year.round(4).to_dict()})  # 输出覆盖与同窗基线检查
display(candidate_wide.round(4))  # 完整展示三十六个候选 MSE

candidate_summary = candidate_table.groupby(['max_depth', 'min_samples_leaf'], as_index=False).agg(mean_validation_mse=('mse', 'mean'), worst_validation_mse=('mse', 'max'))  # 汇总四窗平均与最差误差
assert len(candidate_summary) == 9  # 九种规格必须全部进入汇总
best_mse = candidate_summary['mean_validation_mse'].min()  # 确定最小跨窗均值 MSE
shortlist = candidate_summary.loc[candidate_summary['mean_validation_mse'].le(best_mse + 1e-12)]  # 只保留容差内并列者
selected_spec = shortlist.sort_values(['max_depth', 'min_samples_leaf'], ascending=[True, False]).iloc[0]  # 并列选更浅、再选更大叶节点
final_tree = DecisionTreeRegressor(max_depth=int(selected_spec.max_depth), min_samples_leaf=int(selected_spec.min_samples_leaf), random_state=42).fit(development_pool[feature_names], development_pool['next_eps'])  # 用完整开发期重新拟合
display(candidate_summary.round(4))  # 展示全部九种规格的确定依据
print({'选择结果': {'最大深度': int(selected_spec.max_depth), '叶节点最少样本数': int(selected_spec.min_samples_leaf)}, '并列时的选择': '先选浅树，再选叶节点样本数较大的树'})

sealed_rows, sealed_flow = build_eps_year(2023)  # 确定规格与全开发期重训后才构造最终测试标签
sealed_rows = sealed_rows.replace([np.inf, -np.inf], np.nan).dropna(subset=[*feature_names, 'next_eps'])  # 按确定字段清理最终测试样本
assert sealed_rows['target_year'].eq(2023).all() and len(sealed_rows) > 0  # 核对唯一最终测试目标年与非空样本
assert not sealed_rows.duplicated(['order_book_id', 'target_year']).any()  # 核对公司目标年唯一
sealed_predictions = final_tree.predict(sealed_rows[feature_names])  # 生成最终测试预测
sealed_baseline = np.repeat(development_pool['next_eps'].mean(), len(sealed_rows))  # 建立全开发期均值基线
sealed_results = sealed_rows[['order_book_id', 'next_eps']].assign(tree_prediction=sealed_predictions, mean_baseline=sealed_baseline)  # 对齐最终测试真实值与两组预测
sealed_metrics = {'树MSE': mean_squared_error(sealed_rows['next_eps'], sealed_predictions), '训练均值MSE': mean_squared_error(sealed_rows['next_eps'], sealed_baseline)}  # 只在本块计算一次最终测试指标
tree_has_incremental_value = sealed_metrics['树MSE'] < sealed_metrics['训练均值MSE']  # 以同一最终测试窗基线比较决定是否存在增量预测证据
test_conclusion = '值得用新的时间段继续验证' if tree_has_incremental_value else '暂不采用：没有优于训练均值基线'
print({'目标年': 2023, '最大深度': int(selected_spec.max_depth), '叶节点最少样本数': int(selected_spec.min_samples_leaf), **sealed_metrics, '优于基线': tree_has_incremental_value, '结论': test_conclusion})
display(sealed_results.rename(columns={
    'order_book_id': '证券代码', 'next_eps': '实际每股收益',
    'tree_prediction': '决策树预测', 'mean_baseline': '训练均值预测'
}).round(4))

first_values = sealed_rows.iloc[[0]][feature_names]  # 固定首个最终测试样本
tree_state = final_tree.tree_  # 读取所选决策树结构
node_id = 0  # 从根节点开始追踪规则
rule_card = []  # 收集每个分裂与最终叶节点
while tree_state.children_left[node_id] != tree_state.children_right[node_id]:  # 遍历到叶节点为止
    feature_name = feature_names[tree_state.feature[node_id]]  # 取得当前节点特征名
    threshold = tree_state.threshold[node_id]  # 取得当前节点阈值
    go_left = float(first_values.iloc[0][feature_name]) <= threshold  # 按首个样本决定分支方向
    feature_label = {'log_assets': '总资产对数', 'debt_ratio': '负债率'}[feature_name]
    rule_card.append({'变量': feature_label, '分割值': threshold, '方向': '<=' if go_left else '>', '样本值': float(first_values.iloc[0][feature_name])})
    node_id = tree_state.children_left[node_id] if go_left else tree_state.children_right[node_id]  # 前进到对应子节点
rule_card.append({'变量': '叶节点', '分割值': np.nan, '方向': '预测结果', '样本值': float(tree_state.value[node_id].ravel()[0])})
np.testing.assert_allclose(rule_card[-1]['样本值'], final_tree.predict(first_values)[0])  # 叶节点数值必须等于模型预测
display(pd.DataFrame(rule_card).round(4))  # 输出可检查的逐节点规则卡