在统计和数据分析领域研究中, 一个重要主题就是在响应变量和预测变量之间找到具有更好解释性和预测性的模型. 但是, 随着大数据的到来, 我们所处理的数据和研究的问题越来越复杂, 这种复杂性往往体现在高维性. 例如: 在影像遗传学研究中, 预测变量为基因SNP数据, 响应变量为扫描影像数据, 两者都高达几十万维, 若将全部变量代入模型分析, 模型的解释性和预测性会大大降低, 且没有什么实际意义. 所以, 简化模型对指导研究尤其重要, 较为常用的想法是进行特征选择.

特征选择通过从众多变量中选出重要变量从而获得稀疏的模型, 这些方法多是采用在损失函数后面添加惩罚项的方式来达到变量选择的作用. 针对于单变量选择的惩罚函数有LASSO, bridge^[1], 平滑削边绝对偏离惩罚(SCAD)^[2], 极小极大凹惩罚(MCP)^[3]等. 在处理实际问题时, 经常会遇到变量分组的情况, 比如本文实证部分数据来源于CSMAR数据库, 其就将企业财务数据根据盈利能力、成长能力、偿债能力、分红能力等将变量分成几组, 每组均包含若干个变量. 在这种情况下使用单变量选择方法, 就会忽略分组信息, 从而可能导致变量选择效果大大降低, 因此有些学者开始研究群组变量选择的方法. Yuan等^[4]提出了group LASSO, 其惩罚函数由群组水平上的 $  {L}_{2} $ 范数组成, 从而达到组级别上的稀疏性; Meier等^[5]将这一想法扩展到逻辑回归; 随后受group LASSO的启发, Wang等^[6]提出了Group SCAD, Huang等^[7]提出了Group MCP.

为了实现在选择组变量的同时又可以对组内的变量进行选择, 学者们在已有的方法上进行改进, 提出了双层变量选择方法, Huang 等^[3]对组结构使用 $  {L}_{1}$ 凹惩罚提出了 Group Bridge方法^[8]; Breheny等^[9]2009年提出了CMCP (Composite MCP), 其惩罚项是组间惩罚和组内惩罚的复合函数^[8]; Simon 等^[10]提出了SGL (Sparse Group LASSO), 这里的惩罚项是组内和组间惩罚函数组成的可加惩罚函数; 在此基础上, 王小燕等^[11]提出了adSGL (adaptive SGL), 其通过附权的方式对不同的组系数和单变量系数进行不同程度惩罚, 从而避免了对大系数的过度惩罚. 双层变量选择方法既考虑了分组信息, 也做到了对组间变量和组内变量的同时选择, 灵活性高, 同时变量选择的效果也更好.

针对群组变量下的二分类问题, 已有研究均集中在logistic模型, 即在logistic模型的负对数似然损失函数后面添加惩罚项, 在变量选择的同时直接对二分类问题作出预测, 其预测精度有进一步提升空间. AdaBoost是一种高效的集成学习方法, 其依次生成一系列基学习器, 并结合它们的预测来生成最终结果, 能够显著提升预测效果, 然而, 由于AdaBoost在处理高维数据时使用大量的基学习器来产生最终结果, 因此高内存空间消耗成为一个关键的挑战^[12]. 已有学者将特征选择方法应用于集成剪枝, 通过集成剪枝, 能产生一个规模更小但精度更高的模型. Margineantu等^[13]提出的Kappa剪枝是一种基于顺序的剪枝方法, 通过某种标准对学习器进行排序, 并选择最优的学习器组合组成最终的学习器; Lazarevic等^[14]提出一种基于族的集成剪枝方法, 即先对各学习器进行K-means聚类, 保留每类中表现较优的基学习器组成最终的分类器, 从而提高模型的泛化能力; Azmi等^[15]提出了一种利用LASSO的类关联规则进行集成剪枝的新方法, 该方法对不太重要的基学习器进行剪枝; Jiang等^[16]针对高维单变量选择提出了两阶段集成剪枝方法, 显著地提升了预测性能.

本课题提出并研究了一种基于MCP惩罚的AdaBoost集成剪枝逻辑回归模型(AdaMCPLR), 以产生对高维群组变量下的分类问题的卓越预测精度. 该方法可分为CMCP双层变量选择、集成和集成剪枝3个步骤. 第1步, 使用CMCP惩罚logistic回归来降低数据的维数; 第2步, 经过特征选择后的数据集生成AdaBoost集成; 最后, 再次使用MCP函数对集成系统进行剪枝, 并得到一个简单但更有效的模型. 在此步骤中, 我们还使用LASSO和SCAD对集成进行了修剪, 并通过模拟对比了3种惩罚的分类和选择性能, 显示了基于MCP的剪枝方法的优越性.

1 AdaMCPLR方法

AdaMCPLR可分为CMCP双层变量选择、集成和集成剪枝3个步骤. 首先, 在logistic模型的负对数似然损失函数后添加CMCP函数来降低数据的维数; 接着, 利用特征选择后的数据集生成AdaBoost集成; 最后, 我们再次使用MCP函数对集成系统进行了剪枝, 得到了一个简单但更有效的模型. 以下分别介绍各步骤中涉及的方法.

1.1 CMCP双层变量选择

当处理二分类任务时, logistic回归往往是最常用的模型. 考虑一个标准化的 $n \times d$ 设计矩阵 ${{X}} = ( {{{x}}_1}, $ $ {{{x}}_2}, \cdots ,{{{x}}_n} )^{\text{T}}$ ( $n$ 为观测数, $d$ 为特征数), 已知 ${{X}}$ 的 $ d$ 个特征共分为 $ J$ 组, 且第j组变量包含 ${K}_{j}$ 个变量( $j = 1,2, \cdots , $ $ J$ ), $p$ 维向量 ${{y}} = {\left( {{y_1},{y_2}, \cdots ,{y_n}} \right)^{\text{T}}}$ 为响应变量. 那么, logistic回归模型的形式如下:

$ {y}_{i}=\frac{1}{1+{\text{exp}} \left(-{{{x}}}_{{i}}{{\beta}} \right)}, i=1,\cdots ,n$

(1)

其中, $\;{{\beta}} = \left( {{\beta _1}, \cdots ,{\beta _d}} \right)$ 是模型的系数, logistic模型等价于最小化:

$L\left( {{\beta}} \right) = \mathop \sum \limits_{i = 1}^n \left\{ - {y_i}{{{x}}_i}{{\beta}} + {\rm{ln}}\left[ {1 + {\rm{exp}}\left( {{{{x}}_i}{{\beta}} } \right)} \right] \right\} $

(2)

其中, $L\left( {{\beta}} \right)$ 表示logistic模型的负对数似然损失函数. 通过在损失函数后面添加惩罚项来得到稀疏的模型, 是一种常用的正则化方法. 本文采用CMCP (Composite MCP)作为模型的惩罚项, CMCP是一种复合惩罚项, 其是由组内惩罚和组间惩罚组成的复合函数, 复合惩罚的具体形式为:

$\sum\limits_{j = 1}^J {{f_{\lambda ,b}}} \left( {\sum\limits_{k = 1}^{{K_j}} {{f_{\lambda ,a}}} (|{\beta _{jk}}|) }\right)$

(3)

其中, a、b 分别是内层惩罚函数和外层惩罚函数的调整参数, CMCP中的内层惩罚和外层惩罚都是极小极大凹惩罚(MCP)函数, 其具体形式及其导数形式分别为:

$ {f}_{\lambda ,\mu }\left({\beta }_{j}\right)=\left\{\begin{array}{l}\lambda {\beta }_{j}-\dfrac{{\beta }_{j}^{2}}{2\mu }, 0\le {\beta }_{j}\le \mu \lambda ,\mu \ge 0\\ \dfrac{\mu {\lambda }^{2}}{2},{\beta }_{j}\ge \mu \lambda ,\mu \ge 0\end{array}\right.$

(4)

$f_{\lambda ,u}'\left( {{\beta _j}} \right) = \left\{ {\begin{array}{*{20}{l}} {\lambda - \dfrac{{{\beta _j}}}{\mu },0 \le {\beta _j} \le \mu \lambda ,\mu \ge 0} \\ {0,{\beta _j} \ge \mu \lambda ,\mu \ge 0} \end{array}} \right.$

(5)

由式(5)可以看出, MCP一开始采用了与LASSO相同的处罚率, 然而, 它却不断地放松惩罚, 直到 $  \;{\beta }_{j}\ge $ $ \;\mu \lambda  $ , 惩罚率降为零. 添加惩罚项之后, logistic回归的优化问题变为:

$\begin{array}{*{20}{c}} {{\rm{min}}} \\ {{\beta}} \end{array}\mathop \sum \limits_{i = 1}^n \left\{ { - {y_i}{{{x}}_i}{{\beta}} + {\rm{ln}}\left[ {1 + {\rm{exp}}\left( {{{{x}}_i}{{\beta}} } \right)} \right]} \right\} + {f_{\lambda ,b}}\left( {\mathop \sum \limits_{k = 1}^{{K_j}} {f_{\lambda ,a}}(|{\beta _{jk}}|)} \right)$

(6)

通过最小化目标函数, 求解得到模型系数 $\hat {\;{\beta}}$ , 由于MCP函数具有特征选择的作用(某些变量系数会被压缩至0), 所以当组内惩罚和组间惩罚均为MCP时, 能同时对组间变量和组内变量进行选择, 达到双层变量选择的作用, 经过变量选择之后的设计矩阵记为 ${{{X}}_{\rm{sub}}}$ .

1.2 AdaBoost算法

Boosting是一种高效的集成学习方法, 它将众多基学习器结合在一起, 组成最终的强学习器, AdaBoost就是一个典型的Boosting算法. 本文利用经过特征选择后的数据集 ${{{X}}_{\rm{sub}}}$ 来构造集成模型, 进一步提升预测精度. AdaBoost模型主要由两个主要部分组成: 加法模型和前向分步算法. 加法模型可以表示如下:

$H\left( {{X}} \right) = \mathop \sum \limits_{t = 1}^T {\alpha _t}{h_t}\left( x \right)$

(7)

其中, ${h}_{t}\left(x\right)$ 代表第t步迭代中的弱分类器, 此处采用的基学习器是决策树; ${\alpha _t}$ 是第t个弱分类器在最终的强分类器上所占权重, $H\left( X \right)$ 代表弱分类器的线性组合, 也即最终的强分类器. 在前向分步算法中, 每一步迭代生成的分类器都是在前一次迭代的分类器的基础上产生的, 可以表示为:

${H_m}\left( x \right) = {H_{m - 1}}\left( x \right) + {\alpha _t}{h_t}\left( x \right)$

(8)

其中, ${H_{m - 1}}\left( x \right)$ 是前m−1步迭代中所有弱分类器的线性组合. 算法开始每个样本的权重相同, 在这样的样本分布下训练得到第一个弱分类器以及该分类器在训练集上的误差, 误差决定此分类器的权重. AdaBoost算法使用的损失函数是指数损失函数, 每一轮的弱分类器的权重计算公式如下:

$ {\alpha }_{t}=\frac{1}{2}\rm{ln}\left(\frac{1-{\epsilon}_{t}}{{\epsilon}_{t}}\right) $

(9)

由式(9)可知弱分类器错误率越低, 它在最终的分类器所占的权重就越大, 同时每一步迭代的训练样本的分布也会根据 ${\alpha _t}$ 调整, t轮迭代中第i个样本的权重计算如下:

$\begin{split} {W_{t + 1,i}}& = \frac{{{W_{t,i}}}}{{{Z_t}}}\exp \left( { - {\alpha _t}{y_i}{h_{t,i}}\left( x \right)} \right) \\ & = \left\{ {\begin{array}{*{20}{c}} {\dfrac{{{W_{t,i}}}}{{{Z_t}}}{{\text{e}}^{ - {\alpha _t}}},\;{h_{t,i}}{y_i}\left( x \right) = {y_i}} \\ {\dfrac{{{W_{t,i}}}}{{{Z_t}}}{{\text{e}}^{{\alpha _t}}},\;{h_{t,i}}{y_i}\left( x \right) \ne {y_i}} \end{array}} \right. \end{split} $

(10)

其中, $ {y}_{i }$ 表示第 $ i $ 个样本的分类标签, ${h_{t,i}}\left( x \right)$ 表示第 $ i $ 个样本在第 $ t $ 个弱分类器的预测结果, 由式(10)可知当该样本被误分类时, 样本权重将会增大, 而被正确分类时样本权重将会减少. 假设训练集共有N个样本, $\mathop Z\nolimits_{{t}} = \mathop W\nolimits_{{{t}},i}\Bigg /\displaystyle\sum\nolimits_{i = 1}^N {\mathop W\nolimits_{{{t}},i} }$ 是规范化的因素, 用于确保每轮迭代所有的样本权重总和为1. 最终的强分类器表示为 $f\left({x}\right)={\text{sign}} \left(H\left(X\right)\right)$ , 即最终的分类结果取决于最终组成的强分类器分类结果的符号.

1.3 集成剪枝

AdaBoost通过迭代生成一组基学习器来组成最终的分类器, 可以显著提高分类的准确率. 但是, 它也占用了相当多的内存, 特别是在处理具有许多特征的高维数据时; 此外, 由于错误分类的样本权重较大, 基学习器倾向于适应这些情况, 当数据中存在噪声时, 会导致过拟合. 因此, 必须对集成进行剪枝, 以减少存储空间的需求, 同时提高分类器性能. 受Jiang等^[16]的启发, 本文将集成剪枝扩展应用于群组变量下的分类预测问题.

AdaMCPLR在集成剪枝步骤中的构造方法类似于stacking方法, 在stacking方法中, 称单个基学习器为初级学习器, 将各基学习器结合起来的学习器叫做次级学习器或元学习器. 本文中的AdaMCPLR与stacking方法区别在于stacking的初级学习器是多种模型组成, 而AdaMCPLR的初级学习器均是决策树, 其次级学习器是logistic回归模型, 也即AdaMCPLR以AdaBoost迭代生成的一系列基学习器的预测结果作为logistic回归的输入, 并在logistic回归的负对数似然函数后面加上MCP惩罚, 从而达到集成剪枝的目的.

如上所述, 为了简化模型和提高模型性能, 必须进行集成剪枝. 本文提出的AdaMCPLR将惩罚logistic回归的概念用于Boosting集成剪枝, 考虑以下n×T矩阵:

$Z = \left( {\begin{array}{*{20}{c}} {{h_1}\left( {{{{x}}_1}} \right)}& \ldots &{{h_T}\left( {{{{x}}_1}} \right)} \\ \vdots & \ddots & \vdots \\ {{h_1}\left( {{{{x}}_n}} \right)}& \ldots &{{h_T}\left( {{{{x}}_n}} \right)} \end{array}} \right)$

(11)

其中, $ i\left(i=1,2,\cdots ,n\right)$ 行表示 $T$ 个基学习器在给定样本 ${x}_{{i}}$ 的情况下所做的预测,也即每个基学习器都可以被看作logistic回归模型中的一个特征, 以 $ Z $ 作为logistic回归模型的设计矩阵, 从而将集成剪枝问题转化为一个优化问题, 通过最小化式(12)来寻求稀疏权向量 $ \hat{\theta }: $

$ \sum\limits_{i=1}^{n}\left\{-{y}_{i}{{\textit{z}}}_{i}\theta +\ln\left[1+\exp\left({{\textit{z}}}_{i}\theta \right)\right]\right\}+{f}_{{\lambda }_{1},{\gamma }_{1}}\left({\theta }_{j}\right)$

(12)

其中, ${{{{\textit{z}}}}_{{i}}}$ 表示 ${Z}$ 的第 ${i}$ 行, 这里的惩罚项 ${f_{{{{\lambda }}_1},{{{\gamma }}_1}}}\left( {{{{\alpha }}_{{j}}}} \right)$ 是MCP函数. 将MCP惩罚应用于集成剪枝中, 一些基学习器的相关系数会被压缩至0, 从而达到集成剪枝的目的.

2 算法

Zhao等^[17]提出了PICASSO算法 (PathwIse CalibrAted Sparse Shooting algOrithm), 其分为3个循环: 外循环、中间循环和内循环, 由于在每次迭代中, 中间循环从所有变量中只选定一小部分变量的索引组成活跃集, 并将活跃集传递到内循环中, 利用坐标下降算法求解出对应的解, 从而大大提高了求解效率. 本文对PICASSO算法进行改进, 使其能够应用于群组变量选择.

考虑一个标准化的 $n \times d$ 设计矩阵 $X = {\left( {{x_1},{x_2}, \cdots ,{x_n}} \right)^{\text{T}}}$ ( $n$ 为观测数, $d$ 为特征数) 和二元响应向量 $y \in {\left\{ { - 1,1} \right\}^n}$ , 已知 $X$ 的 $ d$ 个特征共分为 $ J $ 组, 且第 $ j$ 组变量包含 ${K}_{{j}}$ 个变量( $j = 1,2, \cdots ,J$ ). 在生成Boosting 集合之前, 我们筛选变量并通过最小化 ${F_\lambda }\left( \beta \right)$ 来降低维数:

$\begin{split} {F_\lambda }\left( \beta \right) =& \mathop \sum \limits_{i = 1}^n \left\{ { - {y_i}{{{x}}_i}{{\beta}} + {\rm{ln}} \left[ {1 + \exp \left( {{{x}_i}{{\beta}} } \right)} \right]} \right\} \\ & + \mathop \sum \limits_{j = 1}^J {f_{\lambda ,b}}\left( {\mathop \sum \limits_{k = 1}^{{K_j}} {f_{\lambda ,a}}(|{\beta _{jk}}|) } \right) \end{split} $

(13)

其中, 式(13)的第一项表示logistic回归的负对数似然损失函数 $L\left( {{\beta}} \right)$ , 后一项中的 ${f_{\lambda ,a}}$ 、 ${f_{\lambda ,b}}$ 分别表示组内和组间的惩罚函数, 此处都是MCP函数, a和b分别表示组内和组间惩罚的权重, 一般取 $b = {\textit{K}_{j}}a\lambda /2$ .

在外循环中, 生成正则化参数 $ \lambda  $ 的一个递减序列. $\lambda $ 的初始值记为 ${\lambda ^{\left( 0 \right)}}$ , 初始状态下 $\;{\beta }^{\left(0\right)}={\left(0,0,\cdots ,0\right)}^{{\text{T}}}$ , 随着 $ \lambda  $ 的逐渐减小, 解向量出现更多的非零元素. $\lambda $ 递减的速率由参数 $\eta \in \left( {0,1} \right)$ 控制.

在中间循环中, 会从 ${{X}}$ 中选出活跃变量, 并将活跃变量的索引组成活跃集 $A$ , 初始的活跃集 ${A^{\left( 0 \right)}}$ 定义为:

$\begin{split} {A^{\left( 0 \right)}} = \{ j|\beta _j^{\left( 0 \right)} = 0 \; {\text{and}} \; {\nabla _j}L\left( {{\beta ^{\left( 0 \right)}}} \right) \ge \left( {1 - \delta } \right)\lambda \} \cup \{ j{\rm{|}}\beta _j^{\left( 0 \right)} \ne 0{\rm{\} }} \end{split} $

(14)

其中, $\delta $ 是一个很小的值, 此处取0.1, ${\nabla _j}L\left( {{{{\beta}} ^{\left( 0 \right)}}} \right)$ 是梯度向量 $\nabla L\left( {{{{\beta}} ^{\left( 0 \right)}}} \right)$ 的第j个元素. 相应的非活跃集定义为:

${\bar A^{\left( 0 \right)}} = \{ j|\beta _j^{\left( 0 \right)} = 0 \; {\text{and}} \; {\nabla _j}L\left( {{{{\beta}} ^{\left( 0 \right)}}} \right) < \left( {1 - \delta } \right)\lambda \} $

(15)

只有 ${{A}}$ 内的元素会被传递到内循环中, 从而求解出对应的解. 同时中间循环的优化目标变为:

$ \hat{{{\beta}} }=\arg \min{F}_{\lambda }\left({{\beta}} \right) , \; {\text{s.t.}} \; {{{\beta}} }_{\bar {A}}=0$

(16)

在迭代过程中, 一旦活跃集中某元素对应的系数变为0, 如 $\;{\beta _j}$ , 那么 $ { j}$ 将从活跃集中移除. 同时如果非活跃集中的元素 $k$ 满足式(17), 将会被加入到活跃集中.

$k = \arg \max{_{k \in \bar A}}\left| {{\nabla _j}L\left( {{\beta}} \right)} \right|$

(17)

直到 $\left| {{\nabla _j}L\left( {{\beta}} \right)} \right| < \left( {1 - \delta } \right)\lambda $ , 中间循环的迭代停止. 此处的 $\lambda $ 是由外循环给定的.

在内循环中, 利用坐标下降算法对活跃集中元素对应的变量系数进行求解, 在 $t + 1$ 步迭代中( $t = 0,1, $ $ 2, \cdots, T$ ):

$\beta _j^{t + 1} = \arg \min{_{{\beta _j}}}{F_\lambda }\left( {{\beta _j},{{\beta}} _{\backslash j}^{\left( t \right)}} \right)$

(18)

其中, $\beta _j^{\left( t \right)}$ 表示 $\;{{\beta}}$ 的第j个分量在第 $t$ 步迭代中估计值, $\;{{{\beta}} _{\backslash j}}$ 是 $\;{{\beta}}$ 删除了第j个分量后的子向量. 在内循环中, $\;{{{\beta}} ^{t + 1}}$ 每次只更新一个元素而保持其他元素不变, 这也是坐标优化的主要思想.

以上是针对群组变量选择而改进的PICASSO算法. 本文提出的AdaMCPLR分为3个步骤, 在特征选择和集成剪枝中均用到了PICASSO算法, 每个步骤的具体算法如算法1.

算法1. 特征选择

(1)输入设计矩阵 $\scriptstyle X$ ; 响应变量 $\scriptstyle y$ ; 参数 $\scriptstyle\delta$ ; $\scriptstyle \eta \in \left( {0,1} \right)$ ; 算法迭代次数 $\scriptstyle M$ ;

(2)选定初始值 $\scriptstyle {\lambda ^{\left( 0 \right)}} = \max\left\{ {\nabla L\left( 0 \right)} \right\}$ (0是长度为 $\scriptstyle d$ 的零向量), 初始解 $\scriptstyle \;{\beta ^{\left( 0 \right)}} = {\left( {0,0, \cdots ,0} \right)^T}$ ;

(3)定义初始活跃集为:

$\begin{array}{l} \scriptstyle {A}^{\left(0\right)}=\left\{j\right|{\beta }_{j}^{\left(0\right)}=0 \; {\text{and}} \; {\nabla }_{j}L\left({\beta }^{\left(0\right)}\right)\ge \left(1-\delta \right){\lambda }^{\left(0\right)}\}\\ \scriptstyle \cup \{j\rm{|}{\beta }_{j}^{\left(0\right)}\ne 0\},t=0; {\lambda }^{\left(m\right)}={\lambda }^{\left(0\right)}{\eta }^{m}\end{array}$

(4)计算 $\scriptstyle \;{\beta }^{\left(t+0.5\right)}=\arg \min_{\beta }{F}_{{\lambda }^{\left(m\right)}}\left(\rm{β}\right), \; {\text{s.t.}}{\beta }_{\overline{A}}=0$ ;

(5)更新活跃集 $\scriptstyle {A^{\left( {t{\rm{|}}0.5} \right)}} = \left\{ {j{\rm{|}}\beta _j^{\left( {t + 0.5} \right)} \ne 0} \right\}$ , $\scriptstyle {\bar A^{\left( {t + 0.5} \right)}} = $ $ \scriptstyle \left\{ {j{\rm{|}}\beta _j^{\left( {t + 0.5} \right)} = 0} \right\}$ ;

(6) $\scriptstyle {k_t} = \arg \max{_{k \in {{\bar A}^{\left( {t + 0.5} \right)}}}}\left| {{\nabla _k}L\left( {{\beta ^{\left( {t + 0.5} \right)}}} \right)} \right|$ , 将非活跃集中的 $\scriptstyle {k_t}$ 加入到活跃集中, 从而活跃集变为 $\scriptstyle {A}^{\left(t+1\right)}={A}^{\left(t+0.5\right)}\cup $ $\scriptstyle \left\{{k}_{t}\right\},\;{\overline{A}}^{\left(t+1\right)}={\overline{A}}^{\left(t+0.5\right)}\backslash \left\{{k}_{t}\right\}$ , 直到 $\scriptstyle \left| {{\nabla _k}L\left( {{\beta ^{\left( {t + 0.5} \right)}}} \right)} \right| \le \left( {1 + \delta } \right)\lambda$ ;

重复步骤(4)–步骤(6)直至M次迭代结束. 经过特征筛选后, 设计矩阵 ${{X}}$ 包含的特征数量大大减少. 然后将响应向量 ${y}$ 重新编码为 ${\left\{ { - 1,1} \right\}^n}$ , 并将简化后的设计矩阵表示为 ${{{X}}}_{\rm{sub}}$ . 本文第1.2节的算法如算法2所示.

算法2. 集成

(1)初始化样本权重 $\scriptstyle {\omega _i} = \scriptstyle \frac{1}{n}, \; i = 1,2, \cdots ,n,$ $\scriptstyle n$ 为样本容量, 根据权重为 $\scriptstyle {\omega _i}$ 的样本构造决策树 $\scriptstyle {G_t}\left( x \right).\left( {t = 1,2, \cdots ,T} \right)$ ;

(2)计算决策树的加权错误率 $\scriptstyle{ er{r_t} = \frac{{\mathop \sum \nolimits_{i = 1}^n {\omega _i}I\left( {{y_i} \ne {G_t}\left( {{x_i}} \right)} \right)}}{{\mathop \sum \nolimits_{i = 1}^n {\omega _i}}}}$ ;

(3)根据错误率计算决策树的权重 $\scriptstyle {\alpha _t} = {\rm{ln}} ( (1 - er{r_t} )/er{r_t} )$ ;

(4)更新样本权重 $\scriptstyle {\omega _i} \leftarrow {\omega _i} \cdot \exp [ {\alpha _t} \cdot I\left( {{y_i} \ne {G_t}\left( {{x_i}} \right)} \right) ],$ $\scriptstyle i = 1,2, \cdots, n$ ;

迭代 $ T$ 次( $  T$ 为决策树的数量), 输出最终预测 $G\left( x \right) = $ $ {\text{sign}}\left(\displaystyle\sum\nolimits_{t = 1}^T {{\alpha _t}{G_t}\left( x \right)}\right)$ , 将 $ T$ 个基学习器的预测结果组成 $n \times $ $ T$ 的设计矩阵 $ Z$ , 对响应向量 $ { y}$ 重新编码为 ${\left\{ {0,1} \right\}^n}$ , 将集成剪枝转化为特征选择问题, 对特定的冗余项 $ {\theta }_{t}(t=1, $ $ 2,\cdots ,T)$ 压缩到零来达到集成剪枝的目的. 我们再次使用PICASSO和MCP来解决, 集成剪枝算法的目标函数如式(20), 集成剪枝算法除了目标函数与算法1不一致以外, 其余和算法1均相同, 所以此处的集成剪枝算法省略不表.

$\mathop \sum \limits_{i = 1}^n \left\{ { - {y_i}{{{{\textit{z}}}}_i}{{\theta}} + \rm{ln} \left[ {1 + \exp\left( {{{{{\textit{z}}}}_i}{{\theta}} } \right)} \right]} \right\} + {f_{\lambda ,\mu }}\left( {{{{\theta}} _j}} \right)$

(19)

其中,

$ {f}_{\lambda ,\mu }\left({\beta }_{j}\right)=\left(\lambda {\beta }_{j}-\frac{{\beta }_{j}^{2}}{2\mu }\right)I(\left|{\beta }_{j}\right|\le \mu \lambda)+\frac{\mu {\lambda }^{2}}{2}({\rm{|}}{\beta }_{j}|>\mu \lambda)$

(20)

3 模拟分析

为了验证提出的AdaMCPLR方法的优劣, 本节设置了模拟实验, 将不同的双层变量选择方法进行了对比, 同时将不同的惩罚项对集成进行剪枝的性能也作了对比. 采用选择单个变量数(nv)、选择组变量数(ngrp)、错误发现率(FDR)、假阴性率(FNR) 以及交叉验证错误率(cv.error)作为变量选择效果评价指标, 其中nv和ngrp是越接近真实数值越好, 其他指标则是越低模型越优; 同时采用交叉验证错误率(cv.error)、集成规模(Average size)作为集成剪枝效果评价指标, 其中Average size 表示经过集成剪枝之后保留的基学习器的个数.

模拟实验的数据生成过程如下: ① 设置样本数为n, 变量数为p, 变量分为m组, 每组变量数设为5 (n=200, p=50, 则m=10); ② 产生随机变量 ${R_1}, \cdots , $ $ {R_p}, \; {\text{i.i.d}}N\left( {0,I} \right)$ , 以及标准正态分布向量 $Z = {\left( {{Z_1}, \cdots ,{Z_m}} \right)^{\text{T}}}$ , 其中 $E\left( Z \right) = 0$ , $  COV\left({Z}_{{j}_{1}},{Z}_{{j}_{2}}\right)={0.6}^{\left|{j}_{1}-{j}_{2}\right|}$ ; ③ 产生自变量 ${X_{5\left( {j - 1} \right) + k}} = {{{Z_j} + {R_{5\left( {j - 1} \right) + k}}} / {\sqrt 2 }},j = 1, \cdots ,m,k = 1, \cdots ,5\;$ ; ④ 自变量和因变量之间的系数设定为 ${\omega ^{\text{T}}} = \left( {\omega _1^{\text{T}},\omega _2^{\text{T}}, \cdots ,\omega _m^{\text{T}}} \right)$ , 其中 $\omega _1^{\text{T}} = \left( { - 1,1,0,1.5,2} \right)$ , ${\omega }_{2}^{{\text{T}}}=\left(0,0,-0.5,0,1\right)$ , $\omega _3^{\text{T}} = $ $ \left( {0,1,2,0,0} \right)$ , ${\omega }_{4}^{\text{T}}=\left(0,0.5,1,0,0\right)$ , ${\omega }_{5}^{\text{T}}=\left(-0.5,1,0,0,0\right)$ , $\omega _6^{\text{T}} = \ldots = \omega _m^{\text{T}} = 0$ ; ⑤ 因变量 ${y_i}\sim Bernoulli\left( {1,{p_i}} \right)$ , ${p_i} = $ $ prob\left( {{y_i} = 1{\rm{|}}{x_i}} \right) = {1 / {1 + {\text{exp}}\left( { - {\omega ^{\text{T}}}{x_i}} \right)}},i = 1,2, \cdots ,n$ 模拟重复100次, 各评价指标值为100次模拟结果的均值.

由表1可以看出, 群组变量选择方法中, group LASSO总是倾向于选择过多的变量, 导致很多非重要变量被误选; group SCAD恰好相反, 总是倾向于选择过少变量, 导致变量漏选; 在双层变量选择方法中, group bridge也是漏选变量导致FNR指标过大. 综合各指标来看, CMCP在变量选择环节中表现最优, 这也是我们在提出的AdaMCPLR方法中使用CMCP作为变量选择方法的原因.

表2中的第一行表示未集成剪枝之前模型的各项评价指标值, 其余各行分别表示使用LASSO、SCAD和MCP作为集成剪枝惩罚项时模型表现. 将表2中的 cv.error 值与表1中的做对比可以发现, 无论是否进行集成剪枝, 模型的预测精度都有了很大的提升, 同时使用MCP函数作为惩罚项, 能够很大限度地精简模型, 同时模型的预测精度也有所提升, 所以本文提出的方法中使用MCP函数作为集成剪枝的惩罚项.

4 实证检验

当一家公司出现财务困境时, 会对许多外部和内部经济主体产生严重影响, 如股东、公司债权人、客户和供应商、员工以及政府, 因此开发财务困境预测方法一直是金融研究的优先目标之一. 有效的财务困境预测机制不仅可以帮助公司管理者更好地根据自身财务信息提前做好规划, 避免企业遭受财务困境的影响; 也能帮助投资者根据上市公司披露的财务指标信息更好地作出投资选择, 合理规避重大风险; 对于政府监管机构, 有效的财务预警机制有利于政府科学调控经济和优化资源配置, 更好地维护市场稳定.

Beaver^[18]首次提出利用单一指标来构建财务困境预测模型, 但是单一指标并不能涵盖整个财务指标体系对模型的影响, 不能达到理想的预测精度; Altman^[19]提出了Z计分模型, 提出使用多变量判别进行财务困境预测研究; Martin^[20]用logistic回归模型构建财务困境预测模型, 并对银行财务困境预测进行了实证研究; Laitinen等^[21]也使用了logistic回归构建了上市公司财务困境预测模型. 国内方面, 郭斌等^[22]首先构建了完善的预测指标体系, 并在此基础上建立了logistic回归财务困境预测模型; 韩立岩等^[23]利用logistic回归模型对中小公司的财务困境预测进行了实证研究; 在相当长的一段时间内, logistic回归都是财务困境预测研究的主流方法.

随着机器学习的理论和技术的不断发展, 越来越多的学者开始使用机器学习算法建立财务困境预测模型. 殷尹^[24]利用BP神经网络构建了中短期预测模型, Min等^[25]利用神经网络与SVM和 Logistic回归等方法分别构建了财务困境预测模型, 发现SVM取得的预测效果最优; 丁德臣^[26]提出了基于混合全局优化正交遗传算法和支持向量机的财务困境预警模型; 王小燕等^[11]提出了adSGL-logit模型, 并将其应用于信用卡信用评分模型中; 方匡南等^[27]提出了稀疏组 LASSO支持向量机方法(SGL-SVM), 并将其应用到我国制造业上市公司财务困境预测中.

当上市公司出现财务异常, 并且该异常可能会导致其股票终止上市, 或者投资者的投资权益可能会受到重大损害时, 证券交易所会对该上市公司股票交易实行特别处理(Special Treatment, ST). 被特别处理的股票除了会被实行股票报价的日涨跌幅限制为5%以外, 证监会还要求其在股票名称前加上提示性注释“ST”,以区别于一般股票.

本文数据来源于CSMAR经济金融研究数据库, 以全部上市公司为基础样本, 从中选取截止到2020年6月底因财务困境被特别处理(ST)的公司作为正样本, 其余的上市公司作为负样本. 剔除了正样本中缺失值超过30%的指标后, 共选取280个财务指标进行建模. 根据CSMAR数据库对财务指标的分类, 本文选取了偿债能力、披露财务指标、比率结构、经营能力、盈利能力、现金流分析、风险水平、发展能力共 8组280个指标.

在进行数据处理时, 首先剔除指标缺失超过30%的上市公司, 其余的缺失值用指标均值代替; 然后对所有指标进行Z-score标准化处理, 从而消除指标量纲对模型的影响. 最终得到3451个样本, 其中正样本83个, 负样本3368个, 正负样本的数量差异很大, 是非常典型的非平衡数据, 需要采用重抽样的方法使其正负样本均衡. 因此, 本文采用双层重抽样的方法对样本进行处理, 即通过上采样来增加少数类样本的数量, 同时利用下采样来减少多数类样本, 从而来获取相对平衡的数据集. 本文利用SMOTE算法^[28]对少数类进行过抽样, SMOTE算法的基本原理大致描述为: 在每个少数类样本a的最近邻样本中随机选择一个样本b, 接着在样本a和b之间的连线上随机选取一点作为新的少数类样本. Python中的imblearn包提供了更方便的接口, 本文直接调用imblearn中对应的SMOTE方法和下采样中的RandomUnderSampler方法来对样本进行双层重抽样, 从而使正负例比例接近1:1, 最终重抽样之后的总样本个数为830个. 在总样本中随机抽取 70%作为训练集, 剩余的30%的样本作为模型的交叉验证. 接下来, 本文利用AdaMCPLR对处理后的所有样本进行建模, 经过特征选择之后, 保留了39个重要指标. 用这39个重要性指标数据进行集成及集成剪枝. 在此处设置200棵决策树来构建AdaBoost模型, 并用MCP函数来进行集成剪枝, 从而得到最终结果. 将最终AdaMCPLR方法的变量选择效果和预测效果与logistic回归、SVM、SGL-SVM以及adSGL方法的进行了对比, 对比结果如表3.

由表3可知, 传统的logistic回归和SVM方法均不具有变量选择的作用, 保留了全部的变量, 但是预测精度和AUC值均不高. AdaMCPLR的预测精度和AUC值较其他预测方法都是最高的, 另外从选择的变量个数来看, AdaMCPLR仅选择了39个变量, 就达到了很高的样本外预测精度, 同时经过集成剪枝, 模型只保留了47个基学习器. 由此可以看出, 并不是越多变量参与模型, 模型就越优, 反而过多的冗余变量进入模型, 会降低模型的预测精度.

5 总结

针对高维群组变量下的二分类问题, 本文提出的AdaMCPLR分为两阶段, 利用惩罚项既对群组变量进行了特征选择, 又对AdaBoost集成模型的基学习器进行了剪枝, 这是在已有的高维群组变量数据分类预测问题上所没有出现过的方法.

将CMCP用于双层变量筛选, 降低了数据的冗余度, 提高变量选择的精度, 用特征选择后的数据生成AdaBoost集成, 大大提高了模型的预测精度, 由于集成模型由大量基学习器组成, 使用MCP惩罚方法对这些基学习器进行剪枝, 选择的基学习器组合为最终的预测集合, 在提高预测精度的同时, 也大大降低了集成模型的空间消耗. 本文还对PICASSO算法进行了改进, 使其能够应用于群组变量选择问题, 大大提高了模型的求解速度.

本文提出的AdaMCPLR方法适用于生物信息、管理科学、经济学、金融学等领域高维群组变量数据分类预测问题, 在降低空间消耗的同时, 较已有的方法提高了预测精度.