Skip to content

MotifDiscovery and Motif

The public entry point and the result object. See the usage guide for how the pieces fit together.

pystreme.discovery.MotifDiscovery

Fixed-width, GPU-batched de novo motif discovery.

from_bed + scan are the Phase 1 scanner; enrichment the Phase 2 SEA-equivalent; discover the Phase 3 single-motif finder; fit the Phase 4 round loop -- the motifs = disc.fit(idx) call from §3's API sketch. The constructor arguments are the per-object defaults fit and discover fall back to when not given their own.

from_bed classmethod

from_bed(peaks, genome, size: int, bg_order: int = 2, max_n_frac: float = 0.2, device: str = 'cpu', **kwargs) -> 'MotifDiscovery'

Build a store and immediately fit its background (§2.3, order 2 by default). scan() requires a fitted background, and there's no obviously-right implicit default for a self-contained discovery object the way there might be for a bare scanner call -- see pystreme.scanner.scan's docstring for why that stays explicit there. Remaining keyword arguments (widths, control, holdout_frac, revcomp, seed) are the constructor's defaults.

from_fasta classmethod

from_fasta(path, bg_order: int = 2, max_n_frac: float = 0.2, device: str = 'cpu', **kwargs) -> 'MotifDiscovery'

Equal-length sequences already in a FASTA file (e.g. exported by another tool), no genome/BED needed. Same background fitting and constructor defaults as from_bed.

from_sequences classmethod

from_sequences(seqs: list[str], bg_order: int = 2, device: str = 'cpu', **kwargs) -> 'MotifDiscovery'

Equal-length sequences in memory (a list of ACGTN strings).

scan

scan(pwm, idx=None, revcomp: bool | None = None, kernel: str = 'auto')

Standalone PWM scanner (§"Phase 1" deliverable). See pystreme.scanner.scan.

pwm is a single (4, w) probability matrix, or a list of them for a batched multi-motif scan. A single pwm returns squeezed (n_seqs,) tensors instead of scanner.scan's (1, n_seqs); a list returns the full ScanResult with the motif axis intact.

enrichment

enrichment(pwm, control=None, idx=None, revcomp: bool | None = None, kernel: str = 'auto', n_per_seq: int = 1, bg_order: int | None = None, seed: int | None = None) -> ThresholdResult

SEA-equivalent enrichment test for one known motif (§"Phase 2" deliverable): scan pwm against the primary sequences and against a matched control set, then find the score threshold that most significantly enriches primary over control (§2.5).

control is "dinuc_shuffle" (the constructor default), a control.DinucShuffle/control.GCMatched/control.Explicit instance, or anything control.Explicit accepts directly -- a SequenceStore, list[str], or an index array into self.store (see control.resolve_control). idx restricts which of self.store's sequences count as "primary" (the control set is then built to match that subset, not the whole store). bg_order defaults to whatever from_bed already fit (self.store.bg_order).

discover

discover(widths=None, control=None, idx=None, n_per_width: int = 4, n_iter: int = 20, revcomp: bool | None = None, kernel: str = 'auto', n_per_seq: int = 1, bg_order: int | None = None, seed: int | None = None, min_flank_ic: float = 0.3, width_tolerance: float = 2.0) -> Motif

Single-motif de novo discovery (§"Phase 3" deliverable): find seed words per width (§2.1), refine each width's best seeds (§2.1's batched-convolution refinement loop), and return the single best motif across all widths (see _discover_on for how widths are compared: min_flank_ic trims flat flank columns, width_tolerance breaks near-ties by information content).

This is one round of fit without the hold-out split, erasing, or stopping rules -- one motif, one significance number computed on the same data it was found in (Motif.holdout_logp stays None; see Motif's docstring). Positions and central enrichment are filled in the same way fit does.

fit

fit(primary_idx=None, n_motifs: int = 3, pvalue_thresh: float = 0.05, *, widths=None, control=None, holdout_frac: float | None = None, seed: int | None = None, n_per_width: int = 4, n_iter: int = 20, revcomp: bool | None = None, kernel: str = 'auto', n_per_seq: int = 1, bg_order: int | None = None, patience: int = 3, keep_insignificant: bool = False, min_holdout: int = 10, min_flank_ic: float = 0.3, width_tolerance: float = 2.0, verbose: bool = False) -> list[Motif]

The round loop (§"Phase 4" deliverable): motifs = disc.fit(idx).

  1. Resolve primary_idx (all of self.store by default) and build the control set (control, seeded by seed); fit the background on the control (§2.3).
  2. Split primary and control into training / hold-out subsets (holdout_frac, seeded). Hold-out is skipped, with a warning, when either side would have fewer than min_holdout sequences.
  3. Each round: seeds -> refinement -> best motif on the training set (discover's machinery); re-test its threshold on the hold-out set (Motif.holdout_logp); report sites/positions/ central enrichment on the full primary set; then erase the motif's sites from every working set and go again.
  4. Stop after n_motifs significant motifs (holdout_logp <= log(pvalue_thresh), or the training p-value when there is no hold-out), or after patience consecutive insignificant ones, or when no seed words remain.

Insignificant motifs are dropped unless keep_insignificant (they are still erased before continuing, as in STREME). Motifs come back in discovery order, Motif.round numbered from 1. min_flank_ic and width_tolerance control how the per-round width is settled (see _discover_on).

Hold-out caveat (§5): the hold-out p-value is valid for one call. Calling fit repeatedly over overlapping subsets and picking the best result reintroduces the selection the split was meant to remove -- vary seed per call and keep a never-touched set for a final check if that's the workflow.

to_bed

to_bed(motifs, path, passing_only: bool = True) -> int

Write the motifs' sites as BED6 (genomic coordinates, motif|peak name, log-odds score, strand). See pystreme.results.sites_bed.

to_meme

to_meme(motifs, path, background: ndarray | None = None)

Write motifs to MEME minimal format. See pystreme.meme_io.write_meme.

motifs is a list of Motif (as fit returns; named <round>-<consensus>, with nsites and the hold-out p-value -- or the training one when there is none -- in the E= slot), or the {name: (4, w) array} / list-of-MemeMotif forms write_meme accepts directly.

pystreme.discovery.Motif dataclass

A discovered motif: PWM, consensus, significance, sites, positions.

Two significance numbers, deliberately kept apart (§5's hold-out warning): log_pvalue is the enrichment found on the training sequences -- the same sequences the motif was searched for in, so it's an optimistic, selected number, fine for ranking but not for reporting. holdout_logp re-tests the training threshold on the held-out sequences the search never saw; that's the number fit uses to decide whether a motif counts, and the one to quote. It's None when there was no hold-out set (discover() directly, or fit(holdout_frac=0), or too few sequences to split) -- then log_pvalue is all there is.

sites covers every primary sequence (best score/position/strand, with -inf where no valid window exists); n_sites counts those at or above threshold. positions are the centers of those passing sites relative to the sequence center, in bases (§2.7: fixed-width, peak-centered windows make this free), and central /central_logp is the CentriMo-style central-enrichment test over them.

sites_frame

sites_frame()

Per-peak site table -- see pystreme.results.sites_frame.