Análisis de secuencias de ADN con Python
Al finalizar este tema
Podrás calcular el contenido de GC de una secuencia de ADN, analizar archivos FASTA y convertir codones en aminoácidos utilizando Python.
Variables: contenedores de datos
Al igual que se anota en un cuaderno de laboratorio "Concentración de la muestra: 2.5 mg/mL", en Python se almacenan valores en variables.
gene_name = "BRCA1"sequence = "ATGCGATCGATCGATCG"gc_ratio = 0.529
print(f"Gen: {gene_name}")print(f"Longitud de la secuencia: {len(sequence)} bp")print(f"Proporción de GC: {gc_ratio:.1%}")
assert gene_name == "BRCA1"assert len(sequence) == 17Cadenas: Manipulación de secuencias de ADN
Las secuencias de ADN son cadenas (strings) compuestas por cuatro letras: A, T, G y C. Se pueden analizar las secuencias utilizando los métodos de cadenas de Python.
sequence = "ATGCGATCGATCGATCG"
# Contar una base concretag_count = sequence.count("G")c_count = sequence.count("C")print(f"G: {g_count}, C: {c_count}")
# Crear la secuencia complementaria (A↔T, G↔C)complement_table = str.maketrans("ATGC", "TACG")complement = sequence.translate(complement_table)print(f"Original: {sequence}")print(f"Complementaria: {complement}")print(f"Complementaria inversa: {complement[::-1]}")
assert complement == "TACGCTAGCTAGCTAGC"assert complement[::-1] == "CGATCGATCGATCGCAT"Función: Calculadora de contenido GC
Una función es un paso del protocolo experimental empaquetado de forma reutilizable. Una vez creada, se puede aplicar a cualquier secuencia.
def calculate_gc_content(sequence: str) -> float: sequence = sequence.upper() gc_count = sequence.count("G") + sequence.count("C") return (gc_count / len(sequence)) * 100
# Pruebaseq1 = "ATGCGATCGATCGATCG"seq2 = "AAAAAAAAAA"seq3 = "GGGGGGGGGG"
print(f"seq1 GC: {calculate_gc_content(seq1):.1f}%")print(f"seq2 GC: {calculate_gc_content(seq2):.1f}%")print(f"seq3 GC: {calculate_gc_content(seq3):.1f}%")
assert abs(calculate_gc_content(seq1) - 52.9) < 0.1assert calculate_gc_content(seq2) == 0.0assert calculate_gc_content(seq3) == 100.0Listas y bucles: procesamiento de múltiples secuencias a la vez
Al igual que se aplica el mismo tratamiento a todos los pocillos de una placa de 96 pocillos, un bucle for repite la misma operación en múltiples datos.
def calculate_gc_content(sequence: str) -> float: sequence = sequence.upper() gc_count = sequence.count("G") + sequence.count("C") return (gc_count / len(sequence)) * 100
genes = { "BRCA1": "ATGGATTTATCTGCTCTTCGCGTTGAAGAAGTACAAAATGTC", "TP53": "ATGGAGGAGCCGCAGTCAGATCCTAGCGTGAGTTTGCTGTGA", "EGFR": "ATGCGACCCTCCGGGACGGCCGGGGCAGCGCTCCTGGCGCTG",}
results = []for name, seq in genes.items(): gc = calculate_gc_content(seq) results.append(gc) print(f"{name}: GC={gc:.1f} %, longitud={len(seq)} bp")
assert len(results) == 3assert all(0 <= gc <= 100 for gc in results)Diccionario: Procesamiento de archivos FASTA
Un diccionario es como asignar etiquetas a las muestras experimentales. Permite encontrar directamente la secuencia (valor) mediante el nombre del gen (clave).
def parse_fasta(fasta_text: str) -> dict[str, str]: sequences: dict[str, str] = {} current_header = "" for line in fasta_text.strip().split("\n"): if line.startswith(">"): current_header = line[1:].strip() sequences[current_header] = "" else: sequences[current_header] += line.strip() return sequences
sample_fasta = """>BRCA1_humanATGGATTTATCTGCTCTTCGCGTTGAAGAAGTACAAAATGTC>TP53_humanATGGAGGAGCCGCAGTCAGATCCTAGCGTGAGTTTGCTGTGA"""
result = parse_fasta(sample_fasta)print(f"Número de secuencias analizadas: {len(result)}")for name, seq in result.items(): print(f" {name}: {len(seq)}bp")
assert len(result) == 2assert result["BRCA1_human"] == "ATGGATTTATCTGCTCTTCGCGTTGAAGAAGTACAAAATGTC"assert result["TP53_human"] == "ATGGAGGAGCCGCAGTCAGATCCTAGCGTGAGTTTGCTGTGA"Conversión de codón a aminoácido
Una secuencia de ADN de 3 nucleótidos (codón) especifica un aminoácido. Esta tabla de conversión puede representarse como un diccionario.
CODON_TABLE = { "ATG": "M", # Metionina (codón de inicio) "TTT": "F", "TTC": "F", # Phenylalanine "TTA": "L", "TTG": "L", "CTT": "L", "CTC": "L", # Leucine "GAT": "D", "GAC": "D", # Aspartic acid "GAA": "E", "GAG": "E", # Glutamic acid "GCT": "A", "GCC": "A", # Alanine "TAA": "*", "TAG": "*", "TGA": "*", # Codones de terminación}
def translate_sequence(dna: str) -> str: protein = [] for i in range(0, len(dna) - 2, 3): codon = dna[i:i+3] amino_acid = CODON_TABLE.get(codon, "?") if amino_acid == "*": break protein.append(amino_acid) return "".join(protein)
test_seq = "ATGGATTTTGAA"protein = translate_sequence(test_seq)print(f"DNA: {test_seq}")print(f"Protein: {protein}")
assert protein == "MDFE"Inténtalo tú mismo (Ejemplo atenuado)
Completa la función para calcular el contenido de GC rellenando los siguientes espacios en blanco.
def gc_content(seq):seq = seq.upper()gc = seq.count('G') + seq.count('')return / len(seq) * 100
Errores comunes y soluciones
Q: IndentationError: expected an indented block
Python utiliza la sangría (indentation) para distinguir los bloques de código. def, for, if la siguiente línea requiere obligatoriamente una sangría de 4 espacios.
Q: KeyError: 'BRCA1'
Ocurre cuando la clave correspondiente no existe en el diccionario. Al usar dict.get("BRCA1", "no encontrado"), se devuelve un valor predeterminado sin errores, incluso si la clave no existe.
Q: La secuencia contiene letras minúsculas, por lo que el conteo no coincide
Primero, unifique a mayúsculas con .upper() antes de procesar. En un archivo FASTA real, las letras minúsculas pueden representar regiones enmascaradas por repeticiones (repeat-masked region).
En la siguiente publicación, aprenderemos cómo gestionar el control de versiones de este código de análisis con Git.