一、Nadaraya-Waston估计

Nadaraya-Waston估计是非参数模型中的经典核回归（Kernel Regression）模型，如无特别说明，一般情况下核回归就是指Nadaraya-Waston估计。它无需对数据的分布进行假设，运用核函数对数据进行估计。

1.1 思想

在回归分析中，假设因变量 Y Y Y可以由自变量 X X X来解释，即
Y = m ( x ) Y = m(x) Y=m(x)
但我们不知道 m m m的具体形式，也是我们要求的模型。

给定观测数据 { X i , Y i } , i = 1 , … , n \{X_i, Y_i\}, i=1,\dots,n {Xi,Yi},i=1,…,n， Y Y Y的条件期望为：
E ( Y ∣ X = x ) = ∫ y f ( y ∣ x ) d y = ∫ y f ( x , y ) f X ( x ) d y = m ( x ) E(Y|X=x) = \int yf(y|x)dy = \int y \frac{f(x,y)}{f_X(x)} dy = m(x) E(Y∣X=x)=∫yf(y∣x)dy=∫yfX(x)f(x,y)dy=m(x)

运用核密度估计，可以得知 f ( x ) f(x) f(x)的概率密度函数（PDF, Probability Density Function）的估计为： f ^ X ( x ) = 1 n ∑ i = 1 n K h ( x − X i ) \hat{f}_X (x) =\frac{1}{n} \sum_{i=1}^n K_h(x-X_i) f^X(x)=n1i=1∑nKh(x−Xi)

联合概率密度函数 f ( x , y ) f(x, y) f(x,y)的核密度估计为： f ^ h , g ( x , y ) = 1 n ∑ i = 1 n K h ( x − X i h ) K g ( y − Y i g ) \hat{f}_{h,g}(x, y) =\frac{1}{n} \sum_{i=1}^n K_h(\frac{x-X_i}{h}) K_g (\frac{y-Y_i}{g}) f^h,g(x,y)=n1i=1∑nKh(hx−Xi)Kg(gy−Yi)
其中 h , g h,g h,g分别为核函数 K h ( ⋅ ) K_h(·) Kh(⋅)和 K g ( ⋅ ) K_g(·) Kg(⋅)的窗宽参数。

因此，前文条件期望中的分子我们可以化简为：
∫ y f ^ h , g ( x , y ) d y = 1 n ∑ i = 1 n 1 h K ( x − X i h ) y g K ( y − Y i g ) d y = 1 n ∑ i = 1 n K h ( x − X i ) ∫ ( s g + Y i ) K ( s ) d s = 1 n ∑ i = 1 n K h ( x − X i ) Y i \begin{array}{l} \int y\hat{f}_{h,g} (x,y)dy = \frac{1}{n} \sum_{i=1}^n \frac{1}{h}K(\frac{x-X_i}{h}) \frac{y}{g} K (\frac{y-Y_i}{g})dy \\ = \frac{1}{n} \sum_{i=1}^n K_h(x-X_i) \int (sg+Y_i)K(s)ds \\ = \frac{1}{n} \sum_{i=1}^n K_h(x-X_i)Y_i \end{array} ∫yf^h,g(x,y)dy=n1∑i=1nh1K(hx−Xi)gyK(gy−Yi)dy=n1∑i=1nKh(x−Xi)∫(sg+Yi)K(s)ds=n1∑i=1nKh(x−Xi)Yi

于是，我们可以得到NW估计
m ^ h ( x ) = n − 1 ∑ i = 1 n K h ( x − X i ) Y i n − 1 ∑ i = 1 n K h ( x − X i ) = 1 n ∑ i = 1 n ( ∑ i = 1 n K h ( x − X i ) n − 1 ∑ i = 1 n K h ( x − X i ) ) Y i = 1 n ∑ i = 1 n W h i ( x ) Y i \begin{array}{l} \hat{m}_h (x) = \cfrac{n^{-1} \sum_{i=1}^n K_h(x-X_i)Y_i}{n^{-1} \sum_{i=1}^n K_h(x-X_i)} \\ =\frac{1}{n} \sum_{i=1}^n (\cfrac{ \sum_{i=1}^n K_h(x-X_i)}{n^{-1} \sum_{i=1}^n K_h(x-X_i)})Y_i \\ =\frac{1}{n} \sum_{i=1}^n W_{hi}(x)Y_i \end{array} m^h(x)=n−1∑i=1nKh(x−Xi)n−1∑i=1nKh(x−Xi)Yi=n1∑i=1n(n−1∑i=1nKh(x−Xi)∑i=1nKh(x−Xi))Yi=n1∑i=1nWhi(x)Yi

从最后一行式子我们可以看到，NW估计可以看成是因变量 Y i Y_i Yi的局部带权平均，权重大小由核函数确定，越靠近估计点 x x x的 Y Y Y所获得的权重就越大。

补充：核函数性质

1、 K h ( ∙ ) = 1 h K ( ∙ h ) K_h( \bullet ) = \frac{1}{h}K(\frac{\bullet}{h}) Kh(∙)=h1K(h∙)
2、对称且正定： K ( u ) ≥ 0 , K ( u ) = K ( − u ) K(u)\geq 0, K(u) = K(-u) K(u)≥0,K(u)=K(−u)
3、积分为1： ∫ K ( u ) = 1 \int K(u) = 1 ∫K(u)=1
4、期望为0： ∫ u K ( u ) d u = 0 \int uK(u)du=0 ∫uK(u)du=0
5、方差有限： 0 < ∫ u 2 K ( u ) d u < ∞ 0<\int u^2 K(u)du< \infty 0<∫u2K(u)du<∞
6、 ∫ K 2 ( u ) d u < ∞ \int K^2 (u)du < \infty ∫K2(u)du<∞

1.2 优缺点

优点

适用范围广，对数据分布没有要求

缺点

计算量较大，尤其在高维数据中
NW估计在数据两侧（x最大和最小）处的估计效果较差

二、python代码实现

2.1 核函数

首先定义一下核函数，这里我们采用常见的高斯核函数：

K ( t ) = 1 2 π e x p ( 1 2 t 2 ) K(t) = \cfrac{1}{\sqrt{2\pi}} exp(\cfrac{1}{2}t^2) K(t)=2π 1exp(21t2),

K h ( t ) = 1 h K ( t h ) K_h (t) = \cfrac{1}{h} K(\cfrac{t}{h}) Kh(t)=h1K(ht).

代码实现：

import numpy as np
import matplotlib.pyplot as pltdef Gauss_kernel(t, h):'''高斯核函数 $K_h (t)$。input:h: 窗宽，即光滑参数，float.t: 位置参数，float.'''x = t / hout = (1 / h) * 1 / np.sqrt(2*np.pi) * np.exp(-0.5*x*x)return out

2.2 矩阵解法

def NW1(X, Y, h):'''Nadaraya-Waston 估计拟合回归函数。矩阵解法。input:X: 自变量, np.array, [number, 1]Y: 因变量，np.array, [number, 1]h: 窗宽，即光滑参数，float.output:y_pred: Y的拟合值。'''number = X.shape[0]one = np.ones(shape = (number, 1))y_pred = np.zeros(shape=Y.shape)for j in range(number):diag_list = []for i in range(number):diag_list.append(Gauss_kernel(X[i] - X[j], h))d_list = [i[0] for i in diag_list]  # 若报错可将 i[0] 改为 iW = np.diag(d_list)y_pred[j] = np.linalg.inv(one.T @ W @one) @ one.T @ W @ Yreturn y_pred

2.3 for循环解法

def NW2(X, Y, h):'''Nadaraya-Waston 估计拟合回归函数。for循环解法。input:X: 自变量, np.array, [number, 1]Y: 因变量，np.array, [number, 1]h: 窗宽，即光滑参数，float.output:y_pred: Y的拟合值。'''number = X.shape[0]y_pred = np.zeros(shape=Y.shape)for i in range(number):up, low = 0, 0for j in range(number):k = Gauss_kernel(X[j] - X[i], h)up += k * Y[j]low += ky_pred[i] = up / lowreturn y_pred

测试一下：

# 生成仿真数据
number = 201
X = np.linspace(0, 2, number)
X = X[:, np.newaxis]
print('X.shape = ', X.shape)np.random.seed(0)
e = np.random.normal(0, 0.2, number)
e = e[:, np.newaxis]
print('e.shape = ', e.shape)Y = 2*X + np.sin(5*np.pi*X) + e
Y_true = 2*X + np.sin(5*np.pi*X)
print('Y.shape = ', Y.shape)# 测试不同窗宽
h_list = [0.1, 0.05, 0.025, 0.01]y_pre = []
for h in h_list:y_pre.append(NW1(X, Y, h))y_pre = np.array(y_pre)# 画图
plt.figure(figsize=(20, 12), dpi=80)
plt.rc('font', family='Times New Roman', weight = 'medium', size=15)  #设置英文字体
font1={'family' : 'Times New Roman', 'size': 15}# plt.ylabel(r"$x_{a,b}^{m,n}$",fontdict=font1)for i in range(len(h_list)):ax = plt.subplot(2, 2, i+1)ax.scatter(X, Y, c='k', s = 3, label='$Simulation$')ax.plot(X, Y_true, lw=1.5, label='$Truth$')ax.plot(X, y_pre[i], lw=3, label='$\hat{m}_h (x)$')ax.set_xlabel('x')ax.set_ylabel('y')ax.legend()ax.set_title('h = %.2f'%h_list[i])# plt.savefig('compare_h.png')
plt.show()

结果如下图：

可以看到窗宽h越大，NW估计曲线越平滑，此时方差小，偏差大；相反地，窗宽h越小，NW估计曲线越不平整，此时方差大，偏差小。这个时候就要在方差和偏差之间寻找一个平衡，即最优窗宽，使得方差和偏差不至于“偏科”太严重。上图中 h = 0.03 h = 0.03 h=0.03就是一个较优的窗宽。寻找最优窗宽的方法有很多，这里不再赘述，感兴趣的同学可以查阅相关资料。

三、总结

其实说白了，NW估计就是个局部带权均值估计，如果样本量足够，在低维数据上表现还不错，但在高维数据上的就可能比不过其他方法了。当然，这只是我的看法，文章若有错误之处，欢迎各位大佬在评论区畅所欲言~

参考资料

[1] Härdle W, Müller M, Sperlich S, et al. Nonparametric and semiparametric models[M]. Berlin: Springer, 2004.

【非参】python实现Nadaraya-Waston估计相关推荐

[竖立正确的贝叶斯三观] 关于predictive distribution 和非参贝叶斯的理解
[竖立正确的贝叶斯三观] 关于predictive distribution 和非参贝叶斯的理解非参贝叶斯框架:在probabilistic graphical models里,参数往往作为随机变量 ...
ESL第五章基扩张和正则化【自然】三次样条/似然比检验/自然正则化、自由度/光滑矩阵/收缩光滑/局部拟合/等价核、非参逻辑回归、多维张量积/加性样条、RKHS/径向基、小波光滑/自适应滤波、B样条
目录 5.1 导言 5.2 分段多项式和样条 5.2.1 自然三次样条Natural Cubic Splines 5.2.2 例子:南非心脏病 5.2.3 例子:音素识别 5.4 光滑样条 5.4.1 ...
三层Dirichlet 过程(非参贝叶斯模型)-来自Machine Learning
本文作者:合肥工业大学管理学院钱洋 email:1563178220@qq.com 内容可能有不到之处,欢迎交流. 未经本人允许禁止转载. 文章来源 Nguyen V A, Boyd-Graber ...
Python+OpenCV：姿态估计（Pose Estimation）
Python+OpenCV:姿态估计(Pose Estimation) ################################################################ ...
abaqus结构工程分析及实例详解pdf_“结构非线性、材料拟合、冲击碰撞、钣金/金属成形、顺序耦合、多物理场、有/非参优化”专题...
课程背景结构的刚强度.材料参数标定.制造可行性.顺序耦合和多物理场耦合,及其性能优化等,是工程领域的产品设计必须考虑的仿真因素.为了让广大仿真分析人员,更好地掌握结构设计与优化的技术能力,全面理解A ...
t检验（独立样本t检验配对样本t检验非参检验多余两组的比较方差分析非参检验）
文章目录独立样本t检验配对样本t检验非参检验多余两组的比较方差分析非参检验独立样本t检验 library(MASS)t.test(Prob~So, data = UScrime,var. ...
随机森林调参 - python
文章目录 1.一般的模型调参原则 2.随机森林的随机性体现在哪几个方面? 2.1 数据集的随机选取 2.2 待选特征的随机选取 3.为什么使用随机森林? 4.随机森林的构建过程 5.随机森林优缺点总结 ...
随机森林算法及贝叶斯优化调参Python实践
1. 随机森林算法 1.1. 集成模型简介集成学习模型使用一系列弱学习器(也称为基础模型或基模型)进行学习,并将各个弱学习器的结果进行整合,从而获得比单个学习器更好的学习效果. 集成学习模型的常见算 ...
计算机视觉--Python实现人体姿态估计
参考: 在 Python 中使用机器学习进行人体姿势估计姿态检测是计算机视觉领域的一个活跃研究领域.https://mp.weixin.qq.com/s/D_sTpTp_pkLeO2nrcjgpaA ...

【非参】python实现Nadaraya-Waston估计

目录