DOI: 10.12677/aam.2024.1312524, PDF, HTML, XML,  被引量    国家自然科学基金支持
作者: 江春雨：长春理工大学数学与统计学院，吉林 长春；贾小宁：长春理工大学数学与统计学院，吉林 长春；长春理工大学中山研究院遥感技术与大数据分析实验室，广东 中山
关键词: 广义高斯混合分布；低秩矩阵分解；EM算法；图像去噪；Generalized Gaussian Mixed Distribution； Low-Rank Matrix Factorization； EM Algorithm； Image Denoising

文章引用：江春雨, 贾小宁. 基于广义高斯混合模型的图像去噪算法研究[J]. 应用数学进展, 2024, 13(12): 5428-5438. https://doi.org/10.12677/aam.2024.1312524

1. 引言

在图像的拍摄、采集、传输等过程中，不可避免地会引入不必要的信号，这些信号会影响图像的质量，并且在后续图像处理中产生相应的干扰，这些干扰信号称为噪声[1]。图像去噪就是除去图像中的不必要的信号，提高图像的质量，更好地进行图像处理操作。

目前，传统的图像去噪方法[2]包括均值滤波、中值滤波、维纳滤波[3]等，这些方法在去除噪声的同时往往会引入一定的模糊，导致图像细节的损失。为了克服这一问题，研究者们提出了许多先进的图像去噪方法，包括基于稀疏表示、低秩分解[4]、图像先验知识[5] [6]等的方法。

图像去噪在图像处理中占有重要地位，许多去噪方法依赖于对噪声和图像信号的精确建模。高斯分布在建模噪声方面具有重要作用。Meng和La Torre开发了基于MoG的LRMF模型[7] (MoG-LRMF)。MoG的普遍逼近特性使得该模型能够自适应地拟合具有未知统计分布的噪声，这一思想也被应用于低秩张量分解[8]。然而，只有当高斯分量的数量达到无穷大时，MoG才能近似复杂的分布(例如拉普拉斯分布)，因此在实践中是不可行的。针对这个问题，Cao等人提出了一种指数功率分布的混合(PMoEP) [9]来对未知噪声进行建模，并设计了一种惩罚方法来自适应地选择混合分量数。MoEP中的成分非常复杂，可以在一定程度上缓解MoG的不足。高斯混合模型克服了单一高斯分布无法准确描述图像等复杂数据信息的难题，2018年，Chen等人首先观察到不同波段图像的噪声分布在统计上是不同的，但具有一定的相关性，然后提出了非i.i.d. MoG (NMoG)模型[10]来描述这一现象。与高斯模型相比，非对称高斯模型强调了图像的非对称性特征，具有更强的边缘保持能力。

广义高斯混合模型能够适应多种形状的噪声分布，提供比单一高斯模型和传统高斯混合模型更精确的噪声建模。1999年，Molina等人提出了将广义高斯分布用于图像去噪的应用，通过最大后验估计(MAP)方法进行噪声建模和去除[11]。2012年，Yuan等人提出了一种基于广义高斯混合模型的图像去噪方法[12]，通过期望最大化(EM)算法优化模型参数，实现了高效地去噪。接下来，还有将广义高斯混合模型与矩阵补全技术结合，提出了一种新的图像去噪方法，通过对图像矩阵的低秩补全和广义高斯建模，实现了高效的噪声去除。牛津大学采用变分推断方法结合广义高斯混合模型，提出了一种新的图像去噪方法，通过变分推断优化广义高斯模型的参数[13]，有效去除了图像噪声。多伦多大学将贝叶斯推断方法与广义高斯混合模型结合，提出了一种基于贝叶斯框架的图像去噪算法。通过构建复杂的先验模型和后验推断[14]，取得了高质量的去噪效果。2020年至今，出现了许多将GGMM与深度学习等前沿技术结合的方法，并取得了显著的成果。

针对上述问题，本文假设噪声服从广义高斯混合分布，建立基于广义高斯混合模型的低秩矩阵分解模型。模型的参数采用EM算法进行迭代更新，在合成数据和真实图像数据上的大量实验证明了算法的有效性。

2. 相关工作

2.1. 广义高斯混合分布

2.1.1. 广义高斯分布

广义高斯分布在统计学中是一种参数化的连续概率分布，是高斯分布的扩展形式。广义高斯分布是一种灵活的概率模型，通过调整形状参数可以改变分布曲线的形状，从而描述出多种数据特性，如从尖峰到平坦，从轻尾到重尾。

广义高斯分布的模型参数有三个，均值 $\mu$ ，方差 $\sigma$ 和形状参数 $\lambda$ ，它的概率密度函数通常表示为：

$\text{GGD}\left(x\right)=\left[\frac{{\lambda }_n\cdot \eta \left({\lambda }_n,{\sigma }_n\right)}{2\Gamma \left(1/{\lambda }_n\right)}\right]\exp \left\{-\left[\eta \left({\lambda }_n,{\sigma }_n\right)\cdot {\left|x-{\mu }_n\right|}^{{\lambda }_n}\right]\right\}$ (1)

其中， $\eta \left({\lambda }_n,{\sigma }_n\right)=\frac{1}{{\sigma }_n}{\left[\frac{\Gamma \left(3/{\lambda }_n\right)}{\Gamma \left(1/{\lambda }_n\right)}\right]}^{\frac{1}{2}}$ 。

2.1.2. 广义高斯混合分布

广义高斯混合模型是基于广义高斯分布的概率模型，通过将多个广义高斯分布加权组合，能够更灵活地拟合复杂数据的分布。广义高斯混合模型的分布为：

$p\left(x\right)={\displaystyle \underset{n=1}{\overset{N}{\sum }}{\pi }_n\text{GGD}\left({\mu }_n,{\sigma }_n,{\lambda }_n\right)}$ (2)

广义高斯混合模型(GGMM)相比单一广义高斯分布(GGD)，拥有更强的灵活性。GGMM通过多个广义高斯分布的加权组合，能够建模多模态分布，捕捉数据中不同区域的多样化特性。相比之下，单一GGD仅适用于单峰分布，无法刻画复杂分布的局部变化和多模态特性。GGMM在数据拟合、精度提升和多样性描述方面显著优于单一GGD，是复杂分布建模中的更优选择。

2.2. 低秩矩阵分解

低秩矩阵分解[15]是一种将矩阵分解为几个较小矩阵的技术，这些较小矩阵的乘积近似原矩阵，并且具有较低的秩。低秩矩阵分解在数据分析、图像处理和机器学习等领域有广泛应用，尤其在图像去噪中表现出色。常见的低秩矩阵分解方法包括：奇异值分解[16]，主成分分析[17]，非负矩阵分解[18]。

给定一个矩阵 $X\in {ℝ}^{m\times n}$ ，低秩矩阵分解是将它近似表示为两个低维矩阵 $U\in {ℝ}^{m\times r}$ 、 $V\in {ℝ}^{n\times r}$ 的乘积：

$X\approx U{V}^{\text{T}}$ (3)

其中， $r\ll \min \left(m,n\right)$ 是分解矩阵的秩。在低秩矩阵分解中，分解矩阵的秩r决定了所保留的特征信息量。通过选择适当的r，可以在压缩数据维度的同时尽量保持原始矩阵的主要结构特征。

3. 算法模型

3.1. 噪声建模与目标函数

首先，利用GGMM对输入矩阵的噪声元素进行建模，从而获得对数似然优化目标;然后通过假设一个高维的潜变量，在EM框架下迭代求解。

考虑噪声部分(记为 ${\epsilon }_{ij}$ )，每个元素 $x_{ij}$ $\left(i=1,2,\cdots ,I;j=1,2,\cdots ,J\right)$ 在低秩矩阵分解中的表达式可以写成：

$x_{ij}=u_i{v}_j^{\text{T}}+{\epsilon }_{ij}$ (4)

用GGMM对原始数据中的未知噪声进行建模。因此，每个 ${\epsilon }_{ij}$ 服从一个GGMM分布 $p\left(\epsilon \right)$ 定义为：

(5)

则式(5)中的每个 $x_{ij}$ 服从GGMM分布，其均值为 ${\Lambda }_n=u_i{v}_j^{\text{T}}+{\mu }_n$ ，标准差为 ${\sigma }_n$ ，形状参数为 ${\lambda }_n$ 。

因此，输入矩阵 $\chi$ 中每个元素 $x_{ij}$ 的概率可表示为：

$p\left(x_{ij}|\text{Π},\text{Λ},\sigma ,\lambda \right)={\displaystyle \underset{n=1}{\overset{N}{\sum }}{\pi }_n\text{GGD}\left(x_{ij}|{\text{Λ}}_n,{\sigma }_n,{\lambda }_n\right)}$ (6)

然后定义输入矩阵 $\chi$ 的可能性：

(7)

目标是最大化关于参数 $\Pi$ ， $\Lambda$ ， $\sigma$ ， $\lambda$ 的对数似然函数，即：

$\underset{\Lambda ,\Pi ,\sigma ,\lambda }{\max }L\left(\theta \right)={\displaystyle \underset{ij}{\sum }\log {\displaystyle \underset{n=1}{\overset{N}{\sum }}{\pi }_n\text{GGD}\left(x_{ij}|{\Lambda }_n,{\sigma }_n,{\lambda }_n\right)}}$ (8)

EM算法是求解对数似然函数最大化问题的有效算法。因此，在求解式(8)时，我们假设在EM框架下存在一个高维的潜变量。那么原问题可以看作是均值为0的高斯尺度混合(GSM)。定义 $U=\left[u_1,u_2,\cdots ,u_r\right]$ ， $V=\left[v_1,v_2,\cdots ,v_r\right]$ 则公式(8)可以改写为：

${\max }_{U,V,\Pi ,\sigma ,\lambda }L\left(\theta \right)={\displaystyle \underset{ij}{\sum }\log {\displaystyle \underset{n=1}{\overset{N}{\sum }}{\pi }_n\text{GGD}\left(x_{ij}|u_i{v}_j^{\text{T}},{\sigma }_n,{\lambda }_n\right)}}$ (9)

为了最大化误差函数 $L\left(\Theta \right)$ ，我们定义了一个变量 $z_{ijn}$ 。变量 $z_{ijn}$ 定义为：

$z_{ijn}=\frac{{\pi }_n\text{GGD}\left(x_{ij}|u_i{v}_j^{\text{T}},{\sigma }_n,{\lambda }_n\right)}{{\displaystyle \underset{n=1}{\overset{N}{\sum }}{\pi }_n\text{GGD}\left(x_{ij}|u_i{v}_j^{\text{T}},{\sigma }_n,{\lambda }_n\right)}}$ (10)

M步使E步给出的关于每个参数的上界最大化，得到新的误差函数：

$E\left(\Theta \right)={\displaystyle \underset{i,j}{\sum }{\displaystyle \underset{n=1}{\overset{N}{\sum }}{z}_{ijn}\left[\log {\pi }_n+\log \left(\frac{{\lambda }_n\eta \left({\lambda }_n,{\sigma }_n\right)}{2\Gamma \left(1/{\lambda }_n\right)}\right)-{\left[\eta \left({\lambda }_n,{\sigma }_n\right)\left|x_{ij}-u_i{\nu }_j^{\text{T}}\right|\right]}^{{\lambda }_n}\right]}}$ (11)

3.2. 算法求解

M步更新U、V：仅对未知分量U、V重写公式(11)，如下所示：

${\displaystyle \underset{ij}{\sum }{\displaystyle \underset{n=1}{\overset{N}{\sum }}\eta \left({\lambda }_n,{\sigma }_n\right)z_{ijn}{\left|x_{ij}-u_i{v}_j^{\text{T}}\right|}^{{\lambda }_n}}}={\displaystyle \underset{n=1}{\overset{N}{\sum }}{\Vert W\odot \left(X-U{V}^{\text{T}}\right)\Vert }_{{\lambda }_n}^{{\lambda }_n}}$ (12)

这里 $\odot$ 表示Hadamard乘积(组件式乘法)，W的元素 $w_{ij}$ 为

$w_{ij}={\left(\eta \left({\lambda }_n,{\sigma }_n\right)z_{ijn}\right)}^{\frac{1}{{\lambda }_n}}$ (13)

3.2.1. 先验概率更新

构建惩罚PGGMM模型：

${\max }_{\theta }\left\{E_p\left(\theta \right)=E\left(\theta \right)-P\left(\pi ;\delta \right)\right\}$ (14)

其中， $P\left(\pi ;\delta \right)=J\delta {\displaystyle {\sum }_{n=1}^N{D}_n\log \frac{ϵ+{\pi }_n}{ϵ}}$ ，J为数据矩阵的列数， $ϵ$ 是一个非常小的正数， $\delta$ 为调优参数( $\delta >0$ )， $D_n$ 是第n个分量自由参数的数量。

构造Q函数：

$Q\left(\Theta \right)={\displaystyle \underset{i,j}{\sum }{\displaystyle \underset{n=1}{\overset{N}{\sum }}{z}_{ijn}}}\left[\log {\pi }_n+\log \left(\frac{{\lambda }_n\eta \left({\lambda }_n,{\sigma }_n\right)}{2\Gamma \left(1/{\lambda }_n\right)}\right)-{\left[\eta \left({\lambda }_n,{\sigma }_n\right)\left|x_{ij}-u_i{\nu }_j^{\text{T}}\right|\right]}^{{\lambda }_n}\right]-J\delta {\displaystyle \underset{n=1}{\overset{N}{\sum }}{D}_n}\log \left(ϵ+{\pi }_n\right)$ (15)

本文通过最大化Q函数来更新 ${\pi }_n$ ，对Q求一阶偏导数，并令导数为0，即可得到更新方程，

${\pi }_n=\frac{\left|IJ\right|\delta D_n-{\displaystyle \underset{ij}{\sum }{z}_{ijn}}}{\left|IJ\right|\left(\delta \hat{D}-1\right)}=\frac{1}{1-\delta \hat{D}}\left[\frac{{\displaystyle \underset{i,j}{\sum }{z}_{ijn}}}{\left|IJ\right|}-\delta D_n\right]$ (16)

其中 $\hat{D}={\displaystyle \underset{n=1}{\overset{N}{\sum }}{D}_n}$ 。

3.2.2. 标准差更新

在更新参数前，先简单介绍一下稳健统计理论。稳健统计理论是一门研究统计估计在存在异常值或非正常分布数据条件下依然保持可靠性的学科。传统的统计方法，如最小二乘法(LS)，虽然在处理正态分布数据时表现良好，但对异常值极其敏感，可能导致估计结果失真。稳健统计通过引入对异常值不敏感的估计方法，提高估计的稳健性和准确性。

M估计是稳健统计中的一种重要方法，其核心思想是通过选择适当的目标函数 $\rho$ ，使得估计量对异常值的影响降到最低。具体来说，假设有一批观测值 $x_{ij}=U+{\delta }_{ij}$ ，其中U为估计值， ${\delta }_{ij}$ 为误差项。M估

计通过最小化 ${\displaystyle \underset{i,j}{\sum }\rho \left(x_{ij}-U\right)}$ 来估计，即U的M估计量是使得其残差的某个函数值达到最小时的值。目标

函数 $\rho$ 的选择影响估计的稳健性。选定合适的目标函数后，求 $\rho$ 关于U的导数并设其和为零，得到更新后的U，将新的估计值带入 $E\left(\Theta \right)$ 对 ${\sigma }_n$ 的导数，并令其为零，即可得到新的更新公式。具体操作如下。

在迭代步骤设置误差函数 $E\left(\Theta \right)$ 对 ${\sigma }_n$ 的导数，我们有

$\frac{\partial E\left(\Theta \right)}{\partial {\sigma }_n}={\displaystyle \underset{i,j}{\sum }{z}_{ijn}}\left\{\frac{1}{{\sigma }_n}+\frac{1}{{\sigma }_n^2}{\left(\frac{\Gamma \left(3/{\lambda }_n\right)}{\Gamma \left(1/{\lambda }_n\right)}\right)}^{1/2}{\left|x_{ij}-u_i{\nu }_j^{\text{T}}\right|}^{{\lambda }_n}\right\}$ (17)

根据稳健统计理论，任何估计值U是由一个隐式方程定义的：

${\displaystyle \underset{ij}{\sum }\Upsilon }\left(x_{ij}-U\right)=0$ (18)

这就给出了U位置作为加权平均值的数值解：

$U=\frac{{\displaystyle \underset{ij}{\sum }{\omega }_{ij}}{x}_{ij}}{{\displaystyle \underset{ij}{\sum }{\omega }_{ij}}}$ , ${\omega }_{ij}=\frac{\Upsilon \left(x_{ij}-U\right)}{x_{ij}-U}$ (19)

其中， $\Upsilon \left(r\right)$ 是目标函数 $\rho \left(r\right)$ 的导数：

$\rho \left(x_{ij}-U\right)={\left(\frac{\left|x_{ij}-U\right|}{{\sigma }_n}\right)}^{{\lambda }_n}$ (20)

$\Upsilon \left(x_{ij}-U\right)=\frac{\partial \rho \left(x_{ij}-U\right)}{\partial \left(x_{ij}-U\right)}=\frac{{\lambda }_n}{{\sigma }_n^{\lambda }}{\left|x_{ij}-U\right|}^{\lambda -1}\text{sign}\left(x_{ij}-U\right)$ (21)

将其带入权重公式：

${\omega }_{ij}=\frac{{\lambda }_n}{{\sigma }_n^{{\lambda }_n}}{\left|x_{ij}-U\right|}^{{\lambda }_n-1}$ (22)

利用权重更新U，将新的估计值U带入到偏导数公式中，并求解使偏导数为零的 ${\sigma }_n$ 更新公式：

${\displaystyle \underset{i,j}{\sum }{z}_{ijn}}\left\{\frac{1}{{\sigma }_n}+\frac{1}{{\sigma }_n^2}{\left(\frac{\Gamma \left(3/{\lambda }_n\right)}{\Gamma \left(1/{\lambda }_n\right)}\right)}^{1/2}{\left|x_{ij}-U\right|}^{{\lambda }_n}\right\}=0$ (23)

简化并解出 ${\sigma }_n$ 的表达式：

${\sigma }_n^{\left(t+1\right)}={\left[\frac{{\lambda }_n{\left(\frac{\Gamma \left(3/{\lambda }_n\right)}{\Gamma \left(1/{\lambda }_n\right)}\right)}^{{\lambda }_n/2}{\displaystyle \underset{i,j}{\sum }{z}_{ijn}}{\left|x_{ij}-U\right|}^{{\lambda }_n}}{{\displaystyle \underset{i,j}{\sum }{z}_{ijn}}}\right]}^{1/{\lambda }_n}$ (24)

3.2.3. 形状参数更新

下一步是更新 ${\lambda }_n$ 参数。这涉及到保持其他参数不变，并使用Newton Raphson方法改进 ${\lambda }_n$ 的估计。每次迭代都要求目标函数 $E\left(\Theta \right)$ 对参数 ${\lambda }_n$ 求一阶导数和二阶导数：

$\frac{\partial E\left(\theta \right)}{\partial {\lambda }_n}=\frac{1}{{\sigma }_n}\frac{1}{2}\left[\frac{1}{{\lambda }_n^2}\left(\frac{{\Gamma }^{\prime }\left(\frac{3}{{\lambda }_n}\right)}{\Gamma \left(\frac{3}{{\lambda }_n}\right)}-\frac{{\Gamma }^{\prime }\left(\frac{1}{{\lambda }_n}\right)}{\Gamma \left(\frac{1}{{\lambda }_n}\right)}\right)\right]{\displaystyle \underset{ij}{\sum }{z}_{ijn}}{\left|x_{ij}-u_i{v}_j^{\text{T}}\right|}^{{\lambda }_n}$ (25)

${\lambda }_n^{\left(t+1\right)}={\lambda }_n^{\left(t\right)}-\alpha \frac{\partial E\left({\lambda }_n\right)}{\partial {\lambda }_n}/\frac{{\partial }^{2}E\left({\lambda }_n\right)}{\partial {\lambda }_n^2+\gamma }$ (26)

3.2.4. U，V更新

U，V是通过求解 $\underset{U,V}{\min }{\displaystyle \underset{n=1}{\overset{N}{\sum }}{\Vert W\odot \left(X-U{V}^{\text{T}}\right)\Vert }_{{\lambda }_n}^{{\lambda }_n}}$ 来评估的。为了计算公式(12)，本文采用ALMF算法(改

进的梯度下降算法)，U和V矩阵的更新基于Adam优化器，并引入了动量项、学习率调度、梯度裁剪和Lookahead技术。

首先，计算误差矩阵和带权重的误差。误差矩阵 $L=U{V}^{\text{T}}$ 。带权重的误差 $E_w$ ：

$E_w=W\odot \left(L-X\right)$ (27)

其中，W是输入的权重矩阵，X是原始数据矩阵， $\odot$ 是逐元素乘法。

基于误差矩阵，分别U和V进行梯度计算，并引入正则化项 ${\lambda }_1$ 和 ${\lambda }_2$ ，得到U和V的梯度计算结果：

$g_U=E_w\cdot V+{\lambda }_{1}U$ (28)

$g_V=E_w^{\text{T}}\cdot U+{\lambda }_{2}V$ (29)

为了避免梯度爆炸问题，梯度裁剪将梯度限制在一个阈值c之内：

$g_U\text{}=\max \left(\min \left(g_U\text{},c\right),-c\right)$ (30)

$g_V\text{}=\max \left(\min \left(g_V\text{},c\right),-c\right)$ (31)

接下来，通过Adam优化器更新梯度，分别计算一阶和二阶动量估计 $m_u$ ， $v_u$ 以及 $m_V$ ， $v_V$ ，并对动量进行偏差校正，得到修正后的 ${\hat{m}}_U$ ， ${\hat{v}}_U$ 和， ${\hat{m}}_V$ ， ${\hat{v}}_V$ 。利用这些修正后的动量更新矩阵U和V，同时引入权重衰减项，使得优化过程更加稳定。

一阶和二阶动量估计为：

$\{\begin{array}{l}{m}_U={\beta }_1{m}_U+\left(1-{\beta }_1\right)g_U\\ v_U={\beta }_2{v}_U+\left(1-{\beta }_2\right)g_U^2\\ m_V={\beta }_1{m}_V+\left(1-{\beta }_1\right)g_V\\ v_V={\beta }_2{v}_V+\left(1-{\beta }_2\right)g_V^2\end{array}$ (32)

偏差校正结果为：

$\{\begin{array}{l}{\hat{m}}_U=\frac{m_U}{1-{\beta }_{1t}}\\ {\hat{v}}_U=\frac{v_U}{1-{\beta }_{2t}}\\ {\hat{m}}_V=\frac{m_V}{1-{\beta }_{1t}}\\ {\hat{v}}_V=\frac{v_V}{1-{\beta }_{2t}}\end{array}$ (33)

使用校正后的动量来更新U和V，同时引入权重衰减项：

$\{\begin{array}{l}U=U-\eta \left(\frac{{\hat{m}}_U}{\sqrt{{\hat{v}}_U}+ϵ}+{\lambda }_{1}U\right)\\ V=V-\eta \left(\frac{{\hat{m}}_V}{\sqrt{{\hat{v}}_V}+ϵ}+{\lambda }_{2}V\right)\end{array}$ (34)

其中， $\eta$ 是学习率， $ϵ$ 是一个防止分母为零的小常数。

每隔k次迭代时执行Lookahead更新，通过在原始更新基础上慢更新U和V，并将慢更新结果赋值给U和V。

$\{\begin{array}{l}{U}_{slow\text{}}=U_{slow}+0.5\left(U-U_{slow}\right)\\ V_{slow}=V_{slow}+0.5\left(V-V_{slow}\right)\end{array}$ (35)

$U=U_{slow}\text{}$ ， $V=V_{slow}$ 。

每次迭代后，重新计算误差：

$\text{error}={\Vert W\odot \left(X-L\right)\Vert }_F+{\lambda }_1{\Vert U\Vert }_F+{\lambda }_2{\Vert V\Vert }_F$ (36)

如果误差小于设定的阈值，则停止迭代。采用ALMF方法更新U，V的过程由算法1给出，针对GGMM模型提出的算法总结在算法2中。

算法1 ALMF算法求解公式(12)

输入：初始化U，V， ${\lambda }_1$ ， ${\lambda }_2$ ， ${\beta }_1$ ， ${\beta }_2$ ，c。

1. 通过(27)计算 $E_w$ ;

2. 通过(28)(29)计算 $g_U$ ， $g_V$ ;

3. 通过(32)(33)(34)更新 $m_u$ ， $v_u$ ， $m_V$ ， $v_V$ ，U，V；

4. 通过(35)更新 $U_{slow\text{}}$ ， $V_{slow}$

输出U，V。

算法2 模型求解算法

输入：数据X；算法参数：秩r和 $\delta$ 。

初始化 ${\Theta }^{\left(t\right)}=\left\{{\sigma }^{\left(t\right)},{\lambda }^{\left(t\right)},U^{\left(t\right)},V^{\left(t\right)}\right\}$ ，t = 0。

1：通过(16)更新 ${\pi }^{\left(t\right)}$ ；

2：通过(24)更新 ${\sigma }^{\left(t\right)}$ ；

3：通过(26)更新 ${\lambda }^{\left(t\right)}$ ；

4：通过算法1更新 $U^{\left(t\right)}$ 、 $V^{\left(t\right)}$ ;

5：t = t + 1；

输出：参数 $\Theta$

4. 实验

在本节中，为了验证算法的有效性，进行了大量的合成数据实验和真实数据实验，并与MoG WLRMF、GGMM + EPLL进行对比。采用峰值信噪比(PSNR)、结构相似性(SSIM)、特征相似性(FSIM)、相对标准误差(RSE)四种评价指标来评价每种方法的性能。

4.1. 合成数据实验

为了验证算法的有效性，选择了干净图像Facade数据集进行合成数据实验。并设置了4种不同的噪声设置，分别为：

1) 加入零均值高斯噪声，方差为0.06；

2) 加入零均值高斯噪声，方差为30；

3) 加入混合噪声：[−5, 5]均匀噪声，(0, 0.2)和(0, 0.01)的高斯噪声

4) 加入混合噪声：[−35, −35]均匀噪声，(0, 20)和(0, 10)的高斯噪声

表1展示了Facade数据集上4种不同噪声情况下的指标数据，最佳结果由粗体表示。根据表中数据得出，在零均值高斯噪声情况下，本文提出的算法的各项指标均优于另外两个算法；在高斯和均匀噪声混合的情况下，只有FSIM在第三种噪声情况下MoG WLRMF算法表现更好，其他情况均本文提出的算法表现更好。

图1展示了三种去噪方法对Facade数据集的去噪效果图，从图中可以看出，MoG WLRMF对零均值噪声表现效果相对较差，GGMM + EPLL可以对噪声有效去除，相比之下，本文的方法获得的视觉效果更好。

Table 1. Comparison of denoising effects of different algorithms on Facade datasets under different noise conditions

表1. 不同算法在不同噪声情况下对Facade数据集去噪效果对比

情况	评价指标	MoG WLRMF	GGMM + EPLL	Ours
1	PSNR	17.91	19.29	28.48
	SSIM	0.4865	0.8020	0.8526
	FSIM	0.9891	0.9670	0.9978
	RSE	0.2442	0.2120	0.0737
2	PSNR	15.29	18.93	28.42
	SSIM	0.1212	0.6090	0.8521
	FSIM	0.8824	0.9420	0.9966
	RSE	0.2420	0.2210	0.0768
3	PSNR	25.65	5.93	28.49
	SSIM	0.8236	0.0090	0.8528
	FSIM	0.9961	0.8930	0.9898
	RSE	0.1416	0.9880	0.0742
4	PSNR	23.69	5.64	28.45
	SSIM	0.7862	0.0040	0.8528
	FSIM	0.9268	0.8930	0.9878
	RSE	0.1260	1.0210	0.0736

Figure 1. Comparison of Facade datasets by various denoising methods in Case 4

图1. 情况1下各种去噪方法对Facade数据集的对比图

4.2. 真实数据实验

真实高光谱图像的数据实验使用HYDICE城市图像进行演示，此数据集的图像受到大气和水吸收等噪声的严重污染，图2展示了HYDICE城市影像中3个严重污染波段的恢复结果，可以看出，该算法在处理这样的真实粗噪声时仍然有很好的性能。

Figure 2. Real hyperspectral image restoration. (a) Original polluted bands; (b) Corresponding bands recovered by the proposed algorithm

图2. 真实的高光谱图像恢复。(a) 原始污染波段；(b) 本文算法恢复后的相应波段

5. 结论

本文提出了一种新的基于GGMM的图像去噪算法，通过将GGMM引入低秩矩阵分解框架，对图像中的噪声进行建模，从而能够处理不同类型和复杂度的噪声。算法的参数通过EM算法进行求解。本文算法在多种公共数据集上的实验结果表明，与现有经典去噪方法(如MoG WLRMF、GGMM + EPLL等)相比，在PSNR和SSIM指标上具有显著提升，尤其在处理尖峰噪声和重尾噪声时表现出更强的鲁棒性和细节恢复能力。同时主观视觉质量也显著优于其他方法，充分验证了本研究的先进性和实用价值。

基金项目

吉林省自然科学基金，NO. 20240101298JC；国家自然科学基金，NO. 12171054。

参考文献

[1]	翟雪. 基于广义高斯混合模型的数字图像修复算法研究[D]: [硕士学位论文]. 长春: 东北师范大学, 2021.
[2]	刘利平, 乔乐乐, 蒋柳成. 图像去噪方法概述[J]. 计算机科学与探索, 2021, 15(8): 1418-1431.
[3]	邹永宁, 姚功杰. 自适应窗口形状的中值滤波[J]. 光学精密工程, 2018, 26(12): 3028-3039.
[4]	Dong, W., Zhang, L., Shi, G. and Li, X. (2013) Nonlocally Centralized Sparse Representation for Image Restoration. IEEE Transactions on Image Processing, 22, 1620-1630. https://doi.org/10.1109/tip.2012.2235847
[5]	Chen, F., Zhang, L. and Yu, H. (2015) External Patch Prior Guided Internal Clustering for Image Denoising. 2015 IEEE International Conference on Computer Vision (ICCV), Santiago, 7-13 December 2015, 603-611. https://doi.org/10.1109/iccv.2015.76
[6]	Deledalle, C., Parameswaran, S. and Nguyen, T.Q. (2018) Image Denoising with Generalized Gaussian Mixture Model Patch Priors. SIAM Journal on Imaging Sciences, 11, 2568-2609. https://doi.org/10.1137/18m116890x
[7]	Meng, D. and De la Torre, F. (2013) Robust Matrix Factorization with Unknown Noise. 2013 IEEE International Conference on Computer Vision, Sydney, 1-8 December 2013, 1377-1344. https://doi.org/10.1109/iccv.2013.169
[8]	Chen, X., Han, Z., Wang, Y., Zhao, Q., Meng, D. and Tang, Y. (2016) Robust Tensor Factorization with Unknown Noise. 2016 IEEE Conference on Computer Vision and Pattern Recognition (CVPR), Las Vegas, 27-30 June 2016, 5213-5221. https://doi.org/10.1109/cvpr.2016.563
[9]	Cao, X., Zhao, Q., Meng, D., Chen, Y. and Xu, Z. (2016) Robust Low-Rank Matrix Factorization under General Mixture Noise Distributions. IEEE Transactions on Image Processing, 25, 4677-4690. https://doi.org/10.1109/tip.2016.2593343
[10]	Chen, Y., Cao, X., Zhao, Q., Meng, D. and Xu, Z. (2018) Denoising Hyperspectral Image with Non-I.I.D. Noise Structure. IEEE Transactions on Cybernetics, 48, 1054-1066. https://doi.org/10.1109/tcyb.2017.2677944
[11]	Marnissi, Y., Zheng, Y., Chouzenoux, E. and Pesquet, J. (2017) A Variational Bayesian Approach for Image Restoration—Application to Image Deblurring with Poisson-Gaussian Noise. IEEE Transactions on Computational Imaging, 3, 722-737. https://doi.org/10.1109/tci.2017.2700203
[12]	Luo, E., Chan, S.H. and Nguyen, T.Q. (2016) Adaptive Image Denoising by Mixture Adaptation. IEEE Transactions on Image Processing, 25, 4489-4503. https://doi.org/10.1109/tip.2016.2590318
[13]	Amudala, S., Ali, S. and Bouguila, N. (2020) Variational Inference of Infinite Generalized Gaussian Mixture Models with Feature Selection. 2020 IEEE International Conference on Systems, Man, and Cybernetics (SMC), Toronto, 11-14 October 2020, 120-127. https://doi.org/10.1109/smc42975.2020.9283007
[14]	Alroobaea, R., Rubaiee, S., Bourouis, S., Bouguila, N. and Alsufyani, A. (2019) Bayesian Inference Framework for Bounded Generalized Gaussian-Based Mixture Model and Its Application to Biomedical Images Classification. International Journal of Imaging Systems and Technology, 30, 18-30. https://doi.org/10.1002/ima.22391
[15]	Chen, B., Yang, Z. and Yang, Z. (2018) An Algorithm for Low-Rank Matrix Factorization and Its Applications. Neurocomputing, 275, 1012-1020. https://doi.org/10.1016/j.neucom.2017.09.052
[16]	Guo, Q., Zhang, C., Zhang, Y. and Liu, H. (2016) An Efficient SVD-Based Method for Image Denoising. IEEE Transactions on Circuits and Systems for Video Technology, 26, 868-880. https://doi.org/10.1109/tcsvt.2015.2416631
[17]	Bouwmans, T., Javed, S., Zhang, H., Lin, Z. and Otazo, R. (2018) On the Applications of Robust PCA in Image and Video Processing. Proceedings of the IEEE, 106, 1427-1457. https://doi.org/10.1109/jproc.2018.2853589
[18]	Chen, T., Xiang, Q., Zhao, D. and Sun, L. (2023) An Unsupervised Image Denoising Method Using a Nonconvex Low-Rank Model with TV Regularization. Applied Sciences, 13, Article 7184. https://doi.org/10.3390/app13127184