Predicting Alzheimer's Disease Cognitive Assessment via Robust Low-Rank Structured Sparse Model

1 Introduction

Alzheimer's disease (AD), a common form of dementia, affects nerve cells in areas of the brain responsible for memory, cognition, language, and motor activity [Dailey, 2017]. By linear extrapolation of estimates from 2006, the population worldwide who have AD will increase to over 100 million by 2050 [Thompson et al., 2003; Moradi et al., 2015]. In fact, researchers believe that early detection will be key to preventing, slowing and stopping Alzheimer's disease. Neuroimaging as a powerful tool for accurate identification and understanding informative feature is necessary for early Alzheimer's disease prognosis and diagnosis [Liu et al., 2015; Nie et al., 2016]. Therefore, many machine learning methods have been proposed to study neuroimaging measures to detect pathology associated with AD and to predict cognitive scores [Wang et al., 2011b; Huo et al., 2016; Wang et al., 2016]. Among them, structural magnetic resonance imaging (MRI) scans are one of the most extensively used imaging modality in tracking AD progression.

In the association study of selecting effective longitudinal phenotypic markers to predict cognitive scores from imaging features, the input usually consists of two matrices: the imaging feature matrix X = [x1, · · · , xn] ∈ Rd×n and the corresponding cognitive score matrix Y = [y1 , · · · , yn ]T ∈ Rn×m, where n is the number of samples, d is the number of features and m is the number of different measures of a certain cognitive performance.

A forthright method to identify informative imaging markers is to perform feature selection [Chang and Yang, 2016], which has been demonstrated as a useful way to reflect the correlation between cognitive measures after removing the indistinctive neuroimaging markers. More recently, sparse regularization model has been extensively utilized to learn the structure of data and obtain effective feature in different applications. The sparsity-inducing norm based feature selection methods solve the convex optimization problems of the form:

where L is a convex function and Ω can include one or more non-smooth sparsity-inducing norms.

When Ω is the l1-norm, the l1 shrinkage methods such as LASSO identify informative longitudinal phenotypic markers in the brain that are related to pathological changes of AD by imposing flat sparsity [Liu et al., 2014]. However, the selected features distribute randomly among the whole brain that can not be well explained. In predicting cognitive scores task, we expect to select the most informative markers, which are important to all participants including AD, mild cognitive impairment (MCI) and heathy control (HC). To address this issue, group LASSO with a l2,1-norm is used to impose the structured sparsity on parameter matrix for feature selection [Obozinski et al., 2010; Jie et al., 2015; Yang et al., 2017; Chang et al., 2017]. It enforces the important features to have non-zero weights cross all participants, however many important features often are only discriminant to partial classes, i.e. having large weights on these participants. Thus, such important features may be ignored by the above methods.

On this account, Lee et al. and Wang et al. proposed to add one more l1,1-norm regularization term to achieve both structured and flat sparsity [Lee et al., 2010; Wang et al., 2011a]. However, because the l1,1-norm regularization term enforces the flat sparsity and is prone to shrink the non-large values to zeros, the non-zero weights of some important features may also be forced to be zeros, i.e. some features selected by the l2,1-norm regularization term can be totally suppressed by the l1,1-norm regularization term. As a result, many informative longitudinal phenotypic markers are neglected during the feature selection procedure. Thus, more properly designed structured sparsity-inducing norms are desired in feature selection research.

In this paper, we propose a robust low-rank structured sparse regression method (RLSR) to simultaneously select the important neuroimaging markers and learn the underlying structure between cognitive measures. Our main contributions are three-fold: (1) The new mixed structured sparsityinducing norms are introduced to overcome the above overshrinkage drawback in the existing sparse learning based feature selection models. (2) The explicit rank-k low-rank matrix fitting approach is used to extract the underlying interrelation structures between cognitive measures. (3) Because our method leads to a highly non-smooth objective, we derive an efficient algorithm to solve the new objective with proved convergence. We validate our method on cognitive data of the ADNI cohort and obtain promising results.

Notations. We summarize the notations used in this paper. Matrices are written as uppercase letters and vectors are written as bold lowercase letters. For matrix W = {wij}, its i-th row, j-th column are denoted as wi , wj respectively. The lp-norm of the vector v ∈ Rn is defined as ||v||p = (∑ni=1 |vi|p ) 1p for p > 0. The l2,1-norm of matrix W is defined as ||W||2,1 = ∑di=1 ||wi ||2 (in some related papers, people also used the notation l1/l2-norm). l1,2-norm of matrix W is defined as ||W||1,2 = √∑di=1 ||wi ||21 and l1,1-norm of matrix W is defined as ||W||1,1 = ∑di=1 ∑mj=1 |wij |.

2 Robust Low-Rank Structured Sparse Learning

The l2,1-norm based objectives select the informative imaging markers across all the cognitive scores with joint sparsity, i.e. each imaging marker has either small score or large score for all the cognitive measures. However, for accurate identification of effective imaging markers, we utilize participants including AD, MCI and HC during deferent period. Consequently, many features may be irrelevant to each other, which could deteriorate the feature selection performance if we consider all the imaging markers as one group to do feature selection. On the other hand, as we discussed in the introduction section, an extra l1,1-norm regularizer will suppress too many non-zero values and lead unstable feature selection results. For example, when we target to select top 20 features and adjust the trade-off parameter to let l2,1-norm makes about 20 features with relatively large weights, the added l1,1norm will dramatically shrink the weights such that only few important features (much less than 20) can be selected. To tackle this over-shrinkage problem, we add a convex squared l1,2-norm regularizer instead of l1,1-norm and solve:

The typical loss functions are the least square loss and logistic loss. To improve the computational efficiency, we utilize the least square loss in this paper to select informative markers for ADNI data. Thus, our method can be applied to both classification tasks (e.g. AD/MCI versus Normal Controls (NC)) and regression tasks (e.g., estimation of clinical cognitive scores). We perform the latter in this paper.

In Eq. (2), we proposed the novel mixed structured sparsity norms. The standard l2,1-norm enforces the joint sparsity across all cognitive measures to select imaging markers. The new l1,2-norm uses l2-norm between markers, such that at least one non-zero element in the rows of B selected by l2,1-norm regularizer will be kept. Thus, we won't lose the discriminative imaging markers selected by the l2,1-norm regularization. At the same time, the l1,2-norm imposes the l1-norm between cognitive score weights of each marker to shrink the weight values of uncorrelated or irrelevant cognitive measures. For illustration purpose, in Fig. 1, we plot the sparse shrinkage patterns of the matrix B using: (a) l2,1-norm regularizer only, (b) l2,1- norm and l1,1-norm regularizers, (c) l2,1-norm and l1,2-norm regularizers. The l1,1-norm overshrinks the weights and removes the first feature selected by l2,1-norm. The l1,2-norm suppress some weights with supporting the results of l2,1-norm, e.g. the first feature is still kept in the list.

More important, with regard to this specific task that predicting AD progression, we hope to extract and utilize the underlying interrelations between cognitive measures to enhance the accuracy of feature selection. In recent research [Ji and Ye, 2009; Deng et al., 2015], the trace norm regularization has been used to seek the low-rank structured shared common representations:

However, there are two deficiencies: 1) the optimal low-rank approximation is resulted by tediously tuning the parameter λ, which has no direct connection to the rank value; 2) the feature selection isn't enforced in this model.

To address the above problems, we consider the following low-rank regression:

where m is the number of classes and d is the dimension of features. Compared to the parameter λ in Eq. (3), the parameter s is more feasible to be decided by users. Moreover, in order to select features and simultaneously keep the low-rank matrix fitting, we propose the following model,

where W ∈ Rd×s , P ∈ Rs×m and s ≤ min(m, d). The product B = WP is a low-rank matrix with rank(B) ≤ s. Our new objective simultaneously learns the underlying interrelation between cognitive measures by low-rank matrix fitting and selects the informative neuroimaging markers by mixed structured sparsity-inducing norms. We also consider the real world data often have outliers and hence replace the least square loss by the l2,1-norm based loss function, which imposes the l1-norm between data points to reduce the effect of outliers. Our final objective is to solve:

The resulted objective has three non-smooth terms, such that the optimization becomes difficult. To solve our new objective, we will derive an efficient algorithm in next section with proved convergence.

3 Optimization

In this section, an efficient algorithm is proposed to tackle Eq. (6), followed by the proof of its convergence.

3.1 AlgorithmDerivation

We will alternatively and iteratively solve Eq. (6). To begin with, we rewrite Eq. (6) as:

then H ∈ Rn×n is a diagonal matrix and defined as:

where ei (∀i = 1, 2, ..., n) is the i-th row of matrix E in Eq. (8). ˆD ∈ Rd×d is a diagonal matrix and defined as

and Dk ∈ Rd×d is also a diagonal matrix and defined as:

The first step is to fix P,H, ˆD,Dk, and to solve W. Thus, we need solve the following subproblem:

where zk is the k-th column of the matrix XHY P T , ∀k = 1, 2, ..., s. Then taking derivative of Eq. (12) w.r.t. wk and setting it to zero, we get,

The second step is to fix W,H, ˆD,Dk, and to solve P. Because

Then the subproblem becomes:

The solution to Eq. (16) can be obtained by the Theorem 1. The third step is to fix W, P , and to solve H by Eq. (8) and Eq. (9), solve ˆD by Eq. (10), and solve Dk by Eq. (11).

We repeat the above three steps iteratively, until the predefined stopping criterion is satisfied. We summarize the whole algorithm in Alg. 1.

The Step 2 in the iteration of Alg. 1 can be calculated by linear equation system, which can be efficiently solved. Thus, our algorithm can be applied in large-scale datasets.

Algorithm 1 The algorithm to solve Eq. (6)

Input:

1. The training data X ∈ Rd×n with label matrix Y ∈ Rn×m
2. The regularization parameters λ1 and λ2, and the rank s.

Output:

1. The matrices W ∈ Rd×s and P ∈ Rs×m.

Initialization:

1. Set t = 0, initialize H(t) = In×n, ˆD(t) = Id×d, D(t) k = Id×d, ∀k = 1, ..., s, randomly initialize P (t) ∈ Rs×m with P (t)P (t)T = I ∈ Rs×s .

Repeat:

1. Calculate z(t) k , which is the k-th column of the matrix XH(t)Y P (t)T .
2. Calculate W (t) column by column by w(t) k = (XH(t)XT + λ1 ˆD(t) + λ2D(t) k )-1 z(t) k .
3. Calculate M (t) = Y T H(t)XT W (t) .
4. Do SVD of M (t) , M (t) = U (t)Σ(t)V (t)T .
5. Update P (t+1) = V (t) [I, 0]U (t)T .
6. Update E(t+1) = Y -XT W (t)P (t+1) .

1. Update t = t + 1. Until Converge

3.2 Convergence Analysis

Theorem 1. The solution to the optimization problem maxP Tr(PM) s.t. PP T = I, where p ∈ Rs×m, M ∈

Rm×s and s < m is

where U and V are the SVD of M, M = UΛV T and I is the identify matrix, I ∈ Rs×s . 0 ∈ Rs×(m-s) is the matrix with all zeros entries. And [Φ,Ψ] is the matrix operation to horizontally concatenate two matrices Φ and Ψ who have the same number of rows.

Proof. We do SVD of M, M = UΛV T , where U ∈ Rm×m, Λ ∈ Rm×s , V ∈ Rs×s . Then, we have

where Q = V T PU, Q ∈ Rs×m. Note that s < m. λkk and qkk are the k-th element in the diagonal of matrix Λ and Q respectively. Moreover,

where I ∈ Rs×s is the identity matrix. Thus, qkk ≤ 1, ∀k, k = 1, 2, ..., s. Therefore,

and when qkk = 1, ∀k, k = 1, 2, ..., s, the equality holds. In other words, Tr(PM) reaches the maximum when Q = [I, 0]. Recall that Q = V T PU, thus the optimal solution to Eq. (16) is Eq. (17).

The convergence of the Alg. 1 is summarized in the following theorem:

Theorem 2. The Alg. 1 will monotonically decrease the objective of the problem in Eq. (6) in each iteration and converge to the local optimum solution to the problem.

Proof. On one hand, denote the updated P by ˜ P. Because of Theorem 1, when we fixW, H, ˆD and Dk, we get:

After plugging the definition of H by Eq. (9), we can obtain,

At the same time, beginning with (||˜ei ||2 - ||ei ||2)2 ≥ 0, we can get

Therefore, adding Eq. (22) and Eq. (23) together, we get

On the other hand, denote the updated W by ˜W, when we fix P, H. According to the step 2 in the Repeat part in Alg. 1, we have:

where zk is the k-th column of XHY P T .

We plug in the definition of ˆD and Dk by Eq. (10) and Eq. (11) respectively into Eq. (25), and have:

Similarly, beginning with (|| ˜wi ||2 - || ˜wi ||2)2 ≥ 0, we get

Therefore,

Meanwhile, according to am-gm inequality, we have

Thus,

Sum up Eq. (26), λ1×Eq. (28), λ2×Eq. (30), we can arrive at the following conclusion:

We update W and P alternatively, then we arrive at our goal:

In other words, using Alg. 1, we can monotonically decrease the objective function Eq. (6) in each iteration and finally it will converge.

4 Experimental Results

In this section, we evaluate prediction performance of the proposed method by applying it to Alzheimer's Disease Neuroimaging Initiative (ADNI) cohort (adni.loni.usc.edu), where a wide range of imaging markers measured over a period of 2 years are examined and associated to cognitive scores that are relevant to AD.

4.1 Data Descriptions

We apply the proposed method to the ADNI cohort to predict the cognitive scores of the participants from each of their two types of imaging phenotypes, i.e. FreeSurfer markers and voxel-based morphometry (VBM) markers. The detailed information are shown in Table 1. Mean modulated gray matter measures obtained from 90 target regions of interest, normalized by the total intracranial volume, were extracted as features.

4.2 Performance Comparison on the ADNI Cohort

First, we intend to identify a certain set of informative markers that are closely relate to pathological change due to AD. We compared our method against three most related algorithms including multivariate ridge regression (RR), joint l2,1-norm minimization (l2,1) on both loss function and regularization [Nie et al., 2010], linear regression with trace norm. These comparing methods are all widely used in statistical learning and brain image analysis.

In all experiments, we automatically tune the regularization parameters by selecting among the values {10r : r ∈ {-5, ..., 5}} with standard 5-fold cross-validation strategy.After the algorithm converges, we sort the row index of matrix W by the summation of the absolute values in each row, and features are selected by the top ranked indices. To measure prediction performance, we compute the root mean square error (RMSE) between the predicted score and the ground truth.

The prediction experiment evaluated by ridge regression is repeated for 100 times and average results are reported in Fig. 2. As shown in Fig. 2, we can clearly see that the prediction results of our method consistently outperforms other competing methods in nearly all the test cases for all the cognitive tasks except some outlier part in Fig. 2f. But it doesn't really matter since the proposed method catches up soon. The reasons why the proposed method performs best go as follows: RR assumed the cognitive measures to be independent at different time point which neglects the correlations along the time. And for l2,1, since the pathological change of brain structures due to AD usually do not occur in the pre-identified regions with certain shapes, thus it is difficult to define meaningful feature groups. This makes l2,1 performworse. Trace norm is a good way to seek the underlying interrelations between cognitive scores, but it ignores the fact that the informative markers relate to AD among all the imaging measures only occupy a small part. As for the proposed method, we not only detect group structure within longitudinal phenotypic neuroimaging markers, but also capture the correlations among cognitive measures. In addition, for ease of comparison, we also list the RMSE using top 10 and top 20 selected features evaluated by ridge regression in Table 2.

4.3 Identification of Informative Markers

The primary goal of the proposed method is to identify the informative markers which is important for AD diagnosis and prediction. Therefore, we examine the imaging markers selected by our method and show it in Fig. 3. Due to the limit size of display, we only provide one tenth of feature names for both FreeSurfer and VBM markers. As shown in Fig. 3, we observe that hippocampal measures (LHippocampus, RHippocampus, LHippVol and RHippVol) are among the top selected features. These findings are in accordance with the known knowledge that in the pathological pathway of AD, hippocampal is one of sections that can recognize Alzheimer-related changes [Braak and Braak, 1991; Delacourte et al., 1999].

In summary, the identified neuroimaging markers are highly suggestive and effective for tracking the progression of AD, since it strongly agrees with the existing research findings. It also illustrates the necessity and correctness of the selected imaging cognitive associations to reveal the relationships between MRI measures and cognitive scores.

5 Conclusion

To reveal the relationship between cognitive measures and neuroimaging markers, we proposed a novel robust low-rank structured sparse regression model, which selects the most informative imaging markers to predict the cognitive scores for complex brain disorders. Using the new mixed structured sparsity inducing norms and the low-rank approximation function, the proposed method can efficiently identify the effective neuroimaging markers with utilizing the underlying interrelation structures between different cognitive measures.

In addition, we provide an efficient algorithm with proved convergence. Validation experiments conducted on multiple data demonstrate the promise of the proposed method.

Acknowledgements

This work was partially supported by the National Natural Science Foundation of China 61572388, and U.S. NIH R01 AG049371, NSF IIS 1302675, IIS 1344152, DBI 1356628, IIS 1619308, IIS 1633753.

40