基于Tikhonov正则化方法并使用类Nestrov动量法的滤波降噪算法的MATLAB实现

GraviSal

于 2024-08-17 19:13:11 发布

阅读量190

点赞数 5

文章标签： matlab 算法

本文链接：https://blog.csdn.net/2301_80251865/article/details/141283588

版权

在处理时域信号时，常常会遇到数据含有噪音的情况，这在坐标图上的表现就是图像是有很多微x小的的“突刺”；对于一些精度要求偏高的情形，会影响后续信号特征的提取以及频域信号的转换。

我们的任务是尽量让图像尽可能地光滑，而它的一些阶跃点我们仍然保留。基于这样的思想，可以用Tikhonov正则化方法来实现这一任务。

对于一般的拟合任务，我们假定x为去噪信号，y为原始信号，那么我们既要希望去噪信号与原始信号之间的差异尽可能地小，同时去噪信号x还要尽可能地光滑一些。因此我们让正则化项为x的二阶导数，这里用二阶的梯度算符表示，而拟合项为x与y的差（K取单位矩阵）。这样一来，通过选取合适的正则化参数α，我们就可以兼顾去噪信号x的光滑性，和基本特征的保留程度。

$\quad J_0(x)=\|Kx-y\|^2+\alpha\|\nabla^2 x\|^2, \quad x\in X$

如此一来，就有：

$\quad x^{\alpha}=\arg\min_{x\in X}J_{\alpha}(x)$

事实上，当我们选定了合适的正则化参数以后，泛函 $J_{\alpha}(x)$ 就是已经确定了的，这样的话，我们只需要找到使得该函数最小化的x，这就是得到的去噪信号x，具体的求解 $J_{\alpha}(x)$ 的最小值一般有两种方法，一种是令 $\nabla J_{\alpha}(x)=0$ ，将 $\nabla^2 x$ 用二阶中心差分近似代替（步长的系数由正则化系数α决定），这样得到的三对角矩阵记为L，即 $\nabla J(x) = \nabla\left(\frac{1}{2}\|Ax-y\|^2\right) + \nabla\left(\frac{\alpha}{2}\|Lx\|^2\right)$ .

又因为 $\nabla\left(\frac{1}{2}\|Kx-y\|^2\right) = K^T(Kx-y)$ ， $\nabla\left(\frac{\alpha}{2}\|Lx\|^2\right) = \alpha L^TLx$ ，

于是 $\nabla J(x) = K^T(Kx-y) + \alpha L^TLx=0$ ，用最小二乘法解出这个方程中的x，就得到了去噪的信号x。实际求解过程中，使用默认的高斯消元法求解会消耗MATLAB较多的内存和算力，并且求解速度较慢于第二种做法；所以如果要使用本方法，需要配合雅可比迭代法，或者高斯-赛德尔方法快速实现矩阵运算。

而另一种方式是用梯度下降（Gradient Descent）算法，具体的原理这里略去。总而言之，使用这种算法，是不断地迭代去噪信号x，直到迭代后的x最终满足一定的条件，此时即认为得到的x为去噪信号。本文使用的matlab代码即是使用方法二，并在此基础上，使用了变种Nestrov动量法来加速梯度下降的过程。即在梯度下降的过程中，受前一步的影响越来越大，改良算法可以更快地达到收敛。
具体的Nestrov算法是这样的：

$u_k = y_{k-1} - \alpha \nabla f(y_{k-1})$ ，

$y_k = u_k + \frac{k-1}{k+2}(u_k - u_{k-1}) \quad \text{given some } u_0, u_1 \in U$
以下是Matlab代码

function [dinoised_signal,raw_signal] = Tikhonov_Nestrov(signal,eps,h,alpha)
%%INPUT

%signal       去噪信号向量
%eps          精度，当精度达到某要求时，返回所得到的去噪信号
%h            改良的梯度优化算法的步长值，可以设置，以更好地更快地实现去噪
%alpha        正则化系数，显示比重

%%OUTPUT

    if size(signal,1) == 1
        signal = signal';%转置为列向量
    end
    n = size(signal,1);
    L = zeros(n-2,n);
    K = eye(n);
    for idx = 1:n-2
        L(idx,idx) = 1;
        L(idx,idx+1) = -2;
        L(idx,idx+2) = 1;
    end
    f = @(x).5*norm(K*x-signal)+.5*alpha*norm(L*x);
    df = @(x)K'*(K*x-signal)+alpha*(L'*L)*x;
    iter = 1;
    x0 = signal;
    grad = df(x0);
    x = x0;
    U = zeros(n,2);
    while (any(grad>eps))
        U(:,1) = x - h*grad;
        x = U(:,1) + (iter-1)/(iter+2)*(U(:,1)-U(:,2));
        U(:,2) = U(:,1);
        % 上面三行是改进的Nestrov动量法
        % x = x-h*grad;%如果把上面三行代码换成这一行，就是最基本的梯度下降法
        grad = df(x);
        iter = iter+1;
        % if iter>=1000
        %     x = x0;
        %     break;
        % end
    end
    %这里的最终的x即为去噪信号，而x0是原始信号，可以用绘图的方式比较两者的差异
    dinoised_signal = x;raw_signal = x0;
end

具体的效果可以参考下面一张图，这是笔者在一次数模论文中的插图，红色是原始信号，黑色是去噪信号。选取不同的正则化系数有不同的效果，一般为0.001，0.01，0.1，不会太大

具体的滤波降噪算法有很多种，还有比较常见的就是中值滤波算法，卡尔曼滤波算法。对于不同的实际问题，选用不同的滤波算法收获的效果是完全不同的，其他的滤波算法，未来有空我也会继续分享。以上代码都是可以轻松移植并使用的，如果本文对你有帮助的话，请点一个免费的赞，您的支持是我继续更新的动力！

GraviSal

关注

5
点赞
踩
5

收藏

觉得还不错? 一键收藏
0
评论
基于Tikhonov正则化方法并使用类Nestrov动量法的滤波降噪算法的MATLAB实现

本文分享介绍使用Tikhonov正则化方法的滤波降噪算法，并使用Nestrov动量法来实现求解，对于一般时频分析的问题能起到不错的滤波降噪作用。本文的MATLAB代码可以直接使用
复制链接

扫一扫