An Tutorial for Regularized Multi-task Learning using the package RMTL
Abstract
This package provides an efficient implementation of regularized multi-task learning (MTL) comprising 10 algorithms applicable for regression, classification, joint feature selection, task clustering, low-rank learning, sparse learning and network incorporation (Cao, Zhou, and Schwarz 2018). All algorithms are implemented using the accelerated proximal algorithm and feature a complexity of O(1/k^2). In this tutorial, we show the theory and modeling of MTL with examples using the RMTL.
Introduction
This package provides 10 multi-task learning algorithms (5 classification and 5 regression), which incorporate five regularization strategies for knowledge transferring among tasks. All algorithms share the same framework:
$$\min\limits_{W, C} \sum_{i=1}^{t}{L(W_i, C_i |X_i, Y_i)} + \lambda_1\Omega(W) + \lambda_2{||W||}_F^2$$ where L(∘) is the loss function (logistic loss for classification or least square loss for linear regression). Ω(∘) is the cross-task regularization for knowledge transfer, and ||W||F2 is used for improving the generalization. X = {Xi = ni × p|i ∈ {1, ..., t}} and Y = {Yi = ni × 1|i ∈ {1, ..., t}} are predictors matrices and responses of t tasks respectively, while each task i contains ni subjects and p predictors. W = p × t is the coefficient matrix, where Wi is the ith column of W refers to the coefficient vector of task i. Ω(W) jointly modulates multi-task models({W1, W2, ..., Wt}) according to the specific prior structure. In this package, 5 common cross-task regularization methods are implemented to incorporate different priors, i.e. sparse structure, joint predictor selection, low-rank structure, network-based relatedness across tasks and task clustering. The mathematical formulation of each prior structure is demonstrated in the first row of Table 1. For all algorithms, we implemented an solver based on the accelerated gradient descent method, which takes advantage of information from the previous two iterations to calculate the current gradient and then achieves an improved convergent rate (Beck and Teboulle 2009). To solve the non-smooth and convex regularizer, the proximal operator is applied. Moreover, backward line search is used to determine the appropriate step-size in each iteration. Overall, the solver achieves a complexity of ($\frac{1}{k^2}$) and is optimal among first-order gradient descent methods.
All algorithms are summarized in Table 1. Each column corresponds to a MTL algorithm with an specific prior structure.To run an algorithm correctly, users need to tell RMTL the regularization type and problem type (regression or classification), which are summarized in the second and third row of Table 1. λ1 aims to control the effect of cross-task regularization and could be estimated by the cross-validation (CV) procedure, and λ2 is set by users in advance. Their valid ranges are summarized in the table. Please note, for λ1 = 0, the effect of cross-task regularization was canceled. For the algorithm regularized by network and clustering prior (Graph and CMTL), extra information need to provide by users, i.e. G encodes the network information and k assumes the complexity of clustering structure. To train a multi-task model regularized by Lasso, Trace and L21, the sparsity is introduced by optimizing a non-smooth convex term. Therefore, the warm-start function is provided as part of the training procedure for these three algorithms. For this part, We will explain more details in the next section.
| sparse structure | joint feature learning | low-rank structure | network incorperation | task clustering | |
|---|---|---|---|---|---|
| Ω(W) | ||W||1 | ||W||2, 1 | ||W||* | ||WG||F2 | tr(WTW) − tr(FTWTWF) |
| Regularization Type | Lasso | L21 | Trace | Graph | CMTL |
| Problem Type | R/C | R/C | R/C | R/C | R/C |
| λ1 | λ1 > 0 | λ1 > 0 | λ1 > 0 | λ1 ≥ 0 | λ1 > 0 |
| λ2 | λ2 ≥ 0 | λ2 ≥ 0 | λ2 ≥ 0 | λ2 ≥ 0 | λ2 > 0 |
| Extra Input | None | None | None | G | k |
| Warm Start | Yes | Yes | Yes | No | No |
| Reference | (Tibshirani 1996) | (Liu, Ji, and Ye 2009) | (Pong et al. 2010) | (Widmer et al. 2014) | (Jiayu Zhou, Chen, and Ye 2011) |
Cross-validation, Training and Prediction
For all algorithms in RMTL package, λ1 illustrates the strength of relatedness of tasks and high value of λ1 would result in highly similar models. For example, for MTL with low-rank structure (Trace), a large enough λ1 will compress the task space to 1-dimension, thus coefficient vectors of all tasks are proportional. In RMTL, an appropriate λ1 could be estimated using CV based on training data. λ2 is suggested tuning manually by users with the default value of 0 (except for MTL with CMTL). The purpose of λ2 is to introduce the penalty of quadratic form of W leading to several benefits, i.e. promoting the grouping effect of predictors for selecting correlated predictors (Zou and Hastie 2005), stabilizing the numerical results and improving the generalization performance.
Cross-validation
To run the cross-validation and training procedure successfully, the problem type and regularization type have to be given as the arguments to the function cvMTL(X, Y, Regularization=Regulrization, type=problem, …) for cross-validation and function MTL(X, Y, Regularization=Regulrization, type=problem …) for training, while the valid value of problem can be only “Regression” or “Classification” and the value of Regularization has to be selected from Table 1.
In the following example, we show how to use the function cvMTL(X, Y, …) for cross validation and demonstrate the selected parameter as well as CV accuracy curve.
#create simulated data
library(RMTL)
datar <- Create_simulated_data(Regularization="L21", type="Regression")
#perform the cross validation
cvfitr <- cvMTL(datar$X, datar$Y, type="Regression", Regularization="L21", Lam1_seq=10^seq(1,-4, -1), Lam2=0, opts=list(init=0, tol=10^-6, maxIter=1500), nfolds=5, stratify=FALSE, parallel=FALSE)
# meta-information and results of CV
#sequence of lam1
cvfitr$Lam1_seq
#> [1] 1e+01 1e+00 1e-01 1e-02 1e-03 1e-04
#value of lam2
cvfitr$Lam2
#> [1] 0
#the output lam1 value with minimum CV error
print (paste0("estimated lam1: ", cvfitr$Lam1.min))
#> [1] "estimated lam1: 0.1"
#plot CV errors across lam1 sequence in the log space
plot(cvfitr)
Figure 1: CV errors across the sequence of λ1
Stratified Cross-validation procedure is provided specific to the classification problem. The positive and negative subjects are uniformly distributed across folds such that the the performance of parameter selection is not biased by the data imbalance. To turn on the Stratified CV, users need to set cvfit(…,type=“Classification”, stratify=TRUE).
Parallel Computing is allowed to speed up the CV procedure. To run the procedures of k folds simultaneously, the training data has to be replicated k times leading to the dramatically increased memory use. Therefore, users are suggested selecting less cores if the data size is large enough. To trigger on the parallel computing and select, i.e. 3 cores, one has to set cvfit(…,parallel=TRUE, ncores=3).
In the following example, we show how to use Stratified CV and Parallel Computation in practice. We will compare the time consumption between the CV procedure with and without parallelization. Please note, in the Vignette, we are only allowed to use up to 2 cores for parallelization, thus the time comparison is less significant, users are suggested to use k cores in k-folds CV in practice.
datac <- Create_simulated_data(Regularization="L21", type="Classification", n=100)
# CV without parallel computing
start_time <- Sys.time()
cvfitc<-cvMTL(datac$X, datac$Y, type="Classification", Regularization="L21", stratify=TRUE, parallel=FALSE)
Sys.time()-start_time
#> Time difference of 0.1968799 secs
# CV with parallel computing
start_time <- Sys.time()
cvfitc<-cvMTL(datac$X, datac$Y, type="Classification", Regularization="L21", stratify=TRUE, parallel=TRUE, ncores=2)
Sys.time()-start_time
#> Time difference of 0.1846242 secsTraining
The result of CV (cvfit$Lam1.min) is sent to function MTL(X, Y, …) for training. After this, the coefficient matrices of all tasks {W, C} are obtained and ready for predicting new individuals. In the following example, we demonstrate the procedure of model training, the learnt parameters and the historical values of objective function.
#train a MTL model
model<-MTL(datar$X, datar$Y, type="Regression", Regularization="L21",
Lam1=cvfitr$Lam1.min, Lam2=0, opts=list(init=0, tol=10^-6,
maxIter=1500), Lam1_seq=cvfitr$Lam1_seq)
#demo model
model
#>
#> Head Coefficients:
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] -0.05385028 0.003018368 -0.0590724314 -0.08045033 -0.02646358
#> [2,] 0.00000000 0.000000000 0.0000000000 0.00000000 0.00000000
#> [3,] 0.00000000 0.000000000 0.0000000000 0.00000000 0.00000000
#> [4,] -0.07421000 -0.029701551 -0.1418335332 -0.03545683 -0.08477055
#> [5,] -0.14316117 -0.048637249 0.0768739870 0.01546224 -0.27557140
#> [6,] 0.00768479 -0.006718799 0.0009178763 -0.01154268 -0.01698838
#> Call:
#> MTL(X = datar$X, Y = datar$Y, type = "Regression", Regularization = "L21",
#> Lam1 = cvfitr$Lam1.min, Lam1_seq = cvfitr$Lam1_seq, Lam2 = 0,
#> opts = list(init = 0, tol = 10^-6, maxIter = 1500))
#> type:
#> [1] "Regression"
#> Formulation:
#> [1] "SUM_i Loss_i(W) + Lam1*||W||_{2,1} + Lam2*||W||{_2}{^2}"
# learnt models {W, C}
head(model$W)
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] -0.05385028 0.003018368 -0.0590724314 -0.08045033 -0.02646358
#> [2,] 0.00000000 0.000000000 0.0000000000 0.00000000 0.00000000
#> [3,] 0.00000000 0.000000000 0.0000000000 0.00000000 0.00000000
#> [4,] -0.07421000 -0.029701551 -0.1418335332 -0.03545683 -0.08477055
#> [5,] -0.14316117 -0.048637249 0.0768739870 0.01546224 -0.27557140
#> [6,] 0.00768479 -0.006718799 0.0009178763 -0.01154268 -0.01698838
head(model$C)
#> [1] -0.04777385 -0.44886808 0.24526715 -0.12522320 0.14917324
# Historical objective values
str(model$Obj)
#> num [1:98] 2.24 1.9 1.79 1.74 1.71 ...
# other meta infomration
model$Regularization
#> [1] "L21"
model$type
#> [1] "Regression"
model$dim
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] 20 20 20 20 20
#> [2,] 50 50 50 50 50
str(model$opts)
#> List of 5
#> $ init : num 1
#> $ tol : num 1e-06
#> $ maxIter: num 1500
#> $ W0 : num [1:50, 1:5] -0.0539 0 0 -0.0742 -0.1432 ...
#> $ C0 : num [1:5] -0.0478 -0.4489 0.2453 -0.1252 0.1492
#plot the historical objective values in the optimization
plotObj(model)
Figure 2: The convergence of objective values across iterations
Prediction and Error Calculation
To predict the new individuals, two functions are relevant, predict() and calcError(). predict() is used to calculate the outcome given the data of new individuals. Especially, for classification, the outcome is represented as the probability of the individual being assigned to the positive label (P(Y = = 1)). calcError() is used to test if the model works well by calculating the prediction error rate. According to the problem type, specific metric is used, for example, the miss-classification rate ($\frac{1}{n} \sum_i^n 1_{\hat{y}_i=y_i}$) is used for classification problem and the mean square error (MSE) ($\sum_i^n (\hat{y}_i-y_i)^2$) is used for the regression problem. In the following examples, we will show you how to calculate the training and test error, and then make predictions for both classification and regression problem. To achieve it, we have to create the simulated data and train the models in advance just like we did before.
# create simulated data for regression and classification problem
datar <- Create_simulated_data(Regularization="L21", type="Regression")
datac <- Create_simulated_data(Regularization="L21", type="Classification")
# perform CV
cvfitr<-cvMTL(datar$X, datar$Y, type="Regression", Regularization="L21")
cvfitc<-cvMTL(datac$X, datac$Y, type="Classification", Regularization="L21")
# train
modelr<-MTL(datar$X, datar$Y, type="Regression", Regularization="L21",
Lam1=cvfitr$Lam1.min, Lam1_seq=cvfitr$Lam1_seq)
modelc<-MTL(datac$X, datac$Y, type="Classification", Regularization="L21",
Lam1=cvfitc$Lam1.min, Lam1_seq=cvfitc$Lam1_seq)
# test
# for regression problem
calcError(modelr, newX=datar$X, newY=datar$Y) # training error
#> [1] 0.04721407
calcError(modelr, newX=datar$tX, newY=datar$tY) # test error
#> [1] 0.7448988
# for calssification problem
calcError(modelc, newX=datac$X, newY=datac$Y) # training error
#> [1] 0
calcError(modelc, newX=datac$tX, newY=datac$tY) # test error
#> [1] 0.18
# predict
str(predict(modelr, datar$tX)) # for regression
#> List of 5
#> $ : num [1:20, 1] -1.703 0.945 -0.594 1.936 0.379 ...
#> $ : num [1:20, 1] -0.858 2.639 -3.131 1.453 -2.094 ...
#> $ : num [1:20, 1] -0.7481 -0.8172 0.0988 -0.5381 1.2912 ...
#> $ : num [1:20, 1] -2.8349 0.648 -3.7804 4.2249 -0.0662 ...
#> $ : num [1:20, 1] -3.369 2.276 1.408 1.347 -0.799 ...
str(predict(modelc, datac$tX)) # for classification
#> List of 5
#> $ : num [1:20, 1] 0.9592 0.296 0.504 0.0546 0.619 ...
#> $ : num [1:20, 1] 0.25185 0.88845 0.9925 0.93418 0.00439 ...
#> $ : num [1:20, 1] 0.975 0.2872 0.6106 0.0447 0.7386 ...
#> $ : num [1:20, 1] 0.00482 0.96372 0.01087 0.86321 0.99955 ...
#> $ : num [1:20, 1] 0.8668 0.2086 0.9027 0.0167 0.7219 ...Control of Optimization
Options are provided to control the optimization procedure using the argument opts=list(init=0, tol=10^-6, maxIter=1500). opts$init specifies the starting point of the gradient descent algorithm, opts$tol controls tolerance of the acceptable precision of solution to terminate the algorithm, and opts$maxIter is the maximum number of iterations. These options can be modified by users for adapting to the scale of problem and computing resource. For opts$init, two options are provided. opts$init=0 refers to 0 matrix. And opts$init=1 refers to the user-specific starting point. If specified, the algorithm will attempt to access opts$W0 and opts$C0 as the starting point, which has to be given by users in advance. Otherwise, errors are reported. Particularly, the setting opts$init=1 is key to warm-start technique for sparse model training.
In the following example, we aim to calculate a highly precise solution by restricting the tolerance and increasing the number of iterations. In the left figure, the precision of the solution is 0.02 and the algorithm stops after <20 iterations. In the right figure, the precision is set to 10^-8 and the algorithm takes more iterations to run.
par(mfrow=c(1,2))
model<-MTL(datar$X, datar$Y, type="Regression", Regularization="L21",
Lam1=cvfitr$Lam1.min, Lam2=0, opts=list(init=0, tol=10^-2,
maxIter=10))
plotObj(model)
model<-MTL(datar$X, datar$Y, type="Regression", Regularization="L21",
Lam1=cvfitr$Lam1.min, Lam2=0, opts=list(init=0, tol=10^-8,
maxIter=100))
plotObj(model)
Figure 3: Calculate rough (Left) and precise (Right) by controling the optimization
warm-start and cold-start are both provided to train a sparse model. The reason for warm-start is due to the non_strict convexity of the objective function, more than one (usually unlimited) number of optimal solutions exist, which made the model difficult to interpret. In the warm-start, λ1 sequence ({..., λ1(i), λ1(i − 1), ...}) is feed to the training procedure in the decreasing order. Meanwhile, the solution of λ1(i), $\{\overset{*}{W}^{(i)}, \overset{*}{C}^{(i)}\}$, is set as the initial point of the λ1(i − 1), which leads to a unique solution path. In the example, we investigate the algorithm of multi-task feature selection and draw a regularization tree to show the model interpretation. In Figure 4, each line represents the change of significance of one predictor across the λ1 sequence. By relaxing λ1, more predictors are allowed to contribute to the prediction. Please note, significance per predictor is calculated as the euclidean norm of parameters across tasks.
Lam1_seq=10^seq(1,-4, -0.1)
opts=list(init=0, tol=10^-6,maxIter=100)
mat=vector()
for (i in Lam1_seq){
m=MTL(datar$X, datar$Y, type="Regression", Regularization="L21", Lam1=i,opts=opts)
opts$W0=m$W
opts$C0=m$C
opts$init=1
mat=rbind(mat, sqrt(rowSums(m$W^2)))
}
matplot(mat, type="l", xlab="Lambda1", ylab="Significance of Each Predictor")
Figure 4: The solution path using warm-start
The warm-start is triggered by sending two arguments (λ1, λ1 sequence) to the function MTL(…, Lam1=Lam1_value, Lam1_seq=Lam1_sequence), and they can be calculated by the CV procedure or set manually. However, if only Lam1=Lam1_value is given, the sparsity is introduced using cold-start (the initial point is set only via opts$init). See examples below, different setting leads to different solution.
#warm-start
model1<-MTL(datar$X, datar$Y, type="Regression", Regularization="L21", Lam1=0.01, Lam1_seq=10^seq(0,-4, -0.1))
str(model1$W)
#> num [1:50, 1:5] 0.01347 0.11564 0.10777 0.16674 -0.00411 ...
#cold-start
model2<-MTL(datar$X, datar$Y, type="Regression", Regularization="L21", Lam1=0.01)
str(model2$W)
#> num [1:50, 1:5] 0.0578 0.0394 0.0621 0.1085 -0.1761 ...Data Simulation Methods
Create_simulated_data(…) is used to generate examples for testing every MTL algorithm. To create datasets of variant number of tasks, subjects and predictors, the function accepts several arguments: t, p and n illustrating the number of tasks, predictors and subjects. Here, for convenience, we only allow same number of subjects for all tasks. To test the specific functionality of each MTL algorithm, the simulation datasets are generated with the specific prior structure. For details of this part, users are suggested to check out the supplementary of the original paper (Cao, Zhou, and Schwarz 2018). The problem type and regularization type need to be specific as the arguments to create an algorithm-specific dataset. Finally, 5 objects are outputted: (X, Y, tX, tY, W), where (X, Y) are used for training, (tX, tY) are used for testing, and W is the ground truth model. In addition, for multi-task algorithms regularized by Graph and CMTL, extra information is also generated and outputted i.e. matrix G and k.
In the following example, we demonstrate how to create datasets to test the MTL regression method with low-rank structure (Trace), including 100 tasks, as well as 10 subjects, 20 predictors for each task.
Multi-task Learning Algorithms
$$\min\limits_{W, C} \sum_{i=1}^{t}{L(W_i, C_i|X_i, Y_i)} + \lambda_1\Omega(W) + \lambda_2{||W||}_F^2$$
Recall the shared formulation of all algorithms in RMTL, loss function L(∘) could be logistic loss for classification problem $$L(W_i, C_i)= \frac{1}{n_i} \sum_{j=1}^{n_i} log(1+e^{-Y_{i,j}(X_{i,j}W_i^T+C_i)})$$ or least square loss for regression problem. $$L(W_i, C_i)= \frac{1}{n_i} \sum_{j=1}^{n_i} ||Y_{i,j}-X_{i,j}W_i^T-C_i||^2_2$$ where i indexes tasks and and j indexes subjects in each task, therefore Yi, j and Xi, j = 1 × p refer to the outcome and predictors of subject j in task i, while ni refer to the number of subjects in task i.
In the following sections, we will show you how to train a multi-task model with specific cross-task regularization and how to interpret it. Please note, in the Vignettes of R, we are only allowed to demonstrate algorithms on small datasets leading to the less significant results (especially limiting the performance of cross-validation procedure). For the complete analysis, users are directed to the original paper (Cao, Zhou, and Schwarz 2018).
MTL with sparse structure
$$\min\limits_{W, C} \sum_{i=1}^{t}{L(W_i, C_i|X_i, Y_i)} + \lambda_1||W||_1 + \lambda_2{||W||}_F^2$$ This formulation extends the lasso (Tibshirani 1996) to the multi-task scenario such that all models are penalized according to the same L1 strength, and all unimportant coefficients are shrunken to 0, to achieve the global optimum. || ∘ ||1 is the L1 norm.
In the example, we train and test this algorithm and demonstrate the model in the regression problem.
#create data
data <- Create_simulated_data(Regularization="Lasso", type="Regression")
#CV
cvfit<-cvMTL(data$X, data$Y, type="Regression", Regularization="Lasso")
cvfit$Lam1.min
#> [1] 1e-04
#Train
m=MTL(data$X, data$Y, type="Regression", Regularization="Lasso", Lam1=cvfit$Lam1.min, Lam1_seq=cvfit$Lam1_seq)
#Test
paste0("test error: ", calcError(m, data$tX, data$tY))
#> [1] "test error: 17.9476086317892"
#Show models
par(mfrow=c(1,2))
image(t(m$W), xlab="Task Space", ylab="Predictor Space")
title("The Learnt Model")
image(t(data$W), xlab="Task Space", ylab="Predictor Space")
title("The Ground Truth")
Figure 5: Comparision of Ŵ and W
In Figure 5, the learnt model is more sparse than the ground truth due to the fact that L1 penalty only selects one from correlated predictors such that most correlated predictors are penalized to 0 (Zou and Hastie 2005). In this case, users are suggested to use λ2 to turn on the multi-task learning algorithm with elastic net. More analysis could be found in the original paper (Cao, Zhou, and Schwarz 2018).
MTL with joint predictor selection
$$\min\limits_{W, C} \sum_{i=1}^{t}{L(W_i, C_i|X_i, Y_i)} + \lambda_1||W||_{2,1} + \lambda_2{||W||}_F^2$$ This formulation constraints all models to select or reject the same set of features simultaneously. Therefore, the solution only contains predictors which are consistently important to all tasks. ||W||2, 1 = ∑k||Wk, :||2 aims to create the group sparse structure in the feature space. In the multi-task learning scenario, the same feature of all tasks form a group (each row of W), such that features in the same group are equally penalized, while results in the group-wise sparsity. This approach has been used in the transcriptomic markers mining across multiple platforms, where the common set of biomarkers are shared (Xu, Xue, and Yang 2011). The approach observed the improved prediction performance compared to the single-task learning methods.
#create datasets
data <- Create_simulated_data(Regularization="L21", type="Regression")
#CV
cvfit<-cvMTL(data$X, data$Y, type="Regression", Regularization="L21")
#Train
m=MTL(data$X, data$Y, type="Regression", Regularization="L21", Lam1=cvfit$Lam1.min, Lam1_seq=cvfit$Lam1_seq)
#Test
paste0("test error: ", calcError(m, data$tX, data$tY))
#> [1] "test error: 0.434926903996178"
#Show models
par(mfrow=c(1,2))
image(t(m$W), xlab="Task Space", ylab="Predictor Space")
title("The Learnt Model")
image(t(data$W), xlab="Task Space", ylab="Predictor Space")
title("The Ground Truth")
Figure 6: Comparision of Ŵ and W
In the Figure 6, the ground truth model consists of top 5 active predictors while all others are 0. By assuming all tasks sharing the same set of predictors, MTL with L21 successfully locates most true predictors. The complete analysis could be found here (Cao, Zhou, and Schwarz 2018).
MTL with low-rank structure
$$\min\limits_{W, C} \sum_{i=1}^{t}{L(W_i, C_i|X_i, Y_i)} + \lambda_1||W||_* + \lambda_2{||W||}_F^2$$ This formulation constraints all models to a low-rank subspace. With increasing penalty(λ1), the correlation between models increases. ||W||* is the trace norm of W. This method has been applied to improve the drug sensitivity prediction by combining multi-omic data, which successfully captures the information of drug mechanism of action using the low-rank structure (Yuan et al. 2016).
#create data
data <- Create_simulated_data(Regularization="Trace", type="Classification")
#CV
cvfit<-cvMTL(data$X, data$Y, type="Classification", Regularization="Trace")
#Train
m=MTL(data$X, data$Y, type="Classification", Regularization="Trace", Lam1=cvfit$Lam1.min, Lam1_seq=cvfit$Lam1_seq)
#Test
paste0("test error: ", calcError(m, data$tX, data$tY))
#> [1] "test error: 0.28"
#Show task relatedness
par(mfrow=c(1,2))
image(cor(m$W), xlab="Task Space", ylab="Task Space")
title("The Learnt Model")
image(cor(data$W), xlab="Task Space", ylab="Task Space")
title("The Ground Truth")
Figure 7: Comparision of the learnt task relatedness and ground truth
To interpret this model, we show the correlation coefficient of pairwise models as shown in the Figure 7. Through this comparison, users could check if the task relatedness is well captured.
MTL with network structure
$$\min\limits_{W, C} \sum_{i=1}^{t}{L(W_i, C_i|X_i, Y_i)} + \lambda_1||WG||_F^2 + \lambda_2{||W||}_F^2$$ This formulation constraints the models’ relatedness according to the pre-defined graph . If the penalty is heavy enough, the difference of connected tasks is 0. Network matrix G has to be modeled by users. Intuitively, ||WG||F2 equals to an accumulation of differences between related tasks, i.e. ||WG||F2 = ∑||Wτ − Wϕ||F2, where τ and ϕ are closely connected tasks over an network. In general, given an undirected graph representing the network, ||WG||F2 = tr(WLWT), where L is the graph Laplacian. Therefore, penalizing such term improves the task relatedness.
However, it is nontrivial to design G. Here, we give three common examples to demonstrate the construction of G.
1), Assume your tasks are subject to orders i.e. temporal or spatial order. Such order forces the order-oriented smoothness across tasks such that the adjacent models(tasks) are related, then the penalty can be designed as:
$$||WG||_F^2 = \sum_{i=1}^{t-1}||W_i-W_{i+1}||_F^2$$ where G = (t + 1) × t
$$ G_{i,j} = \begin{cases} 1 & i=j \\ -1 & i=j+1 \\ 0 & others \end{cases} $$ 2), The so-called mean-regularized multi-task learning (Evgeniou and Pontil 2004). In this formulation, each model is forced to approximate the mean of all models.
$$||WG||_F^2 = \sum_{i=1}^{t}||W_i - \frac{1}{t} \sum_{j=1}^t W_j||_F^2$$ where G = t × t
$$ G_{i,j} = \begin{cases} \frac{t-1}{t} & i=j \\ \frac{1}{t} & others \end{cases} $$ 3), Assume your tasks are related according to a given graph g = (N, E) where E is the edge set, and N is the node set. Then the penalty is
$$||WG||_F^2 = \sum_{\alpha \in N, \beta \in N, (\alpha, \beta) \in E}^{||E||} ||W_{:, \alpha} - W_{:,\beta}||_F^2$$ where α and β are connected tasks in the network. G = t × ||E||.
$$ G_{i, j} = \begin{cases} 1 & i=\alpha^{(j)} \\ -1 & j=\beta^{(j)} \\ 0 & others \end{cases} $$ where (α(j), β(j)) is the jth term of set E.
In the bio-medical applications, all three constructions are beneficial. The first construction was used for predicting the disease progression of Alzheimer by holding the assumption of temporal smoothness, which successfully identified the progression-related biomarkers (J. Zhou et al. 2013). The second the construction is applied to explore the transcriptomic biomarkers and risk prediction for schizophrenia by integrating heterogeneous transcriptomic datasets, which demonstrates the strong interpretability of MTL method in multi-modal data analysis (Cao, Meyer-Lindenberg, and Schwarz 2018). The third construction has been used in the recognition of splice sites in genome leading to an improved the prediction performance by considering the task-task relatedness (Widmer et al. 2010).
In the following example, we create 5 tasks forming 2 groups(first group: first two tasks; second group: last three tasks). The tasks are highly correlated within the group and show no group-wise correlation. We will show you how to embed the known task-task relatedness in MTL.
#create datasets
data <- Create_simulated_data(Regularization="Graph", type="Classification")
#CV
cvfit<-cvMTL(data$X, data$Y, type="Classification", Regularization="Graph", G=data$G)
#Train
m=MTL(data$X, data$Y, type="Classification", Regularization="Graph", Lam1=cvfit$Lam1.min, Lam1_seq=cvfit$Lam1_seq, G=data$G)
#Test
print(paste0("the test error is: ", calcError(m, newX=data$tX, newY=data$tY)))
#> [1] "the test error is: 0.28"
#Show task relatedness
par(mfrow=c(1,2))
image(cor(m$W), xlab="Task Space", ylab="Task Space")
title("The Learnt Model")
image(cor(data$W), xlab="Task Space", ylab="Task Space")
title("The Ground Truth")
Figure 8: Compare the Learnt Task Relatedness with the Ground Truth
MTL with clustering structure
$$\min\limits_{W, C, F: F^TF=I_k} \sum_{i=1}^{t}{L(W_i, C_i|X_i, Y_i)} + \lambda_1 (tr(W^TW)-tr(F^TW^TWF)) + \lambda_2{||W||}_F^2$$ The formulation combines the data fitting term (L(∘)) and the loss function of k-means clustering ((tr(WTW) − tr(FTWTWF)), and therefore is able to detect a clustered structure among tasks. λ1 is used to balance the clustering and data fitting effects. However the above formulation is non-convex form leading to the challenge for reaching the global optimum, therefore we implements the relaxed convex form of the formulation (Jiayu Zhou, Chen, and Ye 2011).
$$\min\limits_{W, C, F: F^TF=I_k} \sum_{i=1}^{t}{L(W_i, C_i|X_i, Y_i)} + \lambda_1 \eta (1+\eta)||W(\eta I+M)^{-1}W^T||_* $$ $$S.T. \enspace tr(M)=k, S \preceq I, M \in S^t_+, \eta=\frac{\lambda_2}{\lambda_1}$$ where S+t = t × t denotes the symmetric positive semi-definite matrices, S ≼ I restricts I − M to be positive semi-definite. M = t × t reflects the learnt clustered structure, k refers to the complexity of the structure. In the relaxed form, λ2 = 0 leads to the Ω(W) = 0 which canceled the effect of cross-task regularization, therefore λ2 > 0 is suggested.
In the following example, 5 tasks forms two clusters, and we will show you how to train a MTL model while capturing the clustering structure. To demonstrate the result, the task-task correlation matrix is shown in the Figure 9
#Create datasets
data <- Create_simulated_data(Regularization="CMTL", type="Regression")
cvfit<-cvMTL(data$X, data$Y, type="Regression", Regularization="CMTL", k=data$k)
#> Loading required namespace: MASS
#> Loading required namespace: psych
m=MTL(data$X, data$Y, type="Regression", Regularization="CMTL", Lam1=cvfit$Lam1.min, Lam1_seq=cvfit$Lam1_seq, k=data$k)
#Test
paste0("the test error is: ", calcError(m, newX=data$tX, newY=data$tY))
#> [1] "the test error is: 45.3832811686424"
#Show task relatedness
par(mfrow=c(1,2))
image(cor(m$W), xlab="Task Space", ylab="Task Space")
title("The Learnt Model")
image(cor(data$W), xlab="Task Space", ylab="Task Space")
title("Ground Truth")
Figure 9: Compare the Learnt Task Relatedness with the Ground Truth