电子质量
2021年第04期(总第409期)
作者简介院肖哲坤(1999-)，男，湖北武汉，本科，主要研究方向为数学；朱洒洒(1999-)，女，河南周口，本科，研究方向为经济学。
变系数模型的估计方法及应用
Estimation Method and Application of Variable Coefficient Model
肖哲坤1，朱洒洒2(1.兰州大学数学与统计学院，甘肃兰州730107;2.兰州大学经济学院，甘肃兰州730000)
Xiao Zhe-kun 1,Zhu Sa-sa 2(1.School of mathematics and statistics,Lanzhou University,Gansu Lanzhou 730107;2.School of Economics,Lanzhou University,Gansu Lanzhou 730000)
摘要：该文提供了对连续型变系数模型的一种估计方法和求解思路。通过对连续型变系数模型的系数函数进行泰勒展开，把变系数模型问题转化为高维常系数模型来研究。为了防止异常点的影响，选择稳定的损失函数，并加入惩罚函数对高维数据进行变量选择，然后使用模拟退火算法和AIC 准则反复进行求解。使
用LAE+SCAD 估计方法对汽车mpg 数据集进行实证分析，体现出模型及算法的实用价值。关键词：变系数模型；稳健估计；惩罚函数；模拟退火算法；AIC 准则中图分类号：O212
文献标识码：A
文章编号：1003-0107(2021)04-0012-07
Abstract:This paper provides an estimation method and solution idea for continuous variable coefficient model.The problem of variable coefficient model can be transformed into high dimensional constant coefficient model by Taylor expansion of coefficient function of continuous variable coefficient model.In order to prevent the influence of outliers,a stable loss function is selected,and a penalty function is added to select variables for high-dimensional data.Then simulated annealing algorithm and AIC criteria are used to solve the problem repeatedly.This paper uses the estimation method of LAE+SCAD to make an empirical analysis on the vehicle mpg data set,which reflects the practical value of the model and algorithm.
Key words:Variable coefficient model;Robust estimation;Penalty function;Simulated annealing algorithm;AIC Criterion CLC number:O212
Document code:A
Article ID ：1003-0107(2021)04-0012-07
0引言
在高维数据的回归分析中，众多学者致力于寻求一种结构简单、便于估计、容易解释的模型。变系数模型就是这样一类非常重要的非参数回归模型，
不仅可以避免"维数祸根"问题，而且兼具解释性和适应性强的特点[1-2]。
目前，变系数模型的参数估计方法主要有核光滑估计、多项式样条估计和光滑样条估计三种[3-4]，但是这些方法都建立在最小二乘法之上，
而最小二乘法对异常值比较敏感从而会导致误差，需要对变系数模型进行稳健
估计。赵培信、
薛留根(2011)将B 样条和taut string [5]结合得出了变系数模型的一个稳健估计过程。Feng 等(2015)提出了变系数模型中基于稳健样条估计的选择方法；黄绿斓等(2016)在自然立方样条的基础上对变系数模型稳健估计的7种不同方法进行比较判断。李晓亮、陈艳(2019)研究了因变量随机缺失的部分线性变系数模型在解释变量和误差项含有异常点时的稳健估计方法。Tibshirani(1996)提出Lasso 方法之后，学者们就
热衷于研究惩罚的变量选择问题；而高维数据的产生也使得学者更加关注变系数模型研究时的变量选择问题[6-7]。
12
本文将连续型变系数模型的可变系数通过泰勒展开转化为常系数模型，在这一基础上同时考虑其稳健估计和变量选择问题，
采用"损失函数+惩罚函数"的模式，利用模拟退火算法和AIC 准则的方法进行求解。最后对用汽车mpg 数据集进行实证分析，检验模型和估计方法的可行性和准确率。
1模型建立
变系数模型由Hastie 和Tibshirani 于1993年提出[8]。假定线性回归模型中的系数是以指标变量U 为自变量的未知函数，变系数模型的一般形式为：
Y=β0(U )+β1(U )X 1+…+βp (U )X p +ε
(1)
其中，Y =(y 1,y 2,…,y n )T 是因变量，X j =(x j 1,x j 2,…,x jn )T (j =1,2,…,p )是p 个变量，X =(1,X 1,X 2,…,X p )是设计矩阵，β(U )=(β0(U ),β1(U ),…,βp (U ))T 是函数项系数，ε=(ε1,ε2,…,εn )T 是模型误差，独立同分布且满足E (ε|U ,X 1,X 2,…,X p )=0，Var(ε|U ,X 1,X 2,…,X p )=σ2。
当U 是连续变量，即U 的取值范围为连续点集T 时，式(1)中的β0(U ),β1(U ),…,βp (U )是关于U 的连续函数，此时为连续型变系数模型。我们通过对连续型变系数模型的系数函数进行泰勒展开，把变系数模型问题转化为高维常系数模型来研究。
设每一项系数βj (U )都关于U 光滑，对βj (U )(j =0,1,2,…,p )作泰勒展开，可得：
βj (U )=βj0+βj 1U +βj 2U 2
+o (U 2
)
(2)
其中，βj0,βj 1,βj 2是待定常数。于是，泰勒展开后(1)式变为如下形式：Y =(β00+β01U +β02U 2
)+(β10+β11U +β12U 2
)X 1+…
+(βp 0+βp 1U +βp 2U 2
)X p +ε
(3)
进一步，因为X j ,UX j ,U 2
X j (j =1,2,…,p )线性无关，所以自变量为X 1,X 2,…,X p 的变系数模型转换为了自变量为U,U 2,X 1,UX 1,U 2X 1,…,X p ,UX p ,U 2X p 的常系数线性回归模型：
Y =β00+β01U +β02U 2+β10X 1+β11UX 1+β12U 2X 1+…+βp 0X p +βp 1UX p +βp 2U 2
X p +ε
(4)
仍然可以利用常系数线性回归模型的理论来进行求解。
2模型求解
2.1传统求解方法的局限性
为了求解系数β(U )，传统的方法是最小二乘估计(OLS)，即求解使得模型平均残差平方和最小的系数：
β^
(U )=arg min β(U )
Q (β(U ))
=arg min β(U )[1
n
n i =1
∑(y i
-β0
(U )-p
j =1
∑βj
(U )x
ij )2
]
(5)
通过求导，得到β的无偏估计β^
=(X T X )-1X T Y 。但是当p 很大时，会导致两个问题：一是rank(X T X )≤n<p +1，从而X T X 不满秩，也就不可逆，无法代入公式求解；二是rank(X )≤n<p +1，从而X 不列满秩，即X 的列向量线性相关，变量X 1,X 2,…,X p 之间存在多重共线性，可能造成回归结果不稳定。另一方面，式(1)容易形成过拟合，即模型复杂度过高，在训练集上拟合效果好，在测试集上拟合效果差。
解决该问题需要对数据进行降维，主要有两种方法。一种方法是主成分分析(PCA)，通过正交变换将一
组可能存在线性相关的自变量转化为一组线性无关的自变量；另一种方法是正则化估计，通过添加惩罚函数，控制模型的复杂度，将变量选择和对回归系数的估计同时进行。
除此之外，由于OLS 的基本思想是使残差平方和最小，当回归数据中存在异常点时，残差会被"平方"放大，从而对回归结果产生影响，造成回归结果的不稳定。因此，我们需要寻一些对异常点的敏感性小、抗干扰能力强的稳健的估计方法[9-11]。2.2改进的求解过程
第一步：
将OLS 换为稳健的损失函数损失函数是衡量模型真实值和预测值不一样程度的一个函数。稳健回归的基本思想是在OLS 的基础上，将残差平方和这一损失函数替换为更加稳健的其它形式的损失函数，对不同的点施以不同的权重，残差小的点权重大，残差大的点权重小。将损失函数记为L (Y,X ,β(U ))。
第二步：加入惩罚函数
为了防止过拟合现象的发生，我们需要对自变量X 1,X 2,…,X p 进行变量选择。在损失函数后加入惩罚函数，记为p (λ,β(U ))，其中，λ为非负的惩罚调节参数。λ越大，为使目标函数最小，损失函数就要越小，从而越多的系数βj (U )会被压缩为0，这说明相应变量X j (j =1,2,…,p )在模型中不显著，可以被剔除，从而起到变量选择的效果。
至此，我们得到了求解式(1)的模式："损失函数+惩罚函数"，也就是说，我们需要优化目标函数：
min Q (λ,β(U ))=L (Y,X ,β(U ))+p (λ,β(U ))(6)
第三步：利用模拟退火算法求解系数
现在的问题是如何估计式(6)中的λ和β(U )。在给定λ的情况下，对于训练集，我们考虑使用模拟退火算法来求解β(U )。模拟退火算法是一种贪心算法，
但是它在13
电子质量
2021年第04期(总第409期)
搜索过程中引入了随机因素。该算法以一定的概率来接受一个比当前解要差的解，因此有可能跳出局部最优解，达到全局最优解。这里"一定的概率"来源于固体退火原理。根据Metropolis 准则，假设在温度t ，系统由状态x old 变为状态x new ，相应的系统能量由E (x old )变为E (x new )，
则接受新状态的概率为：
p =
1,
若E (x new )≤E (x old )
e
-E (x )-E (x )
t
，若E (x new )>E (x old )
⎧⎩
⏐
⏐⏐⏐⎨⏐⏐⏐⏐(7)
通过MATLAB 运行程序得到β(U )，即可得到该λ下的模型系数解。算法流程图如图1
所示。
图1模拟退火算法流程图
伪代码如下：
步骤1：t :=10,βj :=0(j =0,1,2,…,p ),β'j :=0(j =0,1,2,…,p ),res :=Q (β)；
步骤2：Δj :=random[-0.5,0.5](j =0,1,2,…,p ),βj :=βj +Δj (j =0,1,2,…,p )；
步骤3：判断Q (β)≤res ?成立则转步骤6，否则转步骤4；
步骤4：rd=random[0,1]；步骤5：判断e 100000(res-Q (β))
t
≥rd ?成立则转步骤7，否
则转步骤8；
14
步骤6：res:=Q(β)；
步骤7：β'j:=βj(j=0,1,2,…,p)；
步骤8：t:=0.9996t；
步骤9：判断t≥10-15?成立则转步骤2，否则输出β'j(j =0,1,2,…,p)并结束。
第四步：寻合适的λ
在给定λ的情况下，我们可以利用模拟退火算法求解β(U)。式(6)中λ增大时，模型方差减小，但误差增大，因此我们需要设法到合适的λ，平衡模型的方差和误差。我们需要引入一种衡量λ优劣的标准，这里考虑使用AIC准则。AIC准则由日本统计学家赤池弘次于1974年提出，是衡量统计模型拟合优良性的一种标准。它建立在熵的概念上，可以权衡所估计模型的复杂度和此模型拟合数据的误差。
在第三步中，我们已经求得给定λ后，训练集样本数据的模型参数β(U)。为了进行变量选择，可以设定合适的阈值，令小于给定阈值的βj(j=0,1,2,…,p)等于0。对于测试集样本数据，计算AIC的值：
AIC(λ)=1n
n
i=1
∑(y i-β0(U)-p j=1∑βj(U)x ij)2+2d(8)
其中，d是βj(j=0,1,2,…,p)中非零的个数，βj是用训练集样本数据求解得到的模型参数。比较不同λ下的AIC值来判断λ的优劣，AIC越小代表λ越好。我们使用网格化的方法寻最优的λ。先取定充分大的λmax，使得所有的βj(j=0,1,2,…,p)都小于给定的阈值，此时最优的λ必然落在(0,λmax)中。然后取定合适的较小步长，在(0,λmax)中以该步长均匀取点，选取其中AIC值最小的λ作为近似的最优λ，该λ及其对应的β(U)即为式(6)的求解结果。
3实证分析
下面我们将LAE+SCAD的方法应用于车辆mpg数据集。该数据集共有398个样本，8个变量，来自卡内基梅隆大学维护的StatLib库，用于研究城市汽车的循环油耗。数据可以从网站archive.ics.uci.edu/ml/data sets获得。目前，已经有一些学者对该数据集进行了研究。本文我们使用变系数模型式(1)研究该数据集，分析mpg与其它协变量的关系，并通过预测mpg验证模型的实用性。
mpg是汽车的油耗，即每加仑燃油可以行驶多少英里。本文我们将mpg作为因变量Y，气缸、位移、马力、重量、加速度和产量分别作为自变量X1,X2,X3,X4,X5,X6，车型年作为指标变量U。
(1)描述性统计分析
我们先做描述性统计分析，对数据进行初步探究。使用R软件绘制全部数据，即U=1,2,3时的散点图矩阵，如图2所示，分析各变量之间的相关关系。可以看出mpg与位移、马力、重量这些变量都有较强的负相关性，与加速度和车型年呈现一定的正相关性。同时，位移、马力、重量三个自变量之间相关性较强，说明可能存在多重共线性，需要在变量选择时剔除。选择车型年作为指标变量U，将所有数据分成三组之后，分别对U=1,U=2, U=3
三种情况做散点图矩阵。
(a)U=1,2,3(b)U=1
15
电子质量
2021年第04期(总第409期)
(2)回归分析
正则化的回归分析可以避免表1所示为分组前后所有变量之间的相关系数矩阵。分组前，mpg与气缸、位移、马力、重量有较强的负相关性，与加速度、车型年、产地有一定的正相关性。可以明显看出分组之后mpg与各自变量之间的相关系数有所增加，说明我们选取车型年作为U 是正确的。
为简化计算和提高变系数模型估计的准确率，先来判断βj (U )(j =1,2,…,6)是否是关于指标变量U 的函数。计算Y 与X j 的相关系数ρ(Y,X j )(j =1,2,…,6)，若|ρ(Y,X j )|大于给定阈值，说明Y 与X j 为线性关系，此时βj (U )为常数，此时可以不对其做泰勒展开；若否，说明βj (U )确实是关于U 的函数。
假定|ρ(Y,X j )|(j =1,2,…,6)的阈值为0.8，由表1可知，X 2,X 4与Y 的相关系数大于0.8，所以将β2(U ),β4(U )看做常数，不做泰勒展开，其余系数均进行泰勒展开。然后用LAE+SCAD 的估计方法对训练集进行估计，
得到表2所示的估计结果。若将系数的阈值设为0.0005，则得到回归结果：
Y =23.12263+0.00118X 2+(-0.00097U )X 3-0.00127X 4+(-0.00060+0.00156U )X 5+(-0.00059U +0.00050U 2)X 6
(9)
用测试集计算评价指标，结果为MSE =104.6279,RMSE =10.2288,MAE =7.9949，同时做出拟合效果图，如图3
所示。拟合效果总体来看还是比较好的。
(c)U =2
(d)U =3
图2散点图矩阵
图3连续情况拟合效果图
16

版权声明：本站内容均来自互联网，仅供演示用，请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198，我们将在24小时内删除。