CN101908138A - 基于噪声独立成分分析的合成孔径雷达图像目标识别方法 - Google Patents

Abstract

本发明涉及一种基于噪声独立成分分析的合成孔径雷达图像目标识别方法，属于合成孔径雷达图像目标识别的技术领域。首先对对输入的合成孔径雷达的原始训练样本图像进行预处理，使图像符合白化和零均值化，图像噪声概率密度服从高斯分布；再对预处理之后的训练样本以及将待识别的实时测量样本进行对数噪声独立成分分析，提取出待识别图像的独立成分特征；最后利用独立成分分析法对待识别样本进行识别。本发明方法提取出的独立成分特征更适于分类；省去了理想ICA算法中在预处理阶段必须进行的图像去噪过程，提高了算法对异常数据的鲁棒性，提高了SAR图像目标识别的实时性与可靠性，可用于对SAR图像的特征提取及自动目标识别。

Description

基于噪声独立成分分析的合成孔径雷达图像目标识别方法

技术领域

本发明涉及一种基于噪声独立成分分析的合成孔径雷达图像目标识别方法，属于合成孔径雷达图像目标识别的技术领域。

背景技术

合成孔径雷达(Synthetic Aperture Radar，以下简称SAR)根据目标的后向电磁散射形成SAR图像，能够克服光学成像受距离、气候等条件限制的缺点，在航空测量、航空遥感、卫星海洋观测、航天侦察、图像匹配制导等领域发挥着至关重要的作用。但在SAR成像的过程中，同一分辨率单元内通常分布着较多的散射体，这些散射体回波信号的相位随意分布，相互间的干涉也就导致了SAR图像固有的乘性斑点噪声，使一般的SAR图像与普通的电子光学图像相比，具有分辨率低、对操作条件敏感、高乘性噪声等特点。

SAR图像目标识别技术是利用单个目标或目标群在雷达远区所产生的散射场的特征(如回波信号的幅值、相位、频率和极化信息等)，确定出目标的形状、类型、位置等属性的技术，在军事和民用领域得到了广泛的应用；尤其是在复杂战场环境下，准确、快捷的SAR图像目标识别技术能为提高战场感知能力和反应能力提供强有力地保障。

而在SAR图像目标识别中，SAR图像的噪声往往还包括了各种自然杂波(草地、河流、森林等)和人造杂波(建筑物等)的干扰，并且这些噪声通常都是非高斯性的。因此，针对如何有效克服SAR图像中存在的大量噪声对目标识别的影响，尽可能地提高SAR图像目标识别的正确率，人们提出了很多SAR图像噪声去除以及特征提取的方法，应用较广的如小波滤波，方向扩散滤波，基于马尔可夫随机场的目标检测方法，基于多元统计分析与核准则的主成分分析(Principle Component Analysis：以下简称PCA)、核PCA(Kernel PCA：KPCA)、Fisher线性判别分析(FLDA)、核Fisher判别分析(KFDA)、独立成分分析(Independent Component Analysis：以下简称ICA)等等，涉及到图像处理和模式识别的诸多领域。

在SAR图像目标识别***中，希望通过特征提取达到如下目的：1)使得处于不同方位、差异很大的同类SAR图像具有相同的描述方式；2)消除或降低SAR图像的高噪声特点对分类结果的影响；3)压缩SAR图像庞大的数据量；进而将有效且易于识别的二次特征提供给分类器，由分类器判别出目标的结构、类型等属性，完成SAR图像目标识别的过程。而在针对SAR图像特征提取的方法当中，应用广泛且较为有效的是基于多元统计分析的特征提取方法，主要涉及到PCA和ICA两类。其中，PCA是在最小均方误差的意义上，寻找一组可以最大化样本方差的正交投影方向来实现去除冗余信息、保留含最大能量的成 分成分的特征提取方法，这种方法实际提取的是数据的二阶统计量即方差特征。而ICA作为PCA的推广，提取的是数据的高阶统计量特征，即从观测数据中寻求一组在高阶意义下独立的源数据，由这组源数据线性叠加成观测数据。

SAR图像的很多信息包含于图像的高阶统计量之中，而噪声又通常是非高斯的，因此对于SAR图像而言，ICA特征提取的方法比PCA更具鲁棒性，而包含高阶统计信息的ICA特征也更具可分性。

目前一般的ICA技术通常都将观测数据理想化为零噪声，但这往往与实际情况有较大出入，特别是在具有高乘性非高斯噪声的SAR图像目标识别领域，如何消除或减少噪声对ICA特征的影响是一个亟待解决的问题。目前处理噪声的方法一般可分为三种：①把噪声看作一个独立成分；②在预处理阶段进行去噪；③噪声ICA模型。几乎所有的无噪声ICA模型都默认采用了第一种方法，即把噪声看作原始独立成分中的一个，但在SAR图像目标识别中，较大的噪声成分将削弱正常成分对分类的贡献作用，导致分类精度降低。而预处理阶段的图像去噪能够使去噪过程更具针对性，从而使得ICA提取的特征达到更高的分类精度。例如，上海交通大学的专利ZL 200710046928.3(基于独立成分分析基图像的合成孔径雷达图像消噪方法)对SAR图像利用小波指数对图像进行平滑增强，从而达到去噪效果；美国专利US 7508334(Method and apparatus for processing SAR images based onan anisotropic diffusion filtering algorithm，March 24，2009)提出了针对SAR图像斑点噪声的方向扩散滤波算法来实现SAR图像去噪。然而在预处理阶段进行图像去噪自然会增加不少的计算量，导致图像特征提取的时间过长，因此无法满足强调实时性的SAR图像目标识别的战场需求。

噪声ICA模型在前述理想的无噪声ICA模型基础上，将噪声干扰作为一个不可忽略的影响因素，在此基础上，芬兰赫尔辛基大学的Aapo Hyvarinen等人提出了一种结合噪声移除技术的ICA算法——噪声快速ICA(Noisy FastICA：NFastICA)算法。这种算法引入了高斯矩的概念，成功地从被高斯噪声污染了的观测数据中直接估计出了潜在的随机变量，并且对采样的异常值具有较好的鲁棒性。然而，NFastICA算法的两个重要前提是数据噪声为高斯分布(也就是正态分布)并且噪声协方差矩阵是已知的，但SAR图像中的噪声多为乘性噪声，且除了与待识别目标相关的区域中的非目标回波之外，还包括了自然杂波和人造杂波，几乎不可能用高斯分布加以拟合。因此，NFastICA算法并不能直接应用到SAR图像目标识别领域当中。

发明内容

本发明的目的是提出一种基于噪声独立成分分析的合成孔径雷达图像目标识别方法，以噪声ICA模型为基础，对现有的基于ICA的SAR图像目标识别方法进行改进，以减缓或消除使用ICA方法对SAR图像进行分类识别时大量非高斯乘性噪声对特征提取的影响，并且减少预处理阶段的计算量，提高识别算法的鲁棒性，最终达到提高SAR图像目标识别技 术的识别正确率和识别效率的目的。

本发明提出的基于噪声独立成分分析的合成孔径雷达图像目标识别方法，包括以下步骤：

(1)对输入的合成孔径雷达的原始训练样本图像X_train进行预处理，具体过程如下：

(1-1)对训练样本图像X_train进行对数变换，得到X_ln＝20ln(1+X_Orig)，其中X_Orig表示输入的原始训练样本单幅图像，X_ln表示经对数变换后噪声分布符合高斯分布的图像；

(1-2)对上述X_ln进行去均值化处理，得到零均值图像X_t，X_t＝X_ln-E(X_ln)，其中期望函数E表示对图像求均值；

(1-3)将上述零均值图像X_t划分为目标区域X_o和噪声区域n_o，并满足{X_t}＝{n_o}∪{X_o}，将目标区域和噪声区域的二维图像数据分别拉直成一维行向量，其中的目标区域X_o包含与被识别目标相关的目标、阴影和杂波的回波，噪声区域n_o为背景杂波；

(1-4)重复步骤(1-1)、(1-2)和(1-3)，得到训练样本中的N幅训练图像的一维目标区域数据X_o，构成N行训练矩阵X_O；

(1-5)根据上述得到的N幅训练图像的所有噪声区域n_o，构造训练样本的噪声协方差对角矩阵 

其中的对角元素 

为每幅训练图像的噪声方差估计值， 

表示第i幅训练图像的噪声区域；

(1-6)根据上述训练样本的噪声协方差矩阵∑_O，对上述训练矩阵X_O进行降维和白化处理，得到降维和白化后的合成孔径雷达训练样本图像子空间矩阵X；

(2)获取合成孔径雷达实时测量图像X_test，使用独立成分分析法分别提取上述合成孔径雷达训练样本图像X_train的独立成分特征和合成孔径雷达实时测量图像X_test的独立成分特征，具体过程如下：

(2-1)采用噪声快速独立成分分析法处理上述合成孔径雷达训练样本图像子空间矩阵X，得到一组由基图像估计向量组成的基图像估计矩阵S_e，表示为：

X＝c_train1·s₁+c_train2·s₂+...+c_trainm·s_L＝S_e ^T·c_train.

其中，s₁，s₂，…，s_L表示基图像估计矩阵S_e中L个基图像的估计列向量，c_train＝(c_train1，c_train2，...，c_trainL)^T表示上述训练样本图像子空间矩阵X在由上述基图像估计矩阵S_e构成的基图像子空间中的投影系数，c_train即为合成孔径雷达训练样本X_train的独立成分特征；

(2-2)将合成孔径雷达实时测量图像X_test中待识别的图像数据投影到由上述基图像估计矩阵S_e构成的基图像子空间中，使实时测量图像X_test用基图像估计向量的线性组合表示，为：

X_test＝c_test1·s₁+c_test2·s₂+...+c_testL·s_L＝S_e ^T·c_test

其中，c_test＝(c_test1，c_test2，...，c_testL)^T表示上述实时测量图像X_test在由上述基图像估计矩阵S_e构成的基图像子空间中的投影系数，c_test＝piv(S_e ^T)·X_test，其中的piv(S_e ^T)表示矩阵S_e转置的广义逆，c_test即为待识别的合成孔径雷达实时测量图像X_test的独立成分特征；

(3)根据上述步骤(2)得到的合成孔径雷达训练样本图像X_train的独立成分特征c_train和待识别的合成孔径雷达实时测量图像X_test的独立成分特征c_test，对合成孔径雷达实时测量图像X_test进行识别分类，判断出被测目标的类别。

本发明提出的基于噪声独立成分分析的合成孔径雷达图像目标识别方法，其优点是：

1、本发明提出的基于噪声独立成分分析的合成孔径雷达图像目标识别方法，涉及新的对数正态噪声独立成分分析(Log-normal noise ICA：LnnICA)方法，省掉了已有的无噪声独立成分分析算法在预处理阶段必须进行的图像去噪过程，有效克服了合成孔径雷达图像中大量乘性非高斯噪声对特征提取的影响，使得提取出的独立成分特征更适于分类，从而减少了计算量，提高了合成孔径雷达图像目标识别的可靠性以及对异常数据的鲁棒 性，能较大程度地改善合成孔径雷达目标识别尤其是实时自动目标识别的识别正确率和识别效率。

2、本发明的合成孔径雷达图像目标识别方法，在降维和白化的过程中，改进了PCA算法中主成分特征的选取标准，从降噪的角度提出了最可分成分(Most DiscriminableComponent，以下简称MDC)选取准则，即利用类别可分性来衡量某个特征成分能够表示原模式类别属性的程度，并从中选取可分性最强的多个MDC成分来代表原始图像。采用特征成分的可分性原则在预处理阶段进行降维和白化，可以同时改善以分类为目的的SAR图像去噪效果，提高基于独立成分分析法的SAR图像目标识别的正确率和可靠性。

附图说明

图1为已有技术中基于ICA的SAR图像目标识别技术的一般流程。

图2为本发明方法中，对SAR图像进行对数变换前后噪声分布特性对比。

图3为本发明方法中从单幅训练样本图像提取噪声区域和目标区域的示意图。

图4为本发明提出的基于噪声独立成分分析的合成孔径雷达图像目标识别方法流程。

具体实施方案

本发明提出的基于噪声独立成分分析的合成孔径雷达图像目标识别方法，包括以下步骤：

(1)对输入的合成孔径雷达的原始训练样本图像X_train进行预处理，具体过程如下：

(1-1)对训练样本图像X_train进行对数变换，得到X_ln＝20ln(1+X_Orig)，其中X_Orig表示输入的原始训练样本单幅图像，X_ln表示经对数变换后噪声分布符合高斯分布的图像；

(1-2)对上述X_ln进行去均值化处理，得到零均值图像X_t，X_t＝X_ln-E(X_ln)，其中期望函数E表示对图像求均值；

(1-3)将上述零均值图像X_t划分为目标区域X_o和噪声区域n_o，并满足{X_t}＝{n_o}∪{X_o}，将目标区域和噪声区域的二维图像数据分别拉直成一维行向量，其中的目标区域X_o包含与被识别目标相关的目标、阴影和杂波的回波，噪声区域n_o为背景杂波；

(1-4)重复步骤(1-1)、(1-2)和(1-3)，得到训练样本中的N幅训练图像的 一维目标区域数据X_o，构成N行训练矩阵X_O；

(1-5)根据上述得到的N幅训练图像的所有噪声区域n_o，构造训练样本的噪声协

方差对角矩阵 

其中的对角元素 

为每幅训练图像的噪声方差估计值， 

表示第i幅训练图像的噪声区域；

(1-6)根据上述训练样本的噪声协方差矩阵∑_O，对上述训练矩阵X_O进行降维和白化处理，得到降维和白化后的合成孔径雷达训练样本图像子空间矩阵X；

(2)获取合成孔径雷达实时测量图像X_test，使用独立成分分析法分别提取上述合成孔径雷达训练样本图像X_train的独立成分特征和合成孔径雷达实时测量图像X_test的独立成分特征，具体过程如下：

(2-1)采用噪声快速独立成分分析法处理上述合成孔径雷达训练样本图像子空间矩阵X，得到一组由基图像估计向量组成的基图像估计矩阵S_e，表示为：

X＝c_train1·s₁+c_train2·s₂+...+c_trainm·s_L＝S_e ^T·c_train.

其中，s₁，s₂，…，s_L表示基图像估计矩阵S_e中L个基图像的估计列向量，c_train＝(c_train1，c_train2，...，c_trainL)^T表示上述训练样本图像子空间矩阵X在由上述基图像估计矩阵S_e构成的基图像子空间中的投影系数，c_train即为合成孔径雷达训练样本X_train的独立成分特征；

(2-2)将合成孔径雷达实时测量图像X_test中待识别的图像数据投影到由上述基图像估计矩阵S_e构成的基图像子空间中，使实时测量图像X_test用基图像估计向量的线性组合表示，为：

X_test＝c_test1·s₁+c_test2·s₂+...+c_testL·s_L＝S_e ^T·c_test

其中，c_test＝(c_test1，c_test2，...，c_testL)^T表示上述实时测量图像X_test在由上述基图 像估计矩阵S_e构成的基图像子空间中的投影系数，c_test＝piv(S_e ^T)·X_test，其中的piv(S_e ^T)表示矩阵S_e转置的广义逆，c_test即为待识别的合成孔径雷达实时测量图像X_test的独立成分特征；

(3)根据上述步骤(2)得到的合成孔径雷达训练样本图像X_train的独立成分特征c_train和待识别的合成孔径雷达实时测量图像X_test的独立成分特征c_test，对合成孔径雷达实时测量图像X_test进行识别分类，判断出被测目标的类别。

下面结合附图，详细介绍本发明的内容。

首先简要说明利用独立成分分析(ICA)方法进行合成孔径雷达(SAR)图像目标识别的相关概念和基本原理。

ICA的基本原理可以简单地用下面两个公式表示：

X＝AS

S＝WX

已有技术中利用ICA对SAR图像进行目标识别的基本流程如图1所示。

为了提高SAR图像目标识别的精度，在使用ICA对SAR图像进行特征提取之前，需要对图像进行预处理。预处理方法主要包括：①对原始SAR图像进行某种变换或使用另一种度量方式来表示源数据；②对变换后的图像进行归一化；③从原始图像中剪切出与识别目标相关的目标区域，也称为感兴趣区域(Region of Interest：ROI)，这些区域中可能包含一个或多个目标；其中，目标区域基本上包含所有与目标识别相关的回波，包括目标回波、阴影回波以及其他杂波等等，而噪声区域则几乎完全由背景杂波组成，故通常使用目标区域来对待识别目标进行ICA特征提取，用于目标的分类识别，而噪声区域只用作对噪声进行方差估计。

ICA适用于对一维数据的处理，在预处理之后，还需将二维图像样本拉直成一维行向量。而每个行向量可以看作是由多个向量线性叠加组成的混合信号，即，可以认为观测到的每一幅SAR图像都是由一组潜在的相互独立的基图像线性组合而成的。设x为一幅SAR图像，设S为由基图像向量构成的基图像矩阵，列向量s_i为其中的第i个基图像向量，则x可用S表示为

x＝b₁·s₁+b₂·s₂+...+b_n·s_n＝b·S

其中，系数向量b＝(b₁，b₂，...，b_n)即为图像样本x的ICA特征。

大部分ICA处理的数据要求为零均值白化数据，即满足E{X}＝0，E{XX^T}＝I。令X_O为训练数据，将X_O减去其均值向量m＝E{X_O}即可实现中心化。在使用中心化后的数据求出混淆矩阵A后，将源信号s的均值向量加回到中心化之后的s的估计值上以完成对源信号s的估计。其中，s的均值向量可通过A^-1m求得。

这里假设X_O已经是零均值。而对X_O的白化过程是让X_O的各成分互不相关，等价于X_O的协方差矩阵为单位矩阵(即 

)。可采用如下方法对X_O进行白化。

X＝W_ZX_O

$W_Z=2\·{\left(E\right\{X_O{X_O}^T\left\}\right)}^{-1/2}$

ICA方法要求多元数据的所有独立成分必须是非高斯分布的，或者最多可以有一个独立成分为高斯分布。而评估一个ICA模型的关键是数据独立成分的非高斯性。中心极限定理表明了一组相互独立的随机变量之和的分布通常趋近于高斯分布，由此可以得到下面推论：

推论1：n个独立随机变量的和比原始的任一变量更接近高斯分布。

根据无噪声的ICA模型X＝AS，设y为S中的某一独立成分，即

$y=w^{T}X=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{w}_i{x}_i$

其中w为去混淆矩阵W中对应y的待定列向量，w_i为其第i个元素，x_i为第i个观测信号。采用推论1的结论使w与矩阵A的逆矩阵的对应列相等。令列向量z＝A^Tw，则

$y=w^{T}x=w^T\mathrm{AS}=z^{T}S=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{z}_i{s}_i$

根据推论1，y比原来的独立变量s_i更具高斯性。如果令y的高斯性最小，则y将等于S中的某一成分s_k。y的非高斯性越大时，y与其他成分的独立性就越大，反之则越小；如果y的非高斯性达到最大，便可得到一个独立成分。因此，在ICA理论中，非高斯性等价于独立性。

本发明基于噪声独立成分分析的合成孔径雷达图像目标识别方法的总体流程如图4所示，其具体内容如下：

(1)本发明方法首先对输入的合成孔径雷达的原始训练样本图像X_train进行预处理，其过程是：

(1-1)对训练样本图像X_train进行对数变换，得到X_ln＝20ln(1+X_Orig)，其中X_Orig表示输入的原始训练样本单幅图像，X_ln表示经对数变换后噪声分布符合高斯分布的图像。

此步骤的目的是让合成孔径雷达噪声概率密度服从高斯分布，也就是通常所说的正态分布。通过背景技术部分可知，SAR图像的一个基本特点就是含有大量乘性斑点噪声，而在目标分类识别领域中，SAR图像的噪声还包括大量自然或人造杂波的干扰，因此几乎不可能简单地采用高斯分布对这些噪声进行拟合。

对数正态分布模型是S.F.George提出的一种适用于地面场景的SAR图像杂波统计模型，它是常用的描述非瑞利包络数据的一种统计模型，其主要思想是采用同态滤波器将SAR图像乘性噪声转化为加性高斯白噪声。对数正态分布的概率密度表达式为

$p\left(x\right)=\frac{1}{\sqrt{2\mathrm{\π}}\mathrm{\σ}}\exp [-\frac{{(1\mathrm{nx}-\mathrm{\μ})}^2}{{2\mathrm{\σ}}^2}]$

其中x是像素的灰度值，μ是lnx的均值(尺度参数)，σ是lnx的标准差(即形状参数)。而参数矩估计的表达式为

$\widehat{\mathrm{\μ}}=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\ln x_i,$ ${\widehat{\mathrm{\σ}}}^2=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{(\ln x_i-\widehat{\mathrm{\μ}})}^2$

因此，SAR图像数据经过对数变换之后，其噪声概率密度近似服从高斯分布，如图2所示。

(1-2)对上述X_ln进行去均值化处理，得到零均值图像X_t，X_t＝X_ln-E(X_ln)，其中期望函数E表示对图像求均值。

(1-3)将上述零均值图像X_t划分为目标区域X_o和噪声区域n_o，并满足{X_l}＝{n_o}∪{X_o}，将目标区域和噪声区域的二维图像数据分别拉直成一维行向量，其中的目标区域X_o包含与被识别目标相关的目标、阴影和杂波的回波，噪声区域n_o为背景杂波；

由于本发明采用的ICA算法只能处理一维信号，故在此步骤中还需要将每一幅图像数据的目标区域和噪声区域矩阵分别拉展成一维行向量。

(1-4)重复步骤(1-1)、(1-2)和(1-3)，得到训练样本中的N幅训练图像目标区域的一维拉直向量X_o，构成N行训练矩阵X_O，也即ICA算法中的输入矩阵X；

(1-5)一般情况下，每幅SAR图像都是独立成像，因此其包含的噪声也应是相互统计独立的，从而观测数据(即上步中的X_O)的噪声协方差阵为一个对角矩阵。根据上述得到的N幅训练图像的所有噪声区域n_o，可构造训练样本的噪声协方差对角矩阵 

并且其对角元素值分别对应为每幅图像各自的噪声方差。

从图2可观察到，对数变换后的SAR图像，噪声的概率分布形状近似服从高斯模型，又由于噪声区域n_o已经中心化，故可通过下式直接估计噪声协方差阵的对角元素 

${\mathrm{\Σ}}_o^i=E\left({\left(n_o^i\right)}^2\right)$

其中， 

表示第i幅训练图像的噪声区域。

(1-6)以上步骤中数据X_O已经中心化，其协方差矩阵为 并且，已经估计出了图像噪声n_O的噪声协方差矩阵∑_O。对上述训练矩阵X_O进行降维和白化处理，得到降维和白化后的合成孔径雷达训练样本图像子空间矩阵X。

本步骤可以借助主成分分析(PCA)或最可分成分(MDC)准则同时实现数据的白化及降维，详细说明如下。

首先以PCA方法为例：首先对非噪声数据协方差矩阵C-∑_O进行特征值分解(Eigenvalue Decomposition：EVD)，有

C-∑_O＝EDE^T

于是X_O可以通过如下公式实现白化：

X＝ED^-1/2E^TX_O

很容易证明这样得到的数据X满足E{XX^T}＝I。

由于白化过程并未将X_O中的噪声去除，因此重新估计X中噪声n的协方差矩阵为

∑＝E{nn^T}＝(ED^-1/2E^T)∑_O(ED^-1/2E^T)

通常情况下，原始数据X_O的维数很高，而导致巨大的计算量和冗余量，甚至出现过学习现象。因此，按PDC方法去掉D中太小的成分，选取L(L≤M)个能量最大的成分，最终得到图像新的特征值矩阵D_L，此过程也称为子空间选择。其中，D_L为从D中选取出的前L个最大特征值构成的对角矩阵，E_L为该L个特征值对应的特征向量组成的矩阵。这样，数据X将被降到L维，而∑为L×L大小，并且用E_LD_L ^-1/2E_L ^T取代了原式中的ED^-1/2E^T。从而采用PCA方法同时实现了输入数据的白化和降噪过程，为下一步的ICA特征提取过程输入了一个更有意义的子空间图像数据

X＝E_LD_L ^-1/2E_L ^TX_O

∑＝E{nn^T}＝(E_LD_L ^-1/2E_L ^T)∑_O(E_LD_L ^-1/2E_L ^T)

如背景技术中所述，采用PCA对SAR图像进行降维压缩的方法存在一定的缺陷。PCA选取前L个能量最大的主成分来代表原始图像数据，但成分(投影方向)的能量大小与其能表征图像类别属性的程度并没有直接关系；因此，在极端情况下，当前L个主成分均代表噪声的能量，且不同类别图像的噪声又基本相似时，将完全失去可用于图像分类的特征量，更无法对SAR图像进行正确地识别。

本发明对PCA方法进行了改进，提出了利用类别可分性(Class Discriminability)来衡量一个成分对图像类别属性的显示程度的最可分成分(Most discriminableComponent：MDC)分析法，从中选取前L个最可分成分来代表原始图像。MDC的具体选取方法如下：

设N个图像样本中共包含了C类目标，令属于第i类目标的样本数为N_i，则 

设p_k为PCA变换矩阵P中的第k个投影方向(即第k个成分)，b_ij表示第i类第j个样本在该方向上的投影系数(1≤i≤C，1≤j≤N_i)，所有投影系数b_ij的平均为 第i类样本投影系数b_i·j的平均为 定义类间方差 

类内方差 

MDC方法采用两者之比作为样本成分可分性的度量，表示为

r＝σ_between/σ_within

r的值越大，表示样本该成分的可分性越高。

用PCA或MDC方法也可以同时抑制图像噪声。将提取出的主成分采用前述的反变换得到原始SAR图像的近似，而近似图像与原图像相比，噪声得到了明显地减少和抑制。MDC方法去除了主成分中不具可分性或可分性很小的相似成分(包括主要噪声成分)，利用图像更少的能量获取了更高的可分性。尽管MDC方法去除的成分并非都是噪声成分，但这些成分对分类结果没有太大的积极作用，也皆可以视作噪声成分。

若采用MDC方法实现此步骤，只需在上述子空间选择过程中，按前述的MDC指标选取前L个最可分成分作为选取的子空间，而输入ICA的数据依然为

X＝E_LD_L ^-1/2E_L ^TX_O

∑＝E{nn^T}＝(E_LD_L ^-1/2E_L ^T)∑_O(E_LD_L ^-1/2E_L ^T)

在实现本步骤时，可首先选取主成分中能量占前95％的主成分(图像处理理论中，前95％的能量就可基本代表图像)，再利用MDC指标，分别计算这些主成分的可分性，将这些成分按可分性大小排列，选取前L个投影方向P′_L进行下一步的ICA特征提取。

(2)获取合成孔径雷达实时测量图像X_test，使用独立成分分析法分别提取上述合成孔径雷达训练样本图像X_train的独立成分特征和合成孔径雷达实时测量图像X_test的独立成分特征，具体过程如下：

(2-1)采用噪声快速独立成分分析法处理上述合成孔径雷达训练样本图像子空间矩阵X，得到一组由基图像估计向量组成的基图像估计矩阵S_e，表示为：

X＝c_train1·s₁+c_train2·s₂+...+c_trainm·s_L＝S_e ^T·c_train.

其中，s₁，s₂，…，s_L表示基图像估计矩阵S_e中L个基图像的估计列向量，c_train＝(c_train1，c_train2，...，c_trainL)^T表示上述训练样本图像子空间矩阵X在由上述基图像估计矩阵S_e构成的基图像子空间中的投影系数，c_train即为合成孔径雷达训练样本X_train的独立成分特征；

本发明采用的噪声ICA模型为Noisy FastICA(NFastICA)模型，NFastICA模型结合上述的预处理步骤，形成了一种新的对数正态噪声独立成分分析(Log-normal noise ICA：LnnICA)方法。NFastICA算法通过高斯矩(Gaussian Moments)的特征能直接从被高斯噪声污染了的观测图像数据中估计出潜在的随机变量，下面详述使用基于FastICA方法的NFastICA进行噪声图像独立成分特征提取的具体实现。

在前述理想ICA模型基础上，噪声ICA模型可表示为X＝AS+n.训练样本数据X_O已零均值化，协方差矩阵 

且噪声n_o满足高斯分布，其协方差矩阵为∑_O。对SAR图像数据进行白化处理，并用去除了噪声的协方差矩阵C-∑_O代替C，则白化操作应为

X＝(C-∑_O)^-1/2X_O

其中，X同样为噪声ICA模型的输入，而其噪声的协方差矩阵为

∑＝E{nn^T}＝(C-∑_O)^-1/2∑_O(C-∑_O)^-1/2

前面已经提到，在ICA理论中，非高斯性等价于独立性。在非高斯性的多种度量方式中，一种对随机变量的异常值具有鲁棒性的方法是用近似负熵来度量非高斯性。由信息论可知，无序性越高的随机变量其熵值越大，而对于高斯型变量，一个基本性质是其熵值在所有具相同方差的随机变量中为最大。定义负熵的度量

J(y)＝H(y_Gauss)-H(y)

其中，y_Gauss表示与变量y方差相同的高斯型随机变量，H(·)为求熵函数。这样，负熵J(y)的取值均为非负，而高斯型变量具有值为零的负熵。

为简化计算，通常采用近似的方法对负熵进行估计，即

J(y)≈c[E{G(y)}-E{G(v)}]²

$G_1\left(u\right)=\frac{1}{a_1}\log \cosh a_{1}u,\mathrm{where}1\≤a_1\≤2.$ 或 $G_2\left(u\right)=-\exp \left(\frac{{-u}^2}{2}\right)$

其中v为具有零均值和单位方差的标准正态变量，G为非二次函数，常取这种近似负熵的方法易于理解，计算快捷，具有很强的鲁棒性，因此FastICA中也采用了 这种非高斯性的度量方法，并采用牛顿下降法搜索最优解(各分量相互最独立)。

FastICA的学习规则是寻找一个方向向量w，使得输入数据在该方向上投影y＝w^TX的非高斯性最大，也即使采用负熵度量的近似负熵函数

J(y)≈[E{G(y)}-E{G(v)}]²

的值最大。同时，还要保证每一个独立成分y_i不重复，即W为正交单位阵。

标准正态分布变量v限制了y的方差必须也为1。若原始数据X为白化数据，那么对y单位方差的限制等价于将w的二阶范数归一化，即

E{(w^Tx)}²}＝||w||²＝1.

因此，FastICA算法可描述为以下步骤：

1)单变量搜索最优w_i。

step 1：任意选择一个初始权重向量w_i；

step 2：更新wi，令 

其中，其中g为非二次函数G的二阶导函数，η为牛顿法参数；

step 3：对step2得到的w_i进行归一化，即令

$w_i=w_i^{+}/\left|\right|w_i^{+}\left|\right|$

step 4：判断w_i是否收敛，如果收敛，则本步骤结束，否则返回step2。

2)对去混淆矩阵W进行全局的去相关搜索，保证W正交。

这里使用了直接在W上去相关的方法：将步骤1)得到的所有w_i组成去混淆矩阵W，迭代计算下面的搜索式，直到W收敛。

W＝3W/2-WW^TW/2

从而基图像矩阵可表示为

S＝WX＝WW_ZX_O＝W_IX_O

$X_O=W_I^{-1}\·S$

其中W_I＝WW_Z，包含两层意思，即W_I既去除了二阶相关性(W变换)，也去除了 ICA的高阶相关性(W_z变换)，从而达到了向量S各成分间的相互独立性。因此，不难得到X_O的ICA特征为 

而测试数据X_test的ICA特征为X_testS⁺(+表示广义逆)。

而降维之前的训练样本数据X_O进行PCA分解后，表示如下：

X_O＝RP^T

其中P是PCA的变换矩阵，每一列代表一个投影方向，且P^TP＝I；R为X_O在P上的投影。设P_L为P的前L个投影轴(前L个主成分)，R_L为X_O在其上的投影，即R_L＝X_OP_L；X的最小均方误差估计便为

${\hat{X}}_O=R_L{P}_L^T$

这样，在P_L上进行ICA，取L小于N，就可达到降维减少计算量的作用，并且类似地得到一个去混淆矩阵 

和 

且对应有

${\hat{W}}_I{P}_L^T=\hat{S}$

$P_L^T={\hat{W}}_I^{-1}\hat{S}$

从而得到

${\hat{X}}_O=R_L{\hat{W}}_I^{-1}\·\hat{S}$

即为X_O的ICA特征。同样的，X_test的ICA特征向量为

$b=R_{\mathrm{test}}{\hat{W}}_I^{-1}=X_{\mathrm{test}}{P}_L{\hat{W}}_I^{-1}$

而实际上，由于噪声的存在，E{G(y)}＝E{G(w^TX)}在近似负熵函数J_G(y)≈[E{G(y)}-E{G(v)}]²中已不再代表独立成分的统计量，而是独立成分与噪声成分之和的统计量。基于高斯矩的NFastICA的一个基本思想是选择G使其成为零均值高斯随机变量密度函数或与之相关的函数形式，从而使J_G能连贯地从一系列观测数据中简单地计算出来。若z为零均值非高斯变量，n是方差为σ²的高斯噪声，则可以简单地用代数式将E{G(z)}和E{G(z+n)}的关系表示出来。同样的，当G为一零均值高斯变量的密度函数或相关函数时，E{G(w^TX)}就可以根据含噪声的数据X直接估计出不含噪声的数据的J_G(y)。

根据假设，SAR图像的噪声服从高斯分布(实际中这种情况往往是不成立的，而这也正是本发明所解决的主要问题)。定义方差为c²的高斯密度函数

令 

表示 

的1(1＞0)阶导函数， 

表示 

本身， 

表示 

的1次积函数 

设 表示任意一个分布函数为 

的独立成分(非高斯)， 

表示方差为σ²的独立高斯噪声变量。由于 

函数由高斯函数衍生而来，因此称 

为 

的高斯矩，而对于任意c＞σ²，令 

得到

这说明了 

与其观测 

的高斯矩是等价的。另外还可以证明，用 

替代 

上式同样成立。因此，令 

从而可以使用随机变量的高斯矩从含噪声的观测中直接估计出无噪声的独立成分。令

其中 

应用牛顿法，在w归一化的条件下求上式的极大值，优化w得

$w_i^{+}=E\left\{\mathrm{Xg}\left(w_i^{T}X\right)\right\}-(I+\mathrm{\Σ})E\{\mathrm{Xg}\′\left(w_i^{T}X\right)\}{w}_i$

其中，g为函数G的导数，可取g₁(x)＝tanh(x)，g₂(x)＝x·exp(-x²/2)或g₃(x)＝x³.

其中， 

为高斯累计分布函数。

利用高斯矩的概念，本步骤在噪声呈高斯分布的情况下使用NFastICA，能够简便地从带噪声的观察数据中提取合成孔径雷达训练样本X_train的独立成分特征，用C_train表示。

(2-2)将合成孔径雷达实时测量图像(记为X_test)中待识别的图像数据投影到由上述基图像估计矩阵S_e张成的基图像子空间中，使实时测量图像X_test用基图像估计向量的线性组合表示，为：

X_test＝c_test1·s₁+c_test2·s₂+...+c_testL·s_L＝S_e ^T·c_test

其中，c_test＝(c_test1，c_test2，...，c_testL)^T表示上述实时测量图像X_test在由上述基图 像估计矩阵S_e张成的基图像子空间中的投影系数，通过c_test＝piv(S_e ^T)·X_test求得，其中的piv(S_e ^T)表示矩阵的广义逆，c_test即为合成孔径雷达实时测量图像X_test的独立成分特征；(3)根据上述步骤(2)得到的合成孔径雷达训练样本图像X_train的独立成分特征c_train和待识别的合成孔径雷达实时测量图像X_test的独立成分特征c_test，对合成孔径雷达实时测量图像X_test进行识别分类，判断出被测目标的类别。

(3)根据上述步骤(2)得到的合成孔径雷达训练样本图像X_train的独立成分特征c_train和待识别的合成孔径雷达实时测量图像X_test的独立成分特征c_test，对合成孔径雷达实时测量图像X_test进行识别分类，采用合适的分类器，判断出被测目标的类别。

分类器的选择需要综合考虑实际需求和分类效果等多种因素，例如，可以采用简单易行的最小均方误差(Mean Square Error：MSE)分类器或者实现高维核映射具有较优分类效果的支持向量机(Support Vector Mechine：SVM)分类器。

下面以MSE分类器为例说明对合成孔径雷达进行目标分类识别的具体步骤，为：

1)以向量2-范数的值作为待测目标与训练样本的ICA特征系数之间距离的度量函数(也可定义其它的距离度量，如马尔科夫距离等)，即单个测试样本与某一类训练样本之间的ICA特征距离可表示为

D_i＝||C_traini-C_tost||₂

其中，下标i表示已知的训练样本的类别，C_testi表示第i类训练样本的ICA特征，C_test表示单个待识别样本的ICA特征。

2)按步骤1)中定义的距离度量，对每个测试样本，分别求出其与各类型训练样本之间的距离度量(上述的ICA特征距离)。若有N类训练样本，则通过此步骤可以获得关于每个测试样本的N个距离{D₁，D₂...D_N}。

对每个测试样本，求出其与N类训练样本距离的最小值，将测试样本归入到与其距离最小的训练样本所属的类别当中。即对某测试样本，若其到第i类训练样本的距离D_i＝min{D₁，D₂...D_N}，那么MSE分类器将把该测试样本划分到第i类中，即将此目标识别为第i类。

Claims (1)

1.一种基于噪声独立成分分析的合成孔径雷达图像目标识别方法，其特征在于该方法包括以下步骤：

(1)对输入的合成孔径雷达的原始训练样本图像X_train进行预处理，具体过程如下：

(1-1)对训练样本图像X_train进行对数变换，得到X_ln＝20ln(1+X_Orig)，其中X_Orig表示输入的原始训练样本单幅图像，X_ln表示经对数变换后噪声分布符合高斯分布的图像；

(1-2)对上述X_ln进行去均值化处理，得到零均值图像X_t，X_t＝X_ln-E(X_ln)，其中期望函数E表示对图像求均值；

(1-3)将上述零均值图像X_t划分为目标区域X_o和噪声区域n_o，并满足{X_t}＝{n_o}∪{X_o}，将目标区域和噪声区域的二维图像数据分别拉直成一维行向量，其中的目标区域X_o包含与被识别目标相关的目标、阴影和杂波的回波，噪声区域n_o为背景杂波；

(1-4)重复步骤(1-1)、(1-2)和(1-3)，得到训练样本中的N幅训练图像的一维目标区域数据X_o，构成N行训练矩阵X_O；

(1-5)根据上述得到的N幅训练图像的所有噪声区域n_o，构造训练样本的噪声协方差对角矩阵

其中的对角元素

为每幅训练图像的噪声方差估计值，

表示第i幅训练图像的噪声区域；

(1-6)根据上述训练样本的噪声协方差矩阵∑_O，对上述训练矩阵X_O进行降维和白化处理，得到降维和白化后的合成孔径雷达训练样本图像子空间矩阵X；

(2)获取合成孔径雷达实时测量图像X_test，使用独立成分分析法分别提取上述合成孔径雷达训练样本图像X_train的独立成分特征和合成孔径雷达实时测量图像X_test的独立成分特征，具体过程如下：

(2-1)采用噪声快速独立成分分析法处理上述合成孔径雷达训练样本图像子空间矩阵X，得到一组由基图像估计向量组成的基图像估计矩阵S_e，表示为：

X＝c_train1·s₁+c_train2·s₂+...+c_trainm·s_L＝S_e ^T·c_train.

其中，s₁，s₂，…，s_L表示基图像估计矩阵S_e中L个基图像的估计列向量，c_train＝(c_train1，c_train2，...，c_trainL)^T表示上述训练样本图像子空间矩阵X在由上述基图像估计矩阵S_e构成的基图像子空间中的投影系数，c_train即为合成孔径雷达训练样本X_train的独立成分特征；

(2-2)将合成孔径雷达实时测量图像X_test中待识别的图像数据投影到由上述基图像估计矩阵S_e构成的基图像子空间中，使实时测量图像X_test用基图像估计向量的线性组合表示，为：

X_test＝c_test1·s₁+c_test2·s₂+...+c_testL·s_L＝S_e ^T·c_test

其中，c_test＝(c_test1，c_test2，...，c_testL)^T表示上述实时测量图像X_test在由上述基图像估计矩阵S_e构成的基图像子空间中的投影系数，c_test＝piv(S_e ^T)·X_test，其中的piv(S_e ^T)表示矩阵S_e转置的广义逆，c_test即为待识别的合成孔径雷达实时测量图像X_test的独立成分特征；

(3)根据上述步骤(2)得到的合成孔径雷达训练样本图像X_train的独立成分特征c_train和待识别的合成孔径雷达实时测量图像X_test的独立成分特征c_test，对合成孔径雷达实时测量图像X_test进行识别分类，判断出被测目标的类别。

Images