白永昕 ，田茂再 ,2a，2b

（1.兰州财经大学 统计学院，兰州 730020；2.中国人民大学a.应用统计科学研究中心；b.统计学院，北京 100972）

0 引言

函数型方差分析是函数型数据分析的对象之一。方差分析的目的在于研究诸多因素在不同状态下对因变量的取值是否有显著影响。因此，对函数型的研究对象在某方面的差异进行统计分析时往往需要用到函数型方差分析。函数型方差分析模型拟合和参数估计的研究相对较多[1-3]。在方差分析中，若控制变量确实对观测变量的显著影响存在，那么本文需要进一步研究控制变量的不同水平对观测变量的影响程度。单纯的通过两独立样本的t检验会增大犯弃真的概率，对此可通过多重比较检验解决[4-7]。

甘肃是一个地理分布极其特殊的省份。李栋梁[8]根据气候区划分指标将甘肃省分为8个气候区。从东南到西北包括了北亚热带湿润区到高寒区、半干旱区和干旱区的各种气候类型。年平均气温在0～14℃之间,但是各个气候带的气温差异很大。为此，甘肃省不同气候带对气温的影响很有研究的必要。本文对函数型方差分析模型的参数估计进行介绍，并将多元中方差分析的多重比较检验应用到函数型方差分析中，分析了甘肃省不同气候带对气温的影响。

1 模型及方法

1.1 模型

一般地，假定控制变量A有m个水平，每个水平均有n个样本。那么在水平Aj下的第i次连续的观测结果为：

其中，μ(t)为观测变量总的期望均值函数；aj(t)'j=1'…'m是控制变量水平 j对观测变量产生的附加影响,称为水平 j对观测变量产生的效应,且是残差函数，假定每个εij(t)都是均值为零的独立的高斯过程。

定义一个(mn)×(m+1)的设计阵Z,每一行代表单独的一个观测。符号(ij)表示第 j个水平下第i个观测相应的行,此行的第一列是1,第 j+1列是1,剩余全为0。z(ij)k是设计阵Z中对应第(ij)行第k列的值。同时，定义一个相应的函数型向量 β1=μ'β2=a1,以此类推 βm+1=am,因此得到函数型向量β=(μ'…'am),则模型（1）可以表达为：

假设响应变量是离散数据平滑的结果,因此在实际应用中,一般直接用原始观测数据Y来估计函数型参数。为了简化,其模型（2）的向量表达形式为：

其中，Y是包含mn个气温的函数型向量,β(t)是m+1维的函数型参数向量,ε(t)是m+1维的函数型残差向量。

1.2 参数估计

1.2.1 参数B的估计

下面开始讨论模型中未知函数型参数β(t)的估计。对上述函数型线性模型中的β(t)用基函数展开,则有：

其中,B是5×Kβ的系数矩阵,θ(t)是长度为Kβ的基函数向量,Kβ是基函数展开的个数。对于B的惩罚函数,引入一个线性算子L来定义：

关于微分算子的计算将在下文给出。

如果模型是标准的一般线性模型，β(t)则可以通过标准的最小二乘估计得到。将最小二乘方法拓展到函数型数据，只需以适当的方式重新定义残差平方和。此时Y-ZBθ是一个函数型的向量，由方程（3）、方程（4）和方程（5）得到惩罚最小二乘为：

运用Kronecker积最小化上式可以得到B的估计值为：

1.2.2 微分算子的估计

定义线性微分算子为：

与之对应的微分方程为：

其中，C=(C0'…'CP-1)表示的是微分算子系数；θ=(1'sin(wt)'…)T是傅里叶基函数。

令W=(θ'…'Dp-1θ)，则方程（6）可以表达：

求解以上方程得到向量C的估计值Ĉ为：

由上式则可得到：

进一步，R的估计值R̂为：

2 方差分析

为了研究水平Aj对观测变量产生的效应，上文中估计了不同水平下的观测变量曲线。接下来需要进一步分析控制变量在固定的时间t上是否对观测变量函数产生了显著影响。为此开始进行函数型方差分析。与普通的方差相同，函数型方差分析需要考虑残差平方和函数LMSE和F比值函数FRATIO。不同于标准多元情况的是在特定时间这些量是相互独立的。在上述模型中,残差平方和函数为：

同时，令：

类似与单因素方差分析，可以进一步计算残差的均方函数为：

MSE=SSE/df(error)

其中,df(error)是误差的自由度。

同样,回归的均方函数是SSY(t)与SSE(t)的差值除以回归的自由度。因此：

其中,df(regression)表示两个模型中误差自由度的差。最后,可以得到F检验函数为：

3 假设检验

方差分析只能判断控制变量是否对观测变量产生了显著影响。如果控制变量确实对观测变量产生了显著影响，那么感兴趣的是：控制变量的不同水平对观测变量的影响程度如何；其中哪个水平的作用明显异于其他水平；哪个水平的作用不显著等。尽管这些问题可以通过两独立样本的t检验解决，但是进行多次比较会增大犯第一类错误的概率。对此，本文通过多元方差分析中的多重比较检验来解决。

Tukey的HSD检验是较为常用的多重比较检验。它将控制变量k个水平下观测变量的总体均值做两两对比检验。原假设为第i个总体与第 j个总体的均值函数无显著差异。即：

H0:μ1=…=μm

H1:μ1≠…≠μm

检验统计量为q统计量。第i个总体与第 j个总体均值对比检验的检验统计量定义为：

其中，LMSE为残差平方和函数；ni和nj分别为第i和第 j水平下的样本量。不同于LSD检验的t统计量，q统计量近似服从有n-k个自由度的学生化极差q分布。假设检验的显著性通过q统计量的p值判断。

4 实证分析

4.1 数据来源

数据来自中国气象局提供的天水（站号57006）、宕昌（站号56095）、兰州（站号52889）、靖远（站号52895）、平凉（站号53915）、正宁（站号53935）、张掖（站号52625）、金昌（站号52675）、武威（站号52679）、敦煌（站号52418）、碌曲（站号56071）、卓尼（站号56082）和舟曲（站号56094）共13个地面气候观测站常规4次观测得到的2015年逐日日最高温度、日最低温度。日平均气温由每日四个时点气温值求均值得到。每个站所有缺失数据均少于7天,所缺失的数据由相邻前后5天气温的平均值代替。

4.2 模型

根据甘肃省特殊的地理特征,本文把甘肃省划分为陇南暖温带湿润区、陇中温带半干旱区、河西暖温带干旱区和甘南高寒湿润区四个不同的气候带。同时,在这些气候带中选出13个可靠性较高的气象观测站,它们分别是陇南暖温带湿润区的天水和宕昌站、陇中温带半干旱区的兰州、靖远、平凉和正宁站、河西暖温带干旱区的张掖、金昌、武威和敦煌站以及甘南高寒湿润区的碌曲、卓尼和舟曲站。

每个气候带有不同数量的观测站，假设第g个气候带中第m个气温函数的模型为：

Ymg(t)=μ(t)+ag(t)+εmg

其中,Ymg表示第g个气候带中第m个气温观测站的温度函数。均值函数μ(t)表示甘肃省平均气温，即所有气温观测站的平均温度。ag(t)表示第g个气候带对气温的影响函数，并且满足残差函数εmg表示第g个气候带中第m个气温观测站不可解释的变化。

4.3 方差分析

从上述模型出发，使用R软件进行数据处理，得到了甘肃省4个气候带的气温曲线图（见图1）以及不同气候带中总体均值的95%置信带（见图2）。

图1 实线表示不同气候带的气温曲线，虚线表示甘肃省平均气温

图2 不同气候带中总体均值的95%置信带

从图1中可以看到：陇南地区的全年气温均高于甘肃省平均气温。春季和秋季的气温与均值的差异较小，在夏季的七月份则差异较为明显。冬季气温普遍高出全省均值4℃左右，冬季的差异最明显。陇中地区全年的温度最接近甘肃省平均气温。只有冬季气温略高于平均气温，但差异较小。河西走廊地区相比于全省平均气温则变化较为多样性。春季和秋季气温略低于平均气温，夏季则高于平均气温，但是与平均气温的差异都不太明显。到了冬季，气温则明显下降，比平均气温要低3～4℃。甘南地区的气温变化情况与河西走廊地区恰恰相反。甘南地区在春、夏、秋三个季节的气温都低于全省平均气温。其中，春季和秋季与均值的差异较小，在夏季达到最大，比平均气温要高出3～4℃。冬季气温则要高于平均气温，并且差异相对较小。

方差分析的F比值函数如图3所示。

图3 函数型方差分析的F比值函数

从图3方差分析的F比值函数可以看出，F比值均高于5%显著水平下的2.92，表明4个气候带对气温的影响有显著的差异。

以上分析结果可以表明4个气候带对气温的影响有显著的差异，但是并不能表明每两个气候带之间的差异都显著，也不能表明哪些气候带的影响异于其他气候带的影响，哪些的影响不显著等。因此，本文用单因素方差分析中的HSD多重比较检验进行两两气候带之间的比较。同时，将整个函数型观测变量的值域分成春季（3月至5月）、夏季（6月至8月）、秋季（9月至11月）和冬季（12月至次年2月）四个区间。多重比较检验的结果见表1所示。在显著性数值上用“*”表示均值之间的对比在0.05的置信水平上有显著差异。

表1 HSD多重比较表

多重比较检验的结果分析如下：在春季，陇南和陇中（p=0.3578）、陇南和河西走廊（p=0.1018）以及陇中和河西走廊（p=0.9147）的对比都没有显著差异。陇南和甘南（p=0）、陇中和甘南（p=0.0069）以及河西走廊（p=0.0480）的对比出现了显著差异。在夏季，除了陇中和河西走廊（p=0.9986）的对比无显著差异。其余气候带的两两对比都存在显著差异。同样在秋季，显著差异出现在陇南和河西走廊（p=0.3895）、陇南和甘南（p=0.0051）以及河西走廊和甘南（p=0.0002）的对比中。冬季同夏季相似，除了陇中和甘南（p=0.8677）的对比无显著差异，其余气候带的两两对比均存在显著差异。

最后，感兴趣的是在整个区间上气候带对气温影响的成对比较（见图4）。从图中可以看到：陇南和陇中对气温的影响只有在夏季（p=0.0161）和冬季（p=0）有显著差异；陇南和河西走廊对比，除了春季的其他季节均有显著差异；陇南和甘南则在四个季节均出现了显著性的差异；而陇中与河西走廊对比，只有冬季有显著差异；相比之下，陇中和甘南的显著性差异出现在春季（p=0.069）和夏季（p=0）；类似于陇南和甘南，河西走廊和甘南在四个季节也都出现了显著差异。

图4 不同气候带显著性水平的两两对比

4.4 分析结果

分析结果可以概括如下：

从季节来看,秋季气温受气候带的影响较小,一些地区秋季气温几乎等于平均气温；春季气温受气候带的影响次之；夏季和冬季气温受气候带影响最大。在夏季，由于受到气候带的影响，甘南的气温要比全省平均气温高5℃左右，陇南的气温也要高于全省平均气温2～4℃。在冬季，除了甘南的其他地区气温均受到气候带的影响。陇南冬季的气温比平均气温高5℃左右，而河西走廊冬季气温则相反，比平均气温低5℃左右。总体来看，气温在秋季和春季的波动较小，在夏季和冬季的波动相对较大。

从空间变异来看,陇中地区全年气温受气候带的影响最小。除了冬季气温稍有波动，全年气温曲线几乎与全省平均气温曲线重合。河西走廊地区次之，除了冬季气温有较大波动，全年气温与平均气温相似。陇南地区和甘南地区最大。由于受到气候带的影响,陇南地区和甘南全年气温较平均气温有较大的波动，而且全年气温差异性较大。总体来看,甘肃省西北部受气候带影响较大,东南部则较小。

5 结论

本文详细介绍了函数型方差分析的参数估计和多重比较检验方法，并通过对实例进行分析，给出了函数型方差分析的方法和结果。将多元方差分析中的多重比较检验应用到函数型方差分析中，得到不同气候带两两之间显著性差异的情况，从而确定甘肃省不同气候带对气温的具体影响，为甘肃省气温的研究提供参考。

[1]Wahba G,Klein B.Smoothing Spline ANOVA for Exponential Families,With Application to the Wisconsin Epidemiological Study of Diabetic Retinopathy[J].Annals of Statistics,1995,23(6).

[2]Stone C J,Hansen M H,Kooperberg C,et al.Polynomial Splines and their Tensor Products in Extended Linear Modeling[J].Annals of Statistics,1997,25(4).

[3]Huang J Z.Projection Estimation in Multiple Regression With Application to Functional Anova Models[J].Annals of Statistics,1998,26(1).

[4]Tallarida R J,Murray R B.Newman-Keuls Test[M].Manual of Pharmacologic Calculations,1987.

[5]Keselman H J,Burt H,Cribbie R A.Multiple Comparison Procedures[M].New York:Wiley,1987.

[6]Tuncer D.The Influence of Blood and/or Hemostatic Agent Contamination on Micro-TBS to Dentin[J].Journal of Adhesion Science&Technology,2015,29(13).

[7]Wang T L,Tseng Y K.Do Thinking Styles Matter for Science Achievement and Attitudes Toward Science Class in Male and Female Elementary School Students in Taiwan?[J].International Journal of Science&Mathematics Education,2015,(13).

[8]李栋梁.中国夏季月平均气温异常研究[J].高原气象,1995,14(2).