零均值高斯分布变量的幂的数学期望以及噪声图像块的匹配误差的方差

1 零均值高斯分布变量的幂的数学期望

       现实中我们经常会遇到服从高斯分布的随机变量以及该变量的幂,因为

D ( X ) = E ( X 2 ) − [ E ( X ) 2 ] , (1.1) D(X) = E({X^2}) - [E{(X)^2}], \tag{1.1} D(X)=E(X2)[E(X)2],(1.1)

所以,求解该变量的幂的数学期望对我们求解该变量的幂的方差十分有帮助。但是,显式地给出所有符合高斯分布的变量的幂的数学期望并不是一件易事,这里我们主要讨论均值为零的高斯分布。一方面均值为零可以避免多项式中低阶项的出现,简化了后续的运算,同时概率密度函数关于 y y y 轴对称,这样奇数幂的数学期望很明显为零,我们只需讨论那些偶数幂的情况;另一方面,均值为零的高斯分布本身在现实中也具有十分广泛的应用,例如在信号处理中,高斯噪声的均值一般会被假设为零。因为任意的高斯分布都可以通过简单的线性变换转换为标准正态分布,即

X ∼ N ( 0 ,    σ 2 )    ⇒    Y = X σ ∼ N ( 0 ,    1 ) . (1.2) X \sim \mathcal{N}\left( {0,\;{\sigma ^2}} \right)\; \Rightarrow \;Y = \frac{X}{\sigma } \sim \mathcal{N}\left( {0,\;1} \right). \tag{1.2} XN(0,σ2)Y=σXN(0,1).(1.2)

所以我们可以进一步进行简化,只讨论标准正态分布的情况,而对于其他均值为零的高斯分布,则有

E ( X n ) = E ( σ n ⋅ Y n ) = σ n ⋅ E ( Y n ) . (1.3) E({X^n}) = E({\sigma ^n} \cdot {Y^n}) = {\sigma ^n} \cdot E({Y^n}). \tag{1.3} E(Xn)=E(σnYn)=σnE(Yn).(1.3)

       因为标准正态分布变量的奇数幂的数学期望为零,我们这里只讨论偶数幂的情况。根据数学期望的定义,有

E ( Y n ) = ∫ − ∞ + ∞ y n ⋅ 1 2 π e − y 2 2 d y = 2 n / 2 − 1 π ∫ − ∞ + ∞ ( y 2 2 ) ( n − 1 ) / 2 ⋅ e − y 2 2 d y 2 2 = 2 n / 2 π ∫ 0 + ∞ t ( n + 1 ) / 2 − 1 ⋅ e − t d t = 2 n / 2 π Γ ( n + 1 2 ) . (1.4) \begin{aligned} E\left( {{Y^n}} \right) &= \int_{ - \infty }^{ + \infty } {{y^n} \cdot \tfrac{1}{{\sqrt {2\pi } }}{e^{ - \frac{{{y^2}}}{2}}}dy} = \frac{{{2^{n/2 - 1}}}}{{\sqrt \pi }}\int_{ - \infty }^{ + \infty } {{{\left( {\tfrac{{{y^2}}}{2}} \right)}^{(n - 1)/2}} \cdot {e^{ - \frac{{{y^2}}}{2}}}d\tfrac{{{y^2}}}{2}} \\ &= \frac{{{2^{n/2}}}}{{\sqrt \pi }}\int_0^{ + \infty } {{t^{(n + 1)/2 - 1}} \cdot {e^{ - t}}dt} = \frac{{{2^{n/2}}}}{{\sqrt \pi }}\Gamma \left( {\frac{{n + 1}}{2}} \right). \\ \end{aligned} \tag{1.4} E(Yn)=+yn2π 1e2y2dy=π 2n/21+(2y2)(n1)/2e2y2d2y2=π 2n/20+t(n+1)/21etdt=π 2n/2Γ(2n+1).(1.4)

这里涉及到伽玛函数的知识。伽玛函数的定义为

Γ ( x ) = ∫ 0 + ∞ t x − 1 e − t d t ,    x > 0. (1.5) \Gamma (x) = \int_0^{ + \infty } {{t^{x - 1}}{e^{ - t}}dt,\;x > 0.} \tag{1.5} Γ(x)=0+tx1etdt,x>0.(1.5)

并且有以下性质,

Γ ( x + 1 ) = x Γ ( x ) ,    Γ ( 1 2 ) = π . (1.6) \Gamma (x + 1) = x\Gamma (x),\;\Gamma (\tfrac{1}{2}) = \sqrt \pi . \tag{1.6} Γ(x+1)=xΓ(x),Γ(21)=π .(1.6)

那么可得

E ( Y n ) = 2 n / 2 π ( n − 1 2 ) ( n − 3 2 ) … ( 1 2 ) Γ ( 1 2 ) = ( n − 1 ) ( n − 3 ) … ( 1 ) , n = 2 , 4 , 6 , . . . (1.7) \begin{aligned} E({Y^n}) &= \frac{{{2^{n/2}}}}{{\sqrt \pi }}\left( {\frac{{n - 1}}{2}} \right)\left( {\frac{{n - 3}}{2}} \right) \ldots \left( {\frac{1}{2}} \right)\Gamma \left( {\frac{1}{2}} \right) \\ &= (n - 1)(n - 3) \ldots (1),\quad n = 2,4,6,... \\ \end{aligned} \tag{1.7} E(Yn)=π 2n/2(2n1)(2n3)(21)Γ(21)=(n1)(n3)(1),n=2,4,6,...(1.7)

因此,零均值高斯分布变量的偶数幂的数学期望为

E ( X n ) = σ n ( n − 1 ) ( n − 3 ) … ( 1 ) . (1.8) E({X^n}) = {\sigma ^n}(n - 1)(n - 3) \ldots (1). \tag{1.8} E(Xn)=σn(n1)(n3)(1).(1.8)

2 噪声图像块的匹配误差的方差

       在数字图像处理中,加性高斯白噪声(Additive White Gaussian Noise, AWGN)对各个像素的影响一般可认为是独立同分布的(i.i.d.),即可表示为

v ( x ) = u ( x ) + n ( x ) ,    n ( x ) ∼ N ( 0 ,    σ 2 ) . (2.1) v(x) = u(x) + n(x),\;n(x) \sim \mathcal{N}(0,\;{\sigma ^2}). \tag{2.1} v(x)=u(x)+n(x),n(x)N(0,σ2).(2.1)

其中为 u u u 原始无噪声图像, n n n 为AWGN, v v v 为有噪声图像, x x x 代表像素位置。在这里,噪声 n n n 为随机变量,而无噪声图像 u u u 为常量,所以有

v ( x ) ∼ N ( u ( x ) ,    σ 2 ) . (2.2) v(x) \sim \mathcal{N}\left( {u(x),\;{\sigma ^2}} \right). \tag{2.2} v(x)N(u(x),σ2).(2.2)

对于有限个不同的像素,把它们进行线性叠加,则有

∑ i a i v ( x i ) ∼ N ( ∑ i a i u ( x i ) ,    ∑ i a i 2 σ 2 ) . (2.3) \sum\limits_i {{a_i}v({x_i})} \sim \mathcal{N}\left( {\sum\limits_i {{a_i}u({x_i})} ,\;\sum\limits_i {a_i^2} {\sigma ^2}} \right). \tag{2.3} iaiv(xi)N(iaiu(xi),iai2σ2).(2.3)

因此,如果这些像素的真实值 u u u 相近,那么对它们的和求平均,得到的新的像素值将是一个方差为原始噪声方差的 1 / N 1/N 1/N 的随机变量,从而极大地压缩噪声的量级,提升图像的质量。

       然而,由于我们通常只能得到添加了噪声以后的图像 v v v,要寻找具有相似真实值 u u u 的像素不是一件易事,特别是在噪声水平较高的情况下。如果我们把那些真实值 u u u 相差较大的像素叠加在一起,就可能导致原始图像细节的破坏,造成较大的失真。那么,如果我们直接在有噪声图像 v v v 上对两个像素的相似性进行判断,这种方法的可行性又如何呢?因为有噪声图像 v v v 是一个随机变量,相应地两个像素之间的误差也同样是一个随机变量,所以我们其实可以通过求解像素误差的数学期望以及方差来判断直接在有噪声图像 v v v 上寻找相似像素点可行性。

       有噪声图像 v v v 上两个像素点的误差可表示为

d v ( x 1 ,    x 2 ) = ( v ( x 1 ) − v ( x 2 ) ) 2 = [ ( u ( x 1 ) − u ( x 2 ) ) + ( n ( x 1 ) − n ( x 2 ) ) ] 2 . (2.4) {d_v}\left( {{x_1},\;{x_2}} \right) = {\left( {v\left( {{x_1}} \right) - v\left( {{x_2}} \right)} \right)^2} = {\left[ {\left( {u\left( {{x_1}} \right) - u\left( {{x_2}} \right)} \right) + \left( {n\left( {{x_1}} \right) - n\left( {{x_2}} \right)} \right)} \right]^2}.\tag{2.4} dv(x1,x2)=(v(x1)v(x2))2=[(u(x1)u(x2))+(n(x1)n(x2))]2.(2.4)

这其实属于非中心的卡方分布,即均值不为零的高斯分布变量的平方和的分布。为了简便,后面的推导中我们将 x x x 省去。注意,图像的真实像素值 u u u 为常量,而各个像素上的噪声 n n n 是独立同分布的,那么误差的数学期望为

E ( d v ) = E [ ( u 1 − u 2 ) 2 + 2 ( u 1 − u 2 ) ( n 1 − n 2 ) + ( n 1 − n 2 ) 2 ] = ( u 1 − u 2 ) 2 + 2 ( u 1 − u 2 ) E ( n 1 − n 2 ) + E [ ( n 1 − n 2 ) 2 ] = d u + D ( n 1 − n 2 ) + [ E ( n 1 − n 2 ) ] 2 = d u + D ( n 1 − n 2 ) = d u + 2 σ 2 . (2.5) \begin{aligned} E\left( {{d_v}} \right) &= E\left[ {{{({u_1} - {u_2})}^2} + 2({u_1} - {u_2})({n_1} - {n_2}) + {{({n_1} - {n_2})}^2}} \right] \\ &= {({u_1} - {u_2})^2} + 2({u_1} - {u_2})E\left( {{n_1} - {n_2}} \right) + E\left[ {{{({n_1} - {n_2})}^2}} \right] \\ &= {d_u} + D\left( {{n_1} - {n_2}} \right) + {\left[ {E\left( {{n_1} - {n_2}} \right)} \right]^2} \\ &= {d_u} + D\left( {{n_1} - {n_2}} \right) = {d_u} + 2{\sigma ^2}. \\ \end{aligned} \tag{2.5} E(dv)=E[(u1u2)2+2(u1u2)(n1n2)+(n1n2)2]=(u1u2)2+2(u1u2)E(n1n2)+E[(n1n2)2]=du+D(n1n2)+[E(n1n2)]2=du+D(n1n2)=du+2σ2.(2.5)

因此,从数学期望上看, d v {d_v} dv d u {d_u} du 的大小关系是一样的,在真实无噪声图像 u u u 上两个像素的误差相比于其他像素较小,那么相应地在有噪声图像 v v v 上的误差的数学期望相比于其他像素也会比较小,所以直接在有噪声图像 v v v 上去判断两个像素是否相似确实有一定的可行性。

       然而,虽然数学期望的大小关系一样,如果其方差很大,那么我们还是有比较大的概率匹配出错的,即把真实图像 u u% u 上亮度相差较大的两个像素在有噪声图像 v v v 上判定为相似。所以,接下来我们还得求解 d v {d_v} dv 的方差,即

D ( d v ) = D [ ( u 1 − u 2 ) 2 + 2 ( u 1 − u 2 ) ( n 1 − n 2 ) + ( n 1 − n 2 ) 2 ] = 4 ( u 1 − u 2 ) 2 D ( n 1 − n 2 ) + D [ ( n 1 − n 2 ) 2 ] = 8 σ 2 d u + E [ ( n 1 − n 2 ) 4 ] − { E [ ( n 1 − n 2 ) 2 ] } 2 = 8 σ 2 d u + E [ ( n 1 − n 2 ) 4 ] − 4 σ 4 . (2.6) \begin{aligned} D\left( {{d_v}} \right) &= D\left[ {{{({u_1} - {u_2})}^2} + 2({u_1} - {u_2})({n_1} - {n_2}) + {{({n_1} - {n_2})}^2}} \right] \\ &= 4{\left( {{u_1} - {u_2}} \right)^2}D\left( {{n_1} - {n_2}} \right) + D\left[ {{{({n_1} - {n_2})}^2}} \right] \\ &= 8{\sigma ^2}{d_u} + E\left[ {{{({n_1} - {n_2})}^4}} \right] - {\left\{ {E\left[ {{{({n_1} - {n_2})}^2}} \right]} \right\}^2} \\ &= 8{\sigma ^2}{d_u} + E\left[ {{{({n_1} - {n_2})}^4}} \right] - 4{\sigma ^4}. \\ \end{aligned} \tag{2.6} D(dv)=D[(u1u2)2+2(u1u2)(n1n2)+(n1n2)2]=4(u1u2)2D(n1n2)+D[(n1n2)2]=8σ2du+E[(n1n2)4]{E[(n1n2)2]}2=8σ2du+E[(n1n2)4]4σ4.(2.6)

因为

( n 1 − n 2 ) ∼ N ( 0 ,    2 σ 2 ) . (2.7) \left( {{n_1} - {n_2}} \right) \sim \mathcal{N}\left( {0,\;2{\sigma ^2}} \right).\tag{2.7} (n1n2)N(0,2σ2).(2.7)

根据式(1.8)可得

E [ ( n 1 − n 2 ) 4 ] = ( 2 σ ) 4 ⋅ 3 ⋅ 1 = 12 σ 4 . (2.8) E\left[ {{{\left( {{n_1} - {n_2}} \right)}^4}} \right] = {\left( {\sqrt 2 \sigma } \right)^4} \cdot 3 \cdot 1 = 12{\sigma ^4}.\tag{2.8} E[(n1n2)4]=(2 σ)431=12σ4.(2.8)

所以有

D ( d v ) = 8 σ 2 d u + 8 σ 4 . (2.9) D\left( {{d_v}} \right) = 8{\sigma ^2}{d_u} + 8{\sigma ^4}.\tag{2.9} D(dv)=8σ2du+8σ4.(2.9)

可以看到,有噪声图像 v v v 上两个像素的误差的方差的量级大约为噪声方差的平方,当噪声水平较高时,这种计算方法带来的出错可能性将非常的大,所以对于单个像素来说这种方法是不太可行的。

       既然对于单个像素的计算的方差太大,一种直观的解决办法就是增加像素的个数,例如使用以当前像素为中心的 N × N N \times N N×N 的块进行计算。假设图像中两个块不重叠,那么块中每个像素的误差都是独立的,于是有

E ( 1 N 2 ∑ i d v , i ) = 1 N 2 ∑ i d u , i + 2 σ 2 . (2.10) E\left( {\frac{1}{{{N^2}}}\sum\limits_i {{d_{v,i}}} } \right) = \frac{1}{{{N^2}}}\sum\limits_i {{d_{u,i}}} + 2{\sigma ^2}. \tag{2.10} E(N21idv,i)=N21idu,i+2σ2.(2.10)

D ( 1 N 2 ∑ i d v , i ) = 8 σ 2 N 4 ∑ i d u , i + 8 σ 4 N 2 . (2.11) D\left( {\frac{1}{{{N^2}}}\sum\limits_i {{d_{v,i}}} } \right) = \frac{{8{\sigma ^2}}}{{{N^4}}}\sum\limits_i {{d_{u,i}}} + \frac{{8{\sigma ^4}}}{{{N^2}}}.\tag{2.11} D(N21idv,i)=N48σ2idu,i+N28σ4.(2.11)

因为真实无噪声图像 u u u 中邻域内的像素值一般相差不大,式(2.11)中的方差相当于缩小了 N 2 N^2 N2 倍,所以使用一个邻域块来判断有噪声图像 v v v 中两个像素是否相似相比于单个像素要准确很多,这也是包括非局部平均(Non-Local Means)[1] 以及 BM3D(Block-Matching 3D)[2] 等图像降噪算法使用块进行相似像素匹配的原因。不过,由于计算复杂度的限制,邻域的大小一般不会设置得太大,当噪声的方差很大时,空域块匹配的误差依然会比较大,所以这时 BM3D 算法会先对有噪声图像进行频域的转换,去除那些能量较低的分量,然后直接在频域进行块的匹配,从而降低块匹配的错误率。另外,虽然式(2.10)和(2.11)是基于不重叠块分析的,但实际应用中一般也会对有重叠块进行相似度计算。

参考资料

[1] Antoni Buades, Bartomeu Coll, Jean-Michel Morel. A review of image denoising algorithms, with a new one. Multiscale Modeling and Simulation: A SIAM Interdisciplinary Journal, Society for Industrial and Applied Mathematics, 2005, 4 (2), pp.490-530. ffhal-00271141. https://hal.archives-ouvertes.fr/hal-00271141/document.
[2] Dabov, Kostadin & Foi, Alessandro & Katkovnik, Vladimir & Egiazarian, Karen. (2007). Image Denoising by Sparse 3-D Transform-Domain Collaborative Filtering. IEEE transactions on image processing : a publication of the IEEE Signal Processing Society. 16. 2080-95. 10.1109/TIP.2007.901238. https://www.researchgate.net/publication/6151802_Image_Denoising_by_Sparse_3-D_Transform-Domain_Collaborative_Filtering.

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值