基于Logistic回归模型的高维数据变量选择方法比较研究

A Comparative Study of Variable Selection Methods for High Dimensional Data Based on Logistic Regression Model

摘要: 高维数据已成为现代大数据分析中的热点研究领域。变量选择是一种被广泛用于高维数据分析问题的方法。文献中已出现大量高维变量选择方法，为研究其中有影响的几种方法的适用范围和利弊，本文考虑了lasso、自适应lasso等变量选择方法来研究logistic回归模型中的变量选择问题。首先，通过随机模拟实验研究，分别在低维和高维的情况下比较不同变量选择方法的预测和变量选择效果。然后，在实际数据集中做进一步地实证比较研究。研究结果表明：在同等条件下，自适应lasso在模型预测和可解释性方面均比lasso更具优势。

Abstract: High-dimensional data has become a hot research field in modern large data analysis. Variable selection is a widely-used method for high-dimensional data analysis. A large number of high-dimensional variable selection methods have appeared in the literatures. In order to compare the scope of application, advantages and disadvantages of several influential methods, in this paper, we consider the variable selection methods such as lasso and adaptive lasso to study the variable selection problem in logistic regression model. Firstly, by random simulation experiments, we compare the prediction and selection effects of different variable selection methods in low and high dimensions respectively. Then, we do further empirical analysis in the real data. The results show that under the same conditions, adaptive lasso has more advantages than lasso in model prediction and interpretability.

1. 引言

随着科学技术的发展，很多研究领域中的大规模的高维数据越来越多。高维数据分析逐渐成为现代大数据分析中的一个重要的研究课题。对大规模高维数据分析的现实需求推动了统计思维、数据分析方法及理论研究的发展。

在机器学习领域，降维或变量选择技术被广泛运用于高维数据的分析问题中。降维或变量选择技术的一个关键思想是“稀疏原理”，即假设真实模型仅由预测变量的一个稀疏子集来决定。在统计学习过程中，当真实模型存在稀疏表示时，变量选择不仅能确保拟合模型的预测精度，也能提高模型的可解释性。现有文献中提出了很多变量选择的标准。传统的变量选择方法可追溯到 $C_p$ 、AIC、BIC等，但它们都存在计算开销会随着数据维度的增加而增加的问题。目前，常用变量选择方法主要有子集选择、系数压缩(或正则化)等方法。子集选择同样存在计算开销很大的问题，而且子集选择是一个离散的过程，这导致模型拟合变异性较高。系数压缩或正则化方法是一个连续的变量选择过程，即可同时进行系数估计和变量选择。它对系数施加一定程度的惩罚，使系数被压缩甚至压缩至0。虽然正则化方法得到的估计是有偏差的，但却能在很大程度上改进模型的预测效果，即牺牲较小的偏差而在很大程度上降低方差，在提高预测精度的同时到达变量选择的目的，实现模型预测精度和可解释性的平衡 [1] 。在正则化方法中，具有里程碑意义lasso (Tibshirani, 1996) [2] 是一种基于“ ${\mathcal{l}}_1$ ”惩罚的正则化方法。在此基础上，又发展出了adaptive lasso (Hui ZOU, 2006) [3] 等不同形式的变量选择方法。

本文首先介绍logistic回归模型、lasso、自适应lasso的理论研究，然后通过随机模拟实验，在低维数据和高维数据的情况下，从系数估计效果、模型的可解释性、分类正确率等方面对不同的变量选择方法进行比较，并进一步在实际数据集进行比较。

2. 理论模型介绍

2.1. Logistic回归模型

本文主要考虑分类问题中的一个重要的模型：logistic回归模型。给定数据 $\left(x_i,y_i\right),i=1,2,\cdots ,n$ ，其中， $x_i\in {ℝ}^p$ 是预测变量， $y_i\in \left\{0,1\right\}$ 是服从二项分布的响应变量，logistic回归模型的定义为

$\text{P}\left(y_i=1|x_i\right)=p_i\left(\beta \right)=\frac{\text{exp}\left(x_i^{\text{T}}\beta \right)}{1+\text{exp}\left(x_i^{\text{T}}\beta \right)}$ (1)

其中， $\beta \in {ℝ}^p$ 是 $p\times 1$ 的系数向量。通常用极大似然估计法求解未知系数向量 $\beta$ ，即求解(2)式的优化问题。

${\hat{\beta }}_{MLE}=\arg {\max }_{\beta }\mathcal{l}\left(\beta \right)=\arg {\max }_{\beta }\underset{i=1}{\overset{n}{{\displaystyle \sum }}}\left[y_i\log p_i\left(\beta \right)+\left(1-y_i\right)\log \left\{1-p_i\left(\beta \right)\right\}\right]$ (2)

从分析上来说，由于(2)不存在闭合解，因此用牛顿迭代法来求解(2)式。牛顿迭代法的迭代过程公式如下：

${\hat{\beta }}^{\left(t+1\right)}={\hat{\beta }}^{\left(t\right)}-{\left\{\underset{i=1}{\overset{n}{{\displaystyle \sum }}}{w}_i\left({\hat{\beta }}^{\left(t\right)}\right)x_i{x}_i^{\text{T}}\right\}}^{-1}\frac{\partial \mathcal{l}\left({\hat{\beta }}^{\left(t\right)}\right)}{\partial \beta }$ (3)

其中， $w_i\left(\beta \right)=p_i\left(\beta \right)\left\{1-p_i\left(\beta \right)\right\}$ ， $\mathcal{l}\left(\beta \right)$ 为(2)式中的损失函数。若(3)收敛，迭代过程终止。

2.2. 变量选择方法

Lasso是一种可以同时进行变量选择和系数估计的正则化方法(Tibshirani, 1996) [3] 。Logistic回归模型的lasso估计可定义为

${\hat{\beta }}_{lasso}=\arg {\min }_{\beta }\underset{i=1}{\overset{n}{{\displaystyle \sum }}}\left[-y_i\left(x_i^{\text{T}}\beta \right)+\log \left(1+{\text{e}}^{x_i^{\text{T}}\beta }\right)\right]+\lambda \underset{j=1}{\overset{p}{{\displaystyle \sum }}}\left|{\beta }_j\right|$ (4)

其中， $\lambda$ 是非负的正则参数。(4)式中的第二项也称为“ ${\mathcal{l}}_1$ ”惩罚项。当 $\lambda$ 逐渐增大时，lasso方法可连续地将系数向0压缩，当 $\lambda$ 足够大时，一些系数可以被压缩为0，从而达到变量选择的目的。从(4)显然可见，lasso方法对所有系数所施加的惩罚都是相同的，因此在lasso方法下，系数的值越大，估计结果的偏差越大，而且lasso估计不具有神谕性(oracle properties)，神谕性指正确选择真实模型中的变量的概率收敛到1，且非零系数的估计是渐近正态的。ZOU (2006)提出了自适应lasso (adaptive lasso)，其关键思想是在“ ${\mathcal{l}}_1$ ”惩罚项中，对不同的系数 ${\beta }_j$ 使用不同的自适应权重。Adaptive lasso (以下简称alasso)估计定义为

${\hat{\beta }}_{alasso}=\arg {\min }_{\beta }\underset{i=1}{\overset{n}{{\displaystyle \sum }}}\left[-y_i\left(x_i^{\text{T}}\beta \right)+\text{log}\left(1+{\text{e}}^{x_i^{\text{T}}\beta }\right)\right]+\lambda \underset{j=1}{\overset{p}{{\displaystyle \sum }}}{\omega }_j\left|{\beta }_j\right|$ (5)

其中， $\omega ={\left({\omega }_1,{\omega }_2,\cdots ,{\omega }_P\right)}^{\text{T}}$ 为已知的基于数据集得到的权重向量。不同于lasso，自适应lasso估计具有神谕性。Logistic回归模型中，通常用极大似然估计来设计自适应权重，即 ${\omega }_j={\left(\left|{\left({\hat{\beta }}_{MLE}\right)}_j\right|\right)}^{-\gamma }$ ， $\gamma$ 是一个大于0的常数，该方法记为APLR。然而，在高维数据中，由(2)式定义的极大似然估计可能无法求解，此时基于极大似然估计的自适应权重的变量选择方法就不再适用于高维数据的情况。为解决alasso在高维数据中面临的问题，有学者用lasso估计来设计自适应权重(Bielza, C., Robles, V., & Larrañaga, P. (2011)) [4] ，即 ${\omega }_j={\left(\left|{\left({\hat{\beta }}_{lasso}\right)}_j\right|\right)}^{-\gamma }$ ，该方法记为LAPLR。又有学者提出了基于预测变量的相关系数来设计自适应权重的方法 [5] ，即 ${\omega }_j={\left(\left|{\left({\hat{\beta }}_{cb}\right)}_j\right|\right)}^{-\gamma }$ ，其中 ${\hat{\beta }}_{cb}$ 可通过求解(6)式的优化问题得到，该方法记为CAPLR。

${\hat{\beta }}_{cb}=\arg {\min }_{\beta }\underset{i=1}{\overset{n}{{\displaystyle \sum }}}\left[-y_i\left(x_i^{\text{T}}\beta \right)+\text{log}\left(1+{\text{e}}^{x_i^{\text{T}}\beta }\right)\right]+\lambda \underset{i=1}{\overset{p-1}{{\displaystyle \sum }}}\underset{j>i}{{\displaystyle \sum }}\left\{\frac{{\left({\beta }_i-{\beta }_j\right)}^2}{1-{\rho }_{ij}}+\frac{{\left({\beta }_i+{\beta }_j\right)}^2}{1+{\rho }_{ij}}\right\}$ (6)

${\rho }_{ij}$ 为第i，j个预测变量之间的相关系数，注意 ${\rho }_{ij}\ne 1$ 。

3. 随机模拟研究

本文通过随机模拟实验来比较2.2中提到的lasso及基于不同权重设计方案的自适应lasso在logistic回归模型中的变量选择和预测效果。我们主要通过系数的均方误差 $MSE\left(\hat{\beta }\right)$ ，变量选择的特异性 $Specificity\left(\hat{\beta },\beta \right)$ 和敏感性 $Sensitivity\left(\hat{\beta },\beta \right)$ 、分类正确率率这几个指标来进行比较 [6] 。它们的定义如(7)，(8)，(9)所示。

$MSE\left(\hat{\beta }\right)=\frac{1}{B}\underset{b=1}{\overset{B}{{\displaystyle \sum }}}{\Vert {\hat{\beta }}_b-\beta \Vert }^2$ (7)

$Sensitivity\left(\hat{\beta },\beta \right)=\frac{\#\left\{\left(b,j\right):{\hat{\beta }}_{bj}\ne 0,{\beta }_{bj}\ne 0\right\}}{\#\left\{\left(b,j\right):{\beta }_{bj}\ne 0\right\}}$ (8)

$Specificity\left(\hat{\beta },\beta \right)=\frac{\#\left\{\left(b,j\right):{\hat{\beta }}_{bj}=0,{\beta }_{bj}=0\right\}}{\#\left\{\left(b,j\right):{\beta }_{bj}=0\right\}}$ (9)

其中，B表示随机模拟的次数，“#”表示计数。从(8)式和(9)式可以看出，敏感性指标反映的是真实模型中系数非零的变量的选择情况，敏感性越接近于1，说明变量选择方法越能识别模型的非零系数，若等于1，则说明真实模型的非零参数能完全被识别出来。特异性指标的值越接近于1，说明变量选择的方法能正确识别真实模型中的零系数变量，模型的可解释性越好，若其值等于1，则说明真实模型的零参数能完全被识别出来。两个指标的值越接近于1，变量选择的效果越好。

3.1. 低维数据情形

首先，我们先研究低维情况下的变量选择问题。本文模拟实验的参数设置如下： $n=500,p=20$ ，预测变量的稀疏水平 $s=5$ ，系数真值为 ${\beta }_0={\left(1,1,1,1,1,0,0,0,\cdots ,0\right)}^{\text{T}}$ ， ${\beta }_0$ 的前s个系数为1，其余为0。预测变量之间的自相关系数 $\rho =0.5$ ，均值向量 $\mu =1_p^{\text{T}}$ 。设计矩阵 $X_{n\times p}$ 中每一个观测分别从以下几个多元分布中独立生成：1) 多元正态分布： $X~N\left(\mu ,\Sigma \right)$ ， ${\Sigma }_{ij}={\rho }^{\left|i-j\right|}$ (下同)，记为G；2) 多元t分布： $X~t_d\left(\mu ,\Sigma \right)$ ， $d=1,3,10$ 为自由度，对应数据集分别记为T1，T3，T10。然后按照模型(1)生成响应变量，随机模拟实验次数为 $B=100$ 。

由表1可以看出，在不同的数据分布中，基于不同权重的自适应lasso的预测效果及变量选择结果均优于lasso。从模型预测的角度来说，自适应lasso的系数估计的MSE均小于lasso的MSE，且自适应lasso方法的分类正确率均大于lasso方法。从变量选择的角度来说，四种方法的敏感性均为1，表明四种变量选择方法均能正确识别出真实模型中的非零系数(或重要变量)，但在四个数据集中，lasso的特异性均低于自适应lasso方法，说明lasso方法更容易将系数为零0的变量选入模型中，模型的可解释性没有自适应lasso的好，这也反映了自适应lasso方法的神谕性。值得注意的是，基于不同权重的自适应lasso方法适用于不同数据分布。正态分布中，APLR效果最好，在T1和T3中，LAPLR最优，而在T10中则是CAPLR最好。

3.2. 高维数据情形

本节研究高维数据变量选择的问题，即 $p>n$ 。由于极大似然估计不适用于高维数据的情形，因此在一部分我们只比较lasso、LAPLR、CAPLR的变量选择和预测效果。

Table 1. The variable selection and prediction results in low dimension setting

表1. 低维情形下的变量选择及预测结果

随机模型试验的参数设置如下： $n=500,p=1000,s=5$ ， ${\beta }_0$ 的前五个元素等于1，其余为0。其余参数与低维情形相同。

由表2可知，在高维数据中，两种自适应lasso方法的预测及变量选择效果均优于lasso方法。从模型预测的角度来说，两种自适应lasso系数估计的均方误差MSE更小，分类正确率更高。从变量选择的角度来说，四种方法的敏感性都等于1，表明这两种方法均均能正确识别真实模型中的重要变量，而自LAPLR、CAPLR的特异性在四个数据分布中都等于1，表明自适应lasso能完全识别出真实模型，具有很好的神谕性，模型的可解释性要优于lasso。

Table 2. The variable selection and prediction results in high dimension setting

表2. 高维情形下的变量选择及预测结果

4. 实证分析

我们用数据集“Pima Indians Diabetes”对上述几种方法的预测和变量选择效果进行比较。Pima Indians Diabetes数据集包含392个医疗记录，每个记录包含8个医学检测指标以及糖尿病检测结果。数据集的变量描述如表3所示，“diabetes”为响应变量，取值为 $\left\{0,1\right\}$ ，其余变量为预测变量。

Table 3. The information of Pima Indians Diabetes

表3. 皮马印第安人糖尿病数据集信息

*注：Pima Indians Diabetes来源于UCI机器学习数据库。

在实证分析中，本文按80%、20%的比例将Pima Indians Diabetes数据集划分为训练集，测试集。在训练集中，我们使用10重交叉验证，选择出最优的模型，然后在测试集中用该模型进行预测。由表4可知，lasso、APLR、CAPLR选出的变量个数均为6个，而APLR选出了5个变量。四种方法共同选出的变量有“pregnant”、“glucose”、“triceps”、“mass”、“pedigree”，从预测角度来说，自适应lasso方法的预测准确率均优于lasso方法，其中LAPLR的预测准确率最高，为79.21%。

Table 4. Prediction and variable selection of the four methods

表4. 四种方法预测及变量选择结果

^*选择的变量个数不包含截距项。

5. 研究结论

本文通过随机模拟实验及实际数据分别比较了lasso、APLR、LAPLR、CAPLR这四种变量选择方法的预测及变量选择效果。随机模拟实验的结果表明：在低维和高维的数据中，APLR、LAPLR、CAPLR的预测及变量选择能力均优于lasso，但是APLR、LAPLR、CAPLR方法之间并不存在绝对最优的方法，在不同的数据集中，这三种方法各有优势。实证分析结果与随机模拟的结果基本一致，即：相较与lasso而言，三种自适应lasso方法的预测准确率更高，APLR方法得到的模型更稀疏，可解释性更好。值得强调的是，不同的变量选择方法适用于不同的数据集，在实际分析中，要对多种方法进行综合比较，从而选择出最适合的方法。

基金项目

北京市属高校基本科研业务费NO.110052971921/103资助项目。