Bio Restriction Fragment Analysis

Version Compatibility

Reference examples tested with: BioPython 1.83+ (API verified on 1.86)

Before using code patterns, verify installed versions match. If versions differ:

• Python: pip show biopython, then check the return shape with from Bio.Restriction import EcoRI; from Bio.Seq import Seq; EcoRI.catalyze(Seq('GAATTCGAATTC')) (a tuple of fragment Seqs)

If code throws ImportError, AttributeError, or TypeError, introspect the installed package and adapt the example to match the actual API rather than retrying.

Restriction Fragment Analysis

"What fragments will this digest produce?" -> Simulate the cut and return fragment lengths (and sequences), for one enzyme or a double digest, on linear or circular DNA.

• Python: enzyme.catalyze(seq, linear=...) returns a tuple of fragment Seq objects directly.

Two facts decide whether the prediction is right. First, catalyze() returns the tuple of fragments itself -- EcoRI.catalyze(seq) is (Seq(...), Seq(...), Seq(...)). Do NOT write catalyze(seq)[0] to "get the fragments": that returns only the first fragment, and iterating it then counts its bases, producing nonsense sizes like [1, 1, 1, ...]. Second, topology sets the fragment count: a linear molecule with n cuts yields n+1 fragments; a circular molecule with n cuts yields n. Passing the default linear=True to a plasmid invents one extra fragment that does not exist on the bench.

Fragment Count By Topology

A plasmid showing three bands on a single-enzyme digest has three sites, not two. Counting plasmid bands as if the molecule were linear is the most common digest-interpretation error.

Predict Fragment Sizes

Goal: Turn a digest into the band sizes a gel would show.

Approach: Call catalyze() with the correct topology, take the lengths of the returned fragments, sort descending. The summed sizes must equal the molecule length -- a built-in correctness check.

from Bio import SeqIO
from Bio.Restriction import EcoRI

record = SeqIO.read('sequence.fasta', 'fasta')
seq = record.seq

fragments = EcoRI.catalyze(seq, linear=True)   # tuple of Seq; NO [0]
sizes = sorted((len(f) for f in fragments), reverse=True)
assert sum(sizes) == len(seq)                  # fragments must account for the whole molecule
print(f'{len(sizes)} fragments: {sizes}')

Linear vs Circular Digestion

from Bio.Restriction import EcoRI

linear_frags   = EcoRI.catalyze(seq, linear=True)
circular_frags = EcoRI.catalyze(seq, linear=False)   # plasmid: one fewer fragment
print(f'Linear: {len(linear_frags)} fragments; Circular: {len(circular_frags)} fragments')

Double Digest

Goal: Predict fragments when two enzymes cut the same molecule.

Approach: A double digest cuts at the union of both site sets. The robust way is to build a RestrictionBatch and let Analysis collect every position, then compute fragments from the sorted positions (this generalizes to any number of enzymes and to circular DNA). Sequential catalyze calls also work but are easy to get wrong on circular DNA.

from Bio.Restriction import EcoRI, BamHI, RestrictionBatch, Analysis

def fragments_from_positions(seq_len, positions, linear=True):
    '''Fragment sizes (bp) from a set of 1-based cut positions.'''
    cuts = sorted(set(positions))
    if not cuts:
        return [seq_len]
    spans = [cuts[i + 1] - cuts[i] for i in range(len(cuts) - 1)]
    if linear:
        return [cuts[0]] + spans + [seq_len - cuts[-1]]
    return spans + [(seq_len - cuts[-1]) + cuts[0]]    # circular: wrap-around fragment

batch = RestrictionBatch([EcoRI, BamHI])
positions = [p for sites in Analysis(batch, seq).with_sites().values() for p in sites]
sizes = sorted(fragments_from_positions(len(seq), positions, linear=True), reverse=True)
print(f'Double digest: {len(sizes)} fragments: {sizes}')
assert sum(sizes) == len(seq)

Interpreting The Gel: Size, Resolution, And Topology

Fragment sizes are not what a gel directly reports; migration distance is. The decisions below come from gel physics, not from BioPython, and they determine whether predicted bands are actually resolvable.

Simulate A Gel Pattern (Text)

Goal: Lay predicted digests next to a ladder to see which bands resolve and which co-migrate.

Approach: Pool all band sizes and the ladder, sort descending, and mark each lane. Co-migration shows up as multiple marks at one size.

def simulate_gel(digests, ladder=None):
    '''digests: {lane_name: [sizes]} -> print a text gel against a ladder.'''
    ladder = ladder or [10000, 8000, 6000, 5000, 4000, 3000, 2000, 1500, 1000, 750, 500, 250]
    bands = sorted(set(ladder).union(*[set(s) for s in digests.values()]), reverse=True)
    header = f'{"size":>6} | {"ladder":^6} | ' + ' | '.join(f'{n:^8}' for n in digests)
    print(header); print('-' * len(header))
    for b in bands:
        row = f'{b:>6} | {"---" if b in ladder else "":^6} | '
        row += ' | '.join(f'{"=" * 4 * lane.count(b):^8}' for lane in digests.values())
        print(row)

simulate_gel({'EcoRI': sizes})

Diagnostic Digest Report

Goal: Document an expected digest to plan a clone-verification gel.

Approach: Report site count, positions, fragment sizes, and the total; flag any pair of fragments too close to resolve.

def fragment_report(seq, enzyme, linear=True, resolution=0.10):
    sites = enzyme.search(seq, linear=linear)
    sizes = sorted((len(f) for f in enzyme.catalyze(seq, linear=linear)), reverse=True)
    print(f'{enzyme} ({enzyme.site}): {len(sites)} site(s) at {sites}')
    print(f'  {len(sizes)} fragments, total {sum(sizes)} bp: {sizes}')
    close = [(a, b) for a, b in zip(sizes, sizes[1:]) if a and (a - b) / a < resolution]
    if close:
        print(f'  WARNING likely co-migrating (<{resolution:.0%} apart): {close}')
    return sizes

fragment_report(seq, EcoRI)

Compare Expected vs Observed Bands

def compare_fragments(expected, observed, tolerance=50):
    '''Match predicted sizes to gel-measured sizes within a tolerance (bp).'''
    obs = list(observed)
    matched, missing = [], []
    for exp in expected:
        hit = next((o for o in obs if abs(exp - o) <= tolerance), None)
        if hit is None:
            missing.append(exp)
        else:
            matched.append((exp, hit)); obs.remove(hit)
    print('matched:', matched)
    if missing: print('missing (predicted, not seen -- co-migration or partial digest?):', missing)
    if obs:     print('extra (seen, not predicted -- star activity, contaminant, or wrong map?):', obs)

compare_fragments(expected=[3000, 2000, 1500, 500], observed=[3050, 2000, 1480, 510, 200])

Extra bands and a smeary ladder of intermediate sizes often signal a partial digest (not every site cut), which produces sums of adjacent complete-digest fragments. Drive the reaction to completion (more enzyme or time) before concluding the map is wrong.

Common Errors

Related Skills

• restriction-sites - Find the cut positions that drive fragment calculation
• restriction-mapping - Order sites and compute inter-site distances
• enzyme-selection - Choose enzymes that give a resolvable, diagnostic band pattern
• sequence-manipulation/seq-objects - Work with the fragment Seq objects

References

• Smith HO, Birnstiel ML. A simple method for DNA restriction site mapping. Nucleic Acids Res. 1976;3(9):2387-2398. doi:10.1093/nar/3.9.2387
• Schwartz DC, Cantor CR. Separation of yeast chromosome-sized DNAs by pulsed field gradient gel electrophoresis. Cell. 1984;37(1):67-75. doi:10.1016/0092-8674(84)90301-5