Skip to content

crispr_analysis_utils.guide_alignment.pipeline

crispr_analysis_utils.guide_alignment.pipeline

The guide alignment: read a guide table, scan a reference, write the tables.

One pass of the compiled core finds every exact site of every guide, on both strands, and reads the PAM 3' of each and the base 5' of it. There is no aligner, no index and no error budget: a guide is designed, not sequenced, so its own site is an exact match, and everything the tables report follows from that.

outdir gets three folders and nothing else: alignments/ with the results (guide_genome_annotation.tsv, guide_mapping.bed, guide_genome_annotation_columns.yaml, qc.tsv), run_metadata/ with this invocation's record, and logs/.

RESULTS_DIR module-attribute 

RESULTS_DIR = 'alignments'

Where the step's results go, so nothing of its own sits beside the inputs.

METADATA_DIR module-attribute 

METADATA_DIR = 'run_metadata'

Where each invocation leaves its metadata, one file per invocation.

STEP module-attribute 

STEP = 'guide_alignment'

The step those metadata files are named for.

ORIENTATION_CHECKS module-attribute 

ORIENTATION_CHECKS = ('warn', 'error', 'off')

The values of run's orientation_check.

OrientationError

Bases: RuntimeError

The library reads 3' to 5' and orientation_check is error.

Source code in src/crispr_analysis_utils/guide_alignment/pipeline.py 
60
61
class OrientationError(RuntimeError):
    """The library reads 3' to 5' and `orientation_check` is ``error``."""

RunResult dataclass

What run produced.

Attributes:

  • outdir (Path) –

    The output folder.

  • guides (list[Guide]) –

    The library, as returned by read_guides.

  • summaries (list[GuideSummary]) –

    One per guide, in library order, as returned by summarize.

  • scan (ScanResult) –

    The genome scan: each guide's sites, its orientation, and the counts.

  • metadata (dict) –

    The record written under run_metadata/.

Source code in src/crispr_analysis_utils/guide_alignment/pipeline.py 
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
@dataclass
class RunResult:
    """What `run` produced.

    Attributes
    ----------
    outdir
        The output folder.
    guides
        The library, as returned by `read_guides`.
    summaries
        One per guide, in library order, as returned by `summarize`.
    scan
        The genome scan: each guide's sites, its orientation, and the counts.
    metadata
        The record written under ``run_metadata/``.
    """

    outdir: Path
    guides: list[Guide]
    summaries: list[GuideSummary]
    scan: ScanResult
    metadata: dict = field(default_factory=dict)

run 

run(guides: str | PathLike[str] | DataFrame, reference: str | PathLike[str], outdir: str | PathLike[str], *, pam: str = 'NGG', spacer_length: int = 20, add_leading_g: bool = False, orientation_check: str = 'error', fix_orientation: bool = False, contigs: Collection[str] | None = None, id_col: int | str = 0, spacer_col: int | str = 1, sep: str = '\t', header: bool = False) -> RunResult

Find every guide's genomic site and write the tables.

Parameters:

  • guides (str | PathLike[str] | DataFrame) –

    The guide table, as read_guides takes it.

  • reference (str | PathLike[str]) –

    The reference FASTA, uncompressed.

  • outdir (str | PathLike[str]) –

    Where the tables go. Created if it does not exist; nothing is ever written outside it.

  • pam (str, default: 'NGG' ) –

    The PAM, in IUPAC codes, 3' of the protospacer.

  • spacer_length (int, default: 20 ) –

    See read_guides.

  • add_leading_g (bool, default: False ) –

    See read_guides.

  • orientation_check (str, default: 'error' ) –

    error: stop when entries reading 3' to 5' outnumber those reading 5' to 3', which is a table given reverse-complemented throughout. A handful is reported either way, since one guide can look reversed for want of nothing more than a PAM that differs between assemblies. warn: only report. off: skip the check.

  • fix_orientation (bool, default: False ) –

    Align the reverse complement of the entries that read 3' to 5', instead of leaving them as given.

  • contigs (Collection[str] | None, default: None ) –

    If given, read only these contigs.

  • id_col (int | str, default: 0 ) –

    See read_guides.

  • spacer_col (int | str, default: 0 ) –

    See read_guides.

  • sep (int | str, default: 0 ) –

    See read_guides.

  • header (int | str, default: 0 ) –

    See read_guides.

Returns:

  • RunResult –

    The library, the summaries, the scan and this invocation's metadata.

Raises:

  • FileNotFoundError –

    If reference does not exist.

  • ValueError –

    If reference is compressed, pam is not IUPAC codes, orientation_check is not one of ORIENTATION_CHECKS, or the guide table fails read_guides.

  • OrientationError –

    If the library reads 3' to 5' and orientation_check is error.

Source code in src/crispr_analysis_utils/guide_alignment/pipeline.py 
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
def run(
    guides: str | os.PathLike[str] | pd.DataFrame,
    reference: str | os.PathLike[str],
    outdir: str | os.PathLike[str],
    *,
    pam: str = "NGG",
    spacer_length: int = 20,
    add_leading_g: bool = False,
    orientation_check: str = "error",
    fix_orientation: bool = False,
    contigs: Collection[str] | None = None,
    id_col: int | str = 0,
    spacer_col: int | str = 1,
    sep: str = "\t",
    header: bool = False,
) -> RunResult:
    """Find every guide's genomic site and write the tables.

    Parameters
    ----------
    guides
        The guide table, as `read_guides` takes it.
    reference
        The reference FASTA, uncompressed.
    outdir
        Where the tables go. Created if it does not exist; nothing is ever
        written outside it.
    pam
        The PAM, in IUPAC codes, 3' of the protospacer.
    spacer_length
        See `read_guides`.
    add_leading_g
        See `read_guides`.
    orientation_check
        ``error``: stop when entries reading 3' to 5' outnumber those reading
        5' to 3', which is a table given reverse-complemented throughout. A
        handful is reported either way, since one guide can look reversed for
        want of nothing more than a PAM that differs between assemblies.
        ``warn``: only report. ``off``: skip the check.
    fix_orientation
        Align the reverse complement of the entries that read 3' to 5',
        instead of leaving them as given.
    contigs
        If given, read only these contigs.
    id_col, spacer_col, sep, header
        See `read_guides`.

    Returns
    -------
    RunResult
        The library, the summaries, the scan and this invocation's metadata.

    Raises
    ------
    FileNotFoundError
        If `reference` does not exist.
    ValueError
        If `reference` is compressed, `pam` is not IUPAC codes,
        `orientation_check` is not one of `ORIENTATION_CHECKS`, or the guide
        table fails `read_guides`.
    OrientationError
        If the library reads 3' to 5' and `orientation_check` is ``error``.
    """
    if isinstance(contigs, str):
        raise TypeError("contigs must be a collection of contig names, not a string.")
    if orientation_check not in ORIENTATION_CHECKS:
        raise ValueError(
            f"orientation_check must be one of {', '.join(ORIENTATION_CHECKS)}, "
            f"not {orientation_check!r}."
        )
    patterns = validate_pams([pam])
    reference = check_plain_fasta(reference, "reference")

    outdir = Path(outdir)
    outdir.mkdir(parents=True, exist_ok=True)
    (outdir / "logs").mkdir(exist_ok=True)
    (outdir / METADATA_DIR).mkdir(exist_ok=True)
    (outdir / RESULTS_DIR).mkdir(exist_ok=True)
    record = outdir / METADATA_DIR / f"{STEP}_{_stamp()}.json"

    parameters = {
        "pam": patterns[0],
        "spacer_length": spacer_length,
        "add_leading_g": add_leading_g,
        "orientation_check": orientation_check,
        "fix_orientation": fix_orientation,
        "contigs": None if contigs is None else list(contigs),
        "id_col": id_col,
        "spacer_col": spacer_col,
        "sep": sep,
        "header": header,
    }
    state: dict = {
        "step": STEP,
        "complete": False,
        "error": None,
        "started": _now(),
        "finished": None,
        "cau_version": __version__,
        "inputs": {},
        "parameters": parameters,
        "counts": {},
        "classes": {},
        "orientation": {},
        "orientation_fixed": 0,
        "seconds": {},
    }
    started = time.perf_counter()

    with _log_to(outdir / "logs" / LOG_FILE):
        logger.info("Guide alignment into %s, with %s", outdir, parameters)
        try:
            result = _run(guides, reference, outdir, parameters, state, record, started)
        except Exception as error:
            state["error"] = f"{type(error).__name__}: {error}"
            state["finished"] = _now()
            state["seconds"]["total"] = round(time.perf_counter() - started, 3)
            write_json(record, state)
            logger.error("Guide alignment failed: %s", state["error"])
            raise
    return result