基于狄利克雷MRF混合模型的极化SAR图像分类方法与流程

本发明属于图像处理
技术领域：
，特别涉及一种极化SAR图像的分类方法，可用于极化SAR图像的目标检测与识别。
背景技术：
：极化合成孔径雷达SAR是一种高分辨率成像雷达。其在民用和军事领域的广泛使用需要极化SAR图像解译技术作为支撑，而极化SAR图像分类是机器学习和数据挖掘领域的一项重要技术，也是图像解译的重要内容之一，它可以提供极化SAR图像的整体结构信息，揭示极化SAR图像的本质。近年来，极化SAR图像分类方法一直是该领域前沿研究的热点，其中随机场模型被认为是处理极化SAR图像分类问题的强有力工具，其优势在于能够在极化SAR图像分类过程中充分考虑各像素点间的空域相关性，从而提高算法的抗噪性，目前，基于随机场模型的极化SAR图像分类方法应用较多的是基于马尔科夫随机场MRF的极化SAR图像分类方法。WanyingSong等人发表的论文“TheWGΓDistributionforMulti-LookPolarimetricSARDataandItsApplication”(IEEEGeosci.RemoteSens.Lett.,2014,12(10):2056-2060)提出了将极化SAR数据的乘性模型引入MRF框架中实现极化SAR图像分类的方法。该论文利用极化SAR数据的乘性模型构建MRF框架中的似然概率，能够更精确地描述极化SAR图像的统计特征，从而克服了Wishart分布不具备描述极化SAR图像中复杂场景统计特征的能力。但是该方法的不足是，此类MRF模型需要训练数据学习模型参数，无法在缺乏极化SAR图像先验信息的前提下精确地确定极化SAR图像中的类别数目，只限于处理监督的极化SAR图像分类。PeterOrbanz等人发表的论文“NonparametricBayesianimagesegmentation”(Int.J.Comput.Vision,2008,77:25-45)提出了利用狄利克雷MRF混合模型实现SAR图像的直方图聚类算法。该论文有效地结合了MRF模型和狄利克雷过程混合模型DPMM，并应用到SAR图像的直方图聚类算法中，视SAR图像中的类别数目为模型中的一个随机变量，根据模型和数据不断地对该随机变量进行更新，从而使得SAR图像的直方图聚类算法能够自动地确定SAR图像中的类别数目；且该论文中SAR图像空域结构信息的引入大大提高了模型的抗噪性能。但是该方法的不足是，此类模型是通过SAR图像特征向量的欧式距离来引入空域结构信息的，对极化SAR数据并不适用。技术实现要素：本发明的目的在于克服上述已有问题的缺点，提出一种基于狄利克雷MRF混合模型的极化SAR图像分类方法，以在缺乏图像先验信息的前提下自动地确定极化SAR图像中的类别数目，并实现利用狄利克雷MRF混合模型对极化SAR图像的分类。为达到上述目的，本发明的技术方案包括如下：(1)输入极化SAR图像；(2)从极化SAR图像中提取N个极化散射特征，并对其进行归一化，得到N个归一化的极化散射特征，建立归一化的极化散射特征空间F1；(3)以极化SAR图像中的每个像素点为中心像素点，选取大小为3×3的窗口，使用广义均值方法，对归一化的极化散射特征空间F1逐点进行降噪处理，得到极化散射特征空间F2；(4)初始化模型参数将极化SAR图像的每个像素点的类别标记的集合设定为原始的标记场X＝{Xi}i＝1,2,…,n，其中Xi表示极化SAR图像中像素点i的类别标记，n表示极化SAR图像的像素点个数；初始化极化SAR图像的每个像素点的类别标记为1，即Xi＝1；用极化散射特征空间F2作为观测场Y＝{Yi}i＝1,2,…,n，其中Yi表示极化SAR图像中像素点i的N维降噪的极化散射特征向量；初始化形状参数a＝1，尺度参数b＝5，核宽参数σ＝1，平衡参数λ＝1，维度参数q＝3，第一正态参数μ0初始化为N维零向量，第二正态参数Σ0初始化为N维单位阵；初始化中心参数α0，即对Flevoland地区的极化SAR图像设定中心参数α0＝2；对Oberpfaffenhofen地区的极化SAR图像设定中心参数α0＝10；将当前的迭代次数设为1，设定最大迭代次数为20；(5)根据原始的标记场X、观测场Y和初始化的模型参数，估计极化SAR图像的先验参数sij和似然参数{μk,Σ1}；(6)估计极化SAR图像新的标记场6a)根据先验参数sij、原始的标记场X和初始化的模型参数，得到类别标记的先验概率p3：p3=Σl=1,l≠inδ(Xi,Xl)α0+n-1exp(-λΣj∈∂(i)sij(1-2δ(Xi,Xj))),Xi∈X-iα0α0+n-1,Xi∉X-i,]]>其中Xj表示极化SAR图像中像素点j的类别标记；Xl表示极化SAR图像中像素点l的类别标记；表示极化SAR图像中像素点i的邻域系统；X-i表示极化SAR图像中除了像素点i的类别标记以外的其他像素点的类别标记的集合；δ(Xi,Xl)表示Delta函数，当Xi与Xl相等时，δ(Xi,Xl)的值为1，否则为0；δ(Xi,Xj)表示Delta函数，当Xi与Xj相等时，δ(Xi,Xj)的值为1，否则为0；sij表示极化SAR图像的先验参数；exp(·)表示求指数操作；6b)根据似然参数{μk,Σ1}、原始的标记场X和观测场Y，得到类别标记的似然概率p4：p4=1(2π)N/2|Σ1|1/2exp(-12(Yi-μk)TΣ1-1(Yi-μk)),Xi=k,]]>其中μk表示极化SAR图像中类别标记为k的像素点的均值向量参数；Σ1表示极化SAR图像的协方差矩阵参数；(·)-1表示求矩阵逆的操作；(·)T表示求矩阵转置的操作；|·|表示求矩阵行列式的操作；Xi＝k表示极化SAR图像中像素点i的类别标记为k；6c)根据类别标记的先验概率p3和似然概率p4，得到类别标记的后验概率p5；p5=z-1p4·nk-i·exp(-λΣj∈∂(i)sij(1-2δ(Xi,Xj))),k=1,2,...,Cz-1α01(2π)N|Σ1Σ0|1/2exp(-12(Yi-μ0)T(Σ1Σ0)-1(Yi-μ0)),k=C+1,]]>其中指的是极化SAR图像中除了像素点i以外的类别标记为k的像素点的个数；z为归一化参数；C表示原始的标记场X中的类别数目；6d)通过极大化类别标记的后验概率p5，估计得到极化SAR图像中像素点i的新的类别标记X‾i=maxk=1,2,...,C,C+1{p5};]]>6e)将新的类别标记的集合确定为极化SAR图像的新的标记场(7)判断是否达到最大迭代次数，输出分类结果：将当前的迭代次数加1作为新的迭代次数，判定新的迭代次数是否达到步骤(4)所设定的最大迭代次数20，如果没有达到，则设定新的标记场为原始的标记场X，返回步骤(5)；如果达到，则将新的标记场确定为极化SAR图像的最终的分类结果并输出。本发明与现有技术相比，具有以下优点：1)相比于基于MRF的极化SAR图像的监督分类方法，即WanyingSong等人提出的方法，本发明方法将狄利克雷过程混合模型DPMM引入到MRF框架中，可以克服在缺乏极化SAR图像先验信息的情况下难以确定极化SAR图像中的类别数目这一弊端，因此本发明更适合处理复杂场景的极化SAR图像的分类。2)相比于PeterOrbanz等人提出的方法，本发明利用极化协方差矩阵之间的相似度估计极化SAR图像的先验参数，能增强模型与极化SAR图像的数据相关性，能使同质区域分类结果更平滑且能更好地保持边缘信息，进而获得更好的分类结果，因此本发明可以显著地提高极化SAR图像分类结果的精度。附图说明图1是本发明的实现流程图；图2是用本发明和现有方法对Flevoland地区的极化SAR图像的分类对比图；图3是用本发明和现有方法对Oberpfaffenhofen地区的极化SAR图像的分类对比图。具体实施方式下面结合附图对本发明做进一步的描述。参照图1，本发明的具体实施步骤如下：步骤1，输入极化SAR图像。本发明选取如下两幅极化SAR图像：Flevoland地区的极化SAR图像：图像大小为320×326像素；等效视数为4；分辨率为2.5652m×2m；雷达系统为AIRSAR；Oberpfaffenhofen地区的极化SAR图像：图像大小为500×450像素；等效视数为2；分辨率为3m×0.89m；雷达系统为ESAR。步骤2，提取并归一化极化散射特征，建立归一化的极化散射特征空间F1。2a)从极化SAR图像中提取N个极化散射特征Fr，r＝1,2,...,N，N＝17表示极化散射特征的个数，本实例所选取的N个极化散射特征如表1所示：表1通过极化SAR数据实验验证了层次后向策略应用于极化目标分解后，所提取的极化散射特征有助于提高分类精度，因此本发明提取用结合层次后向策略的Yamaguchi-Sato分解方法对极化SAR图像进行分解所得到的螺旋散射功率、体散射功率、二面角散射功率和表面散射功率作为极化散射特征F1、F2、F3、F4；极化SAR图像的相干矩阵T有利于解译目标的散射机制，例如相干矩阵T的对角元素分别表示Pauli分解的系数，于是极化SAR图像的相干矩阵T的各个元素T11、T12、T13、T22、T23、T33有助于提高分类精度，因此本实例提取T11、|T12|、|T13|、T22、|T23|、T33作为极化散射特征F5、F6、F7、F8、F9、F10；极化SAR图像的相干矩阵T的特征值与极化总功率同样是极化SAR图像分类的重要特征，因此需要提取极化SAR图像的相干矩阵T的3个特征值与极化总功率作为极化散射特征F11、F12、F13、F14；用Cloude分解方法对极化SAR图像进行分解得到的3个极化散射特征：极化熵、平均散射角和反熵是最近十年内极化SAR图像分类中最常用的3个极化散射特征，其中极化熵对植被、树林、非常粗糙地面等目标是非常敏感的，平均散射角用于确定目标的主导散射，反熵用于区分植被与裸地等散射机制差异较小的类别，因此本发明提取极化熵、平均散射角和反熵作为极化散射特征F15、F16、F17；2b)通过下式分别对极化SAR图像的N个极化散射特征Fr，r＝1,2,...,N，进行归一化，得到N个归一化的极化散射特征Fr‾=Fr-min(Fr)max(Fr)-min(Fr);]]>2c)将归一化的极化散射特征中像素点i的特征值的集合确定为极化SAR图像中像素点i的归一化的极化散射特征向量Fi，i＝1,2,...,n：Fi={Fr‾(i)}r=1,2,...,N;]]>2d)将归一化的极化散射特征向量Fi的集合确定为归一化的极化散射特征空间F1：F1＝{Fi}i＝1,2,...,n。步骤3，使用广义均值方法对归一化的极化散射特征空间F1逐点进行降噪处理，得到极化散射特征空间F2。3a)以极化SAR图像中的像素点i为中心像素点，选取大小为3×3的窗口N3×3(i)，使用广义均值方法，根据归一化的极化散射特征在3×3的窗口N3×3(i)内提取极化SAR图像中像素点i的第r个降噪的极化散射特征Gr(i)：Gr(i)=(Πt∈N3×3(i)Fr‾(t))19,]]>其中表示第r个归一化的极化散射特征中像素点t的特征值；N3×3(i)表示以像素点i为中心像素点的3×3的窗口；3b)将N个降噪的极化散射特征Gr(i)的集合确定为极化SAR图像中像素点i的N维降噪的极化散射特征向量Gi：Gi＝{Gr(i)}r＝1,2,...,N；3c)将N维降噪的极化散射特征向量Gi的集合确定为极化散射特征空间F2：F2＝{Gi}i＝1,2,...,n。步骤4，初始化模型参数。4a)初始化MRF模型(X,Y)：将极化SAR图像的每个像素点的类别标记的集合设定为原始的标记场X＝{Xi}i＝1,2,…,n，其中Xi表示极化SAR图像中像素点i的类别标记；初始化极化SAR图像的每个像素点的类别标记为1，即Xi＝1；用极化散射特征空间F2作为观测场Y＝{Yi}i＝1,2,…,n，其中Yi表示极化SAR图像中像素点i的N维降噪的极化散射特征向量；4b)初始化形状参数a＝1，尺度参数b＝5，核宽参数σ＝1，平衡参数λ＝1，维度参数q＝3，第一正态参数μ0初始化为N维零向量，第二正态参数Σ0初始化为N维单位阵；4c)初始化中心参数α0，即对Flevoland地区的极化SAR图像设定中心参数α0＝2；对Oberpfaffenhofen地区的极化SAR图像设定中心参数α0＝10；4d)将当前的迭代次数设为1，设定最大迭代次数为20。步骤5，估计极化SAR图像的先验参数sij和似然参数{μk,Σ1}。先验参数sij的估计方法一般有Derin方法、极化散射特征向量的欧氏距离方法等，但是Derin方法运算量较大，欧氏距离方法不能精确地描述极化SAR数据，因此本发明利用极化协方差矩阵的相似度估计极化SAR图像的先验参数sij，既能降低算法运算量，又能增强与极化SAR图像的数据相关性；似然参数{μk,Σ1}常采用最大似然方法进行估计，但是最大似然方法不能体现狄利克雷模型的随机性优势，因此本发明利用Gibbs采样算法对似然参数{μk,Σ1}的后验分布进行采样，即可估计得到似然参数{μk,Σ1}。5a)从极化SAR图像中选取像素点i的极化协方差矩阵Ci与邻域像素点的极化协方差矩阵Cj，根据初始化的模型参数，计算Ci与Cj之间的对称的改进型Wishart距离d(Ci,Cj)：d(Ci,Cj)=12(Tr(Cj-1Ci-q)+Tr(Ci-1Cj-q)),]]>其中Tr(·)表示求矩阵迹的操作；5b)根据对称的改进型Wishart距离d(Ci,Cj)和初始化的模型参数，计算极化SAR图像的先验参数sij：sij=exp(-d2(Ci,Cj)2σ2),j∈∂(i);]]>5c)根据原始的标记场X、观测场Y和初始化的模型参数，得到均值向量参数μk的后验概率p1：p1∝ΠXi=kN(Yi|μk,Σ1)N(μk|μ0,Σ0),]]>其中N(·|·,·)表示正态分布；Π表示求各项乘积的操作；5d)根据均值向量参数μk的后验概率p1，利用Gibbs采样算法，估计得到极化SAR图像的均值向量参数μk；5e)根据原始的标记场X、观测场Y和初始化的模型参数，得到方差参数的后验概率p2：p2∝Γ(a+n2,b+12Σi=1,Xi=kn(Yir-μkr)2),]]>其中Γ(·,·)表示伽马分布；Yir表示极化SAR图像中像素点i的N维降噪的极化散射特征向量Yi的第r个分量；μkr表示极化SAR图像中类别标记为k的像素点的均值向量参数μk的第r个分量；r＝1,2,…,N；5f)根据方差参数的后验概率p2，利用Gibbs采样算法，估计极化SAR图像的方差参数即可得到协方差矩阵参数5g)由估计得到的均值向量参数μk和协方差矩阵参数Σ1，得到似然参数{μk,Σ1}。步骤6，估计极化SAR图像新的标记场马尔科夫随机场MRF模型在分类问题中根据极大化类别标记的后验概率不断地估计新的标记场直到达到最大迭代次数。因此本发明首先根据先验参数sij和似然参数{μk,Σ1}得到类别标记的先验概率p3和似然概率p4，然后根据先验概率p3和似然概率p4得到后验概率p5，通过极大化类别标记的后验概率p5即可估计得到极化SAR图像新的标记场6a)根据先验参数sij、原始的标记场X和初始化的模型参数，得到类别标记的先验概率p3：p3=Σl=1,l≠inδ(Xi,Xl)α0+n-1exp(-λΣj∈∂(i)sij(1-2δ(Xi,Xj))),Xi∈X-iα0α0+n-1,Xi∉X-i,]]>其中Xj表示极化SAR图像中像素点j的类别标记；Xl表示极化SAR图像中像素点l的类别标记；表示极化SAR图像中像素点i的邻域系统；X-i表示极化SAR图像中除了像素点i的类别标记以外的其他像素点的类别标记的集合；δ(Xi,Xl)表示Delta函数，当Xi与Xl相等时，δ(Xi,Xl)的值为1，否则为0；δ(Xi,Xj)表示Delta函数，当Xi与Xj相等时，δ(Xi,Xj)的值为1，否则为0；sij表示极化SAR图像的先验参数；exp(·)表示求指数操作；6b)根据似然参数{μk,Σ1}、原始的标记场X和观测场Y，得到类别标记的似然概率p4：p4=1(2π)N/2|Σ1|1/2exp(-12(Yi-μk)TΣ1-1(Yi-μk)),Xi=k,]]>其中μk表示极化SAR图像中类别标记为k的像素点的均值向量参数；Σ1表示极化SAR图像的方差矩阵参数；(·)-1表示求矩阵逆的操作；(·)T表示求矩阵转置的操作；|·|表示求矩阵行列式的操作；Xi＝k表示极化SAR图像中像素点i的类别标记为k；6c)基于Baysian准则，根据类别标记的先验概率p3和似然概率p4，得到类别标记的后验概率p5；p5=z-1p4·nk-i·exp(-λΣj∈∂(i)sij(1-2δ(Xi,Xj))),k=1,2,...,Cz-1α01(2π)N|Σ1Σ0|1/2exp(-12(Yi-μ0)T(Σ1Σ0)-1(Yi-μ0)),k=C+1,]]>其中指的是极化SAR图像中除了像素点i以外的类别标记为k的像素点的个数；z为归一化参数；C表示原始的标记场X中的类别数目；6d)通过极大化类别标记的后验概率p5，估计得到极化SAR图像中像素点i的新的类别标记X‾i=maxk=1,2,...,C,C+1{p5};]]>6e)将新的类别标记的集合确定为极化SAR图像的新的标记场(7)判断是否达到最大迭代次数，输出分类结果：将当前的迭代次数加1作为新的迭代次数，判定新的迭代次数是否达到步骤(4)所设定的最大迭代次数20，如果没有达到，则设定新的标记场为原始的标记场X，返回步骤(5)；如果达到，则将新的标记场确定为极化SAR图像的最终的分类结果并输出。下面结合仿真实验对本发明的效果做进一步说明：1)实验条件实验仿真环境为：MATLABR2010b,Intel(R)Core(TM)i7-4770CPU3.4GHz，Window7旗舰版。2)实验内容实验1、分别应用基于MRF的极化SAR图像的监督分类方法，即WanyingSong等人提出的方法、PeterOrbanz等人提出的方法及本发明方法对Flevoland地区的极化SAR图像进行分类，结果如图2所示，其中：图2(a)是Flevoland地区的极化SAR图像的Pauli分解图，图2(b)是WanyingSong等人提出的方法得到的分类结果，图2(c)是PeterOrbanz等人提出的方法得到的分类结果，图2(d)是本发明方法得到的分类结果，图2(e)是Flevoland地区的极化SAR图像的分类结果参考图，图2(f)是Flevoland地区的极化SAR图像各个类别对应的地物类型。从图2可以看到，本发明在抗噪性能和边界信息保持上均优于WanyingSong等人提出的方法和PeterOrbanz等人提出的方法，而且相比于PeterOrbanz等人提出的方法，本发明得到的分类结果中误分类数目明显减少。为了验证分类结果的质量，以分类精度为性能指标参数，评价实验1中本发明和现有方法得到的分类结果，以图2(e)为分类结果参考图。本发明和现有方法得到的分类结果的分类精度(％)如表2所示。表2从表2可以看出，本发明方法得到的分类结果中各个类别的分类精度和总精度都比较高，虽然WanyingSong等人提出的方法也能获得较高的分类精度，但是该方法最大的弊端是不能自动地确定极化SAR图像中的类别数目，而PeterOrbanz等人提出的方法由于缺乏与极化SAR数据的相关性，误分类数目较多，因此分类精度较低。实验2、分别应用基于MRF的极化SAR图像的监督分类方法，即WanyingSong等人提出的方法、PeterOrbanz等人提出的方法及本发明方法对Oberpfaffenhofen地区的极化SAR图像进行分类，结果如图3所示，其中：图3(a)是Oberpfaffenhofen地区的极化SAR图像的Pauli分解图，图3(b)是WanyingSong等人提出的方法得到的分类结果，图3(c)是PeterOrbanz等人提出的方法得到的分类结果，图3(d)是本发明方法得到的分类结果。从图3可以看到，相比于WanyingSong等人提出的方法和PeterOrbanz等人提出的方法，本发明方法得到的分类结果中受噪声影响较小，而且细节信息的保持更加完整，最重要的是能够对相似的地物目标进行区分，误分类区域较少。综上，本发明能使同质区域分类结果更平滑且能更好地保持边缘信息，分类精度明显提高。当前第1页1 2 3