加速试验缺失数据的样条插值分析

徐志昆

（贵州商学院经济与金融学院，贵阳 550014）

0 引言

数据缺失现象普遍存在于各应用领域的统计抽样中。机械原因、人的主观失误、历史局限、有意隐瞒是数据缺失的主导原因。样本数据的缺失在很大程度上影响分析建模和研究结果的准确性。

长寿命、高可靠性是科技进步带来的产品发展趋势，而传统寿命试验周期过长，时间成本巨大，加速试验便成为解决这个问题的极佳选择。但过去的大部分加速试验研究都是在没有数据缺失的前提下进行，若仅因为数据存在缺失就重作试验是极不现实的，代价非常高昂。

目前针对数据缺失主要有加权法、插补法和构造特殊统计模型三类处理方法。对于单元无回答的情况多采用加权法对缺失值进行补救，而对于项目无回答的处理多采用插补的方法，包括单一插补和多重插补，及根据工程背景构造特殊统计模型。Politz 等（1949）提出了经典的Politz-Simmons调整法［1］。Dempster等［2］首次提出一种使得不完全数据得到有效估计的EM 算法。Rubin 等［3］基于EM 算法，率先提出了多重插补方法。Brick 等［4］提出了最近邻插补方法，也即树枝分类的距离函数匹配方法。Liu 等［5］在1994年进一步提出了ECME 等。1998 年金勇进［6］探讨了处理缺失数据中对辅助信息的利用问题。2002 年Schafer 等［7］提出的极大似然估计和多重填补法具有较好的处理精度和较广的应用范围。王乃生等［8］给出恒定应力加速寿命试验中数据缺失时的统计方法。2009 年刘宝慧［9］利用回归插补给出方差分析。杨贵军等［10］对高相关辅助变量择优回归插补法进行了研究。但从20 世纪90 年代初至今，缺乏处理缺失数据的全新思想［11-12］。

针对加速试验数据一般随时间呈现单调变化趋势且精度高的特点，立足于单一插补的角度，提出插值填补法，分别运用Lagrange 插值、三次样条插值和B 样条插值，给出缺失数据合理的替补值，达到减小估计量偏差和数据集完整的目的。

1 基于样条插值原理的缺失值处理

三次样条插值曲线具有良好的性质，在实用中最为普遍。设Δ是［a,b］的一个划分，则：

若函数S（x）满足：

（1）S(x) ∈C2[a,b]（注：Ck[a,b]表示区间［a,b］上具有k阶连续导数的函数集）；

（2）S(xi)=f(xi),i= 0,1,…,n；

（3）S（x）在每个子区间[xi,xi+1](i=0,1,…,n- 1)上都是次数不超过三次的多项式，且至少在一个子区间上为三次多项式。则称S（x）为关于划分Δ的三次样条函数。

提出使用三转角法、三弯矩法和B样条基函数法完成石英摆片加速退化试验的缺失值处理。

1.1 基于三转角法的插值处理

设S（x）在节点xi(i= 0,1,…,n)处的一阶导数值为S'(xi)=mi其中mi是待定参数。记

则有方程组

第二边界条件S"(a) =f"(a),S"(b) =f"(b),有方程组：

第三边界条件m0=mn,m1=mn+1，得方程组：

其中：

1.2 基于三弯矩法的插值处理

选择二阶导数作为待定参数：

三弯矩法基本方程［13］：

其中：

在实际应用中，若三次样条插值没有边界条件，最常用的方法就是采用非扭结条件，即：

再由三弯矩基本方程，可得：

利用上述两种方法均可解出mi(i=0,1,…,n)后，分别代入

即得插值函数S（x），用于插值计算。

1.3 基于B样条基函数法的插值处理

B 样条曲线具有局部性，控制顶点只影响部分曲线的形状，对其余部分不产生影响，比较具有稳健性，且其造型灵活，还可进行统计数据的光滑化处理。

设有控制顶点P0,P1,…,Pn，则p阶（p-1次）B样条曲线的数学表达式为

其中：Ni,p(μ)是p-1 次B 样条曲线的基函数。B样条基函数是一个称为节点矢量的非递减的参数μ的序列所决定的p阶分段多项式，也即为p阶（p-1次）多项式样条。

B样条de Boor-Cox递推定义：

给定一组数据{Qk}(k= 0,1,…,n)，找一条p次B 样条曲线顺序通过这组数值点［14］，即是根据数据点分布情况选定一组合适的节点矢量U=和控制顶点确定p次B样条曲线

2 加速试验的退化数据及结果分析

2.1 2010-11-15号石英摆片在加速应力85℃下试验的等效挠度数据

石英摆片（2010-11-15）在加速应力85℃下试验的等效挠度数据（单位：10-1mm）如表1 所示，为验证方法可行性，设定空白处为缺失数据。实际观测到T1、T2、T6 时刻正面等效挠度分别为7.3856、7.3949、7.3670，T3、T8、T9 时刻反面等效挠度分别为7.3960、7.3900、7.3930。

在实验中，分别用Lagrange 线性和三次插值、三转角和三弯矩插值法、均匀和非均匀B样条基函数法插值进行内推，得到插补结果如表2 所示，并计算误差平方和（SSE）如表3 和4所示。

插补效果如图1、图2所示。

从石英摆片加速退化试验的两组数据来看，通过比较图像和误差平方和，Lagrange 插值、三转角法和非均匀B 样条在两组缺失数据的插补中均取得了较高的精度，效果理想。

表2 插补结果

图1 正面等效挠度插补效果

表4 反面等效挠度误差平方和

2.2 2010-01-09号石英摆片在加速应力85℃下试验的等效挠度数据

石英摆片（2010-01-09）在加速应力85℃下试验的等效挠度数据（单位：10-1mm）如表5 所示，仍设定空白处为缺失数据。实际观测到T9、T10 时刻正面等效挠度分别为7.6390、7.6479，反面等效挠度分别为7.6291、7.6359。

表5 石英摆片（2010-01-09）等效挠度

由于插值法一般在外推时精度不高，甚至可能会发生龙格现象，造成巨大偏差。所以在外推缺失数据时，借鉴均值插补方法，将外推转化为内推处理。步骤如下：

（1）利用已知数据均值来代缺失值相邻的下一时刻数据；

（2）利用插值法计算缺失数据；

（3）重复步骤（1）、（2）直到所有缺失值计算完成。

得到插补结果如表6所示，并计算误差平方和（SSE）如表7和表8所示。

表6 插补结果

外推插补效果如图3、图4所示。

图3 正面等效挠度插补效果

图4 反面等效挠度插补效果

表7 正面等效挠度误差平方和

表8 反面等效挠度误差平方和

从加速试验数据端点缺失的插补结果来看，借助均值插补方法把外推转化为内推可以避免端点处的巨大波动，降低外推风险，且又在一定程度上反映了数据自身的变化趋势。通过比较图像和误差平方和，Lagrange 插值、三转角法和非均匀B 样条在两组缺失数据的外推插补中得到了较高的精度，效果较好。

3 结语

数据缺失是统计工作中普遍存在的现象。掌握数据缺失的处理方法，有助于在进行数据采样、统计分析等环节减少、规避重要信息的丢失，达到提高分析精度的目的。

通过采用Lagrange 插值、样条插值的单一插补方法来研究加速试验缺失数据，发现可用于填补缺失数据的中间点，且相比之下算法简单、容易实现，特别当数据点呈较强规律变化时效果更好。在石英摆片加速退化试验的缺失数据插补中，插值法从数据点自身变化趋势出发，并在外推过程中，借鉴均值插值思想把外推转化为内推，得到了较高精度的缺失数据。

但每一种插值填补方法都不是普遍适用的，都只是对缺失数据分析的一种尝试。尤其对端点的缺失数据填补应进一步研究。在分析具体问题时，应该综合权衡考虑使用一种或者几种方法的综合结果。