De HMM de 2 estados a HMM de 100 estados

Los HMM tratados en M18～M20 eran modelos de juguete con dos estados (H, L). El punto donde los HMM se vuelven realmente potentes en la práctica es cuando el número de estados explota y comienzan a expresar con sofisticación la estructura de las secuencias. Los Profile HMM son la obra maestra de esta expansión, y la búsqueda de genes es la corona de sus aplicaciones prácticas.

En este capítulo repasamos los tres estados del Profile HMM (Match, Insert, Delete) y su estructura topológica, el uso real de las herramientas HMMER, y cómo se utilizan los HMM para la detección de estructuras génicas. Es el último capítulo de la serie sobre HMM y también el más práctico.

Profile HMM — De MSA a HMM

Supongamos que hemos completado un Alineamiento Múltiple de Secuencias (MSA). Al alinear 100 proteínas homólogas, obtenemos una matriz de 200 columnas. Necesitamos un método para resumir esta información y buscar si una nueva secuencia es miembro de esta familia.

El enfoque más simple consiste en almacenar la frecuencia de caracteres en cada columna: una PSSM (Position-Specific Scoring Matrix). Sin embargo, la PSSM no maneja naturalmente los gaps (inserciones y deleciones). El Profile HMM supera esta limitación.

Se definen tres estados para cada columna del MSA:

• Estado Match (M_i): Emite un carácter en la columna i. La probabilidad de emisión se aprende de la distribución de aminoácidos de esa columna.
• Estado Insert (I_i): Inserción entre las columnas i e i+1. La probabilidad de emisión es cercana a la distribución de fondo.
• Estado Delete (D_i): Salta la columna i. No hay emisión (estado silencioso).

Estos tres estados están presentes en cada columna y se conectan únicamente mediante una topología fija.

M_1 → M_2 → M_3 → ... → M_L
      ↕     ↕     ↕            ↕
    D_1 → D_2 → D_3 → ... → D_L
      ↕     ↕     ↕            ↕
    I_0   I_1   I_2   ...    I_L
     ↑     ↑     ↑             ↑
    (loop) (loop) ...       (loop)

Match es Match·Delete·Insert. Delete es Match·Delete. Insert es un bucle consigo mismo o Match·Delete. El número de estados es aproximadamente 3L + 2 (donde L es el número de columnas del AAM); para L=200, hay 602 estados. Es un HMM grande, pero la estructura básica del llenado de la cuadrícula es idéntica a la de M18.

Búsqueda con Profile HMM mediante Viterbi

Cuando llega una nueva secuencia de proteínas X, se ejecuta el algoritmo de Viterbi sobre el Profile HMM.

$V_j(k) = \max_{s} \left[ V_{j-1}(s) \cdot a(s \to k) \right] \cdot b(k, x_j)$

k es ahora un estado (uno entre M_i, I_i, D_i), y j es la posición en la secuencia. El tamaño de la cuadrícula es O(L × T). La ruta óptima resultante es una matriz que asigna cada carácter de X a un estado específico (Match/Insert/Delete) de una columna determinada.

La puntuación total de esta ruta es el bit score (puntuación en bits), que se utiliza para determinar si hay un acierto de dominio. Cada dominio de Pfam se representa mediante este único Profile HMM, y la nueva secuencia de proteínas somete a búsqueda repetida todo Pfam (~20 000 dominios). El principio es el mismo que en M18, pero la escala es completamente diferente.

HMMER: la herramienta estándar para Profile HMM

HMMER (Sean Eddy et al.) es la implementación de facto estándar para Profile HMM. Se utiliza como motor subyacente en Pfam e InterPro; con unos pocos comandos se realiza la búsqueda práctica de dominios.

# 1) Instalación (funciona tal cual en Colab o SageMaker Studio Lab)
wget http://eddylab.org/software/hmmer/hmmer-3.4.tar.gz
tar xzf hmmer-3.4.tar.gz
cd hmmer-3.4 && ./configure && make && sudo make install
hmmscan -h

# 2) Descargar la base de datos Pfam
wget https://ftp.ebi.ac.uk/pub/databases/Pfam/current_release/Pfam-A.hmm.gz
gunzip Pfam-A.hmm.gz
hmmpress Pfam-A.hmm

# 3) Preparar la secuencia de consulta (por ejemplo, p53)
cat > p53.fasta <<'EOF'
>P04637 Human p53
MEEPQSDPSVEPPLSQETFSDLWKLLPENNVLSPLPSQAMDDLMLSPDDIEQWFTEDPGP
DEAPRMPEAAPPVAPAPAAPTPAAPAPAPSWPLSSSVPSQKTYQGSYGFRLGFLHSGTAK
SVTCTYSPALNKMFCQLAKTCPVQLWVDSTPPPGTRVRAMAIYKQSQHMTEVVRRCPHHE
EOF

# 4) Buscar dominios
hmmscan --domtblout p53.domtbl Pfam-A.hmm p53.fasta
head -20 p53.domtbl

--domtblout La salida muestra las columnas con el nombre del dominio, las coordenadas, i-Evalue (E-value individual del dominio), c-Evalue (E-value condicional) y la puntuación. En p53, el dominio P53 (PF00870) se identifica con las coordenadas correctas. PSI-BLAST requería varias iteraciones para obtener la homología distante; ahora se obtiene con una sola búsqueda.

También es importante conocer dos herramientas variantes importantes.

• hmmsearch: Busca múltiples secuencias con un único HMM (perfil).
• hmmscan: Busca una sola secuencia con múltiples HMM (todo Pfam).

Son enfoques opuestos. En la práctica, hmmscan se utiliza para el perfilado de dominios y hmmsearch para la detección de familias específicas.

Crear un nuevo dominio — hmmbuild

Para crear un perfil de dominio personalizado que no esté en Pfam, se utiliza hmmbuild con un MSA como material de partida.

# Archivo MSA (FASTA alineado o formato Stockholm) → archivo HMM
hmmbuild my_domain.hmm my_msa.sto

# Estimación probabilística del entrenamiento previo (Baum-Welch se ejecuta aquí dentro)
head -30 my_domain.hmm

Al abrir el archivo HMM, se observa que las tablas de probabilidad de emisión y transición de cada columna están almacenadas en el dominio logarítmico. Este único archivo constituye la huella estadística de un dominio. Al compartirlo con compañeros del equipo, se reproduce exactamente igual.

Búsqueda de genes mediante HMM — gene finding

La predicción de genes en procariotas es relativamente sencilla: inicio del gen (ATG), repetición de tripletes de codones y terminación (TAA/TAG/TGA). Herramientas como GeneMark o Prodigal expresan esta estructura mediante un HMM con entre 6 y 10 estados.

• Estado intergénico
• Estado de inicio de ORF (codón de inicio)
• Estados de posición del codón (posición 1, 2 y 3 del codón)
• Estado de terminación (codón de parada)

Las transiciones entre estos estados determinan la probabilidad de inicio y terminación del gen. Las probabilidades de emisión reflejan el contenido de GC y la sesgo en el uso de codones, entre otros factores.

En eucariotas, la situación es mucho más compleja: exones, intrones, UTR 5', UTR 3', sitios de empalme (splicing) y promotores. Herramientas como AUGUSTUS o GENSCAN utilizan HMM con cientos de estados. Como adelanto para el próximo capítulo, se presenta un esquema de las fases de los estados:

5'UTR → Start → Exon → Donor → Intron → Acceptor → Exon → ... → Stop → 3'UTR
                          ↓                          ↑
                          └─────── (repetición) ─────┘

La probabilidad de emisión se entrena por especie. En AUGUSTUS, la precisión alcanza el 90 % en especies como humano y ratón con conjuntos de entrenamiento grandes, mientras que es del 60-70 % en especies no modelo con menos datos de entrenamiento. La precisión de los predictores de genes basados en HMM depende de la calidad y cantidad de los datos de entrenamiento.

Práctica — Búsqueda de genes bacterianos con Prodigal

Realicemos una breve práctica centrada exclusivamente en procariotas.

# Instalar Prodigal (ligero y mucho más sencillo que HMMER)
conda install -c bioconda prodigal
prodigal -v

# Descargar el genoma de E. coli K-12
wget https://ftp.ncbi.nlm.nih.gov/genomes/all/GCA/000/005/845/GCA_000005845.2_ASM584v2/GCA_000005845.2_ASM584v2_genomic.fna.gz
gunzip GCA_000005845.2_ASM584v2_genomic.fna.gz

# Ejecutar la predicción génica (modo single = un organismo)
prodigal -i GCA_000005845.2_ASM584v2_genomic.fna \
         -o ecoli.gff -a ecoli.faa -f gff -p single

# Número de genes predichos
grep -c "^>" ecoli.faa

Se predice que E. coli K-12 tiene aproximadamente 4.300 genes (coincide con el valor real dentro del 5%). Es sorprendente que este resultado provenga de tan pocos estados en un HMM. Aquí se reconfirma la potencia de los Profile HMM.

La sensación de complejidad y práctica

• Búsqueda con Profile HMM: O(L × T). L es la longitud del perfil, T la longitud de la consulta. Incluso repasando Pfam completo (20.000 perfiles), se tarda minutos en un portátil.
• Predicción de genes: O(T × N²) respecto al tamaño del genoma T y al número de estados N. Incluso para el genoma humano (3 × 10⁹), toma unos días en un solo servidor.

Estas dos escalas definen la sensación práctica. HMMER es una herramienta para portátil, AUGUSTUS es una herramienta para servidor. La aceleración por GPU apenas ha comenzado a incorporarse (el trasfondo es que la vectorización SIMD de HMMER ya estaba bien optimizada, por lo que la optimización solo con CPU era suficiente).

Mapeo informático — Tablero probabilístico de máquinas de estado finito

El concepto de máquina de estado finito tratado en DryBench se expande aquí directamente.

• FSM determinista: el siguiente estado está determinado desde cada estado.
• FSM no determinista: posibles varios estados siguientes (método de selección indeterminado).
• HMM: asigna probabilidades a cada estado siguiente.
• Profile HMM: añade estados de inserción y eliminación al HMM para manejar huecos.

Si se asimila este proceso de expansión gradual, en el futuro podrá interpretar la atención del Transformer como una "expansión completamente conectada de transiciones de estado con pesos". Todos los modelos de secuencia no son más que recombinaciones de estados, transiciones y emisiones.

Ramas hacia el siguiente capítulo — Cierre del Micro Tier

M21 es el §M.4 final del Micro Tier (M01～M30). A partir del próximo capítulo, pasaremos a los algoritmos de ensamblaje.

• Próximo capítulo (M22): Ensamblaje OLC — Overlap-Layout-Consensus. Gráficos de superposición y camino hamiltoniano.
• Dos capítulos después (M23): Grafo de De Bruijn — Ensamblaje por camino euleriano. El estándar de la era de lecturas cortas.
• Tres capítulos después (M24): hifiasm — Ensamblaje diploide PacBio HiFi. El estándar actual de lecturas largas.

Si la serie HMM fue el modelo probabilístico representativo que podía calificarse en Rosalind, la serie de ensamblaje es la aplicación biológica de los algoritmos de grafos. Prepárese para cambiar su cerebro de la rejilla al grafo.

Si desea profundizar más

• Krogh, Brown, Mian, Sjölander, Haussler (1994), Hidden Markov Models in Computational Biology: Applications to Protein Modeling, JMB — El artículo que estableció por primera vez los Profile HMM. Sigue siendo válido incluso hoy.
• Documentación oficial de HMMER: http://eddylab.org/software/hmmer/ — Uso, derivación probabilística y tutoriales completos.
• Formación de EMBL-EBI — Lista de reproducción práctica de Pfam e InterPro. Permite captar la intuición sobre por qué la búsqueda de dominios es tan útil.
• Stanford CS262 — Clase de Gill Bejerano sobre detección de genes (con subtítulos completos). Se observa el diagrama de fases HMM de GENSCAN mediante pizarra.
• Tutorial de AUGUSTUS: https://github.com/Gaius-Augustus/Augustus/wiki — Pipeline práctico de predicción de genes eucariotas.

Pegue cualquier secuencia proteica en el sitio web de Pfam y ejecute la búsqueda. Al ver cómo se traza el mapa de dominios, comprenderá intuitivamente por qué HMM es un marco algorítmico que ha perdurado durante más de medio siglo.