{"id":"849322e2-7b7a-418e-b4cd-e00117c8566f","arxiv_id":"2501.07964","paper_version":4,"verdict":"REJECT","confidence":"HIGH","novelty_score":2.0,"correctness_risk":"high","formal_verification":"none","parameter_count":0,"one_line_summary":"A tutorial that re-derives the known EM and gradient formulas for multi-task Gaussian processes, with two mathematical errors in the presented derivations.","lead":"This paper re-derives the math behind multi-task Gaussian processes, covering the EM algorithm and gradient-based hyperparameter updates. It aims to help engineers and researchers understand and verify MTGP implementations used in Bayesian optimization.","discovery_kind":"review","skeptic_critique":{"model":"deepseek-v4-flash","headline":"Equation 16 has a sign error in the data-fit gradient term, so the paper's central gradient derivation is not reliable as written.","rationale":"I agree with the reader's REJECT verdict, but the most load-bearing concern is not the stated 'weakest assumption' about separable covariance. The separable model K_x ⊗ K_f is explicitly the Bonilla et al. model that the paper declares it is deriving, so importing it is a legitimate tutorial choice. The decisive defect is internal: Eq. 16 contradicts Eq. 13 by an elementary chain-rule calculation, and the paper's central claim is to provide friendly, reliable derivations of MTGP formulations and gradients. Since Section 4.2 is the recommended practical gradient approach, a wrong sign in its central formula means practitioners cannot use the paper to verify or implement the method. The additional Monte Carlo EM objective error in Section 4.1.3 reinforces that the derivations have not been checked against the formulas being presented. No verdict change is needed; the rejection stands.","tokens_in":12497,"tokens_out":7862,"duration_ms":74972,"concrete_test":"Set up a small MTGP example (e.g., N=2, M=2) with fixed positive-definite K_x, K_f, diagonal Sigma, and arbitrary y. Let A = K_x ⊗ K_f + I_N ⊗ Sigma. Choose one hyperparameter z (e.g., a lengthscale in K_x or an off-diagonal Cholesky entry L_ij). Compute L_mll from Eq. 13 and its gradient by central finite differences. Compare this numerical gradient with Eq. 16 as printed and with the corrected expression -(1/2)tr(QG) + (1/2)y^T Q G Q y. If the printed formula disagrees in the sign of the quadratic term and the corrected formula matches to numerical precision, the sign error is confirmed.","verdict_should_be":"UNCHANGED","load_bearing_attack":"Equation 16 is the main deliverable of the gradient approach in Section 4.2.1: it states dL_mll/dz = -(1/2)(tr(QG) + y^T Q G Q y), where Q = (K_xf + S)^-1 and G = d(K_xf + S)/dz. Re-deriving from Eq. 13, L_mll = -(1/2)log|A| - (1/2)y^T A^-1 y with A = K_xf + S, gives dL_mll/dz = -(1/2)tr(A^-1 A') - (1/2)y^T d(A^-1)/dz y = -(1/2)tr(QG) + (1/2)y^T Q G Q y, using d(A^-1)/dz = -A^-1 A' A^-1. Thus Eq. 16 has the wrong sign on the quadratic, data-fit term. A practitioner following Eq. 16 under gradient ascent would move that term in the wrong direction. This is an internal inconsistency, not a modeling choice, and it directly undermines the paper's stated purpose of providing reliable derivations of MTGP gradients. A related issue appears in Section 4.1.3 Step 2, where the Monte Carlo EM objective is written as N log( (1/K) sum_k |A_k| ) rather than N (1/K) sum_k log|A_k|, which is the correct MC estimate of the expected complete-data log-likelihood term.","agreement_with_reader":"partial"},"referee_report":{"model":"deepseek-v4-flash","summary":"The paper presents detailed derivations of the multi-output (multi-task) Gaussian process model of Bonilla et al. (2007), covering the EM algorithm for the output-correlation kernel and the gradient-based maximization of the marginal log-likelihood. It also discusses Kronecker-structured efficient computations, partial-observation extensions, and practical considerations for Bayesian optimization. The stated purpose is to provide friendly, verifiable derivations of the MTGP formulations and their gradients.","tokens_in":12800,"tokens_out":10049,"duration_ms":81650,"significance":"If the derivations were correct, the paper would be a useful reference for practitioners implementing MTGP, particularly because it fills derivation gaps in the original treatment and connects to modern automatic differentiation and scalable linear algebra. I credit the paper for including explicit proofs of the linear algebra identities and for discussing partial-observation extensions. However, the value of the paper is entirely contingent on mathematical accuracy, and the presence of several incorrect key formulas currently undermines that value.","major_comments":[{"comment":"The derivative of the marginal log-likelihood is stated as ∂L_mll/∂z = -(1/2)(tr(QG) + y^T Q G Q y). Re-deriving from Eq. (13) using d(A^{-1})/dz = -A^{-1} A' A^{-1} gives ∂L_mll/∂z = -(1/2) tr(QG) + (1/2) y^T Q G Q y. The sign of the data-fit term is wrong. Because Eq. (16) is the main deliverable of the gradient derivation in Section 4.2, a practitioner following the paper would move the data-fit gradient in the wrong direction during hyperparameter optimization.","section":"Section 4.2.1, Eq. (16)"},{"comment":"The listed derivatives of K_xf + S are incorrect. Since K_xf = K_x ⊗ K_f, the correct derivatives are ∂(K_xf+S)/∂θ_i = (∂K_x/∂θ_i) ⊗ K_f and ∂(K_xf+S)/∂L_{i,j} = K_x ⊗ (E_{i,j} L^T + L E_{j,i}). The paper instead writes ∂K_x/∂θ_i ⊗ I_M and I_N ⊗ (E_{i,j} L^T + L E_{j,i}), which are not equivalent in general and are only valid under a special structure. Since Eq. (22) for the partial-observation case relies on these formulas, the error propagates to that extension.","section":"Section 4.2.1, Eq. (17)"},{"comment":"The Monte Carlo EM objective is written as N log( (1/K) Σ_k |(F~(k))^T K_x^{-1} F~(k)| ), but the correct Monte Carlo estimate of the expected complete-data log-likelihood term is N (1/K) Σ_k log|(F~(k))^T K_x^{-1} F~(k)|. The logarithm must be inside the average; the stated form is the log of an average determinant, which is not an unbiased estimate of the expected log-determinant and will bias the EM update.","section":"Section 4.1.3, Step 2"}],"minor_comments":[{"comment":"The intermediate expression '-N/M Σ_m ∂log σ_m/∂σ_i' appears to be a typo; the final result agrees with the correct derivative -N/σ_i, but the displayed formula should be corrected to avoid confusion.","section":"Section 4.1.2, Eq. (7)"},{"comment":"The line 'log|F^T K_x^{-1} F| |N^{-M}|' is confusing; since |F^T K_x^{-1} F / N| = |F^T K_x^{-1} F| N^{-M}, the constant should be written additively as -M log N. The final expression is unaffected because constants are discarded, but the presentation needs fixing.","section":"Eq. (11)"},{"comment":"The phrase 'By taking the derivative of zero' should read 'By setting the derivative to zero'.","section":"Section 4.1.2"},{"comment":"The symbol Q is used for K_f^{-1} in Section 4.1.2 and for (K_xf+S)^{-1} in Section 4.2.1; this overloading is confusing and should be replaced with distinct symbols.","section":"Notation"},{"comment":"Several self-citations to the author's TPE papers appear tangential to the derivation content, and Watanabe et al. (2022) and (2023) share the same arXiv identifier; the references should be cleaned up.","section":"Section 5"}],"recommendation":"major_revision","confidential_remarks":"This is a derivation/tutorial manuscript with limited novelty but potentially useful reference value. The main obstacle is correctness: the sign error in Eq. (16) and the incorrect derivative expressions in Eq. (17) are in the central deliverables, and the EM step in Section 4.1.3 contains a further error. I would be willing to reconsider a revised version that corrects these issues and carefully re-verifies all displayed equations."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Two things you should know. First, this is a tutorial, not a research contribution: it re-derives known MTGP results from Bonilla et al. (2007) and Lin et al. (2024), and it introduces no new method, theorem, or experiment. Second, it has two real math errors, one of them in the main deliverable.\n\nWhat it does well: the exposition is clear and well organized. The derivations of the complete-data log-likelihood, the M-step, the Cholesky parameterization, and the partial-observation extension via the Khatri-Rao product are mostly solid and do fill in gaps left by the original papers. The NumPy-style identities (Eqs. 23-24) are genuinely useful for practitioners translating the linear algebra into code. The citation pattern is honest; the self-citations to TPE papers support alternative BO methods and are not used to prop up the MTGP derivation.\n\nThe soft spots are load-bearing. Eq. (16) gives dL/dz = -(1/2)(tr(QG) + y^T Q G Q y). Starting from Eq. (13) and the derivative of a matrix inverse, the quadratic term should be positive: dL/dz = -(1/2)tr(QG) + (1/2)y^T Q G Q y. I re-derived this; the stress-test note is correct. The sign error would send the data-fit term in the wrong direction under gradient ascent, and since the paper's stated purpose is to help practitioners verify their code, this is exactly the kind of mistake that defeats the purpose. The second error is in Section 4.1.3, Step 2: the Monte Carlo EM objective should be argmin [M log|Kx| + N * (1/K) sum_k log|(F~^k)^T Kx^{-1} F~^k|], not N log((1/K) sum_k |...|). Both errors are fixable, but they are not typographical noise; they contradict equations earlier in the same paper.\n\nWho this is for: practitioners implementing MTGP in Bayesian optimization and students who want a step-by-step derivation. A corrected version would be a serviceable reference, but as written I would not hand it to someone as a trustworthy source for gradients. It deserves a serious referee because the topic is useful, the presentation is good, and the errors are localized and repairable. My recommendation for peer review: send it out, but ask reviewers to check the gradient formulas line by line.","headline":"A well-organized tutorial re-derivation of MTGP, but the central gradient formula in Eq. (16) has a sign error and the Monte Carlo EM objective in Step 2 is misspecified; both need fixing before this can be a reliable reference.","tokens_in":13331,"tokens_out":3901,"would_cite":false,"duration_ms":38907,"reading_group":"maybe","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":[],"pacs":[],"model":"deepseek-v4-flash","headline":"A full step-by-step derivation shows how to build and differentiate the multi-task Gaussian process, giving practitioners a way to check their implementations.","keywords":["multi-task Gaussian process","multi-output Gaussian process","Kronecker structure","EM algorithm","marginal likelihood gradients","Bayesian optimization","coregionalization"],"falsifier":"Take a small synthetic two-output, two-input MTGP with chosen hyperparameters, compute the marginal log-likelihood of Eq. (13) directly, and compare its finite-difference gradient to the closed-form gradient formulas in Eqs. (16)–(17); any mismatch beyond numerical tolerance would disprove the derivation.","tokens_in":12278,"feed_emoji":"📐","tokens_out":4357,"duration_ms":39143,"temperature":0.7,"pith_summary":"This paper aims to make the multi-task Gaussian process (MTGP) derivable from scratch. Its central claim is that the EM updates and gradient-based hyperparameter optimization used in the standard MTGP model can be written down explicitly, step by step, with all intermediate algebra shown. If correct, practitioners can verify their MTGP code against these formulas instead of trusting an opaque implementation. The value is practical: multi-task GP is a common tool in Bayesian optimization, and the paper's derivations remove a barrier to using it correctly.","feed_headline":"Step-by-step math closes multi-task GP derivation gaps","feed_subtitle":"Practitioners can now verify MTGP implementations from EM updates to gradient-based hyperparameter learning.","key_machinery":"The central object is the Kronecker-structured covariance $K_{xf}=K_x\\otimes K_f$ combined with the noise term $S=I_N\\otimes\\Sigma$. The argument is carried by standard linear-algebra identities: $(A\\otimes B)^{-1}=A^{-1}\\otimes B^{-1}$, $|A\\otimes B|=|A|^M|B|^N$, and the trace identity $\\operatorname{vec}(C^\\top)^\\top(A\\otimes B)\\operatorname{vec}(C^\\top)=\\operatorname{tr}(A C B^\\top C^\\top)$. These convert high-dimensional vectorized expressions into matrix products of size $N\\times N$ and $M\\times M$, which is what makes both the EM updates and the gradient formulas derivable and computationally efficient.","core_discovery":"On Bonilla et al.'s MTGP model, a zero-mean GP prior over latent functions with covariance $K_x \\otimes K_f$ and observation noise $I_N\\otimes\\Sigma$, the complete-data log-likelihood can be transformed using Kronecker identities and trace identities into a tractable objective. The paper derives the EM updates $\\hat{\\sigma}_i^2 = \\frac{1}{N}\\sum_n (y_{n,i}-f_{n,i})^2$, $\\hat{K}_f = F^\\top K_x^{-1}F/N$, and the reduced objective for $\\boldsymbol{\\theta}$ as $M\\log|K_x| + N\\log|F^\\top K_x^{-1}F|$. It likewise derives the marginal log-likelihood gradients for the gradient approach, including a Cholesky parameterization of $K_f$ and a projection matrix $P$ that extends the same formulas to settings where some outputs are missing at some inputs, preserving $O(N^3+M^3)$ complexity with iterative solvers.","pith_inferences":["The derivations suggest that any implementation of MTGP that fails to reproduce these closed-form updates on a small synthetic dataset very likely has a coding error, so checking against the formulas is a fast smoke test.","The same trace-and-Kronecker machinery presumably extends to other separable multi-output models, such as linear coregionalization models, and the projection-matrix trick may carry over to those settings as a way to handle missing outputs.","The paper's reliance on the separable prior $K_x\\otimes K_f$ implies that for data where output correlations change with input location, the formulas should not be expected to hold, and a non-separable prior would need its own derivation."],"forward_implications":["The derived formulas give practitioners a concrete reference for checking MTGP implementations, particularly the EM updates and the marginal-likelihood gradients.","Because the gradient approach can be written with explicit derivatives, automatic differentiation packages can be used directly, and the closed-form derivatives document what the autodiff is computing.","The partial-observation extension via the projection matrix $P$ covers realistic datasets where some outputs are missing at some inputs, without changing the asymptotic complexity when combined with iterative linear solvers.","The diagonalization of the Kronecker structure recovers the $O(N^3+M^3)$ training cost for the full-observation case, the regime on which practical MTGP implementations rely."],"supporting_citations":[{"why":"Supplies the MTGP model, the EM algorithm, and the gradient approach whose derivation gaps this paper fills.","marker":"Bonilla et al. (2007)"},{"why":"Provides the Kronecker-diagonalization trick that reduces time and memory complexity to $O(N^3+M^3)$ and $O(N^2+M^2)$.","marker":"Stegle et al. (2011)"},{"why":"Extends the efficient computations to the non-Kronecker, partial-observation setting used in Section 4.2.2.","marker":"Lin et al. (2024)"},{"why":"Automatic differentiation in PyTorch is the practical route the paper recommends for computing the derived gradients.","marker":"Paszke et al. (2019)"},{"why":"Provides the posterior-sampling method that leverages the Kronecker structure for efficient sampling.","marker":"Maddox et al. (2021)"},{"why":"Lanczos quadrature is used for computing the log-determinant term in the partial-observation marginal likelihood.","marker":"Dong et al. (2017)"},{"why":"Matheron's rule is used for efficient posterior sampling on top of the Kronecker structure.","marker":"Wilson et al. (2020)"}],"fun_headline_variants":["Multi-task GP math, fully derived","Step-by-step derivations for MTGP","Clear derivations for multi-task GPs","Deriving MTGP gradients and EM updates","MTGP formulas and gradients explained"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The paper's formulas hold for the specific MTGP model where the prior covariance between any two outputs is separable as a product of an input kernel and a fixed output kernel; if output correlations change across the input space, or are otherwise not captured by this product structure, the derived updates and gradients no longer apply.","fun_headline_variants_meta":{"raw":{"variants":["Multi-task GP math, fully derived","Step-by-step derivations for MTGP","Clear derivations for multi-task GPs","Deriving MTGP gradients and EM updates","MTGP formulas and gradients explained"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000211,"raw_usage":{"total_tokens":1362,"prompt_tokens":838,"completion_tokens":524,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":454,"completion_tokens_details":{"reasoning_tokens":461}},"tokens_in":454,"tokens_out":524,"duration_ms":5997,"temperature":1.0,"reasoning_tokens":461,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-10T20:30:20.143111+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Take a small synthetic two-output, two-input MTGP with chosen hyperparameters, compute the marginal log-likelihood of Eq. (13) directly, and compare its finite-difference gradient to the closed-form gradient formulas in Eqs. (16)–(17); any mismatch beyond numerical tolerance would disprove the derivation.","supporting_citations":[{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Supplies the MTGP model, the EM algorithm, and the gradient approach whose derivation gaps this paper fills."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Provides the Kronecker-diagonalization trick that reduces time and memory complexity to $O(N^3+M^3)$ and $O(N^2+M^2)$."},{"cited_title":"Scaling Gaussian Processes for Learning Curve Prediction via Latent Kronecker Structure","cited_arxiv_id":"2410.09239","evidence_quote":"Extends the efficient computations to the non-Kronecker, partial-observation setting used in Section 4.2.2."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Automatic differentiation in PyTorch is the practical route the paper recommends for computing the derived gradients."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Provides the posterior-sampling method that leverages the Kronecker structure for efficient sampling."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Lanczos quadrature is used for computing the log-determinant term in the partial-observation marginal likelihood."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Matheron's rule is used for efficient posterior sampling on top of the Kronecker structure."}],"review_version":1}