大家好，又见面了，我是你们的朋友全栈君。

论文地址

A Multi-Task Learning Formulation for Predicting Disease Progression
Author: Zhou J, Yuan L, Liu J, et al.

摘要

临床判断阿尔茨海默症(AD)的两个指标：Mini Mental State Examination(MMSE)和Alzheimer’s Disease Assessment Scale cognitive subscale(ADAS-Cog)。本文基于岭回归提出了一种多任务学习的算法来通过以上两个标准(criteria)预测疾病进展，同时也可以选出有代表性的特征。本名提出一个多任务学习模型。本文认为，对于AD，每一个时间点的预测算作一个独立任务，通过给岭回归加上时间项约束(temporal group LASSO regularizer)，将独立任务整合为多任务，从而捕获(capture)不同任务中的内在联系(intrinsic relatedness)，前提是认为每个时间点的回归任务近似。时间项约束包括两部分，一个是${\ell }_{2,1}$范数惩罚(penalty)，目的是保证对于所有时间点只挑选一小部分特征。另一个是时间平滑项，目的是为了在连续(successive)的时间中，前后的两个模型的偏差(deviation)程度要小。大量(extensive)实验证明本文提出的模型相较其他的模型更具有优势，同时挑选的特征也同之前他人研究中的结果一致。作者使用的数据集是ADNI。

1. 介绍

目前AD的确诊需要对于脑部进行活检或者解剖，MMSE和ADAS-Cog可以作为评估病人意识状态和诊断潜在AD病人的两个指标。MMSE与潜在(underlying)AD病理学(pathology)特征以及脑功能进一步恶化(progressive deterioration)相关，ADAS-Cog是进行AD药物试验(drug trial)的金标准。作者此时提出两个问题：1）怎样通过意识评分(cognitive score)如MMSE和ADAS-Cog来预测疾病的发展？2）在疾病发展过程中哪些特征是最具有代表性的特征，这些特征的最小集合是什么？这些可以追溯疾病发展的特征包含在脑部MRI成像，脑脊液(CSF)和一些临床诊断评估(baseline clinical assessments)中。
患者的年龄，性别，受教育程度等潜在风险因子与意识评分之间的关系已经被研究过。通过将一些特征(如基于MRI的脑灰质体积、密度，脑血管的形状，海马区)与MMSE关联起来进行研究也已经存在。并且一些研究表明内侧额叶(medial temporal lobe)的密度和大小以及其他风险因素与6个月AD病人的MMSE有关，所以可以使用这些特征来预测未来一段时间的意识水平。Ashford与Schmitt通过使用”time-index”提出时间函数(horologic function)来测定痴呆(dementia)的速度。基于脑萎缩空间模式(spatial patterns of brain atrophy)提出的SPARE-AD的方法也被提出。通过能量函数进行相关研究的也存在。在这些研究中通常都存在“维数灾难(curse of dimensionality)”的问题。所以一些降维的方法也被使用，PCA降维的显著缺点是所有特诊都被聚合，所以降维后的模型不可解释。相关向量回归(relevance vector regression, RVR)方法将特征选择整合到了训练过程中。但是这些方法都只是在单一时间点对于意识评分进行预测，在AD预测中并没有太大用。
我们提出的多任务学习的方法中，每一个任务都只关心单一时间点的预测。多任务学习旨在通过同时学习多重相关任务，提高模型的泛化能力。多任务学习的关键点是发现任务之间的内在联系。对于本文中的疾病，认为一小部分特征具有代表性是合理的。并且，两个连续时间点所对应的意识评分的差距应该小，所以有必要进行时间平滑。所以L2,1范数惩罚用来选一小部分特征，时间平滑项用来减小连续时刻对应的意识评分之间的差距。本文使用ADNI数据集中的三项：MRI影像，CSF和临床诊断评估分数。

2. 多任务回归方程

在纵向(时间向)AD研究中，本文从多个时间点获取患者意识评分。多任务学习中，在不同任务间的时间平滑信息可以被并入模型中作为先验信息(prior knowledge)。假设在多任务回归问题中，有t个时间点，n个训练样本，每个样本有d个特征。一个简单的线性多任务学习模型可以表示如下：
$$mi{n}_{(W)}||XW-Y|{|}_F^2+{\theta }_1||W|{|}_F^2$$
其中X的维度为{n × d} × t，Y的维度为{n × 1} × t， W的维度为d × t。第一项用来衡量训练集上的经验损失(empirical error)，第二项用来控制泛化误差。${\theta }_1>0$为正则化参数。$||.|{|}_F$是Frobenius范数。上式通常被被称为岭回归(ridge regression)并且具有解析解：
$$W=(X^{T}X+{\theta }_{1}I{)}^{-1}{X}^{T}Y$$
其中I是d × d的单位矩阵。上式的缺点为它将不同时间点认为是相互独立的，所以对于本文所考虑的问题其回归效果并不好(没有考虑时间项约束)。

2.1 时间平滑先验

带有时间平滑项的线性回归模型：
$$mi{n}_{(W)}||XW-Y|{|}_F^2+{\theta }_1||W|{|}_F^2+{\theta }_2\sum _{i=1}^{t-1}||w^i-w^{i+1}|{|}_F^2$$
${\theta }_2\ge 0$是用来控制时间平滑的正则化参数。时间平滑项又可以表示为：
$$\sum _{i=1}^{t-1}||w^i-w^{i-1}|{|}_F^2=||WH|{|}_F^2$$
上式中，$H$的维度为t × t-1，定义如下:
$$H_{ij}=\{\begin{array}{rcl}1& & i=j\\ -1& & i=j+1\\ 0& & otherwise\end{array}$$
所以带有时间平滑项的线性回归模型变为：
$$mi{n}_{(W)}||XW-Y|{|}_F^2+{\theta }_1||W|{|}_F^2+{\theta }_2||WH|{|}_F^2$$
上式中也存在解析解，首先上式对于W求偏导并令它等于0：
$$\begin{array}{rl}{X}^{T}XW-X^T+{\theta }_{1}W+{\theta }_{2}WH{H}^T& =0,\\ (X^{T}X+{\theta }_1{I}_d)W+W({\theta }_{2}H{H}^T)& =X^{T}Y,\end{array}$$
其中$I_d$的维度为d × d。因为$(X^T+{\theta }_1{I}_d)$和${\theta }_{2}H{H}^T$是对称矩阵，所以它们的特征分解为$Q_1{\Lambda }_1{Q}_1^T$和$Q_2{\Lambda }_2{Q}_2^T$，其中${\Lambda }_1=diag({\lambda }_1^{(1)},{\lambda }_1^{(2)},...,{\lambda }_1^{(d)})$和${\Lambda }_2=diag({\lambda }_2^{(1)},{\lambda }_2^{(2)},...,{\lambda }_2^{(d)})$是它们的特征值，并且$Q_1$和$Q_2$是正交的。带入以上两个方程，得到：
$$\begin{array}{rl}{Q}_1{\Lambda }_1{Q}_1^{T}W+W{Q}_2{\Lambda }_2{Q}_2^T& =X^{T}Y,\\ {\Lambda }_1{Q}_1^{T}W{Q}_2+Q_1^{T}W{Q}_2{\Lambda }_2& =Q_1^T{X}^{T}Y{Q}_2,\end{array}$$
记$\hat{W}=Q_1^{T}W{Q}_2$；$D=Q_1^T{X}^{T}Y{Q}_2$。则上式变为${\Lambda }_1\hat{W}+\hat{W}{\Lambda }_2=D$，所以$\hat{W}$为：
$${\hat{W}}_{ij}=\frac{D_{i,j}}{{\lambda }_1^{(i)}+{\lambda }_2^{(j)}}$$
所以最优权重矩阵为$W^{\ast }=Q_1\hat{W}{Q}_2^T$。

2.2 解决数据不完整问题

本文尝试将不完整的数据包含进训练数据集中，不完整训练数据指的是部分患者的意识评分在不同时间点可能不完整，此时，带有时间平滑约束的模型不存在解析解。解决这个问题，首先定义一个矩阵$S_{i,j}$，维度为n × t，$i$为样本序号，$j$为时间点，具体形式为：
$$S_{i,j}=\{\begin{array}{rcl}1& & valueexists\\ 0& & otherwise\end{array}$$
则改进后的公式为：
$$mi{n}_{(W)}||S\odot (XW-Y)|{|}_F^2+{\theta }_1||W|{|}_F^2+{\theta }_2||WH|{|}_F^2$$
其中$\odot$为矩阵点乘。定义$P_r(.)$为矩阵行数选择。$P_r(A)$只包括$r_i̸=0$的$A_i$，其中$A_i$是$A$的第$i$行。再令$S^i$为$S$的第$i$列。则可以得到：
$$\begin{array}{rl}{X}_{(i)}& =P_{S^i}(X),\\ y_{(i)}& =P_{S^i}(X^i)\end{array}$$
将改进后的公式对于$w^i$求导并置为0，得到：
$$A{w}^{i-1}+M_i{w}^i+A{w}^{i+1}=T_i$$
其中：
$$\begin{array}{rl}A& =-{\theta }_2{I}_d,\\ M_i& =X_{(i)}^T{X}_{(i)}+{\theta }_1{I}_d+2{\theta }_2{I}_d,\\ T_i& =X_{(i)}^T{y}_{(i)}\end{array}$$

2.3 LASSO时间群正则

ANDI数据集中包含的数据维度太大，会导致维度灾难。传统降维方法不适合，如PCA会导致模型的不可解释性。传统的特征挑选方法也不适合存在缺失数据的多任务回归。本文提出基于${\ell }_{2,1}$范数的惩罚项用于特征选择，目的是希望一小部分具有代表性的特征可以被选择出来。Lasso群正则可以保证所有的回归模型可以共享相同的一部分特征。在之前的公式中加入群正则后，得到：
$$mi{n}_{(W)}||S\odot (XW-Y)|{|}_F^2+{\theta }_1||W|{|}_F^2+{\theta }_2||WH|{|}_F^2+\delta ||W|{|}_{2,1}$$
此时$||W|{|}_{2,1}={\sum }_{i=1}^d\sqrt{{\sum }_{j=1}^t{W}_{i,j}^2}$。注意这边的2范数并不是矩阵范数(L2谱范数)，因为对于$||W|{|}_{2,1}$是逐行来求，所以这边依旧是向量范数。$\delta$是正则化参数。并且特征先被L2约束，再被L1约束。由于上式同时包含平滑约束和非平滑约束，所以可以用梯度下降进行优化求解。

2.3.1 纵向(时间项)稳定特征选择

首先一个需要解决一个问题，选择一定数量的特征用于模型训练。交叉验证可以用于特征选择，但是通常它会需要超过需求数目的特征。本文使用了稳定性选择的方法(Stability Selection)。实验发现所选出的特征对于${\theta }_1$和${\theta }_2$不敏感，所以主要精力放在可以控制模型稀疏性的参数$\delta$上。
具体做法如下：首先记$K$为特征的序号，$\Delta$为正则化参数$\delta$的集合，$\gamma$为迭代次数。首先无放回随机挑选占总体样本数目一半的子样本$B_{(j)}$， B ( j ) = { B ( j ) X , B ( j ) Y } B_{(j)}=\{B_{(j)}^X, B_{(j)}^Y \} B(j)​={
B(j)X​,B(j)Y​}，对于任意给定的$\delta >0$，记${\tilde{W}}^{(j)}$为模型在$B_{(j)}$上的最优解。记 U δ ( B ( j ) ) = { k : W ~ ( j ) ̸ = 0 } U^{\delta}(B_{(j)})=\{k:\tilde{W}^{(j)} \not =0 \} Uδ(B(j)​)={
k:W~(j)̸​=0}为${\tilde{W}}^{(j)}$挑选的特征。重复$\gamma$次，对于每一个特征都会得到特征稳定性${\hat{\prod }}_k^{\delta }$。
$${\hat{\prod }}_k^{\delta }=\sum _{j=1}^{\gamma }\frac{I(k\in U^{\delta }(B_{(j)}))}{\gamma }$$
其中：
$$I(c)=\{\begin{array}{rcl}1& & cistrue\\ 0& & otherwise\end{array}$$
这里${\hat{\prod }}_k^{\delta }$其实计算的是特征k被选中的分数。对于所有的$\delta \in \Delta$，重复以上过程。则对于每一特征都会有一个稳定性分数：$S(k)={\max }_{\delta \in \Delta }({\hat{\prod }}_k^{\delta })$。对于稳定特征集合可以定义为： U ^ s t a b l e = { k : S ( k ) 在 K 中 的 前 η 个 } \hat U^{stable}=\{k: S(k)在K中的前\eta个\} U^stable={
k:S(k)在K中的前η个}，或者可以设定阈值${\pi }_{thr}$，则稳定特征集合可以表示为 U ^ s t a b l e = { k : S ( k ) ≥ π t h r } \hat U^{stable}=\{k: S(k) \ge \pi_{thr}\} U^stable={
k:S(k)≥πthr​}。本文实验证明稳定性分数排前20的特征可以满足本文回归模型的需求。

2.4 本文算法

传统Lasso可能会导致局部最优解。本文提出的算法分为两个阶段。在第一阶段中，使用纵向稳定性特征选择得到${\hat{U}}^{stable}$。在第二阶段中，基于${\hat{U}}^{stable}$使用带有时间平滑正则回归。

3. 实验

实验数据集为ADNI，对比算法为岭回归。由于患者数目并不是很多，所以使用留一法。5折交叉验证用于挑选参数($10^{-3}$到$10^3$，本文模型参数为${\theta }_1$ 和 ${\theta }_2$)，岭回归参数为${\theta }_1$。P值(correlation significance)和R值(correlation coefficient)作为评价指标。好的模型R值高，P值低。
数据集共计7种，分别为脑脊液特征(CSF，记为C)，MRI影像特征(MRI，记为M)，患者个人信息及之前意识水平特征(META，记为E)以及它们的集合：CE，EM，CM和CEM。MRI特征主要分为5类：平均皮质厚度(cortical thickness average, CTA)，皮质厚度标准差(cortical thickness standard deviation, CTStd)，分割后的皮质体积(volume of cortical parcellation, Vol. Cort.)，分割后的白质体积(volume of white matter parcellation, Vol. WM.)和表面面积(surface area, Surf. A.)，共计313维。

4. 结论

实验结果证明多任务学习效果较好，未来考虑使用非线性模型。