Analysing Multi-Task Regression via Random Matrix Theory with Application to Time Series Forecasting

Authors: romain  ilbert
malik  tiomoko
cosme  louart
ambroise  odonnat
vasilii  feofanov
themis  palpanas
ievgen  redko
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-15 16:23:03

Abstract: In this paper, we introduce a novel theoretical framework for multi-task regression, applying random matrix theory to provide precise performance estimations, under high-dimensional, non-Gaussian data distributions. We formulate a multi-task optimization problem as a regularization technique to enable single-task models to leverage multi-task learning information. We derive a closed-form solution for multi-task optimization in the context of linear models. Our analysis provides valuable insights by linking the multi-task learning performance to various model statistics such as raw data covariances, signal-generating hyperplanes, noise levels, as well as the size and number of datasets. We finally propose a consistent estimation of training and testing errors, thereby offering a robust foundation for hyperparameter optimization in multi-task regression scenarios. Experimental validations on both synthetic and real-world datasets in regression and multivariate time series forecasting demonstrate improvements on univariate models, incorporating our method into the training loss and thus leveraging multivariate information.

Main paper: 
 introduction: 
 The principle of multi-task learning, which encompasses concepts such as "learning to learn" and "knowledge transfer," offers an efficient paradigm that mirrors human intelligence. Unlike traditional machine learning, which tackles problems on a task-specific basis with tailored algorithms, multi-task learning leverages shared information across various tasks to enhance overall performance. This approach not only improves accuracy but also addresses challenges related to insufficient data and data cleaning. By analyzing diverse and multimodal data and structuring representations multi-task learning can significantly boost general intelligence, similar to how humans integrate skills. This concept is well-established [33] and has been successfully applied to a wide range of domains, such as computer vision [41] , natural language processing [35, 37] and biology [17, 28, 42] .Our Method. We consider the multi-task regression framework [5] with T tasks with the input space X (t) ⊂ R d and the output space Y (t) ⊂ R q for t ∈ {1, . . . , T }. For each task t, we assume that we are given n t training examples organized into the feature matrix X (t) = [x i ∈ X (t) represents the i-th feature vector of the t-th task and y (t) i ∈ Y (t) is the associated response. In order to have a more tractable and insightful setup, we follow Romera-Paredes et al. [36] and consider the linear multi-task regression model. In particular, we study a straightforward linear signal-plus-noise model that evaluates the response y (t) i for the i-th sample of the t-th task as follows:EQUATIONwhere ε (t) ∈ R nt×q is a matrix of noise vectors with each ε (t)i ∼ N (0, Σ N ), Σ N ∈ R q×q denoting the covariance matrix. The matrix W t ∈ R d×q denotes the signal-generating hyperplane for task t. We denote the concatenation of all task-specific hyperplanes by W = [W ⊤ 1 , . . . , W ⊤ T ] ⊤ ∈ R T d×q . We assume that W t can be decomposed as a sum of a common matrix W 0 ∈ R d×q , which captures the shared information across all the tasks, and a task-specific matrix V t ∈ R d×q , which captures deviations specific to the task t:EQUATIONGiven the multitask regression framework and the linear signal-plus-noise model, we now want to retrieve the common and specific hyperplanes, W 0 and V t , respectively. To achieve this, we study the following minimization problem governed by a parameter λ that controls the balance between the common and specific components of W:EQUATION1. We provide exact train and test risks for the solutions of (3) using random matrix theory. We then decompose the test risk into a signal and noise term responsible for the effectiveness of multi-task learning and the negative transfer, respectively.2. We show how the signal and noise terms compete with each other depending on λ for which we obtain an optimal value optimizing the test risk.3. In a particular case, we derive a closed-form data-dependent solution for the optimal λ ⋆ , involving raw data covariances, signal-generating hyperplane, noise levels, and the size and number of data sets.4. We offer empirical estimates of these model variables to develop a ready-to-use method for hyperparameter optimization in multitask regression problems within our framework.Applications. As an application, we view multivariate time series forecasting (MTSF) as a multitask problem and show how the proposed regularization method can be used to efficiently employ univariate models in the multivariate case. Firstly, we demonstrate that our method improves PatchTST [32] and DLinear [53] compared to the case when these approaches are applied independently to each variable. Secondly, we show that the univariate models enhanced by the proposed regularization achieve performance similar to the current state-of-the-art multivariate forecasting models such as SAMformer [19] and iTransformer [24] .The expected value of the MSE of the test data x ∈ R T ×T d concatenating the feature vector of all the tasks with the corresponding response variable y ∈ R T ×T q reads asR ∞ test = 1 T E[∥y -g(x)∥ 2 2 ] = 1 T E ∥ x ⊤ W √ T d + ε - x ⊤ AZQY T d ∥ 2 2 = 1 T E   ∥ x ⊤ W √ T d + ε - x ⊤ AZQ( Z ⊤ W √ T d + ε) T d ∥ 2 2   = 1 T E ∥ x ⊤ W √ T d + ε - x ⊤ AZQZ ⊤ W T d √ T d - x ⊤ AZQε T d ∥ 2 2 = 1 T E tr W ⊤ ΣW T d - 2tr W ⊤ ΣAZQZ ⊤ W (T d) 2 + tr ε ⊤ ε + tr W ⊤ ZQZ ⊤ AΣAZQZ ⊤ W (T d) 3 + tr ε ⊤ QZ ⊤ AΣAZQε (T d) 2 = 1 T E   tr W ⊤ ΣW T d - 2tr W ⊤ ΣA 1 2 (I T d -Q)A -1 2 W T d + tr ε ⊤ ε + tr W ⊤ ZQZ ⊤ AΣAZQZ ⊤ W (T d) 3 + tr ε ⊤ QZ ⊤ AΣAZQε (T d) 2 = tr W ⊤ ΣW T d -2 tr W ⊤ ΣA 1 2 (I T d -Q)A -1 2 W T d + tr ε ⊤ ε + tr W ⊤ ΣW T d -2 tr W ⊤ ΣA 1 2 QA -1 2 W T d + W ⊤ A -1 2 Q2 (A)A -1 2 W T d + 1 T d tr(Σ N Q2 ) + O 1 √ dThe test risk can be further simplified asR ∞ test = tr (Σ N ) + W ⊤ A -1 2 Q2 (A)A -1 2 W T d + tr Σ N Q2 T d + O 1 √ d D.The different theorems depends on the ground truth W that needs to be estimated through Ŵ.To estimate the test risk, one needs to estimate functionals of the form W ⊤ M Ŵ and ε ⊤ Mε for any matrix M. Using the expression of W = AZQY, we start computingŴ⊤ M Ŵ Ŵ⊤ MW = Y ⊤ QZ ⊤ AMAZQY Using the generative model for Y = Z ⊤ W √ T d + ε, we obtain E Ŵ⊤ MW = E Z ⊤ W √ T d + ε ⊤ QZ ⊤ AMAZQ Z ⊤ W √ T d + ε = 1 T d E W ⊤ ZQZ ⊤ AMAZQZ ⊤ W + E ε ⊤ QZ ⊤ AMAZQεEmploying the elementary relation A1 2 ZQ = QA 1 2Z, one obtains:E Ŵ⊤ MW = 1 T d E W ⊤ A -1 2 QA 1 2 Z ⊤ ZA 1 2 A 1 2 MA 1 2 QA 1 2 ZZ ⊤ W + E ε ⊤ QZ ⊤ AMAZQε = E W ⊤ A -1 2 I -Q A 1 2 MA 1 2 I -Q A -1 2 W + E ε ⊤ QZ ⊤ AMAZQε = E W ⊤ MW -2E W ⊤ MA 1 2 QA -1 2 W + E W ⊤ A -1 2 QA 1 2 MA 1 2 QA -1 2 W + E ε ⊤ QZ ⊤ AMAZQεUsing the deterministic equivalent of Lemma 1, we obtainŴ⊤ M Ŵ ↔ W ⊤ MW -2W ⊤ A -1 2 QA 1 2 MW + trΣ n M(A 1 2 MA 1 2 ) + W ⊤ A -1 2 Q2 (A 1 2 MA 1 2 )A -1 2 W ↔ W ⊤ M -2A -1 2 QA 1 2 M + A -1 2 Q2 (A 1 2 MA 1 2 )A -1 2 W + trΣ n Q2 (A 1 2 MA 1 2 ) ↔ W ⊤ κ(M)W + trΣ n Q2 (A 1 2 MA 1 2 )where We define the mapping κ : R T d×T d → R q×q as followsκ(M) = M -2A -1 2 QA 1 2 M + A -1 2 Q2 (A 1 2 MA 1 2 )A -1 2 .Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?Answer: [Yes]Justification:Guidelines:• The answer NA means that the abstract and introduction do not include the claims made in the paper. • The abstract and/or introduction should clearly state the claims made, including the contributions made in the paper and important assumptions and limitations. A No or NA answer to this question will not be perceived well by the reviewers. • The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings. • It is fine to include aspirational goals as motivation as long as it is clear that these goals are not attained by the paper. 
related work: 
 Usual Bounds in Multi-Task Learning. Since its introduction [1, 5, 26] , multi-task learning (MTL) has demonstrated that transfer learning concepts can effectively leverage shared information to solve related tasks simultaneously. Recent research has shifted focus from specific algorithms or minimization problems to deriving the lowest achievable risk given an MTL data model. The benefits of transfer learning are highlighted through risk bounds, using concepts such as contrasts [23] , which study high-dimensional linear regression with auxiliary samples characterized by the sparsity of their contrast vector or transfer distances [29] , which provide lower bounds for target generalization error based on the number of labeled source and target data and their similarities. Additionally, task-correlation matrices [30] have been used to capture the relationships inherent in the data model.However, these studies remain largely theoretical, lacking practical algorithms to achieve the proposed bounds. They provide broad estimates and orders of magnitude that, while useful for gauging the feasibility of certain objectives, do not offer precise performance evaluations. Moreover, these estimates do not facilitate the optimal tuning of hyperparameters within the MTL framework, a critical aspect of our study. We demonstrate the importance of precise performance evaluations and effective hyperparameter tuning to enhance the practical application of MTL.Random Matrix Theory in Multi-Task Learning. In the high-dimensional regime, Random Matrix Theory allows obtaining precise theoretical estimates on functionals (e.g., train or test risk) of data matrices with the number of samples comparable to data dimension [2, 13, 31, 46] . In the multi-task classification case, large-dimensional analysis has been performed for Least Squares SVM [48] and supervised PCA [47] showing counterproductivity of tackling tasks independently. We inspire our optimization strategy from [48] and conduct theoretical analysis for the multi-task regression, which introduces unique theoretical challenges and with a different intuition since we need to consider another data generation model and a regression learning objective. Our data modeling shares similarity with [52] which theoretically studied conditions of a negative transfer under covariate and model shifts. Our work focuses on the selection of hyperparameters within a general optimization framework, considering both a hyperparameter λ to relate all tasks and specific parameters γ t to balance overfitting and underfitting within each task, accommodating non-Gaussian distributions and employing recent developments such as deterministic equivalents and concentration inequalities.High dimensional analysis for regression. Regression is an important problem that has been extensively explored in the context of single-task learning for high-dimensional data, employing tools such as Random Matrix Theory [12] , physical statistics [9, 15] , and the Convex Gaussian MinMax Theorem (CGMT) [43] , among others. Typically, authors employ a linear signal plus noise model to calculate the asymptotic test risk as a function of the signal-generating parameter and the covariance of the noise. Our research builds upon these studies, extending them to a multi-task learning framework. In doing so, we derive several insights unique to multi-task learning. However, none of these previous studies offer a practical method for estimating the asymptotic test risk, which is dependent on the ground truth signal-generating parameter. Therefore, our work not only extends previous studies into the multi-task case within a generic data distribution (under the assumption of a concentrated random vector), but also provides a practical method for estimating these quantities.In this section, we derive the asymptotic risk for the training data.Question: Does the paper discuss the limitations of the work performed by the authors?Answer: [Yes]Justification:Guidelines:• The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper. • The authors are encouraged to create a separate "Limitations" section in their paper.• The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be. • The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on a few datasets or with a few runs. In general, empirical results often depend on implicit assumptions, which should be articulated. • The authors should reflect on the factors that influence the performance of the approach.For example, a facial recognition algorithm may perform poorly when image resolution is low or images are taken in low lighting. Or a speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon. • The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size. • If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness. • While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an important role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations. 
framework: 
 Notation. Throughout this study, matrices are represented by bold uppercase letters (e.g., matrix A), vectors by bold lowercase letters (e.g., vector v), and scalars by regular, non-bold typeface (e.g., scalar a). The notation A ⊗ B for matrices or vectors A, B is the Kronecker product. D x or Diag(x) stands for a diagonal matrix containing on its diagonal the elements of the vector x. The superscripts t and i are used to denote the task and the sample number, respectively, e.g., xi writes the i-th sample of the t-th task. The canonical vector of R T is denoted by e[T ] t with [e [T ] t ] i = δ ti . Given a matrix M ∈ R p×n , the Frobenius norm of M is denoted ∥M∥ F ≡ tr(M ⊤ M).For our theoretical analysis, we introduce the following notation of training data:Y = [Y (1) , . . . , Y (T ) ] ∈ R q×n , Z = T t=1 e [T ] t e [T ] t ⊤ ⊗ X (t) ∈ R T d×nwhere n = T t=1 n t is the total number of samples in all the tasks.We define V = [V ⊤ 1 , . . . , V ⊤ T ] ⊤ ∈ R dT ×q and propose to solve the linear multi-task problem by finding Ŵ = [ Ŵ⊤ 1 , . . . , Ŵ⊤ T ] ⊤ ∈ R dT ×q which solves the following optimization problem under the assumption of relatedness between the tasks, i.e., W t = W 0 + V t for all tasks t : minEQUATIONwhereJ (W 0 , V) ≡ 1 2λ ∥W 0 ∥ 2 F + 1 2 T t=1 ∥V t ∥ 2 F γ t + 1 2 T t=1 Y (t) - X (t) ⊤ W t √ T d 2 F .The objective function consists of three components: a regularization term for W 0 to mitigate overfitting, a task-specific regularization term that controls deviations V t from the shared weight matrix W 0 , and a loss term quantifying the error between the predicted outputs and the actual responses for each task.The optimization problem is convex for any positive values of λ, γ 1 , . . . , γ T , and it possesses a unique solution. The details of the calculus are given in Appendix B. The formula for Ŵt is given as follows:Ŵt = e [T ] t ⊤ ⊗ I d AZQY √ T d , Ŵ0 = 1 ⊤ T ⊗ λI d ZQY √ T d ,where withQ = Z ⊤ AZ T d + I n -1∈ R n×n , andA = D γ + λ1 T 1 ⊤ T ⊗ I d ∈ R T d×T d .In order to use Random Matrix Theory (RMT) tools, we make two assumptions on the data distribution and the asymptotic regime. Following [50] , we adopt a concentration hypothesis on the feature vectorsx (t)i , which was shown to be highly effective for analyzing machine learning problems [11, 14] . The justification and implications of these assumptions can be found in appendix A.2 and A.3. Assumption 1 (Concentrated Random Vector). We assume that there exists two constants C, c > 0 (independent of dimension d) such that, for any task t, for any 1-Lipschitz function f : R d → R, any feature vector x (t) ∈ X (t) verifies:∀t > 0 : P(|f (x (t) ) -E[f (x (t) )]| ≥ t) ≤ Ce -(t/c) 2 , E[x (t) ] = 0 and Cov[x (t) ] = Σ (t) .In particular, we distinguish the following scenarios:x (t)i ∈ R d are concentrated when they are (i) independent Gaussian random vectors with covariance of bounded norm, (ii) independent random vectors uniformly distributed on the R d sphere of radius √ d, and most importantly (iii) any Lipschitz transformation ϕ(x (t) i ) of the above two cases, with bounded Lipschitz norm. Scenario (iii) is especially pertinent for modeling data in realistic settings. Recent research [39] has demonstrated that images produced by generative adversarial networks (GANs) are inherently qualified as concentrated random vectors.Next, we present a classical RMT assumption that establishes a commensurable relationship between the number of samples and dimension. Assumption 2 (High-dimensional asymptotics). As d → ∞, n t = O(d) and T = O(1). More specifically, we assume that n/d a.s.--→ c 0 < ∞ with n = T t=1 n t .Although different from classical asymptotic where the number of samples is implicitly assumed to be exponentially larger than the dimension, the high-dimensional asymptotic finds many applications including telecommunications [10] , finance [34] and machine learning [11, 14, 48] .Question: For each theoretical result, does the paper provide the full set of assumptions and a complete (and correct) proof?Answer: [Yes]Answer: [Yes] Justification: The datasets used are open-source and all the implementation details are given to reproduce the experimental results. The code will be made available in case of acceptance.Guidelines:• The answer NA means that paper does not include experiments requiring code. • Providing as much information as possible in supplemental material (appended to the paper) is recommended, but including URLs to data and code is permitted. 
estimation of the performances: 
 Given the training dataset X ∈ R n×d with the corresponding response variable Y ∈ R n×q and for any test dataset x (t) ∈ R d and the corresponding response variable y (t) ∈ R q , we aim to derive the theoretical training and test risks defined as follows:R ∞ train = 1 T n E ∥Y -g(X)∥ 2 2 , R ∞ test = 1 T T t=1 E[∥y (t) -g(x (t) )∥ 2 2 ]The output score g(x (t) ) ∈ R q for data x (t) ∈ R d of task t is given by:EQUATIONIn particular, the output for the training score is given by:EQUATIONTo understand the statistical properties of R ∞ train and R ∞ test for the linear generative model, we study the statistical behavior of the resolvent matrix Q defined in (6) . The notion of a deterministic equivalent from random matrix theory [16] is particularly useful here as it allows us to design a deterministic matrix equivalent to Q in the sense of linear forms. Specifically, a deterministic equivalent M of a random matrix M is a deterministic matrix that satisfies u(M -M) → 0 almost surely for any bounded linear form u : R d×d → R. We denote this as M ↔ M. This concept is crucial to our analysis because we need to estimate quantities such as 1 d tr (AQ) or a ⊤ Qb, which are bounded linear forms of Q with A, a and b all having bounded norms. For convenience, we further express some contributions of the performances with the so-called "coresolvent" defined asQ ≡ ( A 1 2 ZZ T A 1 2 T d + I T d ) -1 .Using Lemma 1 provided in the Appendix C.1, whose proofs are included in Appendices C.2, C.3 and C.4, we establish the deterministic equivalents that allow us to introduce our Theorems 5 and 1, characterizing the asymptotic behavior of both training and testing risks.Theorem 1 (Asymptotic train and test risk). Assuming that the training data vectors x (t) i and the test data vectors x (t) are concentrated random vectors, and given the growth rate assumption (Assumption 2), it follows that:R ∞ test = tr W ⊤ A -1 2 Q2 (A)A -1 2 W T d signal term + tr(Σ n Q2 ) T d + tr (Σ n ) noise terms . (ATR)In addition, the asymptotic risk on the training data is given byR ∞ train ↔ tr W ⊤ A -1/2 QA -1/2 W T n - tr W ⊤ A -1/2 Q2 (I T d )A -1/2 W T n + tr Σ n Q2 T n ,where Q, Q2 (I T d ) and Q2 are respectively deterministic equivalents for Q, Q2 and Q 2 .We defer the full proof of this theorem to Appendix D and provide an outline of the test risk proof below. Note that derivations of the asymptotic risk for the training data are more complex and involve computing deterministic equivalents (Lemma 1 of Appendix C.1), which is a powerful tool commonly used in Random Matrix Theory. In the derived theorem, we achieve a stronger convergence result compared to the almost surely convergence. Specifically, we prove a concentration result for the training and test risk with an exponential convergence rate for sufficiently large values of n and d.R ∞ test = 1 T E[∥y -g(x)∥ 2 2 ] = 1 T E ∥ x ⊤ W √ T d + ε - x ⊤ AZQZ ⊤ W T d √ T d - x ⊤ AZQε T d ∥ 2 2 = 1 T E tr W ⊤ ΣW T d - 2tr W ⊤ ΣAZQZ ⊤ W (T d) 2 + tr ε ⊤ ε + tr W ⊤ ZQZ ⊤ AΣAZQZ ⊤ W (T d) 3 + tr ε ⊤ QZ ⊤ AΣAZQε (T d) 2 = tr W ⊤ ΣW T d -2 tr W ⊤ ΣA 1 2 (I T d -Q)A -1 2 W T d + tr ε ⊤ ε + tr W ⊤ ΣW T d -2 tr W ⊤ ΣA 1 2 QA -1 2 W T d + W ⊤ A -1 2 Q2 (A)A -1 2 W T d + 1 T d tr(Σ N Q2 ) + O 1 √ dTo gain theoretical insights, we analyze (ATR) consisting of the signal and the noise components.Signal Term. The signal term can be further approximated, up to some constants as tr(W t) . The matrix (AΣ + I) -2 plays a crucial role in amplifying the signal term tr(W ⊤ W), which in turn allows the test risk to decrease. The off-diagonal elements of AΣ + I) -2 amplify the cross terms (tr(W ⊤ v W t ) for t ̸ = v), enhancing the multi-task aspect, while the diagonal elements amplify the independent terms (∥W t ∥ 2⊤ (AΣ + I) -2 W) with Σ = T t=1 nt d Σ (2 ). This structure is significant in determining the effectiveness of multi-task learning.Furthermore, both terms decrease with an increasing number of samples n t , smaller values of γ t , and a larger value of λ. The cross term, which is crucial for multi-task learning, depends on the matrix Σ -1 t Σ v . This matrix represents the shift in covariates between tasks. When the features are aligned (i.e., Σ -1 t Σ v = I d ), the cross term is maximized, enhancing multi-task learning. However, a larger Fisher distance between the covariances of the tasks results in less favorable correlations for multi-task learning.Noise term. Similar to the signal term, the noise term can be approximated, up to some constants, as tr(Σ N (A -1 + Σ) -1 ). However, there is a major difference between the way both terms are expressed in the test risk. The noise term does not include any cross terms because the noise across different tasks is independent. In this context, only the diagonal elements of the matrix are significant. This diagonal term increases with the sample size and the value of λ. It is responsible for what is known as negative transfer. As the diagonal term increases, it negatively affects the transfer of learning from one task to another. This is a critical aspect to consider in multi-task learning scenarios.In this section, we specialize the theoretical analysis to the simple case of two tasks (T = 2). We assume that the tasks share the same identity covariance and that γ 1 = γ 2 ≡ γ. Under these conditions, the test risk can be approximated, up to some constants, asR ∞ test = D IL ∥W 1 ∥ 2 2 + ∥W 2 ∥ 2 2 + C M T L W ⊤ 1 W 2 + N N TtrΣ n where the diagonal term (independent learning) D IL , the cross term (multi-task learning) C M T L , and the noise term (negative transfer) N N T have closed-form expressions depending on γ and λ: We recall that c 0 has been defined in the Assumption 2. As previously mentioned, the test risk is primarily composed of two terms: the signal term and the noise term, which are in competition with each other. The more similar the tasks are, the stronger the signal term becomes. In the following plot, we illustrate how this competition can influence the risk. Depending on the value of the parameter λ and the sample sizes, the risk can either decrease monotonically, increase monotonically, or exhibit a convex behavior. This competition can lead to an optimal value for λ, which interestingly has a simple closed-form expression that can be obtained by deriving the R ∞ test w.r.t. λ as follows (see details in Appendix E.4):D IL = (c 0 (λ + γ) + 1) 2 + c 2 0 λ 2 (c 0 (λ + γ) + 1) 2 -c 2 0 λ 2 , C M T L = -2c 0 λ(c 0 (λ + γ) + 1) (c 0 (λ + γ) + 1) 2 -c 2 0 λ 2 N N T = (c 0 (λ + γ) 2 + (λ + γ) -c 0 λ 2 ) 2 + λ 2 ((c 0 (λ + γ) + 1) 2 -c 2 0 λ 2 ) (a) n d = 0.5 (b) n d = 1.5 (c) n d = 2.5λ ⋆ = n d SNR - γ 2 , with SNR = ∥W 1 ∥ 2 2 + ∥W 2 ∥ 2 2 trΣ N + W ⊤ 1 W 2 trΣ N .In this section, we compare the theoretical predictions with the empirical results. Our experiment is based on a two-task setting (T = 2) defined as W 1 ∼ N (0, I p ) withW 2 = αW 1 + √ 1 -α 2 W ⊥ 1 . W ⊥1 represents any vector orthogonal to W 1 and α ∈ [0, 1]. This setting allows us to adjust the similarity between tasks through α.Figure 2 shows a comparison of the theoretical and empirical classification errors for different values of λ, highlighting the error-minimizing value of λ. Despite the relatively small values of n and p, there is a very precise match between the asymptotic theory and the practical experiment. This is particularly evident in the accurate estimation of the optimal value for λ. T d tr (Σ N M). These quantities can be consistently estimated under Assumptions 2 as follows (see details in Section F of the Appendix) :1 T d tr W ⊤ MW -ζ(M) a.s. --→ 0, 1 T d tr (Σ N M) -σtrM a.s.--→ 0whereEQUATION2 ).We define the estimate of the noise as σ = lim λ→0 γt→∞ R ∞ train and the function κ -1 is the functional inverse of the mapping κ : R T d×T d → R q×q defined as followsκ(M) = M -2A -1 2 QA 1 2 M + A -1 2 Q2 (A 1 2 MA 1 2 )A -1 2 - n (T d) 2 A -1 2 ΣA -1 2 Q2 (A 1 2 MA 1 2 ).The test risk in the particular of identity covariance matrix can be rewritten asR ∞ test = D IL ∥W 1 ∥ 2 2 + ∥W 2 ∥ 2 2 + C M T L W ⊤ 1 W 2 + N N T trΣ n . Deriving R ∞test with respect to λ leads after some algebraic calculus toλ ⋆ = n d SNR - γ 2where the signal noise ratio is composed of the independent signal to noise ratio and the cross signal to noise ratio SNR =∥W1∥ 2 2 +W2∥ 2 2 trΣn + W ⊤ 1 W2 trΣn F Theoretical Estimations F. 
application to multi-task regression: 
 In our study, we apply the theoretical framework presented in our paper to a real-world regression problem, specifically, the Appliance Energy dataset which aims to predict the total usage of a house. This dataset is a multivariate regression dataset containing 138 time series each of dimension 24. We select two of these features as 2 tasks to cast the problem as a multi-task learning regression problem.Figure 3 presents a comparison between the theoretical predictions and empirical outcomes of our proposed linear framework. Despite the assumptions made in the main body of the paper, the theoretical predictions closely align with the experimental results. This demonstrates that our estimates are effective in practice and provide a reliable approximation of the optimal regularization.In essence, our theoretical framework, when applied to real-world multi-task learning regression problems, yields practical and accurate estimates, thereby validating its effectiveness and applicability.While non-linear models are widely used, establishing their theoretical foundations is challenging. Therefore, we focused on linear models, which, despite their simplicity, provide valuable insights into more complex models.Our results show that test risk curves for non-linear models follow patterns predicted by our theory. This is expected because non-linear models in time series forecasting typically use a linear output layer for prediction. Thus, we can apply our theory to the inputs of this final linear layer. This approach is valid due to data concentration and the Lipschitz nature of neural networks, ensuring outputs of the non-linear part don't deviate significantly from the inputs.Moreover, multivariate time series models often treat channels separately using univariate methods, missing cross-channel information. Our results in Section 5.4 show that our method surpasses univariate baselines by optimally regularizing with λ and γ, supporting our theory's applicability to non-linear models as the final linear layer effectively leverages concentrated inputs.Finally, our regularization approach differs from traditional cross-task regularizations that use one task per dataset. We consider each prediction as a task and introduce γ t parameters alongside λ. These parameters enforce multivariate regularization and control underfitting or overfitting per task. This method is tractable since it's applied at the model's final layer.The similarity between curves for non-linear and linear models indicates our findings are robust; nonlinear models also exhibit optimal regularization parameters, enhancing performance in multivariate forecasting.Our theoretical framework is applied in the context of Multivariate Time Series Forecasting, with related work detailed in Appendix G.1. We previously applied this framework in a linear setting, and now aim to evaluate its empirical validity in the non-linear setting of neural networks. The results presented in this section represent the best test MSE, assuming the ability to find the optimal lambda value, which can be considered as an oracle scenario. A study of these limitations can be found in Appendix H.Motivation. Our approach is applied to the MTSF setting for several reasons. Firstly, many models currently used are essentially univariate, where predictions for individual series are simply concatenated without exploiting the multivariate information inherent in traditional benchmarks.Given that these benchmarks are designed for multivariate forecasting, leveraging multivariate information should yield better results. Secondly, our theoretical framework can benefit this domain, as most predictive models use a linear layer on top of the model to project historical data of length d for predicting the output of length q. This characteristic aligns well with our method, making it a promising fit for enhancing forecasting accuracy.Our approach. We propose a novel method for MTSF by modifying the loss function to incorporate both individual feature transformations f t and a shared transformation f 0 . Each univariate-specific transformation f t is designed to capture the unique dynamics of its respective feature, while f 0 serves as a common transformation applied across all features to capture underlying patterns shared among them. We consider a neural network f with inputs X = [X (1) , X (2) , . . . , X (T ) ], where T is the number of channels and X (t) ∈ R n×d . For a univariate model without MTL regularization, we predict Y = [Y (1) , Y (2) , . . . , Y (T ) ] = [f 1 (X (1) ), . . . , f T (X (T ) )] and Y (t) ∈ R n×q . We compare these models with their corresponding versions that include MTL regularization, formulated as:f M T L t (X (t) ) = f t (X (t) ) + f 0 (Y)with f t : R n×d → R n×q and f 0 : R n×qT → R n×q . We define our regularized loss as follows:L(X, Y) = T t=1 ∥Y (t) -f M T L t (X (t) )∥ 2 F + λ∥f 0 (X)∥ 2 F + T t=1 γ t ∥f t (X (t) ∥ 2 F , ∀t ∈ {1, . . . , T }.where Y (t) are the true predictions, f t represents the univariate model for each channel t, and λ is our regularization parameter, for which we have established a closed form in the case of linear f t . f 0 serves a role equivalent to W 0 , which was defined in our theoretical study and allows for the regularization of the common part. This component can be added at the top of a univariate model. The parameters γ t enable the regularization of the specialized parts f t .In our setup, f t is computed in a similar way as in the model without regularization and f 0 is computed by first flattening the concatenation of the predictions of X (t) , then applying a linear projection leveraging common multivariate information before reshaping. The loss function is specifically designed to balance fitting the multivariate series using f 0 and the specific channels using f t . This Finally, MTL regularization enhances the performance of univariate models, making them often competitive with state-of-the-art multivariate methods like SAMformer and iTransformer. This approach seems to better captures shared dynamics among tasks, leading to more accurate forecasts. 
conclusions and future works: 
 In this article, we have explored linear multi-task learning by deriving a closed-form solution for an optimization problem that capitalizes on the wealth of information available across multiple tasks.Leveraging Random Matrix Theory, we have been able to obtain the asymptotic training and testing risks, and have proposed several insights into high-dimensional multi-task learning regression. Our theoretical analysis, though based on a simplified model, has been effectively applied to multi-task regression and multivariate forecasting, using both synthetic and real-world datasets. We believe that our work lays a solid foundation for future research, paving the way for using random matrix theory with more complex models, such as deep neural networks, within the multi-task learning framework. We would like to note that we are assuming that both the noise ε and the feature x(t)i have zero mean. This is a common assumption in many statistical models and it simplifies the analysis. However, this assumption is not restrictive. In practice, if the data or the response variable are not centered, we can always preprocess the data by subtracting the mean. This preprocessing step brings us back to the zero-mean setting that we consider in our theoretical analysis.Question: Does the paper specify all the training and test details (e.g., data splits, hyperparameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?Answer: [Yes]Justification:Guidelines:• The answer NA means that the paper does not include experiments.• The experimental setting should be presented in the core of the paper to a level of detail that is necessary to appreciate the results and make sense of them. • The full details can be provided either with the code, in appendix, or as supplemental material. 
experiment statistical significance: 
 Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?Answer: [Yes]Justification:Guidelines:• The answer NA means that the paper does not include experiments.• The authors should answer "Yes" if the results are accompanied by error bars, confidence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper. • The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions). • The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.) • The assumptions made should be given (e.g., Normally distributed errors).• It should be clear whether the error bar is the standard deviation or the standard error of the mean. • It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have a 96% CI, if the hypothesis of Normality of errors is not verified. • For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out of range (e.g. negative error rates). Learning. There are many potential societal consequences of our work, none which we feel must be specifically highlighted here. Guidelines:• The answer NA means that there is no societal impact of the work performed.• If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact. • Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.• The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster. • The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology. • If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML). 
safeguards: 
 Question: Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?Answer: [NA]Justification:Guidelines:• The answer NA means that the paper poses no such risks.• Released models that have a high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters. • Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images. • We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make a best faith effort.12. Licenses for existing assets Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?Answer: [Yes] Justification:Guidelines:• The answer NA means that the paper does not use existing assets.• The authors should cite the original paper that produced the code package or dataset.• The authors should state which version of the asset is used and, if possible, include a URL. • The name of the license (e.g., CC-BY 4.0) should be included for each asset.• For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided. • If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswithcode.com/datasets has curated licenses for some datasets. Their licensing guide can help determine the license of a dataset. • For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided. 
 
 Appendix: Data are concentrated random vectors, meaning high-dimensional data maintain stable under complex (Lipschitz) transformations. The strong performance of neural networks on tasks like image recognition and NLP suggests that these models produce stable predictions. As Lipschitz transformations, they maintain controlled distances between inputs, ensuring stability. Recent studies have demonstrated and experimentally confirmed that both real-world data and synthetically generated data using GANs exhibit concentration properties, supporting this assumption. This makes our assumption more realistic than traditional Gaussian assumptions, as it does not rely on specific hypotheses about the shape of the data distribution, but rather on the stability of statistical properties after transformation. Consequently, analyzing a framework of concentrated random vectors is more theoretically challenging than using Gaussian assumptions and represents a key novelty of our theory.The dimension d is of the same order of magnitude as the sample size n. This joint growth captures data complexity better than assuming a fixed feature size with increasing samples, which can oversimplify models. Our theory works for fixed d and n, unbiased by specific parameter choices. The accuracy of empirical predictions depends on both d and n, with variance scaling as O 1 √ dn . Larger d and n reduce variance, making empirical results more reliable and closer to theoretical values. Conversely, smaller d and n increase variance, affecting single predictions. However, this scaling is still better than n growing indefinitely with fixed d, where variance scales as O 1 √ n , leading to increased bias. Implication of the assumptions. In Section 5.2, Figure 3 , the Appliance Energy dataset illustrates our assumptions' realism. Despite the moderate sample size of 42 samples with 142 dimensions and non-synthetic data, the theoretical curve fits the empirical predictions well. Additionally, with synthetic data (d = 100, n = 100), the theory matches the empirical curve across various hyperparameters (Figure 2 ), showcasing the predictive power of Random Matrix Theory.The proposed multi task regression finds Ŵ = [ Ŵ⊤ 1 , . . . , Ŵ⊤ k ] ⊤ ∈ R dT ×q which solves the following optimization problem using the additional assumption of relatedness between the tasks (W t = W 0 + V t for all tasks t):EQUATIONwhereJ (W 0 , V) ≡ 1 2λ tr W ⊤ 0 W 0 + 1 2 T t=1 tr V ⊤ t V t γ t + 1 2 T t=1 tr ξ ⊤ t ξ t ξ t = Y (t) - X (t) ⊤ W t √ T d , ∀t ∈ {1, . . . , T }.The Lagrangian introducing the lagrangian parameters for each task t, α t ∈ R nt×q reads asL(W 0 , V t , ξ t , α t ) = 1 2λ tr W ⊤ 0 W 0 + 1 2 T t=1 tr V ⊤ t V t γ t + 1 2 t=1 tr ξ ⊤ t ξ t + T t=1 tr α ⊤ t Y (t) - X (t) ⊤ (W 0 + V t ) √ T d -ξ tDifferentiating with respect to the unknown variables Ŵ0 , Vt , ξ t , α t and b t , we get the following system of equation1 λ Ŵ0 - T t=1 X (t) α t √ T d = 0 1 γ t Vt - X (t) α t √ T d = 0 ξ t -α t = 0 Y (t) - X (t) ⊤ Ŵ0 √ T d - X (t) ⊤ Vt √ T d -ξ t = 0Plugging the expression of Ŵ0 , Vt and ξ t into the expression of Y (t) givesY (t) = λ T t=1 X (t) ⊤ X (t) T d α t + γ t X (t) ⊤ X (t) T d α t + α twhich can be rewritten asY (t) = (λ + γ t ) X (t) ⊤ X (t) T d α t + λ v̸ =t X (t) ⊤ X (v) T d α v + α t With Y = [Y (1) ⊤ , . . . , Y (T ) ⊤ ] ⊤ ∈ R n×q , α = [α ⊤ 1 , . . . , α ⊤ k ] ⊤ ∈ R n×q , Z = T t=1 e [T ]t e[T ] t ⊤ ⊗ X (t) ∈ R T d×n , this system of equations can be written under the following compact matrix form:Q -1 α = Y with Q = Z ⊤ AZ T d + I n -1∈ R n×n , andA = D γ + λ1 T 1 ⊤ T ⊗ I d ∈ R T d×T d .Solving for α then gives:α = QYMoreover, using Ŵt = Ŵ0 + Vt , the expression of W t becomes:Ŵt = e [T ] t ⊤ ⊗ I d AZα √ T d , Ŵ0 = 1 ⊤ T ⊗ λI d Zα √ T d .C Lemma 1 and proof with Random Matrix Theory C.1 Lemma 1 Lemma 1 (Deterministic equivalents for Q, QM Q and Q 2 for any M ∈ R n×n ). Under the concentrated random vector assumption for each feature vector x (t)i and under the growth rate assumption (Assumption 2), for any deterministic M ∈ R n×n , we have the following convergence:Q ↔ Q, QM Q ↔ Q2 (M), Q 2 ↔ Q2where Q2 , Q and Q2 are defined as followsQ = T t=1 c 0 C (t) 1 + δ t + I T d -1 , δ t = 1 T d tr Σ (t) Q , C (t) = A 1 2 e [T ] t ⊗ Σ (t) A 1 2 Q2 (M) = QM Q + 1 T d T t=1 d t 1 + δ t QC (t) Q, d = I T - 1 T d Ψ -1 Ψ(M) ∈ R T Q2 = I n -Diag t∈[T ] (v t I nt ), v t = 1 T d tr(C (t) Q) (1 + δ t ) 2 + 1 T d tr C (t) Q2 (I n ) (1 + δ t ) 2whereΨ(M ) =   n t T d tr C (t) QM Q 1 + δ t   t∈[T ] ∈ R T , Ψ =   n t T d tr C (t) QC (t ′ ) Q (1 + δ t )(1 + δ t ′ )   t,t ′ ∈[T ] ∈ R T ×T , C.2 Deterministic equivalent of the resolvent QThe evaluation of the expectation of linear forms on Q and Q2 can be found in the literature. To find a result that meets exactly our setting, we will cite [25] that is a bit more general since it treats cases where E[xi ] ̸ = 0 for t ∈ [T ] and i ∈ [n t ]. Unlike the main paper, and to be more general, the study presented below is "quasi asymptotic" meaning that the results are true for finite value of d, n. Let us first rewrite the general required hypotheses, adapting them to our setting. For that purpose, we consider in the rest of this paper a certain asymptoticI ⊂ {(d, n), d ∈ N, n ∈ N} = N 2 satisfying: {d, ∃ n ∈ N : (d, n) ∈ I} = N and {n, ∃ d ∈ N : (d, n) ∈ I} = N.such that n and d can tend to ∞ but with some constraint that is given in the first item of Assumption 3 below. Given two sequences (a d,n) d,n∈I , (b d,n ) d,n∈I > 0, the notation a d,n ≤ O(b d,n ) (or a ≤ O(b)) means that there exists a constant C > 0 such that for all (d, n) ∈ I, a d,n ≤ Cb d,n .Assumption 3. There exists some constants C, c > 0 independent such that:• n ≤ O(d) • Z = (z 1 , . . . , z n ) ∈ R T d×n has independent columns• for any (d, n) ∈ I, and any f : R T d×n → R 1-Lipschitz for the euclidean norm:P (|f (Z) -E[f (Z)]| ≥ t) ≤ Ce -ct 2 . • ∀i ∈ {n, ∃d ∈ N, (d, n) ∈ I}: ∥E[z i ]∥ ≤ O(1).Theorem 2 ([25], Theorem 0.9.). Given T ∈ N, Z ∈ R T d×n and two deterministicA ∈ R T d×T d , we note Q ≡ ( 1 T d A 1 2 ZZ ⊤ A 1 2 + I T d ) -1 .If Z satisfies Assumption 3 and M ∈ R T d×T d is a deterministic matrix satisfying ∥M∥ F ≤ 1, one has the concentration:P tr(M Q) -tr(M Qδ(S) (S)) ≥ t ≤ Ce -ct 2 ,where S = (S 1 , . . . ,S n ) = (E[z 1 z ⊤ 1 ], . . . , E[z n z ⊤ n ]), for δ ∈ R n , Qδ is defined as:Qδ (S) =   1 T d i∈[n] A 1 2 S i A 1 2 1 + δ i + I T d   -1, and δ(S) is the unique solution to the system of equations:∀i ∈ [n] : δ(S) i = 1 n tr A 1 2 S i A 1 2 Qδ(S) .We end this subsection with some results that will be useful for next subsection on the estimation of bilinear forms on Q. Lemma 2 ([25], Lemmas 4.2, 4.6). Under the setting of Theorem 2, given a deterministic vector u ∈ R T d such that ∥u∥ ≤ O(1) and two deterministic matrices U, V such that ∥U∥, ∥V∥ ≤ O(1) and a power r > 0, r ≤ O(1):• E u ⊤ U Q-i Vz i r ≤ O(1) • E 1 T d z ⊤ i U Q-i Vz i -E 1 T d tr Σ i U QB r ≤ O 1 d r 2.To simplify the expression of the following theorem, we take A = I T d . One can replace Z with A 1 2 Z to retrieve the result necessary for the main paper. Theorem 3. Under the setting of Theorem 2, with A = I T d , one can estimate for any deterministic matrices U, V ∈ R T d such that ∥U∥, ∥V∥ ≤ O(1) and any deterministic vectorEQUATION)where we noted:• Q ≡ Qδ (S), δ = δ(S), • Ψ ≡ 1 T d tr Si QSj Q (1+δi)(1+δj ) i,j∈[n] ∈ R n,n • ∀U ∈ R n×n : Ψ(U) ≡ 1 T d tr U QSi Q 1+δi i∈[n] ∈ R n • ∀U, V ∈ R n×n : Ψ(U, V) ≡ 1 T d tr U QV Q ∈ R If there exist T < n dinstinct matrices C 1 , . . . , C T such that: {S 1 , . . . , S n } = {C 1 , . . . , C T } , and if we denote ∀t ∈ [T ] n t = #{i ∈ [n] | S i = C t } and: P ≡   IT -   n t n v (T d) 2 tr S t QS v Q (1 + δ t )(1 + δ v )   t,v∈[T ]    -1 ∈ R T,T ∀U ∈ R T d×T d : Q2 (U) ≡ QU Q + 1 (T d) 2 T t,v=1 tr(S t QU Q)P t,v QS v Q (1 + δ t )(1 + δ v ) ,the result of Theorem 3 rewrites:E QU Q -Q2 (U) ≤ O 1 √ d (9) Proof. Given i ∈ [n], let us note Z -i = (z 1 , . . . , z i-1 , 0, z i+1 , . . . , z n ) and Q-i = ( 1 T d Z -i Z ⊤ -i + I T d ) -1, then we have the identity:EQUATIONGiven u, v ∈ R T d , such that ∥u∥, ∥v∥ ≤ 1, let us express:EQUATION(12) First, given i ∈ [n], let us estimate thanks to (10) :E u ⊤ QS i QU Qv = E u ⊤ Q-i S i QU Qv - 1 T d E u ⊤ Qz i z ⊤ i Q-i S i QU QvHölder inequality combined with Lemma 2 allows us to bound:1 T d E u ⊤ Qz i z ⊤ i Q-i S i QU Qv ≤ 1 T d E u ⊤ Qz i 2 1 2 E z ⊤ i Q-i S i QU Qv 2 1 2 ≤ O 1 d ,one can thus deduce:EQUATIONSecond, one can also estimate thanks to Lemma 10:E u ⊤ Qz i z ⊤ i QU Qv = E u ⊤ Q-i z i z ⊤ i QU Qv 1 + 1 T d z ⊤ i Q -i z i = E u ⊤ Q-i z i z ⊤ i QU Qv 1 + δ i + O 1 √ d ,again thanks to Hölder inequality combined with Lemma 2 that allow us to bound:E δ i -1 T d z ⊤ i Q -i z i (1 + δ i ) 1 + 1 T d z ⊤ i Q -i z i |u ⊤ Q-i z i z ⊤ i QU Qv| ≤ E δ i - 1 T d z ⊤ i Q -i z i 2 1 2 E |u ⊤ Q-i z i z ⊤ i QU Qv| 2 1 2 ≤ O 1 √ d ,The independence between z i and Q-i (and Q) then allow us to deduce (again with formula (10) ):EQUATIONLet us inject ( 13) and ( 14) in (11) to obtain (again with an application of Hölder inequality and Lemma 2 that we do not detail this time):E u ⊤ Q S i 1 + δ i -z i z ⊤ i QU Qv = 1 T d E    u ⊤ Q-i z i z ⊤ i QU Qz i z ⊤ i Q-i v 1 + 1 n z ⊤ i Q-i z i 2    + O 1 √ d , = 1 T d E u ⊤ Q-i S i Q-i v (1 + δ i ) 2 tr S i QU Q + O 1 √ d ,Putting all the estimations together, one finally obtains:EQUATIONOne then see that if we introduce for any V ∈ R n×n the block matrices:• θ = 1 T d ( E[tr(Sj QSi QY )] (1+δi)(1+δj ) ) i,j∈[n] ∈ R n×n • θ(V) = 1 T d ( E[tr(V QSi QY )] 1+δi ) i∈[n] ∈ R n , • θ(U, V) = 1 T d E tr(V QU QY ) ∈ R, then, if ∥V∥ ≤ O(1), multiplying (15) with V and taking the trace leads to:EQUATIONNow, taking U = S1 1+δ1 , . . . , Sn 1+δn , one gets the vectorial equation:θ(V) = Ψ(V) + 1 T d Ψθ(V) + O 1 √ d , When (I T d -1 T d Ψ) is invertible, one gets θ(V) = (I T d -1 T d Ψ) -1 Ψ(V) + O 1 √ d, and combining with ( 16), one finally obtains:θ(U, V) = Ψ(U, V) + 1 T d Ψ(U) ⊤ (I T d - 1 T d Ψ) -1 Ψ(V) + O 1 √ d .C.4 Estimation of the deterministic equivalent of Q 2 Theorem 4. Under the setting of Theorem 3, one can estimate:EQUATIONwith, ∀i ∈ [n]:v i ≡ 1 T d tr S i Q (1 + δ i ) 2 + 1 T d tr S i Q2 (I n ) (1 + δ i ) 2Proof. The justifications are generally the same as in the proof of Theorem 3, we will thus allow ourselves to be quicker in this proof.Using the definition ofQ = Z ⊤ AZ T d + I n -1, we have thatEQUATIONand one can then let appear Q thanks to the relation:EQUATIONthat finally gives us:Q = I n - 1 T d Z ⊤ ZQ = I n - 1 T d Z ⊤ QZOne can then express:Q 2 = I n - 2 T d Z ⊤ QZ + 1 (T d) 2 Z ⊤ QZZ ⊤ QZ = I n - 1 T d Z ⊤ QZ - 1 T d Z ⊤ Q2 Z.Given i, j ∈ [n], i ̸ = j, let us first estimate (thanks to Hölder inequality and Lemma 2):1 T d E z ⊤ i Qz j = 1 T d E z ⊤ i Q-i,j z j (1 + δ i )(1 + δ j ) + O 1 √ d ≤ O 1 √ d , since E[z i ] = E[z j ] = 0. Now, we consider the case j = i to get:1 T d E z ⊤ i Qz i = 1 T d E z ⊤ i Q-i z i (1 + δ i ) 2 + O 1 √ d = 1 T d tr S i Q (1 + δ i ) 2 + O 1 √ d .As before, we know that 1T d E z ⊤ i Qz j ≤ O 1 √ d if i ̸ = j. Considering i ∈ [n],we thus are left to estimate:1 T d E z ⊤ i Q2 z j = 1 T d tr S i Q2 (I n ) (1 + δ i ) 2 + O 1 √ d D Risk Estimation (Proof of Theorem 1) D.i and the test data vectors x (t) are concentrated random vectors, and given the growth rate assumption (Assumption 2), it follows that:R ∞ train ↔ 1 T n tr W ⊤ A -1/2 QA -1/2 W - 1 T n tr W ⊤ A -1/2 Q2 (I T d )A -1/2 W + 1 T n tr Σ N Q2Proof. We aim in this setting of regression, to compute the asymptotic theoretical training risk given by:R ∞ train = 1 T n E Y - Z ⊤ AZ T d QY 2 2Using the definition ofQ = Z ⊤ AZ T d + I T d -1, we have thatZ ⊤ AZ T d Q = Z ⊤ AZ T d + I T d -I T d Z ⊤ AZ T d + I T d -1 = I T d -QPlugging back into the expression of the training risk then leads toR ∞ train = 1 T n E tr Y ⊤ Q 2 YUsing the definition of the linear generative model and in particularY = Z ⊤ W √ T d + ε, we get R ∞ train = 1 T n E tr 1 √ T d Z ⊤ W + ε ⊤ Q 2 1 √ T d Z ⊤ W + ε = 1 T n 1 T d E tr W ⊤ ZQ 2 Z ⊤ W + 1 T n E tr ε ⊤ Q 2 εTo simplify this expression, we will introduced the so-called "coresolvent" defined as:Q = A 1 2 ZZ ⊤ A 1 2 T d + I T d -1, Employing the elementary relation A1 2 ZQ = QA 1 2 Z, one obtains: 1 T d ZQ 2 Z ⊤ = 1 T d A -1 2 QA 1 2 ZQZ ⊤ = A -1 2 Q2 A 1 2 ZZ ⊤ A 1 2 T d A -1 2 = A -1 2 QA -1 2 -A -1 2 Q2 A -1 2 ,Therefore we further getR ∞ train = 1 T n E tr W ⊤ A -1/2 QA -1/2 W - 1 T n E tr W ⊤ A -1/2 Q2 A -1/2 W + 1 T n E tr ε ⊤ Q 2 εUsing deterministic equivalents in Lemma 1, the training risk then leads toR ∞ train = 1 T n tr W ⊤ A -1/2 QA -1/2 W - 1 T n tr W ⊤ A -1/2 Q2 (I T d )A -1/2 W + 1 T n tr Σ N Q2 + O 1 √ dWe recall the test risk asR ∞ test = tr (Σ N ) + W ⊤ A -1 2 Q2 (A)A -1 2 W T d + tr Σ N Q2 T d + O 1 √ dThe test risk is composed of a signal term of a signal term S = W ⊤ A -1 2 Q2(A)A t) . The signal term reads asS = W ⊤ A -1 2 Q2 (A)A -1 2 W.A -1 2 Q2 (A)A -1 2 = A -1 2 QA 1 2 I + Σ A 1 2 QA -1 2 = A Σ + I -1 I + Σ A Σ + I -1This finally leads toS = W ⊤ A Σ + I -1 I + Σ A Σ + I -1The matrixH = A Σ + I -1 I + Σ A Σ + I -1is responsible to amplifying the signal W ⊤ W in order to let the test risk to decrease more or less. It is is decreasing as function of the number of samples in the tasks n t . Furthermore it is composed of two terms (from the independent training W ⊤ t W) and the cross term W ⊤ t W v for t ̸ = v. Both terms decreases as function of the number of samples n t , smaller values of γ t and increasing value of λ. The cross term depends on the matrix Σ -1 t Σ v which materializes the covariate shift between the tasks. More specifically, if the features are aligned Σ -1 t Σ v = I and the cross term is maximal while for bigger Fisher distance between the covariance of the tasks, the correlation is not favorable for multi task learning. To be more specific the off-diagonal term of H are responsible for the cross term therefore for the multi tasks and the diagonal elements are responsible for the independent terms.To analyze more the element of H, let's consider the case where Σ (t) = I and γ t = γ. In this case the diagonal and non diagonal elements D IL and C M T L are respectively given byD IL = (c 0 (λ + γ) + 1) 2 + c 2 0 λ 2 (c 0 (λ + γ) + 1) 2 -c 2 0 λ 2 , C M T L = -2c 0 λ(c 0 (λ + γ) + 1) (c 0 (λ + γ) + 1) 2 -c 20 λ 2 Both function are decreasing function of λ, 1/γ and c 0 .We recall the definition of the noise term N asN = tr Σ N A -1 + Σ -1Now at the difference of the signal term there are no cross terms due to the independence between the noise of the different tasks. In this case on the diagonal elements of A -1 + Σ -1 matters. This diagonal term is increasing for an increasing value of the sample size, the value of λ. Therefore this term is responsible for the negative transfer. In the specific case where Σ (t) = I d and γ t = γ for all task t, the diagonal terms read asN N T = (c 0 (λ + γ) 2 + (λ + γ) -c 0 λ 2 ) 2 + λ 2 ((c 0 (λ + γ) + 1) 2 -c 2 0 λ 2 ) 2 E.The estimation of the noise covariance remains a technical challenge in this process. However, when the noise covariance is isotropic, it is sufficient to estimate only the noise variance. By observing thatlim λ→0,γ→∞ R ∞ train = σ 2 trQ 2 kn ,we can estimate the noise level from the training risk evaluated at large γ and λ = 0.A version of this work with a stronger focus on time series applications can be found at [20] .Multivariate Time Series Forecasting. Previous work has explored multivariate time series methods to enrich task representations and enhance diversity across time series applications [18] .Despite these advancements, random matrix theory, widely applied in multi-task learning, remains underutilized in multivariate time series contexts, particularly for multivariate time series forecasting (MTSF) [52, 54] . MTSF is common in applications like medical data [6] , electricity consumption [49] , temperatures [27] , and stock prices [44] . Various methods, from classical tools [7, 45] and statistical approaches like ARIMA [3, 4] to deep learning techniques [8, 19, 32, 38, 40, 51, 53, 55, 56] , have been developed for this task. Some studies prefer univariate models for multivariate forecasting [32] , while others introduce channel-wise attention mechanisms [19] . We aim to enhance a univariate model by integrating a regularization technique from random matrix theory and multi-task learning.Architectures without MTL regularization. We follow Chen et al. [8] , Nie et al. [32] , and to ensure a fair comparison of baselines, we apply the reversible instance normalization (RevIN) of Kim et al. [21] . All the baselines presented here are univariate i.e. no operation is performed by the network along the channel dimension. For the PatchTST baseline [32] , we used the official implementation than can be found on Github. The network used in Transformer follows the one in [19] , using a single layer Transformer with one head of attention, while RevIN normalization and denormalization are applied respectively before and after the neural network function. The dimension of the model is d m = 16 and remains the same in all our experiments. DLinearU is a single linear layer applied for each channel to directly project the subsequence of historical length into the forecasted subsequence of prediction length. It is the univariate extension of the multivariate DLinear, DLinearM, used in Zeng et al. [53] . The implementation of SAMformer can be found here, and for the iTransformer architecture here. These two multivariate models serve as baseline comparisons. We reported the results found in [19] and [24] . For all of our experiments, we train our baselines PatchTST, DLinearU and Transformer with the Adam optimizer [22] , a batch size of 32 for the ETT datasets and 256 for the Weather dataset , and the learning rates summarized in Table 2 .Architectures with MTL Regularization. We implemented the univariate PatchTST, DLinearU, and Transformer baselines with MTL regularization. Initially, we scale the inputs twice using RevIN normalization. The first scaling is applied to the univariate components, and the second scaling is applied to the multivariate components. For each channel, we then apply our model without MTL regularization. The outputs are concatenated along the channel dimension, and this concatenation is flattened to form a matrix of shape (batch size, q × T ), where q is the prediction horizon and T is the number of channels. We then learn a square matrix W of shape (q × T ) × (q × T ) for projection and reshape the result to obtain an output of shape (batch size, q, T ). This method can be applied on top of any univariate model. Our regularized loss has been introduced in 5. We conduct our experiments on 3 publicly available datasets of real-world time series, widely used for multivariate long-term forecasting [8, 32, 51] . The 2 Electricity Transformer Temperature datasets ETTh1, and ETTh2 [55] contain the time series collected by electricity transformers from July 2016 to July 2018. Whenever possible, we refer to this set of 2 datasets as ETT. Weather [27] contains the time series of meteorological information recorded by 21 weather indicators in 2020. It should be noted Weather is large-scale datasets. The ETT datasets can be downloaded here while the Weather dataset can be downloaded here. Table 3 sums up the characteristics of the datasets used in our experiments. While the study provides valuable insights through its theoretical analysis within a linear framework, it is important to acknowledge its limitations. The linear approach serves as a solid foundation for understanding more complex models, but its practical applications may be constrained. Linear models, though mathematically tractable and often easier to interpret, might not fully capture the intricacies and nonlinear relationships present in real-world data, especially in the context of multivariate time series forecasting.To address this limitation, we decided to extend our algorithm's application to more complex models, specifically within the nonlinear setting of neural networks. This transition aims to evaluate whether the theoretical insights derived from the linear framework hold true empirically when applied to neural networks. As part of this endeavor, an optimal parameter lambda was selected by an oracle, leading to promising results, as detailed in Section 5.4. This oracle-based selection underscores the potential efficacy of our approach when appropriately tuned, even in more complex, nonlinear contexts.It is important to note that the limitations are not related to the real-world data itself, as our setting performs well in the context of multi-task regression for real-world data, as shown in Section 5.2. The difficulty arises from transitioning from a linear to a nonlinear model. The results in Section 5.4 are particularly encouraging, demonstrating that our method can improve upon univariate baselines by regularizing with an optimal lambda, as indicated by our oracle. While the oracle provides an upper bound on performance, actual implementation would require robust methods for hyperparameter optimization in non-linear scenarios, which remains an open area for further research.By expanding the scope of our theoretical framework to encompass nonlinear models, we pave the way for future work that could focus on the theoretical analysis of increasingly complex architecturesNeurIPS Paper Checklist 


Decision: The paper presents a framework for multitask regression and precisely characterizes train and test risk incurred by the resulting estimator, leveraging random matrix theory. The approach is demonstrated on both multitask regression and multivariate time series forecasting. 

The AC and majority of reviewers find the proposed framework and analysis to be significant and exciting. 

We strongly encourage the authors to incorporate the discussion materials in their manuscript, and in particular their excellent points to justify their assumptions. We also urge them to add a small paragraph at the beginning of section 5.3, to reflect their rebuttal comments on the relevance of their theoretical insights beyond the case of linear models.

Reviews: review: paper_summary: paper_summary: The authors derived theoretical insights into the train and test risks of the multi task regression loss using Random Matrix Theory.
strength_weakness: strength_weakness: 1. The paper lacks clarity for instance on line 89, I am not sure what is meant by general optimization and other mathematical tools?
questions: questions: 1. What is $\gamma$ in (3)? 
2. Have you studied the implications of the assumptions made to use random matrix theory to derive the expressions of asymptotic train and test errors?
3. Shouldn't it be $\mathbf{Y}\in \mathbb{R}^{Tq\times n}$ in the equations after line 110?
limitations: limitations: Limitations are not discussed. My guess is that the assumptions made to be able to employ Random Matrix Theory need to be discussed.
paper_summary: paper_summary: The authors characterise the train and test risk of multi-task regression using random matrix theory. Assessment via Figures 2 and 3 shows a good match between theory and empirical. This motivates a regularised objective for learning mutivariate time-series.
strength_weakness: strength_weakness: *Three major weaknesses*

A. There are at least two existing works on the theoretical error bounds for multi-task learning [Chai, NeurIPS 2009; Aston & Sollich, NeurIPS 2012]. While these are using Gaussian processes as a basis, comparison or remarks on the general effect of multi-task learning should be discussed with reference to two existing works. This is especially when all these are essentially linear-in-regressors models.

B. The proof or at least outline of the proof of Theorem 1 *should* be in the main paper, given that this is *the* major contribution of the paper.

C. It is totally unclear the relevance and contribution of the theory itself to section 5.3. Adding a cross-task regularization is already known to help (mostly).

*Minor, but affects clarity*
- The introduction to the paper is on multi-task multi-response/regression. The multi-response/regression part is not made clear.
- Why is there a need to divide by $\sqrt{Td}$ in (1)? Is it to make analysis easier.
- There are two different uses of the letter $t$ in Assumption 1.
- The setup in section 5.3 needs to relate back to the model in section 3.1 to make clearer their connections.
questions: questions: 1. The analysis of the noise term in section 4.2 says that increase in noise negatively affects transfer. However, it is precisely in the presence of  noise that borrowing statistical strength from other tasks is suppose to help. How do we reconcile this?

2. In section 4.3, is it clear that $C_{MTL}$ is always positive, negative or no such conclusion can be made?
limitations: limitations: Authors should include some technical limitations of their current work.
paper_summary: paper_summary: The authors analyse the problem of multi-task regression (for a linear model) under the assumption of concentrated random vectors.  The results obtained provide a comprehensive description of the performance of the model (including critically its generalisation performance).  The authors use this result for hyperparameter optimisation and show on real-world datasets that using the insights from their theoretical model they obtain state-of-the-art performance on multivariate time series forecasting.
strength_weakness: strength_weakness: The paper will not be everyones taste.  The mathematics is not easy to follow, but that is often the nature of making a non-trivial technical contribution.

In the equation for $\lambda^*$ immediately preceding section 4.4 there is a missing norm symbol.
questions: questions: Is there an intuition on why the concentrated random vector assumption justified?  Is it that real data satisfies this assumption accurate, or is that the deviations from this assumption don't seem to be relevant?

You are tackling real world data using a relatively simple linear model and beating what I assume are sophisticated non-linear models.  I find this surprising.  Is there an explanation of why you do so well?  (Are the datasets you are using very simple--linear?,  are the number of examples small compared to the number of features so that more complex models will overfif?).
limitations: limitations: This is fine.
paper_summary: paper_summary: This paper analyzes a linear multi-task regression framework. It derives formulas for asymptotic train and test risk. The formulas provide insights in how the raw data covariances, singal-generating hyperplanes, noise levels, and size of data sets affect the risk. Motivated by the analysis on the linear framework, experiments on multivariate time series forecasting are performed to investigate the effectiveness of MTL regularization.
strength_weakness: strength_weakness: - The analyses are performed on a simple linear model. They do not apply to general nonlinear models, which are much more common in real practice.
- There seems to be a disconnection between the theoretical analysis and the experimental results. The experiments seem to only highlight the effectiveness of MTL regularization and have no connection with the theoretical analysis.
questions: questions: - Is there any connections between the theoretical analysis and the experiments? E.g. making use insights from the analysis to improve experimental results.
- Is splitting the linear operator into a local and global term a novelty of this paper?
- Similarly, is this paper the first work to apply this MTL framework (decompose prediction as a sum of global and local term) to nonlinear models?
- How are the regularization hyperparameter chosen in the experiments?
- In some of the results, incorporating MTL regularization worsen the performance. Why is that?
limitations: limitations: N/A

Comments: title: Official Comment by Authors
Acknowledgement of feedback comments.
comment: Thank you for your advice! We will revise the paper with respect to this point.
Thanks. For the test risks plot, either put it in the main paper, or point the author to the appendix from the main paper. You cannot assume we read the appendix.
We sincerely appreciate your insightful feedback, which has significantly contributed to the improvement of our paper. We believe that we have carefully addressed your comments and concerns, and we hope our answers meet your expectations.

As the reviewer-author discussion period is nearing its end and our window for responding will close soon, we kindly ask you to let us know if you have any additional points that need to be clarified, so that we could  engage in further discussion.

Thank you again for your valuable review. We look forward to hearing from you and hope our efforts align with your suggestions.
Thank you for your encouraging feedback and for recognizing the technical challenges of our work. We greatly appreciate your support and belief in the value of our approach.
Thank you for your constructive feedback and for increasing the score. We have incorporated your suggested changes into the main paper and believe they greatly enhance its clarity. In addition, please note that we have plotted the test risks in Appendix (Section G.4). Please let us know if you have any other question. Your insights are invaluable to improving the paper.
Thank you for your thoughtful feedback and positive evaluation of the paper's theoretical contributions.
Thank you for addressing my concerns. I understand the paper better now. I have no doubt that the proposed method certainly has novel theoretical contributions. I am raising my score.
C. You wrote "Our results show that the test risk curves for non-linear models follow similar patterns to those predicted by our theory." It is important of show this test risk curves to proof your point.

I will increase the score. Content once, it seems that everything is there based on your answers. However, this does suggest some rewrite of the paper which is not insignificant.
Thank your for addressing questions that I appreciate don't have easy answers.  I maintain my believe that your paper represents a strong technical contribution and that using theory as a guide to improving machine learning algorithms is an achievement.  Good luck with convincing the other reviewers and don't be discouraged if the outcome is negative.  The type of work you are doing is technically challenging and difficult to communicate, but I believe it is beneficial to the field.