Day 2: Nucleotide Composition, GC Content Analysis & Six-Frame Sequence Translation

1. Biological Significance of Nucleotide Composition & GC Content

The nucleotide composition of a genome or gene reveals critical insights into its evolutionary history, thermodynamic stability, codon bias, and structural organization. The GC content is defined as the molar fraction of guanine and cytosine nucleotides within a DNA sequence:

$$\text{GC}\% = \frac{G + C}{A + T + G + C} \times 100\%$$

Why GC Content Matters in Molecular Biology

• Thermal Stability: G:C base pairs form three hydrogen bonds, compared to two in A:T base pairs. High GC sequences require higher temperatures for thermal denaturation ($T_m$), influencing PCR optimization, hybridization assays, and thermophilic bacterial adaptations.
• Mammalian Isochore Structure: The mammalian genome is organized into giant genomic mosaics (isochores, > 300 kb) categorized by distinct GC levels (e.g. L1, L2, H1, H2, H3). Gene density, SINE/LINE retrotransposon distribution, and replication timing strongly correlate with isochore GC content.
• Bacterial Phylogenetics: Bacterial genomes exhibit extreme GC variation, ranging from ~13% in endosymbionts (Candidatus Carsonella ruddii) to > 75% in actinobacteria (Streptomyces coelicolor).

2. Replication Origin Prediction via GC Skew

During bacterial DNA replication, the leading and lagging strands undergo asymmetric mutational pressures. The leading strand remains single-stranded for longer durations during replication fork progression, predisposing cytosine bases to spontaneous hydrolytic deamination into uracil (which subsequently replicates as thymine, $C \to T$). Consequently, the leading strand becomes enriched in Guanine over Cytosine. This asymmetry is quantified using GC Skew:

$$\text{GC Skew} = \frac{G - C}{G + C}$$

3. CpG Islands and Epigenetic Architecture

In vertebrate genomes, cytosine in the dinucleotide context 5'-CpG-3' is frequently methylated to 5-methylcytosine ($5\text{mC}$). Spontaneous deamination of $5\text{mC}$ yields Thymine ($T$), which is poorly repaired compared to deamination of unmethylated Cytosine into Uracil ($U$). Over evolutionary timescales, this has caused severe depletion of CpG dinucleotides across the genome (the CpG deficiency phenomenon).

However, regions near gene promoters—termed CpG Islands—remain unmethylated and escape deamination. The standard Takai-Jones criteria define a CpG island as:

1. Sequence window length $\ge 200\text{ bp}$ (typically $\ge 500\text{ bp}$).
2. Calculated $\text{GC content} \ge 50\%$.
3. Observed-to-Expected CpG ratio ($\text{Obs/Exp}$) $\ge 0.60$, calculated as:

$$\text{Obs/Exp CpG} = \frac{N_{\text{CpG}} \times L}{N_{\text{C}} \times N_{\text{G}}}$$

where $N_{\text{CpG}}$ is the count of CpG dinucleotides, $N_{\text{C}}$ and $N_{\text{G}}$ are individual base counts, and $L$ is the total length of the sequence window.

4. Universal Six-Frame Sequence Translation

Because double-stranded DNA contains two antiparallel complementary strands, and a ribosome reads codons as non-overlapping nucleotide triplets, any double-stranded sequence possesses six distinct potential reading frames:

Production Python Script: GC Content, Sliding Skew & 6-Frame Translation

#!/usr/bin/env python3
"""
Day 2 Practical: GC Profiler, Skew Calculator, and 6-Frame Translator
"""

from Bio.Seq import Seq
from Bio.SeqUtils import gc_fraction

def calculate_nucleotide_stats(sequence: str):
    """Compute base counts, GC%, AT%, and molecular weight."""
    seq = sequence.upper()
    total = len(seq)
    if total == 0:
        return {}
    counts = {base: seq.count(base) for base in ["A", "T", "G", "C"]}
    gc_pct = ((counts["G"] + counts["C"]) / total) * 100
    at_pct = ((counts["A"] + counts["T"]) / total) * 100
    return {
        "length": total,
        "counts": counts,
        "gc_percent": round(gc_pct, 2),
        "at_percent": round(at_pct, 2)
    }

def sliding_window_gc_skew(sequence: str, window_size=1000, step=200):
    """Compute GC skew across sliding windows."""
    seq = sequence.upper()
    positions = []
    skews = []
    for i in range(0, len(seq) - window_size + 1, step):
        subseq = seq[i : i + window_size]
        g = subseq.count("G")
        c = subseq.count("C")
        skew = (g - c) / (g + c) if (g + c) > 0 else 0.0
        positions.append(i + window_size // 2)
        skews.append(round(skew, 4))
    return positions, skews

def six_frame_translation(dna_sequence: str, table=1):
    """Translate DNA sequence across all 6 reading frames."""
    dna = Seq(dna_sequence.upper())
    rev_dna = dna.reverse_complement()
    
    translations = {}
    # Forward frames +1, +2, +3
    for frame in range(3):
        # Trim length to multiple of 3
        subseq = dna[frame:]
        trimmed = subseq[:len(subseq) - (len(subseq) % 3)]
        translations[f"+{frame+1}"] = str(trimmed.translate(table=table))
        
    # Reverse frames -1, -2, -3
    for frame in range(3):
        subseq = rev_dna[frame:]
        trimmed = subseq[:len(subseq) - (len(subseq) % 3)]
        translations[f"-{frame+1}"] = str(trimmed.translate(table=table))
        
    return translations

if __name__ == "__main__":
    sample_dna = (
        "ATGGCCCTGTGGATGCGCCTCCTGCCCCTGCTGGCGCTGCTGGCCCTCTGGGGACCTGACCCAGCCGCAGCC"
        "TTTGTGAACCAACACCTGTGCGGCTCACACCTGGTGGAAGCTCTCTACCTAGTGTGCGGGGAACGAGGCTTC"
        "TTCTACACACCCAAGACCCGCCGGGAGGCAGAGGACCTGCAGGTGGGGCAGGTGGAGCTGGGCGGCGGCCCT"
        "GGTGCAGGCAGCCTGCAGCCCTTGGCCCTGGAGGGGTCCCTGCAGAAGCGTGGCATTGTGGAACAATGCTGT"
        "ACCAGCATCTGCTCCCTCTACCAGCTGGAGAACTACTGCAACTAG"
    )

    stats = calculate_nucleotide_stats(sample_dna)
    print("=== NUCLEOTIDE PROFILE ===")
    print(f"Length: {stats['length']} bp | GC%: {stats['gc_percent']}% | AT%: {stats['at_percent']}%")
    print(f"Base Counts: A={stats['counts']['A']}, T={stats['counts']['T']}, G={stats['counts']['G']}, C={stats['counts']['C']}")

    print("
=== SIX-FRAME TRANSLATION ===")
    frames = six_frame_translation(sample_dna)
    for frame_name, protein in frames.items():
        print(f"[{frame_name}]: {protein[:45]}... (Length: {len(protein)} aa)")