网络连接数据即为包含数据点之间关联信息的网状结构数据. 1991年, CERFNet、PSINet及AlterNet商业化, 世界步入互联网时代. 随着计算机技术的发展和物联网的广泛应用, 数据的增长速度和应用范围已达到前所未有的高度. 根据思科(Cisco)的2023年全球网络趋势报告, 全球数据流量在过去5年中增长了近4倍. IDC预测到2025年将达到每年163 ZB^[1]. 伴随数据量的爆炸式扩张, 数据的种类和结构也日益多元化. 不止单个人或物的属性信息被详细记录, 人与人、人与物的关联信息也被以网络结构记录下来. 常见的数据网络包括社交网络、电力网络、生物网络、病毒传播网络、道路交通网络等. 网络数据的记录为其统计分析提供有效的增量信息^[2]. 例如, 在社交网络中, 对用户之间互动模式的分析, 可以深入理解社会趋势和信息传播机制^[3]. 在生物信息学领域, 对网络数据的分析有助于揭示生物分子间的相互作用, 从而推动对疾病机理的理解和新药的开发^[4]. 在金融领域, 对资产之间网络关系的分析, 可以更好地理解市场风险和指导投资决策^[5]. 在物联网技术中, 网络分析对于提高系统效率和设备安全也至关重要.

面对如此快速增长和多样化的数据环境, 各领域的模型和方法也在迅速演变. 传统的统计模型, 如基于独立同分布假设的模型, 面临着新的挑战. 这些模型往往假设数据点之间是独立的, 但在网络连接数据中数据点之间的相互关联和依赖是核心特征. 因此需要新的分析方法来解读这些复杂的关系, 尤其是在面对网络数据的多维特征和动态关联时. 这种转变要求各领域不断创新, 开发出能够有效处理网络连接数据的新方法和技术. 网络连接数据的建模和有效利用已成为当前计算机科学、社会科学和统计学等多个领域的重点研究方向之一.

现有的网络连接模型的研究, 根据研究方向的不同, 可分为两类. 一类聚焦于识别关系数据中隐含的网络结构, 例如机器学习领域的社区检测^[6], 研究了如何在复杂网络中识别出具有紧密联系的节点群组. 另一类旨在发展一般的数据分析框架以利用网络数据进行回归和预测. 网络连接数据的回归问题被关注初期, 一些学者针对某些特定情况下的网络结构加以利用, 例如自回归模型^[7]及它的变体组间交互效应和组间固定效应模型^[8]. 这些模型先假定网络对数据分析的影响的具体形式, 比如内生效应、外生效应或是相关效应等, 之后将研究的重点放在效应的识别上. 融合惩罚的提出为网络连接数据的统计处理提供了新的视角. 这一思想最初由Stephanie等^[9]和Tibshirani等^[10]引入, 旨在缩减网络分析中相互关联节点之间回归系数的差异, 被视为点估计的回归版本, 或是在贝叶斯框架中被视作马尔可夫随机场的先验回归. 基于这一理论基础, Li等^[11]提出网络一致性模型, 为网络连接数据回归分析提供了一个创新而实用的框架.

针对网络连接数据的复杂性, 一个关键的任务是如何从庞大的网络数据中提取出最有意义的特征. 这不仅涉及识别数据点间的关联关系, 还要理解这些关系如何影响整体网络的行为和性质. 在这个过程中, 变量选择^[12]成为一个重要的环节, 它有助于简化模型, 增强其解释性和预测能力. 例如, 在计算机科学领域许多实际应用如自然语言处理^[13]、图像识别^[14,15]和推荐系统^[16]等都涉及大量特征, 并非所有特征都对最终任务有实质性贡献. 有效的变量选择方法能显著提高模型效率和性能, 有助于更深入地理解和利用网络结构数据. 近些年涌现了许多创新的变量选择方法, 包括了经典的前向选择、后向选择、逐步回归^[17]以及基于正则化的方法, 如Lasso、岭回归和弹性网^[18–20]等. 机器学习技术^[21,22]也为变量选择提供了广泛思路.

在现有的网络连接数据变量选择研究中, 通常采用基于p值的方法, 但面临着明显的局限性. 这主要是因为网络连接数据通常具有高度的结构依赖性和复杂性, 传统的基于p值的变量选择方法大多建立在数据点相互独立的假设之上. 在这种情况下, p值方法无法准确反映变量之间的真实关系, 尤其是在面对高维网络数据时, 这种方法很容易导致错误的发现, 从而影响模型的准确性和可靠性.

Benjamini等^[23]最早引入错误发现率(false discovery rate, FDR)的概念, 用于描述在所有选出特征中错选比例的期望. 2015年Barber等^[24]提出Knockoff仿变量方法以有效地控制多重假设检验中的错误发现率(FDR), 对于避免过度发现与伪发现相关问题具有重要意义. Knockoff仿变量方法通过构建用作控制变量的独立于原始变量的仿冒变量, 模仿原始协变量的依赖结构. 通过比较原始变量与仿冒变量的表现情况, 挑选出真正与响应变量相关的特征. 为了适应因变量对协变量的任意依赖结构, Candes等^[25]在一般高维非线性模型中引入Knockoff框架, 提出Model-X Knockoff方法. Fan等^[26]在Model-X Knockoff模型的基础上做了改进, 将Knockoff方法拓展到协变量的联合分布未知的情况, 提出Knockoff变量选择框架下的RANK方法.

为了克服相互依赖数据无法精确计算p值这一局限性, 本文提出一种基于错误发现率控制的网络连接数据变量选择方法(NLKF). NLKF方法通过模拟原始数据的结构创建原始特征的“Knockoff特征”, 并利用网络内聚性惩罚机制, 使相互连接的节点在独立效应上趋于一致, 进而挑选出与目标变量相关的重要特征. 文章从理论和实验两方面说明本文提出的NLKF方法能够将网络连接数据变量选择的错误发现率控制在预设水平. 该方法不依赖于网络结构的先验信息, 广泛适用于网络数据分析的实际应用场景. NLKF方法在处理复杂和高维数据方面表现出色, 能够更有效地满足现代数据分析的需求.

1 具有FDR控制的网络连接数据变量选择方法

在实际回归及预测问题中, 通常只有少数几个样本属性与目标变量相关. 在线性模型问题中, 意味着回归系数向量${\mathbf{\beta }}$中只有少数几个分量是非零的. 本文研究的核心问题在于设计有可靠性保证的网络连接数据的变量选择方法——网络连接数据仿变量方法(NLKF).

1.1 模型设定

考虑下面一个网络内聚性回归模型:

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

(1)

其中, 截距项${\boldsymbol{\alpha }} = {({\alpha _1}, \cdots , {\alpha _n})^{\mathrm{T}}} \in {R^n}$是未知的$n$维独立节点效应向量, ${\alpha _i}, {\text{ }}i \in [1, n]$为节点间各不相同的点效应. $ {\boldsymbol{Y}} = {({y_1}, \cdots, {y_n})^{\mathrm{T}}} \in {R^n} $是n维响应变量, $ {\boldsymbol{X}} = {({x_1}, \cdots, {x_n})^{\mathrm{T}}} \in {R^{n \times p}} $是$n \times p$维设计矩阵, ${x_i}$是第$i$个$p$维协变量, ${\boldsymbol{\beta }} = {({\beta _1}, \cdots, {\beta _p})^{\mathrm{T}}} \in {R^p}$是未知的$p$维系数向量, ${\boldsymbol{\varepsilon}} $是$n$维误差向量, 且$E({\boldsymbol{\varepsilon }}) = 0$, $Var({\boldsymbol{\varepsilon }}) = {\sigma ^2}{\boldsymbol{I}}$.

无向连接图$\mathcal{G}$记录数据间的网络连接信息. $\mathcal{G} = (V, E)$, $V = 1, 2, \cdots, n$为$\mathcal{G}$的顶点集, $E \subset V \times V$为$\mathcal{G}$的边集. 连接矩阵${\boldsymbol{A}} \in {R^{n \times n}}$记录连接图$\mathcal{G}$的节点关系信息. 当节点$u$、$v$相连, 即$(u, v) \in E$时, ${A_{uv}} = 1$; 当节点$u$、$v$不相连, 即$(u, v)\not \in E$时, ${A_{uv}} = 0$. 对于无向连接图$\mathcal{G}$, 连接矩阵满足${A_{uv}} = {A_{vu}}$. 同时假定在连接图中没有环, 即对于任意$v \in V$有 ${A_{vv}} = 0$. ${\boldsymbol{D}}$为连接图$\mathcal{G}$的度矩阵, ${\boldsymbol{D}} = {\mathrm{diag}}({d_1}, {d_2}, \cdots, {d_n})$. 其中${d_u}$为顶点$u$的度, $ {d_u} = \displaystyle\sum\limits_{v \in V} {{A_{uv}}} $. ${\boldsymbol{L}}$为连接图$\mathcal{G}$的拉普拉斯矩阵, ${\boldsymbol{L}} = {\boldsymbol{D}} - {\boldsymbol{A}}$.

上述模型中, 网络数据节点${x_i}, {\text{ }}i \in [1, n]$的节点效应${\alpha _i}, {\text{ }}i \in [1, n]$为模型增加了$n$个未知参数. 在不对${\boldsymbol{\alpha}} $的结构做任何假设的情况下, 引入网络内聚性惩罚, 认为相互连接的网络节点应具有相似的节点效应. 即引入包含网络内聚性惩罚${{\boldsymbol{\alpha }}^{\mathrm{T}}}{\boldsymbol{L\alpha }}$的损失函数:

$ Loss({\boldsymbol{\alpha , \beta}} )=\Vert {\boldsymbol{Y}}-[{\boldsymbol{X}}\tilde{{\boldsymbol{X}}}]{\boldsymbol{b}}-{\boldsymbol{\alpha}} {\Vert }^{2}+\mu {{\boldsymbol{\alpha}} }^{{\mathrm{T}}}{\boldsymbol{L}}{\boldsymbol{\alpha}} +\lambda \left|\right|{\boldsymbol{b}}|{|}_{1} $

(2)

其中, $\lambda $和$\mu $为惩罚系数, $|| \cdot |{|_1}$为${L_1}$范数, ${\boldsymbol{L}}$为数据网络结构的拉普拉斯矩阵, ${\boldsymbol{\alpha }} = {\{ {\alpha _1}, \cdots , {\alpha _n}\} ^{\mathrm{T}}}$为节点效应向量, $ {\boldsymbol{\beta }} = {\{ {\beta _1}, \cdots , {\beta _p}\} ^{\mathrm{T}}} $为模型回归系数.

网络内聚性惩罚项${{\boldsymbol{\alpha }}^{\mathrm{T}}}{\boldsymbol{L\alpha }}$可以写成更直观的形式:

$ {{\boldsymbol{\alpha }}^{\mathrm{T}}}{\boldsymbol{L\alpha }} = \sum\limits_{(u, v) \in E} {{{({\alpha _u} - {\alpha _v})}^2}} $

由上式可知, 网络内聚性惩罚思想, 是通过惩罚相互连接的节点之间的差异, 使得相互连接节点的独立效应${\boldsymbol{\alpha}} $趋于相同.

1.2 网络连接数据的变量选择方法(NLKF)

为协变量${\boldsymbol{X}}$构建模仿其内部依赖结构且与目标变量${\boldsymbol{Y}}$无关的Knockoff特征$\tilde {\boldsymbol{X}} \in {R^{n \times p}}$如下.

定义1. Network-Knockoff特征. 对于网络连接数据$ {\boldsymbol{X}} = {({x_1}, \cdots, {x_n})^{\mathrm{T}}} \in {R^{n \times p}} $, 定义其Network-Knockoff变量为:

$ \tilde {\boldsymbol{X}} = ( {{\tilde x_1}}, \cdots , {\text{ }}{\tilde x_n}) \in {R^{n \times p}} $

(3)

其中, ${\tilde x_i}, i \in [1, n]$独立产生于条件分布:

$ {\tilde x_i}|{x_i}\sim N({\boldsymbol{C}}{x_i}, {{\boldsymbol{B}}^2}) $

$ {\boldsymbol{C}} = {{\boldsymbol{I}}_p} - {\mathrm{diag}}\{ s\} {\boldsymbol{\Omega }} $

$ {\boldsymbol{B}} = {(2{\mathrm{diag}}\{ s\} - {\mathrm{diag}}\{ s\} {\boldsymbol{\Omega }}{\mathrm{diag}}\{ s\} )^{1/2}}$

其中, $ {\boldsymbol{\Omega }} $为网络连接数据$ {\boldsymbol{X}} $的精度矩阵.

由上述定义可知2p维随机向量${(x_i^{\mathrm{T}}, {({\tilde x_i})^{\mathrm{T}}})^{\mathrm{T}}}$独立同分布于均值0, 协方差如下的高斯分布:

$ \left\{\begin{array}{l} cov({x_i}) = {{\mathbf{\Sigma }}_0}\;\\ cov({x_i}, {\tilde x_i}) = {{\mathbf{\Sigma }}_0}{\mathbf{C}}\; \\ cov({\tilde x_i}) = {{\boldsymbol{B}}^2} + {\boldsymbol{C}}{{\boldsymbol{\Sigma }}_0}{{\mathbf{C}}^{\mathrm{T}}} \end{array}\right. $

则$\tilde {\boldsymbol{X}} = {\boldsymbol{X}}{{\boldsymbol{C}}^{\mathrm{T}}} + {\mathbf{Z}}{{\boldsymbol{B}}^2}, $${\boldsymbol{Z}}\mathop \sim \limits^{{\mathrm{i.i.d.}}} N(0, 1)$. s取${{\boldsymbol{\Omega }}^{ - 1}} - {2^{ - 1}}{\mathrm{diag}}\{ s\} $为正定的尽可能大的值.

实际应用中, 一种更符合实际的情况是网络连接数据的协方差未知, 需要提前估计数据的协方差${{\mathbf{\Omega }}^{ - 1}}$. 本文方法对于估计的协方差矩阵同样适用.

求解式(1)中${\boldsymbol{\beta }}$为:

$ \hat{{\boldsymbol{\beta}} }(\lambda )=\mathrm{arg}\mathrm{min}\{\Vert ({\boldsymbol{Y}}-{\boldsymbol{\alpha}} )-{[{\boldsymbol{X}}\tilde{{\boldsymbol{X}}}]}^{{\mathrm{T}}}\hat{{\boldsymbol{b}}}{\Vert }_{2}^{2}+\mu {{\boldsymbol{\alpha}} }^{{\mathrm{T}}}{\boldsymbol{L\alpha}} +\lambda \Vert \hat{{\boldsymbol{b}}}{\Vert }_{1}\} $

(4)

其中, $ \mu $和$ \lambda $为惩罚系数, $ \|\cdot\|_{2} $为$ L_{2} $范数, $ \|\cdot\|_{1} $为$ L_{1} $ 范数. $ \hat{\boldsymbol{b}}=\left[\hat{\boldsymbol{\beta}}^{{\mathrm{T}}}\; \hat{\tilde{\boldsymbol{\beta}}}^{{\mathrm{T}}}\right]^ \mathrm{T} \in R^{2 p} $, $ \hat{\boldsymbol{\beta}}^{{\mathrm{T}}} $和$ \hat{\tilde{\boldsymbol{\beta}}}^{{\mathrm{T}}} $分别对应$ \boldsymbol{X} $和$ \tilde{\boldsymbol{X}} $的回归系数. 基于图拉普拉斯矩阵$ \boldsymbol{L} $的基本性质, 此时, 式(4)可写为:

$ \begin{split} {{\hat {\boldsymbol{\beta}} }}(\lambda ) =& {\mathrm{argmin}}\{ ||({\boldsymbol{Y}} - {\boldsymbol{\alpha }}) - {[{\boldsymbol{X}}\;{{\tilde {\boldsymbol{X}}}}]^{\mathrm{T}}}\hat {\boldsymbol{b}}||_2^2 \\ &+\mu \sum\limits_{(u,v) \in E} {{{({\alpha _u} - {\alpha _v})}^2}} + \lambda ||\hat {\boldsymbol{b}}||{_1}\} \end{split} $

(5)

变量选择的可靠性由错误发现率和统计功效(Power)二者度量, 计算方法如下:

$ {\textit{FDR}}:=E\left[\frac{\#\left\{j: \beta_j=0 \text { and } {j} \in \hat{S}\right\}}{\#\{j: j \in \hat{S}\} \vee 1}\right] $

(6)

$ { Power }:=E\left[\frac{\#\left\{j: \beta_j \neq 0 \text { and } {j} \in \hat{S}\right\}}{\#\left\{j: \beta_j \neq 0\right\}}\right] $

(7)

其中, “#{j:条件}”表示满足条件的j的数量, $ \hat{S} \subset\{1, \cdots, p\} $ 为被选出特征构成的集合, $ a \vee b=\max \{a, b\} $. 由式(6)可知, FDR定义为所有被选出特征中错选特征所占比例的期望. 由式(7)可知, Power定义为所有应选特征中, 被模型正确选出的比例的期望. 本文方法的目标是充分利用网络内聚性现象, 发现并选择与响应变量真实相关的特征, 在正选比例Power足够高的前提下, 将错选比例FDR控制在目标水平.

1.3 构造NLKF统计量

构造一个服从充分性和反对称性质的一般统计量${W_j}$. 对于任意 $S \subseteq \{ 1, \cdots , p\} $, 构造统计量${W_j}$满足:

$ {W_j}({[{\boldsymbol{X}}\;{{\tilde {\boldsymbol{X}}}}]_{swap\left( S \right)}}, {\boldsymbol{y}}) = {W_j}([{\boldsymbol{X}}\;{{\tilde {\boldsymbol{X}}}}], {\boldsymbol{y}}) \cdot \left\{ {\begin{array}{*{20}{l}} { + 1},&{j \notin S } \\ { - 1},&{j \in S} \end{array}} \right.{\text{ }} $

(8)

其中, ${[{\boldsymbol{X}}\;{{\tilde {\boldsymbol{X}}}}]_{swap\left( S \right)}}$表示交换矩阵$[{\boldsymbol{X}}\;{{\tilde {\boldsymbol{X}}}}]$的${{\boldsymbol{X}}_j}, j \in S$列和$ {{{\tilde {\boldsymbol{X}}}}_j} $, $j \in S$列.

式(1)的回归系数$\hat {\boldsymbol{b}} = {[{{{\hat {\boldsymbol{\beta}} }}^{\mathrm{T}}}\;\hat{\tilde{\boldsymbol{\beta}}}^{{\mathrm{T}}}]^{\mathrm{T}}} \in {R^{2p}}$, 可构造满足上述条件的网络连接数据回归模型的${\text{NLKF}}$方法统计量为:

$ {W_j} = |{{{\hat {\boldsymbol{b}}}}_j}| - |{{{\hat {\boldsymbol{b}}}}_{j + p}}|, {\text{ }}j \in \{ 1, \cdots , p\} $

(9)

1.4 计算统计量的阈值

将上述统计量$ {W_j}, {\text{ }}j \in [1, p] $按照$|{W_{{j}}}|$从大到小依次排列, 寻找能够识别相关特征的特征统计量的最小值. 预先指定错误发现率$q$, 定义NLKF方法的阈值$T$如下:

$ T = \min\left\{ t \in \mathcal{W}:\frac{{\# \{ j:{W_j} \leqslant - t\} }}{{\# \{ j:{W_j} \geqslant t\} \vee 1}} \leqslant q\right\} $

(10)

其中, $\mathcal{W} = \{ |{W_j}|:j = 1, \cdots , p\} $. 因统计量$ {W_j} $具有反对称性, 式(10)将$T$定义为使错误发现率控制在指定值$q$以内的最小统计量的值. 当上述集合是空集时, 取$ T = + \infty $.

为使FDR被更好地控制, 定义一个修正的阈值$ {T_ + } $:

$ T_ + = \min\left\{ t \in \mathcal{W}:\frac{{1 + \# \{ j:{W_j} \leqslant - t\} }}{{\# \{ j:{W_j} \geqslant t\} \vee 1}} \leqslant q\right\} $

(11)

当上述集合是空集时, 取$ {T_ + } = + \infty $. 将使用阈值$ {T_ + } $进行特征选择的方法为${\text{NLK}}{{\text{F}}_{\text{+}}}$方法.

由式(10)和式(11)知, ${\text{NLK}}{{\text{F}}_{\text{+}}}$方法比${\text{NLKF}}$方法选出的集合更保守, 因为${\text{NLK}}{{\text{F}}_{\text{+}}}$方法的边界值${T_ + }$总是高于${\text{NLKF}}$方法的边界值$T$.

${\text{NLKF}}$方法可选出特征集合:

$ \hat S = \{ j:{W_j} \geqslant T\} $

${{{\mathrm{NLKF}}} _ + }$方法可选出特征集合:

$ {\hat S_ + } = \{ j:{W_j} \geqslant {T_ + }\} $

设置目标FDR水平为$q$, ${\text{NLKF}}$ (${\text{NLK}}{{\text{F}}_{\text{+}}}$)方法通过模仿原始特征的内部依赖构造与目标变量无关的Network-Knockoff特征. 使用式(2)中提出的网络内聚性惩罚函数, 求解网络内聚性回归模型式(1), 得到模型参数的估计. 借助参数估计值依次计算NLKF统计量$W$及其对应的阈值$T$ (${T_ + }$), 进行网络连接数据的变量选择. 下文从理论和实验两个方面, 证明${\text{NLKF}}$方法能够将网络连接数据变量选择方法的FDR控制在预设水平$q$.

2 理论性质

本节从理论上证明NLKF方法对变量选择FDR的控制. 这一证明过程说明了本文方法的统计可行性, 为网络数据特征筛选的可信度提供保障.

2.1 错误发现率控制在目标水平

首先, 介绍两个关键引理.

引理1. 对于任意集合 $S \subset \{ 1, \cdots , p\} $:

$ [{\boldsymbol{X}}\;\tilde {\boldsymbol{X}}]_{swap(S)}^{\mathrm{T}}{\boldsymbol{A}}{[{\boldsymbol{X}}\;\tilde {\boldsymbol{X}}]_{swap(S)}} = {[{\boldsymbol{X}}\;\tilde {\boldsymbol{X}}]^{\mathrm{T}}}{\boldsymbol{A}}[{\boldsymbol{X}}\;\tilde {\boldsymbol{X}}] $

即为, 对于任意$j \in S$交换${X_j}$和$\tilde {\boldsymbol{X}}$, Gram 矩阵 ${\boldsymbol{G}} = {[{\boldsymbol{X}}\;\tilde {\boldsymbol{X}}]^{\mathrm{T}}} {\boldsymbol{A}}[{\boldsymbol{X}}\;\tilde {\boldsymbol{X}}]$不变. 其中${\boldsymbol{A}} = {{\boldsymbol{I}}_n} - {({{\boldsymbol{I}}_n} + \lambda {\boldsymbol{L}})^{ - 1}}$, L为网络连接数据的拉普拉斯矩阵.

引理2. 记$ {W_j} = |{\hat b_j}| - |{\hat b_{j + p}}|, {\text{ }}j \in \{ 1, \cdots , p\} $为式(9)中定义的${\text{NLKF}}$统计量, 在已知$(|{W_1}|, \cdots , |{W_p}|)$时, 与目标变量无关的变量对应的NLKF统计量${W_j}{\text{, }}j \in S_0^c$的正负性服从标准二项分布. 其中$\;S_0^c$为无关变量下标集.

引理2的证明见附录A.

基于上述引理, 可给出本文关键的错误发现率控制定理及其证明如下.

定理1. FDR控制定理. 对于统计量$ {W_j}, j \in [1, p] $, 给定$q$为FDR的目标水平. 网络连接数据的变量选择方法的阈值$ T = \min\left\{ t \in \mathcal{W}:\dfrac{{\# \{ j:{W_j} \leqslant - t\} }}{{\# \{ j:{W_j} \geqslant t\} \vee 1}} \leqslant q\right\} , {\text{ NLKF}} $挑选出的特征集$\hat S = \{ j:{W_j} \geqslant T\} , $可以将FDR控制在目标水平, 即:

$ E\left[\frac{{\# \{ j:{\beta _j} = 0\;{\text{and}}\;{{j}} \in \hat S\} }}{{\# \{ j:j \in \hat S\} + 1/q}}\right] < q $

其中, $\mathcal{W} = \{ |{W_j}|:j = 1, \cdots , p\} $. 对于一个保守的阈值, $ T _{+} = \min \left\{ t \in \mathcal{W}:\dfrac{{1 + \# \{ j:{W_j} \leqslant - t\} }}{{\# \{ j:{W_j} \geqslant t\} \vee 1}} \leqslant q\right\} ,\; {{ {\mathrm{NLKF}}_{+}}} $ 挑选出的特征集${\hat S_ + } = \{ j:{W_j} \geqslant {T_ + }\} $, 可以将FDR控制在目标水平, 即:

$ E\left[\frac{{\# \{ j:{\beta _j} = 0\;{\text{and}}\;{{j}} \in \hat S\} }}{{\# \{ j:j \in \hat S\} \vee 1}}\right] < q $

其中, $\mathcal{W} = \{ |{W_j}|:j = 1, \cdots , p\} $.

证明: 由引理2可知, 统计量${W_j}$的符号是独立同分布的. 计算统计量${W_j}$之后, 必须选择一个数据依赖的边界值: $ T = \min \left\{ t \in \mathcal{W}:\dfrac{{\# \{ j:{W_j} \leqslant - t\} }}{{\# \{ j:{W_j} \geqslant t\} \vee 1}} \leqslant q\right\} , $不失一般性地假定 $|{W_1}| \geqslant |{W_2}|, \cdots , |{W_{p-1}}|\geqslant |{W_p}|$, 寻找一个使得$FDR \leqslant q$的最小边界值$T$. 从统计量的最小值$t = |{W_p}|$开始测试, 然后逐步测试$t = |{W_{p - 1}}|$, 当测试到一个$t$值满足$FDR(t) \leqslant q$时停止. 通过可选停时定理, 上鞅在随机时间$t = T$的期望值受其在时间$t = 0$的期望值限制: 令${p_0}$为空特征的数量并写出 $Y = \# \{ {\beta _j} = 0\;{\mathrm{and}}\;{W_j} \leqslant 0\} $, 其中“#”表示计数, 有:

$ \begin{gathered} E\left[\frac{{\# \{ j:{\beta _j} = 0\;{\mathrm{and}}\;{W_j} \leqslant - T\} }}{{1 + \# \{ j:{\beta _j} = 0\;and\;{W_j} \geqslant T\} }}\right] \\ \leqslant E\left[\frac{{\# \{ j:{\beta _j} = 0\;{\mathrm{and}}\;{W_j} \leqslant 0\} }}{{1 + \# \{ j:{\beta _j} = 0\;{\mathrm{and}}\;{W_j} \geqslant 0\} }}\right] \\ = E\left[\frac{Y}{{1 + {p_0} - Y}}\right] \leqslant 1 \\ \end{gathered} $

由于$sign({W_j})\mathop \sim\limits^{{\mathrm{i.i.d.}}} \{ \pm 1\} $中$j$为空特征, 因此 $Y = \# \{ j: {\text{ }}{\beta _j} = 0{\text{ }}{\mathrm{and}}{\text{ }}{W_j} \leqslant 0\} $为服从二项分布$({p_0}, \dfrac{1}{2})$的随机变量.

$ \begin{gathered} {\textit{FDR}} = \mathbb{E}\left[\frac{{\# \{ j:{\beta _j} = 0\;{\mathrm{and}}\;{W_j} \geqslant T\} }}{{\# \{ j:{W_j} \geqslant T\} \vee 1}}\right] \\ = E\left[\frac{{\# \{ j:{\beta _j} = 0\;{\mathrm{and}}\;{W_j} \geqslant T\} }}{{1 + \# j:{\beta _j} = 0\;{\mathrm{and}}\;{W_j} \leqslant - T}} \cdot \frac{{1 + \# \{ j:{W_j} \leqslant T\} }}{{\# \{ j:{W_j} \geqslant T\} \vee 1}}\right] \\ \leqslant E\left[\frac{{\# \{ j:{\beta _j} = 0\;{\mathrm{and}}\;{W_j} \geqslant T\} }}{{1 + \# j:{\beta _j} = 0\;{\mathrm{and}}\;{W_j} \leqslant - T}} \cdot q\right] \leqslant q \\ \end{gathered} $

定理1说明了本文方法能确保所选特征集合的FDR维持在预定的目标之内, 为变量选择的可靠性和网络数据分析的可重复性提供了保证. 由第1节可知, ${\text{NLKF}}$ (${\text{NLK}}{{\text{F}}_{\text{+}}}$)方法中统计量的计算涉及原始变量回归系数与Network-Knockoff变量的回归系数共计$2 \times p$个参数, 网络连接数据各节点独立效应的引入又为模型新增了n个未知参数. 为说明模型式(1)中解的存在性, 引入如下定理.

定理2. 解的存在性. 对于网络内聚性回归模型, 见式(1), 网络一致性惩罚的回归系数可求解为:

$ \hat {\boldsymbol{b}} = {({\hat {\boldsymbol{\beta}} ^{\mathrm{T}}}, {\hat {\tilde {\boldsymbol{\beta}} }^{\mathrm{T}}})^{\mathrm{T}}} = {[{\boldsymbol{X}}, \tilde {\boldsymbol{X}}]^{\mathrm{T}}}{\boldsymbol{A}}{[{\boldsymbol{X}}, \tilde {\boldsymbol{X}}]^{ - 1}}{[{\boldsymbol{X}}, \tilde {\boldsymbol{X}}]^{\mathrm{T}}}{\boldsymbol{AY}} $

其中, ${\boldsymbol{A}} = {{\boldsymbol{I}}_n} - {({{\boldsymbol{I}}_n} + \lambda {\boldsymbol{L}})^{ - 1}}$, 是一个对称正定阵.

证明: 求解模型式(1), 最小化下方损失函数: $ Loss({\boldsymbol{\alpha}} , {\boldsymbol{\beta}} )=\Vert {\boldsymbol{Y}}-[{\boldsymbol{X}}\tilde{{\boldsymbol{X}}}]{\boldsymbol{\beta}} -{\boldsymbol{\alpha}} {\Vert }^{2}+{\displaystyle \sum _{(u, v)\in E}{({\alpha }_{u}-{\alpha }_{v})}^{2}}+ \lambda \left|\right|{\boldsymbol{\beta}} |{|}_{1} $, 可得: $ \hat {\boldsymbol{\theta }} = (\hat {\boldsymbol{\alpha }}{{, }}\hat {\boldsymbol{\beta }}) = {({{\boldsymbol{Z}}^{{{\mathrm{T}}}}}{\boldsymbol{Z}} + \lambda {\boldsymbol{M}})^{ - 1}}{{\boldsymbol{Z}}^{{{\mathrm{T}}}}}{\boldsymbol{Y}} $, 其中, $ {\boldsymbol{Z}} = ({{\boldsymbol{I}}_n}, {\boldsymbol{X}}, \tilde {\boldsymbol{X}}) $,$ {\boldsymbol{M}} = \left[ {\begin{array}{*{20}{c}} {\boldsymbol{L}}&{{{\boldsymbol{0}}_{n \times 2p}}} \\ {{{\boldsymbol{0}}_{2p \times n}}}&{{{\boldsymbol{0}}_{2p \times 2p}}} \end{array}} \right] $.

为保证解一定存在, 使用$ {\boldsymbol{L}} + \gamma {\boldsymbol{I}} $代替$ {\boldsymbol{L}} $, $ \gamma $是一个很小的正的常数.

3 算法

在实际应用中, 我们经常遇到具有内部依赖性的网络数据, 这类数据的复杂性使得传统的变量选择方法难以有效处理. 为应对这一挑战, 本文第1节提出了一种新的方法, 即网络连接数据变量选择方法(NLKF), 下文将详细阐述这一方法的具体实现流程.

为网络连接数据变量选择的${\text{NLKF}}$ (${\text{NLK}}{{\text{F}}_{\text{+}}}$)方法设计算法如算法1.

算法1. NLKF算法

输入: $\scriptstyle ({\boldsymbol{X}}, {\boldsymbol{y}}), {\text{ }}{\boldsymbol{L}}, {\text{ }}\lambda , {\text{ }}{q_{FDR}} $.

输出: $\scriptstyle \hat S \subset \{ 1, \cdots, p\} $.

1) 给定协方差矩阵$ \scriptstyle {\boldsymbol{\Sigma }} $求精度矩阵$ \scriptstyle {\boldsymbol{\Omega }} = {{\boldsymbol{\Sigma }}^{ - 1}} $. 当协方差矩阵未知时, 用GLasso或iSEE估计协方差阵, 求$ \scriptstyle \hat {\boldsymbol{\Omega }} $.

2) 根据定义1, 构造Network-Knockoff特征$\scriptstyle \tilde {\boldsymbol{X}}$:

$\left\{\begin{array}{l} \scriptstyle {\boldsymbol{C}} = {{\boldsymbol{I}}_p} - {\mathrm{diag}}\{ s\} {\boldsymbol{\hat \Omega }}, \;{\boldsymbol{B}} = {(2{\mathrm{diag}}\{ s\} - {\mathrm{diag}}\{ s\} {\boldsymbol{\hat \Omega }}{\mathrm{diag}}\{ s\} )^{1/2}}\\\scriptstyle {\tilde x_i}|{x_i}\sim N({\boldsymbol{C}}{x_i}, {{\boldsymbol{B}}^2}),\; \tilde {\boldsymbol{X}} = (\tilde x_1, \cdots , {\text{ }}{\tilde x_n}) \in {R^{n \times p}}\end{array}\right. $

3) 引入网络内聚性惩罚函数, 即式(2):

$\scriptstyle Loss({\boldsymbol{\alpha}} , {\boldsymbol{\beta}} )=\Vert {\boldsymbol{Y}}-[{\boldsymbol{X}}\tilde{{\boldsymbol{X}}}]{\boldsymbol{b}}-{\boldsymbol{\alpha}} {\Vert }^{2}+\mu {{\boldsymbol{\alpha}} }^{{\mathrm{T}}}{\boldsymbol{L}}{\boldsymbol{\alpha}} +\lambda \left|\right|{\boldsymbol{b}}|{|}_{1} $

求解模型$\scriptstyle {\boldsymbol{Y}} = {\boldsymbol{\alpha }} + {\boldsymbol{X\beta }} + \tilde {\boldsymbol{X}}\tilde {\boldsymbol{\beta }} + {\boldsymbol{\varepsilon }} $得:

$\scriptstyle {\boldsymbol{\hat b}} = {({{\boldsymbol{\hat \beta }}^{\mathrm{T}}}, {\hat {\tilde {\boldsymbol{\beta }}}^{\mathrm{T}}})^{\mathrm{T}}} = {[{\boldsymbol{X}}, \tilde {\boldsymbol{X}}]^{\mathrm{T}}}{\boldsymbol{A}}{[{\boldsymbol{X}}, \tilde {\boldsymbol{X}}]^{ - 1}}{[{\boldsymbol{X}}, \tilde {\boldsymbol{X}}]^{\mathrm{T}}}{\boldsymbol{AY}} $

其中, $\scriptstyle {\boldsymbol{A}} = {{\boldsymbol{I}}_n} - {({{\boldsymbol{I}}_n} + \lambda {\boldsymbol{L}})^{ - 1}}$.

4) for j=1 to p do:

$\scriptstyle {W_j} = |{{\boldsymbol{\hat b}}_j}| - |{{\boldsymbol{\hat b}}_{j + p}}| $

end

5) 求阈值$\scriptstyle T$和$\scriptstyle {T_ + }$

6) return $\scriptstyle \hat S = \{ j:{W_j} \geqslant T\} $和$\scriptstyle {\hat S_ + } = \{ j:{W_j} \geqslant {T_ + }\} .$

当原始变量协方差结构已知时, 使用定义1构造原始特征的Network-Knockoff特征$\tilde {\boldsymbol{X}}$. 原始变量协方差未知时, 先估计原始变量的协方差, 再构造$\tilde {\boldsymbol{X}}$. 引入原始变量的Network-Knockoff特征后, 模型(1)可更直观地表达为$ {\boldsymbol{Y}} = {\boldsymbol{\alpha }} + {\boldsymbol{X\beta }} + \tilde {\boldsymbol{X}}\tilde {\boldsymbol{\beta }} + {\boldsymbol{\varepsilon }} $, 其中$ {\boldsymbol{\alpha}} , \text{ }{\boldsymbol{\beta}} 和\text{ }\tilde{{\boldsymbol{\beta}} } $分别对应$ {{n, }}\;p, {\text{ }}p $个未知参数. 定理2说明了模型式(1)的解的存在性, 计算方法如算法1中的步骤3). 基于模型求解结果, 可依次计算特征统计量$W$及阈值$T$ (${T_ + }$) , 进而挑选出目标变量的相关特征.

4 数值模拟

本节主要考察使用网络连接数据仿变量方法(${\text{NLKF}}$)在模拟数据上的方法表现. 为了更直观地评估模型效果, 将本文方法与一般变量选择法(network cohesion Lasso, NCL)进行比较实验.

模拟数据的网络结构由随机图模型生成, 共有$n$个节点, $p$维属性, $K$个随机生成块. 各节点样本由参数为$({\pi _1}, \cdots , {\pi _K})$的多元正态分布生成, 并独立分配给$K$个图块. 各个图块的标签为${c_i}$, $i = 1, \cdots , n$. 图形的边${E_{ij}}$, $1 \leqslant i \leqslant j \leqslant n$, 为$P({E_{ij}} = 1) = {B_{{c_i}{c_j}}}$的Bernoulli随机变量. 设置$K = 2$, ${\pi _1} = {\pi _2} = 1/2$, 同一图块内部节点相互连接的概率$ {B_{kk}} = 0.8 $, 不同图块之间节点相互连接的概率${B_{kl}} = 0.2$.

各个节点${X_i}$独立取样自正态分布$N(0, {\Sigma _X})$, 构造一个系数向量${\mathbf{\beta }} \in {\mathbb{R}^p}$, 其中随机选出15个的位置放置均值为$ {\beta _{{\text{signal}}}} $的15个信号值, 其他$p - k$项的值为0. 独立节点效应${\alpha _i}$独立取样自均值与该节点所在图块相关的正态分布$N({\eta _{{c_i}}}, {s^2})$. 其中, ${\eta _{{c_i}}}$表示该节点所在图块${c_i}$对应的均值, 取${\eta _{{c_1}}} = - 1,\; {\eta _{{c_2}}} = 1$. $s$表示每个图块中个节点的内聚性强度, 取$s = 1$.

设置目标错误发现率(FDR)水平为0.2. 两类方法NLKF和NCL分别借助R程序包“knockoff”和 “glmnet”实现. 用10折交叉验证的方法调优参数.

4.1 数值模拟1: 信号强度比较实验

为观察信号强度的大小对实验结果的影响, 设置$ {\beta _{{\text{signal}}}} $的均值从0.5逐渐增加至5, 分别在高维和低维进行比较试验, 实验结果如图1所示.

实验结果可知, NLKF方法明显优于经典的Model-X Knockoff和Lasso变量选择方法. 低维弱信号强度情况下, NLKF方法的Power比Model-X Knockoff提高0.1, 且随着信号强度的增强, 二者的Power有明显提升, 且相比于Lasso方法, NLKF方法能够将FDR控制在目标水平. 高维弱信号强度情况下, NLKF方法表现与经典的Model-X Knockoff方法差别不大, 随着信号强度的增强, NLKF方法的Power明显优于Model-X Knockoff和Lasso变量选择方法.

4.2 数值模拟2: ${\text{NLKF}}$、${\text{NLK}}{{\text{F}}_{\text{ + }}}$比较实验

设置信号强度均值为3.5, ${\Sigma _{X, ij}} = {0.5^{|i - j|}}$, 数值模拟结果见表1. 其中, $N_\epsilon(0,1)$表示误差项$\epsilon$服从均值为0, 方差为1的正态分布. 表1显示, 使用NLKF方法可以将FDR控制在目标水平(20%)以下, 并且Power值大于NCL方法. 相比于${\text{NLKF}}$方法, 改进的${\text{NLK}}{{\text{F}}_{\text{+}}}$对FDR的控制效果更优.

4.3 数值模拟3: 泛化能力测试

设置信号强度均值3.5, ${\Sigma _{X, ij}} = 0.5 + I(i = j)\times 0.5$, 数值模拟结果见表2. 表2显示, 在原始数据协方差为不同结构时, NLKF方法仍能够在Power足够高的情况下, 将FDR控制在目标水平. 说明了本文方法具有一定的泛化能力.

3个数值模拟分别从不同角度验证本文方法在数据计算上的可行性、实用性和可靠性. 从模拟结果来看, 本文提出的${\text{NLKF}}$ (${\text{NLK}}{{\text{F}}_{\text{+}}}$)方法能够在较高的统计功效(Power)下, 将错误发现率(FDR)控制在目标水平, 显著优于不使用仿变量的传统Lasso变量选择方法(NCL).

5 实际数据分析

使用沪深两市股票的因子数据及股票所属申万一级行业构建的网络关系做实证分析. 因子是指在股票市场中, 能够解释股票收益率的一些特定因素. 常见的因子包括市场风险溢酬因子、市值因子、账面市值比因子等. 使用因子来构建投资组合的方法有很多种, 其中最常见的方法是基于因子收益率的方法. 具体来说, 首先根据投资需求制定投资目标. 然后从因子池中选择合适因子、计算单因子对个股收益率贡献、进而基于各个股票的因子表现来构建股票投资组合. 最后根据自己的需求和风险偏好, 适时调整每个因子在投资组合中的权重以适应市场变化. 其中, 能否从因子库中选出对股票收益率具有更好解释性的因子直接决定了后续投资的表现.

本文使用2022年1–12月沪深两市的4 000只A股($n = 4000$)的200个因子($p = 200$)作为实证数据集. 旨在从中选出对股票收益率有影响的重要因子. 计算过程中的协方差矩阵通过${{\mathrm{GLasso}}} $的方法估计, 参数$\lambda $用5折交叉验证调节. 因子频率有日度、周度和月度, 此处选取月底日度数值, 并将空缺值用当月日均值填补. 根据申万二级行业构建股票标的之间的网络关系. 设置目标错误发现率水平$q = 0.2$. 本文参数$\lambda $和$\mu $使用交叉验证使均方误差(MSE)最小的方法调整. 最终特征选择结果记录见表3.

对相同数据集使用Lasso变量选择方法可选出54个对当前收益率具有较强解释性的因子, 详情见附录B.

对比Lasso变量选择方法选出的54个特征, ${\text{NLKF}}$ (${\text{NLK}}{{\text{F}}_{\text{+}}}$)方法选出的特征数量有明显降低. 参考第4节数值模拟可知, Lasso变量选择方法具有较高的统计功效, 即Lasso方法能够有效地挑出与目标变量相关的特征. 但Lasso方法在选出相关特征的同时也选出较多无关特征, 表现为较高的错误发现率. 实证结果表明, 与传统的变量选择方法Lasso相比, 本文提出的${\text{NLKF}}$和${\text{NLK}}{{\text{F}}_{\text{+}}}$方法在高统计功效的前提下, 错误发现率明显降低.

2022年12月, 食品饮料、酒店航空等消费行业呈现了较为明显的复苏态势. 从分析结果中可以看到, 短期反转、季节反转因子、波动性因子等都是影响A股股票收益率的关键因素. 并且, 公司本身的经营状况和资产质地也是影响公司价值的重要指标, 因此流动资产比总资产、资产周转率等指标也被视为重要变量被挑选出来.

6 结论与展望

2023年2月, 《数字中国建设整体布局规划》明确提出, 到2025年基本形成横向打通、纵向贯通、协调有力的一体化格局, 数字中国建设取得重要进展. 到2035年, 数字化发展水平进入世界前列, 数字中国建设取得重大成就, 同时要把中国数字化建设方面的技术和资源大量输出到国外, 为推动全球数字化经济发展做出重要贡献. “数据”作为新时代的“关键能源”, 将成为新时代的重要生产力和发展引擎.

本文将网络连接数据的融合惩罚和Knockoff思想结合, 提出了针对网络结构数据变量选择的${\text{NLKF}}$和${\text{NLKF}}_{\text{+}}$方法. 从理论水平证明了这个方法对FDR的控制能力, 并从模拟和实证两个方面验证了模型的有效性. NLKF (${\text{NLKF}}{_{\text{+}}}$)方法适用于协方差矩阵未知的场景, 在处理复杂、高维且具有关联结构的数据方面表现出色, 能够满足当前网络连接数据分析的需求. 从FDR角度出发的数据分析, 在高维的、相互关联的数据建模上也将有更多、更广泛的应用.