The interface is a sequence file
What does it mean for an application to work with DNA? A useful starting point is a file of bases, an accompanying record of measurement quality and an API that retrieves a reference sequence.
A sequencer measures a sample; base-calling software turns those measurements into sequence data. The programs below begin at that digital boundary. They read, transform and report on data. Calling a method named transcribe() does not manufacture RNA, and an HTTP request to a sequence database does not control a laboratory instrument.
We'll build five small programs: a FASTA reader, a sequence transformation example, a FASTQ filter, a reference downloader and a JSON report. The local inputs are deliberately tiny, made-up fixtures, so every result is inspectable.
Set up the examples
Download the runnable examples, unzip them and open a terminal in that directory. The examples were checked with Python 3.9 and Biopython 1.85. The dependency is pinned so you can reproduce those results.
python3 -m venv .venv
source .venv/bin/activate
python -m pip install -r requirements.txt
python 01_fasta.py
On Windows, activate the environment with .venv\Scripts\activate instead. Four scripts run entirely offline. Only 04_ncbi.py makes a network request. The scripts write their small output files into the current directory and replace those files when rerun.
Read and validate FASTA
FASTA is a useful first boundary: a header beginning with > followed by sequence text. Biopython's sequence I/O parses it into a SeqRecord with an identifier, description and sequence.
Save this as 01_fasta.py:
from io import StringIO
from Bio import SeqIO
from Bio.SeqUtils import gc_fraction
# A made-up sequence, not an experimental construct.
fasta = ">toy_example\nATGAAAGGTTAA\n"
record = SeqIO.read(StringIO(fasta), "fasta")
sequence = str(record.seq).upper()
# This small importer deliberately accepts only A/C/G/T and unknown bases (N).
invalid = set(sequence) - set("ACGTN")
if not sequence or invalid:
raise ValueError(f"Invalid or empty DNA sequence: {sorted(invalid)}")
known_bases = sum(base in "ACGT" for base in sequence)
gc = gc_fraction(sequence, ambiguous="remove") if known_bases else None
print("id:", record.id)
print("length:", len(sequence))
print("unknown bases:", sequence.count("N"))
print("GC among known bases:", f"{gc:.1%}" if gc is not None else "undefined")
SeqIO.write(record, "sample.fasta", "fasta")
id: toy_example
length: 12
unknown bases: 0
GC among known bases: 25.0%
The script also writes sample.fasta. For an existing file, replace the StringIO(fasta) input with a file path. Use SeqIO.parse() when the file contains multiple records; read() expects exactly one.
The validator deliberately accepts only A, C, G, T and N. Other IUPAC ambiguity codes exist; a general-purpose importer needs an explicit policy for them. Here, N means an unknown base. The GC calculation excludes unknown bases from its denominator, and an all-N input reports an undefined fraction instead of misleadingly reporting zero.
Strands, coordinates and translation
A DNA sequence is conventionally written in its 5′→3′ direction. Its opposite strand is represented with a reverse complement, not just a reversed string. Python slicing uses zero-based, half-open coordinates: [3:6] selects bases 4 through 6.
02_sequence.py uses a deliberately chosen toy coding sequence and the standard genetic code:
from Bio.Seq import Seq
# Explicitly a toy coding sequence in the standard genetic code.
dna = Seq("ATGAAAGGTTAA")
print("bases 4-6:", dna[3:6])
print("reverse complement:", dna.reverse_complement())
print("RNA representation:", dna.transcribe())
print("protein representation:", dna.translate(table=1, cds=True))
assert str(dna.reverse_complement().reverse_complement()) == str(dna)
assert str(dna.translate(table=1, cds=True)) == "MKG"
bases 4-6: AAA
reverse complement: TTAACCTTTCAT
RNA representation: AUGAAAGGUUAA
protein representation: MKG
These sequence operations return data representations. cds=True checks the coding-sequence assumptions, including a valid start, length divisible by three and a terminal stop. It removes the terminal stop from the returned protein.
That does not mean every genomic interval translates into a meaningful protein. Strand, reading frame, genetic code and annotated features matter. This example already knows its coding sequence; it does not discover a gene or predict biological function.
Read sequencer output with FASTQ
FASTA stores bases. FASTQ also stores per-base quality scores. A Phred score Q corresponds to an estimated error probability of 10 ** (-Q / 10). Adding these probabilities estimates the number of errors in a read; it does not identify which bases are wrong.
03_fastq.py creates two synthetic reads with identical bases and different qualities, then retains the cleaner one:
from io import StringIO
from Bio import SeqIO
# Synthetic reads with Sanger/Phred+33 quality characters.
fastq = """@good_read
ATGAAAGGTTAA
+
IIIIIIIIIIII
@noisy_read
ATGAAAGGTTAA
+
!!!!IIIIIIII
"""
kept = []
for read in SeqIO.parse(StringIO(fastq), "fastq"):
quality = read.letter_annotations["phred_quality"]
expected_errors = sum(10 ** (-q / 10) for q in quality)
# An illustrative filter; real thresholds depend on the assay and platform.
accepted = expected_errors < 1 and "N" not in str(read.seq).upper()
print(read.id, f"expected errors={expected_errors:.4f}", "keep" if accepted else "drop")
if accepted:
kept.append(read)
count = SeqIO.write(kept, "filtered.fastq", "fastq")
assert count == 1
print("saved reads:", count)
good_read expected errors=0.0012 keep
noisy_read expected errors=4.0008 drop
saved reads: 1
Open filtered.fastq to inspect the retained read and its quality string. This example uses Sanger/Phred+33 encoding, where I represents Q40 and ! represents Q0. The threshold is illustrative: real quality criteria depend on the sequencing platform and analysis.
For large files, stream accepted records into the writer instead of collecting them in a list. Biopython's FASTQ cookbook covers filtering and preserving quality information.
Fetch a public reference from NCBI
For an application that needs a known sequence, use a versioned accession and a documented database API. Here is a single request for bases 1–120 of the Arabidopsis thaliana chloroplast reference, NC_000932.1.
NCBI's EFetch API uses one-based, inclusive coordinates for seq_start and seq_stop. That differs from Python slicing. strand=1 requests the forward strand.
Save this as 04_ncbi.py:
from datetime import datetime, timezone
from io import StringIO
import json
import os
from pathlib import Path
from urllib.parse import urlencode
from urllib.request import Request, urlopen
from Bio import SeqIO
accession = "NC_000932.1" # Arabidopsis thaliana chloroplast reference.
start, stop = 1, 120 # NCBI uses 1-based, inclusive coordinates here.
params = {
"db": "nuccore", "id": accession,
"rettype": "fasta", "retmode": "text",
"seq_start": start, "seq_stop": stop,
"strand": 1, "tool": "dna_sequence_notes",
}
if os.environ.get("NCBI_EMAIL"):
params["email"] = os.environ["NCBI_EMAIL"]
url = "https://eutils.ncbi.nlm.nih.gov/entrez/eutils/efetch.fcgi?" + urlencode(params)
request = Request(url, headers={"User-Agent": "dna-sequence-notes/1.0"})
with urlopen(request, timeout=20) as response:
text = response.read().decode("utf-8")
if not text.lstrip().startswith(">"):
raise ValueError("NCBI did not return FASTA; do not treat an error page as DNA")
record = SeqIO.read(StringIO(text), "fasta")
if len(record.seq) != stop - start + 1:
raise ValueError("Unexpected sequence length")
SeqIO.write(record, "reference.fasta", "fasta")
metadata = {
"accession_version": accession,
"source": "NCBI EFetch", "start_1based": start, "end_inclusive": stop,
"strand": "+", "retrieved_at": datetime.now(timezone.utc).isoformat(),
}
Path("reference.metadata.json").write_text(json.dumps(metadata, indent=2), encoding="utf-8")
print(record.description)
print("downloaded bases:", len(record.seq))
NC_000932.1:1-120 Arabidopsis thaliana chloroplast, complete genome
downloaded bases: 120
The live request returned that record when this article was updated. The script saves both reference.fasta and reference.metadata.json; the latter records the accession version, interval, strand and retrieval time. Keeping those details is more useful than saving an unlabelled string of bases.
For a recurring application, supply a real contact email through NCBI_EMAIL, follow NCBI's usage guidelines, cache successful responses and handle network failures explicitly. This small example stops on a failed request or malformed response rather than silently substituting data.
Expose a JSON report to an app
Once a sequence is parsed, an ordinary application can display a summary. 05_report.py produces a JSON object from the same known fixture:
from io import StringIO
import json
from Bio import SeqIO
from Bio.SeqUtils import gc_fraction
record = SeqIO.read(StringIO(">toy_example\nATGAAAGGTTAA\n"), "fasta")
sequence = str(record.seq).upper()
report = {
"id": record.id,
"source": "synthetic tutorial fixture",
"length_bp": len(sequence),
"unknown_bases": sequence.count("N"),
"gc_fraction_known_bases": gc_fraction(sequence, ambiguous="remove"),
"interval": {"start": 0, "end": len(sequence), "coordinates": "0-based half-open"},
}
print(json.dumps(report, indent=2))
assert report["length_bp"] == 12
assert report["gc_fraction_known_bases"] == 0.25
The report contains length_bp: 12, unknown_bases: 0 and gc_fraction_known_bases: 0.25, plus the source and coordinate convention. A frontend can format that fraction as 25% and use the interval to label a sequence selection.
This script intentionally uses the validated toy input. If you turn it into an upload endpoint, reuse the first example's alphabet validation and undefined-GC handling, limit upload size, and process records incrementally. Return validation errors as errors; don't quietly remove bases to make an input pass.
A useful API response distinguishes the original measurement, the reference and the calculation performed. None of those fields alone establishes what a sequence does in an organism.
Keep the pipeline reproducible
Run the offline programs first, then the reference download:
python 01_fasta.py
python 02_sequence.py
python 03_fastq.py
python 05_report.py
python 04_ncbi.py
The sequence and filtering examples include assertions against known answers. Keep those fixtures when expanding the code. A reverse complement applied twice should restore the original sequence; the sample filter should retain exactly one record; the JSON report should agree with the parsed input.
Store downloaded reference data with its accession version and interval. Preserve the input file and dependency versions alongside a report. When a result changes, this makes it possible to distinguish a code change from a different reference or input.
For a web service, test your parser and reports with local fixtures. Keep live API checks separate so a database outage doesn't masquerade as a regression in sequence processing.
Build something small and inspectable
These pieces are enough for a sequence-file viewer with base counts, a read-quality summary that explains why a read was filtered, or a reference browser that retrieves a named interval and preserves its provenance.
Start with one input format, a documented coordinate convention and a result you can verify by hand. Then add real datasets without losing the ability to inspect how each output was produced.
The connection to physical DNA is through measurements and laboratory workflows. The software work shown here is a practical part of that connection: parsing those measurements, retrieving references and turning them into useful, testable application data.