Guide alignment

cau guide-alignment finds where a CRISPR guide library maps in a genome. A guide is designed, not sequenced, so its site is an exact match of the protospacer: the tool scans the reference once and answers three questions per guide — where does the protospacer match, is there a PAM 3' of it, and what sits 5' of it. There is no aligner, no index and no error budget.

A 44,000-guide library against hg38 takes about half a minute.

In Claude Code the cau-guide-alignment skill teaches Claude to run it and read the outputs (Agents and skills).

Quick start

cau guide-alignment run \
    --guides guides.tsv \
    --reference GRCh38.fa \
    --outdir results/

The guide table is two columns by default, a name and a sequence, with no header:

240911_EC_TAPseq_Enhancers_1    GAAAAGCCAACATGAATGCAG
240911_EC_TAPseq_Enhancers_2    GTTTGAGCACATAGCTCTAAG

(the separator above is a tab)

The reference must be an uncompressed FASTA. Nothing is written outside --outdir, and nothing is written beside your inputs.

What the sequences must be

The guide as ordered, 5' to 3', without the PAM. The tool appends --pam (default NGG) itself.

--spacer-length (default 20) says how long the protospacer is. Anything 5' of it — typically the G added for the U6 promoter — is aligned too, but a mismatch there is tolerated and reported rather than counted. That is what makes the two sequence columns of the output worth comparing.

If your table holds bare protospacers and the construct adds a G that is not in the file, pass --add-leading-g.

One spacer in sixteen ends in GG

So a table that silently includes the PAM will run without complaint and be wrong by three bases. Check the column before the first run.

Options

Option	Default	Meaning
--guides, --reference, --outdir	required
--pam	NGG	IUPAC codes, 3' of the protospacer
--spacer-length	20	The protospacer is the 3'-most this many bases
--add-leading-g	off	Put a G in front of every sequence
--orientation-check	error	error, warn or off
--fix-orientation	off	Align the reverse complement of flagged entries
--contigs	every contig	Comma-separated names
--id-col, --spacer-col	0, 1	A 0-based position, or a column name
--sep	tab	tab, comma, or one character
--header	off	The first line holds column names

For hg38, give --contigs the primary chromosomes: without it the alt and random scaffolds are read too, and a guide will match its own alt copy.

Coordinates

0-based and half-open, on the reference's forward strand, and the protospacer alone — the PAM and any base 5' of it are outside [start, end).

target_site is the boundary the blunt SpCas9 cut falls on, three bases from the PAM: the cut lies between target_site - 1 and target_site. On + that is end - 3, on - it is start + 3.

Orientation

Both readings of every guide go into the same scan. An entry is reverse_complement only when its reverse complement has a PAM-bearing site and the sequence as given does not. Reading a CCN 5' of a match is not the test — that falls there by chance one time in sixteen.

A table reversed throughout stops the run. A handful of entries only warns, because a single guide can look reversed when its PAM simply differs between assemblies. orientation_evidence records whether the test was distinct or ambiguous; ambiguous means the guide has nothing 5' of its spacer, so a flip and a chance CCN cannot be told apart.

Outputs

--outdir gets three folders and nothing loose, so a run's output never sits beside its inputs:

alignments/      guide_genome_annotation.tsv  guide_mapping.bed
                 guide_genome_annotation_columns.yaml  qc.tsv
run_metadata/    guide_alignment_<UTC timestamp>.json
logs/            guide_alignment.log

guide_genome_annotation_columns.yaml is written from the same definitions the tables are, so it is always an accurate data dictionary of the folder it sits in. Read it first.

guide_genome_annotation.tsv

One row per library entry per site. A guide matching three places gets three rows; one matching nowhere gets a single row whose site columns are empty. The guide's own columns repeat on each of its rows.

Column	Holds
guide_id	The sequence that was aligned, ending in the PAM pattern. Identifies the guide.
alias	The name the input table gives it.
spacer	The protospacer: guide_id without any base 5' of it or the PAM.
oligo	The sequence the library ordered, without the PAM: the spacer plus any base 5' of it.
n_oligos	How many entries share this spacer. Above one the counts are a sum, because one feature is all the reads can resolve.
site_class	unique, multi or no_site: how many places this guide can cut.
n_matches	How many places the protospacer occurs in the genome.
n_target_sites	How many of those have a PAM.
orientation	forward, reverse_complement, no_pam, no_match or not_checked.
orientation_evidence	For reverse_complement only: distinct or ambiguous.
chr, start, end, strand	This row's match. 0-based, half-open, the protospacer alone.
target_site	The boundary the cut falls on: between target_site - 1 and target_site.
pam_matches	TRUE when the genomic PAM fits the pattern.
mapped_sequence	What the genome holds here, to compare against guide_id.

The two sequence columns

They line up base for base, so one diff answers three questions:

guide_id         GAAAAGCCAACATGAATGCAGNGG     what was looked for
mapped_sequence   AAAAGCCAACATGAATGCAGTGG     what the genome holds
                 ^                     ^^^
                 the 5' G did not pair        the real PAM

mapped_sequence keeps the 5' base only where the genome has that same base, so its length tells you whether the G is genomic or was added for the U6 promoter.

n_matches is not site_class

A guide can occur many times and still be cuttable once. This is the most common surprise:

n_matches 89   n_target_sites 1   site_class unique

The protospacer occurs 89 times; one of those has a PAM. Cas9 cuts there and nowhere else. Every one of the 89 is a row in the table, 88 of them with pam_matches FALSE.

guide_mapping.bed

BED6, no header, one row per cuttable site, so a guide that cuts twice gets two rows. start/end are the protospacer, without the PAM. The name field is the spacer, which is the same on every row of one guide, so rows group by name. Matches with no PAM are not here; they are in the TSV.

qc.tsv

Three columns, section, name, value. The class and orientation sections each sum to the number of entries. scan holds what the pass saw: n_matches, n_off_contig_end, n_no_pam, n_sites.

run_metadata/guide_alignment_<timestamp>.json

One file per invocation, named for the step and when it ran: the parameters, the input checksums, the version, the counts and the timings. Nothing per-guide lives here.

What this does not do

Off-target search. Matching is exact, so there are no mismatch counts. Where a guide can also cut with one or two mismatches is a different question from where it comes from, and it is not answered here.

Bulges, scores and other nucleases. No DNA or RNA bulges, no CFD or MIT scores, and 3' PAMs only.

Coming from the GEM version

The tool used to align guides with GEM3. It now scans the reference directly: a guide is designed, not sequenced, so its site is an exact match of the protospacer, and the only other questions are what sits 3' of it (the PAM) and what sits 5' of it. There is no aligner, no index and no error budget.

What went away

Gone	Why
cau guide-alignment index, --index	Nothing is indexed. Point --reference at the FASTA.
--max-mismatches	Matching is exact. Off-target search is separate work.
--alt-pams	Only --pam is searched.
--threads	One pass, no aligner to parallelise.
alignments/*.bam, *.fastq	No reads are written and nothing is aligned.
GEM3 itself	No external tool, so every environment is native on Apple silicon.

A whole library now takes seconds rather than an hour, and the four GEM defects the old version worked around — the -M 5 search cap, the N of NGG eating the error budget, non-exhaustive search modes and garbage coordinates at contig ends — cannot arise.

Renamed files

Was	Now
guides.tsv and sites.tsv	one guide_genome_annotation.tsv
guides.bed	guide_mapping.bed
target_sites.bed	gone; target_site is a column
orientation.tsv	gone; its verdict and site are columns
summary.tsv	qc.tsv, counts only
run.json	run_metadata/guide_alignment_<timestamp>.json

Everything now lives under alignments/, run_metadata/ and logs/.

Renamed columns

Was	Now
guide_id (the lab's name)	alias
—	guide_id is now the aligned sequence
targeting	gone; use n_target_sites
alignment_class	site_class, and off_target_only no longer exists
n_ngg_<m>mm	n_target_sites, exact matches only
min_ngg_mismatches	gone
leading_g	gone; compare the lengths of guide_id and spacer
pam_class	pam_matches, TRUE/FALSE
site_pam, guide_chr/start/end	pam folded into mapped_sequence; chr, start, end
no_perfect_hit (orientation)	no_match

input_sequence and genomic_protospacer are gone: the first is in your input file and in guide_id, and the second is always identical to spacer when matching is exact.

What this changes in a result

• Guides are keyed by sequence. Two entries differing only by a 5' G are two guides, not duplicates. The old version truncated them to the same spacer and collapsed them.
• n_matches is new and often surprising. The old multi counted only PAM-bearing sites and still does; n_matches counts every occurrence of the protospacer, which is usually larger.
• Orientation is decided from both readings, not from a CCN upstream. A library reversed throughout stops the run; a handful of entries warns.