[0073] 以下结合附图对本发明作进一步说明。
[0074] 本发明方法包括以下各步骤：
[0075] (1)利用GJ-5型轨道检测车分别安装在车轴与车厢上的垂直振动加速度计获得车轴和车厢位置的时域振动加速度信号a1(t)和a2(t)，其幅值单位为G(重力加速度，9.8m/s2)，其中a1(t)∈[-0.2,0.2]，a2(t)∈[-15.8,15.5]，GJ-5型轨道检测车运行时速为100千米/小时～150千米/小时，两个加速度计的振动信号均为每隔h米同时采样一次，满足0.2m≤h≤0.3m，共计采集T次，1000≤T<∞，则采样时刻t＝1,2,…,T。
[0076] (2)将步骤(1)中获取的时域振动信号a1(t)和a2(t)进行短时傅里叶变换，获取每个采样时刻的频域频谱，其中设置短时傅里叶变换中窗函数的窗口宽度为τ，且满足20≤τ≤25，通过短时傅里叶变换后得到每个时刻窗口各个频率的幅值，并求每个频率幅值平方的平均值，该平均值即为对应频谱的平均功率，将其作为每个采样时刻所获取振动时域信号a1(t)和a2(t)对应的振动频域特征f1(t)和f2(t)。
[0077] (3)从GJ-5型轨道检测车上获取每个采样时刻轨道高低不平顺的幅值Y(t)，其单位为毫米：
[0078] 在GJ-5型轨道检测车上获取每个采样时刻转向架垂直振动加速度时域信号、惯性基准测量值，以及列车倾角信号之后，利用GJ-5型轨道检测车所携带的数据处理系统中的惯性参考测量算法，可以从这些信号数据中计算出轨道高低不平顺的幅值Y(t)，将f1(t)、f2(t)和Y(t)表示成向量p(t)＝[f1(t),f2(t),Y(t)]，共计可以得到T个向量，它们组成的向量集记为P＝{p(t)|t＝1,2,…,T}。
[0079] (4)建立置信规则库(缩写为BRB)，用其反映车轴及车厢处振动频率特性变量f1和f2(BRB的二维输入)与高低不平顺的幅值变量Y(BRB的输出)之间的非线性关系，其中，BRB的第k条规则记为Rk，其表示形式如下：
[0080] Rk:If f1is AND f2is THEN Y is
[0081]
[0082] Rk的规则权重为θk，满足0≤θk≤1；输入变量f1和f2对应的属性权重分别为δ1,δ2，且0≤δ1,δ2≤1；
[0083] 式(1)中， 和 分别为BRB的输入变量的f1和f2参考值，且有其中j＝1,2，Qj为 的取值空间，其中的元素满足
mj表示对应第j
个输入变量参考值的取值个数，mj≥1；分别在Q1,Q2中抽取一个元素作为f1、f2的参考值，由此组合成规则，共计可以产生L＝m1×m2条规则，L≥1，k＝1,2,3,…,L为规则的编号；
[0084] 式(1)中，Rk后项属性分别为D1,D2,…,DN，并有LY≤D1＜D2＜…＜DN≤RY，N≥2，β1,k,β2,k,…,βN,k分别为D1,D2,…,DN的信度值，并满足0≤βi,k≤1， i＝1,2,…,N；
[0085] 其中，式(1)中，设定初始规则权重为θk＝1，初始属性权重δj＝1。
[0086] 为便于理解，举例说明，假设 中各参考值的取值为A1,1＝8，A1,2＝9，A1,3＝10，中各参考值的取值为A2,1＝0.1，A2,2＝0.2，A3,3＝0.3，假设满足Y＝f1×f2，D1＝0.5,D2＝1.4,D3＝2.2,D3＝3，共计将会产生L＝9条规则，设初始θk＝1,δ1＝δ2＝1，其中部分规则形势如下：
[0087] R1:若f1＝8且f2＝0.1,则[(D1,0.667),(D2,0.333),(D3,0),(D4,0)]；
[0088] R2:若f1＝8且f2＝0.2,则[(D1,0),(D2,0.75),(D3,0.25),(D4,0)]；
[0089] R3:若f1＝8且f2＝0.3,则[(D1,0),(D2,0),(D3,0.75),(D4,0.25)]；
[0090] ……
[0091] R7:若f1＝10且f2＝0.1,则[(D1,0.444),(D2,0.556),(D3,0),(D4,0)]；
[0092] R8:若f1＝10且f2＝0.2,则[(D1,0),(D2,0.25),(D3,0.75),(D4,0)]；
[0093] R9:若f1＝10且f2＝0.3,则[(D1,0),(D2,0),(D3,0),(D4,1)]；
[0094] 这里组成了9条规则，其中βi,k为满足约束根据历史数据分析得到的结果。
[0095] (5)给定振动频率特性f1和f2后，通过置信规则库推理获取它们对应的轨道高低不平顺幅值估计结果 具体步骤如下：
[0096] (5-1)设定f1和f2的取值分别为 和 标I表示BRB系统的输入，并有将它们带入BRB模型，计算它们激活各个规则的权重：
[0097]
[0098] 其中，wk∈[0,1]； 为相对属性权重，表达式为：
[0099]
[0100] 式(2)中， 表示为第k条规则中第j个输入变量相对于参考值 的匹配度(c＝1,2,…,mj)，匹配度的求解方法如下：
[0101] (a)当 和 时， 对于Aj,1和 的匹配度 取值均为1，对于其他参考值的匹配度均为0；
[0102] (b)当 时， 对于Aj,c和Aj,c+1的匹配度 取值分别由式(4)和式(5)给出
[0103]
[0104]
[0105] 此时，输入变量对应的其他参考值的匹配度均为0；
[0106] 步骤(5-2)获得输入变量为 和 时，模型推理后的不同后项输出的信度融合值[0107]
[0108] 其中，
[0109]
[0110] 步骤(5-3)获得输入变量为 和 时轨道高低不平顺幅值估计结果
[0111]
[0112] 为便于理解，举例说明，以步骤(4)中的模型为例，假设模型输入 和此时的真实值Y＝1.118，带入公式(2)-(5)，可得激活规则R1、R2、R4与R5，且：
[0113] 分别抽取两个变量的匹配度组合成计算得到：
[0114] w1＝0.28，w2＝0.12，w4＝0.42，w5＝0.18，其余均为0；
[0115] 所以，将结果带入式(6)-(7)可以获得 带入式(8)，可得 与真实值偏差0.046，这
里为提高精度可提高规则数和后项输出的个数。
[0116] (6)选择训练优化样本集PT
[0117] 随机从步骤(3)中给出的向量集P＝{p(t)|t＝1,2,…,T}，从中抽选出TN个向量组成训练样本集PT＝{p(t)|t＝1,2,…,TN}，TN≥500，这里选择的训练样本集能够尽量激活所有规则，并将样本集中向量的前两维 的取值作为BRB模型的输入，按照步骤(5)得出它们的估计值
[0118] (7)确定BRB非线性优化模型
[0119] 步骤(7-1)确定优化参数向量
[0120] V＝(θk,δj,βi,kk＝1,2,…,L,j＝1,2,i＝1,2,…,N) (9)
[0121] 将各优化参数组成向量V＝[v1,v2,…,vTn]，向量是由规则库中规则权重、属性权重和后项信度组合而成，Tn为优化参数的个数，Tn>2这里，Tn＝L+2+L×N，故满足0≤va≤1，a＝1,2,…,Tn，且对应va为βi,k时满足，
[0122] 步骤(7-2)建立优化目标函数为：
[0123]
[0124] (8)获得基于序列线性规划的BRB模型，具体步骤如下：
[0125] 步骤(8-1)将非线性模型中的目标函数进行一阶泰勒展开
[0126]
[0127] 其中，V0为初始的BRB模型对应的变量值，ξ(V0)表示将V0带入式(10)对应的函数值，并令ξ(V0)＝obj0， 表示给定模型初始值V0时目标函数对应各优化参数va的一阶偏导数，ξ(0)(V)表示非线性规划模型近似的线性规划模型的目标函数，其中，向量V为未知参数向量。由此，可将非线性模型简化为线性规划问题：
[0128] minξ(0)(V)
[0129] s.t.0≤va≤1,a＝1,2,…,Tn (12)
[0130] 且对应va为βi,k时，满足
[0131] 步骤(8-2)确定优化参数向量V中各参数va的寻优区间
[0132] (a)确定各参数va的上限向量UB，其中参数βi,k对应的上限为：
[0133]
[0134] 规则权重及属性权重的上限UB2均为1，由向量UB2与UB1组成向量UB，其中Tk表示在T样本数据中，同时激活第k条规则的样本数据集W对应时刻t组成的集合；
[0135] 为便于理解，举例说明，假设激活第3条规则的Y组成的数据集Q＝{2.01,2.015,2.01,2.01,2.985,2.045,3.015,3.01,3.015,3.09,3.065}，此时t的序列号组成的集合为Tk＝{92,322,367,387,409,411,444,445,446,447,456}，选择对应下标最大的Y_max＝
3.09，取D＝[0,2,4,6,8,10,12]，得到k＝3时的7条对应的上限值UB1(βi,3)＝{1,1,1,1,
0.89,0.712,0.5933}，同理获得所有优化参数对应的上限UB。
[0136] (b)确定线性规划中的移动限move_lim：
[0137] 设置上限的10％作为最初选取的移动限范围，move_lim＝[lower,upper]，其中，lower表示优化参数的移动下限，upper表示优化参数的移动上限，
[0138]
[0139] 这里，tx用于缩小移动限范围，初始tx＝0.5，0.5≤tx≤20；
[0140] 将移动限确定的取值范围与式(12)中给定的各参数的取值范围取交，确定出最终优化参数的寻优区间；
[0141] 步骤(8-3)获得近似线性规划的局部最优解V_yh
[0142] 根据线性规划理论，结合步骤(8-2)确定的各优化参数的寻优区间，在寻优区间内寻找参数的局部最优解V_yh，最常用的两种线性规划算法为单纯形法与内点法；
[0143] 步骤(8-4)判定优化参数的结果是否满足设计要求
[0144] 将步骤(8-3)中线性优化结果V_yh带入式(10)，得到优化参数后的模型对应的目标函数值obj1；
[0145] 如果obj1≥obj0，说明线性规划的结果不如初始模型的结果，此时，tx的值加1，重新带入步骤(8-2-b)，通过缩小移动限的方式，缩小寻优区间继续寻找最优值，tx>20时说明变量的移动限变化不大，停止搜索，重新赋值tx＝0.5，并输出此时模型参数V_yh；
[0146] 如果obj1
[0147] 如果|obj1-obj0|>err，将优化后的V_yh赋给V0，obj1的值赋给obj0，带入步骤(8)重新进行循环迭代，直到|obj1-obj0|≤err，停止迭代，并输出结果V_yh；
[0148] 得到的训练优化结果V_yh组成的模型，即为训练优化后的轨道高低不平顺装置的BRB模型。
[0149] 以下结合附图，详细介绍本发明方法的实施例：
[0150] 本发明方法的流程图如图1所示，核心部分是：首先选取有一定特征的历史数据进行短时傅里叶变化，得到安放在轨检车车轴及车厢上加速度计的时域振动信号对应的频域特征，构建反映输入频域特征信号f1和f2与输出轨道高低不平顺幅值Y之间非线性关系的置信规则库，然后通过选取的训练样本，将BRB模型中各参数组成的非线性优化问题，通过对目标函数的泰勒展开将其转化为更为简单易求的线性规划问题，最后利用逐层线性化的方式近似出非线性优化的解。
[0151] 以下结合我国某既有干线下行区段(1584.5103km～1586.86735km)为例，轨检车以100千米/小时的速度运行，每隔0.25米采集一次轨道的相关参数信号。详细介绍本发明方法的各个步骤，并通过实验结果验证本发明提出的一种基于专家规则推理与序列线性规划的轨道高低不平顺幅值估计方法与利用MATLAB工具箱中fmincon函数相比，在计算效率上的优势，避免fmincon函数中可能存在优化不出结果的情况。
[0152] 1.利用GJ-5型轨检车中分别安装在车轴与车厢上的垂直振动加速度计获得车轴和车厢位置的时域振动加速度信号a1(t)和a2(t)，选择GJ-5型轨检车运行时速为100千米/小时，两个加速度计的振动信号均为每隔0.25米同时采样一次，共计采集T＝(1586.86735-1584.5103)÷(0.25*10-3)＝9428次。
[0153] 2.根据步骤(2)将时域振动信号a1(t)和a2(t)进行短时傅里叶变换，获取每个采样时刻的频域频谱，其中设置短时傅里叶变换中窗函数的窗口宽度为τ＝21，频率特征f1(t)和f2(t)，如图2所示。
[0154] 3.从GJ-5型轨检车上计算出输入f1(t)和f2(t)对应各时刻的轨道高低不平顺的幅值Y(t)，并由此组成向量集P＝{p(t)|t＝1,2,…,9428}。
[0155] 4.建立反映车轴及车厢处振动频率特性变量f1和f2(BRB的二维输入)与高低不平顺的幅值变量Y(BRB的输出)之间的非线性关系的置信规则库(BRB)。
[0156] 选取输入输出变量的语义值，f1的模糊语义值描述为、很小(very small,VS1)、正小(positive small,PS1)、正中(positive medium,PM1)、正大(positive large,PL1)、很大(very large,VL1)；f2的模糊语义值描述为：极小(exceeding small,ES2)、很小(very small,VS2)、正小(positive small,PS2)、正中(positive medium,PM2)、大(large,Z2)、中大(medium large,ML2)、很大(very large,VL2)，极大(exceeding large,EL2)。其参考值如表1-表3所示：
[0157] 表1 f1的语义值和参考值
[0158]
[0159] 表2 f2的语义值和参考值
[0160]
[0161] 表3 Y的语义值和参考值
[0162]
[0163] 建立的初始置信规则库如表4所示，其中后项输出的信度值根据历史数据按照要求给定：
[0164] 表4初始置信规则库
[0165]
[0166]
[0167] 5.结合步骤(6)得到训练样本给定频率特性 与 通过置信规则库推理获取它们对应的轨道高低不平顺幅值估计结果
[0168] 在9428组向量中挑选TN＝500组训练样本，其对应输入的真实值为Y(t)，训练样本输入经过BRB推理机制对应的估计值为
[0169] 首先，确定每个输入变量相对于各参考值的匹配度。对于500个训练样本计算出它们对于各自参考值(语义值)的匹配度。例如，当t＝7500时，输入量为f1＝0.3369、f2＝0.0062，则f1对于PS1和PM1的匹配度分别为0.8215和0.1785，f2对于PS2和PM2的匹配度分别为0.8418和0.1582，对其余参考值的匹配度均为0。
[0170] 然后，计算规则的激活权重。获得输入量对于每个规则中参考值的匹配度 后，即可利用式(2)计算BRB中每条规则的激活权重wk。同样以t＝7500时为例，对于输入变量f1、f2，可以得到其对规则R11～R12和R19～R20的激活权重分别为w11＝0.6915,w12＝0.13,w19＝0.1503,w20＝0.0282，而其他规则的激活权重均为0，亦即激活了BRB中的4条规则。
[0171] 最后，融合激活的规则。利用ER算法融合被激活的置信规则后项的置信结构，得到关于Y的输出置信结构。根据式(3)与式(4)可以求出 例如，将上步中关于f1＝0.3369、f2＝0.0062的wk及初始BRB后项置信度βj,k带入式(6)与式(7)中，可得O＝{(D1,0),(D2,0.4369),(D3,0.5631),(D4,0),(D5,0),(D6,0),(D7,0)}。
[0172] 即，
[0173] 6.确定BRB非线性优化模型
[0174] 得到非线性模型的优化参数指标向量V：
[0175] V＝(θk,δj,βi,kk＝1,2,…,L,j＝1,2,i＝1,2,…,N)
[0176] 根据上面的实例，规则库中规则权重共有L＝40个、属性权重2个和后项信度280个参数，由此组成共计Tn＝322个优化参数的非线性规划模型为：
[0177]
[0178] 且 k＝1,2,…,40,j＝1,2,i＝1,2,…,7
[0179] 7.获得基于序列线性规划的BRB模型，其中，序列线性规划的算法流程如图3所示，具体步骤如下：
[0180] 7.1将非线性模型中的目标函数进行一阶泰勒展开
[0181] 这里，以表4的中的各模型参数组成的向量V0为非线性模型的初始值，V0对应的目标函数值为obj0＝0.2599，将目标函数值代入式(11)，进行泰勒展开由此得到线性规划模型：
[0182] minξ(0)(V)
[0183] s.t.0≤θk≤1,0≤δj≤1,0≤βi,k≤1
[0184] 且 k＝1,2,…,40,j＝1,2,i＝1,2,…,7
[0185] 7.2确定待训练的参数集合V中各参数变量va的寻优区间
[0186] 利用式(13)得到βi,k对应的参数上限UB1，这里以同时激活第3条规则为例，在输入的500组样本中，激活第3条规则的Y组成的数据集Q＝{2.01,2.015,2.01,2.01,2.985,2.045,3.015,3.01,3.015,3.09,3.065}，此时t的序列号组成的集合为Tk＝{92,322,367,
387,409,411,444,445,446,447,456}，选择对应下标最大的Y_max＝3.09，得到k＝3时的7条对应的上限值UB1(βi,3)＝{1,1,1,1,0.89,0.712,0.5933}，同理获得所有优化参数对应的上限UB。
[0187] 根据式(14)以上限为UB1(β6,3,β7,3)为例，再根据初始值V0(β6,3,β7,3)，确定2个优化参数的初始移动限的范围为： 结合变量0≤βi,k≤1的要求，所以，变量搜索的寻优区间变为：
[0188] 7.3获得近似线性规划的局部最优解V_yh
[0189] 利用内点法，对线性规划模型求解局部最优，并利用MATLAB工具箱中的linprog函数实现线性寻优：
[0190] V_yh＝linprog(Der,A,b,Aeq,beq,LB,UB,V0,OPTIONS)；
[0191] Der表示初始模型各参数向量V0相对于目标函数的一阶偏导数向量，A,b为满足移动限要求的不等式约束A·V_yh≤b的矩阵，Aeq,beq为满足对应βi,k参数中， 的等式约束，LB,UB为满足每个参数0≤va≤1最大最小值要求组成的向量，V0为线性规划的初始值由此在可行域内搜索出满足线性要求的V_yh值。
[0192] 7.4判定结果是否满足设计要求
[0193] 将步骤(8-3)中线性优化结果V_yh带入式(10)，得到优化参数后的模型对应的目标函数值obj1＝0.1652；如果obj1≥obj0，则tx自增1，带入式(14)重新计算移动限，在改变后的约束范围内继续寻找局部最优V_yh；
[0194] 而此时obj1
[0195] 表5基于SLP方法优化后的置信规则库参数
[0196]
[0197]
[0198] 8.确定优化后的轨道高低不平顺幅值检测的BRB系统，进行结果的验证测试[0199] 根据以上步骤获得优化后的轨道高低不平顺检测系统中的参数，利用随机抽取的样本数据进行验证。这里，可以获得初始BRB对轨道高低不平顺幅值的估计曲线与优化后的BRB对轨道高低不平顺幅值估计曲线如图4所示，图4(a)表示初始BRB对轨道高低不平顺幅值的估计曲线，图4(b)表示优化后BRB的模型对相同输入变量时轨道高低不平顺幅值的估计曲线，其中由轨检车获得的数据作为真实值用“—”表示，本发明方法给出的轨道高低不平顺幅值估计值用“.”表示，图5表示初始和优化后BRB对轨道高低不平顺幅值的估计分别与真值之间的绝对误差，其中“—”代表初始BRB对轨道高低不平顺幅值估计与真值之间的绝对误差曲线，“--*--”代表优化后的BRB对轨道高低不平顺幅值估计与真值之间的绝对误差曲线。
