CN102662167A - 一种水下目标辐射噪声信号的特征提取方法 - Google Patents

Abstract

本发明提出的一种水下目标辐射噪声信号的特征提取方法，首先对声信号的可视化方法进行了研究，得到了时频可视化图像。接着，以时频可视化图像中声谱图为研究对象，对其进行可视化特征提取方法研究，包括：基于图像纹理特征提取和图像中心矩特征提取。在此基础上，分别利用PCA和LDA算法对上述方法得到的可视化特征进行了优化选择和降维，最后将其应用于水下目标识别的仿真中，实验结果表明，该混合特征下三类水下噪声信号的平均识别率可达91.97％。经过选择后的可视化特征参数的识别率与传统方法提取的特征参数相比，有了显著的提高。

Description

一种水下目标辐射噪声信号的特征提取方法

技术领域

本发明涉及一种水下目标辐射噪声信号的特征提取方法，

背景技术

水下目标的自动识别是现代声纳与水声对抗系统的重要组成部分。目标识别中最为关键的技术是特征提取。特征参数是否有效，在于其中包含的类别信息是否足够多，而干扰信息是否足够少，如果不能采用有效的特征，即便有出色的分类器也无法得到理想的分类识别结果。几十年来，已出现了一系列针对水下噪声信号的特征提取方法，主要有以下几种：基于谱分析的特征提取、基于舰船噪声的非线性特征提取、基于小波变换的特征提取法、神经网络特征提取法以及基于人耳听觉特征的特征提取方法。已有大量研究表明这些方法在实验室条件下具有良好的识别性能，但也都有一定的不足或受到一些因素的限制。

水下目标识别技术分为主动识别和被动识别两类。前者利用声纳发射脉冲声信号，根据接收到的回波信号特征对目标属性作出判决，其优点在于接收到的回波信号中携带着大量有利于分类识别，且反映目标本质特性的信息，但存在一个明显缺点，即易于暴露，不利于自我保护。被动目标识别可以弥补这一不足，它利用被动声纳接收到的水下目标辐射噪声特性进行分类判决，具有隐蔽性好的优点，但是也存在接收信号信噪比低等不利于识别的问题。

特征提取是水下目标识别中最为关键的技术之一，特征参数的种类及维数直接决定着分类识别的正确率。多年来已发展了大量的声学特征提取方法，在一定条件下都具有较好的分类识别效果。但对于复杂环境，提高识别系统鲁棒性仍是一项难题，除了解决通道失配等问题外，仍需要不断探索新的稳健特征提取方法。

到目前为止，还没有哪种特征参数能适用于不同环境下的水下目标识别，因此针对水下噪声信号的特征提取仍就是一个开放性的研究内容。本专利拓展了水下噪声信号的特征提取方法。

发明内容

要解决的技术问题

为了避免现有技术的不足之处，本发明提出一种水下目标辐射噪声信号的特征提取方法，首先将声信号转化为图像即可视化，然后提取其图像特征。这主要是考虑到图像是二维信号，相比于一维时间序列可能包含更多潜在的特征信息，且图像特征提取与识别已有大量成功的应用。

技术方案

一种水下目标辐射噪声信号的特征提取方法，其特征在于步骤如下：

步骤1：对水下目标的离散时域采样信号x(p)进行分帧处理得到分帧处理后信号x_n(m)，其中：p为时域采样点序号，p＝0，1，2，…，P-1，P为信号长度；n是帧序号，m是帧同步的时间序号，N为帧长，m＝0，1，2，…，N-1；

步骤2：对分帧处理后信号x_n(m)进行短时傅里叶变换：

X(n，e^jw)＝∑x_n(m)w(n-m)e^-jwm

其中w(m)为窗序列；

步骤3：计算X(n，e^jw)的短时幅度谱|X(n，k)|，其中0≤k≤N-1，得到功率谱函数Y(n，k)＝|X(n，k)|²＝(X(n，k))×(conj(X(n，k)))；

步骤4：以n为横坐标，以k为纵坐标，将Y(n，k)的值以图像表示得到声谱图；

步骤5：采用纹理特征提取方法和图像中心矩特征提取方法对声谱图进行图像的特征提取，得到30维可视化特征参数；

步骤6：对得到的30维可视化特征参数，利用主成分分析PCA及线性判别分析LDA的方法进行特征的优化选择及降维处理，得到20维混合特征。

所述步骤5中的纹理特征提取方法和图像中心矩特征提取方法的先后次序不受限制。

所述纹理特征提取方法为基于灰度共生矩阵的纹理特征提取或基于灰度-梯度共生矩阵的纹理特征提取。

有益效果

本发明提出的一种水下目标辐射噪声信号的特征提取方法，首先对声信号的可视化方法进行了研究，得到了时频可视化图像。接着，以时频可视化图像中声谱图为研究对象，对其进行可视化特征提取方法研究，包括：基于图像纹理特征提取和图像中心矩特征提取。在此基础上，分别利用PCA和LDA算法对上述方法得到的可视化特征进行了优化选择和降维，最后将其应用于水下目标识别的仿真中，实验结果表明，该混合特征下三类水下噪声信号的平均识别率可达91.97％。经过选择后的可视化特征参数的识别率与传统方法提取的特征参数相比，有了显著的提高。

附图说明

图1：A类中两个声样本的声谱图；

图2：B类中两个声样本的声谱图；

图3：C类中两个声样本的声谱图；

图4：三种方式不同维度下的总平均识别率；

图5：本发明方法的流程图。

具体实施方式

现结合实施例、附图对本发明作进一步描述：

在本实施例中的步骤中：

步骤1：对水下目标的离散时域采样信号x(p)进行分帧处理得到分帧处理后信号x_n(m)，其中：p为时域采样点序号，p＝0，1，2，…，P-1，P为信号长度；n是帧序号，m是帧同步的时间序号，N为帧长，m＝0，1，2，…，N-1；

步骤2：对分帧处理后信号x_n(m)进行短时傅里叶变换：

X(n，e^jw)＝∑x_n(m)w(n-m)e^-jwm

其中w(m)为窗序列；

步骤3：计算X(n，e^jw)的短时幅度谱|X(n，k)|，其中0≤k≤N-1，得到功率谱函数Y(n，k)＝|X(n，k)|²＝(X(n，k))×(conj(X(n，k)))；

步骤4：以n为横坐标，以k为纵坐标，将Y(n，k)的值以图像表示得到声谱图；

步骤5：采用纹理特征提取方法和图像中心矩特征提取方法对声谱图进行图像的特征提取，得到30维可视化特征参数；

步骤6：对得到的30维可视化特征参数，利用主成分分析PCA及线性判别分析LDA的方法进行特征的优化选择及降维处理，得到20维混合特征。

所述基于灰度共生矩阵的纹理特征提取方法是在图像灰度共生矩阵的基础上提取出能反映图像关于方向、相邻间隔、变化幅度的综合信息的8维特征参数，具体如下：

设图像的水平方向的像素总数为N_x，垂直方向的像素总数为N_y，为了降低计算量，将图像灰度作归一化，令其最高灰度级为N_g，记：

L_x＝{1，2，…，N_x}，L_y＝{1，2，…，N_y}           (4)

G＝{1，2，…，N_g}                                (5)

可以把待分析的图像理解为从L_x×L_y到G的一个映射，即L_x×L_y中的每个点，对应属于G中的一个灰度。

灰度共生矩阵描述的是在θ方向上，相隔d像元距离的一对像元，分别具有灰度层i和j出现的概率，其元素可记为P(i，j，d，θ)，通常θ取0°，45°，90°，135°。易知灰度共生矩阵是一个对称矩阵，其阶数是由图像中的灰度层数来决定。这个矩阵是距离和方向的函数，在规定的计算窗口或图像区域内统计符合条件的像元对数。四个方向的灰度共生矩阵的定义如下：

P(i，j，0°，d)＝#{[k，l)，(m，n)]∈(L_x×L_y)×(L_x×L_y)，           (6)

k-m＝0，|l-n|＝d，f(k，l)＝i，f(m，n)＝j}

P(i，j，45°，d)＝#{[(k，l)，(m，n)]∈(L_x×L_y)×(L_x×L_y)，

                                                                   (7)

|k-m|＝d，|l-n|＝d，f(k，l)＝i，f(m，n)＝j}

P(i，j，90°，d)＝#{[(k，l)，(m，n)]∈(L_x×L_y)×(L_x×L_y)，

                                                                   (8)

|k-m|＝d，l-n＝0，f(k，l)＝i，f(m，n)＝j}

P(i，j，135°，d)＝#{[(k，l)，(m，n)]∈(L_x×L_y)×(L_x×L_y)，k-m＝d，

l-n＝-d或k-m＝-d，l-n＝d，                                         (9)

f(k，l)＝i，f(m，n)＝j}

一般来讲较小的距离d会取得较好的效果，这里取d为1。

为了方便地的对共生矩阵的统计量进行描述，首先对共生矩阵进行正规化处理：

P(i，j)＝P(i，j)/R                                               (10)

这里R代表正规化常数，当θ＝0°或90°时，R＝2N_y(N_x-1)，当θ＝45°或135°时，R＝2(N_y-1)(N_x-1)。记：

$P_x\left(i\right)=\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}P(i,j)(i=\mathrm{1,2},...,N_g),$

$P_y\left(j\right)=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}P(i,j)\left(j=\mathrm{1,2},...,N_g\right),$

$P_{x+y}\left(k\right)=\underset{\stackrel{i=1}{i+j=k}}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}P(i,j)(k=\mathrm{2,3},...,{2N}_g),---\left(11\right)$

$P_{x-y}\left(k\right)=\underset{\stackrel{i=1}{|i-j|=k}}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}P(i,j)(k=\mathrm{0,1},...,N_g-1),$

像元间隔d取值为1，每个灰度图的灰度级量化为16，提取了4个方向(0°，45°，90°，135°)的角二阶矩、熵、惯量和相关性四个特征参数，针对每个参数，求取4个方向上的均值及均方差，从而形成了一个8维的特征矢量。各特征参数具体的定义如下：

1、角二阶距(能量)

$f_1=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}{\left[P(i,j)\right]}^2---\left(12\right)$

2、熵

$f_2=-\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}P(i,j)\log \left[P(i,j)\right]---\left(13\right)$

3、对比度

$f_3=\underset{n=0}{\overset{N_{g-1}}{\mathrm{\Σ}}}{n}^2\left\{\underset{\stackrel{i=1}{|i-j|=1}}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}P(i,j)\right\}---\left(14\right)$

4、相关性

$f_4=\{\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}i\×j\×P(i,j)={\mathrm{\μ}}_x{\mathrm{\μ}}_y\}/{\mathrm{\δ}}_x{\mathrm{\δ}}_y---\left(15\right)$

其中μ_x，δ_x分别为{P_x(i)；i＝1，2，…，N_g}的均值和均方差μ_y，δ_y分别为{P_y(j)；j＝1，2，…，N_g}的均值和均方差；

角二阶矩是灰度共生矩阵阵元分布不均匀性的度量，与它相关的另一个度量是能量，两者性质相似；熵主要用来检测图像空间的复杂度和内部的均匀性，细纹理的熵较大，粗纹理的熵较小；对比度，又称主对角惯性矩，是图像小区域内灰度变化总量，反映了高值阵元远离对角线的程度；相关性是用来描述灰度共生矩阵中行或列元素之间相似程度的，它反映某种灰度值沿某些方向的延伸长度，若延伸得越长，则相关系数越大。

所述基于灰度-梯度共生矩阵的纹理特征方法是利用图像的灰度和梯度的综合信息提取出15维特征参数。图像的灰度-梯度共生矩阵集中反映了图像中像素的灰度和梯度的相互关系，各像素的灰度是构成一幅图像的基础，而梯度则是构成图像边缘轮廓的要素。灰度-梯度空间可以清晰地表现图像内像素灰度与梯度的分布规律，同时也体现了各像素与其领域像素的空间关系，对图像的纹理能很好地描绘。具体如下：

灰度-梯度共生矩阵的元素H(x，y)定义为在归一化的灰度图像F(i，j)及其归一化的梯度图像G(i，j)中具有灰度值x和梯度值y的像素数，即在集合{(i，j)|F(i，j)＝x∩G(i，j)＝y，i，j＝0，1，…，N-1}中元素的个数。其中，F(i，j)∈[0，L-1]，G(i，j)∈[0，L_g-1]。

对灰度-梯度共生矩阵进行归一化处理，使其各元素之和为1。如式16所示：

${\hat{H}}_{(x,y)}=\frac{H(x,y)}{\underset{x=0}{\overset{L-1}{\mathrm{\Σ}}}\underset{y=0}{\overset{L_g-1}{\mathrm{\Σ}}}H(x,y)}---\left(16\right)$

而 $\underset{x=0}{\overset{L-1}{\mathrm{\Σ}}}\underset{y=0}{\overset{L_g-1}{\mathrm{\Σ}}}H(x,y)=N\×N=N^2,$ 所以上式可以表示为：

${\hat{H}}_{(x,y)}=\frac{H(x,y)}{N^2}---\left(17\right)$

该灰度-梯度共生矩阵的原点在左上角，向右梯度值增加，向下灰度值增加。对于粗纹理的图像，在灰度轴附近集中分布，仅少量边界点远离灰度轴分布。对于细纹理图像，

离开灰度轴而沿梯度轴散开分布。灰度-梯度共生矩阵常用的数字特征见表1。

表1灰度-梯度数字特征

因梯度表达了图像灰度的最大变化率，图像各点梯度的大小则表达了该点处等灰度线的密集程度，所以大梯度优势和小梯度优势可从一定程度上反映图像的灰度变化剧烈程度。当图像的灰度变化平缓时，小梯度像素数多，小梯度优势大。而当图像的灰度变化剧烈时，大梯度像素数多，大梯度优势大。对于每一幅声谱图，最终可以得到15维大小的特征矢量。

所述图像中心矩特征提取是对图像的灰度值进行基于灰度统计量的计算过程，具体如下：

归一化中心矩特征：矩是描述图像属性的一种重要方法，下面来对其进行定义。假设图像的灰度值用f(x，y)表示，于是图像的p+q阶矩定义为：

$m_{\mathrm{pq}}=\underset{x}{\mathrm{\Σ}}\underset{y}{\mathrm{\Σ}}{x}^p{y}^{q}f(x,y)(p,q=\mathrm{0,1,2},...)---\left(18\right)$

其对应的中心矩定义为：

${\mathrm{\μ}}_{\mathrm{pq}}=\underset{x}{\mathrm{\Σ}}\underset{y}{\mathrm{\Σ}}{(x-\stackrel{\‾}{x})}^p{(y-\stackrel{\‾}{y})}^{q}f(x,y)(p,q=\mathrm{0,1,2},...)---\left(19\right)$

其中，

为质心坐标。图像的归一化中心矩定义为：

$M_{\mathrm{pq}}=\frac{{\mathrm{\μ}}_{\mathrm{pq}}}{{\mathrm{\μ}}_{00}^r}(p,q=\mathrm{0,1,2},...)---\left(20\right)$

其中，

(p+q＝2，3，…)。取M₀₂，M₀₃，M₁₁，M₁₂，M₂₀，M₂₁，M₃₀这7个不同阶数的归一化中心矩作为特征参数。

由上可以看出，对于一个噪声信号的声谱图，经特征提取后最终可以得到一个30维的特征矢量。包括：8维基于灰度共生矩阵的纹理特征、15维基于灰度-梯度共生矩阵的纹理特征以及7维归一化中心矩特征。

所述主成分分析法是在数据空间中找到一组向量以尽可能地解释数据的方差，从而用降维后的低维向量保存原数据中的主要信息，使不同类别的样本更易于区分。主成分分析(Principal Component Analysis，PCA)主要原理就是基于对原始数据进行统计分析，利用线性变换，对高维数据进行分析和压缩。其求解过程如下：

首先需要将特征矩阵转化成一维行向量，这里将转换后的p个训练样本记为{X₁，X₂，…，X_p}，基于PCA的特征提取的主要步骤如下：

(1)计算训练样本的均值

$\stackrel{\‾}{X}=\frac{1}{p}\underset{i=1}{\overset{p}{\mathrm{\Σ}}}{X}_i---\left(21\right)$

(2)计算训练样本的协方差矩阵

$Q=\underset{i=1}{\overset{p}{\mathrm{\Σ}}}{(X_i-\stackrel{\‾}{X})}^T(X_i-\stackrel{\‾}{X})---\left(22\right)$

(3)计算协方差矩阵的特征值λ及对应的特征向量。取最大的k个正特征值对应的特征向量组合得到转换矩阵V，V为N×k维矩阵；而k的取值通常是根据贡献率来确定的，定义前k个主分量的累积贡献率为：

当

足够大(例如0.7，0.8，0.9等)时，认为由这k个特征值对应的特征向量组成的矩阵可以很好的表征原训练样本矩阵。

(4)将所有训练样本组合得到训练样本矩阵X，利用转换矩阵V对X进行特征提取，得到一个p×k维训练样本特征矩阵：

C＝XV                     (24)

所述线性判别分析法是将具有标签(类别)的数据通过投影的方法，投影到新的维度空间中，使得投影后的点，形成按类别区分、一簇一簇的情况。最终相同类别的点，将会在投影后的空间中更接近。线性判别分析方法(Linear Discriminant Analysis，LDA)是一种监督性学习算法。这种方法的最终目的就是找到一些特征使得类间离散度与类内离散度的比值最大。算法的实施过程如下：

假设用来区分二分类的直线(投影函数)为：

y＝w^Tx                       (25)

类别i的原始中心点为：

$m_i=\frac{1}{n_i}\underset{x\∈C_i}{\mathrm{\Σ}}x---\left(26\right)$

其中C_i表示类别为i的点的集合。

类别i投影后的中心点为：

${\tilde{m}}_i=w^T{m}_i---\left(27\right)$

投影后类别点之间的离散度(方差)为：

${\tilde{s}}_i=\underset{y\∈Y_i}{\mathrm{\Σ}}{(y-{\tilde{m}}_i)}^2---\left(28\right)$

于是可以得到投影到w后的损失函数：

$J\left(w\right)=\frac{{|{\tilde{m}}_1-{\tilde{m}}_2|}^2}{{\tilde{s}}_1^2+{\tilde{s}}_2^2}---\left(29\right)$

上式中分母表示每一个类别内的方差之和，方差越大表示一个类别内的点越分散，分子为两个类别各自的中心点的距离的平方，最大化J(w)就可以求出最优的w。

样本的内类离散度矩阵可以定义为：

$s_w=\underset{i=1}{\overset{c}{\mathrm{\Σ}}}{p}_i{s}_i---\left(30\right)$

S_i定义为：

$s_i=\underset{x\∈C_i}{\mathrm{\Σ}}(x-m_i){(x-m_i)}^T(i=\mathrm{1,2},...,C)---\left(31\right)$

其中P_i为先验概率，m_i是C_i类的均值，m是所有样本的均值。

将s_i代入J(w)表达式中，可以得到Fisher准则函数如下：

$J\left(w\right)=\frac{w^T{s}_{b}w}{w^T{s}_{w}w}---\left(32\right)$

求出J(w)最大值时对应的特征向量，可通过数学变换得到，令J(W_opt)＝λ，即要求最大的λ。

$\mathrm{\λ}=\frac{w^T{s}_{b}w}{w^T{s}_{w}w}---\left(33\right)$

w^Ts_bw＝λw^Ts_ww                                         (34)

s_bw＝λs_ww                                             (35)

对于一个m分类问题，w就是满足如下等式的解：

s_bw_i＝λs_ww_i(i＝1，2，…，m)                           (36)

s_w ^-1s_bw_i＝λw_i(i＝1，2，…，m)

即求s_w ^-1s_b得特征值及其特征向量问题。

PCA的目的在于找到对分类贡献大的一组成分，最终实现特征的降维；而LDA可以使同一类样本间距离变小，而不同类样本间距离变大，这一过程并没有对特征进行降维，而是将样本空间映射到另一个空间中，使得同类样本更集中，不同类样本间间距变大。

本实施例中用来论证的仿真实验中对象是A、B、C三类水下目标辐射噪声信号，采样频率为8KHz，每类中有8个样本，每个样本的时长为5s。为了增加实验的样本数，得到更为有参考价值的信息，对每一幅噪声信号的可视化图像进行旋转，在0°～180°的范围内，每次旋转10°得到一幅图像，将每次旋转后的图像保存下来作为一个新的样本，于是每一幅原始声谱图经选择后可以得到18幅图像，对应着18个图像样本。因此三类噪声信号就有3*8*18＝432，在分类识别实验中训练集与测试集的样本数各取216个。表2给出的三类样本的30维可视化特征下各自及总的平均识别率。

表2 30维可视化特征下样本的识别率(％)

在实际处理中，除了可以单独使用LDA和PCA算法进行处理，还可以将两类算法结合起来使用，一种方式是先进行PCA处理，去除样本空间中冗余的信息，在此基础上对降维后的空间进行LDA处理，增大各类样本间距离；另一种方式则是相反，先进行LDA处理，然后由PCA进行降维。PCA、LDA+PCA以及PCA+LDA处理后的可视化特征的平均识别率见表3。

表3三种处理方法下的平均识别率(％)

由表3可以看出：

(1)经过降维后(维数大于等于8时)的样本的识别率均有所提高，这表明经过PCA处理后的数据中的对分类不利的信息被去除了很多，且样本间的相关性变小，有利于区分各类样本。

(2)与表2相比，LDA+PCA这种处理方式效果不理想。

(3)PCA+LDA这两种方法效果均比较好，在维度不低于8的情况下，A类的识别率都在95％以上，而B类的识别率也都在72％以上，均远高于处理前的识别率。尤其是在维度为20的情况下，识别率最为理想。

图4所示的是PCA、LDA+PCA以及PCA+LDA这三种处理方式下的三类样本总平均识别率随着维度降低而变化的趋势图。可以看出：在维度不小于8的时候，PCA+LDA的效果最佳，尤其是在维度为20的时候，总平均识别率最高；而在维度小于24的时候，PCA的效果好于LDA+PCA。

综上所述，本实施例选择了先PCA后LDA的特征选择算法；经特征选择后的20维特征参数较之于原始30维特征来讲，不仅维度大幅度下降，而且三类目标各自以及总的平均识别率均有显著提高。

Claims (3)

1.一种水下目标辐射噪声信号的特征提取方法，其特征在于步骤如下：

步骤1：对水下目标的离散时域采样信号x(p)进行分帧处理得到分帧处理后信号x_n(m)，其中：p为时域采样点序号，p＝0，1，2，…，P-1，P为信号长度；n是帧序号，m是帧同步的时间序号，N为帧长，m＝0，1，2，…，N-1；

步骤2：对分帧处理后信号x_n(m)进行短时傅里叶变换：

X(n，e^jw)＝∑x_n(m)w(n-m)e^-jwm

其中w(m)为窗序列；

步骤3：计算X(n，e^jw)的短时幅度谱|X(n，k)|，其中0≤k≤N-1，得到功率谱函数Y(n，k)＝|X(n，k)|²＝(X(n，k))×(conj(X(n，k)))；

步骤4：以n为横坐标，以k为纵坐标，将Y(n，k)的值以图像表示得到声谱图；

步骤5：采用纹理特征提取方法和图像中心矩特征提取方法对声谱图进行图像的特征提取，得到30维可视化特征参数；

步骤6：对得到的30维可视化特征参数，利用主成分分析PCA及线性判别分析LDA的方法进行特征的优化选择及降维处理，得到20维混合特征。

2.根据权利要求1所述水下目标辐射噪声信号的特征提取方法，其特征在于：所述步骤5中的纹理特征提取方法和图像中心矩特征提取方法的先后次序不受限制。

3.根据权利要求1或2所述水下目标辐射噪声信号的特征提取方法，其特征在于：所述纹理特征提取方法为基于灰度共生矩阵的纹理特征提取或基于灰度-梯度共生矩阵的纹理特征提取。

Images