Deep learning¶

Motifs as convolution filters, sequences as one-hot arrays, and TF-MoDISco results and contribution scores, for working with sequence models. Arrays are numpy arrays in the layouts PyTorch and most sequence models use, so they can be passed to torch.from_numpy(), and PyTorch tensors can be passed back as they are. Reading TF-MoDISco results and chromBPNet contribution scores needs the h5py extra:

pip install universalmotif[h5py]

Attribution and contribution scores from a model, for each position or for each letter at each position, can be drawn along their sequence with view_attributions(). Scores saved by chromBPNet or in the numpy files TF-MoDISco takes are read with read_contributions(), and scores kept as bigWig tracks with signal_values(). The motifs and deep learning guide puts these together on the outputs of a published model.

Reading and drawing scores¶

read_modisco ¶

read_modisco(
    path: str | PathLike[str],
    type: str = "CWM",
    trim_threshold: float | None = 0.3,
    subpatterns: bool = False,
) -> list[Motif]

Read the patterns of a TF-MoDISco results file, as written by modisco-lite's modisco motifs.

type picks which matrix of each pattern to read: "CWM" for its contribution scores, "hCWM" for its hypothetical contribution scores, also as a CWM, or "PPM" for the letter frequencies of its seqlets. Motifs are named after the pattern, such as pos_patterns.pattern_0, with the number of seqlets as nsites. Patterns span a fixed window with weak flanks, so by default they are trimmed as modisco's own reports do: the columns kept run from the first to the last whose absolute contribution reaches trim_threshold times the largest, and the other matrix types are trimmed to the same columns. With subpatterns, the subpatterns of each pattern are read as well.

Needs the h5py extra.

Examples:

>>> import universalmotif as um
>>> cwms = um.read_modisco("modisco_results.h5")
>>> ppms = um.read_modisco("modisco_results.h5", type="PPM")

read_contributions ¶

read_contributions(
    path: str | PathLike[str],
    sequences: Any = None,
    scores: str = "hypothetical",
    regions: Any = None,
    indices: Sequence[int] | None = None,
    alphabet: str = "ACGT",
) -> ContributionScores

Read the contribution scores a sequence model gives its inputs.

Two layouts are read:

The HDF5 files chromBPNet writes, <prefix>.counts_scores.h5 and <prefix>.profile_scores.h5, which hold the one-hot sequences and their scores. scores picks the hypothetical scores of every letter at every position ("hypothetical"), or those times the one-hot sequences ("projected"). The genomic range of each sequence comes from <prefix>.interpreted_regions.bed next to the file, or from regions, a narrowPeak file or a table with chrom, start and summit columns.
numpy .npy or .npz files of scores, as TF-MoDISco takes, with the one-hot sequences given as sequences, a file or an array, or the sequences themselves as strings.

Arrays may have shape (sequences, letters, length) or (sequences, length, letters), with the letters in the order of alphabet. Scores are read as 32-bit floats, and indices reads just some of the sequences of a large HDF5 or .npy file, which is read a part at a time; a .npz file is compressed, so it is read whole.

Needs the h5py extra for HDF5 files.

Examples:

>>> import universalmotif as um
>>> c = um.read_contributions("model.counts_scores.h5")
>>> c.view(0, kind="all")
>>> c = um.read_contributions("shap.npz", sequences="ohe.npz")

ContributionScores `dataclass` ¶

ContributionScores(
    scores: ndarray,
    sequences: list[str],
    names: list[str],
    ranges: Any,
)

Contribution scores of a set of sequences, from read_contributions().

scores is an array of shape (sequences, letters, length) holding a score for every letter at every position, sequences the sequences themselves, names their names, and ranges a DataFrame of the genomic range of each sequence, with positions counted from 0 and the end left out, when the file gives them. Indexing gives the scores and sequence of one sequence, and view() draws them.

Examples:

>>> import numpy as np
>>> import universalmotif as um
>>> c = um.ContributionScores(np.zeros((2, 4, 5)), ["ACGTA", "TTGAC"], ["a", "b"], None)
>>> c
ContributionScores(2 sequences of 5 letters)
>>> scores, seq = c[1]
>>> seq
'TTGAC'

view ¶

view(i: int, kind: str = 'observed', **kwargs: Any)

Draw the scores of sequence i with view_attributions(), with positions on the genome when the ranges are known, so that region is given in genome positions too. Other arguments are passed on.

Sequences and filters¶

one_hot ¶

one_hot(
    sequences: str | Sequence[str] | Mapping[str, str],
    alphabet: str = "ACGT",
    dtype: Any = np.float32,
) -> np.ndarray

Sequences one-hot encoded, as an array of shape (sequences, letters, length).

This is the layout most sequence models take, with one channel per letter of alphabet in the order given. Letters are matched in either case, and letters outside the alphabet, such as N, become columns of zeros. The sequences must all have the same length.

Examples:

>>> import universalmotif as um
>>> x = um.one_hot(["ACGTN", "TTTTT"])
>>> x.shape
(2, 4, 5)
>>> x[0].astype(int).tolist()
[[1, 0, 0, 0, 0], [0, 1, 0, 0, 0], [0, 0, 1, 0, 0], [0, 0, 0, 1, 0]]

from_one_hot ¶

from_one_hot(
    x: Any, alphabet: str = "ACGT", unknown: str = "N"
) -> list[str]

Sequences from one-hot arrays of shape (sequences, letters, length), as made by one_hot(). Each position takes the letter of its largest value, so soft encodings work too, and positions with no positive value become unknown.

Examples:

>>> import universalmotif as um
>>> um.from_one_hot(um.one_hot(["ACGTN", "TTTTT"]))
['ACGTN', 'TTTTT']

one_hot_batches ¶

one_hot_batches(
    genome: Any,
    ranges: Any,
    batch_size: int = 256,
    width: int | None = None,
    alphabet: str = "ACGT",
    dtype: Any = np.float32,
    nthreads: int = 1,
)

Genomic ranges as one-hot batches for a model, read a batch at a time.

genome and ranges are as for extract_ranges(): a Genome read from disk, PackedSequences or a mapping from chromosome names to sequences, and a DataFrame or tuples of ranges counted from 0 with the end left out. Ranges on the - strand are reverse complemented. Models take sequences of one length, so the ranges must all have the same width, or with width each is replaced by a window of that width around its middle, the lower one for even widths. Every range is checked before any is read, and must fit within its chromosome.

Returns an iterable of arrays of shape (sequences, letters, width), up to batch_size sequences each and in the order of the ranges, with a length so that progress can be shown. Only one batch is held in memory at a time, so large sets of regions can be fed to a model, for example as model(torch.from_numpy(x)) for each batch x.

Examples:

>>> import universalmotif as um
>>> genome = {"chr1": "AACCGGTTACGTACGTAAAA"}
>>> batches = um.one_hot_batches(genome, [("chr1", 0, 9), ("chr1", 4, 20, "-")],
...                              width=4, batch_size=1)
>>> len(batches)
2
>>> [um.from_one_hot(x) for x in batches]
[['CCGG'], ['GTAC']]

to_filters ¶

to_filters(
    motifs: Motif | Iterable[Motif],
    use_type: str = "PWM",
    align: str = "left",
    rc: bool = False,
) -> np.ndarray

Motifs as convolution filters, an array of shape (motifs, letters, width).

This is the layout of the weights of a one-dimensional convolution, such as PyTorch's Conv1d with the letters as input channels, so that torch.from_numpy(um.to_filters(motifs)) can be set as a layer's weights, or used to scan one-hot sequences from one_hot(). Filters hold the motifs as use_type, log2 odds by default, and motifs narrower than the widest are padded with zeros, on the right with align="left" or on both sides with "centre". With rc, each motif's reverse complement follows it, for DNA and RNA motifs. All motifs must share an alphabet.

Examples:

>>> import universalmotif as um
>>> filters = um.to_filters(um.example_motifs()[:3])
>>> filters.shape, filters.dtype
((3, 4, 10), dtype('float32'))
>>> bool((filters[2, :, 8:] == 0).all())
True

from_filters ¶

from_filters(
    filters: Any,
    names: Sequence[str] | None = None,
    type: str = "PWM",
    alphabet: str = "DNA",
    bkg: Sequence[float] | None = None,
    trim: bool = True,
) -> list[Motif]

Motifs from convolution filters or matrices from other tools.

filters is an array of shape (motifs, letters, width), such as the weights of a convolution layer, or a mapping from names to (letters, width) matrices, as tangermeme and memelite keep motifs. Tensors are turned into numpy arrays first, so PyTorch tensors can be passed as they are, once on the CPU. type says what the values are, and alphabet gives the letters in row order. With trim, all-zero columns at either end, such as the padding of to_filters(), are removed.

Examples:

>>> import universalmotif as um
>>> motifs = um.example_motifs()[:3]
>>> back = um.from_filters(um.to_filters(motifs), names=[m.name for m in motifs])
>>> [(m.name, m.width, m.type) for m in back]
[('TATA-box', 10, 'PWM'), ('AP-1', 9, 'PWM'), ('CRE', 8, 'PWM')]

Probing a model¶

These make batches of sequences for a model of your own, without depending on any deep-learning library. mutagenesis() makes every single-letter substitution of a sequence, and mutagenesis_effects() turns the model's predictions for them into an array that view_attributions() can draw. insert_sites() writes the same sites into background sequences, and site_spacings() places two sites at a range of distances apart. view_variant() compares the scores of a sequence with those of a variant of it.

mutagenesis ¶

mutagenesis(
    sequence: str,
    alphabet: str = "ACGT",
    region: tuple[int, int] | None = None,
    dtype: Any = np.float32,
) -> np.ndarray

Every single-letter substitution of a sequence, one-hot encoded for a model.

Returns an array of shape (variants, letters, length) whose first entry is the sequence itself. It is followed, for each position of region in turn, by the sequence with each other letter of alphabet at that position, in alphabet order. region is a (start, end) pair counted from 0 with the end left out, and is the whole sequence by default. Positions holding a letter outside the alphabet, such as N, get every letter. A model's predictions for the batch, in the same order, are turned into the effect of each substitution by mutagenesis_effects(). The batch has three variants per position of a DNA sequence, so for long sequences, region keeps it to the part of interest.

Examples:

>>> import universalmotif as um
>>> batch = um.mutagenesis("ACGT")
>>> batch.shape
(13, 4, 4)
>>> um.from_one_hot(batch[:4])
['ACGT', 'CCGT', 'GCGT', 'TCGT']

mutagenesis_effects ¶

mutagenesis_effects(
    predictions: Any,
    sequence: str,
    alphabet: str = "ACGT",
    region: tuple[int, int] | None = None,
    centre: bool = False,
) -> np.ndarray

The effect of each single-letter substitution, from a model's predictions for the batch made by mutagenesis().

predictions holds one prediction per entry of the batch, in order, with the sequence's own first; sequence, alphabet and region must be those the batch was made with. The effect of a substitution is its prediction minus the sequence's, and the letters of the sequence get 0. Returns an array of shape (letters, positions of region), or with several outputs per prediction, given as an array of shape (variants, outputs), one of shape (outputs, letters, positions). With centre, the mean effect at each position is subtracted, as is common before drawing the effects as attributions with view_attributions().

Examples:

>>> import numpy as np
>>> import universalmotif as um
>>> batch = um.mutagenesis("ACGT")
>>> preds = batch[:, 0, :].sum(axis=1)  # a model that counts As
>>> um.mutagenesis_effects(preds, "ACGT").astype(int)
array([[ 0,  1,  1,  1],
       [-1,  0,  0,  0],
       [-1,  0,  0,  0],
       [-1,  0,  0,  0]])

insert_sites ¶

insert_sites(
    sequences: str | Sequence[str] | Mapping[str, str],
    sites: Any,
    positions: int | Sequence[int] | None = None,
)

Write the same sites into every sequence, for measuring a model's response to them.

sites is a site given as a string, a motif, whose most probable letter at each position is used, or a list of these to insert together. positions gives the start of each site, counted from 0; by default a single site is centred in each sequence, rounding down. Sites overwrite the letters where they land, later ones over earlier ones, and must fit within every sequence.

Unlike implant_motifs(), which draws a new site from a motif for every placement, the same letters go into every sequence, as when comparing a model's predictions for background sequences, such as from shuffle_sequences(), with and without a motif. Returns the new sequences in the same form as the input.

Examples:

>>> import universalmotif as um
>>> bkg = ["AAAAAAAAAA", "CCCCCCCCCC"]
>>> um.insert_sites(bkg, "TGACTCA")
['ATGACTCAAA', 'CTGACTCACC']
>>> um.insert_sites(bkg, ["GG", "TT"], positions=[0, 8])
['GGAAAAAATT', 'GGCCCCCCTT']

site_spacings ¶

site_spacings(
    sequences: str | Sequence[str] | Mapping[str, str],
    first: Any,
    second: Any,
    spacings: Iterable[int],
    start: int | None = None,
) -> dict[int, Any]

Two sites at a range of spacings, for measuring how a model responds to the distance between them.

For each spacing, first is written into every sequence at start and second begins that many letters after first ends. Sites are given as for insert_sites(). By default start places the pair in the middle of each sequence at the widest spacing, so that first stays put as second moves. Returns a dict from each spacing to the new sequences, in the same form as the input.

Examples:

>>> import universalmotif as um
>>> res = um.site_spacings("AAAAAAAAAAAA", "GG", "TT", spacings=[0, 2, 4])
>>> res
{0: 'AAGGTTAAAAAA', 2: 'AAGGAATTAAAA', 4: 'AAGGAAAATTAA'}

Deep learning¶

Reading and drawing scores¶

read_modisco ¶

read_contributions ¶

ContributionScores dataclass ¶

view ¶

Sequences and filters¶

one_hot ¶

from_one_hot ¶

one_hot_batches ¶

to_filters ¶

from_filters ¶

Probing a model¶

mutagenesis ¶

mutagenesis_effects ¶

insert_sites ¶

site_spacings ¶

ContributionScores `dataclass` ¶