局部近似

2024年12月1日

训练深度神经网络本质上是一个压缩任务。我们希望将训练数据分布表示为由一组矩阵参数化的函数。分布越复杂，所需的参数就越多。近似整个分布的理由是，我们可以在推理时使用相同的模型和权重，对任何有效点进行前向传播。

但是，如果我们的模型在推理时即时训练呢？那么，在传播 $x$ 时，我们只需要对 $x$ 周围的局部分布进行建模。由于局部区域的维度应该比整个训练集低得多，因此一个更简单的模型就足够了！

这就是局部近似或局部回归背后的思想。让我们考虑一个简单的回归任务。

任务

我们得到了以下数据的 $100$ 个样本：

Y = sin (4 X) + ϵ

其中

ϵ \sim N (0, \frac{1}{3})

绘图代码

import numpy as np
import plotly.graph_objects as go

# 生成数据
np.random.seed(42)
n_points = 100
X = np.random.uniform(0, 1, n_points)
epsilon = np.random.normal(0, 1 / 3, n_points)
Y = np.sin(4 * X) + epsilon

# 真实函数
x_true = np.linspace(0, 1, 500)
y_true = np.sin(4 * x_true)

# 创建绘图
fig = go.Figure()

# 添加噪声数据的散点
fig.add_trace(
    go.Scatter(
        x=X,
        y=Y,
        mode="markers",
        name="噪声数据",
        marker=dict(color="gray"),
    )
)

# 添加真实函数
fig.add_trace(
    go.Scatter(
        x=x_true,
        y=y_true,
        mode="lines",
        name="真实函数",
        line=dict(color="red"),
    )
)

# 更新布局
fig.update_layout(
    title="数据",
    xaxis_title="X",
    yaxis_title="Y",
    template="plotly_dark",
    height=400,
    width=730,
)

# 将绘图保存为 HTML 文件
filename = "local_approximation_data.html"
fig.write_html(filename)
print(f"已保存绘图至 {filename}")

# 显示绘图
fig.show()

我们将数据集记为 $D$ ，其中包含样本 $(x_{i}, y_{i}) \in D$ 。

我们的任务是通过数据拟合一条合理的曲线，使其近似匹配真实函数。我们将这条曲线记为 $\hat{f}$ 。

K近邻算法

给定某个 $x$ ，一种方法是取 $x$ 的 $k$ 个最近值 $x_{i}$ ，并将它们的 $y_{i}$ 值平均作为估计值。即，

\hat{f} (x) = Ave (y_{i} ∣ x_{i} \in N_{k} (x))

其中 $N_{k} (x)$ 表示 $x$ 的 $k$ 个最近点。

绘图代码

import plotly.graph_objects as go
import numpy as np

# 生成数据
np.random.seed(42)
n_points = 100
X = np.random.uniform(0, 1, n_points)
epsilon = np.random.normal(0, 1 / 3, n_points)
Y = np.sin(4 * X) + epsilon

# 真实函数
x_true = np.linspace(0, 1, 500)
y_true = np.sin(4 * x_true)

# 对一系列 k 值进行 k-NN
x_curve = np.arange(0, 1, 0.01)
k_range = range(1, 21)
y_curves_knn = {}

for k in k_range:
    y_curve = []
    for x in x_curve:
        distances = np.square(X - x)
        nearest_indices = np.argsort(distances)[:k]
        y_curve.append(np.mean(Y[nearest_indices]))
    y_curves_knn[k] = y_curve

# 创建 Plotly 图形
fig = go.Figure()

# 添加静态轨迹
fig.add_trace(
    go.Scatter(x=X, y=Y, mode="markers", name="噪声数据", marker=dict(color="gray"))
)

fig.add_trace(
    go.Scatter(
        x=x_true, y=y_true, mode="lines", name="真实函数", line=dict(color="red")
    )
)

# 添加第一个 k-NN 曲线（k=13，默认滑块位置）
initial_k = 13
fig.add_trace(
    go.Scatter(
        x=x_curve,
        y=y_curves_knn[initial_k],
        mode="lines",
        name="k-NN 曲线",
        line=dict(color="yellow"),
    )
)

# 定义滑块步骤
steps = []
for k in k_range:
    step = dict(
        method="update",
        args=[
            {"y": [Y, y_true, y_curves_knn[k]]},  # 更新轨迹的 y 数据
            {
                "title": f"交互式 k-NN 曲线，k = {k} 的滑块"
            },  # 动态更新标题
        ],
        label=f"{k}",
    )
    steps.append(step)

# 将滑块添加到布局中
sliders = [
    dict(
        active=initial_k - 1,
        currentvalue={"prefix": "k = "},
        pad={"t": 50},
        steps=steps,
    )
]

fig.update_layout(
    sliders=sliders,
    title=f"交互式 k-NN 曲线，k = {initial_k} 的滑块",
    xaxis_title="X",
    yaxis_title="Y",
    template="plotly_dark",
    height=400,
    width=730,
)

# 显示并保存图形
fig.show()
html_path = "./knn_slider.html"
fig.write_html(html_path)
print(f"已将交互式图形保存到 {html_path}")

通过使用滑块，你可以看到较大的 $k$ 值会导致曲线更平滑，但较小的 $k$ 值曲线会包含一些噪声。在极端情况下， $k = 1$ 会完全追踪训练数据，而 $k = 100$ 会给出一个平坦的全局平均值。

纳达拉亚-沃森核回归

与其将数据子集限制为 $k$ 个点，不如考虑集合中的所有点，但根据每个点与 $x$ 的接近程度来加权其贡献。考虑以下模型

\hat{f} (x) = \frac{\sum _{i = 1}^{∣ D ∣} K _{λ} ( x , x _{i} ) y _{i}}{\sum _{i = 1}^{∣ D ∣} K _{λ} ( x , x _{i} )}

其中 $K_{λ}$ 是一个核函数，我们将用它作为接近度的度量。

K_{λ} (x_{0}, x) = D (\frac{∣ x - x _{0} ∣}{λ})

该函数由 $λ$ 参数化，称为带宽，它控制数据中哪些范围的 $x$ 值会影响 $\hat{f}$ 的输出。如果我们绘制这些函数，这一点就会变得清晰。

核函数

下图绘制的是

f (x) = α K_{λ, D} (0, x)

其中 $α$ 使得 $f$ 在其支持域上积分为 $1$ 。

D (u) = ⎩ ⎨ ⎧ \frac{3}{4} (1 - u^{2}) 0 如果 ∣ u ∣ \leq 1, 如果 ∣ u ∣ > 1.

D (u) = ⎩ ⎨ ⎧ (1 - ∣ u ∣^{3})^{3} 0 如果 ∣ u ∣ \leq 1, 如果 ∣ u ∣ > 1.

D (u) = \frac{1}{2 π} e^{- \frac{1}{2} u^{2}} .

绘图代码

import numpy as np
import plotly.graph_objects as go
from scipy.integrate import quad

# 定义核函数
def epanechnikov_kernel(u):
    return np.maximum(0, 0.75 * (1 - u**2))

def tricube_kernel(u):
    return np.maximum(0, (1 - np.abs(u) ** 3) ** 3)

def gaussian_kernel(u):
    return np.exp(-0.5 * u**2) / np.sqrt(2 * np.pi)

def renormalized_kernel(kernel_func, u_range, bandwidth):
    def kernel_with_lambda(u):
        scaled_u = u / bandwidth
        normalization_factor, _ = quad(lambda v: kernel_func(v / bandwidth), *u_range)
        return kernel_func(scaled_u) / normalization_factor

    return kernel_with_lambda

# 核函数绘图生成器
def generate_kernel_plot(
    kernel_name, kernel_func, x_range, u_range, lambda_values, y_range
):
    fig = go.Figure()

    # 初始 lambda
    initial_lambda = lambda_values[len(lambda_values) // 2]

    # 生成初始核函数曲线
    x = np.linspace(*x_range, 500)
    kernel_with_lambda = renormalized_kernel(kernel_func, u_range, initial_lambda)
    y = kernel_with_lambda(x)
    fig.add_trace(
        go.Scatter(
            x=x,
            y=y,
            mode="lines",
            name=f"{kernel_name} 核函数 (λ={initial_lambda:.2f})",
            line=dict(color="green"),
        )
    )

    # 创建滑块的帧
    frames = []
    for bandwidth in lambda_values:
        kernel_with_lambda = renormalized_kernel(kernel_func, u_range, bandwidth)
        y = kernel_with_lambda(x)
        frames.append(
            go.Frame(
                data=[
                    go.Scatter(
                        x=x,
                        y=y,
                        mode="lines",
                        name=f"{kernel_name} 核函数 (λ={bandwidth:.2f})",
                        line=dict(color="green"),
                    )
                ],
                name=f"{bandwidth:.2f}",
            )
        )

    # 将帧添加到图中
    fig.frames = frames

    # 添加滑块
    sliders = [
        {
            "active": len(lambda_values) // 2,
            "currentvalue": {"prefix": "带宽 λ: "},
            "steps": [
                {
                    "args": [
                        [f"{bandwidth:.2f}"],
                        {"frame": {"duration": 0, "redraw": True}, "mode": "immediate"},
                    ],
                    "label": f"{bandwidth:.2f}",
                    "method": "animate",
                }
                for bandwidth in lambda_values
            ],
        }
    ]

    # 更新布局
    fig.update_layout(
        title=f"{kernel_name} 核函数",
        xaxis_title="u",
        yaxis_title="K(u)",
        yaxis_range=y_range,
        template="plotly_dark",
        sliders=sliders,
        height=400,  # 调整为与之前大小匹配
        width=730,  # 调整为与之前大小匹配
        updatemenus=[
            {
                "direction": "left",
                "pad": {"r": 10, "t": 87},
                "showactive": False,
                "type": "buttons",
                "x": 0.1,
                "xanchor": "right",
                "y": 0,
                "yanchor": "top",
            }
        ],
    )

    return fig

# 核函数
kernels = {
    "Epanechnikov": epanechnikov_kernel,
    "Tricube": tricube_kernel,
    "Gaussian": gaussian_kernel,
}

# 参数
x_range_plot = (-3, 3)  # u 值的绘图范围
u_range_integration = (-3, 3)  # 归一化范围
lambda_values = np.linspace(0.01, 2, 20)  # 从 0.01 到 2 的线性 lambda 值
y_range_plot = (0, 1.5)  # 调整范围以适应归一化函数

# 为每个核函数生成并显示图表
for kernel_name, kernel_func in kernels.items():
    fig = generate_kernel_plot(
        kernel_name,
        kernel_func,
        x_range_plot,
        u_range_integration,
        lambda_values,
        y_range_plot,
    )

    # 将图表保存为 HTML 文件
    filename = f"{kernel_name}_dynamic_normalization_kernel_function.html"
    fig.write_html(filename, auto_play=False)
    print(f"已保存 {kernel_name} 核函数图表到 {filename}")

    # 显示图表
    fig.show()

结果

我们现在绘制每个核函数的结果。每个图都有一个 $λ$ 滑块，可以实时控制输出。

绘图代码

import numpy as np
import plotly.graph_objects as go

# 定义核函数
def epanechnikov_kernel(u):
    return np.maximum(0, 0.75 * (1 - u**2))

def tricube_kernel(u):
    return np.maximum(0, (1 - np.abs(u) ** 3) ** 3)

def gaussian_kernel(u):
    return np.exp(-0.5 * u**2) / np.sqrt(2 * np.pi)

# 核回归函数
def kernel_regression(X, Y, x_curve, kernel_func, bandwidth):
    y_curve = []
    for x in x_curve:
        distances = np.abs(X - x) / bandwidth
        weights = kernel_func(distances)
        weighted_average = (
            np.sum(weights * Y) / np.sum(weights) if np.sum(weights) > 0 else 0
        )
        y_curve.append(weighted_average)
    return y_curve

# 生成数据
np.random.seed(42)
n_points = 100
X = np.random.uniform(0, 1, n_points)
epsilon = np.random.normal(0, 1 / 3, n_points)
Y = np.sin(4 * X) + epsilon

# 真实曲线
x_true = np.linspace(0, 1, 500)
y_true = np.sin(4 * x_true)

# 核估计的点
x_curve = x_true

# 核函数
kernels = {
    "Epanechnikov": epanechnikov_kernel,
    "Tricube": tricube_kernel,
    "Gaussian": gaussian_kernel,
}

# 滑块的对数空间带宽范围
lambda_values = np.logspace(-2, 0, 20)  # 从 0.01 到 1

# 为每个核生成单独的图
for kernel_name, kernel_func in kernels.items():
    fig = go.Figure()

    # 添加噪声数据的散点图
    fig.add_trace(
        go.Scatter(
            x=X, y=Y, mode="markers", name="噪声数据", marker=dict(color="gray")
        )
    )

    # 添加真实函数
    fig.add_trace(
        go.Scatter(
            x=x_true,
            y=y_true,
            mode="lines",
            name="真实函数",
            line=dict(color="red"),
        )
    )

    # 添加初始核曲线
    initial_bandwidth = lambda_values[0]
    y_curve = kernel_regression(X, Y, x_curve, kernel_func, initial_bandwidth)
    fig.add_trace(
        go.Scatter(
            x=x_curve,
            y=y_curve,
            mode="lines",
            name=f"Nadaraya-Watson ({kernel_name})",
            line=dict(color="green"),
        )
    )

    # 创建滑块的帧
    frames = []
    for bandwidth in lambda_values:
        y_curve = kernel_regression(X, Y, x_curve, kernel_func, bandwidth)
        frames.append(
            go.Frame(
                data=[
                    go.Scatter(
                        x=X,
                        y=Y,
                        mode="markers",
                        name="噪声数据",
                        marker=dict(color="gray"),
                    ),
                    go.Scatter(
                        x=x_true,
                        y=y_true,
                        mode="lines",
                        name="真实函数",
                        line=dict(color="red"),
                    ),
                    go.Scatter(
                        x=x_curve,
                        y=y_curve,
                        mode="lines",
                        name=f"Nadaraya-Watson ({kernel_name})",
                        line=dict(color="green"),
                    ),
                ],
                name=f"{bandwidth:.2f}",
            )
        )

    # 将帧添加到图中
    fig.frames = frames

    # 添加滑块
    sliders = [
        {
            "active": 0,
            "currentvalue": {"prefix": "带宽 λ: "},
            "steps": [
                {
                    "args": [
                        [f"{bandwidth:.2f}"],
                        {"frame": {"duration": 0, "redraw": True}, "mode": "immediate"},
                    ],
                    "label": f"{bandwidth:.2f}",
                    "method": "animate",
                }
                for bandwidth in lambda_values
            ],
        }
    ]

    # 更新布局
    fig.update_layout(
        title=f"Nadaraya-Watson 核回归 ({kernel_name} 核)",
        xaxis_title="X",
        yaxis_title="Y",
        template="plotly_dark",
        sliders=sliders,
        height=400,
        width=730,
        updatemenus=[
            {
                "buttons": [
                    {
                        "args": [
                            None,
                            {
                                "frame": {"duration": 500, "redraw": True},
                                "fromcurrent": True,
                            },
                        ],
                        "label": "播放",
                        "method": "animate",
                    },
                    {
                        "args": [
                            [None],
                            {
                                "frame": {"duration": 0, "redraw": True},
                                "mode": "immediate",
                            },
                        ],
                        "label": "暂停",
                        "method": "animate",
                    },
                ],
                "direction": "left",
                "pad": {"r": 10, "t": 87},
                "showactive": False,
                "type": "buttons",
                "x": 0.1,
                "xanchor": "right",
                "y": 0,
                "yanchor": "top",
            }
        ],
    )

    # 将图保存为 HTML 文件
    filename = f"{kernel_name}_kernel_regression.html"
    fig.write_html(filename, auto_play=False)
    print(f"已保存 {kernel_name} 核图到 {filename}")

    # 显示图
    fig.show()

我们可以看到，数据的简单加权平均值能够很好地模拟正弦曲线。

局部线性回归

在Nadaraya-Watson核回归中，我们通过核函数 $K_{λ}$ 定义的邻域内进行加权平均。这种方法的一个潜在问题是局部邻域内的平滑插值，因为我们实际上并没有假设该区域遵循任何模型。

如果我们假设每个区域都是局部线性的呢？那么，我们可以求解最小二乘拟合并自由插值！

区域：$k$-近邻

让我们将局部区域定义为输入点的 $k$ 个最近邻。设 $X = [N_{k} (x_{0}) 1]$ ， $Y$ 为对应的 $y$ 值。最小二乘拟合系数为

β = (X^{⊤} X)^{- 1} XY

绘图代码

import plotly.graph_objects as go
import numpy as np

# 生成数据
np.random.seed(42)
n_points = 100
X = np.random.uniform(0, 1, n_points)
epsilon = np.random.normal(0, 1 / 3, n_points)
Y = np.sin(4 * X) + epsilon

# 真实函数
x_true = np.linspace(0, 1, 500)
y_true = np.sin(4 * x_true)

# k-近邻局部线性回归
def knn_linear_regression(X, Y, x_curve, k_range):
    y_curves = {}
    for k in k_range:
        y_curve = []
        for x in x_curve:
            # 找到 k 个最近邻
            distances = np.abs(X - x)
            nearest_indices = np.argsort(distances)[:k]

            # 选择 k 个最近邻
            X_knn = X[nearest_indices]
            Y_knn = Y[nearest_indices]

            # 为 k-近邻创建设计矩阵
            X_design = np.vstack((np.ones_like(X_knn), X_knn)).T

            # 使用普通最小二乘法求解 beta
            beta = np.linalg.pinv(X_design.T @ X_design) @ X_design.T @ Y_knn

            # 预测 y 值
            y_curve.append(beta[0] + beta[1] * x)
        y_curves[k] = y_curve
    return y_curves

# 公共变量
x_curve = np.arange(0, 1, 0.01)
k_range = range(1, 21)  # k 的取值范围从 1 到 20
initial_k = 10  # k 的默认值

# 使用 k-近邻计算局部线性回归
y_curves_knn = knn_linear_regression(X, Y, x_curve, k_range)

# 创建 Plotly 图表
fig = go.Figure()

# 添加静态轨迹
fig.add_trace(
    go.Scatter(x=X, y=Y, mode="markers", name="噪声数据", marker=dict(color="gray"))
)

fig.add_trace(
    go.Scatter(
        x=x_true, y=y_true, mode="lines", name="真实函数", line=dict(color="red"))
)

# 添加第一个 k-近邻曲线 (k=initial_k)
fig.add_trace(
    go.Scatter(
        x=x_curve,
        y=y_curves_knn[initial_k],
        mode="lines",
        name="k-近邻曲线",
        line=dict(color="yellow"),
    )
)

# 定义滑块步骤
steps = []
for k in k_range:
    step = dict(
        method="update",
        args=[
            {"y": [Y, y_true, y_curves_knn[k]]},  # 更新轨迹的 y 数据
            {
                "title": f"k-近邻局部线性回归曲线，k = {k}"
            },  # 动态更新标题
        ],
        label=f"{k}",
    )
    steps.append(step)

# 将滑块添加到布局中
sliders = [
    dict(
        active=k_range.index(initial_k),  # 使用 initial_k 的索引
        currentvalue={"prefix": "k = "},
        pad={"t": 50},
        steps=steps,
    )
]

fig.update_layout(
    sliders=sliders,
    title=f"k-近邻局部线性回归曲线，k = {initial_k}",
    xaxis_title="X",
    yaxis_title="Y",
    template="plotly_dark",
    height=400,
    width=730,
)

# 显示并保存图表
fig.show()
html_path = "./knn_slider_llr.html"
fig.write_html(html_path)
print(f"已将交互式 k-近邻图表保存至 {html_path}")

我们可以看到，当 $k$ 较小时，输出可能会显得相当粗糙。

区域：核函数

或许我们可以借鉴 Nadaraya-Watson 核函数的一些思想。我们希望不同程度地考虑训练集中的所有点，局部区域内的点赋予较高权重，区域外的点赋予较低权重。

为此，我们可以使用加权最小二乘目标函数，权重为 $W (x_{0}) = diag (K_{λ} (x_{0}, x_{i}))$ 。其解为

β = (X^{⊤} WX)^{- 1} X^{⊤} WY

绘制不同核函数 $D$ 的结果：

绘图代码

import plotly.graph_objects as go
import numpy as np

# 生成数据
np.random.seed(42)
n_points = 100
X = np.random.uniform(0, 1, n_points)
epsilon = np.random.normal(0, 1 / 3, n_points)
Y = np.sin(4 * X) + epsilon

# 真实函数
x_true = np.linspace(0, 1, 500)
y_true = np.sin(4 * x_true)

# 核函数
def gaussian_kernel(u):
    return np.exp(-0.5 * u**2)

def epanechnikov_kernel(u):
    return np.maximum(0, 1 - u**2)

def tricube_kernel(u):
    return np.maximum(0, (1 - np.abs(u) ** 3) ** 3)

# 局部线性回归，针对特定核函数
def local_linear_regression(X, Y, x_curve, bandwidths, kernel):
    y_curves = {}
    for λ in bandwidths:
        λ_rounded = round(λ, 2)
        y_curve = []
        for x in x_curve:
            # 使用指定核函数计算权重
            distances = (X - x) / λ
            weights = kernel(distances)
            W = np.diag(weights)

            # 创建设计矩阵
            X_design = np.vstack((np.ones_like(X), X)).T

            # 使用加权最小二乘法求解 beta
            beta = np.linalg.pinv(X_design.T @ W @ X_design) @ X_design.T @ W @ Y

            # 预测 y 值
            y_curve.append(beta[0] + beta[1] * x)
        y_curves[λ_rounded] = y_curve
    return y_curves

# 公共变量
x_curve = np.arange(0, 1, 0.01)
bandwidths = np.linspace(0.05, 0.5, 20)
initial_λ = bandwidths[len(bandwidths) // 2]

# 为每个核函数生成绘图
kernels = {
    "Gaussian Kernel": gaussian_kernel,
    "Epanechnikov Kernel": epanechnikov_kernel,
    "Tricube Kernel": tricube_kernel,
}
plots = []

for kernel_name, kernel_func in kernels.items():
    # 使用指定核函数计算 LLR
    y_curves = local_linear_regression(X, Y, x_curve, bandwidths, kernel_func)

    # 创建 Plotly 图形
    fig = go.Figure()

    # 添加静态轨迹
    fig.add_trace(
        go.Scatter(
            x=X, y=Y, mode="markers", name="噪声数据", marker=dict(color="gray")
        )
    )

    fig.add_trace(
        go.Scatter(
            x=x_true,
            y=y_true,
            mode="lines",
            name="真实函数",
            line=dict(color="red"),
        )
    )

    # 添加第一个 LLR 曲线（使用带宽的中间值）
    fig.add_trace(
        go.Scatter(
            x=x_curve,
            y=y_curves[round(initial_λ, 2)],
            mode="lines",
            name=f"{kernel_name} 曲线",
            line=dict(color="yellow"),
        )
    )

    # 定义滑块步骤
    steps = []
    for λ in bandwidths:
        λ_rounded = round(λ, 2)
        step = dict(
            method="update",
            args=[
                {"y": [Y, y_true, y_curves[λ_rounded]]},  # 更新轨迹的 y 数据
                {
                    "title": f"LLR: {kernel_name} 带宽 λ = {λ_rounded}"
                },  # 动态更新标题
            ],
            label=f"{λ_rounded}",
        )
        steps.append(step)

    # 将滑块添加到布局中
    sliders = [
        dict(
            active=len(bandwidths) // 2,  # 使用中间带宽的索引
            currentvalue={"prefix": "λ = "},
            pad={"t": 50},
            steps=steps,
        )
    ]

    fig.update_layout(
        sliders=sliders,
        title=f"LLR: {kernel_name} 带宽 λ = {round(initial_λ, 2)}",
        xaxis_title="X",
        yaxis_title="Y",
        template="plotly_dark",
        height=400,
        width=730,
    )

    plots.append(fig)

# 显示并保存绘图
for i, (kernel_name, fig) in enumerate(zip(kernels.keys(), plots)):
    fig.show()
    html_path = f"./llr_{kernel_name.lower().replace(' ', '_')}.html"
    fig.write_html(html_path)
    print(f"已保存 {kernel_name} 的交互式绘图到 {html_path}")

我觉得结果看起来平滑多了！

参考文献

统计学习基础 - Hastie, Tibshirani, 和 Friedman (2009). 一本关于数据挖掘、推断和预测的全面指南。了解更多.

←

因果自注意力机制的工作原理

交互式高斯混合模型

→