L曲线正则化与Tikhonov参数选择：从原理到Python实现

简介：MATLAB环境下的Tikhonov正则化完整工具包，面向需要借助岭回归解决过拟合问题、处理不适定反问题的研究者和学生。包内12个文件均为.m脚本，涵盖核心算法、L曲线绘制、广义交叉验证（GCV）、奇异值分解（SVD）及最小二乘求解等关键环节，并附有phillips、shaw等经典测试问题，方便直接运行验证。作者将Tikhonov正则化的理论实现、正则化参数λ的L曲线选择策略与实际算例整合在一起，使用者可结合l_curve.m、plot_lc.m和gcv.m系统对比不同选参方式的差异。已有1330人浏览学习。整套代码包体积仅14KB，却完整覆盖了从数学模型推导到Matlab代码实现的正则化落地流程，无论用于课堂实验、论文复现还是工程快速原型验证，都能帮助读者显著缩短从原理到应用的转化时间，尤其适合机器学习、信号处理和数值计算方向的入门与进阶实践。

1. L曲线正则化：Tikhonov 正则化里最该被重视的选参方法

手里拿到 tikhonov.zip 这类打包好的 Tikhonov 正则化代码包，不等于知道 λ 怎么定。最常见的情况是：方程 Ax=b 病态到最小二乘解范数爆炸、符号乱跳，加 Tikhonov 正则化项把解稳住，但正则化系数 λ 选大了解被压成一条直线，选小了等于没加。L曲线正则化就是专门解决这个「选 λ」问题的——把不同 λ 下的残差范数和解范数画在双对数坐标里，得到一条 L 形曲线，拐点对应的就是推荐 λ。这篇笔记面向做反演、图像恢复、数值微分这类病态逆问题的人，从原理讲到可复现的最小实现，再给你五个真实的翻车现场，读完可以直接照着调自己的数据。

2. Tikhonov 正则化先立住：病态问题、解的唯一性与正则化系数

L 曲线不是独立工具，它是给 Tikhonov 正则化挑 λ 的配套方法。先把 Tikhonov 为什么存在、λ 在干什么讲清楚，后面看曲线才不是一头雾水。这一章不涉及复杂推导，只讲「什么时候必须用」「λ 大了小了各会发生什么」「它和岭回归、弹性网正则化到底有什么区别」。

2.1 病态问题：满秩矩阵为什么给不出可信解

第一类积分方程离散化、图像去模糊里的点扩散函数、数值微分的差分格式，最后都会落到一个线性方程组 Ax=b 上。这类问题的共同特征是矩阵 A 的条件数极大。拿 Hilbert 矩阵举例：H_{ij}=1/(i+j-1)，到 8 阶时条件数已经超过 1e10，矩阵满秩、可逆、求逆毫无问题，但你往右端项里加一个 1e-6 量级的扰动，解就会完全变样。

它的本质是 A 的奇异值衰减太快。最小二乘解 x_ls=(A^T A)^{-1}A^T b 在代数值上完全合法，但高频方向上的小奇异值会把噪声放大几个数量级。于是你得到一个残差很小、范数巨大、元素符号随机翻转的解——从数值上看它精确拟合了数据，从物理上看它不是你要的那个场。这就是病态问题最反直觉的地方：残差小不等于解可信。

数值微分是这类问题里最常见的例子。观测数据带噪声，你对它做差分，差分算子作用在噪声上会产生高频振荡。把差分格式写成矩阵，它的奇异值从大到小衰减得比指数还快，任何一点噪声都会主导解的形状。此时直接求解没有意义，必须给解加约束。Tikhonov 正则化干的正是这件事：在「拟合数据」和「控制解的形状」之间加一个折中项。

2.2 正则化系数 λ 的三种效应：欠正则化、过正则化与临界点

Tikhonov 正则化的目标函数写成：

J(x)=||Ax-b||²+λ²||Lx||²

其中 λ 是正则化系数，L 是正则化矩阵。λ=0 时退化为最小二乘；λ→∞ 时，解会被压向 L 的零空间——L=I 时就是压向零向量。实际计算中 λ 落在三个区间，对应三种完全不同的结果：

λ 太小，惩罚项形同虚设，小奇异值仍然主导解，||Lx|| 巨大，残差却逼近噪声下限，这是过拟合段。λ 太大，解的范数被压得太狠，残差远大于噪声本身能解释的范围，这是欠拟合段。只有中间某个临界 λ，让残差和解范数达到平衡，解既拟合数据又保持稳定。

正则化系数 λ 是整个 Tikhonov 流程里唯一需要手调的全局参数。它不像矩阵分解里那种「算一次就完事」的参数，而是直接决定解是噪声主导还是过度平滑。调 λ 曾经基本靠经验试错，在不同数量级之间来回扫，非常玄学。L 曲线法就是把这个过程变成一套可视化的几何判断，让你不用盲试。后面所有实现都围绕 λ 展开，先记住它的目标函数形式，后面代码直接对应。

2.3 Tikhonov 与岭回归、弹性网正则化的关系：L 矩阵的引入

岭回归是 Tikhonov 在 L=I 时的特例，目标函数是 min||Ax-b||²+λ²||x||²。它在统计里很常用，作用是压缩线性回归的系数。区别在于岭回归里 A 是设计矩阵，特征已经标准化，λ 让系数整体收缩；而 Tikhonov 里 A 是离散算子，L 可以换成差分矩阵，作用不只是压缩，还包括平滑和形状约束。

弹性网正则化是 L1+L2 的组合，面向高维变量选择，它解决的变量稀疏性问题和 Tikhonov 解决的病态逆问题完全是两件事。如果项目里有人把弹性网正则化直接套到反演问题里，通常会出两个毛病：L1 项强制一部分系数归零，物理上连续的场被切成一块一块的；L2 项的系数又没法像 Tikhonov 的 L 矩阵那样表达先验结构。所以别看到「正则化」三个字就互相替代，目标不同，工具不能乱换。

L 矩阵的选择是 Tikhonov 真正的灵活之处。零阶：L=I，约束解的能量，适合解本身没有平滑先验的情况。一阶：L 取一阶差分矩阵，惩罚相邻元素之差，让解平滑。二阶：L 取二阶差分矩阵，惩罚曲率，让解的斜率变化平稳。L 的维度随 A 列数变化，一阶差分 L∈R^{(n-1)×n}。说白了，L 矩阵里装的是物理先验，这也是 L 曲线纵轴「解范数」到底在衡量什么的关键——它衡量的不是原始解的大小，而是解在 L 定义下的形状代价。

3. L 曲线方法的原理：把选 λ 变成找拐点的几何问题

L 曲线这名字听起来复杂，其实核心就是把二维曲线的一个点拉成一条曲线。这一章讲清楚为什么残差范数和解范数画出来像「L」，拐点为什么是最优 λ，以及数值上怎么把拐点算出来。理解这三个问题，第 4 章的代码就只是翻译。

3.1 L 曲线为什么是 L 形：残差范数与解范数的对抗

对每个候选 λ，用 Tikhonov 求出一个解 x_λ，然后算两个标量：残差范数 ρ(λ)=||Ax_λ-b||，解范数 η(λ)=||Lx_λ||。把 log ρ 做横轴、log η 做纵轴，把所有 λ 对应的点连起来，就得到一条 L 形曲线。

这条曲线为什么是 L 形？回到 λ 的三种效应：λ 极小时，解过拟合，残差趋近于零，解范数巨大，所以点落在曲线左边靠上的竖直段；λ 极大时，解被严重压缩，解范数趋近于零，残差巨大，点落在曲线底部靠右的水平段；中间过渡段把这两段连起来，形成一个明显的拐角。

拐点为什么是好的选择？看几何意义。在拐点左侧，你往 λ 小的方向移动，解范数急剧上升，但残差几乎不变——你在用很大的「解代价」换取极小的「拟合改善」。在拐点右侧，你往 λ 大的方向移动，残差急剧上升，但解范数几乎不变——你在用很大的「拟合恶化」换取极小的「解改善」。拐点处这两个方向的边际变化率达到平衡，相当于折中。这就是 L 曲线法全部直觉所在：最优 λ 让残差和解范数都不想再为对方让步。

3.2 拐点定位：数值曲率计算与重参数化

有了 L 形曲线，下一步是让计算机找拐点。拐点在几何上对应曲率最大的点。平面曲线 (x(t), y(t)) 的曲率公式是：

κ = (x'y'' − y'x'') / (x'² + y'²)^{3/2}

其中 x=log ρ，y=log η，t 是参数。理论上把每个 λ 代入就能算出曲率序列，取最大值对应的 λ* 即可。但实际操作有个坑：λ 通常横跨 6~8 个数量级，如果直接用 λ 作为参数 t，采样点在曲线两端密集、中间稀疏，一阶导和二阶导的差分会严重失真。我一般会先用 log λ 作为参数，再对 log ρ 和 log η 做三次样条拟合，在密集的 log λ 网格上重新采样，然后用样条的一阶、二阶导数值代曲率公式。三次样条保证二阶导数连续，曲率曲线不会出现人为锯齿。如果你用线性插值去做这件事，二阶导恒为零，曲率公式直接失效，这是新手最容易踩的坑。

一个值得注意的细节是：曲线的拐点在 log-log 空间里定义，而不是在线性空间里定义。因为 ρ 和 η 经常跨多个数量级，线性坐标系下曲线会被挤成一条直角的折线，视觉上像 L，但曲率数值不稳定。取 log 之后两个轴都是对数尺度，曲线在各段上的几何特征才均匀，最大曲率点才是稳定的 λ*。

3.3 L 曲线的边界：什么时候拐点会失效

L 曲线法不是万能的，有几种情况拐点会失效，必须在用之前判断。

第一种是问题本身良态或噪声极低。此时曲线没有明显的 L 形，最大曲率的位置对 λ 网格的选择极其敏感，λ* 会在网格端点之间跳来跳去。这种情况说明最小二乘已经够用，Tikhonov 是多余的正则化。

第二种是 ρ 和 η 的尺度相差几个数量级。比如残差范数在 1e-2 量级，解范数在 1e10 量级，画出来是一条近乎竖直的线，拐点肉眼都找不到。解决办法是对 A、b、L 做归一化，或者用广义奇异值分解把两个范数放到可比尺度上，后面代码里会给出具体做法。

第三种是 λ 采样过密且噪声主导。样条会拟合出局部毛刺，产生伪拐点，导致 λ* 落在完全错误的区域。常见做法是把采样点控制在 100~200 个，太多反而坏事。这三个边界直接决定第 4 章代码里 λ 网格的范围、采样密度和数据预处理方式，属于 L 曲线的使用说明书。

4. 用 Python 复现 tikhonov.zip 的核心：L 曲线选 λ 的最小实现

标题里的 tikhonov.zip，看名字就是一个把 Tikhonov 正则化与 L 曲线选参封装好的工具包。这类包不管界面长什么样，核心逻辑都是同一套：生成 λ 网格、对每个 λ 求解、计算残差范数与解范数、定位 L 曲线拐点。下面我用 Python 从零把这条链路实现一遍。跑通之后，你拿到任何封装包都能快速验证它的 λ 选得对不对，而不是把它当黑匣子。

4.1 生成病态测试问题：从 Hilbert 矩阵开始

写代码之前先造一个可控的病态问题。Hilbert 矩阵是经典选择，不需要外部数据，一段代码就能生成，条件数随阶数急剧上升，非常适合做基准测试。

import numpy as np def hilbert_matrix(n): """生成 n 阶 Hilbert 矩阵，条件数随 n 指数增长""" i, j = np.indices((n, n)) return 1.0 / (i + j + 1) n = 8 A = hilbert_matrix(n) cond = np.linalg.cond(A) print(f"A 的条件数: {cond:.3e}") # 构造一个已知的"真实解"，让后续能对比还原效果 x_true = np.sin(np.linspace(1.0, 3.0, n)) b0 = A @ x_true # 加入噪声，模拟观测误差；噪声放大会让病态问题更加明显 rng = np.random.default_rng(42) noise_std = 1e-6 b = b0 + noise_std * rng.standard_normal(n)

这段代码的逻辑是：先建立 A 和已知解 x_true，算出无噪声的右端项 b0，再加一个微小噪声得到实际观测 b。为什么要有 x_true？因为后续验证 L 曲线选出的 λ 到底好不好时，需要和真实解比较。n 选 8 是因为 8 阶 Hilbert 矩阵的条件数足够高，能明显看出最小二乘解爆炸，又不至于让数值计算完全失效。噪声标准差 1e-6 看似很小，但在这个条件数下已经足以把解搅乱，这正是病态问题的典型特征。

4.2 实现 Tikhonov 求解与 L 曲线曲率定位

有了测试问题，接下来是两步核心函数：Tikhonov 求解器和 L 曲线拐点定位。求解器用正规方程实现，目标函数和第 2 章保持一致：||Ax-b||²+λ²||Lx||²。

def tikhonov_solve(A, b, lam, L=None): """求解 Tikhonov 正则化问题 目标函数: ||Ax - b||^2 + lam^2 * ||Lx||^2 法方程: (A^T A + lam^2 * L^T L) x = A^T b """ A = np.asarray(A, dtype=float) b = np.asarray(b, dtype=float) if L is None: L = np.eye(A.shape[1]) lhs = A.T @ A + (lam ** 2) * (L.T @ L) rhs = A.T @ b return np.linalg.solve(lhs, rhs)

正规方程的写法最简单直白，但有个已知缺陷：A^T A 会把条件数平方。所以在 lam 很小时数值上可能有不稳定风险，我会在第 4.4 节给出更稳的 SVD 替代方案。这里先用正规方程，因为逻辑清楚，适合演示主流程。注意目标函数里 lam 的幂次是 2，对应的法方程修正项是 lam²L^T L，这个细节直接决定后面 λ 网格的物理含义，不要和别处用 lam 一次方的写法混用。

然后是 L 曲线的核心部分：对一组 λ 分别求解放置曲线，再用三次样条求最大曲率点。

from scipy.interpolate import CubicSpline def l_curve_corner(A, b, lams, L=None): """计算 L 曲线并返回最大曲率点对应的 lambda""" rho, eta = [], [] for lam in lams: x = tikhonov_solve(A, b, lam, L) rho.append(np.linalg.norm(A @ x - b)) eta.append(np.linalg.norm(L @ x)) # 取对数，进入 log-log 空间 X = np.log(np.array(rho)) Y = np.log(np.array(eta)) # 以 log(lambda) 为自变量，用三次样条拟合两条曲线 T = np.log(np.array(lams)) spline_x = CubicSpline(T, X) spline_y = CubicSpline(T, Y) # 在密集 log-lambda 网格上计算一阶、二阶导数和曲率 T_dense = np.linspace(T.min(), T.max(), 1000) x1 = spline_x(T_dense, 1) x2 = spline_x(T_dense, 2) y1 = spline_y(T_dense, 1) y2 = spline_y(T_dense, 2) kappa = (x1 * y2 - y1 * x2) / np.power(x1**2 + y1**2, 1.5) # 曲率最大值对应的 lambda 就是拐点 lam_star = np.exp(T_dense[np.argmax(kappa)]) return lam_star, np.array(rho), np.array(eta)

逻辑说明：对每个 λ 解一次线性方程组，收集残差范数 ρ 和解范数 η，然后全部取对数。这里选择 log λ 作为样条自变量，而不是 λ 本身，因为 λ 网格通常横跨多个数量级，直接用 λ 会让曲线在两端严重压缩，导数不稳定。CubicSpline 保证二阶导数连续，曲率曲线是光滑的，最大点就是 L 曲线的拐角。返回的 lam_star 就是 L 曲线法推荐的 λ，rho 和 eta 留着画图用。

4.3 主流程：把 λ 网格、求解、拐点串起来

两个核心函数写完后，主流程就只剩生成 λ 网格并调用。λ 网格用对数等距，覆盖范围要足够大。

# 对数等距的 lambda 网格，覆盖多个数量级 lams = np.logspace(-12, 2, 200) L = np.eye(n) # 先用零阶 Tikhonov（岭回归），后续可换差分矩阵 lam_star, rho, eta = l_curve_corner(A, b, lams, L) # 计算最小二乘解作为对照 x_ls = np.linalg.solve(A.T @ A, A.T @ b) x_reg = tikhonov_solve(A, b, lam_star, L) print(f"L 曲线推荐 lambda: {lam_star:.3e}") print(f"最小二乘解范数: {np.linalg.norm(x_ls):.3e}") print(f"Tikhonov 解范数: {np.linalg.norm(x_reg):.3e}") print(f"与真实解误差: {np.linalg.norm(x_reg - x_true):.3e}")

这一段的重点是先打印出 λ*，再对比最小二乘和 Tikhonov 的解范数。你会发现最小二乘解范数可能是 Tikhonov 的几千倍，这就是病态问题不加约束的后果。λ* 落在网格两端时说明范围不对，要按下一节的动态范围规则重新设置，而不是手动硬凑。另外可以把 x_reg 和 x_true 画在一起看形状，Tikhonov 解应该基本还原真实解的趋势，而不是一条高频振荡的折线。

4.4 关键参数与边界条件：λ 网格范围、采样密度、正则化矩阵

L 曲线选 λ 的结果对三个参数敏感，用的时候按下面的原则设。

注意：用正规方程实现时，λ 下界不能取得太小，因为 A^T A 会把条件数从 1e10 变成 1e20，double 精度下已经不可靠。

更稳的做法是用 SVD 替代正规方程。当 L=I 时，Tikhonov 解有解析形式：对 A 做奇异值分解 A=UΣV^T，则 x_λ=Σ(σ_i²/(σ_i²+λ²))(u_i^T b/σ_i)v_i。这个实现不需要构造 A^T A，数值稳定性好得多，λ 网格下界可以进一步缩小。如果 L 不是单位阵，需要用到广义奇异值分解（GSVD），Python 标准库没有现成的 gsvd，可以用scipy.linalg里的qz或自己封装 QR 分解来做，工程上我会先把 L=I 的 SVD 版本跑通，再用 GSVD 处理带差分矩阵的情况。

5. 避坑指南：L 曲线正则化最常见的五个翻车现场

L 曲线法原理不复杂，但落地时我见过也踩过不少坑。这一章挑五个高频问题，按「现象 → 原因 → 解决」写清楚，每条都是可以直接对着排查的清单。

5.1 λ 与 λ² 混用，同一个网格却得到完全不同的拐点

现象：按网上某篇笔记的公式写了目标函数，却用另一篇代码的 λ 网格，结果 L 曲线的拐点位置诡异，λ* 和之前调好的经验值差了 100 倍。

原因：不同实现里目标函数写法不统一。有的写 ||Ax-b||²+λ²||Lx||²，有的写 ||Ax-b||²+λ||Lx||²。用第二种写法时，实际正则化强度是 sqrt(λ)，同一个 λ 数值下解的形状完全不同，L 曲线的横纵坐标对不上，拐点当然偏移。

解决：动手前先确定自己代码里 λ 的幂次。我一般统一用 λ² 形式，因为 SVD 解析解里出现的就是 σ_i²+λ²，这个写法最自然。如果接手别人的代码，第一步先看目标函数里 λ 是从一次方还是二次方进入的，再决定网格范围。

5.2 拐点卡在 λ 网格的两端

现象：运行曲率计算后，λ* 等于 lams[0] 或 lams[-1]，或者非常靠近端点。这时候 L 曲线法给出的不是「最优 λ」，而是「网格边界」。

原因：λ 网格范围没覆盖曲线真正的拐角。常见有几种：网格下界太大，曲线左侧竖直段根本没画出来；网格上界太小，水平段也没画出来；再或者数据没归一化，ρ 和 η 尺度悬殊，拐点被挤到角落。

解决：先用 4.4 节的规则把网格扩到 σ_min 的百分之一到 σ_max 的十倍，再检查曲线形态。如果拐点还是贴边，用 SVD 实现正则化而不是正规方程，把下界再降几个数量级。还有一种快速目检法：打印出 rho 和 eta 的最大最小值，如果横轴或纵轴只变化了一个数量级，说明问题良态，L 曲线法本身就不适用。

5.3 零阶 Tikhonov 处理平滑性问题，解被均匀压缩

现象：数据是地震道或图像行，需要平滑的剖面，用 L=I 做 Tikhonov，L 曲线选出的 λ 让解整体变小，但高频毛刺还在。

原因：零阶 Tikhonov（L=I）惩罚的是解向量本身的大小，它把大值和小值一起往零压，并不惩罚相邻点的剧烈变化。对平滑性起作用的是一阶差分或二阶差分矩阵，它们惩罚的是相邻差异和曲率。L=I 适合压缩幅值，不适合平滑结构。

解决：把 L 换成一阶差分矩阵 D1，目标函数变成 ||Ax-b||²+λ²||D1x||²。这时 η=||D1x|| 衡量的是解的总变差，L 曲线拐点选出的 λ 会让解既拟合数据又保持分段平滑。如果数据有断裂或台阶，用二阶差分 D2 惩罚曲率即可。改 L 之后，L 曲线的纵轴含义变了，λ* 数值不能和零阶版本直接对比。

5.4 GCV 与 L 曲线给出矛盾结论

现象：同一个问题，广义交叉验证（GCV）选出的 λ 比 L 曲线小两三个数量级，两个结果差距太大，不知道该信谁。

原因：两者对噪声的假设不同。GCV 假设残差是独立同分布的高斯噪声，在噪声有相关性或非高斯时严重偏向小 λ；L 曲线法对噪声分布的依赖相对弱，但在曲线没有明显拐角时也不稳定。工程上这俩给出同一数量级的 λ 是正常情况，差一两个数量级就需要检查数据。

解决：我用 L 曲线做候选范围，再用 GCV 在候选范围内做二次筛选。具体操作是：先取 L 曲线拐点的 λ*，在 λ*/10 到 λ*×10 之间做 GCV 曲线，取 GCV 最小值为最终 λ。这样既保留 L 曲线对相关噪声的鲁棒性，又用 GCV 做了精修。千万别只拿其中一个当真理，在反演项目里两个方法打架是常态。

5.5 把一致性正则化机制与 Tikhonov 正则化当成一回事

现象：项目文档里同时出现了 Tikhonov 正则化和一致性正则化机制，有人把半监督学习里的正则化系数直接挪到反演问题上，结果发现目标函数和收敛行为完全对不上。

原因：这俩名字都带「正则化」，但根本不是一个层面的东西。Tikhonov 正则化是逆问题里的解先验约束，作用于解向量的范数或平滑度；一致性正则化机制是半监督学习里的扰动不变性约束，要求模型对输入扰动给出一致的预测。一个约束的是解的形状，一个约束的是模型的行为，作用对象、数学形式、调参方法全都不同。

解决：拿到一个正则化相关的工具包，先看它作用在谁身上。如果约束项里出现的是 ||Lx|| 这类解向量的函数，大概率是 Tikhonov 体系；如果约束项是模型对两个不同增强视图输出的距离，那属于一致性正则化，不要混用。尤其注意机器学习框架里的正则化系数和反演问题里的 λ，单位、含义、搜索方式都不一样，照搬会翻车。

6. 进阶验证：用模拟数据检验 L 曲线选出的 λ 是否可信

L 曲线给出 λ* 后，怎么知道这个推荐值到底有多好？在真实数据里你没法验证，但模拟数据能：构造一个已知真实解的问题，穷举所有 λ 找出真正的最优解，再和 L 曲线的推荐值对比。这个验证流程我每次换数据形态都会跑一遍，十几分钟就能定位实现里的问题。

6.1 构造已知真实解的模拟实验并穷举最优 λ

基于第 4 章的 A 和 x_true，遍历 λ 网格，对每个 λ 计算重建解与真实解的误差，误差最小的 λ 就是理论最优。

def best_lambda_by_error(A, b, x_true, lams, L=None): """穷举 lambda，返回重建误差最小的最优 lambda""" errors = [] for lam in lams: x = tikhonov_solve(A, b, lam, L) errors.append(np.linalg.norm(x - x_true)) return lams[np.argmin(errors)], errors best_lam, errors = best_lambda_by_error(A, b, x_true, lams, L) print(f"理论最优 lambda: {best_lam:.3e}") print(f"L 曲线推荐 lambda: {lam_star:.3e}")

误差序列关于 λ 通常呈 U 形：λ 太小重建误差被噪声主导，λ 太大重建误差被过度平滑主导，中间最低点就是理论最优。如果 L 曲线选出的 λ* 和理论最优在同一个数量级内，说明实现没问题，可以放心用；如果差了两个数量级以上，多半是 λ 网格范围没覆盖拐点，或者 L 矩阵选得不对，回去改。

6.2 用误差比作为验收指标

我一般还会算一个定量指标：把 L 曲线 λ* 对应的误差除以理论最优误差，得到一个误差比。误差比在 1.5 以内说明 L 曲线选参已经接近完美；在 3 以内说明可用，但如果你追求极致可以再用 GCV 精修；超过 3 就要回头排查数据归一化和网格范围。这套验证的成本只是几十次小规模线性方程求解，而收获是确认整个选参链路没有系统性偏差，以后换到真实数据，你至少知道你的 λ 不是瞎猜的。

说句血泪经验：我最早做数值微分恢复时，图省事固定了一个 λ 用到底，换了一组噪声更大的数据后解直接变成高频振荡，项目差点推翻重做。后来所有 Tikhonov 求解一律先用 L 曲线开路，再用模拟数据验证，再也不用拍脑袋定正则化系数。亲眼看着 L 曲线拐点和理论最优落进同一个数量级，比任何理论推导都让人安心。希望帮到你。

本文还有配套的精品资源，点击获取