IV Alineamiento mltiple de secuencias Andrs Moreira Departamento

IV Alineamiento múltiple de secuencias Andrés Moreira Departamento de Informática UTFSM

Alineamiento de secuencias Hemos hablado del alineamiento de a pares : se alinean dos secuencias de la mejor manera posible para maximizar identidad y conservación; se busca evaluar similaridad, y posible homología. Sirve para ver si dos secuencias se parecen, significativamente. Sirve para buscar secuencias parecidas a una secuencia dada. . Y también es la base para poder alinear de a más de 2.

Alineamiento múltiple de secuencias Alineamiento múltiple: se alinea un grupo de secuencias ( no es lo mismo que alinearlas a todas de a pares!). El alineamiento múltiple busca identificar letras homólogas en secuencias homólogas. En general se aplica sobre grupos de secuencias que ya se saben homólogas. Aún en sus variantes más simples, es NP-duro.

Alineamiento múltiple de secuencias El alineamiento múltiple ayuda a: Hacer árboles filogenéticos, y estudiar la evolución a nivel molecular.

Alineamiento múltiple de secuencias El alineamiento múltiple ayuda a: Hacer árboles filogenéticos, y estudiar la evolución a nivel molecular. Preguntas del tipo: • ¿De dónde salió la secuencia X? • ¿Qué filogenia tienen las secuencias X, Y, Z? • ¿Este genoma, contiene algún gen que haya llegado por transferencia horizontal? • ¿Hay genes duplicados? • . . . etc

Alineamiento múltiple de secuencias El alineamiento múltiple ayuda a: Caracterizar familias de proteínas identificación de regiones conservadas, funcionalmente importantes (sitios activos de enzimas, sitios de asociación con otras moléculas. . . ) inferencia de función (cuando alguien entra a la familia) inferencia de estructura (global y de regiones específicas) para miembros nuevos de la familia

Alineamiento múltiple de secuencias El alineamiento múltiple ayuda a: Caracterizar familias de proteínas También sirve para buscar miembros nuevos, buscando homólogos en una BD, que sigan el mismo “perfil”.

Perfiles y logos Alineamiento y perfil: A C A T A T A A A T A C G C A T T C G G A T C G A 0. 8 0 0. 4 0. 2 1 0 C 0 0. 6 0. 2 0 0 G 0 0 0. 4 0 0. 2 T 0. 2 0. 4 0 0. 2 0 0. 8

Perfiles y logos Sequence Logo A 0. 8 0 0. 4 0. 2 1 0 C 0 0. 6 0. 2 0 0 G 0 0 0. 4 0 0. 2 T 0. 2 0. 4 0 0. 2 0 0. 8

Perfiles y logos En cada columna tenemos una densidad de probabilidad; podemos calcular su entropía: Si pia es la frecuencia de la letra a en la columna i, entonces H(i) = −Σa pia log pia • Es una forma de evaluar la calidad del alineamiento en esa columna: H(i)=0 todas iguales (la misma letra) H(i) máximo cero información! • Se suelen tomar los gaps como una letra más, aunque hay variaciones.

Perfiles y logos En el logo, se escala cada columna por su entropía (con respecto al máximo, o mejor dicho, a su diferencia con respecto al máximo: columnas más altas reducen más la incerteza). Más útil que

Perfiles y logos También se hacen de proteínas:

Perfiles y logos En ese caso se colorean los aminoácidos por grupos: • Hidrofóbicos: A, V, I, L, M, F, P, W • Básicos: K, R, H • Acídicos: D, E • Polares: S, T, N, Q, G, Y, C

¿Cómo evaluar? ¿Cómo elegimos el alineamiento? Opción 1: el que construye un experto, a ojo. Se usan. Y a veces son excelentes. Opción 2: el que encontramos con un algoritmo. Pero en ese caso necesitamos una función a maximizar. ¿Qué es un buen alineamiento múltiple?

¿Cómo evaluar? En rigor, quisiéramos buscar el árbol filogenético más probable. . . X V U A W C D E B (nos deshacemos de lo desconocido, sumándolo) . . . pero es impracticable. Falta demasiada información (e. g. , las ramas pueden ser de distintos largos), y además es demasiado cómputo.

¿Cómo evaluar? Lo más común: se dan puntajes por columnas del alineamiento, y se busca optimizar la suma de esos puntajes (es la extensión más obvia de lo que se hace en el alineamiento de a pares). Se usan varias opciones. Entre otras: • Cantidad máxima de letras idénticas (trivial, e inútil) • Entropía • “Suma de a pares”: suma el puntaje de cada par de letras de la columna (calculado como al alinear de a pares).

¿Cómo evaluar? Suma de a pares: . . . A. . . S. . . T. . . AA: 4, AS: 1, AT: 0 ST: 1 Puntaje: 4+1+0+1 = 7 Desventaja de la suma de a pares: • No tiene detrás un modelo probabilista razonable. Desventaja de la entropía (y de la suma de a pares): • No tiene detrás un modelo de la filogenia! Se parcha. Se ponen pesos para distintos niveles de semejanza (en alguna medida, recoge filogenia).

Prog. dinámica multidimensional Ok, escogemos una forma de poner puntaje, por columna. ¿Cómo lo optimizamos? Primera idea: extender el método de programación dinámica usado para alinear de a pares. Así como antes avanzábamos por cada una de las secuencias, moviéndonos en un plano, ahora avanzamos por k secuencias. 0 1 1 2 3 4 A -- T G C A A T -- C -- A T G C (0, 0, 0) (1, 1, 0) (1, 2, 1) (2, 3, 2) (3, 3, 3) (4, 4, 4)

Prog. dinámica multidimensional 2 D V W 3 D

Prog. dinámica multidimensional si, j, k = max si-1, j-1, k-1 + (vi, wj, uk) si-1, j-1, k + (vi, wj, _ ) si-1, j, k-1 + (vi, _, uk) si, j-1, k-1 + (_, wj, uk) si-1, j, k + (vi, _ , _) si, j-1, k + (_, wj, _) si, j, k-1 + (_, _, uk)

Prog. dinámica multidimensional El orden es O(2 knk), para k secuencias de largo n. Para una proteína típica, n=100. Más allá de k=4, 5, se vuelve inservible. Para k=7, n=100, requeriría varios terabytes de RAM. No se usa!!

Otras ideas. . . Otra idea: alinear todos los pares de secuencias, y combinar todos esos alineamientos en uno común. Después de todo, un alineamiento múltiple induce alineamientos de a pares: (se proyecta el camino 3 d sobre los lados del cubo)

Otras ideas. . . induce: x: y: z: x: ACGCGG-C; y: ACGC-GAC; AC-GCGG-C AC-GC-GAG GCCGC-GAG x: AC-GCGG-C; z: GCCGC-GAG; y: AC-GCGAG z: GCCGCGAG . . . así que cabría la esperanza de hacer el proceso inverso.

Otras ideas. . . Pero: • Alineamiento óptimo múltiple no induce óptimos de a pares. • Alineamientos óptimos de a pares son difíciles de combinar en un buen alineamiento múltiple. • Alineamientos de a pares pueden ser inconsistentes (es decir, puede no existir uno múltiple que los induzca). no va por ahí la cosa.

Alineamiento progresivo Algo que sí se usa (y mucho): alineamiento progresivo. Idea: alineo primero dos secuencias. Luego alineo el resultado de eso con una tercera secuencia. Y así sucesivamente. ¿En qué orden? Una opción sería tomar las dos más cercanas (evaluando la similaridad de a pares). Luego, la que se acerque más al resultado de eso. Y así, sucesivamente. No se hace.

Alineamiento progresivo CLUSTAL: D C • Calcula la similaridad de cada par. A E B • Construye un “árbol guía”. • Alinea primero lo más cercano, y va subiendo por el árbol. Ojo: se requiere alinear secuencias con alineamientos (y alineamientos con alineamientos!). A X V U W C B D E

Alineamiento progresivo El arbol guía se construye con neighbour joining (un algoritmo que veremos al hablar de filogenia). v 1 v 2 v 3 v 4. 17. 87. 28. 59. 33. 62 Calcular: v 1, 3 = alineamiento(v 1, v 3) v 1, 3, 4 = alineamiento((v 1, 3), v 4) v 1, 2, 3, 4 = alineamiento((v 1, 3, 4), v 2) v 1 v 3 v 4 v 2

Alineamiento progresivo Pueden aparecer nuevos gaps en el camino. S 1 S 3 S 4 S 2 S 4 S 1 S 3 S 2 S 4 Gap nuevo para optimizar el alineamiento de (S 2, S 4) con (S 1, S 3)

Alineamiento progresivo ¿Cómo asignar puntajes al alinear una letra contra un perfil, o un perfil contra otro? • Lo sencillo: escribir la nueva columna completa (con lo que ya llevo alineado en cada lado para esa posición), y evaluar el puntaje de la columna. • Alternativa: comparar sólo con el perfil: No será A contra AC, sino A contra (0. 5, 0, 0). Esto es útil cuando hay muchas secuencias, y cuando se usa un perfil para buscar en una base de datos.

Alineamiento progresivo Pero ojo: si se compara una letra y contra un perfil p (en que p(x) es la frecuencia de la letra x ), y estamos usando la matriz de sustitución Q, entonces el puntaje NO es X pues ese Q se obtuvo como Q(a, b)=log q (a, b), donde q (a, b) representaba la probabilidad de que el alineamiento fuera no aleatorio. Para que siga teniendo sentido sumar, hay que usar y para dos perfiles,

Alineamiento progresivo Clustal es lejos lo más usado hoy en día por los biólogos para alinear secuencias. http: //www. clustal. org/ Viene en dos sabores: Clustal. W (línea de comando), Clustal. X (+GUI). La GUI no es algo menor: alinear es un arte. La W viene de “weighted”: añade pesos en distintos niveles del alineamiento (según la distancia a la raíz), quitándole importancia excesiva a grupos de secuencias muy similares.

Alineamiento progresivo Otras cosas que agrega: • Usa matrices de sustitución distintas, según nivel de semejanza entre secuencias. • La penalización de gaps también la va variando: • Se penalizan menos allí donde ya hay gaps • Se penalizan más cerca de donde ya hay gaps. • Además los modifica consideraciones bioquímicas. • Hay posibilidades para que el usuario predetermine el árbol guía, o para realinear tramos de la secuencia.

Alineamiento progresivo Problemas del alineamiento progresivo: • Se casa mucho con los alineamientos con que parte; esos gaps y alineación se acarrean después. Si las secuencias no son muy parecidas, puede ser malo. Respuesta: métodos iterativos, y otros. • Es global! La semejanza puede estar en trozos de las secuencias.

Alineamiento progresivo Lo que debiera ser Lo que da Clustal Respuestas a eso: 1. Cortar antes los tramos similares y alinearlos. 2. Usar algoritmos que se preocupan de eso.

Alineamiento iterativo Iterar: para deshacerse de la dependencia de lo inicial. En Clustal es posible iterar de varias formas: • Usar el alineamiento para construir un árbol filogenético, y usar eso como árbol guía en un nuevo alineamiento. [Eso requiere pasar por otro software]. • Realinear las secuencias ya alineadas (con sus gaps). +Hacerlo a partir de varios sets de parámetros distintos. Ir viendo qué zonas son importantes. Realinear a mano. Etc. (El manual recomienda generar al menos 10 alineamientos).

Alineamiento iterativo Otros métodos iterativos: • DIALIGN: pone énfasis en alineamientos locales entre las secuencias. • MUSCLE: va revisando el árbol guía

Alineamiento iterativo: MUSCLE: va revisando el árbol guía. 1) Parte como Clustal 2) Sobre la marcha va revisando: acaso el árbol guía ha cambiado, y acaso alguna parte se puede realinear mejor.

Alineamiento iterativo: MUSCLE ¿Ha cambiado el árbol? Calcula la similaridad de los alineamientos de a pares inducidos por el actual alineamiento múltiple. Si eso da un árbol distinto al que está usando, realínea según eso, y ve si el puntaje general mejora.

Alineamiento iterativo: MUSCLE ¿Hay algo que realinear? Prueba a cortar el árbol. S T U X Z Borra arista S T U X Z Realínea los perfiles resultantes de cada árbol

Alineamiento iterativo: alg. evolucionarios Métodos evolucionarios: Lo que se muta es la ubicación de los gaps (un alineamiento no es otra cosa que una distribución de gaps!). • Algoritmos genéticos: SAGA • Simulated Annealing: MSASA Variante “hill climbing”: muta en direcciones provechosas.

Alineamiento iterativo: alg. evolucionarios Cruce de alineamientos Una mutación . . . etc

Alineamiento : DIALIGN II Énfasis en los alineamientos locales. • Busca buenos (=buen P-value) alineamientos locales (de a pares). • Los va re-evaluando según qué tan compatibles son unos con otros. . . • Hasta que obtiene un buen set de alineamientos locales ordenados y coherentes. Y con eso, alínea.

Alineamiento : DIALIGN II Es especialmente útil para: • similitudes que son sólo locales • secuencias muy lejanas Y es una gracia el no tener que tratar los “gaps” de manera ad-hoc.

Alineamiento: nota al margen • Hay muuuchos algoritmos y software de MSA. • Casi toda técnica de optimización conocida por el hombre, se ha aplicado. Para bien o para mal. • Clasificar los algoritmos es difícil; combinan características de varios lados (estocásticos/deterministas, progresivos/simultaneos, exactos/heurísticos, iterativos o no. . . ), y mutan de año en año. • Aquí daré sólo las ideas principales de algunos; meternos al detalle sería interminable.

Alineamiento: nota al margen • Esa diversidad, y los avances permanentes, son importantes. El alineamiento (de a pares y múltiple) es ubicuo en bioinformática. Dos reviews instructivos, legibles y no muy largos: • Multiple sequence alignment, R. C. Edgar & S. Batzoglou (2006), Curr Op Struct Biol 16(3) 368 -373; doi: 10. 1016/j. sbi. 2006. 04. 004 Disponible en http: //ai. stanford. edu/~serafim/Publications/2006_MSA_COSB. pdf • Recent evolutions of multiple sequence alignment algorithms, C. Notredame (2007), PLo. S Comput Biol 3(8) e 123, doi: 10. 1371/journal. pcbi. 0030123 Disponible ahí.

Alineamiento: T-Coffee Idea: Generar una biblioteca de alineamientos (globales y locales) de a pares entre las secuencias. Tratar de hacer el alineamiento múltiple que sea consistente con esa biblioteca. • Primero genera una biblioteca, haciendo alineamientos de a pares, global y local, entre todas las secuencias.

Alineamiento: T-Coffee • Asigna puntaje a los matches de esos alineamientos. En un proceso llamado extensión de la biblioteca, usa terceras secuencias para modificar los puntajes en un par dado. • Al final tiene, sobre los alineamientos de a pares, puntajes que reflejan la compatibilidad con otros alineamientos. • Y usa eso para hacer la prog. din. , en alineamiento progresivo.

Alineamiento: T-Coffee

Alineamiento: T-Coffee + Presta más atención a homologías locales que Clustal. + Permite incluir información de otras fuentes! A la biblioteca primaria se le pueden agregar alineamientos obtenidos con Clustal, Dialign, información biológica, lo que sea. Es de lo mejorcito que ha aparecido.

Alineamiento: POA: Partial Order Alignment Idea: no resumir los alineamientos en perfiles, sino que dibujar un grafo. Y alinear vía programación dinámica, pero contra ese grafo.

Alineamiento: POA La idea no es mala. . . Y permite conservar la información de un alineamiento, en una representación compacta.

Alineamiento: POA Alínea permitiendo agregar caminos alternativos por sobre la matriz 2 d. Y ayuda a visualizar estructuras compartidas parcialmente:

Alineamiento: POA Problemas: • Se vuelve poco práctico si hay demasiada bifurcación: con sets muy grandes de secuencias, o secuencias evolutivamente distantes. • Con frecuencia lo queremos es generalizar: discernir los motivos que caracterizan una familia de proteínas, o una secuencia de regulación en el DNA, etc.

Alineamiento: benchmarks Como ya dijimos, los benchmarks aquí son humanos: sets curados a mano, a partir de amplia información filogenética, estructural, etc. El más conocido es Bali. Base: • Son “correctos” • Tiene sets “difíciles”, agrupados según naturaleza de la dificultad. • Lo más difícil es donde hay mucha distancia evolutiva: la twilight zone! http: //www. people. virginia. edu/~wrp/papers/ismb 2000. pdf

Alineamiento: benchmarks Bali. Base, Thompson et al, NAR, 1999, PROBLEMA Descripción 1. Grupo de secuencias cercanas 2. Una secuencia lejana 3. Dos grupos distantes 4. Inserción interna 5. Inserción en un extremo

Alineamiento: benchmarks

Alineamiento: benchmarks

Alineamiento: benchmarks

Alineamiento: benchmarks PROBLEMA Estrategia Clustal. W, T-coffee, MSA T-Coffee Dialign T-Coffee

Alineamiento: benchmarks. . . aunque, de otra fuente: Pero parece ser que: Dialign para cosas lejanas, T-Coffee para lo demás. Como regla a primera vista, claro. Toda la fauna que anda afuera tiene razón de ser, según el caso. Y en último término, T-Coffee puede combinar.

Alineamiento: benchmarks Lista de algunos: MAFFT POA MUSCLE Progressive/Simulataneous Progressive/Iterative www. biophys. kyoto-u. jp/katoh www. bioinformatics. ucla. edu/poa www. drive 5. com/muscle/

Alineamiento: benchmarks En todo caso, la mayoría de los biólogos sigue usando Clustal. Volvamos a eso de generalizar. Extracción de motivos patrones.

Motivos Motivo ( motif ): un patrón que puede estar presente en una secuencia. Puede ser más o menos estricto. Puede reflejar un sitio activo de una proteína, o un punto de interacción entre moléculas.

Motivos Con los motivos se necesita: • Un lenguaje para describirlos: para saber que tipos de motivos se pueden encontrar. • Una forma de evaluarlos: una función objetivo, que defina los motivos que son “interesantes”. • Algoritmos para determinarlos, y para localizarlos. Determinarlos: extraer un motivo a partir de un conjunto de secuencias. Localizarlos: dado un motivo, saber reconocerlo.

Motivos Lenguaje para describirlos: Se usan varios. • Expresiones regulares (no necesito explicárselas) • Perfiles (ya vistos, es la matriz de frecuencias x ubicación) • Modelos markovianos ocultos (HMM)

Cadenas de Markov En una cadena de Markov, el sistema puede estar en una serie de posibles estados (en un momento dado), y entre un instante y el siguiente cambia de estado de acuerdo a una cierta matriz de probabilidad. A Componentes: • Set de estados. • Una distribución inicial de probabilidad sobre los estados • La matriz de transición. B C

Cadenas de Markov Estados: día soleado, lluvioso, nublado. Matriz de transición: probabilidades para hoy, dado el estado de ayer Una distribución inicial.

Cadenas de Markov Hay un montón de teoría sobre cadenas de Markov; se usan mucho. Aquí nos interesa un tipo un poco más complicado de proceso Markoviano: 1. El grafo será acíclico: no hay “vuelta atrás” en los estados. 2. En cada paso, no observamos el estado, sino un estado “de salida”, que depende del estado en que estemos, con cierta probabilidad. 3. Como no observamos directamente el estado en que estamos, se habla de modelos ocultos de Markov.

HMM "Hidden Markov Model" (HMM). Ejemplo (no acíclico): un casino tiene un dado bueno, y un dado cargado (con P(6)=½, el resto uniforme). Cambia del bueno al cargado con probabilidad 0. 05, y del cargado al bueno con probabilidad 0. 1. 0. 95 Bueno 0. 05 0. 1 0. 90 Cargado Eso sería como cadena de Markov. Pero supongamos que no sabemos que dado se está usando. Sólo vemos el resultado de las tiradas.

HMM 0. 95 Bueno 0. 05 0. 1 0. 90 Cargado HMM: la variable que 0. 90 0. 95 observamos depende de 1: 1/6 1: 1/10 estados internos que no vemos. 2: 1/6 0. 05 2: 1/10 Si aparece 312453666641 ¿cuál es la secuencia de estados internos más probable? BBBBBBCCCCCC 3: 1/6 3: 1/10 4: 1/6 4: 1/10 5: 1/6 0. 10 5: 1/10 6: 1/6 6: 1/2 Bueno Cargado

HMM Esquema general de un HMM (hidden markov model): Componentes: • Set de estados ocultos (aquí, x 1, x 2, x 3) • Set de estados de salida, o “alfabeto” (aquí, y 1, y 2, y 3) • Distribución inicial • Matriz de transición entre estados (Pij) • Matriz de output (prob. de estados de salida, dados los ocultos) (Qij)

HMM Problemas usuales con HMM: • Dado un HMM y una secuencia observada, ¿cuál es la probabilidad de que la secuencia sea producida por ese HMM? • Dado un HMM y una secuencia observada, ¿cuál es la secuencia más probable de estados internos que la genera? • Dado un conjunto de secuencias observadas, ¿cuál es el mejor HMM? Algoritmos: forward, Viterbi, Baum-Welch

HMM Aplicando a secuencias genéticas: Los estados internos son las columnas “relevantes” de un alineamiento. Si no hay gaps, entonces se reduce a un perfil: M 1 M 2 “Motif HMM” M 3 M 4 M 5

HMM Los “problemas usuales” se convierten en: • Dado un perfil HMM, esta secuencia, ¿es probable que pertenezca al club? • Dado un perfil HMM, esta secuencia, ¿cómo la alíneo con las previas? • Dado este conjunto de secuencias, ¿cómo lo resumo con un perfil tipo HMM? Todas son importantes!

HMM ACA TCA AGA ACC --ACT C---G-- ATG ATC AGC ATC inserción Prob. de transición Prob. de output Sec. de consenso: ACAC - - ATC P (ACACATC) = 0. 8 x 1 x 0. 8 x 0. 6 x 0. 4 x 0. 6 x 1 x 0. 8 = 4. 7 x 10 -2

HMM Hay que distinguir entre la secuencia observada, y el camino que se usa para generarla (que no se ve!). Para una secuencia observada, pueden haber varios caminos alternativos. Si conocemos la secuencia observada a=a 1, . . . a. L, y también el camino, s=s 1, . . . , s. L (supondremos que el carácter vacío es parte del alfabeto, para tener la misma longitud), entonces la probabilidad conjunta es O con la notación de antes,

HMM Pero si no conocemos los estados internos (s), entonces sólo podemos calcular donde sumamos sobre todos los s que pueden generar a. Nota: se suele agregar un estado de inicio y uno de término, ambos “mudos” (o bien, con outputs especiales, “^” y “$”, por ejemplo).

HMM Este es un caso fácil: cada a tiene un único s. P(ACA---ATG) P(TCAACTATC) P(ACAC--AGC) P(AGA---ATC) P(ACCG--ATC) = = = 3. 3 x 10 -2 0. 0075 x 10 -2 1. 2 x 10 -2 3. 3 x 10 -2 0. 59 x 10 -2 4. 9 3. 0 5. 3 4. 9 4. 6 Consenso: 6. 7 P(ACAC--ATC) = 4. 7 x 10 -2 6. 3 P(ACA---ATC) = 13. 1 x 10 -2 Excepcional: P(TGCT--AGG) = 0. 0023 x 10 -2 -0. 97 ACA TCA AGA ACC --ACT C---G-- ATG ATC AGC ATC Usamos log odds para tener cifras más fáciles. Odds: prob. / prob. del caso nulo. Nulo en este caso: secuencia aleatoria de largo L P=(0. 25)L Log Odds(a) = log(P(a)/(0. 25)L) Valor positivo Más probable que modelo nulo!

HMM: Viterbi Dada a y un HMM M, ¿cuál es la s más probable? Respuesta: algoritmo de Viterbi. Programación dinámica, similar a Smith-Waterman. También rellenamos una matriz. En un eje ponemos la secuencia observada, en el otro los estados ocultos. En la posición (x, y) la matriz contiene la máxima probabilidad de un camino que termine en el estado oculto x-ésimo con la observación ay. (pasando por a 1, . . . ay-1).

HMM: Viterbi ^ ACG 0 2 3 Secuencia observada y Estados ocultos . . . y vamos guardando el argmax, para poder reconstruir la ruta al llegar al final, desde (3, $). 1 $ x ^ 0 (inicio) 1 1 0 2 0 3 (fin) 0 A C V(2, c) A C G $

HMM: Forward Dada a, y un HMM M, ¿cuánto vale P( a/M )? Respuesta: algoritmo Forward. Idéntico a Viterbi, pero reemplazamos “max” por “sum”. En la posición (x, y) la matriz contiene la suma de las probabilidades de los caminos que terminan en el estado oculto x-ésimo con la observación ay. (pasando por a 1, . . . ay-1).

HMM: Decodificación posterior Posterior decoding contesta a la pregunta: dada a y un HMM M que generó a, ¿cuál es la probabilidad de que el estado oculto x haya generado la letra ay? Esto sirve, por ejemplo, para saber cómo se corresponden los aminoácidos de una proteína con posiciones estructurales conservadas. Respuesta: 1) Se hace un algoritmo Backward, análogo a Forward. La variable B(x, y) será la suma de probalidades de caminos que parten en el estado x con la observación ay. 2) La respuesta a la pregunta será Demostración: ejercicio

HMM: Construcción Queda la pregunta del millón: ¿cómo construímos el HMM? Hay dos problemas: • Determinar la arquitectura (estados ocultos y sus conexiones) • Determinar los parámetros (probabilidades de transición y de output) En el caso más general, la arquitectura es un arte. Para alineamientos, se puede usar esta forma estándar. ¿Cuántos nodos “M”? La cantidad de columnas mayoritariamente “no gap” o el promedio de longitudes

HMM: Construcción A veces las secuencias ya están alineadas, y queremos un HMM que represente el alineamiento. En ese caso, si las secuencias no son muchas, puede hacerse “a mano”: contamos las frecuencias de letras, vemos los desvíos. . Donde haya gaps, ponemos insert o delete según la cantidad de gaps (>50% insert). La corrección que se suele hacer, para generalizar: agregar probabilidades pequeñas de emisión para letras ausentes. (Pondero las frecuencias observadas con las frecuencias “a priori”, dando un cierto peso a cada parte). ACA TCA AGA ACC --ACT C---G-- ATG ATC AGC ATC

HMM: Baum-Welch Si las secuencias no están alineadas, entonces es más complicado. Recibo un conjunto de secuencias, {a 1, . . a. N}. Si supiera cuál es el camino que genera cada una, podría estimar las probabilidades de transición y de emisión contando la frecuencia con que ocurren. Notemos Aij la cantidad de veces que usan la transición i j, y Bij la cantidad de veces que la letra j se emite en el estado i. Entonces podría estimar

HMM: Baum-Welch Por otro lado, si tuviera los parámetros de la HMM, entonces podría saber la probabilidad de los caminos que generan las secuencias, y por lo tanto la frecuencia con que usarán cada emisión y cada transición. En particular, la probabilidad de que se una transición de x a w sería Y la probabilidad de emitir una letra e desde el estado x sería

HMM: Baum-Welch Entonces, si tuviera los parámetros, sabría las probabilidades de los caminos, y si tuviera los caminos, podría estimar los parámetros. algoritmo EM “Expectation. Maximization” • Se aplica cuando desconocemos parámetros de un modelo que genera ciertos datos, y tampoco conocemos una parte de eso que se generó. E-step: estimar datos faltantes Parámetros Información faltante M-step: estimar parámetros (maximizando verosimilitud)

HMM: Baum-Welch En este caso, los parámetros desconocidos son los del HMM, y de lo que se genera, conocemos la secuencia pero no el camino. E-step: estimar datos faltantes: • Partimos con alguna “adivinanza” (posiblemente uniforme) sobre los parámetros. • Iteramos hasta converger. Parámetros: Py. Q Información faltante: Caminos para {ak} M-step: estimar parámetros:

HMM: Construcción Entonces, para construir un HMM con secuencias no alineadas: • Definir arquitectura (longitud); por lo general se usa la longitud promedio. • Aplicar Baum-Welch. Ojalá con varias condiciones de inicio distintas, para evitar óptimos locales. • Una vez obtenido el HMM, alinear las secuencias encontrando para c/u su alineación con el HMM (vía Viterbi).

HMM presentes localmente En principio, una secuencia se alínea globalmente contra un HMM. Pero podemos permitir alineamientos locales! Q emite letras no alineadas. Viterbi sigue funcionando.

HMM Sirven para: • Representar un perfil, de un alineamiento nuevo, o una familia de secuencias o subsecuencias (lugares en DNA, familias de proteínas, dominios funcionales en proteínas. . . ) Se puede pensar en los HMM como secuencias “generalizadas”, con probabilidades para las variaciones respecto al patrón. • Crear alineamientos (usando Baum-Welch). Pero sólo cuando hay hartas secuencias; tienen demasiados parámetros para ser entrenados con poca información.

HMM Sirven para: • Alinear una secuencia nueva c/r al HMM (usar Viterbi). • Clasificar una secuencia nueva: usando Forward evaluamos la probabilidad de que pertenezca a distintas familias (cada una caracterizada por su HMM). • Buscar nuevos miembros (idem, Forward, pero varío la query: es el HMM, no la secuencia).

HMM: software • SAM (Sequence Alignment and Modeling System): el nombre no menciona HMM, pero se trata de eso. • Meta-MEME: “is a software toolkit for building and using motif-based hidden Markov models of biological sequences”. La parte “MEME” es “Multiple EM for Motif Elicitation”. • HMMER: búsqueda en base de datos mediante HMM.

HMM: datos Pfam: base de datos de dominios y familias de proteínas. • Un MSA curado a mano se usa como “semilla” Pfam-A. • Con eso se crea un HMM. • Se usa el HMM para agregar otros miembros de manera automática (el umbral de puntaje se escoge a mano) Pfam-B. • Contiene más de 9000 familias; el 75% de las proteínas conocidas tienen al menos un match en Pfam.

+ datos Bases de datos relacionadas: hay varias que también agrupan dominios y familias. • Pro. Site: Curada a mano, altamente confiable. Usa expresiones regulares, y ahora también perfiles. • Pro. Dom: Usa secuencias de consenso para la representación. Generada haciendo clustering sobre secuencias de Swiss. Prot y Tr. Embl. • Inter. Pro: integra a varias (Pfam, Pro. Site, Pro. Dom. . . )

Algo más sobre motivos • Hay softwares que se dedican a buscar motivos conocidos en una secuencia (consultan Pro. Site, Pfam, etc. ). • Otro problema es encontrar un motivo nuevo en un conjunto de secuencias dado. Una opción es hacer un alineamiento múltiple, y ver las regiones mejor alineadas. • Eso falla si las secuencias sólo comparten el motivo, y lo tienen en lugares distintos. Entonces podríamos pensar en algoritmos locales, como Dialign. • Hay algoritmos más específicos, vía EM, y Gibbs sampling.

Idea del Gibbs sampling • N secuencias, largos posiblemente distintos. • Una longitud L para el motivo desconocido. • Una adivinanza inicial, con una ubicación para el motivo en cada secuencia. Eso da un MSA (aunque pueda ser pésimo). Itero los siguientes pasos: • elimino una de las ubicaciones • escojo una nueva ubicación en esa secuencia, de manera probabilista

Idea del Gibbs sampling Para escoger la nueva ubicación, le asigno a cada posible nueva ubicación una probabilidad que refleja qué tan bien esa ubicación le vendría al MSA actual de las otras N-1 ubicaciones: es el producto, sobre i=1. . L, de la frecuencia de la letra candidata en esa posición del MSA, dividido por la frecuencia “background”. Puede agregarse un paso que haga variar L, si conviene.

Idea del Gibbs sampling 5’- TCTCCACGGCTAATTAGGTGATCATGAAAAAA TGAAAAATTCATGAGAAAAGAGTCAGACATCGAAACAT …HIS 7 5’- ATGGCAGAATCACTTTAAAACGTGGCCCCACCCGCTGCACCCTGTGCATTTTGTACGTTACTGCGAAATGACTCAACG …ARO 4 5’- CACATCCAACGAATCACCGTTATCGTGACTCACTTTCGCATC GCCGAAGTGCCATAAAAAATATTTTTT …ILV 6 5’- TGCGAACAAAAGAGTCATTACAACGAGGAAATAGAAGAAAATGAAAAATTTTCGACAAAATGTATAGTCATTTCTATC …THR 4 5’- ACAAAGGTACCTTCCTGGCCAATCTCACAGATTTAATATA GTAAATTGTCATGCATATGACTCATCCCGAACATGAAA …ARO 1 5’- ATTGACTCATTTTCCTCTGACTACTACCAGTTCAAAATGTTAGAGAAAAATAGAAAAGCAGAAAAAATAA …HOM 2 5’- GGCGCCACAGTCCGCGTTTGGTTATCCGGCTGACTCATTCTGACTCTTTTTTGGAAAGTGTGGCAT GTGCTTCACACA …PRO 3 TGAAAAATTC GACATCGAAA GCACTTCGGC GAGTCATTAC GTAAATTGTC CCACAGTCCG TGTGAAGCAC ***** (En este ejemplo se le da aún más libertad a las ubicaciones: pueden estar en cualquier hebra, y en cualquier secuencia –de modo que puede haber dos en una secuencia, y ninguna en otra. )

Idea del Gibbs sampling Agrego? 5’- TCTCCACGGCTAATTAGGTGATCATGAAAAATTCATGAGAAAAGAGTCA GACATCGAAACAT …HIS 7 5’- ATGGCAGAATCACTTTAAAACGTGGCCCCACCCGCTGCACCCTGTGCATTTTGTACGTTACTGCGAAATGACTCAACG …ARO 4 5’- CACATCCAACGAATCACCGTTATCGTGACTCACTTTCGCATC GCCGAAGTGCCATAAAAAATATTTTTT …ILV 6 5’- TGCGAACAAAAGAGTCATTACAACGAGGAAATAGAAGAAAATGAAAAATTTTCGACAAAATGTATAGTCATTTCTATC …THR 4 5’- ACAAAGGTACCTTCCTGGCCAATCTCACAGATTTAATATA GTAAATTGTCATGCATATGACTCATCCCGAACATGAAA …ARO 1 5’- ATTGACTCATTTTCCTCTGACTACTACCAGTTCAAAATGTTAGAGAAAAATAGAAAAGCAGAAAAAATAA …HOM 2 5’- GGCGCCACAGTCCGCGTTTGGTTATCCGGCTGACTCATTCTGACTCTTTTTTGGAAAGTGTGGCAT GTGCTTCACACA …PRO 3 TGAAAAATTC GACATCGAAA GCACTTCGGC GAGTCATTAC GTAAATTGTC CCACAGTCCG TGTGAAGCAC Mejora? TCTCCA GACATCGAAA GCACTTCGGC GAGTCATTAC GTAAATTGTC CCACAGTCCG TGTGAAGCAC

Idea del Gibbs sampling 5’- TCTCCACGGCTAATTAGGTGATCATGAAAAATTCATGAG AAAAGAGTCAGACATCGAAACAT …HIS 7 5’- ATGGCAGAATCACTTTAAAACGTGGCCCCACCCGCTGCACCCTGTGCATTTTGTACGTTACTGCG AAATGACTCAACG …ARO 4 5’- CACATCCAACGAATCACCGTTATCG TGACTCACTTTCGCATCGCCGAAGTGCCATAAAAAATATTTTTT …ILV 6 5’- TGCGAACAAAAGAGTCATTACAACGAGGAAATAGAAGAAAATGAAAAATTTTCGACAAAATGTATAGTCATTTCTATC …THR 4 5’- ACAAAGGTACCTTCCTGGCCAATCTCACAGATTTAATATAGTAAATTGTCATGCATA TGACTCATCCCGAACATGAAA …ARO 1 5’- ATTGACTCATTTTCCTCTGACTACTACCAGTTCAAAATGTTAGAGAAAAATAGAAAAGCAGAAAAAATAA …HOM 2 5’- GGCGCCACAGTCCGCGTTTGGTTATCCGGC TGACTCATTCTGACTCTTTTTTGGAAAGTGTGGCATGTGCTTCACACA …PRO 3 AAAAGAGTCA AAATGACTCA AAGTGAGTCA AAAAGAGTCA GGATGAGTCA AAATGAGTCA GAATGAGTCA AAAAGAGTCA (Es el motivo de un sitio de regulación en genes de S. cerevisiae)

Gibbs sampling, EM The Gibbs Motif Sampler Homepage: http: //bayesweb. wadsworth. org/gibbs. html software bajable, y también servicio online Otro servidor: http: //www. cbs. dtu. dk/biotools/Easy. Gibbs/

Idea con EM Los métodos EM no escogen ubicación ni hacen un MSA: • su parámetro a iterar es un motivo • le asignan una probabilidad a cada ubicación, según ese motivo • recalculan el motivo como un perfil de todas las ubicaciones, ponderadas por esas probabilidades Hay diversas variaciones sobre estas ideas (EM y Gibbs).

Dominios, motivos, familias • Motivos: patrones locales. Pueden ser sitios activos de proteínas, dominios funcionales, regiones reconocibles de DNA, etc. • Dominios: los segmentos funcionales en las proteínas (que son relativamente modulares). Un dominio suele estar caracterizado por un motivo.

Dominios, motivos, familias • Familia: una familia de proteínas comparte la mayoría de sus dominios; puede haber variaciones. Y dentro de “superfamilias” de proteínas, más variación aún. Los motivos salen de MSA, pero también son útiles para construir MSA, si detectamos primero motivos conocidos en las secuencias, y los usamos para anclar el alineamiento.

Proteínas con múltiples dominios Alinear bien proteínas multi-dominio es un problema. Especialmente si lo queremos es no sólo identificar dominios, sino relación filogenética.

Proteínas con múltiples dominios En este caso, las proteínas alineadas en A tienen ancestro común, pero no las alineadas en B. Sin embargo las de B resultan similares por culpa de los dominios IG compartidos (IG es un dominio “promiscuo”, que suele insertarse en proteínas). ¿Cómo distinguir?

Proteínas con múltiples dominios Una aproximación: Sequence Similarity Network Reveals Common Ancestry of Multidomain Proteins N. Song et al, PLo. S Comput Biol. 2008 4(5): e 1000063. 10. 1371/journal. pcbi. 1000063 Estudian el grafo de las secuencias similares, dentro de una base de datos, a dos secuencias dadas. Según la forma de ese grafo, se distingue la ascendencia común vs. la homología sólo de un trozo.

Proteínas con múltiples dominios Los vecinos comunes son más que los vecinos únicos, en el caso de verdadera homología de las secuencias. En el caso en que sólo hay homología de algunos dominios, resulta al revés. Definen un coeficiente de “correlación de vecindades”, que discrimina bien entre homólogos verdaderos y los que no.