Skip to contents

Assembles .compute_global_variance() and factorises it once via base::chol(), returning the matrix itself together with its Cholesky factor, log-determinant (via the factor's diagonal, not base::det()), and inverse (via base::chol2inv(), which reuses the factor rather than an independent base::solve()).

Usage

.sigma_p_factorisation(p, Sigma)

Arguments

p

Numeric vector \(\boldsymbol{p}\in\mathbb{R}^{J}\).

Sigma

Array of cell-type covariances in \(\mathcal{M}_{G\times G\times J}\).

Value

A list with elements: matrix (\(\boldsymbol{\Sigma}(\boldsymbol{p})\) itself), chol (upper-triangular Cholesky factor), log_det (\(\log\det\boldsymbol{\Sigma}(\boldsymbol{p})\)) and inverse (\(\boldsymbol{\Sigma}(\boldsymbol{p})^{-1}\)).

Details

stats::optim(), stats::nlminb() and marqLevAlg::marqLevAlg() each treat the log-likelihood, gradient and Hessian as three independent callback functions, but all Newton-type solvers evaluate them at the SAME trial point \(\boldsymbol{p}\) within one iteration (the Hessian chain rule in hessian_loglik_constrained() even re-evaluates the unconstrained gradient a second time internally). Without sharing work, one (log-lik, gradient, Hessian) triple pays for the \(O(G^{3})\) assembly-and-factorisation of \(\boldsymbol{\Sigma}(\boldsymbol{p})\) up to four times over; this single-slot cache, keyed on exact equality of p and Sigma, means only the first of those calls actually factorises, and the rest simply return the cached result in \(O(G^{2})\) (the cost of the equality check). Profiling on a 38-gene / 3-cell-type scenario showed the redundancy