随着科学技术的进步, 不同领域的数据都呈现出网络连接的趋势, 许多科学领域都涉及某种形式的网络研究, 例如人际关系研究、学术论文合著和引用、蛋白质相互作用模式等. 20年前, 关于网络的流行书籍及其研究开始出现^[1], 而像Facebook、MySpace和LinkedIn这样的在线网络社区在近10年间也是蓬勃兴起, 这更加增强了人们对网络数据的研究兴趣. 网络连接数据由节点和边组成, 社交网络是此类网络模型的一个典型代表. 社交网络中, 每个节点代表一个人, 边代表人与人之间的沟通交流, 此外, 还有商业网络、基因网络等.

目前关于网络连接数据的研究主要分为两个方面. 一方面是关于网络结构的研究. 另一方面主要是将网络连接数据中的结构信息与统计学习中常用的经典模型结合起来研究.

在网络结构方面, 最早被应用于社区检测. 社区检测兴起于物理学和计算机科学领域, 而后开始应用于统计领域. 其中一类社区检测算法是通过在节点的所有可能分区上优化启发式全局准则来检测社区^{[2, 3]}. 基于概率模型的方法^{[4, 5]}是另一类社区检测算法. 一些学者从观察到的邻接矩阵中检测社区或潜在结构^[6-8], 从其他节点之间的信息估计特定节点之间的边缘概率^[9]. 社交网络是此类网络模型的代表, 因此针对社交网络的研究也受到了大量的关注^{[10, 11]}.

在与经典模型结合方面, 一般是与常用的模型相结合. 例如, 时间序列模型^[12], 线性模型^[13], 变系数模型^[14], 随机效应模型^[15], 变化点检测问题^[16], 自回归模型^{[17, 18]}等.

线性回归模型是统计学习中的经典模型之一, 应用十分广泛, 关于网络数据的回归模型也开始引起学者的关注. 例如, Asur等^[19]将网络数据应用于预测模型, 通过研究网络结构来预测现实生活中某一现象的结果. Li等^[13]将网络连接数据应用于回归预测模型, Zhu等^[17]和Tang等^[18]将网络连接数据与自回归模型相结合, 都表明网络连接数据在回归模型中的研究价值. 随着科技的发展, 数据的采集变得更加容易, 高维数据也越来越受到研究学者的关注, 但是高维数据中存在大量的冗余信息, 如何选出有研究价值的数据？变量选择领域应运而生. 故将网络连接数据应用到变量选择领域是一个值得研究的课题.

对于线性回归模型, 超高的维度使得传统的普通最小二乘法不再适用. 正则化是稀疏建模和变量选择的有效方法, 通过在目标函数上添加惩罚函数来降低模型的复杂度. 根据惩罚函数的不同, 正则化方法一般可以分为凸正则化和非凸正则化.

凸正则化方法主要包括岭回归、LASSO、弹性网以及Dantzig Selector等. 虽然凸正则化的研究已经很成熟, 但由于惩罚函数的凸性, 使得凸正则化估计量都是有偏的. Fan等提出了一个非凸正则化方法—SCAD(smoothly clipped absolute deviation)^[20], 并证明了其Oracle性质. 非凸惩罚函数回归的渐进无偏估计, 能进一步降低模型的预测总误差. 此后, 非凸惩罚受到了广泛的关注, 例如MCP (minimax concave penalty)^[21]、限制Capped- ${L_1}$ ^[22]、Hard阈值惩罚^[23]等.

关于网络连接数据的变量选择问题近年来也有学者做过相关研究^{[24, 25]}. 例如Li等^[24]和Kim等^[25]考虑样本系数之间的网络凝聚效应, 即网络中连接节点表现出相似的行为, 对系数同时施加了 ${L_1}$ 惩罚和凝聚效应惩罚 $\;{{\boldsymbol{\beta}} ^{\rm {T}}}{\boldsymbol{L}}{\boldsymbol{\beta}}$ , 从而能够解决网络连接数据的变量选择问题, 但他们针对的是同质性网络连接数据, 即假设每个样本的个体效应值 ${\boldsymbol{\alpha}}$ 相同, 并没有考虑到异质性, 异质性是指不同样本的个体效应 ${\boldsymbol{\alpha}}$ 不同. 在现实生活中, 因为网络凝聚效应的存在而使得网络中的样本存在群组效应, 联系密切的样本组成一个群组, 他们之间的行为会相互影响而慢慢趋同. 针对线性回归模型, 这种群组效应的一个直观体现就是群组内样本的个体效应 ${\boldsymbol{\alpha}}$ 相同, 不同群组间个体效应 ${\boldsymbol{\alpha}}$ 不同. 若忽略群组间个体效应的差异性, 将所有样本的个体效应视为相同, 在进行变量选择和预测估计时都会产生较大偏差, 影响模型精度. 故考虑异质性, 能够提高模型精度. 因此, 针对异质性网络连接数据的研究具有重要的价值和实际意义. Li等^[13]考虑到个体效应之间的异质性, 并惩罚相连样本个体效应的差异性, 提高了回归模型中估计和预测的精度, 但他主要关注的是预测问题, 没有涉及到变量选择.

本文的目标是对因网络凝聚效应而产生个体效应的组异质性的网络连接数据进行变量选择, 我们对组内样本间个体效应的差异性 ${\boldsymbol{L\alpha}}$ 和变量系数 $\;{\boldsymbol{\beta}}$ 进行联合惩罚, 从而保证组内样本的个体效应具有相同的估计值. 本文提出的方法不仅能够处理含有组异质性的网络连接数据的变量选择问题, 而且能够改善变量选择、估计和预测的结果. 在本文中, 我们主要使用 ${L_1}$ 、MCP和SCAD罚函数, 并且运用ADMM算法进行求解, 同时证明了算法的收敛性.

1 网络连接数据的变量选择方法 1.1 模型设定

本文中所有的向量都是列向量. 考虑一般的线性回归模型, ${\boldsymbol{Y}}{\text{ = (}}{{{y}}_1}, {{{y}}_2}, \cdots , {{{y}}_{_n}}{)^{\rm {T}}}$ 是 $n$ 维响应变量, ${\boldsymbol{X}}{\text{ = (}}{{{x}}_1}, $ $ {{{x}}_2}, \cdots , {{{x}}_{_n}}{)^{\rm T}}$ 是 $n \times p$ 设计矩阵. 假设 $ {\boldsymbol{X}} $ 是固定的且其列已经标准化. 样本 $ {\boldsymbol{X}} $ 的结构网络为 $G = (V, E)$ , 其中 $V = \{ 1, 2, \cdots , $ $ n\} $ 为样本节点集合, $E \subset V \times V$ 为边的集合. 我们用邻接矩阵 ${\boldsymbol{A}} = {({A_{uv}})_{n \times n}} \in {R_{n \times n}}$ 表示该网络以及样本节点和节点之间的连接关系, 若 $(u, v) \in E$ , 则 ${A_{uv}} = 1$ , 否则为0. ${A_{uu}} = 0, {A_{uv}} = {A_{vu}}$ . 网络 $G$ 的拉普拉斯矩阵 ${\boldsymbol{L}} = {\boldsymbol{D}} - {\boldsymbol{A}}$ , ${\boldsymbol{D}} = diag({d_1}, {d_2}, \cdots , {d_n})$ 为度矩阵, ${\boldsymbol{D}}$ 的对角线元素为每个节点的度 ${d_u} = \displaystyle\sum\limits_v {{A_{uv}}} $ . 建立如下线性回归模型:

$ {\boldsymbol{Y}} = {\boldsymbol{\alpha}} + {\boldsymbol{X\beta}} + {\boldsymbol{\varepsilon}} $

(1)

其中, ${\boldsymbol{\alpha}} = {({{\boldsymbol{\alpha}} _1}, {{\boldsymbol{\alpha}} _2}, \cdots , {{\boldsymbol{\alpha}} _n})^{\rm T}}$ 是节点个体效应向量. 假设相连样本的个体效应相等, 不相连样本的个体效应不等, 即样本之间存在组异质性. $\;{\boldsymbol{\beta}} = {({{\boldsymbol{\beta}} _1}, {{\boldsymbol{\beta}} _2}, \cdots , {{\boldsymbol{\beta}} _p})^{\rm T}}$ 是模型的回归系数向量. ${\boldsymbol{\varepsilon}} = {({{\boldsymbol{\varepsilon}} _1}, {{\boldsymbol{\varepsilon}} _2}, \cdots , {{\boldsymbol{\varepsilon}} _n})^{\rm T}}$ 是 $n$ 维误差向量, $E({\boldsymbol{\varepsilon}} ) = {\boldsymbol{0}}, \; {\boldsymbol{var}}({\boldsymbol{\varepsilon}} ) = {\sigma ^2}{{\boldsymbol{I}}_n}$ .

Li等^[13]提出了网络连接数据的预测方法(the regression with network cohesion, RNC), 其主要思想是最小化如下损失函数:

$ {\boldsymbol{L}}({\boldsymbol{\alpha}} , {\boldsymbol{\beta}} ) = ||{\boldsymbol{Y}} - {\boldsymbol{\alpha}} - {\boldsymbol{X\beta}} |{|^2} + {\boldsymbol{\mu}} {{\boldsymbol{\alpha}} ^{\rm T}}{\boldsymbol{L\alpha}} $

(2)

其中, $\;{{\mu}} > 0$ 是调整参数. RNC主要是惩罚网络中相连节点个体效应的差异性, 该惩罚可以推导出一个等价的、更直观的形式 ${{\boldsymbol{\alpha}} ^{\rm T}}{\boldsymbol{L\alpha}} =\displaystyle \sum\limits_{(u, v) \in E} ( {{\boldsymbol{\alpha}} _u} - {{\boldsymbol{\alpha}} _v}{)^2}$ .

1.2 组异质性网络连接数据的变量选择方法(SNC)

RNC中假设各样本的个体效应不相等, 惩罚项 $\;{{\mu}} {{\boldsymbol{\alpha}} ^{\rm T}}{\boldsymbol{L\alpha}}$ 用来惩罚相连样本个体效应的差异性, 从而能够处理异质性网络连接数据的回归预测问题. 但是由于连接网络中的个体常常存在群组效应, 同一个群组中行为特征存在统一准则而基本相同. 因此, 在本文中我们假设样本之间存在组异质性, 即组内样本(相连样本)的个体效应相等, 组间样本(不相连样本)的个体效应不相等. 通过对 ${\boldsymbol{L\alpha}} $ 施加惩罚, 惩罚组内样本个体效应的差异性并压缩至0, ${\boldsymbol{L\alpha}} $ 中的元素是 ${({{{\alpha}} _u} - {{{\alpha}} _v})_{(u, v) \in E}}$ 或其等价形式. 为了产生 $\;{\boldsymbol{\beta}}$ 的一个稀疏估计, 我们将同时惩罚 $\;{\boldsymbol{\beta}}$ 和 ${\boldsymbol{L\alpha}} $ , 这就是我们提出的方法—网络连接数据的变量选择(variable selection with network cohesion, SNC).

令 $ {\boldsymbol{\theta}} = {({{\boldsymbol{\beta}} ^{\rm T}}, {{\boldsymbol{\alpha}} ^{\rm T}})^{\rm T}},\; {\boldsymbol{H}} = \left( {\begin{array}{*{20}{c}} {\boldsymbol{I}}& {\boldsymbol{0}} \\ {\boldsymbol{0}} & {\boldsymbol{L}} \end{array}} \right) $ , 则 ${\boldsymbol{H}}{\boldsymbol{\theta}} = \left( \begin{array}{c} {\boldsymbol{\beta}} \\ {\boldsymbol{L\alpha}} \end{array} \right)$ . SNC的目标函数为:

$ Q({\boldsymbol{\alpha}} , {\boldsymbol{\beta}} ) = \frac{1}{{2n}}||{\boldsymbol{Y}} - ({\boldsymbol{X}}, {\boldsymbol{I}}){\boldsymbol{\theta}} ||_2^2 + {p_\lambda }(|{\boldsymbol{H}}{\boldsymbol{\theta}} |) $

(3)

在本文中, 对 ${\boldsymbol{H\theta}} $ 的惩罚主要使用 ${L_1}$ 和非凸惩罚, 非凸惩罚包括MCP和SCAD罚函数. MCP罚函数为 ${p_M}(t, \lambda ) = \lambda \int_0^t \left(1 - \dfrac{x}{{a\lambda }}\right) {_ + }dx,\; a > 0$ , SCAD罚函数为 ${p_S}(t, \lambda ) $ $ = \lambda \int_0^t \min \{ 1, \left(a - \dfrac{x}{\lambda } \right){_ + }/(a - 1)\} dx, \;a > 2$ .

将SNC方法的估计结果与没有对节点个体效应的差异进行惩罚的情况下进行对比, 能够提高估计和预测的精度.

2 算法

直接最小化目标函数(3)很难求解出估计量的值, 因为惩罚函数对于每个 ${{{\alpha}} _i}$ 是不可分的. 因此, 我们通过引入一组新的参数 ${\boldsymbol{\gamma}} = {\boldsymbol{H\theta}} $ 来重新参数化准则. 最小化式(3)等价于最小化如下约束优化问题:

$\left\{ \begin{array}{l} Q({\boldsymbol{\theta}} , {\boldsymbol{\gamma}} ) = \dfrac{1}{{2n}}||{\boldsymbol{Y}} - ({\boldsymbol{X}}, {\boldsymbol{I}}){\boldsymbol{\theta}} ||_2^2 + {p_\lambda }(|{\boldsymbol{\gamma}} |)\\ {\rm {s.t.}}\;\;\;\;\;\;{\boldsymbol{H\theta}} - {\boldsymbol{\gamma}} = {\boldsymbol{0}} \end{array} \right.$

基于文献[26]中的思路, 利用增广拉格朗日方法, 通过最小化如下损失函数得到参数的估计:

$ L({\boldsymbol{\theta}} , {\boldsymbol{\gamma}} , {\boldsymbol{\varphi}} ) = Q({\boldsymbol{\theta}} , {\boldsymbol{\gamma}} ) + {{\boldsymbol{\varphi}} ^{\rm T}}({\boldsymbol{H\theta}} - {\boldsymbol{\gamma}} ) + \rho /2||{\boldsymbol{H\theta}} - {\boldsymbol{\gamma}} ||_2^2 $

(4)

其中, 对偶变量 ${\boldsymbol{\varphi}} $ 是拉格朗日乘数, $\;\rho > 0$ 是惩罚因子. 我们通过交替方向乘子法(alternating direction multiplier method, ADMM)来迭代求解 $({\boldsymbol{\theta}} , {\boldsymbol{\gamma}} , {\boldsymbol{\varphi}} )$ 的估计. 对于给定的 $({\boldsymbol{\theta}} , {\boldsymbol{\gamma}} , {\boldsymbol{\varphi}} )$ , $L({\boldsymbol{\theta}} , {\boldsymbol{\gamma}} , {\boldsymbol{\varphi}} )$ 关于 ${\boldsymbol{\gamma}} $ 的最小值是唯一的, 并且在 ${L_1}$ 惩罚或非凸惩罚下有一个近似的形式. 当给定 $({\boldsymbol{\theta}} , {\boldsymbol{\gamma}} , {\boldsymbol{\varphi}} )$ , 上述最小化问题等价于:

$ \rho /2||{\boldsymbol{\tau}} - {\boldsymbol{\gamma}} ||_2^2 + {p_\lambda }(|{\boldsymbol{\gamma}} |) $

其中, ${\boldsymbol{\tau}} = {\boldsymbol{H\theta}} + {\rho ^{ - 1}}{\boldsymbol{\varphi}} $ , 故在 ${L_1}$ 或非凸惩罚下估计量的近似的形式为:

$ \hat {\boldsymbol{\gamma}} = ST({\boldsymbol{\tau}} , \lambda /\rho ), $

其中, $ST(t, \lambda ) = sign(t){(|t| - \lambda )_ + }$ 是soft阈值准则, ${(x)_ + } = $ $ x, x > 0$ , 否则 ${(x)_ + } = 0$ .

对于MCP罚函数 $(a > 1/\rho )$ , 有:

$ {\hat{\gamma} _{ij}} = \left\{ {\begin{array}{*{20}{c}} \begin{array}{l} ST({{{\tau}} _{ij}}, \lambda /\rho ), \\ {{{\tau}} _{ij}}, \end{array} &\begin{array}{l} |{{{\tau}} _{ij}}| \le a\lambda \\ |{{{\tau}} _{ij}}| > a\lambda \end{array} \end{array}} \right. $

对于SCAD罚函数 $(a > 1/\rho + 1)$ , 有:

$ {\hat{\gamma} _{ij}} = \left\{ {\begin{array}{*{20}{c}} \begin{array}{l} ST({{{\tau}} _{ij}}, \lambda /\rho ), \\ \dfrac{{ST({{{\tau}} _{ij}}, a\lambda /(a - 1)\rho )}}{{1 - 1/((a - 1)\rho )}}, \\ {{{\tau}} _{ij}}, \end{array} &\begin{array}{l} |{{{\tau}} _{ij}}| \le \lambda + \lambda /\rho \\ \lambda + \lambda /\rho < |{{{\tau}} _{ij}}| \le a\lambda \\ |{{{\tau}} _{ij}}| > a\lambda \end{array} \end{array}} \right. $

算法步骤如算法1.

算法1. ADMM算法

输入: 预测变量 $\scriptstyle {\boldsymbol{X}}$ , 响应变量 $\scriptstyle {\boldsymbol{Y}}$ , 邻接矩阵 $\scriptstyle {\boldsymbol{A}}$ , 惩罚因子 $\scriptstyle\; \rho $ , 停止准则 $\scriptstyle \eta $ ;

输出: $ \scriptstyle {\widehat{\boldsymbol{\theta}} }, {\widehat{\boldsymbol{\gamma}} }, {\widehat{\boldsymbol{\varphi}} } $ ;

目标: 迭代求解获得 $ \scriptstyle {\boldsymbol{\theta}} , {\boldsymbol{\gamma}} $ 和 $\scriptstyle {\boldsymbol{\varphi}} $ .

初始化 $\scriptstyle {{\boldsymbol{\theta}} }^{(0)}, {{\boldsymbol{\gamma}} }^{(0)}={{\boldsymbol{H\theta}} }^{(0)}, {{\boldsymbol{\varphi}} }^{(0)}={\boldsymbol 0}, m=0, \eta =\mathrm{0.03.}$

While $\scriptstyle m \ge 0$ , do

　　　　 $\scriptstyle {{\boldsymbol{\theta}} ^{(m + 1)}} = {[{{{n}}^{ - 1}}{({\boldsymbol{X}}, {\boldsymbol{I}})^{\rm T}}({\boldsymbol{X}}, {\boldsymbol{I}}) + \rho {{\boldsymbol{H}}^{\rm T}}{\boldsymbol{H}}]^{ - 1}}$

　　　　　　 $\scriptstyle *[{n}^{-1}{({\boldsymbol{X}}, {\boldsymbol{I}})}^{\rm T}{\boldsymbol{Y}}+\rho {{\boldsymbol{H}}}^{\rm T}{{\boldsymbol{\gamma}} }^{(m)}-{{\boldsymbol{H}}}^{\rm T}{{\boldsymbol{\varphi}} }^{(m)}];$

　　　　 $\scriptstyle {{\boldsymbol{\gamma}} ^{(m + 1)}} = ST({{\boldsymbol{\tau}} ^{(m + 1)}}, \lambda /\rho )$ ;

　　　　 $\scriptstyle {{\boldsymbol{\varphi}} ^{(m + 1)}} = {{\boldsymbol{\varphi}} ^{(m)}} + \rho ({{\boldsymbol{H\theta}} ^{(m + 1)}} - {{\boldsymbol{\gamma}} ^{(m + 1)}})$ .

　　　If $\scriptstyle {{\boldsymbol{r}}^{(m + 1)}} = {{\boldsymbol{H\theta}} ^{(m + 1)}} - {{\boldsymbol{\gamma}} ^{(m + 1)}}$ , $\scriptstyle ||{{\boldsymbol{r}}^{(m + 1)}}|| < \eta $

　　　　then

　　　　　　 $ \scriptstyle ({\widehat{\boldsymbol{\theta}} }, {\widehat{\boldsymbol{\gamma}} }, {\widehat{\boldsymbol{\varphi}} })=({{\boldsymbol{\theta}} }^{(m+1)}, {{\boldsymbol{\gamma}} }^{(m+1)}, {{\boldsymbol{\varphi}} }^{(m+1)}) $ ;

　　　　　　Break;

　　　Else

　　　　　　 $\scriptstyle m = m + 1;$

　　　End

End

对ADMM算法过程中的原始变量进行追踪, ${{\boldsymbol{r}}^{(m + 1)}} = {{\boldsymbol{H\theta}} ^{(m + 1)}} - {{\boldsymbol{\gamma}} ^{(m + 1)}}$ . 停止准则为 $||{{\boldsymbol{r}}^{(m + 1)}}|| < \eta $ , 其中 $\eta > 0$ 为一个非常小的常数.

下面考虑ADMM算法的收敛性.

命题1. 对于MCP和SCAD函数, ADMM算法的原始残差 ${{\boldsymbol{r}}^{(m)}} = {{\boldsymbol{H\theta}} ^{(m)}} - {{\boldsymbol{\gamma}} ^{(m)}}$ 和对偶残差 ${{\boldsymbol{s}}^{(m + 1)}} = \rho {{\boldsymbol{H}}^{\rm T}} $ $ ({{\boldsymbol{r}}^{(m + 1)}} - {{\boldsymbol{r}}^{(m)}})$ 满足 $ \mathop {\lim }\limits_{m \to \infty } ||{{\boldsymbol{r}}^{(m)}}|{|^2} = 0 $ , $ \mathop {\lim }\limits_{m \to \infty } ||{{\boldsymbol{s}}^{(m)}}|{|^2} = 0 $ .

命题1表明该算法实现了原可行性和对偶可行性, 证明材料见附录. 因此, 它收敛于一个局部最优点. 当采用非凸惩罚函数, 如MCP和SCAD罚函数时, 此最优点是目标函数的局部最优解. 综上, 算法收敛性和稳定性得到证明. 因为 ${{\boldsymbol{\theta}} ^{(m)}} = {({({{\boldsymbol{\beta}} ^{(m)}})^{\rm T}}, {({{\boldsymbol{\alpha}} ^{(m)}})^{\rm T}})^{\rm T}}$ 是不稀疏的, 但我们已证明 ${{\boldsymbol{H\theta}} ^{(m)}} = {({({{\boldsymbol{\beta}} ^{(m)}})^{\rm T}}, {({{\boldsymbol{L\alpha}} ^{(m)}})^{\rm T}})^{\rm T}}$ 是收敛于 ${{\boldsymbol{\gamma}} ^{(m)}}$ , 故我们令 ${{\boldsymbol{\gamma}} ^{(m)}}$ 的前 ${{p}}$ 项作为 ${\boldsymbol{\beta}} $ 的估计值, 即可得到 $\;{\boldsymbol{\beta}} $ 的稀疏解.

3 数值模拟

在数值模拟中, 主要比较本文提出的SNC方法和没有对个体节点效应的差异性进行惩罚的LASSO、MCP、SCAD方法在变量选择和预测方面的效果. 网络凝聚效应下的变量选择方法就是考虑了样本之间的连接关系网络的方法, 即我们的SNC方法.无网络凝聚效应下的变量选择方法, 就是不考虑样本之间的连接网络的惩罚方法. 在这里, 我们首先定义几个效果评估指标:

(1)预测损失(prediction error, PE): ${\boldsymbol{E}}({{\boldsymbol{X}}^{\rm T}}{{\boldsymbol{\beta}} _0} + {{\boldsymbol{\alpha}} _0} - $ $ {{\boldsymbol{X}}^{\rm T}}\widehat {\boldsymbol{\beta}} - \widehat {\boldsymbol{\alpha}} )^{^2}$ ;

(2) ${L_q}$ 损失: $||\widehat {{{\boldsymbol{\beta}}}} - {\boldsymbol{\beta}}_0 |{|_q}, q = 1, 2, \infty $ ;

(3)均方误差(MSE( $\widehat {\boldsymbol{\alpha}} $ )): $E{(\widehat {\boldsymbol{\alpha}} - {{\boldsymbol{\alpha}} _0})^2}$ ;

(4)假阳性数(false positives, FP): 真实为反例却被预测为正例的个数;

(5)假阴性数(false negatives, FN): 真实为正例却被预测为反例的个数;

(6)真阳性数(true positives, TP): 真实为正例预测也为正例的个数;

(7)真阴性数(true negatives, TN): 真实为反例预测也为反例的个数;

(8) ${F_1}$ -score: 2TP/(2TP+FP+FN).

3.1 模拟1

对于式(1)中的线性回归模型, 我们从该模型中随机生成100个数据集. 训练样本的大小考虑两种情况 $(n, p) = (100, 200)$ 和 $(n, p) = (100, 500)$ , 设计矩阵 ${\boldsymbol{X}}$ 中的每一行从正态分布 ${{N}}({\boldsymbol{0}}, {\boldsymbol{\Sigma}} ), {\boldsymbol{\Sigma}} = {({0.5^{|i - j|}})_{1 \le i, j \le p}}$ 中随机抽样. 真实回归系数 $\;{{\boldsymbol{\beta}} _0} = {(\underbrace {0.6, \cdots , 0.6}_{10}, {{\boldsymbol{0}}^{\rm T}})^{\rm T}}$ , $\sigma = 0.3$ 为随机误差 ${\boldsymbol{\varepsilon}} $ 的标准差. 惩罚因子 $\;\rho = 1$ , $\lambda $ 用交叉验证来选取, 停止条件 $\eta = 0.03$ .

为了生成含有组异质性样本间的邻接矩阵 ${\boldsymbol{A}}$ , 我们用ER随机图模型生成一个包含 $n = 100$ 个节点的样本网络, 样本网络由4个不相连的部分 ${G_1}, {G_2}, {G_3}, {G_4}$ 组成, 每个部分包含25个节点. 每个单独的部分都是一个ER随机图, 节点与节点之间以 ${{{p}}_b}$ 的概率生成边, 即 ${A_{ij}} = 1$ , 否则为0, 令 ${{{p}}_b} = 0.1$ . 4个部分中相连样本的个体节点效应 ${{{\alpha}} _{{i}}}$ 的值分别为1, −1, 0.5, −0.5, 独立样本的个体节点效应为0.3.

表1展示了两种方法在预测评估指标上的结果对比. 与没有利用相连节点的网络凝聚效应对个体效应进行惩罚的LASSO、MCP和SCAD结果相比, SNC-LASSO、SNC-MCP和SNC-SCAD都明显改善了估计和预测误差. 这表明将网络凝聚效应加入变量选择模型中, 可以改善模型变量选择、估计和预测的精度.

表2展示了两种方法在100次模拟实验下变量选择评估指标结果. 我们可以看出各项指标下, SNC方法的变量选择效果都明显优于没有利用网络凝聚效应进行惩罚的方法. 另外, SNC-MCP和SNC-SCAD都要优于SNC-LASSO. 尤其对于假阳性数FP, 100次模拟中, SNC-LASSO的FP平均为15.41 ( $p$ =200)和17.21 ( $p$ =500), 而SNC-MCP分别为0.05 ( $p$ =200)和0.3 ( $p$ =500), SNC-SCAD分别为1.06 ( $p$ =200) 和0.2 ( $p$ =500), MCP和SCAD变量选择的准确性比LASSO显著提高, 主要是由于LASSO的有偏性.

3.2 模拟2

模拟1中的结果表明网络凝聚效应惩罚能够改善变量选择、估计和预测效果, 网络凝聚效应主要与邻接矩阵中个体之间产生联系的概率 ${{{p}}_b}$ 有关, 接下来我们将研究 ${{{p}}_b}$ 对SNC方法的变量选择、估计和预测效果的影响. 模型2中的设定与模型1类似, 不同的是我们取 ${{{p}}_b} = seq(0, 0.02, 0.2)$ , R语言函数seq(a, b, c)用于生成一组从a到b, 间隔为c的序列.

图1和图2分别展示了 ${p_b}$ 对预测和变量选择效果的影响. 从图1可以看出, 随着 ${p_b}$ 的增大, 即网络的凝聚效应增强, SNC方法能够明显降低预测损失, 并在 ${p_b} = 0.08$ 附近趋于稳定. 图2表示 ${p_b}$ 对 ${F_1}$ 分数的影响, ${F_1}$ 分数是查准率和查全率的调和平均数, 当 ${p_b} = 0$ 即样本之间没有连接关系时, ${F_1}$ 分数值很低. 随着 ${p_b}$ 的增大, ${F_1}$ 分数值逐渐增大, 同样地, 在 ${p_b} = 0.08$ 附近达到最大值, 此时SNC方法变量选择的效果较好.

4 实际数据分析

我们研究的真实数据案例来自于Teenagers Friends and Lifestyle Study^[27]. 这项研究主要是青少年友谊网对他们自身某些行为的影响. 该实际数据与本文中的模型设定保持一致, 因青少年时期学生喜爱团体活动, 故凝聚效应使得网络之间存在组异质性.

Teenagers Friends and Lifestyle Study旨在确定在青少年早期到中期不良习性的变化过程. 实验记录了3个时间点 $ \left( {{T_1},{T_2},{T_3}} \right)$ 的数据. 样本来自于160名学生, 通过每个学生及学生之间的朋友(最多6个)关系来建立友谊网络. 研究给出了3个时期的友谊网络, 网络中“1”表示“best friend”, “2”表示“just a friend”, “0”表示“no friend”, “10”表示缺失值, 我们根据学生之间的友谊网络来获取邻接矩阵 ${\boldsymbol{A}}$ .

本文使用的数据集 ${\boldsymbol{X}}$ 包含160个样本, 40个特征变量包括青少年的年龄、性别、生活方式、休闲活动以及家庭成员吸烟等情况, 考虑特征之间的交互作用, 最终特征变量为250个. 我们的目标是利用友谊网络找出影响青少年不良习性的关键因素, 并预测青少年自身不良行为的活动频率. 我们分别选取alcohol、tobacco和cannabis作为响应变量 ${\boldsymbol{Y}}$ , 对于tobacco, 元素1表示从未抽过烟, 2表示偶尔吸烟, 3表示经常吸烟, 故我们将其取对数作为响应变量 ${\boldsymbol{Y}}$ 的值.

时间点 ${T_1}$ 的友谊网络如图3所示. 我们只展示了学生之间的“best friendship”(包括“just a friend”和“best friend”). 根据友谊网络建立邻接矩阵 ${\boldsymbol{A}}$ 时, 当学生 $i$ 和学生 $j$ 为“best friend”, 则 ${A_{ij}} = {A_{ji}} = 1$ , 否则 ${A_{ij}} = 1$ .

分别选取alcohol、tobacco和cannabis作为响应变量来研究影响青少年酗酒、吸烟和吸毒的因素. 将样本随机分成两份: 训练集和测试集, 重复实验100次. 由于不知道真实情况下的参数设定, 无法像模拟实验中那样对比假阴性数、假阳性数等指标. 因此, 主要从预测损失和变量选择两个方面来验证SNC方法的有效性.

表3展示了SNC方法SNC-LASSO、SNC-MCP、SNC-SCAD与无网络凝聚效应下的变量选择方法LASSO、MCP和SCAD对青少年不良习性(酗酒、抽烟以及吸食大麻)的预测损失, 从结果中可以看出SNC方法预测的相对更准确一点. 青少年时期大家都是团体活动, 生活习惯很容易相互影响而慢慢趋同, 而网络凝聚效应正是考虑了这一点, 团体内个体的表现行为更具相似性, 惩罚团体内个体效应的差异性, 提高了个体效应的预测精度, 从而降低了整个模型的预测误差.

为了使挑选出来的变量更具可解释性, 下面我们不考虑特征之间的交互作用, 用SNC方法和无网络凝聚效应下的变量选择方法来挑选变量, 并重复实验100次, 计算100次实验下挑选出来的变量的比例.

表4中我们看到, LASSO、MCP和SCAD挑选出更多的冗余变量. 显然, 两种方法下, 特征变量parent smoking, sibling smoking, “I hang round in the streets”, “I play computer games”和“I go to dance clubs or raves”是最显著的. 青少年时期他们的世界观、人生观和价值观还在形成阶段, 易受他人或团体的影响, 在街上闲逛、经常打电脑游戏、参加俱乐部以及兄弟姐妹抽烟等行为都容易使青少年沾染上不良习性. 通过研究分析, 我们知道了青少年时期朋友以及家人行为的重要性, 家人、朋友以及整个社会需要给青少年营造一个良好健康的成长环境, 给他们树立积极向上的榜样.

针对各种方法挑选出来变量之后的模型进行回归, 我们得到回归后各变量系数的显著性检验以及调整可决系数 ${R^2}$ 和标准误差如表5所示.

由表5可知, SNC方法选取了sex.F、I hang out in the streets、I play computer games、money、parent.smoking和sibling.smoking 6个变量, 根据值可以看出这些变量都通过了显著性检验. 而LASSO、MCP和SCAD方法选出了少许的冗余变量. 另外, 从表中的调整可决系数和标准误差来看, SNC方法的效果也是优于没有网络凝聚效应下的变量选择方法.

5 总结

本文主要对线性回归模型中因网络凝聚效应而产生个体效应的组异质性的网络连接数据进行变量选择, 使用非凸惩罚MCP和SCAD罚函数同时惩罚变量系数 $\;{\boldsymbol{\beta}} $ 和组内样本的个体效应的差异性 ${\boldsymbol{L\alpha}} $ , 使得能够对含有组异质性的网络连接数据筛选出有用变量.

针对本文提出的方法, 我们运用ADMM算法进行求解, 并证明了算法的收敛性. 针对SNC方法, 本文进行了相关模拟, 从变量选择和预测两个方面来衡量该方法的效果. 从实验结果来看, 无论是预测损失还是变量选择的准确性都有明显改善. 实例分析中, 我们将SNC方法应用于青少年友谊网络和生活方式的研究, 分析预测青少年吸烟等不良习性的活动频率以及挑选出影响青少年吸烟等不良习性的特征变量.

本文提出的方法, 为含有组异质性网络连接数据的变量选择问题提供了一种解决思路. 我们将变量选择方法进一步拓展了应用领域, 对于基因网络、交通网络、公司网络等网络连接数据, SNC方法都能适用.