grassp.tl.mgsa

Contents

grassp.tl.mgsa#

mgsa(o, sets, population=None, *, alpha_grid=None, beta_grid=None, p_grid=None, alpha_prior=None, beta_prior=None, p_prior=None, method='auto', max_active=4, exact_max_configs=2000000, n_steps=1000000, n_restarts=1, burn_in=None, thin=100, flip_freq=0.8, species='hsap', deduplicate_terms=True, seed=None)[source]#

Run model-based gene set analysis.

Parameters:
o Iterable[str]

The study set: identifiers of the observed (“active”) genes.

sets Union[Mapping[str, Sequence[str]], str, None]

Gene sets as a {name: [genes]} mapping, or anything accepted by load_gmt() (a GMT path, gseapy library name, or None for the grassp-bundled sets).

population Optional[Iterable[str]] (default: None)

The universe of gene identifiers. If None, the union of all set members is used. Study genes outside the population are dropped (and the count is logged).

alpha_grid Optional[Sequence[float]] (default: None)

Discrete grids for the false-positive rate, false-negative rate and set-activation prior. Defaults: alpha = linspace(0.01, 0.3, 10), beta = linspace(0.1, 0.95, 10) and p = linspace(1, min(20, floor(n_sets/3)), 10) / n_sets. The beta ceiling (0.95) is raised from the Bioconductor default (0.8) so a small cluster can activate a large, sparsely-covered compartment.

beta_grid Optional[Sequence[float]] (default: None)

Discrete grids for the false-positive rate, false-negative rate and set-activation prior. Defaults: alpha = linspace(0.01, 0.3, 10), beta = linspace(0.1, 0.95, 10) and p = linspace(1, min(20, floor(n_sets/3)), 10) / n_sets. The beta ceiling (0.95) is raised from the Bioconductor default (0.8) so a small cluster can activate a large, sparsely-covered compartment.

p_grid Optional[Sequence[float]] (default: None)

Discrete grids for the false-positive rate, false-negative rate and set-activation prior. Defaults: alpha = linspace(0.01, 0.3, 10), beta = linspace(0.1, 0.95, 10) and p = linspace(1, min(20, floor(n_sets/3)), 10) / n_sets. The beta ceiling (0.95) is raised from the Bioconductor default (0.8) so a small cluster can activate a large, sparsely-covered compartment.

alpha_prior Optional[callable] (default: None)

Optional priors over the corresponding grids. None (default) is a uniform prior (matching Bioconductor). Otherwise a callable mapping the grid values to non-negative (unnormalised) weights — e.g. beta_prior=lambda b: (1 - b) ** k to penalise high false-negative rates while still allowing them when the likelihood (enough observed hits) outweighs the penalty. Applied to both the exact and MCMC paths.

beta_prior Optional[callable] (default: None)

Optional priors over the corresponding grids. None (default) is a uniform prior (matching Bioconductor). Otherwise a callable mapping the grid values to non-negative (unnormalised) weights — e.g. beta_prior=lambda b: (1 - b) ** k to penalise high false-negative rates while still allowing them when the likelihood (enough observed hits) outweighs the penalty. Applied to both the exact and MCMC paths.

p_prior Optional[callable] (default: None)

Optional priors over the corresponding grids. None (default) is a uniform prior (matching Bioconductor). Otherwise a callable mapping the grid values to non-negative (unnormalised) weights — e.g. beta_prior=lambda b: (1 - b) ** k to penalise high false-negative rates while still allowing them when the likelihood (enough observed hits) outweighs the penalty. Applied to both the exact and MCMC paths.

method str (default: 'auto')

Inference method. "exact" enumerates all configurations with at most max_active active sets and integrates alpha/beta/p analytically, giving variance-free marginals (see Notes). "mcmc" runs the Metropolis-Hastings sampler mirroring the R package. "auto" (default) uses "exact" when the enumeration is feasible (candidate sets pruned to those overlapping the study set; config count <= exact_max_configs) and falls back to "mcmc" otherwise.

max_active int (default: 4)

Exact method: maximum number of simultaneously active sets to enumerate. The posterior is sparse (the p prior penalises many active sets), so a small cap is effectively exact; raise it until marginals stop changing.

exact_max_configs int (default: 2000000)

Exact method: if the number of configurations would exceed this, auto falls back to MCMC (and method="exact" raises).

n_steps int (default: 1000000)

MCMC steps per restart.

n_restarts int (default: 1)

Number of independent chains; posterior estimates are averaged and the standard error is their spread. Use >= 5 for meaningful standard errors.

burn_in Optional[int] (default: None)

Steps discarded before recording. Defaults to n_steps // 2.

thin int (default: 100)

Record one sample every thin steps after burn-in.

flip_freq float (default: 0.8)

Probability of proposing a set-state move (vs. a parameter move).

species str (default: 'hsap')

Forwarded to load_gmt() when sets needs resolving.

deduplicate_terms bool (default: True)

Forwarded to load_gmt(); if True (default) sets with identical gene membership are collapsed to a single term before analysis, so synonymous/duplicate sets do not split the posterior mass.

seed Optional[int] (default: None)

MCMC base RNG seed. Restart r uses seed + r, making the whole run reproducible. If None, non-deterministic seeds are drawn.

Return type:

MgsaResult

Returns:

MgsaResult Posterior summaries for the sets and for alpha/beta/p, the MAP configuration, and diagnostics.

Notes

The exact method enumerates every configuration with up to max_active active sets (over all sets) and integrates the parameters analytically; it is exact up to the max_active cap, and std_error is 0 (no Monte-Carlo noise). It is feasible for modest set counts (C(n_sets, <= max_active) configurations); for large collections auto falls back to MCMC.

The MCMC RNG stream differs from the reference C implementation (a Mersenne Twister with a different call order), so MCMC results agree with the Bioconductor package only up to Monte-Carlo error, not bit-for-bit.