基于信号时频特征和支持向量机的供水管道泄漏识别方法与流程

本发明涉及漏水检测定位技术领域，特别是涉及一种基于信号时频特征和支持向量机的供水管道泄漏识别方法。

背景技术：

水是人类和一切生物赖以生存的物质基础,是人类社会发展不可缺少的自然资源。联合国发布的《2018年世界水资源开发报告》显示，由于人口增长、经济发展和消费方式转变等因素，全球水资源的需求正在以每年1％的速度增长，而这一速度在未来20年还将大幅加快。随着人口的不断增长和污染的持续恶化，水资源的缺乏日益加剧。在全球水资源严重短缺的情况下，水资源的浪费问题却十分严重。其中，供水管道泄漏造成的水资源浪费十分严重。世界银行的一项研究表明，全世界由于供水管道的泄漏造成的水损失量达到每年486亿立方米，相应的经济损失约为每年146亿美元。因此，研究高效的供水管道泄漏检测定位技术对于保护水资源和促进经济发展具有重要意义。

为了能够检测出地下供水管道的泄漏情况，学术界和工业界开展了大量研究工作，开发了许多有效的检测方法。最早出现的检测方法是音听法。该方法是检测人员借助听音检测设备根据泄漏声的大小与音质特点来判断漏水区域。这种方法虽然操作简便，但依赖检测人员的经验，并且由于供水管网分布广泛，所以该方法工作量大、可靠性低。地面穿透雷达可以通过检测由漏水引起的土壤空隙来确定管道中的泄漏位置，然而，由于不同地区地质结构相对复杂，所以该方法实用性较差，并且成本昂贵。根据供水管道内部压力变化学者们相继提出了压力梯度法、负压波法、流量平衡法。这些方法对管内水流压力值较为敏感，但是，由于供水管网水流持续波动,当波动幅度较大时容易产生误报。研究发现泄漏信号频谱较为集中，且管道的振动频率仅与泄漏情况有关。利用这一特性通过对管道上的压电加速度传感器获取的信号进行谱分析，从而进行泄漏检测。但是，当存在与泄漏信号频谱相似的环境噪声时，该方法容易出现误判。研究人员结合泄漏声信号线性预测编码双谱系数(lpcc)和隐马尔可夫模型(hiddenmarkovmodels,hmm)进一步提高了对泄漏信号和环境干扰信号的辨识能力。然而受算法本身的限制，随着系统运行时间的增加，该方法的错误概率会上升。对于大规模的供水管网系统，学者们尝试使用实时模型法在管网上建模，将管网上的测量数据与流量模型的预测值进行比较，但是在实际应用中该方法建模难度系数高、数据运算量大。

技术实现要素：

为了解决以上问题，本发明提出了一种基于信号时频特征和支持向量机的供水管道泄漏识别方法，该方法利用漏水信号随机性、频谱较集中等特点，根据泄漏信号和非泄漏信号时频特性的差异性基于信号固有模态函数、近似熵、主成分提出了三个时频特征。利用这些特征构造特征矩阵作为支持向量机输入，将支持向量机作为分类器对信号进行识别并输出识别结果，从而解决了现有管道泄漏检测技术存在的建模难度系数大、误判率高等问题。

根据本发明的基于信号时频特征和支持向量机的供水管道泄漏识别方法，所述方法根据泄漏信号和非泄漏信号时频特性的差异性，基于信号固有模态函数、近似熵、主成分提出了三个时频特征，并利用这些特征构造特征矩阵作为支持向量机输入，将支持向量机作为分类器对信号进行识别并输出识别结果，确定信号为泄漏信号还是非泄漏信号。

进一步的，所述方法包括：

s1：输入被检测的信号；

s2：对输入的信号进行特征提取；

s3：将所提取的特征集输入已进行优化的支持向量机，利用支持向量机对特征进行识别；

s4：支持向量机根据输入的信号特征输出识别结果，确定信号为泄漏信号还是非泄漏信号。

进一步的，s2中所提取的特征包括以下三个时频特征：基于固有模态函数的频域特征、基于近似熵的特征和基于主成分分析的特征。

进一步的，提取基于固有模态函数的频域特征步骤包括：

利用经验模态分解法对输入的信号进行处理，分解得到所述输入的信号的多个固有模态函数；

对所得到的固有模态函数进行处理，得到信号的固有模态函数功率谱；

计算固有模态函数功率谱的均值，作为泄漏信号的时频特征。

进一步的，提取基于近似熵的特征步骤包括：

根据采集管道信号的n个样本u(1),u(2),…,u(n)构造出两个长度为m的序列，x(i)＝[u(i),u(i+1),…,u(i+m-1)]，x(j)＝[u(j),u(j+1),…,u(j+m-1)]，其中i,j≤n-m+1，计算序列x(i)和x(j)之间的距离，

d[x(i),x(j)]＝maxk＝1,2,…,m[|u(i+m-1)-u(j+k-1)|]

给定一个阈值r，对于每一个i＜n-m+1统计出d[x(i),x(j)]≤r的数目，并计算出这个数目与矢量个数之间的比值：

对于所有i值，求的平均值为φ^m(r)，

将m增加1，重复以上步骤得到φ^m+1(r)，根据φ^m+1(r)和φ^m(r)得到近似熵的值为：

apen(m,r)＝φ^m(r)-φ^m+1(r)，作为泄漏信号的时频特征。

进一步的，r的取值为被检测信号标准差的0.1至0.2倍。

进一步的，提取基于主成分分析的特征步骤包括：

采集m组管道信号x1，x2，…，xm，每组信号含有n个样本，表示为xi＝(x1i,x2i,…,xni)^t，则由此构成的n×m阶矩阵x＝[x1x2…xm]为

利用x的协方差矩阵前l个特征值，其中，0＜l≤m，所述特征值从大到小排列，对应的特征向量αi＝(α1i,α2i,…,αmi)^t(i＝1,2,···,l)可以求得l个新的向量，

yi＝xαi,(i＝1,2,…,l)

yi为x的主成分；

利用信号的主成分yi构造n×l阶的主成分信号矩阵y＝[y1y2…yl]，根据主成分信号矩阵与原始信号矩阵内积gji＝[yj,xi]构造矩阵g，

g＝y^t(x-e[x])

选取gj＝[gj1gj2···gjm],0＜j≤l作为泄漏信号的时频特征。

进一步的，l的取值由贡献率和累计贡献率来确定。

进一步的，所述支持向量机采用以下方式进行优化：

在不同时间段对不同区域管道发生泄漏和未发生泄漏时进行信号采集构造信号样本库，从信号样本库中随机选取信号作为训练样本信号和测试样本信号；

利用训练样本信号的特征集对支持向量机进行训练，形成支持向量机辨识模型；

利用测试样本信号的特征集对训练完成的支持向量机辨识模型进行测试；

根据测试结果对支持向量机做进一步的优化，直到测试输出的准确率满足要求；

形成用于管道泄漏识别的支持向量机模型。

本发明的有益效果：

(1)本发明基于信号固有模态函数、近似熵、主成分提出了三个时频特征。利用这些特征构造特征矩阵作为支持向量机输入，使得本发明的检测效果更全面，避免单一考虑某一特征时产生误判概率较大的问题，有效提高泄漏检测正确率；

(2)本发明采用经验模态分解(empiricalmodedecomposition，emd)对管道信号做时频分析，把复杂信号分解为有限的固有模态函数(intrinsicmodefunction，imf)之和的形式，从而实现了利用多个imf对泄漏信号的功率谱密度进行多尺度分析，实现高精度漏水定位；

(3)主成分。

本发明利用信号主成分分量构造主成分矩阵，然后利用主成分矩阵与原始信号矩阵的内积gji＝[yj,xi]构造基于主成分的特征矩阵，从而实现了多个可以用于泄漏检测的基于主成分的时域特征。

附图说明

图1示出了根据本发明的基于信号时频特征和支持向量机的供水管道泄漏检测方法流程图；

图2a示出了管道正常运行时的信号功率谱；

图2b示出了管道发生泄漏时的信号功率谱；

图3示出了管道漏水信号的固有模态函数；

图4示出了管道泄漏信号与未泄漏信号的近似熵；

图5示出了支持向量机辨识模型训练与优化原理图；

图6示出了支持向量机参数组合(c,γ)下的识别准确率；

图7示出了管道泄漏与未泄漏检测结果。

具体实施方式

下文将结合具体附图详细描述本发明具体实施例。应当注意的是，下述实施例中描述的技术特征或者技术特征的组合不应当被认为是孤立的，它们可以被相互组合从而达到更好的技术效果。在下述实施例的附图中，各附图所出现的相同标号代表相同的特征或者部件，可应用于不同实施例中。

如图1所示，根据本发明的基于信号时频特征和支持向量机的供水管道泄漏检测方法包括以下步骤：

步骤101：输入被检测的信号；

步骤102：对输入的信号进行特征提取；

步骤103：将所提取的特征集输入已进行优化的支持向量机，利用支持向量机对特征进行识别；

步骤104：支持向量机根据输入的信号特征输出识别结果，确定信号为泄漏信号还是非泄漏信号。

基于固有模态函数的频域特征

大量研究表明，未发生泄漏时管道上的信号与发生泄漏时的信号功率谱具有显著差异，泄漏信号的功率谱成分主要集中在特定频带内。因此，可以将信号的功率谱作为管道泄漏识别的特征。为了提取信号功率谱密度的差异性，我们采用经验模态分解(empiricalmodedecomposition，emd)对管道信号做时频分析。经验模态分解能够选择分解信号的尺度，把复杂信号分解为有限的固有模态函数(intrinsicmodefunction，imf)之和的形式，从而实现了利用多个imf对泄漏信号的功率谱密度进行多尺度分析。

对一个函数s(t)做hilbert变换，得

然后，基于x(t)和y(t)构造解析函数，

解析信号的相位函数为

将相位函数对时间求导，可以得到解析信号的瞬时频率函数为

从瞬时频率函数的定义可知，瞬时频率是时间的函数，在某一特定时刻仅有一个频率值与其相对应。瞬时频率在某些情况下会出现没有意义的负频率，如果全部瞬时频率都为正频率，则将此函数s(t)称为固有模态函数(imf)。因此，固有模态函数必须满足以下两个条件：

(1)该函数极值点的数量ne(包括极小值点和极大值点)与过零点的数量相等或者相差最多为1，即

(ns-1)≤ne≤(ns+1)(5)

(2)在任意时刻ti上，函数局部极大值确定的上包络线和局部极小值确定的下包络线的均值为零，即

[smax(t)+smin(t)]/2＝0,ti∈[ta,tb](6)

式中[ta,tb]为一段时间长度。

限定条件(1)表明s(t)中既不能出现小于零的极大值也不出现大于零的极小值。条件(2)是去除局部波动而产生的波形不对称。通常，一个信号可以包含多个固有模态函数，利用经验模态分解法可以提取一个复杂信号的固有模态函数。为了求得泄漏信号的功率谱密度特征，我们首先利用经验模态分解法对信号进行处理，得到该信号的固有模态函数。然后，根据固有模态函数可以求出不同模态函数的功率谱，进而提取频域特征。

首先，将原始信号x(t)的所有极值点分别用三次样条曲线连接起来，得到x(t)的上、下包络线，使该信号处于两条包络线之间。同时，令两条包络线均值组成的函数为m(t)。从原始信号x(t)中减去其上、下包络线均值m(t)，得到：

h1(t)＝x(t)-m(t)(7)

然后，检测h1(t)是否满足imf的两个条件，如果不满足则对h1(t)重复上述操作，直到满足imf的条件。将此时的h1(t)表示为c1(t)，则c1(t)为信号x(t)的第一个imf，

c1(t)＝h1(t)(8)

进一步，从原始信号x(t)中减去c1(t)，

r1(t)＝x(t)-c1(t)(9)

将r1(t)作为一个新的信号，利用上述方法求得r1(t)的第一个imf，该函数就是x(t)的第二个imf，记为c2(t)。依次类推，可以逐步求得信号x(t)的第n个imfcn(t)和余项rn(t)。

通过上述步骤，原始信号x(t)被分解成为n个imf和一个余项之和，

通常情况下，判断imf的条件(2)是无法满足的，一般设定一个停止准则，当满足停止准则时，即认为条件(2)成立。为此，设两个连续的处理结果之间的标准差sd为，本方法中当标准差sd满足

即认为条件(2)满足。其中，t为信号的观测长度，hk-1(t)和hk(t)是imf求解过程中两个连续处理结果。研究表明，标准差sd的阈值通常可以取0.2-0.3。

图2a示出了管道正常运行时的信号功率谱，图2b示出了管道发生泄漏时的信号功率谱。分析结果表明，管道泄漏声信号的频率主要集中在1.6khz附近。图3为管道泄漏信号的经验模态分解结果，实验中标准差sd的阈值设为0.3。实验结果表明，该信号经过经验模态分解得到了5个固有模态函数。可以利用这5个固有模态函数从多尺度对泄漏信号的频域特征进行分析。

对emd分解后的imf分量ci(n)求其离散傅立叶变换得ci(k)，

对ci(k)求其模平方，得到信号的功率谱

然后，求式(13)的均值

本发明将信号固有模态函数功率谱的均值作为泄漏信号的频域特征。

基于近似熵的特征

管道泄漏是一种局部小概率事件，所以泄漏信号与未泄漏信号在随机性上存在一定差异，可以从分析信号随机性的角度提取泄漏信号的特征。近似熵(approximateentropy，apen)是相似向量由m维增加到m+1维时继续保持其相似性的条件概率，是当维数变化时产生新模式概率的大小。产生新模式的概率越大，信号越复杂，相应的近似熵越大。因此，本发明选取近似熵作为泄漏信号识别的特征之一。

首先，根据采集管道信号的n个样本u(1),u(2),…,u(n)构造出两个长度为m的序列，x(i)＝[u(i),u(i+1),…,u(i+m-1)]，x(j)＝[u(j),u(j+1),…,u(j+m-1)]，其中i,j≤n-m+1。然后，计算序列x(i)和x(j)之间的距离，

d[x(i),x(j)]＝maxk＝1,2,…,m[|u(i+m-1)-u(j+k-1)|](15)

给定一个阈值r，对于每一个i＜n-m+1统计出d[x(i),x(j)]≤r的数目，并计算出这个数目与矢量个数之间的比值：

对于所有i值，求的平均值为φ^m(r)，

将m增加1，重复(15)—(17)得到φ^m+1(r)，根据φ^m+1(r)和φ^m(r)可以得到近似熵的值为，

apen(m,r)＝φ^m(r)-φ^m+1(r)(18)

以上分析表明近似熵是一个无量纲的标量，其值的大小与m和r有关。为了使近似熵的取值具有较为合理的统计特性，根据经验m＝2，r的取值一般为序列标准差(sd,standdeviation)的0.1-0.2倍。图4为铸铁管道发生泄漏与未发生泄漏时信号的近似熵。实验中，两种情况各提取50组数据，每组数据长度为5000,m＝2,r＝0.2sd。从图4结果可以看出，泄漏时信号的apen均值明显高于未泄漏时的均值，这表明泄漏信号的随机性特征要高于无泄漏信号，可以作为泄漏识别的特征。

基于信号主成分的特征

主成分分析(principalcomponentsanalysis，pca)法是经典的特征提取方法。该方法是利用降维的思想，将多个变量转化为少数几个综合变量(即主成分)。其中，每个主成分都是原始变量的线性组合，各主成分间互不相关。主成分能够反映原始变量的绝大部分信息，且所有信息互不重叠。本发明利用主成分分析法来分析管道泄漏和未泄漏信号的差异性。

采集m组管道信号x1，x2，…，xm，每组信号含有n个样本，可以表示为xi＝(x1i,x2i,…,xni)^t，则由此构成的n×m阶矩阵x＝[x1x2…xm]为

根据pca的原理，利用x的协方差矩阵前l(0＜l≤m)个特征值(从大到小排列)对应的特征向量αi＝(α1i,α2i,…,αmi)^t(i＝1,2,···,l)可以求得l个新的向量，

yi＝xαi,(i＝1,2,…,l)(20)

称yi为x的主成分，式(20)中满足yi与yj(i≠j；i,j＝1,2,…,l)不相关。y1为x1，x2，…，xm的一切线性组合中方差最大者，y2为与y1不相关的x1，x2，…，xm所有线性组合中方差最大者。以此类推，yl为与y1，y2，…，yl-1均不相关的x1，x2，…，xm所有线性组合中方差最大者。在实际应用中，l的取值可以由贡献率和累计贡献率来确定。贡献率的定义为，

其中，λ＝[λ1λ2···λm]是x协方差矩阵特征值由大到小排列。然后，通过设定一个阈值，使累积贡献率达到设定的阈值来确定l的取值。

本方法中我们利用信号的主成分构造n×l阶的主成分信号矩阵y＝[y1y2…yl]，然后，根据主成分信号矩阵与原始信号矩阵内积gji＝[yj,xi]构造矩阵g，

g＝y^t(x-e[x])(22)

进一步选取gj＝[gj1gj2···gjm],0＜j≤l作为漏水识别的特征。

基于时频特征和支持向量机的泄漏信号检测

上文提出的信号时频特征在辨识管道泄漏时具有不同的特点，但是，单一考虑某一特征时产生误判的概率较大。例如，如果泄漏信号与未泄漏信号的功率谱分布具有明显差异时，固有模态函数的功率谱均值具有较好的识别效果，但如果存在同频带干扰时则该方法容易出现误判。当管道漏点的泄漏量较小时，泄漏信号的近似熵均值与未泄漏信号的区别会不明显，这容易导致该方法出现误判。

为了提高泄漏检测的准确率，本发明综合利用以上提出的特征组合作为识别特征，采用支持向量机(svm)对信号特征进行分类，判断管道是否发生泄漏。svm是一种以统计学习理论为基础的数据挖掘方法，在解决小样本问题、非线性问题以及高维数据等问题中具有优势，广泛应用于数据预测、数据拟合、模式识别等领域中。假设训练集数据样本为(xi,yi)，其中1≤i≤n，每个样本xi∈r^d，d为输入空间的维度，yi∈{-1,1}为类别标签。如果该训练集可被一个超平面线性划分，则该超平面可以表示为w·x+b＝0，其中w和b是决定了超平面的位置。满足下面条件的样本称为支持向量，

yi(w·xi+b)＝1(23)

样本的最优化分实际上是最优分类超平面的求解问题，

式中，w是特征空间中分类超平面的系数向量；b是分类面的阈值；ξi是考虑分类误差而引入的松弛因子且ξi≥0；c是对于错分样本的惩罚因子。由此得到的最优分类超平面可以表示为，

w0·x+b0＝0(25)

利用拉格朗日乘子法，可以将式(24)的最优化超平面求解问题转化为它的对偶问题，

式中α＝(α1,…,αn)为拉格朗日乘子，满足αi＞0的样本为支持向量。利用式(26)求得的最优超平面可以得到其对应的决策函数，

其中，为式(33)最优解时的αi值。

对于非线性可分的情况，可以通过一个映射函数(在svm称为核函数)将低维的输入空间r^d映射到高维的特征空间h，使训练样本由低维的线性不可分问题转变为高维空间的线性可分问题。此时优化问题的对偶问题为：

其中，k(xi,xj)＝φ(xi)·φ(xj)为核函数。从式(28)可以看出，对于非线性可分问题需要选择合适的核函数k(·)构造svm模型。利用式(28)对应的决策函数为，

为了提高漏水检测的准确率，我们需要利用已知信号对支持向量机进行优化，优化过程如图5所示。考虑到地下供水管道所处环境因素的影响，需要在不同时间段对不同区域管道发生泄漏和未发生泄漏时进行信号采集构造信号样本库，从信号样本库中随机选取信号作为训练样本和测试样本。首先，利用训练样本的特征集对svm进行训练，形成初步的辨识模型。然后，利用测试样本的特征集对训练完成的svm模型进行测试。根据测试结果对svm做进一步的优化，直到测试输出的准确率满足要求，形成用于管道泄漏识别的svm模型。应用中，检测系统上位机将无线传感器网络采集的信号输入svm模型，根据模型输出的标签判断管道是否发生泄漏。

式(26)-(29)的理论分析表明，决定svm模型性能的因素为核函数和惩罚因子c。svm的核函数主要分为线性核、多项式核、sigmoid核和径向基核。根据漏水信号的特点和svm模型的泛化能力，本发明选择了径向基核作为核函数。径向基核函数的表达式为：

k(xi,xj)＝exp(-γ||xi-xj||²)(30)

因此，对于一个基于径向基核函数的svm，其性能由参数(c,γ)决定。为了使svm模型达到较高的识别效果，图5中的优化过程实际上是利用训练样本和测试样本对参数c和γ进行优化调整。大量研究表明，实际应用中c和γ的取指数序列可以达到好的效果，通常为c＝2^-5,2^-4,···,2¹⁵，γ＝2^-15,2^-14,···,2⁵。

本发明基于cross-validation的”grid-search”方法，对svm参数进行优化。

首先，根据参数的取值范围对参数进行取值，c＝2^x,x∈[-5,15]，γ＝2^y,y∈[-15,5]。然后，利用训练样本和测试样本对不同2^x和2^y组合下的svm模型进行测试，并输出测试的准确率。最后，选取最优cross-validationaccuracy的c和γ作为漏水识别svm的参数。

实施例1

实验中选择了pvc供水管道进行信号采集。在不同时段对泄漏和非泄漏情况各采集了100组数据，每组数据长度为5000，用于训练和优化支持向量机。同时，在相对寂静的清晨时段对泄漏和非泄漏情况各采集100数据，通过人为加入噪声的方法验证本发明研究的泄漏检测和时延估计算法的有效性。

首先，从不同时段采集的泄漏信号和非泄漏信号中各提取50组数据组成一个100组的训练样本。然后，利用剩余样本构造一个100组数据的测试样本。支持向量机的参数(c,γ)取值为2的整数次幂，并且，c的取值范围为c∈[2^-5,2¹⁵]，γ的取值范围为γ∈[2^-15,2⁵]。通过网格搜索法，采用21×21＝441个参数组合下的(c,γ)对模型训练，并使用测试集数据对模型性能检测，得到的检测准确率为如图6所示。图6的结果表明，本发明提出算法的最高识别准确率为98％。并且，可以看出当惩罚因子c≥2²和核函数参数γ≤20，且两个参数的乘积2¹≤c×γ≤2⁷时基于径向基核的svm模型对管道泄漏信号具有较好的识别性能。

图7为选择参数组合(c,γ)＝(2⁹,2^-4)时的检测结果，其中管道发生泄漏信号标签设为-1，管道未泄漏情况下的信号标签设为1。从识别结果中可以看出，本发明提出的方法仅将2组泄漏情况下的信号判别为未泄漏信号，其余情况均判断正确。表1为寂静段采集的泄漏信号经人工加入高斯噪声和脉冲噪声后利用本发明提出算法进行泄漏识别的情况。结果表明，本发明提出的基于信号多特征与支持向量机的供水管道泄漏检测方法能够有效检测管道的泄漏情况。

表1.高斯噪声和脉冲噪声环境中泄漏信号识别率

本文虽然已经给出了本发明的实施例，但是本领域的技术人员应当理解，在不脱离本发明精神的情况下，可以对本文的实施例进行改变。上述实施例只是示例性的，不应以本文的实施例作为本发明权利范围的限定。