基于概率矩阵分解的高光谱图像锐化方法与流程

技术特征：

1.基于概率矩阵分解的高光谱图像锐化方法，其特征在于：

所述高光谱图像锐化是指在通信接收端把输入的一幅低分辨的高光谱图像和一幅高分辨的多光谱图像融合生成一幅高分辨的高光谱图像，从而折中地解决遥感领域中图像的空间分辨率和光谱分辨率矛盾的一种方法，所述高光谱图像锐化方法是在通信接收端的一台计算机中依次按一下步骤实现的：

步骤(1)：输入要求参数，实现计算机初始化

用一台机载可见光/红外成像光谱仪对同一目标拍摄一幅低分辨的高光谱图像简称为X，以及一幅高分辨的多光谱图像简称为Y，其中：

L为所述图像X的光谱通道数，l为所述图像Y的光谱通道数，l＜＜L,

n为所述图像X的各光谱波段的成像像素点数，N为所述图像Y的各光谱波段的成像像素点数，n＜＜N；

输入多光谱摄像头对应高光谱图像的频率响应矩阵简称光谱响应矩阵F；

手动设置分解矩阵维数r和算法迭代次数Δ，其中参数r控制分解矩阵的行数；

令所述融合后的高分辨率的高光谱图像为简称为待求图像Z，光学成像原理表明：X＝ZB，Y＝FZ，其中为高分辨图像到低分辨图像的响应矩阵；

步骤(2)：对所述图像X用线性插值运算矩阵C插值，且是一个已知固定常数矩阵；由于X＝ZB，因此，得到的插值图像满足即所述插值图像是待求图像Z经分辨率的降低后插值形成的图像；

步骤(3)：在下述假设和推导下，待求图像Z近似为表示转置运算；

假设：所述图像Z中的像素光谱矢量是由少量隐藏的光谱特征矢量经过线性叠加形成，这些光谱特征矢量为矩阵U的每一行，像素光谱矢量对应的叠加系数为矩阵T的每一行，则有Z＝U'T+R，其中r表示分解矩阵U的行数，也表示隐藏光谱特征矢量的个数，该参数在步骤(1)中人为预先设定；残余项趋近于0，故设R为各向同性的高斯噪声；

令G＝BC,W＝TG,V＝T-W,N＝RG,则若近似认为所述残余项R满足R≈N，则待求图像Z满足其中

步骤(4)：按以下步骤计算所述图像Y的残差变换图像其中

步骤(4.1)：对所述图像Y进行预处理，消去其中只含在所述插值图像中的高光谱成分，得到残差图像矩阵

光学成像原理表明：Y＝FZ，所以

步骤(4.2)：通过残差图像变换矩阵对所述残差图像矩阵E进行变换，以便把所述残差图像E中各向异性高斯噪声转变为和所述残余项同分布的各向同性高斯噪声，其中：分解矩阵Q和分解矩阵D通过对矩阵(FF')^-1进行奇异值分解获得，Q满足QDQ'＝(FF')^-1，D为对角矩阵，表示对对角矩阵D中元素进行算术平方根运算；通过残差图像变换矩阵Φ，得到变换后的光谱响应矩阵以及残差变换图像

步骤(5)：根据计算所得的残差变换图像和线性插值图像用变分贝叶斯方法迭代求取待求图像Z的分解矩阵U,V和W，得到U'V，从而估算出待求图像Z，所用变分贝叶斯方法步骤如下：

步骤(5.1)：根据分解矩阵的维数r、所述图像Y的各光谱波段的成像像素点数N和所述图像X的光谱通道数L，见步骤(1)，确定分解矩阵和的大小,在此基础上初始化下列参数：

1_r×N表示大小为r×N的全1矩阵；

${\stackrel{\‾}{\mathrm{\α}}}_n^{\left(0\right)}={\stackrel{\‾}{\mathrm{\α}}}_u^{\left(0\right)}={\stackrel{\‾}{\mathrm{\α}}}_v^{\left(0\right)}={\stackrel{\‾}{\mathrm{\α}}}_w^{\left(0\right)}=1;$

I_r表示大小为r×r的单位矩阵；

${\[\stackrel{\‾}{{\mathrm{WW}}^{\′}}\]}^{\left(0\right)}={\stackrel{\‾}{W}}^{\left(0\right)}{\left({\stackrel{\‾}{W}}^{\left(0\right)}\right)}^{\′}+N{\Σ}_w^{\left(0\right)};{\[\stackrel{\‾}{{\mathrm{VV}}^{\′}}\]}^{\left(0\right)}={\stackrel{\‾}{V}}^{\left(0\right)}{\left({\stackrel{\‾}{V}}^{\left(0\right)}\right)}^{\′}+N{\Σ}_V^{\left(0\right)};$

${\stackrel{\‾}{U}}^{\left(0\right)}=0;{\Σ}_u^{\left(0\right)}=0;{\[\stackrel{\‾}{U{\tilde{F}}^{\′}\tilde{F}{U}^{\′}}\]}^{\left(0\right)}={\[\stackrel{\‾}{{\mathrm{UU}}^{\′}}\]}^{\left(0\right)}=0;$

步骤(5.2)：总迭代次数用Δ表示，第t次迭代的各分解矩阵的取值表示为(·)^(t)，符号的上划线表示符号的概率意义上的均值，t＝1,2,3,...,Δ；

步骤(5.3)：输入第t-1次迭代，即上一次迭代获得的参数按下式计算第t次迭代参数：和

$\left\{\begin{array}{c}{\Σ}_u^{\left(t\right)}={({\stackrel{\‾}{\mathrm{\α}}}_n^{(t-1)}{I}_L\⊗{\[\stackrel{\‾}{{\mathrm{WW}}^{\′}}\]}^{(t-1)}+{\stackrel{\‾}{\mathrm{\α}}}_n^{(t-1)}{\tilde{F}}^{\′}\tilde{F}\⊗{\[\stackrel{\‾}{{\mathrm{VV}}^{\′}}\]}^{(t-1)}+{\stackrel{\‾}{\mathrm{\α}}}_u^{(t-1)}{I}_{Lr})}^{-1}\\ vec\[{\stackrel{\‾}{U}}^{\left(t\right)}\]={\stackrel{\‾}{\mathrm{\α}}}_n^{(t-1)}{\Σ}_u^{\left(t\right)}vec\[{\stackrel{\‾}{W}}^{(t-1)}{\tilde{X}}^{\′}+{\stackrel{\‾}{V}}^{(t-1)}{\tilde{E}}^{\′}\tilde{F}\]\\ {\[\stackrel{\‾}{{\mathrm{UU}}^{\′}}\]}^{\left(t\right)}={\stackrel{\‾}{U}}^{\left(t\right)}{\left({\stackrel{\‾}{U}}^{\left(t\right)}\right)}^{\′}+\underset{i=1}{\overset{L}{\Σ}}\underset{j=1}{\overset{L}{\Σ}}{\[{\Σ}_u^{\left(t\right)}\]}_{r\×r}^{i,j}\\ {\[\stackrel{\‾}{U{\tilde{F}}^{\′}\tilde{F}{U}^{\′}}\]}^{\left(t\right)}={\stackrel{\‾}{U}}^{\left(t\right)}{\tilde{F}}^{\′}\tilde{F}{\left({\stackrel{\‾}{U}}^{\left(t\right)}\right)}^{\′}+\underset{i=1}{\overset{L}{\Σ}}\underset{j=1}{\overset{L}{\Σ}}{\[{\tilde{F}}^{\′}\tilde{F}\]}_{1\×1}^{i,j}{\[{\Σ}_u^{\left(t\right)}\]}_{r\×r}^{i,j},\end{array}\right.,$

其中I_L表示大小为L×L的单位矩阵；I_Lr表示大小为Lr×Lr的单位矩阵；L为所述图像X的光谱通道数；r为步骤(1)输入的分解矩阵的维数；表示两个矩阵的Kronecker积；vec[·]表示矩阵的矢量化表示形式；为一种运算，定义如下：

其中a_i,j表示矩阵A第i行第j列的元素；当r＝1时，该运算值为a_i,j；否则为一个大小为r×r的矩阵；

步骤(5.4)：输入第t-1次迭代，即上一次迭代的获得的参数以及步骤(5.3)中第t次迭代得到的参数根据下列式子，计算获得和

$\left\{\begin{array}{c}{\Σ}_V^{\left(t\right)}={({\stackrel{\‾}{\mathrm{\α}}}_n^{(t-1)}{\[\stackrel{\‾}{U{\tilde{F}}^{\′}\tilde{F}{U}^{\′}}\]}^{\left(t\right)}+{\stackrel{\‾}{\mathrm{\α}}}_v^{(t-1)}{I}_r)}^{-1}\\ {\stackrel{\‾}{V}}^{\left(t\right)}={\stackrel{\‾}{\mathrm{\α}}}_n^{(t-1)}{\Σ}_V^{\left(t\right)}{\stackrel{\‾}{U}}^{\left(t\right)}{\tilde{F}}^{\′}\tilde{E}\\ {\[\stackrel{\‾}{{\mathrm{VV}}^{\′}}\]}^{\left(t\right)}={\stackrel{\‾}{V}}^{\left(t\right)}{\left({\stackrel{\‾}{V}}^{\left(t\right)}\right)}^{\′}+N{\Σ}_V^{\left(t\right)}\end{array}\right.,$

步骤(5.5)：输入第t-1次迭代，即上一次迭代的获得的参数以及步骤(5.3)中第t次迭代得到的参数按下式计算和

$\left\{\begin{array}{c}{\Σ}_u^{\left(t\right)}={({\stackrel{\‾}{\mathrm{\α}}}_n^{(t-1)}{\[\stackrel{\‾}{{\mathrm{UU}}^{\′}}\]}^{\left(t\right)}+{\stackrel{\‾}{\mathrm{\α}}}_w^{(t-1)}{I}_r)}^{-1}\\ {\stackrel{\‾}{W}}^{\left(t\right)}={\stackrel{\‾}{\mathrm{\α}}}_n^{(t-1)}{\Σ}_w^{\left(t\right)}{\stackrel{\‾}{U}}^{\left(t\right)}\tilde{X}\\ {\[\stackrel{\‾}{{\mathrm{WW}}^{\′}}\]}^{\left(t\right)}={\stackrel{\‾}{W}}^{\left(t\right)}{\left({\stackrel{\‾}{W}}^{\left(t\right)}\right)}^{\′}+h{\Σ}_w^{\left(t\right)},\end{array}\right.,$

步骤(5.6)：输入步骤(5.3)～步骤(5.5)所计算获得的所有参数，按下式计算和

$C_1=tr\{{\[\stackrel{\‾}{U{\tilde{F}}^{\′}\tilde{F}{U}^{\′}}\]}^{\left(t\right)}{\[\stackrel{\‾}{{\mathrm{VV}}^{\′}}\]}^{\left(t\right)}+{\[\stackrel{\‾}{{\mathrm{UU}}^{\′}}\]}^{\left(t\right)}{\[\stackrel{\‾}{{\mathrm{WW}}^{\′}}\]}^{\left(t\right)}\},$

$C_2=tr\{{\stackrel{\‾}{U}}^{\left(t\right)}{\left({\stackrel{\‾}{U}}^{\left(t\right)}\right)}^{\′}{\stackrel{\‾}{W}}^{\left(t\right)}{\left({\stackrel{\‾}{W}}^{\left(t\right)}\right)}^{\′}+{\stackrel{\‾}{U}}^{\left(t\right)}{\tilde{F}}^{\′}\tilde{F}{\left({\stackrel{\‾}{U}}^{\left(t\right)}\right)}^{\′}{\stackrel{\‾}{V}}^{\left(t\right)}{\left({\stackrel{\‾}{V}}^{\left(t\right)}\right)}^{\′}\},$

${\stackrel{\‾}{\mathrm{\α}}}_n^{\left(t\right)}=\frac{2\×{10}^{-6}+NL+Nl}{2\×{10}^{-6}+\left|\right|\tilde{X}-{\left({\stackrel{\‾}{U}}^{\left(t\right)}\right)}^{\′}{\stackrel{\‾}{W}}^{\left(t\right)}|{|}_F^2+||\tilde{E}-\tilde{F}{\left({\stackrel{\‾}{U}}^{\left(t\right)}\right)}^{\′}{\stackrel{\‾}{V}}^{\left(t\right)}|{|}_F^2+C_1-C_2},$

${\stackrel{\‾}{\mathrm{\α}}}_u^{\left(t\right)}=\frac{2\×{10}^{-6}+Lr}{2\×{10}^{-6}+tr\left\{{\[\stackrel{\‾}{{\mathrm{UU}}^{\′}}\]}^{\left(t\right)}\right\}},{\stackrel{\‾}{\mathrm{\α}}}_v^{\left(t\right)}=\frac{2\×{10}^{-6}+Nr}{2\×{10}^{-6}+tr\left\{{\[\stackrel{\‾}{{\mathrm{VV}}^{\′}}\]}^{\left(t\right)}\right\}},{\stackrel{\‾}{\mathrm{\α}}}_w^{\left(t\right)}=\frac{2\×{10}^{-6}+Nr}{2\×{10}^{-6}+tr\left\{{\[\stackrel{\‾}{{\mathrm{WW}}^{\′}}\]}^{\left(t\right)}\right\}},$

其中，||·||_F表示矩阵的Frobenius范数，tr{·}表示矩阵的迹运算，N、L和r的定义见步骤(1)；

步骤(5.7)：根据步骤(5.2)～步骤(5.6)的顺序，重复迭代Δ次，得到和再按下式计算输出锐化后的高分辨率的高光谱图像