半球谐振陀螺(hemispherical resonator gyro, HRG)是高精度的振动陀螺仪，其敏感器件全部用熔融石英加工而成，包括半球振子、激励罩和读出基座3个元件. 当HRG的敏感轴绕惯性空间转动时，HRG的谐振子在绕中心轴旋转时有哥氏效应产生，在哥氏效应下谐振子振型在环向相对壳体进动^[1]. HRG的输出会随着环境温度变化而变化, 进行温度建模补偿分析对提升HRG的输出精度意义重大. 周小刚等^[2]采用独立成分分析方法建立陀螺的温度误差补偿模型，并对与温度有关的确定性漂移进行补偿，同时使用Allan方差对不确定性漂移进行分析补偿. 周强等^[3]认为HRG谐振子频率随温度变化导致的输出变化是线性的，进而设计新的温度补偿电路，精准测试HRG的输出温度特性，将漂移由0.05°/h提高到0.03°/h. 吴宗收等^[4]提出将改进的粒子群优化(particle swarm optimization, PSO)算法和自回归滑动平均(auto-regressive moving average, ARMA)建模方法结合起来，经过HRG升温实验检验，参数寻优能力相比传统方法提高了一倍. 李广胜等^[5]将自回归(auto-regressive, AR)多变量模型引入HRG温度漂移模型的时间序列，借助系统辨识的相关理论，建立并分析了漂移和温度的多变量自回归模型. Li等^[6]提出基于多元回归的HRG温度误差补偿方法，建立了补偿温度误差的多元回归模型，并通过升温和降温实验验证了该模型的补偿效果. 虽然补偿效果明显，但多元回归模型对原始数据要求较高，需要数据具有较好的可重复性.

HRG的输出数据具有非线性和非平稳性，传统的时间序列分析模型只能捕捉线性关系，不适合HRG预测分析. 核主成分分析法(kernel principal component analysis, KPCA)和数据处理组合方法(group method of data handing, GMDH)是处理非线性数据的常用方法. 汪国新等^[7]将KPCA引入火力发电设备状态诊断与分析，提出多元时序分割算法. 实验结果表明，该算法能有效分析模型并准确提取非线性特征，可用于解决火力发电设备的状态诊断问题. 黄泽丰等^[8]为了提升光场成像的空间分辨率，通过多尺度潜在低秩分解图像，并引入KPCA对各层特征系数进行融合，使得基础层和显著层的特征信息更加明显，所生成的图像具有更高的分辨率和更好的视觉效果. 高运广等^[9]将KPCA模型引入非线性微惯性测量单元传感器故障诊断，在缩短参数选择时间和减少工作量的同时，使故障检测的准确率提高了18.44%. 段伟等^[10]将GMDH神经网络引入地震液化场地分析，建立侧向变形的预测模型，该模型的准确度和可靠性都明显高于其他神经网络. 凯立德·艾尔巴兹^[11]以盾构机为研究对象，提出集成GMDH神经网络和遗传算法的智能方法，预测了盾构机滚刀寿命；相比经验模型，该智能方法能够提供更高的准确度. Ahmadi等^[12]针对一类逆变器，提出攻击弹性模型预测自适应控制器，借助GMDH网络估计系统的不确定性使系统免受脉冲、缩放和随机攻击的影响. Xie等^[13]将GMDH神经网络应用于离群值检测模型，在有效消除离群值的同时，提高了训练后分类模型的分类精度. 针对光伏/燃料/电池系统的电压和功率调节问题，Band等^[14]设计了基于GMDH神经网络的控制系统, 考虑了可变未知动态、未知温度和辐照以及输出负荷突然变化等运行条件，该系统在几种场景的仿真中相比其他方法表现出较好的性能.

针对HRG温度建模与补偿问题，本研究1）提出基于KPCA-GMDH神经网络的建模补偿方法，分析HRG全寿命数据的大数据特征，初步选取输入特征向量；2）基于KPCA去除初选特征向量之间的冗余性3）将筛选后特征向量代入GMDH网络建立模型；4）进行单一样本和多样本的训练和预测，对比不同模型的预测效果以验证所提方法的性能.

在输入角速度为零时，陀螺的输出叫做零偏. 区别于传统机械陀螺，HRG的核心部件半球谐振子通过径向振动产生驻波，再由驻波的进动敏感输入角速度. HRG在工作时，环境温度的变化以及陀螺内部产热会影响谐振子的密度、半径、弹性模量和泊松比等参数，使得谐振子的固有频率发生变化，表现为陀螺的输出产生温度漂移，因此温度是引起陀螺漂移的重要因素之一. 零偏的温度漂移严重制约陀螺的精度，有必要对零偏的温度漂移进行抑制或补偿.

半球谐振子四波腹振动的谐振频率约为

(1)$ f = 1.512\;7\frac{h}{{{r^2}}}\sqrt {\frac{E}{{(1+\mu )\rho }}} . $

式中：h为半球谐振子的厚度，r为半球谐振子的半径，E为弹性模量，μ为谐振子材料的泊松比，ρ为谐振子的材料密度. 上述变量中弹性模量受温度影响较大，该参数随温度变化的表达式为

(2)$ E(T) = {E_0}[1 - {k_E}(T - {T_0})]. $

式中：E₀为取常温时谐振子的弹性模量，k_E为弹性模量随温度变化的系数，量级约为10⁻⁵. 式(2)是线性关系式，且E(T)与E₀非常接近，将式(2)带入式(1)后在E₀处进行泰勒展开即可获得较好的线性近似关系：

(3)$ \begin{split} f = &\frac{{1.512\;7h}}{{{r^2}}} \sqrt {\frac{{{E_0}}}{{(1+\mu )\rho }}} \left[1 - \frac{1}{2}{k_E}(T - {T_0}) - \right. \\ & \left.\frac{1}{8}k_E^2{(T - {T_0})^2}+ \cdots \right]. \end{split} $

k_E数值很小，保留一阶项即可获得理想的精度：

(4)$ f = \frac{{1.512\;7h}}{{{r^2}}} \sqrt {\frac{{{E_0}}}{{(1+\mu )\rho }}} \left[1 - \frac{1}{2}{k_E}(T - {T_0})\right]. $

当保留一阶项时，谐振子谐振频率f与温度T为线性关系，因此，可以选取谐振频率作为HRG的温度漂移补偿模型的自变量特征. 式(4)中，温度T表示谐振子材料的实时温度，而非真实的环境温度. 谐振子从环境吸热的过程类似于热工过程，传递函数用典型一阶惯性环节近似表示为

(5)$ {G_0}(s) = \frac{K}{{Fs+1}}. $

式(5)表明，谐振子谐振频率与环境温度存在非线性关系，在建立HRG温度漂移补偿模型时应当考虑；惯性环节表明谐振频率的导数项(离散变量为差分项或前一周期采样值)应当作为自变量特征.

工业大数据是以互联网大数据为基础，结合工业过程背景提出的新概念. 通过在时间上不断存储、积累过程运行数据，在空间上扩展、采集、运输数据，工业过程中遍布着不同尺度的时空间数据，以及散落于各级工业部门的不同来源、不同类别的数据. 工业大数据的主要特性是数据规模大，数据总类多，数据要求处理速度快，数据价值密度低，数据真实性低^[15]. 为了处理工业过程中的大数据，解决复杂的控制、优化和故障评估等问题，大数据建模成为热门的研究方向. 王晓军^[16]建立基于大数据的风洞试验段马赫数预测模型，提出特征子集集成方法、子模型学习算法和集成修建算法，有效提高了马赫数的预测速度及精度. 王龙晖^[17]以调节阀为研究对象，提出基于大数据驱动的调节阀故障诊断方法、调节阀阀后压力预测方法，有效提升了调节阀的智能化水平. 杨小佳^[18]将腐蚀大数据技术应用于低合金结构钢，并结合深度学习模型挖掘低合金结构钢内在腐蚀规律，验证了该技术的可靠性.

通过数据管理、分析和挖掘体系，可以得到HRG全寿命全周期测试的数据^[19]. HRG测试数据数量巨大，符合“数据规模大”的特性；数据的采集和建模须进行实时处理，符合“数据要求处理速度快”的特性；每组数据中均包含冗余数据且数据价值实现程度低，符合“数据价值密度低”的特性；数据里通常包含有噪声数据，符合“数据真实性低”特性. 因此，可以认定HRG测试数据为工业大数据.

在神经网络中，特征选择的好坏直接影响预测效果. 根据特征和因变量的相关程度，特征分为无关、弱相关和强相关3种，特征之间可能存在关系的相互重叠，根据特征间的相关程度，分为冗余和非冗余特征. 特征选择的作用是，消除原有集合中嘈杂、冗余和不相关的特征，优化训练模型，提高算法搜索效率，提高分类精度^[20].

设定测试采集的${\mathrm{HRG}} $谐振频率为$ {{\boldsymbol{X}}_f} = [{f_1},{f_2},\cdots, {f_{\mathrm{N}}}] \in {{{{\bf{{{R}}}}}}^{N \times 1}} $, 其中$ f_i, i = 1, 2,\cdots, N $为${\mathrm{HRG }}$谐振频率，N为样本个数. 一方面，谐振频率${\boldsymbol{X}}_f $数量级较大，而${\mathrm{HRG}} $的温度漂移输出值与谐振频率差分值${\mathrm{d}}{\boldsymbol{X}}_f $数量级都很小，从实际特性分析，各自变量特征阶数不应超过二阶，否则实际系数的量级远小于或大于其他自变量特征的系数，对应的自变量变化值对模型的影响可以忽略. 另一方面，$ X_f \cdot {\mathrm{d}}{\boldsymbol{X}}_f $与温度漂移输出值的数量级相接近，须考虑各自变量特征之间的交叉耦合项. 综合前文的分析，初步考虑以$ {\mathrm{d}}{\boldsymbol{X}}_f $、$ \sqrt{{\boldsymbol{X}}_f} $、${\boldsymbol{X}}_f $、常值偏置项及其各自的耦合项作为模型的自变量特征，形成12维特征输入向量$ {\boldsymbol{X}}_{12} $：

(6)$ \begin{split} {{\boldsymbol{X}}_{12}} = &[1,{\rm{ }}\sqrt {{{\boldsymbol{X}}_f}} ,{\rm{ }}{{\boldsymbol{X}}_f},{\boldsymbol{X}}_f^2{\rm{, }}{\mathrm{d}}{{\boldsymbol{X}}_f},{\rm{ }}\sqrt {{{\boldsymbol{X}}_f}} \cdot {\mathrm{d}}{{\boldsymbol{X}}_f},{\rm{ }}\\&{\rm{ }}{{\boldsymbol{X}}_f} \cdot {\mathrm{d}}{{\boldsymbol{X}}_f},{\rm{ }}{\boldsymbol{X}}_f^2 \cdot {\mathrm{d}}{{\boldsymbol{X}}_f},{\rm{ }}{({\mathrm{d}}{{\boldsymbol{X}}_f})^2},{\rm{ }}\\& {\rm{ }}\sqrt {{{\boldsymbol{X}}_f}} \cdot {({\mathrm{d}}{{\boldsymbol{X}}_f})^2},{\rm{ }}{{\boldsymbol{X}}_f} \cdot {({\mathrm{d}}{{\boldsymbol{X}}_f})^2},{\rm{ }}{{\boldsymbol{X}}_f^2} \cdot {({\mathrm{d}}{{\boldsymbol{X}}_f})^2}].\end{split} $

初选的非线性特征向量之间存在相关性，不利于后续计算，为此，基于KPCA进一步筛选特征向量. 通过非线性映射$ \varPhi (\cdot) $将训练数据${\boldsymbol{X}} \in {{{{{{{{{\bf{R}}}}}}}}}^{N \times M}}$投影到新的高维空间${\boldsymbol{F}}$，其中$N$是样本数量，M为变量数，再对新的数据矩阵$ \varPhi ({\boldsymbol{X}}) \in {{\mathrm{{\bf{R}}}}^N} \times {\boldsymbol{F}} $进行线性主成分分析，得到加载向量${\boldsymbol{p}} \in{\boldsymbol{ F}}$，使得得分向量${\boldsymbol{t}} = \varPhi ({\boldsymbol{X}}){\boldsymbol{p}}$的方差最大化^[21]. 优化过程的表达式为

(7)$ {{\mathrm{max}} _{\boldsymbol{p}}}{\text{ }}\frac{{{{\boldsymbol{t}}^{\mathrm{T}}}{\boldsymbol{t}}}}{{N - 1}} = {{\mathrm{max}} _{\boldsymbol{p}}}\frac{{{{\boldsymbol{p}}^{\mathrm{T}}}{\varPhi ^{\mathrm{T}}}({\boldsymbol{X}})\varPhi ({\boldsymbol{X}}){\boldsymbol{p}}}}{{N - 1}}. $

式中：$ {{\boldsymbol{p}}^{\mathrm{T}}}{\boldsymbol{p}} = 1 $, $ {\boldsymbol{X}} = {\left[ {{{\boldsymbol{x}}_1},{{\boldsymbol{x}}_2},\cdots,{{\boldsymbol{x}}_N}} \right]^{\text{T}}} $，存在系数向量${\boldsymbol{\alpha}} = {\left[ {{{{\alpha}} _1}{\text{,}}{\alpha _2}{\text{,}}\cdots{\text{,}}{\alpha _N}} \right]^{\mathrm{T}}} \in {{\bf{R}}^{{N}}}$满足

(8)$ {\boldsymbol{p}} = \sum\limits_{j = 1}^N \varPhi \left( {{{\boldsymbol{x}}_j}} \right){\alpha _j} = {\varPhi ^{\mathrm{T}}}({\boldsymbol{X}}){\boldsymbol{\alpha}} . $

组合式（7）和式（8），KPCA优化变为

(9)$ {{\mathrm{max}} _{\boldsymbol{\alpha }}}\frac{1}{{N - 1}}{{\boldsymbol{\alpha }}^{\mathrm{T}}}\varPhi ({\boldsymbol{X}}){\varPhi ^{\mathrm{T}}}({\boldsymbol{X}})\varPhi ({\boldsymbol{X}}){\varPhi ^{\mathrm{T}}}({\boldsymbol{X}}){\boldsymbol{\alpha }}. $

式中: $ {{\boldsymbol{\alpha }}^{\mathrm{T}}}\varPhi ({\boldsymbol{X}}){\varPhi ^{\mathrm{T}}}({\boldsymbol{X}}){\boldsymbol{\alpha }} = 1 $. 为了避免定义非线性映射$\varPhi ( \cdot )$，引入核矩阵${\boldsymbol{K}} = \varPhi ({\boldsymbol{X}}){\varPhi ^{\mathrm{T}}}({\boldsymbol{X}})$，${\boldsymbol{K}}$的第$(i,j)$个元素计算式为

(10)$ {\boldsymbol{K}}(i,j) = {\varPhi ^{\mathrm{T}}}\left( {{{\boldsymbol{x}}_i}} \right)\varPhi \left( {{{\boldsymbol{x}}_j}} \right) = \ker \left( {{{\boldsymbol{x}}_i},{{\boldsymbol{x}}_j}} \right). $

典型的核函数是高斯核$ \ker ( {{{\boldsymbol{x}}_i},{{\boldsymbol{x}}_j}} ) = \exp ( { - {{\| {{{\boldsymbol{x}}_i} - {{\boldsymbol{x}}_{{j}}}} \|}^2}/\sigma } ) $，$\sigma > 0$为核宽度；多项式核$ \ker \left( {{{\boldsymbol{x}}_i},{{\boldsymbol{x}}_j}} \right) = = {\left( {{\boldsymbol{x}}_i^{\mathrm{T}}{{\boldsymbol{x}}_j}+{d_0}} \right)^{{d_1}}} $，其中$ {d_0} $、$ {d_1} $为多项式核参数. 求解式（9）可以得到特征值分解为

(11)$ {\boldsymbol{K}}{\boldsymbol{\alpha}} = (N - 1)\lambda {\boldsymbol{\alpha}} . $

当${\text{ }}\tilde N \leqslant N$时，式（11）的解带来非零特征值${\lambda _1} \geqslant {\lambda _2} \geqslant \cdots{\lambda _{\tilde N}}$和对应的特征向量${{\boldsymbol{\alpha}} _i},1 \leqslant i \leqslant \tilde N$. 所有的特征向量${\boldsymbol{A}} = \left[ {{{\boldsymbol{\alpha }}_1}{\text{,}}{{\boldsymbol{\alpha }}_2}{\text{,}}\cdots{\text{,}}{{\boldsymbol{\alpha }}_{\tilde N}}} \right] \in {{\bf{R}}^{N \times \tilde N}}$被分成2个部分，其中前$K$个特征向量（对应的特征值不小于平均特征值）定义核的主成分子空间(principal component subspace, PCS)，其余$\tilde N - K$个特征向量构成核的残差子空间(residual subspace, RS). KPCA模型训练计算内核的PCS，${\boldsymbol{T}} = {\boldsymbol{KA}} \in {{\bf{R}}^{N \times \tilde N}}$为$ {\boldsymbol{X}} $的非线性特征. 测试向量${{\boldsymbol{x}}_t} \in {{\bf{R}}^M}$的第$ i $个核得分计算式为

(12)$ {t_{t,i}} = {\varPhi ^{\mathrm{T}}}\left( {{{\boldsymbol{x}}_t}} \right){{\boldsymbol{p}}_i} = {\varPhi ^{\mathrm{T}}}\left( {{{\boldsymbol{x}}_t}} \right){\varPhi ^{\mathrm{T}}}({\boldsymbol{X}}){{\boldsymbol{\alpha }}_i} = {\boldsymbol{k}}_t^{\mathrm{T}}{{\boldsymbol{\alpha }}_i}. $

式中：${{\boldsymbol{k}}_t} = {\left[ {{k_{t,1}},{k_{t,2}},\cdots,{k_{t,N}}} \right]^{\mathrm{T}}} = \varPhi ({\boldsymbol{X}})\varPhi \left( {{{\boldsymbol{x}}_t}} \right) \in {{\bf{R}}^N}$，当${\text{ }}1 \leqslant j \leqslant N$时，${k_{t,j}} = \ker \left( {{{\boldsymbol{x}}_j},{{\boldsymbol{x}}_t}} \right)$. 累计方差贡献率(cumulative percent variance，CPV)准则根据主元方差的累计和百分比来确定主元个数. 把测量向量${\boldsymbol{x}}$协方差矩阵的前A个特征值的和除以所有特征值的和称为前A个主元的累计贡献率^[22]. 前A个主元的累计贡献率为

(13)$ {\text{CPV}} = \left.{{\displaystyle\sum\limits_{i = 1}^{{A}} {{\lambda _i}} }} \right/{{\displaystyle\sum\limits_{i = 1}^M {{\lambda _i}} }}. $

式中：$ {\lambda _i},{\text{ }}i = 1,2,\cdots,M $为对应的特征值，通常认为当CPV达到85%时，前A个主元包含了原数据足够多的信息. 初步的特征向量用式（6）表示，再根据KPCA和CPV准则计算，即可得到最终的特征向量.

GMDH网络是以KPCA-GMDH多项式为基础的复杂非线性系统辨识方法^[23]. 假设非线性系统以$ {x_1}, {x_2},\cdots,{x_n} $为输入，$y$为输出，关系为

(14)$ y = g\left( {{x_1},{x_2}, \cdots ,{x_n}} \right). $

函数$g$的离散沃尔泰拉(Volterra)级数展开式（KPCA-GMDH多项式）为

(15)$ \begin{split} y = &{a_0}+\sum\limits_{i = 1}^n {{a_i}} {x_i}+\sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n {{a_{ij}}} } {x_i}{x_j}+ \\& \sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n {\sum\limits_{k = 1}^n {{a_{ijk}}} } } {x_i}{x_j}{x_k}+ \cdots .\end{split} $

被广泛用于非线性模型的完全描述. 当$n$较大时，多项式的维度很高，因此完全确定$ {a_0},{a_i},{a_{ij}},\cdots $等参数的数值不现实. GMDH通过局部简单模型不断逼近式（15），得到近似模型如图1所示. 图中，部分多项式$G$为2个输入变量的完全多项式；$ {\theta _0},{\theta _1},\cdots,{\theta _m} $分别为每层的模型评价准则；$ {y_{ij}} $为由部分模型计算得到的输出，部分模型由拟合实测数据辨识得到；中间变量$ {x_{ij}} $是从$ {y_{ij}} $中按每层的模型评价准则值选出来的，将作为下一层的输入. 将最终筛选的特征向量带入GMDH神经网络，训练给出的输出值为$ {\hat y_i} $，即

(16)$ {\hat y_i} = \hat g\left( {{x_{i1}},{x_{i2}}, \cdots ,{x_{in}}} \right);\;\;\; i = 1,2, \cdots ,n. $

$ {\hat y_i} $与实际输出${y_i}$往往存在偏差，残差平方和为

(17)$ e = \sum\limits_{i = 1}^N {{{\left[ {\hat g\left( {{x_{i1}},{x_{i2}}, \cdots ,{x_{in}}} \right) - {y_i}} \right]}^2}} . $

将$ e $最小化便得到GMDH神经网络结果.

如图2所示为基于KPCA-GMDH神经网络模型的建模补偿方法的基本步骤. 1）分析HRG数据的大数据特征，初步考虑以${\mathrm{d}}f$、$\sqrt f $、$f$、常值偏置项及其各自的耦合项作为模型的自变量特征，形成12维特征输入向量${{\boldsymbol{X}}_{12}}$. 2）基于KPCA去除初选特征向量之间的冗余性，借助累计方差贡献率确定筛选后的特征向量${{\boldsymbol{X}}_A}$, A为筛选后的主元个数. 3）将实验数据分为训练集和验证集. 通过训练集确定网络权值和网络结构，完成完整KPCA-GMDH神经网络模型的构建.

以HRG的多次输出数据为分析对象. 数据采集方案：将HRG安装在平台单轴回路中，测量地纬度为${29^ \circ }35'$，如图3所示. 实验输入量为地速竖直分量，输出量为平台框架输出角. 数据采集从陀螺通电瞬间开始，连续采集数据6~10 h，采样频率为0.04 Hz，采集参数包括陀螺转角、半球谐振子谐振频率. 共进行7次数据采集，每次间隔24 h.

实验共测得7组数据，以时间$ t $为横坐标，分别以谐振频率$ f $和陀螺转角$ \theta $为纵坐标，绘制原始数据曲线如图4所示.

在测试过程中，陀螺仪静止置于稳定平台上，因此将转台运动从实际测量输出中扣除. 考虑到转台近似以匀角速度转动，以一阶线性模型对输出值进行拟合，再由陀螺输出值减去拟合值获得实际输出值，扣除转台运动前、后HRG输出如图5所示.

平滑处理也称模糊处理，是图像处理领域常用的操作. 平滑样条曲线拟合方法是常用的平滑方法之一，它基于样条曲线和平滑系数来拟合曲线，能够最佳适应曲线中的数据点^[24]. 如图6所示，原始数据存在“阶梯”和大量震动，推断为量化噪声，可以借助均值滤波函数使数据曲线更加光滑.

扣除转台运动后，将7次陀螺样本输出值依次解析. 查看原始数据发现，样本1存在500 s的未启动状态，其余样本存在40±3 s的未启动状态. 将未启动状态截除后，所有样本的输出统一按0起始时刻依次画出，如图7所示. 此时样本1不仅在输出时长上与其余样本不同，用于分析的初始段数据丢失较多，且整体输出值也明显异于其余样本，存在大量粗糙数据. 为了提高整体数据的处理精度，综合考虑后剔除样本1.

分析陀螺的角速度$ \omega $. 数据量较大不利于后续处理研究，因此进行数据的重采样. 原始数据采样间隔为40 ms，现将每1 000个点取均值融合成1个点，即对原始数据进行重构. 每个样本的数据长度不一，为了方便后续处理，将末端的数据截断，取前480 min的数据，如图8所示.

选取样本后，通过初选特征向量和KPCA确定最终特征向量，随机取数据的80%作为训练集，剩余20%作为验证集，带入GMDH网络训练，并与多项式模型进行对比. 以样本2为例，其初选特征向量、经KPCA筛选后的特征向量、HRG的转动角速度的统计特性分别如表1、表2、表3所示.

GMDH网络训练集的输入为

(18)$ {\bf{Tin}} = [{{\boldsymbol{X}}_1},{{\boldsymbol{X}}_2},{{\boldsymbol{X}}_3},{\boldsymbol{Y}}]. $

式中：${\boldsymbol{Y}}$为样本2陀螺的转动角速度. 通过训练确定神经网络模型后，带入验证集数据. GMDH网络验证集的输入为

(19)$ {\bf{Vin}} = [{{\boldsymbol{X}}_1},{{\boldsymbol{X}}_2},{{\boldsymbol{X}}_3}]. $

GMDH网络验证集的输出为

(20)$ {\bf{Vot}} = [\hat {\boldsymbol{Y}}]. $

预测效果的评价指标选择均方根误差RMSE，计算式为

(21)$ {\mathrm{RMSE}} = \sqrt {\frac{1}{N}\sum\limits_{i = 1}^N {{(\hat Y_i - Y_i)}^2}} . $

式中：$N$为验证集数据个数. 最终预测结果如图9所示，其中$ \omega $为陀螺的角速度，$ \Delta \omega $为陀螺角速度预测误差，其余单一样本的训练结果如表4所示. 由表可知，不同训练样本下，KPCA-GMDH神经网络预测结果的RMSE均小于多项式预测结果的RMSE，说明前者的预测效果优于后者.

进一步验证KPCA-GMDH神经网络的鲁棒性，将多样本进行融合，其中样本A由原始样本2、3、5、6、7融合而成，样本B由原始样本3、5、6、7融合而成，样本C由原始样本2、5、6、7融合而成，样本D由原始样本3、5、6融合而成. 从实验数据预处理结果来看，样本2和样本7数据趋势较其余4个样本不同，因此推断，去掉样本2和样本7可以提升网络精度. 分别就融合样本A、B、C和D参与训练的4种情况进行实验，考虑到样本4在样本2、3、4、5、6和7中最为稳定，将样本4作为测试组，将KPCA-GMDH神经网络与多项式、GMDH网络进行对比，评价指标为误差标准差：

(22)$ \sigma = \sqrt {\frac{1}{{N - 1}}\sum {{{({\varDelta _i} - \overline \varDelta )}^2}} } . $

式中：$\varDelta $为预测方法的预测误差. 实验结果如图10~图13和表5所示. 由表可知，KPCA-GMDH模型相比传统多项式模型精度分别提升了48.5%、54.0%、56.3%、68.4%，相比GMDH模型分别提升了3.6%、5.1%、3.8%、8.8%，说明所提模型适用于HRG建模补偿，且降噪效果良好.

以HRG大数据为研究对象，本研究提出基于KPCA-GMDH神经网络的建模补偿方法，实验并分析不同样本不同条件下HRG的输出精度，验证模型的有效性. 基于HRG大样本数据建立KPCA-GMDH神经网络模型. 实验表明，单一样本预测时，所提方法的预测效果明显好于传统多项式；在多样本预测时，所提方法的预测精度相比多项式模型分别提升了48.5%，54.0%，56.3%，68.4%；相比GMDH模型分别提升了3.6%、5.1%、3.8%、8.8%，说明所建模型适用于HRG建模补偿，且降噪效果良好. KPCA-GMDH神经网络模型会受训练样本的输入自变量(频率)的范围影响，测试样本的频率在训练样本的频率范围内就能够预测出较好的结果. 当新输入的测试数据不在训练集的频率范围内时，预测结果会出现较大的偏差，在后续工作中，计划通过补充实验扩大样本的方式，尽可能扩大训练样本的频率范围. 无论在何种训练模式下，KPCA-GMDH神经网络模型预测的HRG数据的标准差均在千分之一量级，基本满足零偏稳定性和快速启动性的要求.