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 byload_gmt()(a GMT path, gseapy library name, orNonefor 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)andp = 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)andp = 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)andp = 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) ** kto 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) ** kto 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) ** kto 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 mostmax_activeactive sets and integratesalpha/beta/panalytically, 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
pprior 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,
autofalls back to MCMC (andmethod="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
thinsteps 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()whensetsneeds resolving.- deduplicate_terms
bool(default:True) Forwarded to
load_gmt(); ifTrue(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
rusesseed + r, making the whole run reproducible. IfNone, non-deterministic seeds are drawn.
- o
- Return type:
- 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_activeactive sets (over all sets) and integrates the parameters analytically; it is exact up to themax_activecap, andstd_erroris 0 (no Monte-Carlo noise). It is feasible for modest set counts (C(n_sets, <= max_active)configurations); for large collectionsautofalls 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.