pythonforbio.
[WASM idle]
Biopython11/15

Seq

Starting sandbox…
Beginnerlesson

Seq

In one line: an immutable, string-like object of biological letters with biology-aware methods — no alphabet, just letters and your knowledge of what they mean.

Creating and using

from Bio.Seq import Seq, MutableSeq

s = Seq("ATGGCCATTGTAATGGGCCGC")

len(s); s[0]; s[3:9]; s[::-1]
str(s)                      # → a plain Python string
s + Seq("TTT")              # concatenation
"ATG" in s                  # membership
s.count("GC")               # NON-overlapping count, like str.count
s.count_overlap("GC")       # overlapping count
s.find("ATG"); s.index("ATG")
s.upper(); s.lower()
s.startswith("ATG"); s.endswith("TAA")
s.split("N"); s.strip("N")
s.replace("-", "")          # this is how you remove gaps (ungap() was removed)

Seq supports nearly the whole str interface, so most string intuition transfers directly.

Bio.Alphabet is gone

Removed in Biopython 1.78. This is a TypeError now:

from Bio.Alphabet import IUPAC          # ImportError
Seq("ATG", IUPAC.unambiguous_dna)        # TypeError

A Seq is just letters. The library does not know whether ATG is DNA or protein, and mostly does not need to. Where a format requires it (writing GenBank), you supply it as metadata:

record.annotations["molecule_type"] = "DNA"

Every tutorial written before 2020 will show the alphabet API. See Biopython Deprecations.

The biology methods

s.complement()                 # A↔T, G↔C — same 5'→3' orientation
s.reverse_complement()         # the one you actually want
s.transcribe()                 # DNA → RNA (T → U)
s.back_transcribe()            # RNA → DNA
s.translate()                  # → protein
s.complement_rna(); s.reverse_complement_rna()

complement() vs reverse_complement(): complement alone gives you the opposite strand read in the same direction, which is almost never biologically meaningful. For "the other strand's sequence as it would be read", you want reverse_complement().

Translation, carefully

s.translate()                              # standard code; * for stops
s.translate(to_stop=True)                  # stop at the first stop codon
s.translate(table=2)                       # vertebrate mitochondrial
s.translate(table="Bacterial", cds=True)   # validate as a proper CDS
s.translate(stop_symbol="")                # omit stops entirely
s.translate(gap="-")                        # handle gapped alignments

Points to internalise:

  • to_stop=True is usually what you mean. Without it, translation runs to the end of the sequence and you get MAIVMGR*KGAR* with internal stops.
  • cds=True validates: it requires a length divisible by 3, a valid start codon, a single terminal stop, and no internal stops. It raises TranslationError otherwise. Use it when you believe you have a real CDS — the exception is information.
  • Partial codons warn (BiopythonWarning: Partial codon) and are dropped.
  • Table matters. Mitochondrial genomes use different codon tables; TGA is a stop in the standard code and tryptophan in vertebrate mitochondria. Getting this wrong silently truncates your protein.
from Bio.Data import CodonTable
CodonTable.unambiguous_dna_by_id[1]        # the standard table
print(CodonTable.unambiguous_dna_by_name["Vertebrate Mitochondrial"])

Immutability

s[0] = "G"                     # TypeError — Seq is immutable

m = MutableSeq("ATGGCC")
m[0] = "G"
m.remove("C"); m.reverse()
back = Seq(m)                  # freeze it again

Immutability means Seq is hashable and safe to share. Use MutableSeq for in-place editing, then convert back.

For simulating a substitution the string route is usually clearer:

mutated = Seq(str(s)[:pos] + alt + str(s)[pos + len(ref):])

Undefined and partial sequences

from Bio.Seq import Seq
s = Seq(None, length=1_000_000)      # known length, unknown content
len(s)                                # 1000000
str(s)                                # UndefinedSequenceError

This is how Biopython represents a record whose length is known from a header but whose sequence was not loaded (common with GenBank contig records and some SeqIO indexing modes). Guard with a try/except if you parse arbitrary files.

Bioinformatics examples

from Bio.Seq import Seq
from Bio.SeqUtils import gc_fraction, molecular_weight

s = Seq("ATGGCCATTGTAATGGGCCGCTGAAAGGGTGCCCGATAG")

gc_fraction(s)                             # 0.0–1.0 (NOT the removed GC(), 0–100)
molecular_weight(s, seq_type="DNA")
s.translate(to_stop=True)                  # MAIVMGR

# find ORFs in all six frames
def orfs(seq, min_aa=50):
    for strand, nuc in [(+1, seq), (-1, seq.reverse_complement())]:
        for frame in range(3):
            trimmed = nuc[frame: len(nuc) - (len(nuc) - frame) % 3]
            for prot in trimmed.translate().split("*"):
                if len(prot) >= min_aa:
                    yield strand, frame, prot

# reverse-strand variant: get the sequence as the transcript reads it
if strand == "-":
    coding = genomic_seq.reverse_complement()

# codon usage
codons = [str(s[i:i+3]) for i in range(0, len(s) - 2, 3)]
pd.Series(codons).value_counts()

# GC at the third codon position (GC3) — a codon-bias measure
gc3 = sum(1 for i in range(2, len(s), 3) if s[i] in "GC") / (len(s) // 3)

Ambiguity codes

IUPAC ambiguity letters (N, R, Y, W, S, K, M, ...) are handled correctly by complement/reverse_complementR (A or G) complements to Y (C or T). Translation of an ambiguous codon gives X where the amino acid is undetermined, or the specific residue where all possibilities agree.

Real-world sequences contain N (unknown base) and soft-masked lowercase (repeat-masked regions). s.upper() before comparisons, and decide deliberately whether N counts in your GC denominator — gc_fraction has an ambiguous= parameter for this.

Common mistakes

  • Using Bio.Alphabet. Removed.
  • ungap(), tomutable(), toseq(). Removed. Use .replace("-", "") and the constructors.
  • translate() without to_stop=True, getting internal stops.
  • Wrong codon table for mitochondrial or bacterial sequences.
  • complement() where you meant reverse_complement().
  • Forgetting reverse-strand genes need reverse-complementing before translation.
  • Assuming Seq is mutable.
  • Ignoring soft-masking. Lowercase bases fail a naive == "G" test.
  • SeqUtils.GC() — removed in 1.82, and it returned a percentage while gc_fraction returns a fraction. Off by 100×.

See also

SeqRecord · SeqIO · Bio.SeqUtils · Biopython Deprecations · Applied - Sequence Composition Analysis · Applied - FASTA and FASTQ Workflows

lesson example

No output yet — run the code to populate this drawer.