Assembly API

Module for full-length TCR sequence assembly.

Overview

The assembly module builds full-length TCR sequences from CDR3 and V/J gene information. It supports:

• Leader sequences: Per-chain configuration - extract from contig FASTAs or use standard signal peptides
• Constant regions: Fetch from Ensembl (TRAC, TRBC1, TRBC2) or use built-in sequences
• 2A linkers: Join alpha and beta chains with self-cleaving peptides (T2A, P2A, E2A, F2A)
• Single-chain constructs: Generate β-linker-α format for expression

Leader Sequence Options

Each chain (alpha and beta) can have its own leader configuration:

Option	Description
None	No leader sequence
"from_contig"	Extract native leader from CellRanger FASTA
"CD8A", "CD28", etc.	Use a standard signal peptide

Default configuration: CD28 on alpha chain, CD8A on beta chain (distinct leaders for identification).

Available Leader Sequences

Leader	Source	Species	Sequence
CD8A	CD8A signal peptide (UniProt P01732)	Human	MALPVTALLLPLALLLHAARP
CD28	CD28 signal peptide (UniProt P10747)	Human	MLRLLLALNLFPSIQVTG
IgK	IgGκ light chain signal peptide	Mouse	METDTLLLWVLLLWVPGSTG
TRAC	TCR alpha constant signal peptide	Human	MAGTWLLLLLALGCPALPTG
TRBC	TCR beta constant signal peptide	Human	MGTSLLCWMALCLLGADHADG

Available Linkers

Linker	Source	Sequence
T2A	Thosea asigna virus	EGRGSLLTCGDVEENPGP
P2A	Porcine teschovirus-1	GSGATNFSLLKQAGDVEENPGP
E2A	Equine rhinitis A virus	QCTNYALLKLAGDVESNPGP
F2A	Foot-and-mouth disease virus	VKQTLNFDLLKLAGDVESNPGP

Usage Examples

from tcrsift import assemble_full_sequences

# Default: CD28 on alpha, CD8A on beta
assembled = assemble_full_sequences(clonotypes)

# Custom leaders
assembled = assemble_full_sequences(
    clonotypes,
    alpha_leader="CD8A",
    beta_leader="CD28",
)

# Leader only on beta chain (first in 2A construct)
assembled = assemble_full_sequences(
    clonotypes,
    alpha_leader=None,
    beta_leader="CD8A",
)

# Extract native leaders from contigs
assembled = assemble_full_sequences(
    clonotypes,
    contigs_dir="/path/to/contigs",
    alpha_leader="from_contig",
    beta_leader="from_contig",
)

# No leaders at all
assembled = assemble_full_sequences(
    clonotypes,
    alpha_leader=None,
    beta_leader=None,
)

API Reference

assemble

Full-length TCR sequence assembly for TCRsift.

Builds complete TCR sequences including leader peptides and constant regions.

CODON_TABLE module-attribute

CODON_TABLE = {'ATA': 'I', 'ATC': 'I', 'ATT': 'I', 'ATG': 'M', 'ACA': 'T', 'ACC': 'T', 'ACG': 'T', 'ACT': 'T', 'AAC': 'N', 'AAT': 'N', 'AAA': 'K', 'AAG': 'K', 'AGC': 'S', 'AGT': 'S', 'AGA': 'R', 'AGG': 'R', 'CTA': 'L', 'CTC': 'L', 'CTG': 'L', 'CTT': 'L', 'CCA': 'P', 'CCC': 'P', 'CCG': 'P', 'CCT': 'P', 'CAC': 'H', 'CAT': 'H', 'CAA': 'Q', 'CAG': 'Q', 'CGA': 'R', 'CGC': 'R', 'CGG': 'R', 'CGT': 'R', 'GTA': 'V', 'GTC': 'V', 'GTG': 'V', 'GTT': 'V', 'GCA': 'A', 'GCC': 'A', 'GCG': 'A', 'GCT': 'A', 'GAC': 'D', 'GAT': 'D', 'GAA': 'E', 'GAG': 'E', 'GGA': 'G', 'GGC': 'G', 'GGG': 'G', 'GGT': 'G', 'TCA': 'S', 'TCC': 'S', 'TCG': 'S', 'TCT': 'S', 'TTC': 'F', 'TTT': 'F', 'TTA': 'L', 'TTG': 'L', 'TAC': 'Y', 'TAT': 'Y', 'TAA': '*', 'TAG': '*', 'TGC': 'C', 'TGT': 'C', 'TGA': '*', 'TGG': 'W'}

LINKERS module-attribute

LINKERS = {'T2A': {'dna': 'GAGGGCAGAGGAAGTCTGCTAACATGCGGTGACGTCGAGGAGAATCCTGGCCCG', 'aa': 'EGRGSLLTCGDVEENPGP', 'source': 'Thosea asigna virus'}, 'P2A': {'dna': 'GGAAGCGGAGCTACTAACTTCAGCCTGCTGAAGCAGGCTGGAGACGTGGAGGAGAACCCTGGACCT', 'aa': 'GSGATNFSLLKQAGDVEENPGP', 'source': 'Porcine teschovirus-1'}, 'E2A': {'dna': 'CAGTGTACTAATTATGCTCTCTTGAAATTGGCTGGAGATGTTGAGAGCAACCCAGGTCCC', 'aa': 'QCTNYALLKLAGDVESNPGP', 'source': 'Equine rhinitis A virus'}, 'F2A': {'dna': 'GTGAAACAGACTTTGAATTTTGACCTTCTCAAGTTGGCGGGAGACGTGGAGTCCAACCCAGGGCCC', 'aa': 'VKQTLNFDLLKLAGDVESNPGP', 'source': 'Foot-and-mouth disease virus'}}

DEFAULT_LEADERS module-attribute

DEFAULT_LEADERS = {'CD8A': {'aa': 'MALPVTALLLPLALLLHAARP', 'dna': 'ATGGCCCTGCCTGTGACAGCCCTGCTGCTGCCTCTGGCTCTGCTGCTGCATGCCGCTAGACCC', 'source': 'Human CD8A signal peptide (UniProt P01732)', 'species': 'human'}, 'CD28': {'aa': 'MLRLLLALNLFPSIQVTG', 'dna': 'ATGCTCCGCCTGCTGCTGGCCCTGAACCTGTTCCCCAGCATCCAGGTGACCGGC', 'source': 'Human CD28 signal peptide (UniProt P10747)', 'species': 'human'}, 'IgK': {'aa': 'METDTLLLWVLLLWVPGSTG', 'dna': 'ATGGAGACAGACACACTCCTGCTATGGGTACTGCTGCTCTGGGTTCCAGGTTCCACTGGT', 'source': 'Murine IgGκ light chain signal peptide', 'species': 'mouse', 'note': 'Widely used for high secretion efficiency in mammalian expression'}, 'TRAC': {'aa': 'MAGTWLLLLLALGCPALPTG', 'dna': 'ATGGCTGGCACCTGGCTGCTGCTGCTGCTGGCCCTGGGATGCCCAGCACTGCCCACAGGC', 'source': 'Human TRAC native signal peptide', 'species': 'human'}, 'TRBC': {'aa': 'MGTSLLCWMALCLLGADHADG', 'dna': 'ATGGGCACCAGCCTGCTGTGCTGGATGGCCCTGTGCCTGCTGGGAGCAGACCACGCCGATGGC', 'source': 'Human TRBC native signal peptide', 'species': 'human'}}

assemble_full_sequences

assemble_full_sequences(clonotypes: DataFrame, contigs_dir: str | Path | None = None, alpha_leader: str | None = 'CD28', beta_leader: str | None = 'CD8A', include_constant: bool = True, constant_source: str = 'ensembl', linker: str = 'T2A', trac_allele: str = 'auto', trbc1_allele: str = 'auto', trbc2_allele: str = 'auto', stop_codons: tuple[str, ...] = ('TAA', 'TGA'), verbose: bool = True, show_progress: bool = True, sample_name_from: str = 'parent', cellranger_dir: str | Path | None = None, sample_sheet: DataFrame | None = None, leader_fallback: str = 'germline', force_alpha_leader: str | None = None, force_beta_leader: str | None = None, secondary_alpha_leader: str | None = None, secondary_beta_leader: str | None = None, codon_optimization: str = 'focal', max_codon_repeats: int = 1, avoid_enzymes: tuple[str, ...] = ()) -> pd.DataFrame

Assemble full-length TCR sequences.

Parameters:

Name	Type	Description	Default
clonotypes	DataFrame	Clonotype DataFrame with VDJ sequences (from fwr1/cdr1/fwr2/cdr2/fwr3/cdr3/fwr4)	required
contigs_dir	str or Path	Directory with CellRanger contig FASTA files. Required if alpha_leader or beta_leader is set to "from_contig".	None
alpha_leader	str or None	Leader sequence for alpha chain. Options: - None: No leader sequence - "from_contig": Extract native leader from contig FASTA (requires contigs_dir) - Key from DEFAULT_LEADERS: "CD8A", "CD28", "IgK", "TRAC", "TRBC" Default is "CD28" to provide distinct sequences from beta chain.	'CD28'
beta_leader	str or None	Leader sequence for beta chain. Same options as alpha_leader. Default is "CD8A" to provide distinct sequences from alpha chain.	'CD8A'
leader_fallback	str or None	Curated signal peptide (CD8A/CD28/IgK/TRAC/TRBC) to substitute when a from_contig leader is implausible — weak-Kozak over-capture, out-of-window length, missing Met/h-region (#263). Default None keeps the contig-extracted leader and leaves its {chain}_leader_qc flag set (visible, not silent). Only affects from_contig leaders.	'germline'
include_constant	bool	Include constant region sequences.	True
constant_source	str	Source for constant regions: "canonical" (default): splice in the hardcoded canonical human TRAC / TRBC1 / TRBC2 (:data:HUMAN_CONSTANT_REGIONS_AA). When CellRanger contigs are available, the canonical's first ~15 residues are verified against the observed contig start and a qc_warnings entry is emitted on mismatch. "ensembl": back-compat alias for "canonical". The earlier pyensembl-backed path was removed in #66 — it read the full mRNA at the wrong frame offset and silently truncated constants to 2–11 residues. "from-data": read {chain}_constant_aa / {chain}_constant_nt directly from the input frame.	'ensembl'
linker	str	Linker sequence for single-chain constructs: "T2A", "P2A", "E2A", "F2A"	'T2A'
trac_allele	str	Per-gene allele selection (#113). Default "auto" runs the allele picker over packaged alleles in :data:HUMAN_CONSTANT_ALLELES scoring each against the donor's contig translation; pass an explicit allele label (e.g. "01" or "03") to force a specific canonical. Invalid labels emit a QC warning and fall back to the default. constant_source="from-data" ignores these knobs (the constant comes from the input frame, not the FASTA); a warning is logged in that combination.	'auto'
trbc1_allele	str	Per-gene allele selection (#113). Default "auto" runs the allele picker over packaged alleles in :data:HUMAN_CONSTANT_ALLELES scoring each against the donor's contig translation; pass an explicit allele label (e.g. "01" or "03") to force a specific canonical. Invalid labels emit a QC warning and fall back to the default. constant_source="from-data" ignores these knobs (the constant comes from the input frame, not the FASTA); a warning is logged in that combination.	'auto'
trbc2_allele	str	Per-gene allele selection (#113). Default "auto" runs the allele picker over packaged alleles in :data:HUMAN_CONSTANT_ALLELES scoring each against the donor's contig translation; pass an explicit allele label (e.g. "01" or "03") to force a specific canonical. Invalid labels emit a QC warning and fall back to the default. constant_source="from-data" ignores these knobs (the constant comes from the input frame, not the FASTA); a warning is logged in that combination.	'auto'
stop_codons	tuple of str	Stop codons appended to the codon-optimized constant CDS. Default uses two non-redundant stops (recognized by different release factors → reduced read-through in synthesized constructs). Pass ("TAA",) for single-stop (pre-2.4 behavior) or () to omit stops entirely. Each entry must be one of "TAA"/"TAG"/"TGA"; invalid entries raise ValueError. Stops are only appended to the _optimized columns and the final _constant_nt (alias) — not to _contig columns, which preserve the donor's CellRanger contig bytes verbatim (and are typically truncated past the contig coverage edge).	``("TAA", "TGA")``
verbose	bool	Print progress information	True
show_progress	bool	Show progress bar	True

Returns:

Type

Description

DataFrame

Clonotypes with full sequences added. Allele audit columns emitted per clone (#113): {chain}_allele_called — e.g. "TRBC2*01" or None when no contig + non-explicit override. {chain}_allele_score — float ∈ [0, 1] (1.0 = perfect AA agreement; 1.0 also when the user forced an allele). {chain}_allele_alternatives — semicolon-joined string of runner-up alleles + scores, sorted by score desc with FASTA order as tie-break. Example: "TRBC2*01:1.000;TRBC2*03:0.933". None when the picker no-decided or only one allele was packaged. Constant NT triad emitted per chain (#116): {chain}_constant_nt_contig — pure CellRanger contig NT past the J→C junction (no canonical splicing, no codon optimization). Truncated where contig coverage ends. None when no contig is available for the clone. {chain}_constant_nt_optimized — codon-optimized canonical CDS for the picked allele (motif-aware via :func:optimize_codons; avoids common restriction sites and homopolymer runs) plus the requested stop codons. {chain}_constant_nt — the assembly-aware blend: uses donor-real bytes from the contig where they agree with the canonical AA, falls back to _optimized for the rest. This is the column most callers want; kept as the default for back-compat. The same triad is also exposed for full_{chain}_nt (leader + VDJ + constant + stop).

Examples:

>>> # Default: CD28 on alpha, CD8A on beta (distinct leaders)
>>> assembled = assemble_full_sequences(clonotypes)

>>> # No leader sequences
>>> assembled = assemble_full_sequences(clonotypes, alpha_leader=None, beta_leader=None)

>>> # Leader only on beta chain (first in 2A construct)
>>> assembled = assemble_full_sequences(clonotypes, alpha_leader=None, beta_leader="CD8A")

>>> # Extract native leaders from contig FASTAs
>>> assembled = assemble_full_sequences(
...     clonotypes,
...     contigs_dir="/path/to/contigs",
...     alpha_leader="from_contig",
...     beta_leader="from_contig",
... )

Source code in tcrsift/assemble.py

def assemble_full_sequences(
    clonotypes: pd.DataFrame,
    contigs_dir: str | Path | None = None,
    alpha_leader: str | None = "CD28",
    beta_leader: str | None = "CD8A",
    include_constant: bool = True,
    constant_source: str = "ensembl",
    linker: str = "T2A",
    trac_allele: str = "auto",
    trbc1_allele: str = "auto",
    trbc2_allele: str = "auto",
    stop_codons: tuple[str, ...] = ("TAA", "TGA"),
    verbose: bool = True,
    show_progress: bool = True,
    sample_name_from: str = "parent",
    cellranger_dir: str | Path | None = None,
    sample_sheet: pd.DataFrame | None = None,
    leader_fallback: str = "germline",
    force_alpha_leader: str | None = None,
    force_beta_leader: str | None = None,
    secondary_alpha_leader: str | None = None,
    secondary_beta_leader: str | None = None,
    codon_optimization: str = "focal",
    max_codon_repeats: int = 1,
    avoid_enzymes: tuple[str, ...] = (),
) -> pd.DataFrame:
    """
    Assemble full-length TCR sequences.

    Parameters
    ----------
    clonotypes : pd.DataFrame
        Clonotype DataFrame with VDJ sequences (from fwr1/cdr1/fwr2/cdr2/fwr3/cdr3/fwr4)
    contigs_dir : str or Path, optional
        Directory with CellRanger contig FASTA files. Required if alpha_leader or
        beta_leader is set to "from_contig".
    alpha_leader : str or None
        Leader sequence for alpha chain. Options:
        - None: No leader sequence
        - "from_contig": Extract native leader from contig FASTA (requires contigs_dir)
        - Key from DEFAULT_LEADERS: "CD8A", "CD28", "IgK", "TRAC", "TRBC"
        Default is "CD28" to provide distinct sequences from beta chain.
    beta_leader : str or None
        Leader sequence for beta chain. Same options as alpha_leader.
        Default is "CD8A" to provide distinct sequences from alpha chain.
    leader_fallback : str or None
        Curated signal peptide (CD8A/CD28/IgK/TRAC/TRBC) to substitute when a
        ``from_contig`` leader is implausible — weak-Kozak over-capture,
        out-of-window length, missing Met/h-region (#263). Default ``None``
        keeps the contig-extracted leader and leaves its ``{chain}_leader_qc``
        flag set (visible, not silent). Only affects ``from_contig`` leaders.
    include_constant : bool
        Include constant region sequences.
    constant_source : str
        Source for constant regions:

        - ``"canonical"`` (default): splice in the hardcoded canonical
          human TRAC / TRBC1 / TRBC2 (:data:`HUMAN_CONSTANT_REGIONS_AA`).
          When CellRanger contigs are available, the canonical's first
          ~15 residues are verified against the observed contig start
          and a ``qc_warnings`` entry is emitted on mismatch.
        - ``"ensembl"``: back-compat alias for ``"canonical"``. The
          earlier pyensembl-backed path was removed in #66 — it read
          the full mRNA at the wrong frame offset and silently
          truncated constants to 2–11 residues.
        - ``"from-data"``: read ``{chain}_constant_aa`` /
          ``{chain}_constant_nt`` directly from the input frame.
    linker : str
        Linker sequence for single-chain constructs: "T2A", "P2A", "E2A", "F2A"
    trac_allele, trbc1_allele, trbc2_allele : str
        Per-gene allele selection (#113). Default ``"auto"`` runs the
        allele picker over packaged alleles in
        :data:`HUMAN_CONSTANT_ALLELES` scoring each against the
        donor's contig translation; pass an explicit allele label
        (e.g. ``"01"`` or ``"03"``) to force a specific canonical.
        Invalid labels emit a QC warning and fall back to the
        default. ``constant_source="from-data"`` ignores these knobs
        (the constant comes from the input frame, not the FASTA);
        a warning is logged in that combination.
    stop_codons : tuple of str, default ``("TAA", "TGA")``
        Stop codons appended to the codon-optimized constant CDS.
        Default uses two non-redundant stops (recognized by
        different release factors → reduced read-through in
        synthesized constructs). Pass ``("TAA",)`` for single-stop
        (pre-2.4 behavior) or ``()`` to omit stops entirely. Each
        entry must be one of ``"TAA"``/``"TAG"``/``"TGA"``;
        invalid entries raise ``ValueError``. Stops are only
        appended to the ``_optimized`` columns and the final
        ``_constant_nt`` (alias) — not to ``_contig`` columns,
        which preserve the donor's CellRanger contig bytes verbatim
        (and are typically truncated past the contig coverage edge).
    verbose : bool
        Print progress information
    show_progress : bool
        Show progress bar

    Returns
    -------
    pd.DataFrame
        Clonotypes with full sequences added. Allele audit columns
        emitted per clone (#113):

        - ``{chain}_allele_called`` — e.g. ``"TRBC2*01"`` or ``None``
          when no contig + non-explicit override.
        - ``{chain}_allele_score`` — float ∈ [0, 1] (1.0 = perfect AA
          agreement; 1.0 also when the user forced an allele).
        - ``{chain}_allele_alternatives`` — semicolon-joined string
          of runner-up alleles + scores, sorted by score desc with
          FASTA order as tie-break. Example: ``"TRBC2*01:1.000;TRBC2*03:0.933"``.
          ``None`` when the picker no-decided or only one allele was
          packaged.

        Constant NT triad emitted per chain (#116):

        - ``{chain}_constant_nt_contig`` — pure CellRanger contig
          NT past the J→C junction (no canonical splicing, no
          codon optimization). Truncated where contig coverage
          ends. ``None`` when no contig is available for the clone.
        - ``{chain}_constant_nt_optimized`` — codon-optimized
          canonical CDS for the picked allele (motif-aware via
          :func:`optimize_codons`; avoids common restriction sites
          and homopolymer runs) plus the requested stop codons.
        - ``{chain}_constant_nt`` — the assembly-aware blend: uses
          donor-real bytes from the contig where they agree with
          the canonical AA, falls back to ``_optimized`` for the
          rest. This is the column most callers want; kept as the
          default for back-compat.

        The same triad is also exposed for ``full_{chain}_nt``
        (leader + VDJ + constant + stop).

    Examples
    --------
    >>> # Default: CD28 on alpha, CD8A on beta (distinct leaders)
    >>> assembled = assemble_full_sequences(clonotypes)

    >>> # No leader sequences
    >>> assembled = assemble_full_sequences(clonotypes, alpha_leader=None, beta_leader=None)

    >>> # Leader only on beta chain (first in 2A construct)
    >>> assembled = assemble_full_sequences(clonotypes, alpha_leader=None, beta_leader="CD8A")

    >>> # Extract native leaders from contig FASTAs
    >>> assembled = assemble_full_sequences(
    ...     clonotypes,
    ...     contigs_dir="/path/to/contigs",
    ...     alpha_leader="from_contig",
    ...     beta_leader="from_contig",
    ... )
    """
    # Validate inputs
    clonotypes = validate_clonotype_df(clonotypes, for_assembly=True)

    # Validate stop_codons early so users get the error before any
    # heavy lifting. Empty tuple is allowed (caller wants no stops).
    stops_nt = stop_codons_nt(stop_codons)
    recode = CodonOptimization(
        mode=codon_optimization,
        max_codon_repeats=max_codon_repeats,
        avoid_enzymes=tuple(avoid_enzymes or ()),
    )

    valid_constant_sources = ["canonical", "ensembl", "from-data"]
    if constant_source not in valid_constant_sources:
        raise TCRsiftValidationError(
            f"Invalid constant_source: '{constant_source}'",
            hint=f"Valid options are: {valid_constant_sources}",
        )

    # #113: warn if the user combined ``constant_source="from-data"``
    # with an explicit allele override — the from-data path reads the
    # constant directly from the input frame and doesn't consult the
    # packaged FASTA, so allele kwargs are silently ignored. Surfacing
    # this loudly prevents a "I set --trbc2-allele 03 but my output is
    # still *01" confusion.
    if constant_source == "from-data" and any(
        a != "auto" for a in (trac_allele, trbc1_allele, trbc2_allele)
    ):
        logger.warning(
            "constant_source='from-data' ignores trac_allele / "
            "trbc1_allele / trbc2_allele — the constant comes from "
            "the input frame, not the packaged canonical FASTA. "
            "Reset the allele kwargs to 'auto' or switch to "
            "constant_source='canonical' to use them."
        )

    # Validate and resolve leader options for each chain
    leader_config = {}
    for chain, leader_param in [("alpha", alpha_leader), ("beta", beta_leader)]:
        if leader_param is None:
            leader_config[chain] = None
        elif leader_param.lower() == "from_contig":
            if not contigs_dir and not cellranger_dir:
                raise TCRsiftValidationError(
                    f"{chain}_leader='from_contig' requires contigs_dir or "
                    "cellranger_dir to be specified",
                    hint="Provide contigs_dir/cellranger_dir with CellRanger FASTA files, or use a default leader like 'CD8A'",
                )
            leader_config[chain] = "from_contig"
        elif leader_param.upper() in DEFAULT_LEADERS:
            leader_config[chain] = DEFAULT_LEADERS[leader_param.upper()]
        else:
            raise TCRsiftValidationError(
                f"Unknown {chain}_leader: '{leader_param}'",
                hint=f"Valid options are: None, 'from_contig', or one of {list(DEFAULT_LEADERS.keys())}",
            )

    if verbose:
        alpha_desc = _describe_leader(alpha_leader, leader_config["alpha"])
        beta_desc = _describe_leader(beta_leader, leader_config["beta"])
        logger.info(f"Assembling full sequences for {len(clonotypes):,} clonotypes")
        logger.info(f"  Alpha leader: {alpha_desc}")
        logger.info(f"  Beta leader: {beta_desc}")
        logger.info(f"  Constant regions: {include_constant} (source: {constant_source})")
        logger.info(f"  Linker: {linker}")

    df = clonotypes.copy()

    # Load constant regions if needed. "canonical" / "ensembl" both
    # hit the same hardcoded-AA path now (#66, #67); "ensembl" is
    # retained as a back-compat alias and emits an info note.
    constant_seqs = {}
    if include_constant and constant_source in ("canonical", "ensembl"):
        if constant_source == "ensembl" and verbose:
            logger.info(
                "  constant_source='ensembl' is now an alias for 'canonical' "
                "(#66) — splicing canonical TRAC / TRBC1 / TRBC2 from "
                "HUMAN_CONSTANT_REGIONS_AA."
            )
        constant_seqs = get_constant_region_sequences(stop_codons=stop_codons)
        if verbose:
            logger.info(
                f"    Loaded {len(constant_seqs)} canonical constant region sequences"
            )
            if stop_codons:
                logger.info(
                    f"    Appending stop codons {list(stop_codons)} to "
                    "codon-optimized constants."
                )
            else:
                logger.info(
                    "    No stop codons configured (stop_codons=()); "
                    "constants will not be terminated."
                )

    # Warn if from-data constants requested but not present
    if include_constant and constant_source == "from-data":
        constant_cols = [
            "alpha_constant_aa",
            "alpha_constant_nt",
            "beta_constant_aa",
            "beta_constant_nt",
        ]
        if not any(col in df.columns for col in constant_cols):
            logger.warning(
                "  constant_source='from-data' but no constant region columns found in input. "
                "Constants will be omitted."
            )

    # Load contigs whenever a contigs_dir is given. Pre-1.3 we only
    # loaded when a leader was "from_contig"; that gate is too narrow
    # now that the C-region NT blend (#103) also consumes contigs to
    # splice donor-real bytes into the J→C boundary. Loading is cheap
    # relative to assembly and avoids a silent "contigs given but
    # blend silently noop'd" failure mode.
    sample_contigs: dict = {}
    # cellranger_dir is shorthand for contigs_dir + sample_name_from='grandparent'
    # (#124). Resolved here so the rest of the function only sees contigs_dir.
    if cellranger_dir is not None:
        if contigs_dir is not None:
            raise ValueError(
                "assemble_full_sequences: pass exactly one of contigs_dir or "
                "cellranger_dir, not both."
            )
        contigs_dir = cellranger_dir
        sample_name_from = "grandparent"
    if contigs_dir:
        contigs_dir = validate_directory_exists(Path(contigs_dir), "contigs directory")
        if verbose:
            logger.info(
                f"  Loading contigs from {contigs_dir} (sample_name_from={sample_name_from!r})..."
            )
        sample_contigs = load_contigs(
            contigs_dir,
            sample_name_from=sample_name_from,
            sample_sheet=sample_sheet,
        )
        if verbose:
            total_contigs = sum(len(c) for c in sample_contigs.values())
            logger.info(f"    Loaded {total_contigs:,} contigs from {len(sample_contigs)} samples")

    # Process each clonotype
    if verbose:
        logger.info("  Assembling sequences...")

    assembly_results = []

    # Build per-gene allele override dict for the picker (#113).
    allele_overrides = {
        "TRAC": trac_allele,
        "TRBC1": trbc1_allele,
        "TRBC2": trbc2_allele,
    }

    # Create iterator with optional progress bar
    row_iter = df.iterrows()
    if show_progress:
        row_iter = tqdm(
            list(df.iterrows()),
            desc="Assembling sequences",
            unit="clone",
        )

    for idx, row in row_iter:
        result = _assemble_clone(
            row,
            sample_contigs,
            constant_seqs,
            leader_config,
            include_constant,
            constant_source,
            allele_overrides=allele_overrides,
            stops_nt=stops_nt,
            recode=recode,
        )
        assembly_results.append(result)

    # Add assembly columns to dataframe
    result_df = pd.DataFrame(assembly_results)
    for col in result_df.columns:
        df[col] = result_df[col].values

    # Add single-chain construct if requested
    if linker and "full_beta_aa" in df.columns and "full_alpha_aa" in df.columns:
        if verbose:
            logger.info(f"  Creating single-chain constructs with {linker} linker...")
        df = _add_single_chain(df, linker)

    # Leader policy (#270): keep an SP-sound, consistent contig leader; warn +
    # switch a bad/inconsistent one to the germline (or configured) leader;
    # honor per-chain force / secondary-construct directives. The default runs
    # on every from_contig assembly (no-op when nothing needs changing).
    df = apply_leader_policy(
        df, linker,
        leader_fallback=leader_fallback or "germline",
        force_alpha=force_alpha_leader, force_beta=force_beta_leader,
        secondary_alpha=secondary_alpha_leader, secondary_beta=secondary_beta_leader,
        verbose=verbose,
    )

    # Summary
    if verbose:
        n_with_alpha = df["full_alpha_aa"].notna().sum() if "full_alpha_aa" in df.columns else 0
        n_with_beta = df["full_beta_aa"].notna().sum() if "full_beta_aa" in df.columns else 0
        n_single_chain = (
            df["single_chain_aa"].notna().sum() if "single_chain_aa" in df.columns else 0
        )
        logger.info("  Assembly complete:")
        logger.info(f"    With full alpha: {n_with_alpha:,}")
        logger.info(f"    With full beta: {n_with_beta:,}")
        logger.info(f"    Single-chain constructs: {n_single_chain:,}")

    # Per-row override audit column (#90). True where the canonical
    # picked by the J-family rule disagreed with CellRanger's raw
    # `beta_c_gene` call. Lets users filter / inspect post-hoc
    # without going through the log. Written even when no overrides
    # fired so the column is always present (defaults to False).
    if "beta_c_gene" in df.columns:
        df["beta_c_gene_overridden"] = _trbc_override_mask(df)

    # One aggregate line for J-family overrides of CellRanger TRBC
    # (#90). pick_canonical_constant is silent per-call to avoid a
    # log flood on cohorts where every other clone hits the case.
    # Gated by `verbose` so silent runs stay silent — the override
    # itself still happens, this is the audit trail. The per-row
    # `beta_c_gene_overridden` column above is the data-side audit
    # that survives the log.
    if verbose:
        n_overrides = (
            int(df["beta_c_gene_overridden"].sum())
            if "beta_c_gene_overridden" in df.columns
            else 0
        )
        if n_overrides:
            n_beta = df["beta_c_gene"].notna().sum() if "beta_c_gene" in df.columns else 0
            logger.warning(
                "  Overrode CellRanger TRBC call on %d / %d β clones to match "
                "J-family parity (TRBJ1→TRBC1, TRBJ2→TRBC2). CellRanger's "
                "TRBC1/2 discrimination is unreliable; the locus rule is "
                "authoritative.",
                n_overrides, n_beta,
            )

        # Cohort allele audit (#119, #120). Emits the per-chain
        # called/no-called breakdown + any novel-allele candidates
        # detected across the cohort. Silent on small cohorts where
        # the heuristic thresholds don't fire.
        if any(
            f"{c}_allele_called_reason" in df.columns for c in ("alpha", "beta")
        ):
            for line in allele_audit_report(df).splitlines():
                logger.info(line)

    # Constant-allele fidelity (#187). Derive a scalar divergence count per
    # chain from the per-position column the assembler recorded against the
    # donor's sequenced bases, and warn loudly (regardless of verbose) so an
    # allele-level fidelity loss isn't silent. Both are no-ops without contigs
    # (the divergence columns stay None / 0).
    for chain in ("alpha", "beta"):
        div_col = f"{chain}_allele_divergence_positions"
        if div_col in df.columns:
            df[f"{chain}_constant_allele_divergence"] = df[div_col].apply(
                _count_divergence_positions
            )
    _warn_constant_allele_divergence(df)

    # Synthesis-hazard QC (#206). Scores each assembled construct NT for
    # gene-synthesis pitfalls (GC window, homopolymers, long repeats, surviving
    # restriction sites) plus cross-construct hazards (duplicate CDS, α/β swap).
    # No-op when no construct-NT column exists. Warn (regardless of verbose) on
    # any flagged construct so a vendor-rejection risk isn't silent.
    # Correct blanket-N α junctions using a per-J consensus from this cohort's
    # contig-verified clones, before QC/provenance (#242).
    df = _apply_cohort_alpha_junctions(df)

    df = add_synthesis_qc(df)

    # Surface dual-α (allelic-inclusion) clones explicitly (#237) — easy to miss
    # when the second α is packed into ``merged_alpha_partners``. ``alpha_count``
    # is the number of distinct α the clone carries; ``dual_alpha`` flags >1.
    if "merged_alpha_partners" in df.columns and len(df):
        def _alpha_partners(val: object) -> list:
            s = "" if val is None or (isinstance(val, float) and pd.isna(val)) else str(val)
            return [p.strip() for p in s.split(";") if p.strip() and p.strip().lower() != "nan"]
        partners = df["merged_alpha_partners"].map(_alpha_partners)
        df["alpha_count"] = partners.map(lambda p: len(p) if p else 1)
        df["dual_alpha"] = df["alpha_count"] > 1
        # Explicit second-α identity (#237): the partner α that isn't this row's
        # primary CDR3_alpha (else the 2nd listed), so the dual nature reads
        # without unpacking merged_alpha_partners. None for single-α clones.
        prim = df["CDR3_alpha"] if "CDR3_alpha" in df.columns else pd.Series([None] * len(df), index=df.index)
        df["CDR3_alpha_2"] = [
            next((a for a in ps if a != p), (ps[1] if len(ps) > 1 else None))
            for ps, p in zip(partners, prim)
        ]

    # Per-construct contig-verification provenance (#243/#244): a construct is
    # contig-verified only when every present chain's J→C junction residue was
    # read from the donor's contig (not the canonical fallback). Drives the
    # fail-closed fidelity gate.
    if len(df):
        df["construct_contig_verified"] = [
            not _unverified_chains(row) for _, row in df.iterrows()
        ]

    if "synth_gc_ok" in df.columns and len(df):
        n_hazard = int(
            (~df["synth_gc_ok"]).sum()
            + (df.get("synth_max_homopolymer", pd.Series(dtype=int)) >= 9).sum()
            + (df.get("synth_max_repeat", pd.Series(dtype=int)) >= 20).sum()
            + (df.get("synth_restriction_sites", pd.Series(dtype=str)) != "").sum()
        )
        if n_hazard:
            logger.warning(
                "  Synthesis-hazard QC flagged %d construct-level issue(s) "
                "(GC window / homopolymer / repeat / restriction site). See the "
                "synth_* columns and synthesis_qc_report() for details.",
                n_hazard,
            )

    # Framework (FR1–FR3) germline divergence (#286 follow-up). Compares each
    # clone's V-domain framework to the IMGT V-REGION germline and records the
    # divergence for the `framework` rows of collect_germline_variants. No-op
    # (adds no columns) unless the V-REGION reference is cached
    # (`tcrsift data download --db imgt_trav_vregion imgt_trbv_vregion`).
    df = annotate_vregion_divergence(df)

    # Schema stability (#278). `qc_warnings` is stashed per-clone only when a
    # warning fires, so the column is absent on a clean cohort and present on a
    # noisy one — 144 vs 145 columns across donors, which breaks naive
    # `df["qc_warnings"]` consumers and cross-cohort concats. Always emit it,
    # with an empty list per row where no warning was raised (a list, not NaN,
    # so the surfacing loop's `isinstance(qcs, list)` and callers' `... or []`
    # both stay correct).
    df = _ensure_qc_warnings_column(df)

    return df

translate_dna

translate_dna(dna_seq: str) -> tuple[str, str]

Translate DNA sequence to amino acids.

Returns:

Type	Description
tuple	(amino_acid_sequence, ragged_3p_nucleotides)

Source code in tcrsift/assemble.py

def translate_dna(dna_seq: str) -> tuple[str, str]:
    """
    Translate DNA sequence to amino acids.

    Returns
    -------
    tuple
        (amino_acid_sequence, ragged_3p_nucleotides)
    """
    seq_len = len(dna_seq)
    seq_len_trimmed = (seq_len // 3) * 3

    if seq_len != seq_len_trimmed:
        ragged_nt = dna_seq[seq_len_trimmed:]
        dna_seq = dna_seq[:seq_len_trimmed]
    else:
        ragged_nt = ""

    aa_seq = "".join([CODON_TABLE.get(dna_seq[i : i + 3], "X") for i in range(0, len(dna_seq), 3)])

    # Stop at first stop codon
    if "*" in aa_seq:
        ragged_nt = ""
        aa_seq = aa_seq[: aa_seq.index("*")]

    return aa_seq, ragged_nt

find_longest_orf

find_longest_orf(dna_seq: str) -> tuple[str, int, str]

Find and translate the longest open reading frame.

Returns:

Type	Description
tuple	(amino_acid_sequence, start_offset, ragged_3p_nucleotides)

Source code in tcrsift/assemble.py

def find_longest_orf(dna_seq: str) -> tuple[str, int, str]:
    """
    Find and translate the longest open reading frame.

    Returns
    -------
    tuple
        (amino_acid_sequence, start_offset, ragged_3p_nucleotides)
    """
    start_positions = [i for i in range(len(dna_seq)) if dna_seq[i : i + 3] == "ATG"]

    longest_aa = ""
    longest_offset = 0
    longest_ragged = ""

    for start in start_positions:
        subseq = dna_seq[start:]
        aa, ragged = translate_dna(subseq)
        if len(aa) > len(longest_aa):
            longest_aa = aa
            longest_offset = start
            longest_ragged = ragged

    return longest_aa, longest_offset, longest_ragged

parse_fasta

parse_fasta(path: str | Path) -> dict[str, str]

Parse a FASTA file.

Returns:

Type	Description
dict	Mapping from sequence ID to sequence

Source code in tcrsift/assemble.py

def parse_fasta(path: str | Path) -> dict[str, str]:
    """
    Parse a FASTA file.

    Returns
    -------
    dict
        Mapping from sequence ID to sequence
    """
    path = Path(path)
    results = {}
    curr_id = None
    lines = []

    with open(path) as f:
        for line in f:
            line = line.strip()
            if line.startswith(">"):
                if curr_id and lines:
                    results[curr_id] = "".join(lines)
                    lines = []
                curr_id = line[1:].split()[0]  # Take first word after >
            else:
                lines.append(line)

        # Don't forget last entry
        if curr_id and lines:
            results[curr_id] = "".join(lines)

    return results

load_contigs

load_contigs(contig_dir: str | Path | None = None, *, sample_name_from: str = 'parent', cellranger_dir: str | Path | None = None, sample_sheet: DataFrame | None = None) -> dict[str, dict[str, str]]

Load contig sequences from a CellRanger-style output tree (#124).

Parameters:

Name	Type	Description	Default
contig_dir	str or Path	Root directory to scan for *contig*.fasta. Mutually exclusive with cellranger_dir.	None
sample_name_from	('parent', 'grandparent', 'sheet')	How to derive the sample name from each discovered FASTA's path: 'parent' (default, backward-compat) — sample = FASTA's immediate parent directory name. Matches the symlinked layout contig_dir/{sample}/filtered_contig.fasta. 'grandparent' — sample = FASTA's grandparent directory name. Matches CellRanger multi output: per_sample_outs/{sample}/vdj_t/filtered_contig.fasta. 'sheet' — match each FASTA against sample_sheet['vdj_dir'] and take the corresponding sample column. Most explicit; useful for non-standard layouts.	'parent'
cellranger_dir	str or Path	Shorthand for contig_dir=cellranger_dir, sample_name_from='grandparent'. Pass a raw CellRanger per_sample_outs/ directory and the sample names will resolve correctly without a symlink tree. Mutually exclusive with contig_dir.	None
sample_sheet	DataFrame	Required when sample_name_from='sheet'. Must have columns sample and vdj_dir; each FASTA's path must be inside one of the listed vdj_dir paths.	None

Returns:

Type	Description
dict	Nested dict: sample -> contig_id -> sequence.

Source code in tcrsift/assemble.py

def load_contigs(
    contig_dir: str | Path | None = None,
    *,
    sample_name_from: str = "parent",
    cellranger_dir: str | Path | None = None,
    sample_sheet: pd.DataFrame | None = None,
) -> dict[str, dict[str, str]]:
    """
    Load contig sequences from a CellRanger-style output tree (#124).

    Parameters
    ----------
    contig_dir : str or Path, optional
        Root directory to scan for ``*contig*.fasta``. Mutually exclusive
        with ``cellranger_dir``.
    sample_name_from : {'parent', 'grandparent', 'sheet'}, default 'parent'
        How to derive the sample name from each discovered FASTA's path:

        - ``'parent'`` (default, backward-compat) — sample = FASTA's
          immediate parent directory name. Matches the symlinked layout
          ``contig_dir/{sample}/filtered_contig.fasta``.
        - ``'grandparent'`` — sample = FASTA's grandparent directory name.
          Matches CellRanger ``multi`` output:
          ``per_sample_outs/{sample}/vdj_t/filtered_contig.fasta``.
        - ``'sheet'`` — match each FASTA against ``sample_sheet['vdj_dir']``
          and take the corresponding ``sample`` column. Most explicit;
          useful for non-standard layouts.
    cellranger_dir : str or Path, optional
        Shorthand for ``contig_dir=cellranger_dir, sample_name_from='grandparent'``.
        Pass a raw CellRanger ``per_sample_outs/`` directory and the sample
        names will resolve correctly without a symlink tree. Mutually
        exclusive with ``contig_dir``.
    sample_sheet : pandas.DataFrame, optional
        Required when ``sample_name_from='sheet'``. Must have columns
        ``sample`` and ``vdj_dir``; each FASTA's path must be inside one of
        the listed ``vdj_dir`` paths.

    Returns
    -------
    dict
        Nested dict: ``sample -> contig_id -> sequence``.
    """
    if cellranger_dir is not None:
        if contig_dir is not None:
            raise ValueError(
                "load_contigs: pass exactly one of contig_dir or cellranger_dir, not both."
            )
        contig_dir = cellranger_dir
        sample_name_from = "grandparent"
    if contig_dir is None:
        raise ValueError(
            "load_contigs: must pass either contig_dir or cellranger_dir."
        )
    if sample_name_from not in SAMPLE_NAME_FROM_CHOICES:
        raise ValueError(
            f"load_contigs: sample_name_from={sample_name_from!r} not in "
            f"{SAMPLE_NAME_FROM_CHOICES}"
        )
    if sample_name_from == "sheet" and sample_sheet is None:
        raise ValueError(
            "load_contigs: sample_name_from='sheet' requires sample_sheet=..."
        )

    contig_dir = Path(contig_dir)
    sample_contigs: dict[str, dict[str, str]] = {}

    sheet_dir_to_sample: dict[Path, str] = {}
    if sample_name_from == "sheet":
        missing = {"sample", "vdj_dir"} - set(sample_sheet.columns)
        if missing:
            raise ValueError(
                f"load_contigs: sample_sheet is missing required columns {sorted(missing)}"
            )
        sheet_dir_to_sample = {
            Path(row.vdj_dir).resolve(): row.sample
            for row in sample_sheet.itertuples(index=False)
        }

    def _sample_for(fasta_path: Path) -> str | None:
        if sample_name_from == "parent":
            return fasta_path.parent.name
        if sample_name_from == "grandparent":
            return fasta_path.parent.parent.name
        # sheet
        resolved = fasta_path.resolve()
        for sheet_path, sample in sheet_dir_to_sample.items():
            try:
                resolved.relative_to(sheet_path)
            except ValueError:
                continue
            return sample
        return None

    for fasta_path in contig_dir.rglob("*contig*.fasta"):
        sample_name = _sample_for(fasta_path)
        if sample_name is None:
            logger.warning(
                "load_contigs: %s did not match any sample_sheet vdj_dir; skipping",
                fasta_path,
            )
            continue
        sample_contigs.setdefault(sample_name, {}).update(parse_fasta(fasta_path))

    # Flat-layout fallback: bare ``*.fasta`` files directly in contig_dir.
    # Only meaningful for the 'parent' default — the grandparent/sheet
    # modes assume the CellRanger-style nested layout where this wouldn't
    # match anything.
    if sample_name_from == "parent":
        for fasta_path in contig_dir.glob("*.fasta"):
            sample_name = fasta_path.stem.split("_")[0]
            sample_contigs.setdefault(sample_name, {}).update(parse_fasta(fasta_path))

    logger.info(f"Loaded contigs from {len(sample_contigs)} samples")
    return sample_contigs

get_constant_region_sequences

get_constant_region_sequences(*, stop_codons: tuple[str, ...] = ('TAA', 'TGA')) -> dict[str, str]

Return human TCR constant-region CDS (DNA) via codon-aware back- translation.

Sources NT from :data:HUMAN_CONSTANT_REGIONS_AA and :func:optimize_codons (motif-aware — avoids 5+ mononucleotide runs and common Type-II restriction sites; see #116). The earlier pyensembl-backed implementation read the full mRNA at frame offset 2 and silently truncated TRAC / TRBC1 / TRBC2 to 2–11 residues for every assembled clonotype (#66). Hardcoding eliminates the frame bug, drops the pyensembl dependency, and matches the canonical sequences that downstream cloning constructs need (#67).

Parameters:

Name	Type	Description	Default
stop_codons	tuple of str	Stop codons to append to each CDS. Default is the two non-redundant stops (different release factors → reduces read-through in synthesized constructs).	``("TAA", "TGA")``

Returns:

Type	Description
dict	Gene name → CDS (DNA, post-junction-codon-aligned … stop). Each CDS is the codon-optimized translation of the canonical AA plus the requested stop codons.

Source code in tcrsift/assemble.py

def get_constant_region_sequences(
    *,
    stop_codons: tuple[str, ...] = ("TAA", "TGA"),
) -> dict[str, str]:
    """
    Return human TCR constant-region CDS (DNA) via codon-aware back-
    translation.

    Sources NT from :data:`HUMAN_CONSTANT_REGIONS_AA` and
    :func:`optimize_codons` (motif-aware — avoids 5+ mononucleotide
    runs and common Type-II restriction sites; see #116). The earlier
    pyensembl-backed implementation read the full mRNA at frame
    offset 2 and silently truncated TRAC / TRBC1 / TRBC2 to 2–11
    residues for every assembled clonotype (#66). Hardcoding
    eliminates the frame bug, drops the pyensembl dependency, and
    matches the canonical sequences that downstream cloning
    constructs need (#67).

    Parameters
    ----------
    stop_codons : tuple of str, default ``("TAA", "TGA")``
        Stop codons to append to each CDS. Default is the two
        non-redundant stops (different release factors → reduces
        read-through in synthesized constructs).

    Returns
    -------
    dict
        Gene name → CDS (DNA, post-junction-codon-aligned … stop).
        Each CDS is the codon-optimized translation of the canonical
        AA plus the requested stop codons.
    """
    stops_nt = stop_codons_nt(stop_codons)
    out: dict[str, str] = {}
    for name, aa in HUMAN_CONSTANT_REGIONS_AA.items():
        out[name] = optimize_codons(aa) + stops_nt
    return out

validate_sequences

validate_sequences(df: DataFrame, strict: bool = False, fix: bool = False) -> list[ValidationMessage]

Validate assembled sequences end-to-end per row (#67, #68).

Per-chain checks on each full_{chain}_aa row:

• Length window (200–450 aa).
• CDR3 substring present in the assembled sequence.
• Constant region's canonical C-terminus (from :data:CONSTANT_REGION_ENDINGS).
• Constant region's canonical start (first 8+ residues of :data:HUMAN_CONSTANT_REGIONS_AA).
• Constant length floor (:data:CONSTANT_AA_FLOOR).
• No premature stop codon (* mid-chain).
• Methionine start (when a leader was included).
• Standard residue alphabet only (:data:_VALID_AA_CHARS).
• Byte-for-byte equality of leader + vdj + constant and the assembled full_{chain}_aa, when all three parts are present.

Cross-row checks:

• β-chain J→C in-cis parity: TRBJ1* pairs with TRBC1, TRBJ2* with TRBC2.
• Single-chain integrity: linker present in single_chain_aa.
• Per-row qc_warnings (stashed by _add_constant_regions) are surfaced.

Returns a flat list of :class:ValidationMessage (a str subclass with .idx and .severity attributes). When strict is True, raises :class:TCRsiftValidationError if any load-bearing checks fail. Informational notes ("didn't have enough info to check") are returned but never raise — distinguishing those from a silent pass was the gap that hid #66.

When fix is True, autocorrect β J→C parity mismatches in-place via :func:fix_jc_parity before validating (#89). Autocorrection notes are prepended to the returned messages and the input frame is mutated unconditionally — the mutation survives even when a later check fails and strict=True raises.

Source code in tcrsift/assemble.py

def validate_sequences(
    df: pd.DataFrame, strict: bool = False, fix: bool = False
) -> list[ValidationMessage]:
    """Validate assembled sequences end-to-end per row (#67, #68).

    Per-chain checks on each ``full_{chain}_aa`` row:

    - Length window (200–450 aa).
    - CDR3 substring present in the assembled sequence.
    - Constant region's canonical C-terminus (from
      :data:`CONSTANT_REGION_ENDINGS`).
    - Constant region's canonical start (first 8+ residues of
      :data:`HUMAN_CONSTANT_REGIONS_AA`).
    - Constant length floor (:data:`CONSTANT_AA_FLOOR`).
    - No premature stop codon (``*`` mid-chain).
    - Methionine start (when a leader was included).
    - Standard residue alphabet only (:data:`_VALID_AA_CHARS`).
    - Byte-for-byte equality of ``leader + vdj + constant`` and the
      assembled ``full_{chain}_aa``, when all three parts are present.

    Cross-row checks:

    - β-chain J→C in-cis parity: ``TRBJ1*`` pairs with TRBC1,
      ``TRBJ2*`` with TRBC2.
    - Single-chain integrity: linker present in ``single_chain_aa``.
    - Per-row ``qc_warnings`` (stashed by ``_add_constant_regions``)
      are surfaced.

    Returns a flat list of :class:`ValidationMessage` (a ``str``
    subclass with ``.idx`` and ``.severity`` attributes). When
    ``strict`` is True, raises :class:`TCRsiftValidationError` if any
    load-bearing checks fail. Informational notes ("didn't have enough
    info to check") are returned but never raise — distinguishing
    those from a silent pass was the gap that hid #66.

    When ``fix`` is True, autocorrect β J→C parity mismatches in-place
    via :func:`fix_jc_parity` before validating (#89). Autocorrection
    notes are prepended to the returned messages and the input frame
    is mutated *unconditionally* — the mutation survives even when a
    later check fails and ``strict=True`` raises.
    """
    load_bearing: list[ValidationMessage] = []
    informational: list[ValidationMessage] = []
    autocorrect_notes: list[ValidationMessage] = []

    def _lb(idx, text):
        load_bearing.append(
            ValidationMessage(f"Clone {idx}: {text}", idx=idx, severity="load_bearing")
        )

    def _info(idx, text):
        informational.append(
            ValidationMessage(f"Clone {idx}: {text}", idx=idx, severity="informational")
        )

    if fix:
        autocorrect_notes = fix_jc_parity(df)

    # 1. Per-chain, per-row shape and content checks.
    for chain in ["alpha", "beta"]:
        col = f"full_{chain}_aa"
        if col not in df.columns:
            continue
        for idx, row in df.iterrows():
            seq = row.get(col, "")
            if not seq or not isinstance(seq, str):
                continue

            if len(seq) < 200:
                _lb(idx, f"{chain} chain too short ({len(seq)} aa)")
            elif len(seq) > 450:
                _lb(idx, f"{chain} chain too long ({len(seq)} aa)")

            # Standard residue alphabet only.
            invalid = sorted(set(seq) - _VALID_AA_CHARS)
            if invalid:
                _lb(idx,
                    f"full_{chain}_aa has invalid residues {''.join(invalid)!r}")

            # Premature stop codon: ``*`` anywhere except the final position.
            body = seq.rstrip("*")
            if "*" in body:
                pos = body.index("*")
                _lb(idx,
                    f"full_{chain}_aa has premature stop "
                    f"at position {pos} of {len(seq)}")

            # Methionine start (only when leader was included).
            leader = row.get(f"{chain}_leader_aa")
            if isinstance(leader, str) and leader and not seq.startswith("M"):
                _lb(idx,
                    f"full_{chain}_aa doesn't start with M "
                    f"(starts with {seq[:5]!r}; leader present)")

            # Surface a from_contig signal-peptide QC flag (#263) as an
            # informational note — the extracted leader is kept (it's the best
            # available), but a weak-Kozak over-capture / out-of-range length /
            # missing h-region shouldn't pass silently. Never load-bearing here:
            # the assembler already substitutes a curated SP when leader_fallback
            # is set; this just keeps the flag visible in a validation pass.
            lqc = row.get(f"{chain}_leader_qc")
            if isinstance(lqc, str) and lqc not in ("ok", "curated_fallback"):
                support = row.get(f"{chain}_leader_support")
                _info(idx,
                    f"{chain} signal peptide QC={lqc} "
                    f"(len {row.get(f'{chain}_leader_len')}, "
                    f"kozak {row.get(f'{chain}_leader_kozak_score')}, "
                    f"support {support})")

            cdr3_col = f"CDR3_{chain}"
            if cdr3_col in row:
                cdr3 = row[cdr3_col]
                if isinstance(cdr3, str) and cdr3 and cdr3 not in seq:
                    _lb(idx,
                        f"CDR3_{chain}={cdr3!r} not found in full sequence")

            # Constant region length floor (separate from full-chain
            # length check — catches truncations that hid behind a
            # long leader).
            const = row.get(f"{chain}_constant_aa")
            if isinstance(const, str) and const:
                floor = CONSTANT_AA_FLOOR[chain]
                if len(const) < floor:
                    _lb(idx,
                        f"{chain}_constant_aa too short "
                        f"({len(const)} aa, floor {floor})")

            # J→C junction-seam guard (#235). The constant must begin with a
            # valid junction residue (β: E; α: N/Y/D/H) followed by the bare
            # canonical — NOT the bare canonical directly. A direct-canonical
            # start means the junction residue was dropped and the chain is 1 aa
            # short at the seam. This is the check the prior self-consistent
            # canonical-start comparison (observed vs the row's own constant_aa)
            # could not make, which is how #235 shipped silently.
            c_base = _resolve_c_gene(row, chain).split("*")[0]
            bare = HUMAN_CONSTANT_REGIONS_AA.get(c_base)
            if isinstance(const, str) and const and bare:
                bare8 = bare[:8]
                if const.startswith(bare8):
                    _lb(idx,
                        f"{chain}_constant_aa starts at the bare canonical "
                        f"{bare8!r} — J→C junction residue missing (chain 1 aa "
                        f"short at the seam)")
                elif not (
                    const[0] in _VALID_JUNCTION_RESIDUES.get(chain, frozenset())
                    and const[1:].startswith(bare8)
                ):
                    _lb(idx,
                        f"{chain}_constant_aa has an unexpected J→C seam: "
                        f"{const[:9]!r} (expected a valid junction residue "
                        f"{sorted(_VALID_JUNCTION_RESIDUES.get(chain, ()))} "
                        f"+ canonical {bare8!r})")
                # α J-gene cross-check (#242/#236): the seam residue must match
                # the J gene's germline residue. Catches a wrongly-defaulted N
                # (or a contig misread) that the syntactic {N,Y,D,H} check above
                # can't — a defaulted N for a Y/D/H J gene is ~30% of α.
                elif chain == "alpha":
                    expected = _traj_junction_residue(row.get("alpha_j_gene"))
                    if expected and const[0] != expected:
                        _lb(idx,
                            f"alpha J→C junction {const[0]!r} disagrees with the "
                            f"germline residue {expected!r} for "
                            f"{str(row.get('alpha_j_gene'))!r} — likely a "
                            f"defaulted/misread junction (#242)")

            # Byte-for-byte: full == leader + vdj + constant when all
            # parts are available. Catches dropped/added residues
            # during assembly (the failure mode no other check finds).
            vdj = row.get(f"vdj_{chain}_aa")
            if (
                isinstance(leader, str) and leader
                and isinstance(vdj, str) and vdj
                and isinstance(const, str) and const
            ):
                expected = leader + vdj + const
                if expected != seq:
                    _lb(idx,
                        f"full_{chain}_aa != leader+vdj+constant "
                        f"(len {len(seq)} vs {len(expected)})")

            # Canonical-ending / canonical-start checks. Prefer
            # ``{chain}_c_gene_canonical`` (written by
            # ``_add_constant_regions``); fall back to the raw
            # CellRanger ``{chain}_c_gene`` call. If neither resolves
            # to a known canonical gene name, emit an informational
            # "did-not-check" note so the user knows validation was
            # incomplete.
            c_gene = _resolve_c_gene(row, chain)
            # Strip any allele suffix (e.g. "TRBC1*01" → "TRBC1").
            c_gene_base = c_gene.split("*")[0]
            if c_gene_base in CONSTANT_REGION_ENDINGS:
                expected_end = CONSTANT_REGION_ENDINGS[c_gene_base]
                if not seq.endswith(expected_end):
                    _lb(idx,
                        f"{chain} doesn't end with canonical "
                        f"{c_gene_base} C-terminus ({expected_end!r}); "
                        f"got {seq[-len(expected_end):]!r}")
                # Prefer the row's actual ``{chain}_constant_aa`` as the
                # expected start. After 2.0 (#105) this carries the
                # per-clone junction residue prepend (e.g.
                # ``NIQNPDPAVY...`` for α with N junction;
                # ``EDLNKVFP...`` for β with universal-E junction).
                # Falling back to ``HUMAN_CONSTANT_REGIONS_AA[base]``
                # — the bare-mature canonical — was the #107 regression:
                # observed includes the junction; expected didn't.
                expected_const = row.get(f"{chain}_constant_aa")
                if not isinstance(expected_const, str) or not expected_const:
                    expected_const = HUMAN_CONSTANT_REGIONS_AA.get(c_gene_base)
                if expected_const and not verify_canonical_constant_start(
                    _expected_constant_start_from_full(seq, row, chain),
                    expected_const,
                    min_match=8,
                ):
                    _lb(idx,
                        f"{chain} constant start doesn't match canonical "
                        f"{c_gene_base} (expected start "
                        f"{expected_const[:15]!r})")
            else:
                _info(idx,
                    f"{chain} c_gene={c_gene!r} not in "
                    "CONSTANT_REGION_ENDINGS — canonical ending/start unverifiable")

    # 2. β-chain J→C in-cis parity disagreement (#67).
    if "beta_j_gene" in df.columns:
        for idx, row in df.iterrows():
            j_gene = row.get("beta_j_gene", "")
            c_gene = _resolve_c_gene(row, "beta")
            j_family = _beta_j_family(j_gene)
            c_base = _beta_c_base(c_gene)
            expected_c = BETA_JC_PARITY.get(j_family) if j_family else None
            if expected_c and c_base and c_base != expected_c:
                _lb(idx,
                    f"β J→C parity mismatch — "
                    f"{j_gene} (family {j_family}) should pair with "
                    f"{expected_c}, got {c_base}")

    # 3. Single-chain (β-linker-α) construct integrity. Three checks:
    #    a) the linker appears in ``single_chain_aa`` exactly once,
    #    b) ``single_chain_aa == β.rstrip('*') + linker + α`` byte-for-byte
    #       (catches dropped residues or wrong concatenation order),
    #    c) if the linker AA matches a known 2A peptide name, it should
    #       be the canonical sequence from :data:`LINKERS`.
    if "single_chain_aa" in df.columns and "linker" in df.columns:
        for idx, row in df.iterrows():
            sc = row.get("single_chain_aa")
            linker = row.get("linker")
            if not isinstance(sc, str) or not sc:
                continue
            if not isinstance(linker, str) or not linker:
                continue
            count = sc.count(linker)
            if count == 0:
                _lb(idx, f"single_chain_aa missing linker {linker!r}")
            elif count > 1:
                _lb(idx,
                    f"single_chain_aa contains linker "
                    f"{linker!r} {count} times (expected 1)")

            beta = row.get("full_beta_aa")
            alpha = row.get("full_alpha_aa")
            if (
                isinstance(beta, str) and beta
                and isinstance(alpha, str) and alpha
                and count == 1
            ):
                expected_sc = beta.rstrip("*") + linker + alpha
                if sc != expected_sc:
                    _lb(idx,
                        f"single_chain_aa != β+linker+α "
                        f"(len {len(sc)} vs {len(expected_sc)})")

    # 4. NT → AA round-trip invariants (#91). Every NT column the
    # assembler emits must back-translate to its corresponding AA
    # column. This is the integration check that catches splice-time
    # frame bugs — the original #91 (+1 trailing nt in VDJ_*_nt was
    # not trimmed before concatenation, frame-shifting everything
    # past the VDJ→C boundary) would have been caught here on day one.
    _validate_nt_aa_roundtrip(df, _lb)

    # 5. Surface any per-row qc_warnings the assembler stashed. These are
    # *informational, self-corrected* notes (fell back to canonical, allele
    # not called, contig diverged at the C-region boundary, …) — the
    # assembler already handled them and the assembled sequence is valid, so
    # they must NOT be load-bearing (that aborted `tcrsift run` on real
    # CellRanger data, #129). Genuine structural failures are caught by the
    # explicit checks above and remain load-bearing.
    if "qc_warnings" in df.columns:
        for idx, qcs in df["qc_warnings"].items():
            if isinstance(qcs, list):
                for msg in qcs:
                    _info(idx, str(msg))

    all_messages = autocorrect_notes + load_bearing + informational
    if strict and load_bearing:
        raise TCRsiftValidationError(
            f"Sequence validation failed ({len(load_bearing)} load-bearing issues):\n  "
            + "\n  ".join(load_bearing[:10])
            + ("\n  ..." if len(load_bearing) > 10 else "")
        )
    return all_messages

export_fasta

export_fasta(df: DataFrame, output_path: str | Path, sequence_col: str = 'single_chain_aa')

Export sequences to FASTA format.

Parameters:

Name	Type	Description	Default
df	DataFrame	DataFrame with sequences	required
output_path	str or Path	Output file path	required
sequence_col	str	Column containing sequences to export	'single_chain_aa'

Source code in tcrsift/assemble.py

def export_fasta(df: pd.DataFrame, output_path: str | Path, sequence_col: str = "single_chain_aa"):
    """
    Export sequences to FASTA format.

    Parameters
    ----------
    df : pd.DataFrame
        DataFrame with sequences
    output_path : str or Path
        Output file path
    sequence_col : str
        Column containing sequences to export
    """
    with open(output_path, "w") as f:
        for idx, row in df.iterrows():
            seq = row.get(sequence_col, "")
            if not seq:
                continue

            # Build header
            cdr3ab = row.get("CDR3ab", idx)
            cdr3a = row.get("CDR3_alpha", "")
            cdr3b = row.get("CDR3_beta", "")

            header = f">{cdr3ab} CDR3a={cdr3a} CDR3b={cdr3b}"
            f.write(f"{header}\n{seq}\n")

    logger.info(f"Exported {len(df)} sequences to {output_path}")