Skip to contents

Maximises the DeCovarT log-likelihood over the simplex while holding a subset of coordinates at prescribed values. Writing \(A\) for the constrained index set and \(s=\sum_{j\in A}c_j\), the free block is reparametrised as $$ \boldsymbol{p}_{A^{c}} = (1-s)\,\boldsymbol{\psi}(\boldsymbol{z}), \qquad \boldsymbol{z}\in\mathbb{R}^{|A^{c}|-1}, $$ with \(\boldsymbol{\psi}\) the isometric logistic map (isometric_logistic()) on a Helmert basis of the free face. The constrained coordinates are substituted rather than pushed through a logarithm, so a null such as \(p_j=0\) is representable exactly: this is what makes boundary likelihood-ratio tests computable (see lrt_decovart()).

Usage

restricted_mle_decovart(
  y,
  mean_signature_matrix,
  Sigma,
  fixed,
  epsilon = 10^-8,
  itmax = 500L
)

Arguments

y

Numeric vector (or one-column matrix) \(\boldsymbol{y}\in\mathbb{R}^{G}\).

mean_signature_matrix

Numeric matrix \(\boldsymbol{\mu}\in\mathcal{M}_{G\times J}\) (plug-in means).

Sigma

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

fixed

Named or integer-indexed numeric vector of constrained ratios, e.g. c(celltype_1 = 0). Names are matched against colnames(mean_signature_matrix).

epsilon, itmax

Relative convergence tolerance and iteration budget passed to stats::optim().

Value

A list with coefficients (the restricted \(J\) vector), loglik, fixed, and convergence.

Details

Zero mass on \(A\) keeps \(\boldsymbol{\Sigma}(\boldsymbol{p})=\sum_j p_j^{2} \boldsymbol{\Sigma}_j\) positive definite as long as at least one coordinate remains strictly positive. Optimisation uses stats::optim() ("BFGS") with the analytic score obtained by chaining gradient_loglik_unconstrained() through \((1-s)\mathbf{J}_{\boldsymbol{\psi}}\) on the free-face ILR chart.

Examples

mu <- matrix(c(20, 22, 18, 22, 20, 24), nrow = 2)
colnames(mu) <- paste0("ct", 1:3)
Sigma <- array(c(diag(2), diag(2), diag(2)), dim = c(2, 2, 3))
y <- drop(mu %*% c(0.5, 0.5, 0))
restricted_mle_decovart(y, mu, Sigma, fixed = c(ct3 = 0))$coefficients
#> ct1 ct2 ct3 
#> 0.5 0.5 0.0