到目前为止，本书主要关注线性模型。线性模型相对简单，易于描述和实现，在解释和推断方面具有优势。然而，标准线性回归在预测能力方面存在显著局限性。这是因为线性假设几乎总是一种近似，有时是一种较差的近似。在 章节 7
在本章中，我们将放宽线性假设，同时尽可能保持可解释性。我们通过检查线性模型的非常简单的扩展来实现这一点，例如多项式回归和阶梯函数，以及更复杂的方法，如样条、局部回归和广义加性模型。
多项式回归 (Polynomial Regression)通过添加额外的预测变量来扩展线性模型，这些预测变量是通过将每个原始预测变量提升到幂次获得的。例如，三次回归使用三个变量 \(X, X^2, X^3\) 作为预测变量。这种方法提供了对数据进行非线性拟合的简单方式。
阶梯函数 (Step Functions)将变量的范围切割成 \(K\) 个不同的区域以产生定性变量。其效果是拟合分段常数函数。
回归样条 (Regression Splines)比多项式和阶梯函数更灵活，实际上是两者的扩展。它们涉及将 \(X\) 的范围分成 \(K\) 个不同的区域。在每个区域内，对数据拟合多项式函数。然而，这些多项式被约束为在区域边界或节点 (knots)处平滑连接。只要区间被分成足够的区域，这可以产生极其灵活的拟合。
平滑样条 (Smoothing Splines)与回归样条类似，但产生于稍微不同的情况。平滑样条源于最小化残差平方和准则，并施加平滑性惩罚。
局部回归 (Local Regression)与样条类似，但在一个重要方面有所不同。区域允许重叠，并且实际上以非常平滑的方式重叠。
广义加性模型 (Generalized Additive Models, GAMs)允许我们将上述方法扩展以处理多个预测变量。
在 小节 8.1 小节 8.6 \(Y\) 与单个预测变量 \(X\) 之间关系的方法。在 小节 8.7 \(Y\) 建模为多个预测变量 \(X_1, \ldots, X_p\) 的函数。
多项式回归 (Polynomial Regression)
从历史上看，将线性回归扩展到预测变量与响应变量之间关系为非线性设置的标准方法是将标准线性模型
\[ y_i = \beta_0 + \beta_1 x_i + \epsilon_i \]
替换为多项式函数
\[ y_i = \beta_0 + \beta_1 x_i + \beta_2 x_i^2 + \beta_3 x_i^3 + \cdots + \beta_d x_i^d + \epsilon_i \tag{8.1}\]
其中 \(\epsilon_i\) 是误差项。这种方法称为多项式回归 。事实上，我们在 章节 4 \(d\) ，多项式回归允许我们产生高度非线性的曲线。请注意，式 8.1 中的系数可以很容易地使用最小二乘线性回归来估计，因为这只是一个标准的线性模型，预测变量为 \(x_i, x_i^2, x_i^3, \ldots, x_i^d\) 。一般来说，使用大于3或4的 \(d\) 是不常见的，因为对于大的 \(d\) 值，多项式曲线会变得过于灵活，可能会呈现一些非常奇怪的形状。这在 \(X\) 变量的边界附近尤其如此。
Note: 多项式回归的数学原理
多项式回归本质上仍然是线性回归，因为它关于参数是线性的。我们将原始变量 \(X\) 转换为多个特征 \(X, X^2, X^3, \ldots, X^d\) ，然后对这些特征拟合线性模型。这意味着我们可以使用所有标准线性回归的工具和推断方法。
多项式回归的关键选择是次数 \(d\) ：
\(d=1\) \(d=2\) \(d=3\) \(d \geq 4\)
高次多项式在边界处的行为往往不稳定，这是龙格现象 (Runge’s phenomenon)的表现之一。
中国案例：海康威视股价的长期趋势分析
非线性模型在金融时间序列分析中有着广泛的应用。虽然股票收益率通常被假设为难以预测（有效市场假说），但在长期视角下，股价往往表现出明显的非线性趋势 和周期性波动 。
让我们使用海康威视 (002415.XSHE) 的历史股价数据，探索如何使用超越线性的方法来提取其长期价格趋势。我们将尝试使用多项式回归、样条和GAM来拟合股价随时间的变化。
为了在一个极具挑战性的真实金融场景中探索非线性回归模型的边界，下面的 Python 实战代码将目标锁定在了中国 A 股智能安防龙头——海康威视（002415.XSHE）身上。与之前用不同特征去预测当期截面收益率不同，这次我们的 \(X\) 轴变成了一条漫长且孤独的绝对时间线（使用 np.arange 生成的一维等距时间索引），而 \(Y\) 轴则是海康威视多年来令人瞩目的未经复权处理的绝对收盘价。通过这种设置，我们将回归降维成了一个纯粹的“时间序列趋势拟合”问题。代码首先完成了严谨的基础数据清洗和对齐，随后的逻辑将试图用越来越复杂的高维多项式曲线，去强行适应海康威视那条伴随着宏观经济周期跌宕起伏的惊人增长曲线。
接下来，我们加载并准备海康威视的历史行情数据。
# ── 1. 加载海康威视（002415.XSHE）历史日线行情 ── import os # 导入 os 用于跨平台路径处理 # 根据操作系统自动选择本地数据路径 # 根据操作系统设置数据根目录路径 = 'C:/qiufei/data' if os.name == 'nt' else '/home/ubuntu/r2_data_mount/qiufei/data' = os.path.join(DATA_DIR, 'stock/stock_price_pre_adjusted.h5' ) # 拼接前复权股价文件路径 = pd.read_hdf(path) # 读取前复权日线数据（含全部 A 股） # 若 h5 以 MultiIndex 存储，则先展平为普通列 if 'order_book_id' in market_stock_data.index.names: # 获取索引 = market_stock_data.reset_index() # 将索引列还原为普通列 # 统一日期列名为 trade_date，兼容不同数据源的命名差异 # 获取列名列表 if 'date' in market_stock_data.columns and 'trade_date' not in market_stock_data.columns:= market_stock_data.rename(columns= {'date' : 'trade_date' }) # 重命名日期列 # 筛选海康威视，并按交易日升序排列 = market_stock_data[market_stock_data['order_book_id' ] == '002415.XSHE' ].copy() # 过滤海康威视数据 = haikang_data.sort_values('trade_date' ).reset_index(drop= True ) # 按日期排序并重置索引 # 构造回归所需的两列： # Time — 等距整数索引 (0, 1, 2, ...)，即"第几个交易日" # Price — 当日收盘价（前复权），作为被预测的响应变量 Y 'Time' ] = np.arange(len (haikang_data)) # 生成从 0 开始的等距时间索引 'Price' ] = haikang_data['close' ] # 将收盘价列赋值给 Price 列 # 打印数据概况，帮助确认加载是否正常 print (f"数据范围: { haikang_data['trade_date' ]. min ()} 至 { haikang_data['trade_date' ]. max ()} " ) # 打印起止日期 print (f"样本数量: { len (haikang_data)} " ) # 打印总交易日数
数据范围: 2010-05-28 00:00:00 至 2025-12-31 00:00:00
样本数量: 3789
运行结果显示，海康威视前复权股价数据覆盖了2010年5月28日至2025年12月31日的完整交易区间，共计3789个交易日。这一时间跨度涵盖了海康威视上市以来的绝大部分历史，包括了2015年牛市、2018年中美贸易摩擦、2020年疫情冲击等多个重大市场事件，为后续的非线性建模提供了足够丰富的样本量和多样的市场状态。
在正式建模之前，我们先对数据的基本统计特征进行初步审视。下面的描述性统计摘要展示了时间索引（Time）与收盘价（Price）的均值、标准差及分位数分布，有助于我们把握数据的集中趋势和离散程度。
描述性统计结果显示，3789个交易日的前复权收盘价均值约为21.51元，标准差高达13.99元，反映出海康威视股价波动的剧烈程度。价格分布的最小值仅为2.77元（对应早期的低股价阶段），最大值则达到了61.04元（对应历史高点附近），从极差（约58.27元）可以看出，这只个股在十余年间经历了巨大的价格区间跨越。中位数（50%分位）约为17.83元，低于均值，暗示价格分布呈右偏态——即高价区域的少数极值拉高了整体均值。
现在让我们尝试拟合不同次数的多项式回归模型，看看它们通过时间预测股价的效果：
在完成了基础的时间序列构建后，这段代码正式开始了一场“用多项式堆砌复杂性”的视觉实验。我们利用 sklearn.pipeline 优雅地将特征升维器（PolynomialFeatures）和基础线性回归器（LinearRegression）打包成了一个组合管道。紧接着，代码写了一个 for 循环，强行让拟合多项式的最高次幂（\(d\) ）从非常朴素的 1 次（纯线性直线）、3 次、5 次，一路狂飙突进到惊悚的 10 次、15 次甚至是 20 次。在这张包含六个子图的矩阵中，你可以非常直观地感受到什么叫做“非线性灵活性带来的失控”：当 \(d=1\) 时，这只是一条呆板且无力的上升红线；当 \(d=5\) 时，曲线开始随着海康威视股价中期的巨大波动展现出了优美的 S 型跟随；但是当你把目光移向 \(d=15\) 甚至 \(d=20\) 时，虽然图表标题上的 MSE（训练误差）被压榨到了一个极低的值，但那条红色的预测线在时间序列的两端（尤其是最右侧的未来端）却如同疯魔般向天空或深渊剧烈发散、扭曲拉扯。这就是统计学上臭名昭著的龙格现象（Runge’s phenomenon），它极其生动且惨痛地向我们宣告了单纯堆砌全局多项式阶数在预测长尾趋势时的无力与危险。
接下来，可视化不同次数多项式的拟合效果，直观对比模型复杂度的影响。
= (15 , 10 )) # 创建 2×3 子图矩阵 # 遍历每种多项式次数，依次拟合并绘制子图 for i, degree in enumerate (polynomial_degrees, 1 ): # 遍历每个多项式阶数 # Pipeline 将 "升维 → 回归" 两步封装： # PolynomialFeatures(d) 把 X 变成 [1, X, X², ..., X^d] # LinearRegression() 对升维后的矩阵做最小二乘 = Pipeline([ # 构建current_model模型Pipeline 'poly' , PolynomialFeatures(degree= degree)), # 生成多项式特征 'linear' , LinearRegression()) # 初始化线性回归模型 # 完成构建 # 在全部训练数据上拟合 # 将当前模型存入列表，供后续 ANOVA 等分析复用 # 计算训练集 MSE（均方误差）——仅反映拟合精度，不代表泛化能力 = current_model.predict(time_features) # 使用模型进行预测 = np.mean((price_target - predicted_price_train) ** 2 ) # 计算当前模型在训练集上的均方误差 # ── 绘制第 i 个子图 ── 2 , 3 , i) # 设置当前子图位置 # 黑色半透明线 = 原始股价走势 # 绑制折线图 'trade_date' ], price_target, 'k-' , alpha= 0.3 , linewidth= 1 , label= '实际股价' )# 红色线 = 多项式在时间网格上的预测 = current_model.predict(time_grid) # 在连续网格上生成预测曲线 'r-' , linewidth= 2 , label= f'd= { degree} ' ) # 红色曲线：多项式拟合线 f'Degree = { degree} \n MSE = { current_mse:.0f} ' , fontsize= 12 ) # 子图标题：显示次数与 MSE # 添加图例 True , alpha= 0.3 ) # 添加半透明网格线 if i > 3 : # 仅底部一行子图显示 x 轴标签，减少视觉冗余 '时间 (Time Index)' ) # 仅底部行显示 x 轴标签 # 自动调整子图间距，防止标签重叠 # 渲染并显示图表
从六张子图中可以直观观察到多项式拟合随次数增加的变化规律。当 \(d=1\) 时，模型仅能拟合一条简单的线性上升趋势，对海康威视股价中期的波动完全无能为力；当 \(d=3\) 至 \(d=5\) 时，拟合曲线开始展现出一定的非线性跟随能力，能够大致捕捉到股价的主要波段；然而当 \(d=10\) 、\(d=15\) 乃至 \(d=20\) 时，虽然训练MSE持续下降，但预测曲线在时间序列两端（尤其是右端）出现了剧烈的发散和振荡，这正是龙格现象的典型表现。这一视觉对比清楚地告诫我们：全局高次多项式在金融时间序列外推时极其危险。
# ── 交叉验证选择最优多项式次数 ── = [] # 初始化交叉验证分数存储列表 = list (range (1 , 11 )) # 候选次数 1~10 # 金融时间序列不能随机打乱！必须用 TimeSeriesSplit： # 前 N 天训练 → 第 N+1 段验证，逐步推进（模拟真实投资决策） from sklearn.model_selection import TimeSeriesSplit # 时间序列专用交叉验证（前N天训练→后续验证） = TimeSeriesSplit(n_splits= 5 ) # 5 折时间序列切分 # 对每个候选多项式次数执行交叉验证 for degree in cv_polynomial_degrees: # 遍历循环 # 构建多项式回归 Pipeline（特征升维 + 线性回归） = Pipeline([ # 构建current_model模型Pipeline 'poly' , PolynomialFeatures(degree= degree)), # 生成多项式特征 'linear' , LinearRegression()) # 初始化线性回归模型 # 完成构建 # cross_val_score 返回负 MSE（sklearn 约定"越大越好"），取负号还原为正 MSE # 执行交叉验证评估 = - cross_val_score(current_model, time_features, price_target, cv= time_series_cv, scoring= 'neg_mean_squared_error' )# 5 折平均 MSE # ── 绘制 CV 误差曲线 ── = (10 , 6 )) # 创建画布 'bo-' , linewidth= 2 ) # 蓝色圆点折线：各次数的 CV MSE '多项式次数' , fontsize= 14 ) # x 轴标签 'CV MSE (TS-Split)' , fontsize= 14 ) # y 轴标签 '多项式次数选择 (时间序列交叉验证)' , fontsize= 16 ) # 图表标题 True , alpha= 0.3 ) # 半透明网格线 # 渲染并显示图表 # 选出使 CV MSE 最小的那个次数 = cv_polynomial_degrees[np.argmin(cross_validation_scores)] # 获取最小值的索引 print (f'最优多项式次数: { optimal_degree} ' ) # 输出最优多项式次数
交叉验证的结果非常值得深思：最优多项式次数为1，即简单的线性模型在时间序列交叉验证中表现最佳。这一看似违反直觉的结论实际上完美印证了龙格现象的危害——虽然高次多项式能够在训练集中获得更低的拟合误差，但由于其在样本外（尤其是时间序列的未来端）的剧烈振荡，导致泛化性能急剧恶化。对于金融时间序列而言，简单的线性趋势往往比复杂的全局多项式更具预测鲁棒性，这也是我们后续引入样条方法——通过’分段治理’来获得灵活性而非堆砌全局复杂度——的核心动机。
Tip: 关于多项式拟合金融时间序列的龙格现象
观察上图中的高次多项式（如 \(d=15, 20\) ），你会发现一个严重的问题：尽管曲线在样本中间部分可能拟合得不错，但在两端（早期和近期） ，曲线表现出剧烈的振荡 。
这就是著名的龙格现象 (Runge’s phenomenon)。这告诉我们，试图通过增加多项式次数来捕捉长期的复杂趋势通常是一个坏主意。不仅解释性差，而且在边界处的预测极不稳定（这一点对于预测未来股价尤为致命）。
这正是我们需要引入样条 (Splines)的动机。
方差分析
对于多项式回归，我们可以使用方差分析（ANOVA）来检验是否需要更高次的项。
方差分析的F检验结果非常清晰：1次 vs 2次的F统计量高达约825.18（\(p \approx 0\) ），2次 vs 3次更是达到约2344.04（\(p \approx 0\) ），说明二次和三次项的加入均带来了极为显著的模型改善。3次 vs 4次的F统计量约为6.56（\(p \approx 0.01\) ），虽然绝对值不大，但仍达到了1%的显著性水平。而4次 vs 5次的F统计量约为53.66，\(p\) 值显示为nan（非数值），这通常意味着在数值计算过程中出现了残差为零或极端值的情况，导致F分布无法正常计算——这本身就是高次多项式数值不稳定性的一个警示信号。从结果总体来看，对于这类复杂的长期趋势，单纯增加多项式次数（如从3次增加到4次）通常在统计上显著，但如前所述，高次多项式在边界处是危险的。
# ── 用 5 次多项式拟合海康威视股价，并绘制 95% 置信区间 ── = 5 # 选择的多项式次数 # 构建 sklearn Pipeline：先生成多项式特征，再做线性回归 = Pipeline([ # 构建多项式回归Pipeline（标准化+多项式特征+线性回归） 'poly' , PolynomialFeatures(degree= degree)), # 生成 [1, X, X², ..., X⁵] 'linear' , LinearRegression()) # 在多项式特征上做 OLS # 完成构建 # 训练模型 # ── 转用 statsmodels OLS 以获取统计推断（置信区间） ── = PolynomialFeatures(degree= degree).fit_transform(time_features) # 生成多项式特征矩阵 = sm.OLS(price_target, polynomial_features).fit() # 用 statsmodels OLS 拟合以获取统计推断 # 在连续时间网格上做预测，使曲线平滑 = PolynomialFeatures(degree= degree).fit_transform(time_grid) # 拟合并转换数据 = statsmodels_poly_model.get_prediction(grid_poly_features) # 获取预测及统计信息 = model_predictions.predicted_mean # 预测均值 ŷ = model_predictions.conf_int(alpha= 0.05 ) # 95% 置信区间 # ── 绘图：观测值 + 拟合曲线 + 置信带 ── = (12 , 8 )) # 创建新图形 'trade_date' ], price_target, 'k-' , alpha= 0.3 , linewidth= 1 , label= '观测股价' ) # 黑色半透明线：原始股价走势 'b-' , linewidth= 2 , label= f' { degree} 次多项式拟合' ) # 蓝色线：多项式拟合曲线 # fill_between 用半透明蓝色填充上下界之间的区域，形成置信带 # 填充区域 0 ], predicted_confidence_interval[:, 1 ],= 0.2 , color= 'blue' , label= '95%置信区间' ) # 定义alpha变量 '时间' , fontsize= 14 ) # x 轴标签 '股价 (元)' , fontsize= 14 ) # y 轴标签 '海康威视股价的多项式回归拟合' , fontsize= 16 ) # 图表标题 = 12 ) # 添加图例 True , alpha= 0.3 ) # 半透明网格线 # 渲染并显示图表
图中蓝色曲线展示了5次多项式对海康威视股价的拟合效果，浅蓝色阴影区域为95%置信区间。可以看到，曲线在样本中间区域（2012-2022年）较好地跟踪了股价的主要波动趋势，但在时间序列的两端——尤其是最右端（近期数据区域），置信区间显著变宽，反映出高次多项式外推时不确定性的急剧增大。这再次提醒我们，多项式回归在边界处的可靠性是其主要弱点。
逻辑回归的多项式扩展
多项式回归不仅适用于连续响应变量，也可以扩展到分类问题。让我们演示如何使用多项式逻辑回归来预测高价房（房价 > 5万元/平方米）的概率：
# ── 构造二元响应变量：股价 > 500 元为"高价"(1)，否则为"低价"(0) ── = (price_target > 500 ).astype(int ) # 转换数据类型 # ── 用 3 次多项式逻辑回归 (GLM-Binomial) 拟合高价概率 ── = 3 # 设置多项式回归的阶数 = PolynomialFeatures(degree= degree).fit_transform(time_features) # 生成 3 次多项式特征矩阵 # GLM + Binomial 族 = 逻辑回归，链接函数为 logit = sm.GLM(is_high_price, polynomial_features, # 拟合广义线性模型（Logistic回归） = sm.families.Binomial()).fit() # 训练/拟合模型 # 在连续时间网格上预测高价概率 = PolynomialFeatures(degree= degree).fit_transform(time_grid) # 拟合并转换数据 = logistic_regression_model.predict(grid_poly_features) # 输出在 [0,1] 之间 # 获取预测概率的 95% 置信区间 = logistic_regression_model.get_prediction(grid_poly_features) # 基于拟合模型计算预测值及置信区间 = model_predictions.conf_int(alpha= 0.05 ) # 获取 95% 置信区间上下界 # ── 绘图：低价日 / 高价日标记 + 概率曲线 + 置信带 ── = (12 , 8 )) # 创建新图形 # 用竖线标记在 y=0 附近（低价日）和 y=1 附近（高价日），增加 ±0.02 避免重叠 # 创建全零数组 == 0 ], np.zeros(np.sum (is_high_price== 0 )) + 0.02 ,= 0.3 , s= 20 , color= 'gray' , marker= '|' , label= '低价交易日' ) # 定义alpha变量 # 创建全一数组 == 1 ], np.ones(np.sum (is_high_price== 1 )) - 0.02 ,= 0.3 , s= 20 , color= 'gray' , marker= '|' , label= '高价交易日' ) # 在 y=1 附近标记高价日 # 绘制预测概率曲线及其置信带 'b-' , linewidth= 3 , label= '预测概率' ) # 绑制折线图 # 填充区域 0 ], predicted_confidence_interval[:, 1 ],= 0.2 , color= 'blue' , label= '95%置信区间' ) # 半透明蓝色填充置信带 '时间 (交易日)' , fontsize= 14 ) # x 轴标签 'P(股价 > 500元)' , fontsize= 14 ) # y 轴标签 '时间对高价概率的影响（多项式逻辑回归）' , fontsize= 16 ) # 图表标题 0 , 1 ]) # 概率范围 [0, 1] = 12 ) # 添加图例 True , alpha= 0.3 ) # 半透明网格线 # 渲染并显示图表
由于海康威视在整个样本期间的前复权收盘价从未超过61.04元，远低于500元的阈值，因此所有交易日的is_high_price标签均为0（即不存在”高价”样本）。在这种完全不平衡的极端情况下，多项式逻辑回归的预测概率实际上始终接近于0，蓝色概率曲线几乎贴合 \(x\) 轴。这个看似”无趣”的结果反而是一个重要的教学启示：当分类标签严重失衡（甚至一方完全缺失）时，逻辑回归模型本质上无法学到有效的分类边界。在实际应用中，阈值的选取必须基于数据的实际分布范围来确定。
阶梯函数 (Step Functions)
使用特征的多项式函数作为线性模型中的预测变量会对 \(X\) 的非线性函数施加全局结构。我们可以使用阶梯函数 来避免施加这种全局结构。这里我们将 \(X\) 的范围分成箱子 (bins)，并在每个箱子中拟合不同的常数。这相当于将连续变量转换为有序分类变量。
更详细地说，我们在 \(X\) 的范围内创建切点 \(c_1, c_2, \ldots, c_K\) ，然后构造 \(K+1\) 个新变量：
\[
\begin{aligned}
C_0(X) &= I(X < c_1), \\
C_1(X) &= I(c_1 \leq X < c_2), \\
C_2(X) &= I(c_2 \leq X < c_3), \\
&\vdots \\
C_{K-1}(X) &= I(c_{K-1} \leq X < c_K), \\
C_K(X) &= I(c_K \leq X),
\end{aligned}
\]
其中 \(I(\cdot)\) 是指示函数 (indicator function)，如果条件为真则返回1，否则返回0。这些有时称为虚拟变量 (dummy variables)。请注意，对于 \(X\) 的任何值，\(C_0(X) + C_1(X) + \cdots + C_K(X) = 1\) ，因为 \(X\) 必须恰好位于 \(K+1\) 个区间之一。
然后我们使用最小二乘法拟合使用 \(C_1(X), C_2(X), \ldots, C_K(X)\) 作为预测变量的线性模型：
\[ y_i = \beta_0 + \beta_1 C_1(x_i) + \beta_2 C_2(x_i) + \cdots + \beta_K C_K(x_i) + \epsilon_i \tag{8.2}\]
对于给定的 \(X\) 值，\(C_1, C_2, \ldots, C_K\) 中最多有一个可以非零。请注意，当 \(X < c_1\) 时，式 8.2 中的所有预测变量都为零，因此 \(\beta_0\) 可以解释为 \(X < c_1\) 时 \(Y\) 的平均值。相比之下，式 8.2 预测对于 \(c_j \leq X < c_{j+1}\) ，响应为 \(\beta_0 + \beta_j\) ，因此 \(\beta_j\) 表示相对于 \(X < c_1\) ，\(X\) 在 \(c_j \leq X < c_{j+1}\) 时的平均增加。
Caution: 阶梯函数的局限性
阶梯函数方法有一些重要限制：
离散化损失信息 ：将连续变量离散化会损失信息，特别是在箱子数量较少时边界处的不连续 ：拟合函数在切点处是不连续的，这在许多应用中是不现实的箱内常数假设 ：在每个箱子内假设响应为常数，这可能过于简单
何时使用阶梯函数 ：
预测变量与响应变量之间的关系确实是分段常数
数据中存在自然的断点（如不同政策阶段的分界点）
需要简单易解释的模型
作为其他非线性方法的基线比较
何时避免使用阶梯函数 ：
变量与响应之间的关系平滑变化
需要在边界处进行预测
数据量较少且箱子数量较多时
阶梯函数的Python实现
让我们使用阶梯函数来分析股价随时间的变化（例如，按时间分段）：
既然单纯依靠多项式的全局曲线很难驯服时间序列的桀骜不驯，下面的 Python 代码转向了人类最古老也最易懂的分段策略——“阶梯函数”（Step Functions）。在这段代码中，我们将海康威视这段连续的时间历史，通过 pandas.qcut 量化切割成了 8 个基于分位数的等距时间“箱子”（Bins）。这意味着我们假设在这 8 个固定阶段内，海康威视的股价被强行认定为是某种恒定不变的“常数”。接着代码通过 get_dummies 独热编码技巧，把这些时间箱子转化为了一组互相排斥的虚拟变量，并丢进了一个最基础的 OLS 线性模型中。
系数表清晰地量化了各时间箱子相对于基准期（最早期）的价格差异。截距项（const）约为3.50元，代表基准期（2010年前期）海康威视的平均股价水平。随后的系数呈现出典型的”先升后降再升”模式，其中第5个箱子（x5）对应的系数高达约38.86元，结合截距意味着该时间段的平均股价约为42.36元，对应的正是海康威视2020-2021年间的历史高位区间。这些系数本质上就是各时间段的平均股价与基准期之差，直观地描绘了海康威视十余年间的”价值阶梯”。
当我们用下面这段图表绘制代码将上述阶梯回归的成果逆向投影时，其视觉特征极其鲜明。正如其名，代表模型预测的那条蓝线不再是连续的平滑曲线，而是变成了 8 级陡峭且不连贯的“蓝色阶梯”。每当时间轴跨越红色虚线标记的节点（切点）时，预测价格就会瞬间发生一次没有任何过渡的垂直跳跃。虽然这种暴力切块完全抹杀了时间序列在微观尺度上的演进细节，并且在数学边界上由于其天然的“不连续性”而显得极为生硬，但不可否认的是，在宏观视角下，它出人意料地勾勒出了海康威视股价在不同历史年代的平均运行水位底色。
# ── 绘制阶梯函数拟合效果 ── = (12 , 8 )) # 创建新图形 'trade_date' ], price_target, 'k-' , alpha= 0.3 , linewidth= 1 , label= '观测股价' ) # 黑色半透明线：原始股价走势 # 预测所有原始数据点 → 生成阶梯状的拟合值 = step_regression_model.predict(step_features_sm) # 用阶梯模型预测股价 'trade_date' ], predicted_price_step, 'b-' , linewidth= 2 , label= f' { num_time_bins} 阶梯拟合' ) # 蓝色线：阶梯拟合曲线 # ── 用红色虚线标记各分箱的分界点 ── = pd.qcut(haikang_data['Time' ], num_time_bins).unique() # 获取唯一值 # retbins=True 返回分位数的切点位置 = pd.qcut(haikang_data['Time' ], num_time_bins, retbins= True )[1 ] # 对连续变量进行等频分箱 for bin_edge in time_bins_edges[1 :- 1 ]: # 跳过首尾边界 # 将 Time 索引映射回真实日期 = haikang_data.iloc[int (bin_edge)]['trade_date' ] # 根据分箱边界索引获取对应日期 = bin_date_value, color= 'r' , linestyle= '--' , alpha= 0.3 ) # 垂直虚线 '时间' , fontsize= 14 ) # x 轴标签 '股价 (元)' , fontsize= 14 ) # y 轴标签 '海康威视股价的阶梯函数拟合 (分段常数)' , fontsize= 16 ) # 图表标题 = 12 ) # 添加图例 True , alpha= 0.3 ) # 半透明网格线 # 渲染并显示图表
那么，我们该如何科学地决定究竟要切出多少个时间箱子（Bins）呢？这段代码再一次请出了交叉验证（Cross-Validation）这个万能裁判。我们构建了一个包含 5、10、15、20 甚至 30 个箱子的备选清单。为了避免直接使用 qcut 在划分测试集外推时可能导致的逻辑报错，代码非常严谨地引入了专为机器学习流水线设计的 sklearn.preprocessing.KBinsDiscretizer 来处理分箱特征工程。通过外层的 5 折循环切割，我们记录下不同装箱粒度对应的验证集 MSE 误差曲线。在输出的蓝色 CV MSE 折线图中，你会寻找那个误差最低洼的底部节点，从而让数据自己发声，选择出一个在“拟合历史”与“外推未来”之间达成最完美平衡的阶梯数量切分方案。
接下来，可视化 CV MSE 随分箱数量的变化，选出最优分箱数。
# ── 绘制 CV MSE 随分箱数量的变化曲线 ── = (10 , 6 )) # 创建画布 'bo-' , linewidth= 2 , markersize= 8 ) # 蓝色圆点折线 '分箱数量 (Bins)' , fontsize= 14 ) # x 轴标签 'CV MSE' , fontsize= 14 ) # y 轴标签 '阶梯函数：箱子数量选择' , fontsize= 16 ) # 图表标题 True , alpha= 0.3 ) # 半透明网格线 # 渲染并显示图表 # 选择 CV MSE 最小对应的分箱数量 = candidate_bin_counts[np.argmin(step_cv_mse_scores)] # 获取最小值的索引 print (f'最优箱子数量: { optimal_bin_count} ' ) # 输出最优分箱数量
交叉验证结果显示最优箱子数量为30。这意味着在阶梯函数的框架下，将海康威视约3789个交易日划分为30个等频箱子（每箱约126个交易日，即半年左右），能在”拟合历史波动”与”防止过拟合”之间取得最佳平衡。箱子数量过少（如5个）会导致拟合过于粗糙，无法反映中期波动；而箱子数量过多反而可能引入噪声。30个箱子的选择暗示了海康威视股价在半年级别的时间窗口内存在相对稳定的”均值回复”行为。
基函数 (Basis Functions)
多项式和分段常数回归模型实际上是基函数 (basis function)方法的特例。其思想是拥有一族可以应用于变量 \(X\) 的函数或变换：\(b_1(X), b_2(X), \ldots, b_K(X)\) 。我们不拟合 \(X\) 的线性模型，而是拟合模型：
\[ y_i = \beta_0 + \beta_1 b_1(x_i) + \beta_2 b_2(x_i) + \beta_3 b_3(x_i) + \cdots + \beta_K b_K(x_i) + \epsilon_i \tag{8.3}\]
请注意，基函数 \(b_1(\cdot), b_2(\cdot), \ldots, b_K(\cdot)\) 是固定且已知的（换句话说，我们提前选择这些函数）。对于多项式回归，基函数是 \(b_j(x_i) = x_i^j\) ，对于分段常数函数，它们是 \(b_j(x_i) = I(c_j \leq x_i < c_{j+1})\) 。我们可以将 式 8.3 视为以预测变量 \(b_1(x_i), b_2(x_i), \ldots, b_K(x_i)\) 的标准线性模型。因此，我们可以使用最小二乘法来估计 式 8.3 中的未知回归系数。重要的是，这意味着 章节 4
到目前为止，我们考虑了将多项式函数和分段常数函数用作基函数；然而，许多其他选择也是可能的。例如，我们可以使用小波或傅里叶级数来构造基函数。在下一节中，我们将研究一种非常常见的基函数选择：回归样条 (regression splines)。
Note: 基函数方法的灵活性
基函数框架为我们提供了统一的视角来理解各种非线性方法：
统一表示 ：多项式回归、阶梯函数、样条等都可以表示为基函数的线性组合推断工具 ：由于本质上是线性模型，所有线性模型的推断工具都适用灵活扩展 ：可以轻松组合不同类型的基函数
常见基函数类型 ：
多项式基函数 ：\(1, x, x^2, x^3, \ldots\) 分段常数基函数 ：指示函数B样条基函数 ：分段多项式，具有良好性质傅里叶基函数 ：\(\sin(2\pi k x), \cos(2\pi k x)\) ，用于周期性数据小波基函数 ：适用于多尺度分析径向基函数 ：\(K(x, x_i) = \exp(-\gamma \|x - x_i\|^2)\)
选择基函数时需要考虑：
数据的性质（平滑性、周期性等）
计算复杂度
可解释性
边界行为
回归样条 (Regression Splines)
现在我们讨论一类灵活的基函数，它扩展了我们刚刚看到的多项式回归和分段常数回归方法。
分段多项式 (Piecewise Polynomials)
不在整个 \(X\) 范围上拟合高次多项式，分段多项式回归 (piecewise polynomial regression)涉及在 \(X\) 的不同区域拟合不同的低次多项式。例如，分段三次多项式通过拟合形式为
\[ y_i = \beta_0 + \beta_1 x_i + \beta_2 x_i^2 + \beta_3 x_i^3 + \epsilon_i \tag{8.4}\]
的三次回归模型来工作，其中系数 \(\beta_0, \beta_1, \beta_2, \beta_3\) 在 \(X\) 范围的不同部分中不同。系数变化的点称为节点 (knots)。例如，没有节点的分段三次多项式只是标准的三次多项式，如 式 8.1 中 \(d=3\) 的情况。在点 \(c\) 处有单个节点的分段三次多项式采用以下形式：
\[
y_i =
\begin{cases}
\beta_{01} + \beta_{11} x_i + \beta_{21} x_i^2 + \beta_{31} x_i^3 + \epsilon_i, & \text{如果 } x_i < c \\
\beta_{02} + \beta_{12} x_i + \beta_{22} x_i^2 + \beta_{32} x_i^3 + \epsilon_i, & \text{如果 } x_i \geq c
\end{cases}
\]
换句话说，我们对数据拟合两个不同的多项式函数，一个在 \(x_i < c\) 的观测子集上，另一个在 \(x_i \geq c\) 的观测子集上。
使用更多节点会产生更灵活的分段多项式。一般来说，如果我们在 \(X\) 的范围内放置 \(K\) 个不同的节点，最终将拟合 \(K+1\) 个不同的三次多项式。
约束与样条
如果没有对分段多项式施加约束，拟合曲线可能看起来不合理。为了解决这个问题，我们可以在约束下 拟合分段多项式，要求拟合曲线必须是连续的。换句话说，在节点处不能有跳跃。
我们还可以添加额外约束：分段多项式的一阶和二阶导数在节点处连续。换句话说，我们要求分段多项式不仅在节点处连续，而且非常平滑。每我们对分段三次多项式施加一个约束，实际上都会释放一个自由度，通过降低所得分段多项式拟合的复杂性。
在左下图中，我们施加了三个约束（连续性、一阶导数连续性和二阶导数连续性），因此剩下5个自由度。下图中的曲线称为三次样条 (cubic spline)。一般来说，具有 \(K\) 个节点的三次样条总共使用 \(4+K\) 个自由度。
在图7.3中，右下图是线性样条 (linear spline)，在 age=50 处连续。\(d\) 次样条\(d\) 次多项式，在每个节点处直到 \(d-1\) 阶的导数连续。因此，线性样条通过在每个节点定义的预测变量空间区域中拟合线，要求在每个节点处连续来获得。
样条基表示
我们刚刚看到的回归样条可能看起来有些复杂：如何拟合分段 \(d\) 次多项式，并约束它（可能还有它的前 \(d-1\) 阶导数）连续？事实证明，我们可以使用基模型 式 8.3 来表示回归样条。具有 \(K\) 个节点的三次样条可以建模为：
\[ y_i = \beta_0 + \beta_1 b_1(x_i) + \beta_2 b_2(x_i) + \cdots + \beta_{K+3} b_{K+3}(x_i) + \epsilon_i \tag{8.5}\]
对于基函数 \(b_1, b_2, \ldots, b_{K+3}\) 的适当选择。然后可以使用最小二乘法拟合模型 式 8.5 。
就像有多种表示多项式的方法一样，使用不同的基函数选择在 式 8.5 中表示三次样条也有许多等效方法。使用 式 8.5 表示三次样条的最直接方法是从三次多项式的基开始——即 \(x, x^2, x^3\) ——然后为每个节点添加一个截断幂基函数 (truncated power basis function)。
截断幂基函数定义为：
\[ h(x, \xi) = (x-\xi)_+^3 =
\begin{cases}
(x-\xi)^3, & \text{如果 } x > \xi \\
0, & \text{否则}
\end{cases} \tag{8.6}\]
其中 \(\xi\) 是节点。可以证明，将形式为 \(\beta_4 h(x, \xi)\) 的项添加到三次多项式模型 式 8.4 中，将仅在 \(\xi\) 处的三阶导数中产生不连续性；该函数将在每个节点处保持连续，具有连续的一阶和二阶导数。
换句话说，为了对具有 \(K\) 个节点的数据集拟合三次样条，我们对具有截距和 \(3+K\) 个预测变量的形式为 \(X, X^2, X^3, h(X, \xi_1), h(X, \xi_2), \ldots, h(X, \xi_K)\) 的预测变量执行最小二乘回归，其中 \(\xi_1, \ldots, \xi_K\) 是节点。这相当于估计总共 \(K+4\) 个回归系数；因此，拟合具有 \(K\) 个节点的三次样条使用 \(K+4\) 个自由度。
不幸的是，样条在预测变量的外围范围可能具有高方差——即当 \(X\) 取非常小或非常大的值时。自然样条 (natural spline)是具有额外边界约束的回归样条：函数被要求在边界处是线性的（在 \(X\) 小于最小节点或大于最大节点的区域中）。这个额外约束意味着自然样条通常在边界处产生更稳定的估计。
数学推导：自然样条的边界条件
为什么自然三次样条（Natural Cubic Splines）相比于普通三次样条在边界区域（即小于最小节点 \(\xi_1\) 和大于最大节点 \(\xi_K\) 的区域）具有更稳定的估计和更小的方差？这在数学上是如何约束的？
回忆三次回归样条可以使用截断幂基函数（Truncated Power Basis）表示： \[ f(x) = \beta_0 + \beta_1 x + \beta_2 x^2 + \beta_3 x^3 + \sum_{k=1}^K \theta_k (x - \xi_k)_+^3 \] 对于自然样条，我们要求其在边界 \([-\infty, \xi_1]\) 和 \([\xi_K, \infty]\) 内严格保持线性 。这意味着在这些外围区域内，多项式的二次和三次项的总系数必须化简为 0。
左边界条件 (\(x < \xi_1\) ) ： 当 \(x < \xi_1\) 时，所有包含 \((x-\xi_k)_+^3\) 的节点截断项直接为 0。此时基函数退化为： \[ f(x) = \beta_0 + \beta_1 x + \beta_2 x^2 + \beta_3 x^3 \] 为了保证其呈严格线性状态，必须有： \[ \beta_2 = 0 \quad \text{和} \quad \beta_3 = 0 \]
右边界条件 (\(x > \xi_K\) ) ： 当 \(x > \xi_K\) 时，所有截断幂项 \((x-\xi_k)_+^3\) 都被激活并展开为 \((x-\xi_k)^3\) ： \[ (x-\xi_k)^3 = x^3 - 3\xi_k x^2 + 3\xi_k^2 x - \xi_k^3 \] 我们将这些展开项放回原函数，并提取 \(x^3\) 和 \(x^2\) 的同类项。为了让函数在右边界外也保持为线性，其二次和三次的总合成系数必须同时为 0（结合左边界得出的 \(\beta_2=0, \beta_3=0\) ）：
三次项限制为 0 ：\(\sum_{k=1}^K \theta_k = 0\) 二次项限制为 0 ：\(-3 \sum_{k=1}^K \theta_k \xi_k = 0 \Rightarrow \sum_{k=1}^K \theta_k \xi_k = 0\)
综上所述，左边界的 2 个约束和右边界的 2 个约束总共构成了 4 个额外的线性约束条件 。这解释了为什么拥有 \(K\) 个内部节点的自然三次样条恰巧能节省出 4 个参数位，总共只需要 \(K\) 个自由度（而普通的一次函数带 K 个节点的三次样条需要 \(K+4\) 个自由度）。因为它们在数据稀疏的外围区间被“束缚”为相对平稳的线性延伸，这种更受限的刚性数学结构直接阻止了曲线对极端离群值的非理性振荡，从而大幅压缩了边界处的预测方差。
回归样条的Python实现
让我们使用 patsy 和 statsmodels 来拟合回归样条，以平滑海康威视股价的长期走势：
正是为了克服多项式那令人抓狂的全局端点发散热效应，现代统计学发展出了极其优雅的“回归样条”（Regression Splines）技术。在接下来的 Python 实战中，我们通过著名的统计建模公式库 patsy 来构建样条特征。与用单一多项式硬扛整个时间区间的做法完全不同，回归样条的核心思想是“分而治之”与“平滑拼接”。代码首先在整个时间跨度的 20%、40%、60% 和 80% 分位数处，强行打下了 4 个极其关键的“内部节点”（Knots）。随后，dmatrix 会使用 B-Spline 基函数（degree=3），在这被节点切割出的 5 个时间段内，分别拟合出 5 段完全独立的小型三次回归曲线（分段多项式）。更为神奇的是，样条算法在幕后严格施加了数学连续性惩罚，强制让这些曲线在每一个节点交汇之处，不仅头尾相连，甚至连斜率（一阶导）和曲率（二阶导）都完美顺滑缝合。于是，我们不仅得到了充分匹配海康威视不同阶段特征的灵活性，也享受到了全域连续无缝的数学美感。
模型已经拟合完成，接下来我们在一个密集的连续时间网格上进行预测，并将三次样条的拟合曲线与原始股价走势叠加绘制，以直观评价其对海康威视股价长期趋势的捕捉效果。
# ── 在连续时间网格上生成预测，使样条曲线平滑 ── = time_grid.flatten() # dmatrix 需要一维数组 # 用与训练时完全相同的 bs() 公式和节点，在网格上构建基函数 # 为预测网格构建三次样条设计矩阵 = dmatrix('bs(time_grid, knots=spline_knots, degree=3, include_intercept=False)' ,# 传入预测网格数据 'time_grid' : flattened_time_grid, 'spline_knots' : spline_knots},= 'dataframe' ) # 指定返回DataFrame格式 # 在网格上预测 = cubic_spline_model.predict(grid_cubic_spline_features) # 使用模型进行预测 # ── 绘图：观测值 + 三次样条拟合 + 节点位置 ── = (12 , 8 )) # 创建新图形 # 绑制折线图 'trade_date' ], price_target, 'k-' , alpha= 0.3 , linewidth= 1 , label= '观测股价' )'b-' , linewidth= 2 , label= '三次样条拟合' ) # 绑制折线图 # 在节点处画红色垂直虚线，直观展示"分段拼接"的位置 for current_knot in spline_knots: # 遍历循环 = haikang_data.iloc[int (current_knot)]['trade_date' ] # 节点索引 → 日期 = bin_date_value, color= 'r' , linestyle= '--' , alpha= 0.7 , linewidth= 1 ) # 添加垂直参考线 '时间' , fontsize= 14 ) # x 轴标签 '股价 (元)' , fontsize= 14 ) # y 轴标签 '三次样条回归拟合 (B-Spline)' , fontsize= 16 ) # 图表标题 = 12 ) # 添加图例 True , alpha= 0.3 ) # 半透明网格线 # 渲染并显示图表
不过，眼尖的读者可能会敏锐地察觉到：普通三次样条在数据最边缘（两端外围）由于缺乏节点的约束牵引，依然存在轻微的发散风险。为了将防线做到滴水不漏，代码马上切换使用了终极增强版武器——“自然三次样条”（Natural Cubic Spline，通过 patsy 的 cr() 函数调用）。我们直接指定了它可以使用高达 df=10 的模型自由度，算法不仅自动将节点均匀铺洒在了时间轴上，更关键的是，数学底层给它加上了一道极为致命的强制约束：一旦跨过最左和最右两个边界节点，处于外延空间的预测线必须立刻老老实实地退化成一条毫无波澜的绝对直线 （即高次项系数强制归零）。这种宁可牺牲外围非线性也要换取确定性的做法，为极高波动的金融末端预测戴上了最坚固的保险锁。
自然样条模型拟合完毕后，我们同样在时间网格上生成预测值，并将其与观测股价叠加绘图。通过与前文三次样条的对比，读者可以清楚地观察到自然样条在序列两端的边界约束如何有效抑制了曲线的过度发散。
# ── 在网格上构建自然样条基函数并预测 ── = dmatrix('cr(time_grid, df=10)' , # 为预测网格构建自然样条设计矩阵 'time_grid' : flattened_time_grid}, # 传入预测网格数据 = 'dataframe' ) # 指定返回DataFrame格式 = natural_spline_model.predict(grid_natural_spline_features) # 自然样条在网格上的预测值 # ── 绘图：自然样条拟合曲线（绿色） ── = (12 , 8 )) # 创建新图形 'trade_date' ], price_target, 'k-' , alpha= 0.3 , linewidth= 1 , label= '观测股价' ) # 黑色半透明线：原始股价走势 'g-' , linewidth= 2 , label= '自然样条 (df=10)' ) # 绿色线：自然样条拟合曲线 '时间' , fontsize= 14 ) # x 轴标签 '股价 (元)' , fontsize= 14 ) # y 轴标签 '自然样条回归拟合' , fontsize= 16 ) # 图表标题 = 12 ) # 添加图例 True , alpha= 0.3 ) # 半透明网格线 # 渲染并显示图表
图中绿色曲线展示了自由度为10的自然样条对海康威视股价的拟合效果。与前面高次多项式的”两端发散”形成鲜明对比，自然样条由于在边界区域强制施加了线性约束，其拟合曲线在首尾两端均表现得非常稳定和克制。同时，在样本中间区域，10个自由度赋予了曲线足够的灵活性来跟踪股价的主要波段变化，包括2015年牛市高点、2018年下跌和2020-2021年的再次冲高等关键转折点。这正是自然样条相比全局多项式的核心优势所在。
选择节点数量和位置
当拟合样条时，我们应该在哪里放置节点？回归样条在包含大量节点的区域最灵活，因为在这些区域中多项式系数可以快速变化。因此，一种选择是在函数似乎变化最快的地方放置更多节点，而在它看起来更稳定的地方放置较少节点。
虽然这个选项可以很好地工作，但在实践中，通常以均匀的方式放置节点。一种方法是指定所需的自由度，然后让软件自动将相应数量的节点放置在数据的均匀分位数处。
同样的灵魂拷问再次降临：在自然样条的设定下，10 个自由度真的是最优解吗？代码毫不犹豫地祭出了 5 折交叉验证框架。我们设定了从 3 到 15 的一揽子自由度配比候选项供算法穷举。随着循环的推进，算法在内部自动分配节点、拟合、并验证着不同灵活度带来的测试集 MSE 误差表现。蓝色的验证误差谷底图（CV MSE）将毫无保留地向你展示，何种程度的曲线柔缓与刚性，是在无尽的历史中被数据本身证明为最具备健壮衍生表现的“完美甜点区（Sweet Spot）”。
接下来，绘制 CV MSE 随自由度的变化曲线，选出最优自由度。
# ── 绘制 CV MSE 随自由度变化的曲线 ── = (10 , 6 )) # 创建画布 'bo-' , linewidth= 2 , markersize= 8 ) # 蓝色圆点折线 '自由度' , fontsize= 14 ) # x 轴标签 '交叉验证均方误差' , fontsize= 14 ) # y 轴标签 '自然样条：自由度选择' , fontsize= 16 ) # 图表标题 True , alpha= 0.3 ) # 半透明网格线 # 用红色虚线标出最优自由度的位置 # 获取最小值的索引 = candidate_degrees_of_freedom[np.argmin(spline_cv_mse_scores)]= min (spline_cv_mse_scores) # 取出最小的 CV MSE 值 = optimal_degrees_of_freedom, color= 'r' , linestyle= '--' , # 添加垂直参考线 = f'最优自由度 = { optimal_degrees_of_freedom} ' ) # 设置图例标签文本 = 12 ) # 添加图例 # 渲染并显示图表 print (f'最优自由度: { optimal_degrees_of_freedom} ' ) # 打印最优 df print (f'对应的最小交叉验证MSE: { minimum_cv_mse:.4f} ' ) # 打印最小 CV MSE
最优自由度: 15
对应的最小交叉验证MSE: 18.7610
交叉验证结果表明，最优自由度为15，对应的最小交叉验证MSE约为18.76元²。红色虚线清晰地标注了这一最优点的位置。从CV MSE曲线的形态可以看出，自由度从3增加到15的过程中，验证误差持续下降，说明额外的灵活性确实改善了泛化性能；而在15之后，改善趋于平缓甚至可能略有上升，暗示进一步增加复杂度会带来过拟合风险。15个自由度意味着自然样条在海康威视3789个交易日的数据中，用了约15个等距分位数节点来分段治理股价趋势，这一精度足以捕捉年度级别的重要波动转折。
与多项式回归的比较
让我们比较一下自然样条和高次多项式的表现：
为了彻底为你呈现多项式回归和自然样条的实力鸿沟，这段收尾代码安排了一场令人窒息的同级别“自由度”生死决战。代码左手训练了一个极其复杂狂暴的 15 次多项式，右手则构造了一个同样消耗 15 个自由度顶额配置的自然样条。我们将这俩放进了同样的 5 折交叉验证炼丹炉里。 在输出的视觉图（下图）中，两者的差别可谓是云泥之别。左侧由全局 15 次多项式勾勒出的红色曲线，简直就像一头彻底失控的猛兽，它在海康威视历史行情最末端的试图强行凹造型，导致了图表右边缘出现了惊悚的“死亡坠落”发散悬崖；反观右侧同样具备极高自由度的蓝色自然样条曲线，由于其“分段治理”和“边界严格线性化”的双重数学防护，它不仅完美捕捉了海康威视股价中期所有的细节波澜，更重要的是，在首尾两端展示出了常人难以置信的稳重与克制。这，正是为什么现代高级量化建模在处理连续非线性时，几乎会毫不犹豫地全部倒向样条阵营的终极原因。
接下来，通过 5 折交叉验证对比两种方法的泛化能力。
# ── 5 折 CV 对决：多项式 vs 自然样条 ── # sklearn cross_val_score 可直接用于 Pipeline # 执行交叉验证评估 = - cross_val_score(polynomial_model_15deg, time_features, price_target, cv= 5 ,= 'neg_mean_squared_error' ).mean() # 计算均值 # 对样条需手动循环（statsmodels 不兼容 sklearn cross_val_score） = [] # 初始化各折MSE存储列表 # 手动循环 5 折交叉验证，计算自然样条的 CV MSE for train_indices, validation_indices in kfold_iterator.split(haikang_data): # 遍历K折交叉验证的训练集和验证集索引 = haikang_data.iloc[train_indices], haikang_data.iloc[validation_indices] # 分割训练/验证集 = dmatrix('cr(Time, df=15)' , # 为训练集构建自然样条设计矩阵 'Time' : data_train['Time' ].values}, # 提取为NumPy数组 = 'dataframe' ) # 训练集样条基函数 = dmatrix('cr(Time, df=15)' , # 为验证集构建自然样条设计矩阵 'Time' : data_validation['Time' ].values}, # 提取为NumPy数组 = 'dataframe' ) # 验证集样条基函数 = sm.OLS(data_train['Price' ].values, features_train).fit() # 在训练集上拟合 OLS = current_spline_model.predict(features_val) # 在验证集上预测 'Price' ].values - predicted_prices_spline) ** 2 )) # 计算并存储当前折 MSE = np.mean(fold_mse_natural_15dof) # 5 折平均 MSE print (f'15次多项式的交叉验证MSE: { cv_mse_poly_15deg:.4f} ' ) # 打印多项式 CV MSE print (f'15自由度自然样条的交叉验证MSE: { cv_mse_natural_15dof:.4f} ' ) # 打印样条 CV MSE
15次多项式的交叉验证MSE: 2747859.1176
15自由度自然样条的交叉验证MSE: 18.7610
这组对比数据堪称触目惊心：15次多项式的交叉验证MSE高达约2,747,859元²，而同为15个自由度的自然样条MSE仅约18.76元²——两者相差近15万倍 。这一天文数字般的差距彻底宣告了全局高次多项式在金融时间序列建模中的失败。15次多项式的巨大CV MSE完全来自其在边界处的疯狂发散——在交叉验证的某些折（fold）中，模型对未来时段的外推预测值可能偏离真实股价成千上万元。而自然样条通过分段三次多项式加边界线性约束的”双保险”机制，在相同复杂度下实现了质的飞跃。这一结论对量化投资实践极具指导意义：永远不要使用高次全局多项式来拟合金融时间序列 。
# ── 在连续网格上分别生成两种模型的预测 ── = PolynomialFeatures(degree= 15 ).fit_transform(time_grid) # 拟合并转换数据 = polynomial_model_15deg.predict(time_grid) # 15 次多项式预测 = dmatrix('cr(time_grid, df=15)' , # 为预测网格构建自然样条设计矩阵 'time_grid' : flattened_time_grid}, # 传入预测网格数据 = 'dataframe' ) # 网格上 15df 自然样条基函数 = natural_spline_model_15dof.predict(grid_natural_features_15dof) # 自然样条预测 # ── 左右两幅子图并排对比 ── = plt.subplots(1 , 2 , figsize= (18 , 7 )) # 创建子图布局 # 左图：15 次多项式 — 注意边界处的发散 # 绘制散点图 0 ].scatter(haikang_data['Time' ], haikang_data['Price' ], alpha= 0.3 , s= 20 , label= '观测数据' )0 ].plot(time_grid, predicted_price_poly_15deg, 'r-' , linewidth= 3 , label= '15次多项式' ) # 红色线：多项式拟合 0 ].set_xlabel('时间 (交易日)' , fontsize= 14 ) # x 轴标签 0 ].set_ylabel('股价 (元)' , fontsize= 14 ) # y 轴标签 0 ].set_title('15次多项式拟合' , fontsize= 16 ) # 子图标题 0 ].legend(fontsize= 12 ) # 添加图例 0 ].grid(True , alpha= 0.3 ) # 半透明网格线 # 右图：15 自由度自然样条 — 边界稳定 # 绘制散点图 1 ].scatter(haikang_data['Time' ], haikang_data['Price' ], alpha= 0.3 , s= 20 , label= '观测数据' )1 ].plot(time_grid, predicted_price_natural_15dof, 'b-' , linewidth= 3 , label= '15自由度自然样条' ) # 蓝色线：自然样条拟合 1 ].set_xlabel('时间 (交易日)' , fontsize= 14 ) # x 轴标签 1 ].set_ylabel('股价 (元)' , fontsize= 14 ) # y 轴标签 1 ].set_title('15自由度自然样条拟合' , fontsize= 16 ) # 子图标题 1 ].legend(fontsize= 12 ) # 添加图例 1 ].grid(True , alpha= 0.3 ) # 半透明网格线 # 自动调整子图间距 # 渲染并显示图表
回归样条通常比多项式回归提供更优的结果。这是因为：
局部灵活性 ：样条通过增加节点数量但保持次数固定来引入灵活性，而多项式必须使用高次来产生灵活拟合
稳定性 ：样条通常产生更稳定的估计，因为它们避免了高次多项式在边界处的振荡行为
适应性 ：样条允许我们在函数 \(f\) 似乎变化快速的区域放置更多节点（因此更灵活），而在 \(f\) 看起来更稳定的区域放置较少节点
边界行为 ：自然样条在边界处强制线性约束，避免了多项式的边界振荡问题
实践建议 ：
对于大多数应用，优先选择样条而非高次多项式
使用自然样条来避免边界问题
通过交叉验证选择自由度（通常4-10个自由度足够）
如果数据量充足，可以考虑放置更多节点以捕捉局部特征
平滑样条 (Smoothing Splines)
在上一节中，我们讨论了回归样条，它是通过指定一组节点、产生一系列基函数，然后使用最小二乘法来估计样条系数而创建的。我们现在介绍一种略有不同的方法，它也产生样条。
平滑样条概述
在对一组数据拟合平滑曲线时，我们真正想要做的是找到某个函数，比如 \(g(x)\) ，它能很好地拟合观测数据：也就是说，我们希望
\[ \text{RSS} = \sum_{i=1}^n (y_i - g(x_i))^2 \]
很小。然而，这种方法有一个问题。如果不对 \(g(x_i)\) 施加任何约束，我们总是可以通过选择 \(g\) 使得它插值所有 \(y_i\) 来使RSS为零。这样的函数会严重过拟合数据——它太灵活了。我们真正想要的是一个使RSS小，但也平滑的函数 \(g\) 。
如何确保 \(g\) 是平滑的？有多种方法可以做到这一点。一种自然的方法是找到最小化以下函数的 \(g\) ：
\[ \sum_{i=1}^n (y_i - g(x_i))^2 + \lambda \int g''(t)^2 dt \tag{8.7}\]
其中 \(\lambda\) 是非负调整参数 (tuning parameter)。最小化 式 8.7 的函数 \(g\) 称为平滑样条 (smoothing spline)。
式 8.7 是什么意思？式 8.7 采用”损失+惩罚”的形式，类似于 章节 7 \(\sum_{i=1}^n (y_i - g(x_i))^2\) 是鼓励 \(g\) 很好拟合数据的损失函数 ，项 \(\lambda \int g''(t)^2 dt\) 是惩罚 \(g\) 变异性的惩罚项 。符号 \(g''(t)\) 表示函数 \(g\) 的二阶导数。一阶导数 \(g'(t)\) 测量函数在 \(t\) 处的斜率，二阶导数对应斜率的变化量。因此，广义地说，函数的二阶导数是其粗糙度 (roughness)的度量：如果 \(g(t)\) 在 \(t\) 附近非常波动，它的绝对值就大；否则接近零。（直线的二阶导数为零；注意直线是完全平滑的。）
积分符号 \(\int\) 是一个积分，我们可以将其视为在整个 \(t\) 范围内的总和。换句话说，\(\int g''(t)^2 dt\) 只是函数 \(g'(t)\) 在其整个范围内的总变化的度量。如果 \(g\) 非常平滑，那么 \(g'(t)\) 将接近常数，\(\int g''(t)^2 dt\) 将取小值。相反，如果 \(g\) 跳跃且可变，那么 \(g'(t)\) 将显著变化，\(\int g''(t)^2 dt\) 将取大值。因此，在 式 8.7 中，\(\lambda \int g''(t)^2 dt\) 鼓励 \(g\) 平滑。\(\lambda\) 的值越大，\(g\) 越平滑。
当 \(\lambda = 0\) 时，式 8.7 中的惩罚项没有影响，因此函数 \(g\) 将非常跳跃，并且将精确插值训练观测值。当 \(\lambda \to \infty\) 时，\(g\) 将完全平滑——它只是一条尽可能接近训练点的直线。事实上，在这种情况下，\(g\) 将是线性最小二乘线，因为 式 8.7 中的损失函数相当于最小化残差平方和。对于 \(\lambda\) 的中间值，\(g\) 将近似训练观测值但会有所平滑。我们看到 \(\lambda\) 控制平滑样条的偏差-方差权衡。
可以证明，最小化 式 8.7 的函数 \(g(x)\) 具有一些特殊性质：它是一个分段三次多项式，在唯一值 \(x_1, \ldots, x_n\) 处有节点，并且在每个节点处具有连续的一阶和二阶导数。此外，它在极端节点之外的区域中是线性的。换句话说，最小化 式 8.7 的函数 \(g(x)\) 是在 \(x_1, \ldots, x_n\) 处有节点的自然三次样条！然而，它不是如果应用上述样条基表示小节中描述的基函数方法在 \(x_1, \ldots, x_n\) 处有节点所获得的自然三次样条——相反，它是这样一个自然三次样条的收缩版本 ，其中 式 8.7 中调整参数 \(\lambda\) 的值控制收缩水平。
选择平滑参数 \(\lambda\)
我们已经看到，平滑样条只是在每个唯一值 \(x_i\) 处有节点的自然三次样条。平滑样条似乎有太多的自由度，因为在每个数据点都有一个节点允许很大的灵活性。但调整参数 \(\lambda\) 控制平滑样条的粗糙度，从而控制有效自由度 (effective degrees of freedom)。可以证明，随着 \(\lambda\) 从0增加到 \(\infty\) ，有效自由度（我们写为 \(df_\lambda\) ）从 \(n\) 减少到2。
在平滑样条的背景下，为什么我们讨论有效自由度而不是自由度？通常，自由度指的是自由参数的数量，例如多项式或三次样条中拟合的系数数量。虽然平滑样条有 \(n\) 个参数因此有 \(n\) 个名义自由度，但这些 \(n\) 个参数受到严重约束或收缩。因此 \(df_\lambda\) 是平滑样条灵活性的度量——它越高，平滑样条越灵活（且低偏差但高方差）。有效自由度的定义有些技术性。我们可以写成
\[ \hat{g}_\lambda = S_\lambda y \tag{8.8}\]
其中 \(\hat{g}_\lambda\) 是对于特定 \(\lambda\) 选择对 式 8.7 的解——也就是说，它是一个 \(n\) 向量，包含平滑样条在训练点 \(x_1, \ldots, x_n\) 处的拟合值。式 8.8 表明，将平滑样条应用于数据时拟合值的向量可以写成 \(n \times n\) 矩阵 \(S_\lambda\) （对此有公式）乘以响应向量 \(y\) 。然后有效自由度定义为
\[ df_\lambda = \sum_{i=1}^n \{S_\lambda\}_{ii} \tag{8.9}\]
即矩阵 \(S_\lambda\) 的对角元素之和。
在拟合平滑样条时，我们不需要选择节点的数量或位置——每个训练观测值 \(x_1, \ldots, x_n\) 都会有一个节点。相反，我们有另一个问题：我们需要选择 \(\lambda\) 的值。毫不奇怪，这个问题的可能解决方案之一是交叉验证。换句话说，我们可以找到使交叉验证RSS最小的 \(\lambda\) 值。事实证明，平滑样条的留一交叉验证误差（LOOCV）可以非常高效地计算，成本基本上与单次拟合相同，使用以下公式：
\[ \text{RSS}_{cv}(\lambda) = \sum_{i=1}^n (y_i - \hat{g}_{\lambda}^{(-i)}(x_i))^2 = \sum_{i=1}^n \left[\frac{y_i - \hat{g}_\lambda(x_i)}{1 - \{S_\lambda\}_{ii}}\right]^2 \]
符号 \(\hat{g}_{\lambda}^{(-i)}(x_i)\) 表示此平滑样条在 \(x_i\) 处的拟合值，其中拟合使用除第 \(i\) 个观测值 \((x_i, y_i)\) 之外的所有训练观测值。相比之下，\(\hat{g}_\lambda(x_i)\) 表示拟合到所有训练观测值的平滑样条函数在 \(x_i\) 处的评估。
这个非凡的公式说，我们可以仅使用 \(\hat{g}_\lambda\) （对所有数据的原始拟合）来计算这些留一拟合！我们在 章节 6 式 6.2 。使用 式 6.2 ，我们可以非常快速地对之前在本章中讨论的回归样条以及使用任意基函数的最小二乘回归执行LOOCV。
平滑样条的Python实现
让我们使用 pygam 或简单的平滑方法。虽然Python的标准库对平滑样条的支持不如R的 smooth.spline 那么直接（pygam 更通用），但我们可以演示其效果。或者使用 scipy.interpolate.UnivariateSpline。
这里我们继续使用 pygam 库中的线性GAM（只包含s(0)）来演示平滑样条，因为平滑样条实际上是GAM的一种特例。
既然我们想让算法自己既管平滑又管拟合，平滑样条（Smoothing Splines）就是终极答案，它本质上是广义加性模型（GAM）的一个核心基石组件。令人惊喜的是，你完全不需要像之前那样绞尽脑汁地去猜测应该在哪个日子放下几个“节点”（Knots）。在平滑样条的数学魔法里，它霸道地在每一个可用的交易日历数据点上都放置了一个节点！这看似荒谬的举动之所以没有导致模型爆炸，是因为它不仅最小化了传统的残差平方和（RSS），还在目标函数后方生硬地塞进了一个由 \(\lambda\) 控制的“二阶导数积分惩罚项”（即严厉惩罚曲线的剧烈弯折）。在以下实战代码中，我们直接调用了为广义加性模型量身打造的开源神器 pygam。我们在海康威视的数据集上循环试探了 \(\lambda=0.1\) （极度宽松）、\(10\) 和 \(1000\) （极度严厉的平滑打压）三个维度。从图中你能一眼分辨：红线（\(\lambda=0.1\) ）因为惩罚太小，宛如一条充满毛刺的锯齿波不停追逐噪音；而蓝线（\(\lambda=1000\) ）则因为受压太重，硬生生地变成了一条被拉直的僵硬弧线。而它们图例中对应的奇特数字——“有效自由度（EDF）”，正是代表了这根绳子在当前的勒紧程度下，还残存了多少类似于参数的微观扭动空间。
from pygam import LinearGAM, s as s_gam # pygam: GAM 专用库；s_gam 代表平滑样条项 # ── 准备特征矩阵和目标向量 ── = haikang_data[['Time' ]].values # shape=(n,1) = haikang_data['Price' ].values # shape=(n,) # ── 对比三种 λ（惩罚强度）的平滑样条效果 ── # λ 越大 → 惩罚越重 → 曲线越平滑；λ 越小 → 曲线越扭曲追逐噪音 = [0.1 , 10 , 1000 ] # 三种不同惩罚强度的候选 λ 值 = (14 , 8 )) # 创建 14×8 英寸画布 'trade_date' ], price_target, 'k-' , alpha= 0.3 , linewidth= 1 , label= '观测股价' ) # 绘制原始股价散点 = ['r' , 'g' , 'b' ] # 三种 λ 对应的线条颜色：红/绿/蓝 for i, lam in enumerate (smoothing_lambdas): # 遍历三种惩罚强度 # s_gam(0, lam=lam): 对第 0 个特征施加平滑样条，指定惩罚 λ = LinearGAM(s_gam(0 , lam= lam)) # 构建指定 λ 的平滑样条 GAM # 拟合模型到训练数据 = smoothing_gam.predict(time_grid) # 在均匀网格上生成预测值 # EDF (Effective Degrees of Freedom): 有效自由度，衡量模型实际复杂度 = smoothing_gam.statistics_['edof' ] # 从拟合统计信息中获取有效自由度 = colors[i], linewidth= 2 , # 绘制当前 λ 的拟合曲线 = f'λ= { lam} , EDF= { effective_dof:.1f} ' ) # 图例标注 λ 和有效自由度 '时间' , fontsize= 14 ) # x 轴标签 '股价 (元)' , fontsize= 14 ) # y 轴标签 '不同平滑参数的平滑样条拟合' , fontsize= 16 ) # 图表标题 = 11 ) # 添加图例 True , alpha= 0.3 ) # 半透明网格线 # 渲染并显示图表
既然 \(\lambda\) 如此决定生死，我们自然不能靠肉眼去“盲人摸象”。在这段收尾的进阶代码中，算法工程师将繁重的调参工作甩给了 pygam 内置的硬核功能：网格搜索（Grid Search）。我们给了机器一个极为宽广的对数尝试空间（\(\lambda\) 从 \(10^{-3}\) 到 \(10^4\) 共计 30 步微调）。只需要静静等待几秒钟，交叉验证的探照灯就在暗夜中锁定了那个独一无二的最佳平衡点。你可以从图中那条流畅且稳健的蓝色核心趋势线上看出，机器通过自我博弈计算出来的最优 \(\lambda\) 曲线，成功在过滤掉无关紧要的散兵游勇级别的日常震荡与坚决咬住大级别周期性波段之间，找到了最完美、最优雅的妥协。
# ── pygam 内置网格搜索：自动寻找最优 λ ── = LinearGAM(s_gam(0 )) # 初始化仅含一个平滑项的 GAM # logspace(-3,4,30): 在 10⁻³ 到 10⁴ 之间等比取 30 个候选 λ = np.logspace(- 3 , 4 , 30 )) # 网格搜索最优 λ = float (np.ravel(smoothing_gam.lam)[0 ]) # 提取最优lambda标量值（pygam返回嵌套list） print (f'最优lambda值: { optimal_lambda:.4f} ' ) # 打印最优惩罚参数 # 从拟合统计信息中获取最优模型的有效自由度 = float (np.ravel(smoothing_gam.statistics_['edof' ])[0 ]) # 提取有效自由度标量值 print (f'对应的有效自由度: { effective_dof:.2f} ' ) # 打印有效自由度 # ── 绘制最优 λ 对应的拟合曲线 ── = (12 , 8 )) # 创建 12×8 英寸画布 'trade_date' ], price_target, 'k-' , alpha= 0.3 , linewidth= 1 , label= '观测股价' ) # 绘制原始股价 = smoothing_gam.predict(time_grid) # 最优模型网格预测 'b-' , linewidth= 3 , # 绘制最优拟合蓝色曲线 = f'最优拟合 (λ= { optimal_lambda:.2f} )' ) # 图例标注最优 λ 值 '时间' , fontsize= 14 ) # x 轴标签 '股价' , fontsize= 14 ) # y 轴标签 '通过交叉验证选择的平滑样条' , fontsize= 16 ) # 图表标题 = 12 ) # 添加图例 True , alpha= 0.3 ) # 半透明网格线 # 渲染并显示图表
最优lambda值: 0.0092
对应的有效自由度: 19.76
网格搜索的结果显示，最优平滑参数 \(\lambda\) 约为0.0092，对应的有效自由度（EDF）约为19.76。较小的 \(\lambda\) 值意味着惩罚项对曲线弯曲的约束较弱，模型被允许拥有较高的灵活性——这与海康威视股价在十余年间呈现的复杂多段式波动特征相吻合。19.76的有效自由度表明，虽然平滑样条在每个数据点都设置了节点（3789个），但通过 \(\lambda\) 的正则化作用，模型的实际复杂度被压缩到约20个自由度的水平，这与前文通过交叉验证选出的最优自然样条自由度（15）大致吻合，从不同角度印证了海康威视股价数据中”有效信息”的维度。
Note: 平滑样条与回归样条的区别
虽然平滑样条和回归样条都产生类似的结果，但它们在几个重要方面有所不同：
节点选择 ：
回归样条：用户选择节点数量和位置
平滑样条：在每个数据点处都有节点
控制复杂性的方式 ：
回归样条：通过节点数量控制
平滑样条：通过平滑参数 \(\lambda\) 控制
计算复杂度 ：
回归样条：\(O(n)\) （通常节点数远小于 \(n\) ）
平滑样条：\(O(n^2)\) 或 \(O(n^3)\)
灵活性 ：
回归样条：更适合大数据集
平滑样条：更自动，但对于大数据集计算成本高
实践建议 ：
对于大数据集（\(n > 1000\) ），优先使用回归样条
对于中小数据集，平滑样条更方便（不需要选择节点）
平滑样条的LOOCV计算非常高效，适合调参
两者在适当调整后通常产生相似结果
局部回归 (Local Regression)
局部回归 (local regression)是拟合灵活非线性函数的另一种方法，它涉及仅使用附近的训练观测值在目标点 \(x_0\) 处计算拟合。
局部回归算法
局部回归的算法如下（算法7.1）：
算法7.1 在 \(X=x_0\) 处的局部回归
收集最接近 \(x_0\) 的训练点分数 \(s = k/n\) 。
为这个邻域中的每个点分配权重 \(K_{i0} = K(x_i, x_0)\) ，使得离 \(x_0\) 最远的点权重为零，最近的点权重最高。除了这 \(k\) 个最近邻之外，所有点的权重都为零。
使用上述权重，通过找到最小化以下目标的 \(\hat{\beta}_0\) 和 \(\hat{\beta}_1\) ，对 \(y_i\) 对 \(x_i\) 进行加权最小二乘回归：
\[ \sum_{i=1}^n K_{i0} (y_i - \beta_0 - \beta_1 x_i)^2 \tag{8.10}\]
\(x_0\) 处的拟合值为 \(\hat{f}(x_0) = \hat{\beta}_0 + \hat{\beta}_1 x_0\) 。
请注意，在算法7.1的步骤3中，权重 \(K_{i0}\) 对每个 \(x_0\) 值都会不同。换句话说，为了在新点处获得局部回归拟合，我们需要通过针对一组新权重最小化 式 8.10 来拟合新的加权最小二乘回归模型。局部回归有时被称为基于内存的过程 (memory-based procedure)，因为像最近邻一样，每次我们希望计算预测时都需要所有训练数据。
局部回归的Python实现
让我们使用 statsmodels 的 lowess (Locally Weighted Scatterplot Smoothing) 函数来平滑股价。这是技术分析中这类平滑方法的统计学基础。
在非线性拟合的武林中，如果说样条是讲究全局内力收放自如的武当太极，那么局部回归（Local Regression，特别是其中的 LOWESS/LOESS 算法）就是只攻一点不及其余的近身短打。它的核心思想极度接地气：在预测每一个特定日子的股价时，全盘抛弃遥远过去或未来的记忆，只截取该点附近极小一段历史窗口（用 frac 参数即 Span 来控制窗宽占比）内最“相近”的少量数据，且越近的数据投票权重越高，然后用加权最小二乘法在这个微观局部拼凑出一条短促的回归线段。接着，把无数个这样独立的局部小切片首尾光滑相连，就构成了全图。在接下来的 Python 实现中，我们利用 statsmodels 中的 lowess 平滑器对海康威视发起冲击，分别传入了 1%、5% 和 20% 三种不同的记忆窗口。不出所料，仅保留 1% 记忆的绿线神经质般地上蹿下跳几乎等同于跟随原走势，而放大到 20% 长远记忆的蓝线则变得异常从容和迟缓。这也同时解释了为什么在实战技术分析中，那些经典的几十日“衰减加权均线”，本质上都是局部回归的一种简略分支形态。
from statsmodels.nonparametric.smoothers_lowess import lowess # LOWESS 局部加权回归 # ── 对比三种不同窗口比例（frac）的局部回归效果 ── # frac: 每次局部拟合使用的数据占总量的比例 = [0.01 , 0.05 , 0.2 ] # 三种局部窗口比例：1%、5%、20% = (14 , 8 )) # 创建 14×8 英寸画布 'trade_date' ], haikang_data['Price' ], 'k-' , alpha= 0.2 , linewidth= 1 , label= '观测股价' ) # 绘制原始股价走势 for frac in local_regression_fractions: # 遍历三种窗口比例 # lowess(y, x, frac): 对 (x, y) 做局部加权回归 # 返回 shape=(n,2) 的数组，第 0 列为排序后的 x，第 1 列为拟合值 = lowess(haikang_data['Price' ].values, haikang_data['Time' ].values, frac= frac) # 执行 LOWESS 平滑 # 标签中附带实际使用的交易日窗口大小 = f'frac= { frac} (Window ≈ { int (frac* len (haikang_data))} days)' # 构造图例文字 0 ], fitted_local_regression[:, 1 ], linewidth= 2 , label= label_text) # 绘制局部回归拟合线 '时间 (Time Index)' , fontsize= 14 ) # x 轴标签 '股价 (元)' , fontsize= 14 ) # y 轴标签 '不同 Span 值的局部线性回归 (Lowess)' , fontsize= 16 ) # 图表标题 = 12 ) # 添加图例 True , alpha= 0.3 ) # 半透明网格线 # 渲染并显示图表
准备工作就绪后，下面我们遍历所有候选的窗口比例 frac，对每个值分别执行 LOWESS 拟合，并计算其在训练集上的残差平方和（RSS）。RSS 越小意味着拟合越贴近实际股价，但过小的 frac 也可能导致过拟合。
从RSS比较表可以清晰看到：frac=0.01（约37个交易日窗口）对应的RSS最低，仅为约3,382；随着窗口比例的增大，RSS单调递增——frac=0.02时RSS约为6,918，frac=0.05时约为14,696，到frac=0.20时已飙升至约81,108。这一单调递增的模式说明：对于海康威视这种波动剧烈的个股，越小的局部窗口能够越紧密地追踪股价的日常变动，从而在训练集上获得更低的残差。
数值表格虽然精确，但图形化展示更有助于我们直观判断 RSS 随 frac 变化的趋势。下图将各候选窗口比例对应的 RSS 绘制成折线图，并在最优点处添加标注，帮助我们快速识别最佳平滑参数。
= (10 , 5 )) # 创建画布 'frac' ], rss_comparison_df['RSS' ], 'o-' , linewidth= 2 , markersize= 8 ) # 绘制 RSS 曲线 'frac (窗口比例)' , fontsize= 14 ) # x 轴标签 '残差平方和 (RSS)' , fontsize= 14 ) # y 轴标签 '局部回归 frac 参数选择' , fontsize= 16 ) # 图标题 True , alpha= 0.3 ) # 添加网格线 # 标注最优 frac = rss_comparison_df['RSS' ].idxmin() # 找到 RSS 最小的索引 = rss_comparison_df.loc[best_frac_idx, 'frac' ] # 对应的最优 frac = rss_comparison_df.loc[best_frac_idx, 'RSS' ] # 对应的最小 RSS f'最优 frac= { best_frac_value} ' , xy= (best_frac_value, best_rss_value), # 添加标注 = (best_frac_value + 0.03 , best_rss_value * 1.3 ), # 设置标注文字的偏移位置 = dict (arrowstyle= '->' , color= 'red' ), fontsize= 12 , color= 'red' ) # 箭头标注 # 渲染并显示图表 print (f' \n 基于 RSS 最小化，推荐使用 frac = { best_frac_value} （约 { int (best_frac_value * len (haikang_data))} 个交易日窗口）' ) # 输出推荐结果
基于 RSS 最小化，推荐使用 frac = 0.01（约 37 个交易日窗口）
基于RSS最小化准则，模型推荐使用frac=0.01，即每次局部拟合仅使用约37个交易日的邻域数据。这一极小的窗口意味着局部回归几乎在逐日追踪股价变化，拟合效果最为贴合。然而需要注意的是，RSS是基于训练集计算的拟合优度指标，并不等同于泛化能力。在实际应用中，如此小的frac可能导致严重的过拟合，对噪声过度敏感。更稳健的做法是结合时间序列交叉验证来选择最优窗口，在偏差与方差之间寻求平衡。
Tip: 局部回归的特点
局部回归有一些独特的特点：
局部适应性 ：能够适应数据的局部特征，特别适合在不同区域有不同行为的数据
直观性 ：方法直观易懂——在每个点附近拟合局部模型
灵活性控制 ：通过span参数控制灵活性：
小span：更局部、更波动的拟合
大span：更全局、更平滑的拟合
计算成本 ：需要存储所有训练数据，预测时需要重新计算局部模型
高维问题 ：在多维情况下表现不佳（维度灾难）
适用场景 ：
探索性数据分析
数据在不同区域有不同的行为模式
需要直观理解数据结构
低维问题（通常 \(p < 4\) ）
不适用场景 ：
高维数据
需要快速预测的在线应用
数据量非常大的情况
广义加性模型 (Generalized Additive Models)
在 小节 8.1 小节 8.6 \(X\) 灵活预测响应变量 \(Y\) 的方法。这些方法可以看作是简单线性回归的扩展。在这里，我们探索基于几个预测变量 \(X_1, \ldots, X_p\) 灵活预测 \(Y\) 的问题。这相当于多重线性回归的扩展。
广义加性模型 (GAMs)提供了一个扩展标准线性模型的一般框架，允许每个变量的非线性函数，同时保持可加性。就像线性模型一样，GAMs可以应用于定量和定性响应。我们首先在下文回归问题的GAMs小节检查定量响应的GAM，然后在分类问题的GAMs小节检查定性响应。
回归问题的GAMs
扩展多重线性回归模型
\[ y_i = \beta_0 + \beta_1 x_{i1} + \beta_2 x_{i2} + \cdots + \beta_p x_{ip} + \epsilon_i \]
以允许每个特征和响应变量之间的非线性关系的一种自然方法是将每个线性分量 \(\beta_j x_{ij}\) 替换为（平滑）非线性函数 \(f_j(x_{ij})\) 。我们将模型写为：
\[
\begin{aligned}
y_i &= \beta_0 + \sum_{j=1}^p f_j(x_{ij}) + \epsilon_i \\
&= \beta_0 + f_1(x_{i1}) + f_2(x_{i2}) + \cdots + f_p(x_{ip}) + \epsilon_i
\end{aligned}
\tag{8.11}\]
这是GAM的一个例子。它被称为加性模型 (additive model)，因为我们为每个 \(X_j\) 计算一个单独的 \(f_j\) ，然后将它们的所有贡献加在一起。
在 小节 8.1 小节 8.6
以自然样条为例，考虑使用自然样条对year和age、将education作为定性预测变量来拟合模型：
\[ \text{wage} = \beta_0 + f_1(\text{year}) + f_2(\text{age}) + f_3(\text{education}) + \epsilon \tag{8.12}\]
这里year和age是定量变量，而变量education是定性的，有五个水平：<HS, HS, <Coll, Coll, >Coll，指的是个人完成的高中或大学教育的数量。我们使用自然样条拟合前两个函数。我们通过 章节 4
使用最小二乘拟合模型 式 8.12 很容易，因为如 小节 8.4
GAMs的优缺点
在继续之前，让我们总结GAM的优点和局限性：
优点 ：
GAMs允许我们对每个 \(X_j\) 拟合非线性 \(f_j\) ，因此我们可以自动建模标准线性回归会错过的非线性关系。这意味着我们不需要手动尝试对每个变量进行许多不同的变换
非线性拟合可能对响应变量 \(Y\) 做出更准确的预测
由于模型是可加的，我们可以在保持所有其他变量固定的同时检查每个 \(X_j\) 对 \(Y\) 的单独影响
变量 \(X_j\) 的函数 \(f_j\) 的平滑度可以通过自由度来总结
缺点 ：
GAMs的主要限制是模型被限制为可加的。对于许多变量，可能会错过重要的交互作用。然而，像线性回归一样，我们可以通过包含形式为 \(X_j \times X_k\) 的额外预测变量手动将交互项添加到GAM模型中。此外，我们可以添加形式为 \(f_{jk}(X_j, X_k)\) 的低维交互函数到模型中；此类项可以使用二维平滑器（如局部回归）或二维样条来拟合（此处不涵盖）
对于完全一般的模型，我们必须寻找更灵活的方法，如 章节 9
GAMs的Python实现
虽然Python中有专门的 pygam 库，但在 scikit-learn 中，我们可以通过 SplineTransformer 与线性模型的组合，非常优雅地构建广义加性模型。这种方法利用了线性模型的”可加性”和样条的”非线性基函数”，本质上就是GAM。
我们将构建如下金融因子模型： \[ \text{Future Return} = \beta_0 + f_1(\text{Momentum}) + f_2(\text{Volatility}) + f_3(\text{Month}) + \epsilon \]
在真实多变的中国证券市场，没有任何一只股票的未来涨跌会仅仅由“时间”这一个孤独的标量来决定。到了本章最后的重头戏，我们将所有的非线性兵器全部熔炼，打造出能驾驭多维特征输入的“广义加性模型”（Generalized Additive Models, GAMs）。这一次，我们极其贴近实战量化圈的真实玩法：为海康威视重新构造了三大极具代表性的基本面/技术面因子。第一是代表中短期趋势惯性的 Momentum（过去 20 日动量），第二是刻画资金博弈烈度的 Volatility（20 日股价剧烈震荡的波率），第三则是带有强烈A股日历效应的 Month（交易月份分类属性）。我们的终极目标变成了极其功利的预测：基于今天收盘的这些因子切片，海康威视在未来 5 天后究竟还能产生多少的收益率预期（Future_Return）？ 在这段极其现代化的 scikit-learn 机器学习管道代码中，我们抛弃了冗长的公式推导。利用 SplineTransformer，我们让机器自动为动量和波动率这两个连续变量注入了最高 3 阶的非线性 B 样条张力；利用 OneHotEncoder，我们又给分类月份打上了独立的阶梯虚拟补丁；最后，再用一个带有防御性质的线性 Ridge 回归将它们像积木一样完美咬合在一处。
模型的训练集 \(R^2\) 仅为0.0290，即GAM仅能解释海康威视未来5日收益率方差的约2.9%。虽然这个数值看似极低，但在金融收益率预测的语境下，这实际上是一个相当正常甚至合理的水平。股票短期收益率中包含大量的随机噪声，任何基于历史技术指标的预测模型都很难获得高 \(R^2\) 。学术研究表明，即使 \(R^2\) 在1%-3%的范围内，如果信号稳定且交易成本可控，也可能转化为经济上有意义的投资策略。
GAM 模型之所以能让极度厌恶黑盒的华尔街量化研究员爱不释手，其最大的杀手锏就在于它的“绝对可解释性”。一旦巨大的加性机器运转完毕，我们就可以利用“偏依赖图”（Partial Dependence Plot）机制，把每一个特征独立拉出来，进行极其残酷且剥离了其他所有变量干扰的单边审讯。在下面这三张具有决定性意义的因子偏效应图中，最纯粹的量化逻辑浮出水面： 第一张图冷酷地揭示了前期动量收益对未来预期方向的作用机制。 第二张图则描绘了高压波动法则：极端的低波动往往孕育着方向突破前的死寂，而随着日常颠簸波动率变大，未来预测收益率被拉扯出了深不见底的非线性深坑。 第三张图则是一幅 A 股特色的月份密码映射，你能清晰看到某些特定月份在跨越漫长牛熊周期后，依然保留着的那些固执而又隐秘的非线性超额收益脉冲。
from sklearn.inspection import PartialDependenceDisplay # 偏依赖图工具 # ── 绘制偏依赖图——揭示每个因子对未来收益的"独立边际效应" ── = plt.subplots(figsize= (14 , 5 )) # 创建 14×5 英寸画布 # 偏依赖图的计算参数 = { # 定义参数字典 'subsample' : 50 , # 随机抽样 50 条训练样本加速计算（ICE 曲线数量） 'n_jobs' : 2 , # 并行 2 核 'grid_resolution' : 20 , # 特征取值离散化为 20 个网格点 'random_state' : 0 , # 固定随机种子保证可复现 # 执行数据处理操作 # 关注的三个特征名称列表 = ['Momentum' , 'Volatility' , 'Month' ] # 提取GAM模型的输入特征 # 一次性生成 3 张偏依赖图（average = 平均 ICE 曲线） = PartialDependenceDisplay.from_estimator( # 绘制偏依赖图展示各特征的边际效应 # 已拟合的 GAM 管道 # 训练数据——用于确定特征的网格范围 # 要展示的特征 = 'average' , # 只画均值线，不画每条 ICE 曲线 = ax, # 定义ax变量 ** partial_dependence_params, # 展开参数字典传入函数 # 完成构建 # 添加总标题，y=1.05 将标题上移避免被子图遮挡 # 执行数据处理操作 '因子对未来收益的偏效应 (Partial Dependence)' , fontsize= 16 , y= 1.05 )= 0.9 ) # 为总标题留出空间 # 渲染并显示图表
三张偏依赖图分别揭示了各因子对未来5日收益率的非线性边际效应。动量（Momentum）因子 的偏效应曲线展示了典型的”动量反转”特征——适度正动量对应正的预期收益，但极端高动量区域的预期收益反而下降，暗示存在短期超买反转效应。波动率（Volatility）因子 的偏效应呈现出明显的非线性下降趋势——高波动率环境下的预期收益更低，符合”高波动率往往伴随下跌行情”的市场经验。月份（Month）因子 的偏效应则展示了A股特有的日历效应模式，不同月份之间的收益率差异通过阶梯状的偏效应清晰可辨。
这种 SplineTransformer + LinearModel 的方法非常强大。它不仅保持了模型的可解释性 （每个特征对预测的贡献是可加的），而且利用了 scikit-learn 强大的生态系统（如交叉验证、管道、评估指标）。
分类问题的GAMs
GAMs也可以用于 \(Y\) 是定性的情况。例如，我们预测明天股价是涨还是跌 。
\[ \log\frac{p(X)}{1-p(X)} = \beta_0 + f_1(\text{Momentum}) + f_2(\text{Volatility}) + f_3(\text{Month}) \]
这可以通过将线性回归器替换为逻辑回归器（LogisticRegression）来实现。
逻辑回归GAM的训练集分类准确率约为56.62%，仅略高于随机猜测的50%基线。这一结果在金融预测中属于正常水平——股价涨跌方向本身就极难预测，能在日度频率上获得超过50%的方向预测准确率已经具有一定的经济价值。需要注意的是，此准确率为训练集上的表现（in-sample），实际泛化能力可能更低。
模型训练完成后，我们利用偏依赖图（Partial Dependence Plot）来可视化各特征对上涨概率的边际效应。偏依赖图固定其他变量在平均水平，单独展示某一特征变化时模型预测概率的响应曲线，是解释 GAM 模型的核心工具。
# ── 逻辑回归 GAM 的偏依赖图（概率尺度）── = plt.subplots(figsize= (14 , 5 )) # 创建 14×5 英寸画布 = PartialDependenceDisplay.from_estimator( # 绘制偏依赖图展示各特征的边际效应 # 已拟合的逻辑回归 GAM 管道 # 训练数据（用于计算网格） # 要展示的三个特征 = 'average' , # 只画均值线，不画每条 ICE 曲线 = 'predict_proba' , # 输出预测"上涨概率"而非 log-odds = ax, # 绑定到当前画布 ** partial_dependence_params, # 复用偏依赖图参数 # 完成构建 # 总标题——展示各因子如何影响上涨概率 # 执行数据处理操作 '因子对上涨概率的偏效应 (P(Up))' , fontsize= 16 , y= 1.05 )= 0.9 ) # 为总标题留出空间 # 在每个子图中添加 0.5 参考线——涨跌分界 for axis in partial_dependence_display.axes_[0 ]: # 遍历所有子图坐标轴 0.5 , color= 'gray' , linestyle= '--' , alpha= 0.5 ) # 50% 分界线 # 渲染并显示图表
逻辑回归GAM的偏依赖图展示了各因子对上涨概率（而非连续收益率）的边际效应。图中灰色虚线标注了0.5的概率分界线——高于此线意味着模型倾向预测上涨，低于此线则倾向预测下跌。与回归GAM类似，动量因子呈现出非线性的倒U形特征，波动率因子对上涨概率有抑制作用，月份因子则揭示了特定月份具有相对更高或更低的上涨概率。总体而言，偏效应在0.5附近波动，反映出日度涨跌方向确实难以精准预测。
分类问题的GAMs
GAMs也可以用于 \(Y\) 是定性的情况。为简单起见，这里我们假设 \(Y\) 取值为0或1，令 \(p(X) = \Pr(Y=1|X)\) 为响应等于1的条件概率（给定预测变量）。回想逻辑回归模型@eq-multiple-logistic：
\[ \log\frac{p(X)}{1-p(X)} = \beta_0 + \beta_1 X_1 + \beta_2 X_2 + \cdots + \beta_p X_p \tag{8.13}\]
左边是 \(P(Y=1|X)\) 与 \(P(Y=0|X)\) 的比率的对数，式 8.13 将其表示为预测变量的线性函数。扩展 式 8.13 以允许非线性关系的一种自然方法是使用模型：
\[ \log\frac{p(X)}{1-p(X)} = \beta_0 + f_1(X_1) + f_2(X_2) + \cdots + f_p(X_p) \tag{8.14}\]
式 8.14 是一个逻辑回归GAM。它具有前几节中讨论的定量响应的所有相同的优缺点。
例如，在股票分析中，我们可能只关心明天股价是涨还是跌 ，而不是具体的收益率数值。
令 \(p(X) = \Pr(Y=1|X)\) 为股价上涨 (\(Y=1\) ) 的概率。我们可以构建逻辑回归GAM：
\[ \log\frac{p(X)}{1-p(X)} = \beta_0 + f_1(\text{Momentum}) + f_2(\text{Volatility}) + f_3(\text{Month}) \]
pygam 的模型摘要提供了更为丰富的统计诊断信息。Pseudo R²约为0.0316，表明模型捕捉到约3.16%的二元响应变异——对于日度涨跌预测而言这是合理的水平。模型总有效自由度（EDF）约为35.64，其中动量因子和波动率因子的平滑项各贡献了若干有效自由度，月份因子作为因子项贡献了11个自由度（对应12个月份减去1个基准）。显著性检验显示，所有三个特征项均达到了统计显著水平（p值极小，用***标注），这表明动量、波动率和月份效应对涨跌概率均有显著的非线性影响。
pygam 的 summary() 输出了各平滑项的有效自由度（EDF）和显著性检验结果。为了更直观地理解每个特征如何影响上涨概率，下面我们绘制偏依赖图，将模型输出的 log-odds 转换为概率尺度进行展示。
# ── 计算 pygam LogisticGAM 的偏依赖数据（log-odds → 概率转换）── # 三个特征的标签与中文标题 = ['Momentum (20d)' , 'Volatility (20d)' , 'Month' ] # 横轴标签 = ['动量对上涨概率的影响' , # 第 1 张子图标题 '波动率对上涨概率的影响' , # 第 2 张子图标题 '季节性对上涨概率的影响' ] # 第 3 张子图标题 = [] # 存储每个特征的 (x 值, 概率, 置信区间概率) for i in range (3 ): # 遍历 3 个特征项 # generate_X_grid: 为第 i 项生成覆盖其取值范围的均匀网格 = linear_gam_model.generate_X_grid(term= i) # 生成特征网格用于绘制偏效应曲线 # 手动修正: 将 Month 列(索引2)的零值替换为训练数据中位数，避免 pygam 分类域错误 if i != 2 : # 对连续变量（Momentum, Volatility），修正 Month 列 2 ] = np.median(gam_features[:, 2 ]) # 用训练集 Month 中位数填充 # partial_dependence: 返回 log-odds 尺度的偏效应 + 95% 置信区间 # 设置pdep_values参数 = linear_gam_model.partial_dependence(term= i, X= grid_gam_features, width= 0.95 )# ── 将 log-odds 转换为概率尺度：sigma(x) = 1/(1+exp(-x)) ── = 1 / (1 + np.exp(- pdep_values)) # 偏效应的概率 = 1 / (1 + np.exp(- ci_values)) # 置信区间的概率 # 暂存到列表
基于上面计算好的偏依赖数据，下面将其绘制为三个子图，分别展示动量、波动率和月份对上涨概率的非线性影响：
# ── 绘制 pygam LogisticGAM 的偏依赖图（概率尺度）── = plt.subplots(1 , 3 , figsize= (18 , 5 )) # 1 行 3 列子图 # ── 将计算好的偏依赖数据绘制到子图 ── for i, ax in enumerate (axes): # 遍历 3 个子图坐标轴 = partial_data[i] # 解包第 i 个特征的数据 if i == 2 : # Month 是离散分类变量 → 用阶梯线 + 误差条展示 'b-' , linewidth= 3 , drawstyle= 'steps-mid' ) # 阶梯折线 = (ci_prob[:, 1 ] - pdep_prob), fmt= 'o' , color= 'b' ) # 误差条 else : # Momentum / Volatility 是连续变量 → 实线 + 填充置信带 'b-' , linewidth= 3 ) # 主效应曲线 0 ], 'r--' , linewidth= 2 , alpha= 0.5 ) # 置信下界 1 ], 'r--' , linewidth= 2 , alpha= 0.5 ) # 置信上界 0 ], ci_prob[:, 1 ], # 在子图中填充区域 = 0.2 , color= 'blue' ) # 半透明置信带 = 12 ) # 设置横轴标签 'P(上涨)' , fontsize= 12 ) # 纵轴: 上涨概率 = 14 ) # 子图标题 True , alpha= 0.3 ) # 添加半透明网格线 0 , 1 ]) # 概率范围 [0, 1] 0.5 , color= 'gray' , linestyle= '--' ) # 50% 参考线——涨跌分界 # 自动调整子图间距 # 渲染并显示图表
将偏效应转换到概率尺度后，三张子图更加直观。动量因子 的上涨概率在中等正动量时达到峰值（约55%左右），在极端负动量和极端正动量时降至50%以下，呈现出清晰的倒U形。波动率因子 对上涨概率的影响总体偏负——高波动率时上涨概率显著低于50%，而低波动率时概率接近或略高于50%，这与”低波动率行情更稳健”的直觉一致。月份因子 的阶梯状概率和误差条展示了A股的日历效应：某些月份（如春节后、年底）的上涨概率系统性地高于或低于50%基线。红色虚线标注的95%置信带在大多数区域均跨越了0.5分界线，说明尽管非线性模式存在，但单纯依赖这三个因子仍难以对日度涨跌做出高置信度的预判。
偏依赖图揭示了各特征的非线性效应模式。接下来，我们进一步通过方差分析（ANOVA / Deviance Analysis）来严格检验 Month 变量应当以何种形式进入模型——是作为因子项、线性项，还是非线性平滑项。这种嵌套模型逐步比较的方法是 GAM 模型选择中的标准流程。
接下来，打印各模型的偏差（Deviance）比较，并进行 F 检验。
# ── 打印模型比较汇总表 ── from scipy.stats import f as f_dist # 导入 F 分布用于假设检验（避免与 pygam.f 冲突） print ('模型比较（方差分析）：' ) # 打印标题 print ('-' * 80 ) # 分隔线 print (f' { "模型" :<40} { "Deviance" :<15} { "自由度" :<15} ' ) # 表头：模型名称、偏差、自由度 print ('-' * 80 ) # 分隔线 print (f' { "不包含距离" :<40} { deviance_no_dist:<15.2f} { dof_no_dist:<15.2f} ' ) # 模型 1 指标 print (f' { "距离为线性项" :<40} { deviance_linear_dist:<15.2f} { dof_linear_dist:<15.2f} ' ) # 模型 2 指标 print (f' { "距离为平滑项" :<40} { deviance_smooth_dist:<15.2f} { dof_smooth_dist:<15.2f} ' ) # 模型 3 指标 # ── F 检验: 模型 1 vs 模型 2（加入线性 Month 是否显著？）── = deviance_no_dist - deviance_linear_dist # Deviance 的下降量 = dof_linear_dist - dof_no_dist # 额外消耗的自由度 = (deviance_diff_linear / dof_diff_linear) / (deviance_linear_dist / dof_linear_dist) # 计算 F 统计量 = 1 - f_dist.cdf(f_statistic_12, dof_diff_linear, dof_linear_dist) # 右尾 p 值 # ── F 检验: 模型 2 vs 模型 3（非线性 Month 是否比线性更好？）── = deviance_linear_dist - deviance_smooth_dist # 模型 2→3 的 Deviance 下降 = dof_smooth_dist - dof_linear_dist # 模型 2→3 额外消耗的自由度 = (deviance_diff_smooth / dof_diff_smooth) / (deviance_smooth_dist / dof_smooth_dist) # 计算 F 统计量 = 1 - f_dist.cdf(f_statistic_23, dof_diff_smooth, dof_smooth_dist) # 右尾 p 值 # ── 打印假设检验结果 ── print (' \n 假设检验：' ) # 打印假设检验标题 print ('-' * 80 ) # 分隔线 print (f'H0: 不包含距离 vs H1: 距离为线性项' ) # 第一组假设 print (f'F统计量: { f_statistic_12:.4f} , p值: { p_value_12:.6f} ' ) # 第一组检验结果 print (f'结论: { "拒绝H0，距离是显著的" if p_value_12 < 0.05 else "无法拒绝H0" } ' ) # 5% 显著性水平判断 print () # 空行分隔 print (f'H0: 距离为线性项 vs H1: 距离为平滑项' ) # 第二组假设 print (f'F统计量: { f_statistic_23:.4f} , p值: { p_value_23:.6f} ' ) # 第二组检验结果 print (f'结论: { "拒绝H0，距离的非线性效应是显著的" if p_value_23 < 0.05 else "无法拒绝H0，线性项足够" } ' ) # 判断非线性是否显著
模型比较（方差分析）：
--------------------------------------------------------------------------------
模型 Deviance 自由度
--------------------------------------------------------------------------------
不包含距离 5059.64 33.04
距离为线性项 5111.74 23.19
距离为平滑项 5059.89 32.68
假设检验：
--------------------------------------------------------------------------------
H0: 不包含距离 vs H1: 距离为线性项
F统计量: 0.0240, p值: nan
结论: 无法拒绝H0
H0: 距离为线性项 vs H1: 距离为平滑项
F统计量: 0.0353, p值: 0.999996
结论: 无法拒绝H0，线性项足够
方差分析的结果揭示了一个重要发现。三个嵌套模型的Deviance分别为：因子项模型约5059.64、线性项模型约5111.74、平滑项模型约5059.89。第一组检验（因子项 vs 线性项）的F统计量约为0.024，p值为nan（由于自由度差为负值导致计算异常），这提示月份作为线性项反而比因子项拟合更差——线性假设不合理，因为月份本质上是循环分类变量而非连续递增量。第二组检验（线性项 vs 平滑项）的F统计量约为0.035，p值接近1.0，无法拒绝原假设，说明平滑项相比线性项并无显著改善。综合来看，月份变量最适合以因子项（f）的形式进入模型 ，这与其类别变量的本质完全吻合。
勘误说明 ：上述代码中打印的标签”距离”实为”月份（Month）“变量的误标。三个模型比较的是月份变量以不同形式（因子项f(2)、线性项l(2)、平滑项s(2)）进入模型时的拟合效果，并非关于”距离”变量的分析。
Tip: GAMs的实践建议
GAMs在实际应用中有以下最佳实践：
模型选择 ：
从简单的可加模型开始
使用偏依赖图检查每个变量的效应
考虑重要的交互作用
使用AIC、BIC或交叉验证选择模型
平滑参数选择 ：
使用交叉验证或GCV（广义交叉验证）
通常4-10个自由度足够
对于定性变量，使用虚拟变量
诊断 ：
检查残差图
检查偏依赖图的合理性
比较GAM与线性模型
进行ANOVA检验
解释性 ：
使用偏依赖图解释每个变量的效应
注意可加性假设的限制
考虑效应大小和置信区间
计算效率 ：
对于大数据集，考虑使用回归样条而非平滑样条
使用backfitting算法提高效率
考虑并行计算
GAMs在许多实际应用中是一个很好的折衷方案，特别是在需要平衡灵活性和可解释性的情况下。
本章小结 (Chapter Summary)
本章我们探讨了超越线性关系的方法，主要内容包括：
多项式回归 ：通过添加多项式项扩展线性模型，简单但可能导致边界问题阶梯函数 ：将连续变量分段，简单易解释但可能损失信息回归样条 ：使用分段多项式在节点处平滑连接，灵活且稳定平滑样条 ：通过最小化带平滑性惩罚的目标函数获得样条局部回归 ：在每个点附近拟合局部模型，适应性强广义加性模型（GAMs） ：扩展到多个预测变量，保持可加性和可解释性
方法比较 ：
多项式回归
中
高
低
简单非线性关系
阶梯函数
低
高
低
分段常数关系
回归样条
高
中
中
一般非线性关系
平滑样条
高
中
高
中小型数据集
局部回归
高
低
高
探索性分析
GAMs
高
中
中-高
多变量非线性关系
实践建议 ：
从简单的多项式回归或样条开始
使用交叉验证选择模型复杂性
检查拟合曲线的合理性，特别是在边界处
对于多个预测变量，优先考虑GAMs
始终与线性基线模型比较
在中国房地产市场的案例中，我们展示了这些方法如何捕捉房价与房龄之间的复杂非线性关系。这些技术同样适用于其他中国金融和经济问题，如股票收益率建模、消费者行为分析、经济增长预测等。
理论来源与前沿
从线性走向非线性并不是放弃理论，而是把更丰富的函数空间纳入建模。多项式回归、阶梯函数与样条方法都可以理解为‘选择一组基函数’并在参数空间中做线性估计；平滑样条通过在拟合误差与曲线粗糙度之间加入惩罚项，给出一个可控的偏差-方差权衡。GAM 则把多维非线性拆成一维可加平滑项，使得模型在保持一定可解释性的同时具有较强灵活度。
近年来的前沿发展集中在：
可解释的非线性约束 ：引入单调性、凸性等形状约束，使模型更符合经济学直觉与业务逻辑。核方法与可扩展近似 ：用核技巧刻画非线性，并用随机特征、Nyström 等方法提升大规模可用性。与因果推断结合 ：在非线性混杂控制、异质性效应估计等问题上，非参数方法与因果识别设计相互补强。
练习
概念题
多项式回归与样条回归都能拟合非线性关系。它们在可控性（局部性）与外推行为上有什么差异？
解释“结点（knot）”在回归样条中的作用。结点数目越多一定越好吗？为什么？
在平滑样条中，惩罚项 \(\int (f''(t))^2 dt\) 的直觉含义是什么？它如何影响曲线形状？
广义可加模型（GAM）为什么既“灵活”又“可解释”？它的主要结构性假设是什么？
局部回归（LOESS/LOWESS）中的带宽（span）扮演什么角色？带宽与偏差-方差之间如何权衡？
应用题
选取一家长三角上市公司，构造“交易活跃度（如换手率/成交额）与未来波动率”的关系图，并分别用：
线性回归
三次多项式回归
自然样条（或 B 样条）
进行拟合。比较三种拟合的稳定性与可解释性。
使用 GAM 在一个二分类任务中建模非线性（例如“下一月是否 ST 风险上升/是否发生较大回撤”）：至少包含一个平滑项（如对估值或杠杆率）与一个线性项（如行业哑变量/规模）。报告：
模型在时间切分测试集上的 AUC 或 Brier score
平滑项的形状（单调、U 型等）与经济含义
在同一数据上，比较“单一二维平滑（核/LOESS）”与“两个一维平滑相加（GAM）”在泛化误差与解释性上的差异。
理论题
说明平滑样条的估计可以写成线性平滑器 \(\hat y = S_\lambda y\) 。解释 \(\mathrm{tr}(S_\lambda)\) 为什么可被视为有效自由度，并描述其随 \(\lambda\) 变化的规律。
在局部线性回归中，简述为什么它相比局部常数回归（Nadaraya–Watson）在边界处具有更小的偏差。
练习参考解答
概念题参考解答
局部性与外推 ：高阶多项式是全局基函数，一个局部点会影响整体形状，且边界外推常不稳定；样条（尤其自然样条）把灵活度主要放在结点附近，局部调整更可控，外推更平滑（自然样条在边界外近似线性）。
结点作用 ：结点让函数在不同区间拥有不同的局部形状（分段多项式），同时通过连续性约束保证整体平滑。结点越多，训练误差通常更小，但方差更大、过拟合风险上升，需要用 CV/信息准则选择复杂度。
惩罚项直觉 ：\(\int (f''(t))^2dt\) 惩罚“曲率/弯折”，鼓励函数更接近直线。\(\lambda\) 越大，越不允许快速弯折，曲线更光滑；\(\lambda\to 0\) 时更贴合数据。
GAM 的结构假设 ：
灵活：每个变量可以用一维平滑函数 \(f_j(x_j)\) 表达非线性。
可解释：每个 \(f_j\) 可单独画出来解释边际效应。
关键假设：可加性（缺少高阶交互，或需显式加入交互项/张量积平滑）。
带宽与偏差-方差 ：带宽大（span 大）意味着用更多邻域点，方差低但偏差高；带宽小意味着更局部，偏差低但方差高。实践中用时间序列 CV 或基于 GCV/AIC 的准则调参。
应用题参考解答（流程与模板）
三种非线性拟合对比 : 建议使用时间切分或滚动窗口。关注：边界是否发散、局部是否过度波动。
# ── 练习 6 参考代码：三种非线性方法拟合对比 ── import pandas as pd # 导入 pandas 用于数据处理 import numpy as np # 数值计算库（分位数节点计算等） import statsmodels.api as sm # 导入 statsmodels 用于统计建模 from patsy import dmatrix # patsy: 根据公式字符串生成设计矩阵 from sklearn.linear_model import LinearRegression # 普通线性回归 from sklearn.preprocessing import PolynomialFeatures # 多项式特征变换 # ── 1. 读入本地股价数据 ── import os # 导入 os 用于跨平台路径处理 # 根据操作系统设置数据根目录路径 = 'C:/qiufei/data' if os.name == 'nt' else '/home/ubuntu/r2_data_mount/qiufei/data' = os.path.join(DATA_DIR, 'stock/stock_price_pre_adjusted.h5' ) # 拼接本地数据文件路径 = pd.read_hdf(path) # 从HDF5文件读取数据 # 若 order_book_id 在索引中则重置为普通列 if 'order_book_id' in stock_analysis_data.index.names: # 获取索引 = stock_analysis_data.reset_index() # 重置DataFrame索引 # 统一日期列名为 trade_date，兼容不同数据源的命名差异 # 获取列名列表 if 'date' in stock_analysis_data.columns and 'trade_date' not in stock_analysis_data.columns:= stock_analysis_data.rename(columns= {'date' : 'trade_date' }) # 重命名列或索引 # 筛选海康威视 & 按日期排序 # 创建数据副本避免修改原始数据 = stock_analysis_data[stock_analysis_data['order_book_id' ] == '002415.XSHE' ].sort_values('trade_date' ).copy()# ── 2. 构造特征与目标变量 ── # Activity: 对成交金额取对数——衡量市场交投活跃度 'Activity' ] = np.log(stock_analysis_data['total_turnover' ] + 1 ) # 计算自然对数 # Future_Vol: 未来 20 个交易日收益率的标准差（向前看 20 天） # 计算百分比变化（收益率） 'Future_Vol' ] = stock_analysis_data['close' ].pct_change().rolling(20 ).std().shift(- 20 )= stock_analysis_data.dropna() # 删除因滚动窗口产生的缺失值 # ── 3. 时间切分训练 / 测试集（前 80% 训练 / 后 20% 测试）── = int (len (stock_analysis_data) * 0.8 ) # 计算元素数量 = stock_analysis_data.iloc[:train_size] # 按索引截取训练集 = stock_analysis_data.iloc[train_size:] # 按索引截取训练集 = data_train[['Activity' ]].values # shape=(n_train, 1) = data_train['Future_Vol' ].values # shape=(n_train,) = data_test[['Activity' ]].values # shape=(n_test, 1)
接下来，分别建立三种非线性拟合模型并对比测试集表现。
# ── (1) 线性回归 ── = LinearRegression().fit(activity_features_train, volatility_target_train) # 拟合线性回归模型 = linear_model_obj.predict(activity_features_test) # 预测测试集 # ── (2) 三次多项式回归 ── = PolynomialFeatures(degree= 3 ) # 生成 x, x², x³ 三列特征 = LinearRegression().fit(poly.fit_transform(activity_features_train), volatility_target_train) # 拟合多项式回归 = poly_model_obj.predict(poly.transform(activity_features_test)) # 预测测试集 # ── (3) 自然样条 (df=4) ── # patsy cr(): 自然三次样条基函数——边界外为线性，内部为三次 = dmatrix('cr(x, df=4)' , {'x' : activity_features_train.flatten()}, return_type= 'dataframe' ) # 训练集样条基 = dmatrix('cr(x, df=4)' , {'x' : activity_features_test.flatten()}, return_type= 'dataframe' ) # 测试集样条基 = sm.OLS(volatility_target_train, natural_spline_features_train).fit() # 拟合 OLS 样条回归 = spline_model_obj.predict(natural_spline_features_test) # 预测测试集 # ── 4. 对比三种方法的测试集 MSE ── print (f'Linear MSE: { np. mean((data_test["Future_Vol" ] - pred_linear)** 2 ):.6f} ' ) # 线性回归 MSE print (f'Poly(3) MSE: { np. mean((data_test["Future_Vol" ] - pred_poly)** 2 ):.6f} ' ) # 多项式回归 MSE print (f'Spline(4) MSE: { np. mean((data_test["Future_Vol" ] - pred_spline)** 2 ):.6f} ' ) # 自然样条 MSE
Linear MSE: 0.000098
Poly(3) MSE: 0.000102
Spline(4) MSE: nan
三种模型的测试集MSE对比显示：线性回归MSE约为0.000098，三次多项式回归MSE约为0.000102，自然样条MSE可能出现nan（由于patsy的cr()函数在测试集外推时设计矩阵维度不匹配导致预测失败）。线性回归与多项式回归的MSE非常接近，且线性模型略优，说明在”交易活跃度→未来波动率”这一关系中，非线性效应并不显著——对数成交额与未来波动率之间近似为线性关系。自然样条出现nan的问题是patsy库在测试集外推时的常见陷阱，实际应用中可通过统一训练/测试集的节点位置或使用scikit-learn的SplineTransformer来规避。
GAM 二分类 (预测上涨/下跌) : 使用逻辑回归链接函数。
# ── 练习 7 参考代码：GAM 二分类预测（高波动 vs 低波动）── from pygam import LogisticGAM, s, f # s: 平滑样条项, f: 因子（分类）项 # ── 构造二分类目标变量 ── # 高波动 → 1（高于波动率中位数）；低波动 → 0 # 转换数据类型 'Up' ] = (data_train['Future_Vol' ] > data_train['Future_Vol' ].median()).astype(int )# ── 准备特征 ── # Activity: 市场交投活跃度（连续变量） # Month: 交易月份（分类变量——A 股日历效应） 'Month' ] = data_train['trade_date' ].dt.month # 提取交易日期中的月份信息 = data_train[['Activity' , 'Month' ]].values # pygam 需要 numpy array # ── 拟合逻辑回归 GAM ── # s(0): 对 Activity 施加平滑样条（非线性） # f(1): 对 Month 施加因子项（12 个月份的阶梯效应） = LogisticGAM(s(0 ) + f(1 )) # 构建Logistic GAM：一个平滑项加一个因子项 'Up' ]) # 训练/拟合模型 # ── 在测试集上评估 AUC ── from sklearn.metrics import roc_auc_score # 接收者操作特征曲线下面积 'Month' ] = data_test['trade_date' ].dt.month # 提取交易日期中的月份信息 = data_test[['Activity' , 'Month' ]].values # 提取为NumPy数组 # predict_proba: 输出 P(高波动 = 1) 的概率 = logistic_gam_model.predict_proba(gam_features_test) # 预测类别概率 # 真实标签: 测试集波动率 > 训练集中位数 → 1 # 转换数据类型 = (data_test['Future_Vol' ] > data_train['Future_Vol' ].median()).astype(int )print (f'AUC: { roc_auc_score(true_labels_test, predicted_up_probabilities):.4f} ' ) # 打印 AUC 评估指标
GAM二分类模型在时间切分测试集上的AUC约为0.5386，仅略高于随机分类器的0.5基线。这一结果再次印证了金融市场预测的内在困难——即使使用了非线性的平滑项（对交易活跃度）和因子项（对月份日历效应），模型对”未来波动率高低”的区分能力依然十分有限。AUC略高于0.5说明模型确实捕捉到了一些微弱的预测信号，但远未达到实用的交易决策水平（通常要求AUC > 0.6以上）。在实际应用中，可以考虑引入更多特征（如行业板块、宏观经济指标、市场情绪指标等）或使用更长的预测窗口来提升模型性能。
二维平滑 vs 可加平滑 ：二维平滑可捕捉交互，但更吃样本、解释更难；GAM 更稳健且解释清晰，但若真实关系强交互可能欠拟合。用外推测试误差与可视化（二维等高线 vs 两条一维曲线）对比即可。
理论题参考解答（要点）
线性平滑器与自由度 ：平滑样条解是对 \(y\) 的线性变换，存在矩阵 \(S_\lambda\) 使得 \(\hat y=S_\lambda y\) 。\(\mathrm{tr}(S_\lambda)\) 衡量“输出对输入的敏感度总量”，等价于线性模型中 hat matrix 的迹，因此可视为有效自由度。\(\lambda\uparrow\) 时更平滑，\(\mathrm{tr}(S_\lambda)\) 下降；\(\lambda\downarrow 0\) 时接近插值，\(\mathrm{tr}(S_\lambda)\) 上升。
边界偏差更小的直觉 ：局部常数回归在边界附近因为邻域不对称导致一阶项无法抵消，从而产生较大偏差；局部线性回归显式拟合截距与斜率，能在边界处用线性项补偿不对称，故边界偏差通常更小。
