Multi-Task Gaussian Process

• Multi-task Gaussian processes are an extension of Gaussian Process regression that jointly models multiple correlated outputs by factorizing the covariance into input and task kernels.
• They utilize methods such as intrinsic coregionalization and linear models of coregionalization to share statistical strength, enabling efficient information transfer across varied tasks.
• Applications span nuclear physics, multi-robot coverage, and engineering, with tailored approaches addressing scalability, misaligned inputs, and heterogeneous observation domains.

Multi-task Gaussian process (MTGP), also called multi-output Gaussian process, extends Gaussian process regression from a single scalar-valued function to a collection of correlated outputs or tasks. In a standard formulation, $f_\cdot:\{1,\dots,M\}\times \mathcal X\to\mathbb R$ is jointly Gaussian with covariance that factors into an input-space kernel and an output-space kernel, so that $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$; closely related constructions write $$\langle f_\tau(x)f_{\tau'}(x')\rangle=D_{\tau\tau'}C(x,x')$$ or, in linear coregionalization form, $$\mathbf K(x,x')=\sum_{q=1}^Q \mathbf B^{(q)}k_q(x,x')$$ (Watanabe, 14 Jan 2025, Ashton et al., 2012, Comlek et al., 9 Jan 2026). The framework is used to transfer statistical strength across related outputs, fidelities, sensing modalities, or intervention functions while retaining posterior means, posterior covariances, and marginal-likelihood-based hyperparameter learning; later work has broadened the paradigm to common-mean models, clustered mixtures, heterogeneous input domains, temporally misaligned tasks, causal intervention functions, and fast structured solvers (Leroy et al., 2020, Leroy et al., 2020, Liu et al., 2022, Mikheeva et al., 2021, Aglietti et al., 2020, Sorokin et al., 16 Mar 2026).

1. Formal definition and covariance constructions

A recurring MTGP construction is the separable or product covariance model. For $T$ regression tasks $f_\tau(x)$, one writes

$$\langle f_\tau(x)f_{\tau'}(x')\rangle = D_{\tau\tau'}\,C(x,x'),$$

where $C(x,x')$ governs smoothness, length scale, and roughness in the input space, while $D$ is a free-form inter-task covariance matrix (Ashton et al., 2012). In vectorized form, for $N$ input points and $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$0 tasks, the latent covariance can be written as $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$1, with task-specific observation noise $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$2, yielding

$\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$3

in the fully observed case (Watanabe, 14 Jan 2025). This Kronecker factorization is the standard algebraic backbone of many MTGPs.

The intrinsic coregionalization model (ICM) is a direct specialization of this product form. In the nuclear mass and charge-radius setting, the covariance between data point $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$4 from task $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$5 and data point $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$6 from task $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$7 is

$\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$8

with a $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$9 task covariance matrix $$\langle f_\tau(x)f_{\tau'}(x')\rangle=D_{\tau\tau'}C(x,x')$$0; the off-diagonal blocks $$\langle f_\tau(x)f_{\tau'}(x')\rangle=D_{\tau\tau'}C(x,x')$$1 and $$\langle f_\tau(x)f_{\tau'}(x')\rangle=D_{\tau\tau'}C(x,x')$$2 are precisely what permit cross-task information sharing (Ye et al., 23 Jul 2025).

A more flexible family is the linear model of coregionalization (LMC),

$$\langle f_\tau(x)f_{\tau'}(x')\rangle=D_{\tau\tau'}C(x,x')$$3

which represents the matrix-valued kernel as a sum of latent kernel components, each with its own task-correlation matrix (Comlek et al., 9 Jan 2026). The semiparametric latent factor model (SLFM) is a low-rank special case with $$\langle f_\tau(x)f_{\tau'}(x')\rangle=D_{\tau\tau'}C(x,x')$$4, reducing parameter count from $$\langle f_\tau(x)f_{\tau'}(x')\rangle=D_{\tau\tau'}C(x,x')$$5 to $$\langle f_\tau(x)f_{\tau'}(x')\rangle=D_{\tau\tau'}C(x,x')$$6 (Comlek et al., 9 Jan 2026).

The same structural idea appears in graph-based multitask fields. In multi-robot coverage, the stacked demand vector over vertices and tasks is assigned a Gaussian prior

$$\langle f_\tau(x)f_{\tau'}(x')\rangle=D_{\tau\tau'}C(x,x')$$7

explicitly described as consistent with the multitask GP model of Bonilla et al.; here $$\langle f_\tau(x)f_{\tau'}(x')\rangle=D_{\tau\tau'}C(x,x')$$8 encodes spatial covariance across vertices and $$\langle f_\tau(x)f_{\tau'}(x')\rangle=D_{\tau\tau'}C(x,x')$$9 encodes inter-task correlation (Wei et al., 11 Mar 2026). The covariance therefore carries two geometries at once: a geometry over inputs and a geometry over tasks.

2. Alternative information-sharing mechanisms

Classical MTGPs share information through covariance coupling, but the cited literature also shows that this is not the only probabilistic mechanism available. In MAGMA, each task is modeled as

$$\mathbf K(x,x')=\sum_{q=1}^Q \mathbf B^{(q)}k_q(x,x')$$0

where $$\mathbf K(x,x')=\sum_{q=1}^Q \mathbf B^{(q)}k_q(x,x')$$1 is a common mean process shared by all tasks, $$\mathbf K(x,x')=\sum_{q=1}^Q \mathbf B^{(q)}k_q(x,x')$$2 is a task-specific zero-mean latent GP, and $$\mathbf K(x,x')=\sum_{q=1}^Q \mathbf B^{(q)}k_q(x,x')$$3 is Gaussian noise (Leroy et al., 2020). Tasks are conditionally independent given $$\mathbf K(x,x')=\sum_{q=1}^Q \mathbf B^{(q)}k_q(x,x')$$4, so inter-task dependence is induced by a shared random mean rather than by a task covariance matrix in the usual LMC sense. This replaces covariance-only coupling by prior-trend coupling.

MagmaClust pushes this idea further by introducing a mixture of multi-task GPs with cluster-specific mean processes. Conditional on cluster $$\mathbf K(x,x')=\sum_{q=1}^Q \mathbf B^{(q)}k_q(x,x')$$5, the observed function for individual $$\mathbf K(x,x')=\sum_{q=1}^Q \mathbf B^{(q)}k_q(x,x')$$6 is

$$\mathbf K(x,x')=\sum_{q=1}^Q \mathbf B^{(q)}k_q(x,x')$$7

with $$\mathbf K(x,x')=\sum_{q=1}^Q \mathbf B^{(q)}k_q(x,x')$$8 a cluster-specific mean process, $$\mathbf K(x,x')=\sum_{q=1}^Q \mathbf B^{(q)}k_q(x,x')$$9 an individual-specific centered GP, and $T$0 observation noise (Leroy et al., 2020). The predictive distribution for a new partially observed task is then a mixture of cluster-specific GP posteriors rather than a single Gaussian. This architecture is designed for group-structured functional data, where a single global sharing mechanism is too restrictive.

Another departure from additive latent decomposition is the hierarchical-latent-interaction model. Standard LMC and convolution-based MTGPs assemble outputs from additive independent latent functions, and the cited work argues that this misses interactions among latent functions and among task-mixing coefficients (Chen et al., 2018). The proposed kernel

$T$1

combines function interaction through generalized convolution spectral mixture terms $T$2 and coefficient interaction through cross-coregionalization terms $T$3 (Chen et al., 2018). This suggests that MTGP dependence can be organized at several levels: task space, latent-function space, and interactions between the two.

3. Posterior inference and training objectives

In the standard fully observed setting, MTGP learning proceeds by maximizing the marginal log-likelihood associated with the covariance $T$4. The output covariance is commonly parameterized as $T$5 with lower-triangular $T$6, which guarantees positive semidefiniteness and permits gradient-based optimization over unconstrained parameters (Watanabe, 14 Jan 2025). This is the direct analogue of scalar GP evidence maximization, but over a much richer covariance geometry.

A different inference route is the conditional one-output likelihood formulation. Cool-MTGP rewrites the multitask likelihood as a chain of conditioned univariate GPs,

$T$7

with each task $T$8 modeled conditionally on previous tasks as

$T$9

(García-Hinde et al., 2020). The method learns only $f_\tau(x)$0 task-specific parameters, namely $f_\tau(x)$1 and $f_\tau(x)$2, and reconstructs both the multitask noise covariance and the intertask covariance afterward (García-Hinde et al., 2020). Its stated purpose is to avoid low-rank approximation and the need to validate the rank hyperparameter $f_\tau(x)$3.

Scalable MTGPs typically replace exact inference by sparse variational approximations. In NSVLMC, each latent GP $f_\tau(x)$4 is represented through inducing variables $f_\tau(x)$5, and the evidence lower bound combines an expected log-likelihood with analytic KL terms for $f_\tau(x)$6 and $f_\tau(x)$7 (Liu et al., 2021). The neural embedding of coregionalization introduces an input-dependent latent transform $f_\tau(x)$8 so that

$f_\tau(x)$9

retaining GP uncertainty while increasing latent diversity without simply increasing the number of base latent GPs $$\langle f_\tau(x)f_{\tau'}(x')\rangle = D_{\tau\tau'}\,C(x,x'),$$0 (Liu et al., 2021).

Continual MTGPs modify the prior-posterior recursion itself. Continual multi-task Gaussian processes reconstruct a new GP prior from the old variational posterior by integrating the GP conditional $$\langle f_\tau(x)f_{\tau'}(x')\rangle = D_{\tau\tau'}\,C(x,x'),$$1 against the previous $$\langle f_\tau(x)f_{\tau'}(x')\rangle = D_{\tau\tau'}\,C(x,x'),$$2, thereby propagating posterior uncertainty forward in function space (Moreno-Muñoz et al., 2019). The resulting continual lower bound contains a triple-KL structure rather than a single regularizer, and the framework is designed to handle sequential batches, non-Gaussian likelihoods, and asynchronous multi-channel observations (Moreno-Muñoz et al., 2019).

4. Transfer, learning curves, and common misconceptions

A recurrent misconception is that adding related tasks necessarily improves asymptotic prediction. The learning-curve analysis for product-form MTGPs shows a more qualified picture. For arbitrary $$\langle f_\tau(x)f_{\tau'}(x')\rangle = D_{\tau\tau'}\,C(x,x'),$$3, the paper derives a self-consistency approximation for the average Bayes error and uses it to study pure transfer, large-$$\langle f_\tau(x)f_{\tau'}(x')\rangle = D_{\tau\tau'}\,C(x,x'),$$4 asymptotics, and the many-task regime (Ashton et al., 2012). In pure transfer learning, if some tasks have no data while other tasks are perfectly learned, the error on an unobserved task cannot go below a finite floor,

$$\langle f_\tau(x)f_{\tau'}(x')\rangle = D_{\tau\tau'}\,C(x,x'),$$5

so infinitely many examples from other tasks do not eliminate uncertainty altogether (Ashton et al., 2012).

The same analysis shows that, as $$\langle f_\tau(x)f_{\tau'}(x')\rangle = D_{\tau\tau'}\,C(x,x'),$$6, multi-task learning can become asymptotically essentially useless unless the degree of inter-task correlation is near its maximal value $$\langle f_\tau(x)f_{\tau'}(x')\rangle = D_{\tau\tau'}\,C(x,x'),$$7, and that this effect is strongest for smooth target functions such as those associated with squared exponential kernels (Ashton et al., 2012). In the many-task regime, the learning curves separate into an initial collective-learning phase, where the Bayes error drops to a plateau value, and a later decay phase once the number of examples becomes proportional to the number of tasks (Ashton et al., 2012). These results place a theoretical limit on the frequently assumed monotonic benefit of cross-task transfer.

Another oversimplification is the idea that multitask gains require strongly positive task correlation. In the nuclear mass and charge-radius study, the learned correlation depends on the feature set: $$\langle f_\tau(x)f_{\tau'}(x')\rangle = D_{\tau\tau'}\,C(x,x'),$$8 for MTL-M9 and $$\langle f_\tau(x)f_{\tau'}(x')\rangle = D_{\tau\tau'}\,C(x,x'),$$9 for MTL-M12, yet the MTGP still consistently outperforms single-task learning (Ye et al., 23 Jul 2025). The paper explicitly notes that the multi-task benefit does not require the tasks to be trivially positively correlated; rather, it depends on how the chosen features interact with the underlying physics (Ye et al., 23 Jul 2025).

A control-oriented theory appears in multi-robot coverage. There the multitask GP is not only a regressor but also a learning engine driving exploration, and the resulting DSMLC algorithm achieves

$C(x,x')$0

under the stated assumptions (Wei et al., 11 Mar 2026). A remark in that work emphasizes that exploiting task correlation changes the uncertainty-decay rate from $C(x,x')$1 when correlation is ignored to $C(x,x')$2 when multitask correlation is used (Wei et al., 11 Mar 2026). This is a stronger operational notion of transfer than Bayes-error reduction alone.

5. Heterogeneous domains, alignment, and nonstandard observation structures

Standard MTGPs are often formulated under a same-input-domain assumption, but several lines of work generalize the model to settings where this assumption fails. HSVLMC considers tasks with different input domains $C(x,x')$3 and introduces prior domain mappings

$C(x,x')$4

into a common latent domain (Liu et al., 2022). Its Bayesian calibration step treats aligned inputs as latent variables with a variational posterior $C(x,x')$5, and the ELBO includes both $C(x,x')$6 and $C(x,x')$7, thereby correcting for dimensionality-reduction error induced by the prior mappings (Liu et al., 2022).

Temporal misalignment is handled in aligned multi-task Gaussian processes by introducing a separate monotonic warp $C(x,x')$8 for each task and evaluating the input kernel on warped times,

$C(x,x')$9

(Mikheeva et al., 2021). The warp is modeled through an ODE-based monotonic stochastic process with GP prior over the drift function, and full Bayesian inference is carried out by variational methods rather than by MAP alignment (Mikheeva et al., 2021). The purpose is to prevent erroneous correlation estimates and poor uncertainty quantification when tasks are phase-shifted or stretched in time.

DAG-GP extends MTGP learning to intervention functions defined on different input spaces and even different dimensionalities. Each causal task $D$0 is written as an integral transform of a shared base function,

$D$1

with the operator $D$2 determined by the DAG structure (Aglietti et al., 2020). This is a genuinely heterogeneous-input MTGP: task relatedness is induced by causal semantics rather than by a fixed task index set alone.

Multi-resolution multi-task Gaussian processes generalize the observation model itself. In MRGP, an observation is not a pointwise sample but an average over a support $D$3,

$D$4

so the same latent field may be observed at different spatial or temporal resolutions and fidelities (Hamelijnck et al., 2019). The framework is explicitly designed for dependent observation processes, multi-fidelity sensing, and biased aggregation operators.

A further extension addresses function-valued inputs. For mechanical systems with functional covariates, the multitask kernel is defined on $D$5 by the fully separable form

$D$6

where $D$7 is built from weighted $D$8 distances between functional inputs and $D$9 captures inter-task correlations (Gninkou et al., 24 Feb 2026). This preserves exact GP conditioning while moving beyond finite-dimensional vector inputs.

6. Scalability and structured computation

The principal computational obstacle in MTGPs is the $N$0 Gram matrix, with $N$1. Standard MTGP regression therefore incurs $N$2 storage and $N$3 computations for factorization and determinant evaluation (Sorokin et al., 16 Mar 2026). Much of the modern MTGP literature is a response to this bottleneck.

One strategy is to exploit algebraic structure. In the time-varying Markov transition-matrix model, the multitask covariance has the Kronecker form $N$4, and when the temporal grid is regular and the kernel is stationary, Toeplitz matrices can be embedded into circulant matrices to reduce complexity from $N$5 to $N$6 for matrix-vector products (Ugurel, 2023). For a 2D grid of size $N$7, the paper states that the cost of estimating the inverse can be reduced to approximately $N$8 instead of the naive $N$9, and informal experiments reported training-time reductions of about $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$00 to $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$01 when exploiting grid-based Toeplitz structure (Ugurel, 2023).

An exact Kronecker reduction appears in functional-covariate MTGPs. With covariance

$\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$02

naive dense cost $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$03 and memory $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$04 are replaced by

$\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$05

and memory $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$06 through mode-wise triangular solves and tensor algebra (Gninkou et al., 24 Feb 2026).

Fast MTGP regression pushes structure further by pairing low-discrepancy designs with special product kernels $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$07 so that the multitask Gram matrix becomes transform-diagonalizable after FFT or Walsh-Hadamard changes of basis (Sorokin et al., 16 Mar 2026). In the equal-sample-size case $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$08, the paper gives the structured costs: evaluate $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$09 in $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$10, invert and compute the determinant in $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$11, and apply $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$12 in $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$13 (Sorokin et al., 16 Mar 2026). The method is targeted at problems with few tasks and very large per-task sample sizes.

Sparse latent-structure parameterizations offer another route. Grouped GP regression with direct sparse Cholesky functional representations exploits pivot-based conditional independence inside each group so that a group Cholesky factor has only $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$14 nonzero elements, reducing storage from $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$15 to $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$16 and computation from $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$17 to $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$18 (Dahl et al., 2019). Cool-MTGP reaches a different tradeoff by replacing full multitask likelihood optimization by conditioned univariate GPs, and its approximate version is reported to scale around $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$19, with an embarrassingly parallel version giving $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$20 per task (García-Hinde et al., 2020).

7. Applications and empirical record

MTGPs have been used in domains where outputs are physically or operationally coupled. In nuclear structure modeling, an MTGP with an ICM covariance over two tasks—nuclear mass and charge radius—uses 12 physical input features and reports root-mean-square deviations of $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$21 MeV for masses and $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$22 fm for charge radii in the MTL-M12 configuration (Ye et al., 23 Jul 2025). The same study reports learned task correlations of $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$23 for MTL-M9 and $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$24 for MTL-M12, showing that useful cross-task transfer can occur even when the fitted correlation is negative (Ye et al., 23 Jul 2025).

In adaptive multi-robot coverage, the multitask GP models sensory demands over a discrete graph as a jointly Gaussian field over vertices and tasks, uses mutual-information-driven exploration, and feeds its posterior mean into a federated multitask coverage controller (Wei et al., 11 Mar 2026). The resulting DSMLC algorithm is proved to achieve sublinear cumulative regret and, in simulations, has smaller cumulative regret than a randomized multitask learning and coverage baseline; in the two-task case with inter-task correlation $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$25, regret is higher than in the single-task case but DSMLC still outperforms RMLC (Wei et al., 11 Mar 2026).

Engineering applications emphasize sparse high-fidelity data. The engineering MTGP framework validated on the Forrester benchmark, 3D ellipsoidal void modeling, and friction-stir welding reports that MTGP improves the primary-task RMSE relative to single-task GP, with gains that are largest when task correlation is high (Comlek et al., 9 Jan 2026). In the ellipsoidal-void setting, reported RMSE improvements are roughly $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$26 to $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$27 in the all-elastic case and $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$28 to $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$29 in the mixed elastic-plastic case; in friction-stir welding, mean RMSE improvements over single-task GP are $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$30, $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$31, $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$32, and $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$33 for $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$34, $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$35, $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$36, and $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$37 experimental training samples, respectively (Comlek et al., 9 Jan 2026).

Spatiotemporal scientific modeling provides a further benchmark. The physics-augmented multi-task GP for cardiac electrodynamics combines a task kernel, a Laplacian-spectral spatial kernel on a ventricular surface mesh, a temporal Matérn-$\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$38 kernel, and PDE-residual regularization from the FitzHugh–Nagumo reaction-diffusion model (Zhang et al., 15 Oct 2025). On a mesh with $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$39 nodes and $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$40 elements, the method predicts normalized electric potential $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$41 and recovery current $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$42; for Protocol I with $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$43 and $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$44, it reports $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$45 at $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$46 and $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$47 at $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$48, with reductions of $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$49 and $\mathrm{Cov}(f_i(x),f_j(x'))=k_x(x,x')k_f(i,j)$50 versus the geometry-aware M-GP baseline (Zhang et al., 15 Oct 2025).

Across these examples, MTGP functions less as a single model than as a family of probabilistic constructions for correlated outputs. The shared core is Gaussian-process conditioning over vector-valued latent structure; the major differences concern how task dependence is encoded, how mismatched domains are aligned, and how the resulting covariance is made computationally tractable.