Converge an existing fit, gene by gene

A post-hoc adjustment: takes a fitted fitSpiDE() object and converges every gene to its own penalised negative-binomial optimum, re-estimating the dispersion at the converged mean. This is exactly the stage that fitSpiDE(converge = TRUE) (the default) runs after SpaNorm::fitNB returns; as a separate function it can be applied to a fit made with converge = FALSE, to a fit serialised by an older version of the package, or re-applied with different iteration settings, without refitting. The two routes share one implementation and one penalty vector, so polishSpiDE(fitSpiDE(..., converge = FALSE)) and fitSpiDE(..., converge = TRUE) give the same fit.

Usage

polishSpiDE(object, spe, ...)

# S4 method for class 'SpiDEResults'
polishSpiDE(
  object,
  spe,
  assay = "counts",
  lambda.a = 0,
  maxit = 50L,
  tol = 1e-08,
  block.size = NULL,
  BPPARAM = BiocParallel::SerialParam(),
  verbose = TRUE
)

Arguments

object

a SpiDEResults from fitSpiDE() (one GLM fit per bandwidth). A twoStageSpiDE() result has no per-gene GLM fit and is refused.

spe

the SpatialExperiment the fit was made from (for the counts).

...

further arguments passed to the method.

assay

a character, the counts assay.

lambda.a

the ridge the fit was made with, for a fixed-effects fit (random = "none"), which does not record it; ignored for a mixed fit, which carries its own penalty vector.

maxit, tol

iteration cap and relative log-likelihood tolerance per gene (the converge.maxit / converge.tol of fitSpiDE()).

block.size

genes per block; NULL splits one block per BPPARAM worker.

BPPARAM

a BiocParallelParam; the stage is blocked over genes.

verbose

report progress.

Value

the object with converged alpha and psi, per-gene diagnostics in @polish, and inference cleared.

Details

Why it exists: fitNB fits all genes in one IRLS loop with a shared cell-weight vector and an aggregate convergence criterion, and for bright, cell-type-restricted genes it stops one to four standard errors short of the gene's own optimum, with a dispersion estimated off the optimum that is too large. Converging each gene calibrates the null and, on the toy fixture, sharpens the planted signal (see the model vignette).

Every inference slot derived from the old coefficients (t_stat, se, the combined p-values, the results table and the cross-bandwidth weights) is cleared; call testSpiDE() again afterwards. The reference degrees of freedom (@df) and the variance components are kept, as they are under fitSpiDE(converge = TRUE).

See also

fitSpiDE() for the same stage as a default, testSpiDE() to recompute inference afterwards.

Examples

data(toySpiDE)
spe <- buildNiches(toySpiDE, sigma = 20)
fit0 <- fitSpiDE(spe, condition = "condition", sigma = 20, random = "none",
                 converge = FALSE, verbose = FALSE)
fit1 <- polishSpiDE(fit0, spe, verbose = FALSE)
head(fits(fit1)[[1]]@polish)
#>    iterations psi_fitnb restarted capped singular psi_bound polished
#> G1         13  1.098450     FALSE  FALSE    FALSE     FALSE     TRUE
#> G2         16  1.072685     FALSE  FALSE    FALSE     FALSE     TRUE
#> G3         18  1.071431     FALSE  FALSE    FALSE     FALSE     TRUE
#> G4         11  1.071833     FALSE  FALSE    FALSE     FALSE     TRUE
#> G5         16  1.072373     FALSE  FALSE    FALSE     FALSE     TRUE
#> G6          9  1.071476     FALSE  FALSE    FALSE     FALSE     TRUE