Repositorio Institucional de Documentos

Tesis de M´aster M´aster en Ingenier´ıa de Sistemas e Inform´atica Curso 2010/2011 An´alisis filogen´etico molecular: Dise˜no e implementaci´on de algoritmos escalables y fiables y verificaci´on autom´atica de propiedades de una filogenia. Autor: Jorge ´ Alvarez Jarreta Bajo la direcci´on de: Elvira Mayordomo C´amara Gregorio de Miguel Casado Departamento de Inform´atica e Ingenier´ıa de Sistemas ´ Area de Lenguajes y Sistemas Inform´aticos Escuela de Ingenier´ıa y Arquitectura Universidad de Zaragoza Septiembre de 2011 An´alisis filogen´etico molecular: Dise˜no e implementaci´on de algoritmos escalables y fiables y verificaci´on autom´atica de propiedades de una filogenia. RESUMEN La filogen´etica es la ciencia que estudia las relaciones entre organismos bas´andose en lo cercano que est´an unos de otros. La forma m´as visual y conveniente de representar estas relaciones evolutivas entre un grupo de organismos es a trav´es de los ´arboles filogen´eticos. La secuenciaci´on de cadenas biol´ogicas es un proceso mecanizado que permite obtener las secuencias biol´ogicas que se utilizar´an para construir estos ´arboles, que se desarrollan a partir de modelos evolutivos: modelos matem´aticos que intentan explicar de la forma m´as fiel posible la evoluci´on real de dichas secuencias. El an´alisis filogen´etico es un proceso formado por distintas etapas, que pueden variar seg´un los objetivos, pero cuya finalidad es siempre la misma: poder reconstruir el ´arbol filogen´etico. Estas etapas pueden incluir: estudio de modelos evolutivos, an´alisis estad´ıstico, alineamiento de secuencias, . . . Esta tesis de m´aster tiene como objetivo principal el desarrollo del criterio te´orico y herramientas pr´acticas que permitan llevar a cabo un an´alisis filogen´etico completo. Los procesos de secuenciaci´on de cadenas biol´ogicas no est´an exentos de errores, los cuales pueden aparecer en cualquier parte de la secuencia. Hasta la fecha, todo proceso de verificaci´on de la secuencia requiere la actuaci´on manual de un experto en el campo, lo que, sin duda, es un proceso muy costoso. Por otro lado, actualmente el coste computacional limita de forma pr´actica tanto la realizaci´on de filogenias extensivas (tratando miles y decenas de miles de secuencias) como la aplicaci´on de modelos evolutivos m´as generales, interesantes y explicativos que el modelo uniforme (que ´unicamente se utilizan en tama˜nos de problema reducidos). Las aportaciones de este trabajo abordan tres aspectos espec´ıficos: por un lado el desarrollo e implementaci´on de un sistema de inferencia filogen´etica que concentra varios m´etodos sobre an´alisis de secuencias y estudio de filogenias que no se hab´ıan unido hasta el momento; por otro lado el desarrollo e implementaci´on de una aplicaci´on para la detecci´on autom´atica de errores de las cadenas obtenidas en los procesos de secuenciaci´on; y por ´ultimo el estudio te´orico de nuevos algoritmos para la caracterizaci´on de problemas entre aquellos que se consideran como no resolubles en la actualidad. Los resultados de todo este trabajo han concluido en la creaci´on de dos art´ıculos (uno publicado y otro en fase de revisi´on). Adicionalmente, hay un tercero que est´a en fase de desarrollo, lo que refleja el amplio inter´es en estos temas por la comunidad cient´ıfica y su utilidad pr´actica. I Agradecimientos En este rinconcito tan personal, dentro de un documento lleno de tanto formalismo y palabras t´ecnicas, me gustar´ıa destacar a unas cuantas personas que me han ayudado de forma especial a lo largo de este trabajo. En primer lugar, como no pod´ıa ser de otra manera, me gustar´ıa agradecer a Elvira toda su paciencia y dedicaci´on que ha tenido conmigo durante este ´ultimo a˜no (el resto ya se lo agradec´ı en el PFC, tampoco quiero ser pelota), que incluso cuando sus proyectos personales han ocupado la mayor parte de su tiempo, ha conseguido sacar lo posible para ayudarme. En Goyo he descubierto un gran amigo, y me gustar´ıa agradecerle el haberme aceptado en su vida privada cuando ya me sent´ıa parte de la perspectiva laboral (y no ten´ıa por qu´e extenderse). A Jose Manuel me gustar´ıa agradecerle, sobre todo, esas charlas que hemos tenido, en las que he podido absorber algunos de esos conocimientos que solo la edad (¡que no vejez!) nos va dando. A Nacho y Roberto los meto en el mismo saco, ya que son mis camaradas, como Nacho suele decir, y me gustar´ıa agradecerles las charlas, la ayuda y el apoyo que he recibido de ellos. No me gustar´ıa olvidarme de Eduardo, que me ha hecho recordar m´as de una vez que cuando hablo de trabajo con gente de otros campos he de dejar de “hablar en chino” para que pueda existir una comunicaci´on completa. Alej´andome un poco de la reci´en nacida Escuela de Ingenieros y Arquitectos, me gustar´ıa dar las gracias a mis padres por su inestimable apoyo, cari˜no y, c´omo no, sus broncas, que m´as de una vez me han venido bien para poner los pies en la tierra. A mis amigos Antonio y Fronti, quienes han tenido que aguantar m´as de una vez mis “rollos” bioinform´aticos. Y, muy especialmente, a mi novia, Laura, por estar siempre ah´ı, y confesarle mi admiraci´on secreta por su dedicaci´on al trabajo, que me ha servido muchas veces de inspiraci´on para querer superarme a m´ı mismo. No quiero olvidar a quienes me han concedido una ayuda econ´omica para poder continuar mi labor investigadora. Por ello, me gustar´ıa agradecer al Instituto de Investigaci´on e Ingenier´ıa de Arag´on (I3A), al Ministerio de Educaci´on, al Departamento de Ciencia, Tecnolog´ıa y Universidad del Gobierno de Arag´on y al Fondo Social Europeo por ambas becas que he percibido durante este ´ultimo a˜no ([AP2008-03447] y [grupo GISED T27]). Me gustar´ıa terminar citando una frase de An´onimo, el autor inmortal e incesante de tantos escritos y dichos, “la prosperidad hace amistades, y la adversidad las prueba”. V ´ Indice general 1. Introducci´on 1 1.1. Contexto del trabajo . . . . . . . . . . . . . . . . . . . . . . . 1 1.2. Objetivos ............................. 1 1.3. Estructura de la memoria . . . . . . . . . . . . . . . . . . . . 2 2. Biolog´ıa y bioinform´atica 5 2.1. Fundamentos biol´ogicos . . . . . . . . . . . . . . . . . . . . . 5 2.2. Introducci´on a la bioinform´atica . . . . . . . . . . . . . . . . 6 2.2.1. Filogen´etica . . . . . . . . . . . . . . . . . . . . . . . . 6 2.2.2. Super´arboles . . . . . . . . . . . . . . . . . . . . . . . 8 3. Flujos de trabajo con selecci´on de modelos: una aproximaci´on multilocus al an´alisis filogen´etico 9 3.1. Introducci´on............................ 9 3.2. Descomposici´on del problema . . . . . . . . . . . . . . . . . . 11 3.3. Dise˜no del flujo de trabajo . . . . . . . . . . . . . . . . . . . . 13 3.4. Aspectos de implementaci´on . . . . . . . . . . . . . . . . . . . 15 3.5. Resultados y an´alisis de rendimiento . . . . . . . . . . . . . . 16 4. Detecci´on autom´atica de errores de secuenciaci´on a partir de informaci´on filogen´etica 19 4.1. Introducci´on............................ 19 4.2. Detecci´on de errores de secuenciaci´on . . . . . . . . . . . . . . 22 4.3. Pruebas y resultados . . . . . . . . . . . . . . . . . . . . . . . 23 4.3.1. Estudio de comportamiento . . . . . . . . . . . . . . . 24 4.3.2. Estudio de rendimiento . . . . . . . . . . . . . . . . . 28 5. Complejidad param´etrica en bioinform´atica: filogenias casi perfectas 29 5.1. Introducci´on............................ 29 5.2. Preliminares............................ 30 5.2.1. Notaci´on ......................... 31 5.2.2. Algoritmos PP y NPP . . . . . . . . . . . . . . . . . . 31 5.3. An´alisis de complejidad . . . . . . . . . . . . . . . . . . . . . 33 5.4. Propuesta............................. 34 5.4.1. Modelo evolutivo . . . . . . . . . . . . . . . . . . . . . 34 5.4.2. Funci´on de penalizaci´on . . . . . . . . . . . . . . . . . 35 5.4.3. Lemas........................... 36 6. Conclusiones 37 VII 2 Biolog´ıa y bioinform´atica Este cap´ıtulo contiene todos aquellos t´erminos y conceptos relacionados con la biolog´ıa y la bioinform´atica que se han considerado fundamentales para la comprensi´on completa del trabajo. Tanto la estructura como el contenido han sido modificados y ampliados tomando como base el Capitulo 2 y el Ap´endice B incluidos en mi proyecto fin de carrera [3]. 2.1 Fundamentos biol´ogicos Las mitocondrias son unos org´anulos contenidos en algunas c´elulas en los que se produce la oxidaci´on de las mol´eculas de glucosa, proceso por el cual la c´elula obtiene energ´ıa. El ADN mitocondrial humano (ADNmt de ahora en adelante) es un elemento biol´ogico muy interesante desde muchos puntos de vista. Al estar dentro de las mitocondrias y ser independiente del ADN celular, es centro de muchas teor´ıas y controversias respecto al origen o inclusi´on de estos corp´usculos en la c´elula, que es esencial para la vida de muchos seres vivos actualmente. El ADNmt posee tanto una tasa de mutaci´on como de conservaci´on muy elevadas, lo que lo hace id´oneo para el estudio de individuos y su relaci´on evolutiva, como miembros de una misma especie [20, 34, 43]. Adem´as, el ADNmt est´a muy ligado a ciertas enfermedades que provocan una muerte prematura en los individuos que las padecen. El ADNmt ha sido largamente estudiado, y actualmente se sabe que lo forman entre 16557 y 16576 nucle´otidos, divididos en 37 genes distintos. Un gen es un segmento de una secuencia de ADN o ARN que contiene toda la informaci´on necesaria para codificar un elemento funcional de la c´elula. De esos 37 genes, 13 tienen como objetivo la creaci´on de prote´ınas, 22 codifican ARNt (ARN transferente o de transferencia) y los dos restantes corresponden a las dos unidades que conforman el ribosoma (ARNr). Hay que indicar que el ADNmt es circular, a diferencia del ADN celular (muy conocido por su estructura de doble h´elice). Posee una regi´on de control, tambi´en conocida como bucle D, que no tiene como objetivo la codificaci´on sino la uni´on de los extremos para conformar su estructura circular. Esta zona contiene dos regiones hipervariables: HVR1y HVR2, las cuales pueden 5 CAP´ ITULO 2. BIOLOG´ IA Y BIOINFORM´ ATICA determinar el haplogrupo de un organismo [37]. Un haplogrupo es una agrupaci´on de haplotipos, los cuales representan conjuntos de polimorfismos con alguna caracter´ıstica estad´ıstica com´un. La clasificaci´on y determinaci´on de haplogrupos tiene como fundamento el detectar estas caracter´ısticas comunes, las cuales pertenecer´an de forma ´unica a la familia en cuesti´on. Los haplogrupos son muy importantes en el estudio de la relaci´on filogen´etica entre individuos, pudiendo determinar el punto de origen del linaje que se est´e estudiando. En el caso del ser humano, el ancestro matrilineal com´un m´as reciente del ser humano se ha denominado “Eva mitocondrial” y vivi´o hace unos 140000 a˜nos. El ancestro com´un m´as reciente global, incluyendo la herencia patrilineal, es mucho m´as reciente. Es necesario destacar la gran variabilidad del ADNmt, dado que es fundamental para poder realizar una clasificaci´on en haplogrupos de forma fiable. Para dar cuenta de dicha fiabilidad, por ejemplo, se estima que entre dos seres humanos cualesquiera existen ´unicamente entre 50 y 70 nucle´otidos diferentes, es decir, un 0.42 % del total de nucle´otidos del ADNmt, una proporci´on muy elevada. 2.2 Introducci´on a la bioinform´atica La bioinform´atica se puede considerar una rama de reciente aparici´on en la inform´atica cuyo objetivo es el uso de la inform´atica para la resoluci´on de problemas biol´ogicos. Su necesidad se ha visto acentuada en los ´ultimos a˜nos debido a la magnitud y complejidad que est´an adquiriendo ciertas investigaciones referentes a la biolog´ıa, que hacen intratable su resoluci´on de forma manual [28, 24]. En concreto en este trabajo se han tratado temas del ´area de an´alisis de secuencias, m´as en concreto, sobre construcci´on de filogenias. Para poder desarrollar el estudio de filogenias es indispensable disponer de secuencias previamente alineadas y verificadas. El alineamiento de secuencias es una operaci´on que establece equivalencias entre los distintos caracteres de las secuencias introducidas. Haciendo la suposici´on de que entre dichas secuencias existe un parentesco, pretende detectar aquellos hechos que las separan. Por otro lado, es importante que la informaci´on contenida en las secuencias haya sido verificada tras el proceso de secuenciaci´on, dado que los errores introducidos en algunas posiciones podr´ıan alterar la estructura final de la filogenia. 2.2.1. Filogen´etica La filogen´etica la ciencia que tiene como objetivo la clasificaci´on de las especies, tanto las existentes como las extinguidas, en base no a caracter´ısticas fenot´ıpicas sino a su relaci´on evolutiva. Si el lector esta familiarizado con estos temas puede que vea en esta definici´on muchas coincidencias con 6 2.2. INTRODUCCI´ ON A LA BIOINFORM´ ATICA la clad´ıstica, y no se equivoca: a menudo filogen´etica y clad´ıstica son tratadas como sin´onimos por tener el mismo objetivo. La clad´ıstica tiende a definir distintos tipos de grupos de especies. Los m´as conocidos son: Grupo monofil´etico o clado: compuesto por un organismo y todos sus descendientes. Al igual que pasaba con “Eva mitocondrial”, la ra´ız u origen de este grupo es el ancestro com´un m´as reciente del mismo. Grupo parafil´etico: compuesto por aquellos organismos cuya ra´ız u origen no incluye a todos los descendientes. Grupo polifil´etico: compuesto por varios grupos monofil´eticos no solapados. Dada su complejidad suele desaconsejarse su uso. Estos grupos normalmente se representan en ´arboles biol´ogicos donde se expresa la relaci´on de parentesco entre organismos o especies (almacenados en sus hojas). Estos ´arboles se conocen como ´arboles filogen´eticos. Aunque es com´un, sobre todo en la clad´ıstica, que estos ´arboles sean binarios, la aridad o n´umero de sub´arboles que se derivan de cada nodo interno puede ser mayor que 2. Modelos evolutivos Los modelos evolutivos, tambi´en conocidos como modelos de substituci´on, son modelos matem´aticos que se crearon para definir la probabilidad de cambio entre los distintos nucle´otidos en secuencias de ADN o ARN, o amino´acidos, en el caso de prote´ınas. Su objetivo es intentar explicar la evoluci´on que han sufrido dichas secuencias a lo largo del tiempo. El modelo puede contemplar la posibilidad de que exista reversi´on, es decir, que una mutaci´on que se haya producido se deshaga con el paso del tiempo. Un modelo queda completamente definido con la asignaci´on, parcial o total, de una serie de par´ametros. Los t´erminos que quedan libres ser´an determinados en el momento que se eval´ue el modelo, dependiendo su valor de las secuencias de entrada. Por ejemplo, algunos modelos establecen que el cambio entre nucle´otidos en secuencias de ADN o ARN es equiprobable, mientras que otros postulan que las transiciones (cambio entre nucle´otidos del mismo tipo) son m´as frecuentes que las transversiones (cambio entre nucle´otidos de distinto tipo). Recordar que con tipos se hace referencia a la divisi´on entre purinas, que son los nucle´otidos A y G, y las pirimidinas, C y T. No se debe olvidar que los modelos requieren que se asigne o calcule las frecuencias iniciales de los distintos nucle´otidos o amino´acidos. Selecci´on de modelos Actualmente existe una extensa variedad de modelos, y su elecci´on, lejos de ser irrelevante, puede acarrear serios problemas en trabajos de an´alisis o 7 CAP´ ITULO 2. BIOLOG´ IA Y BIOINFORM´ ATICA estudio de secuencias o prote´ınas. Esencialmente, una mala elecci´on en un modelo evolutivo puede establecer relaciones incorrectas entre espec´ımenes o, en el mejor de los casos, derivar en un ´arbol de peor calidad [27, 32, 42]. La complejidad y desconocimiento que se posee actualmente de la evoluci´on tiene como consecuencia inmediata que ni el mejor de los ´arboles sea fiel a la realidad [40]. Por tanto, nunca se debe olvidar que la filogen´etica pretende obtener filogenias pr´oximas a la realidad, pero en ning´un momento se puede estar seguro de que la coincidencia sea total. Aunque la importancia de la elecci´on de un modelo evolutivo levanta opiniones contradictorias entre los investigadores, debido a los inconvenientes expuestos, se intenta demostrar la efectividad de estos modelos comparando los resultados que se obtienen con conjuntos de pruebas creados de forma artificial, de los cuales se conoce su relaci´on evolutiva “real” y, por tanto, su modelo evolutivo correspondiente. Existen varios m´etodos aplicables a un ´arbol filogen´etico para conseguir un valor que permita determinar qu´e modelos se aproximan m´as a la realidad intr´ınseca de los datos. Un m´etodo muy aplicado es el de m´axima verosimilitud (Maximum Likelihood en ingl´es) [15]. Este m´etodo tiene como objetivo deducir los datos de entrada a partir de los resultados observables [40, 27, 41]. Asignando una probabilidad a las distintas mutaciones que ha podido sufrir una secuencia, se estima la distribuci´on de probabilidad del espacio de ´arboles. Es uno de los m´etodos m´as flexibles pero, a su vez, resulta uno de los m´as costosos, computacionalmente hablando. Se sospecha que puede tratarse de un problema NP-completo; pero as´ı como para otros m´etodo ha podido demostrarse, para el m´etodo de m´axima verosimilitud sigue siendo ´unicamente una suposici´on. Adem´as, pese a presentar una aproximaci´on m´as sencilla y pontente que otros m´etodos, la m´axima verosimilitud tiene el inconveniente de valorar m´as positivamente a los modelos m´as complejos (aquellos con un mayor n´umero de par´ametros libres). Para compensar esta deficiencia se han desarrollado los criterios de informaci´on, cuya base se sustenta en penalizar aquellos modelos de mayor complejidad. As´ı, mediante una combinaci´on de ambos par´ametros (complejidad y verosimilitud), se puede escoger un modelo evolutivo equitativo, dejando que sea el investigador el que decida en ´ultimo lugar. Como parte de los criterios de informaci´on se pueden encontrar el de Akaike (conocido como AIC) y el bayesiano (conocido como BIC). 2.2.2. Super´arboles Se obtiene como resultado de aunar varios ´arboles filogen´eticos bajo una misma ra´ız. Normalmente se requiere que, para que el resultado tenga sentido, exista cierto solapamiento entre las hojas que componen los distintos ´arboles filogen´eticos [38]. 8 3 Flujos de trabajo con selecci´on de modelos: una aproximaci´on multilocus al an´alisis filogen´etico En este cap´ıtulo se ha reflejado todo el trabajo de investigaci´on desempe˜nado con el objetivo de desarrollar un sistema para el an´alisis filogen´etico completo con conjuntos de datos muy grandes (miles y decenas de miles de secuencias biol´ogicas), tambi´en conocidas como filogenias extensivas. Como se ha indicado en el Cap´ıtulo 1, este trabajo ha partido de mi proyecto fin de carrera [3]. 3.1 Introducci´on A pesar de que la informaci´on de caracter biol´ogico que se tiene actualmente es imperfecta, una gran cantidad de conocimientos cient´ıficos se ha ido acumulando a lo largo de las ´ultimas d´ecadas. As´ı, los avances t´ecnicos en los computadores expanden continuamente las fronteras de lo que es viable computacionalmente. La inferencia de filogenias pertenece a una clase de problemas que se resisten a los enfoques de resoluci´on convencionales debido a su naturaleza altamente combinatoria. Por otra parte, el uso de modelos evolutivos, que reflejan patrones de cambio espec´ıficos observados en distintos conjuntos de datos, es crucial para la obtenci´on de filogenias lo m´as fieles posibles a la realidad. Sin embargo, este objetivo no ha sido abordado en aquellos casos en que los conjuntos de datos son grandes, debido a su elevado coste computacional asociado. La paralelizaci´on es una t´ecnica muy ´util para reducir el tiempo de ejecuci´on de los algoritmos (siempre que sea posible aplicarla), apoy´andose en la gran cantidad de recursos computacionales que est´an disponibles, por ejemplo, en forma de procesadores independientes interconectados por grandes redes de comunicaci´on. Los esfuerzos en esta direcci´on se han centrado, sobre todo, en la paralelizaci´on denominada “de grano fino” aplicada a los algoritmos est´andar y el uso de t´ecnicas algor´ıtmicas para mejorar las implementaciones existentes [26, 39]. Algunas de estas empresas han obtenido 9 CAP´ ITULO 3. FLUJOS DE TRABAJO CON SELECCI´ ON DE MODELOS: UNA APROXIMACI´ ON MULTILOCUS. . . medidas excepcionalmente buenas, pese a que los algoritmos est´andar no hayan sido dise˜nados teniendo en cuenta la concurrencia. Adem´as, el n´umero de tareas independientes que puede haber en un determinado momento es limitado, as´ı como su carga individual, la cual podr´ıa poner trabas a esquemas de asignaci´on simples o restringir su uso a redes de corta distancia o con baja latencia. Los esquemas tipo maestro-esclavo dominan esta clase de enfoques. Por otra parte, los flujos de trabajo (workflows) son formalismos abstractos que describen tareas complejas formadas por subtareas relacionadas entre s´ı. Esta sistematizaci´on ha sido desarrollada de forma natural en entornos de empresa y manufactura. Cuando son utilizados en la experimentaci´on cient´ıfica, no se encuentran ´unicamente en la estructuraci´on y documentaci´on de experimentos, sino que adem´as, si se dan las especificaciones e implementaciones adecuadas, forman la base de entornos de experimentaci´on de ejecuci´on autom´atica [18]. Una gran cantidad de proyectos software se han desarrollado con esta finalidad, incluyendo muchos enfocados a las aplicaciones en bioinform´atica [25]. Desgraciadamente, los flujos de trabajo no est´an concebidos para expresar y manejar los posibles niveles de paralelizaci´on alcanzables, aunque, con esfuerzo y algunos cambios, se puede lograr una adaptaci´on aproximada. Actualmente su elevada naturaleza interactiva es su mayor debilidad para alcanzar estos prop´ositos. El trabajo previo en el uso de flujos de trabajo en filogen´etica ha seguido esta misma filosof´ıa [11]. Aunque los entornos de ejecuci´on de bajo nivel y los entornos interactivos tengan distintas metas, ambos pueden beneficiarse consiguiendo costes computacionales inferiores. Para logarlo, bajo la propuesta de uso de flujos de trabajo que se va a exponer a lo largo de este cap´ıtulo, se deben incorporar los primeros dentro de los segundos, de forma anidada. Se defiende con la propuesta el uso eficiente de fuentes conocidas de tareas independientes y potencialmente concurrentes, y de informaci´on biol´ogica conocida o inferida que simplifica el caso general en algoritmos no informados y ofrece soluciones de mayor calidad en un menor tiempo. Con estas ideas en mente, el resto del cap´ıtulo aborda dos objetivos principales: en primer lugar, revelar los procesos concurrentes existentes a alto nivel en las aproximaciones tradicionales a los problemas de reconstrucci´on filogen´etica, incluyendo m´etodos de particionado que incrementan la granularidad y mejoran el rendimiento, que, adicionalmente, permiten incorporar criterios adicionales que se manejan en biolog´ıa (por ejemplo, diferentes genes y loci evolucionan de forma diferente seg´un muestran distintos estudios multilocus); y en segundo lugar, dise˜nar e implementar flujos de trabajo autom´aticos y arbitrariamente escalables, haciendo hincapi´e en lo relativo a la modularidad, la facilidad en el mantenimiento y la posibilidad de integrar nuevas etapas (como se ha hecho con la selecci´on de modelos evolutivos). 10 3.2. DESCOMPOSICI´ ON DEL PROBLEMA 3.2 Descomposici´on del problema A continuaci´on se va a centrar el estudio en el problema de la reconstrucci´on can´onica en la filogen´etica computacional. Dado un conjunto de secuencias Salineadas, el objetivo es generar un ´arbol Tque satisface (o se aproxima a) un cierto criterio de optimizaci´on. La conexi´on entre SyTse establece mediante la asignaci´on de etiquetas a las hojas del ´arbol por las secuencias del conjunto inicial. Las dimensiones principales del problema, que gobiernan la complejidad del algoritmo, son el n´umero de secuencias s, compartido por SyT, y la longitud de las secuencias l. Dicha longitud ldesigna el n´umero total de elementos clad´ısticos del conjunto, y, equivalentemente, la longitud del alineamiento. El problema de inferencia filogen´etica puede descomponerse en un conjunto de subproblemas independientes (como la selecci´on de modelos o la robustez estad´ıstica,) para los que es posible encontrar en la literatura cient´ıfica un abanico amplio de m´etodos algor´ıtmicos. Esto hace posible incorporar mejoras parciales e incluirlas directamente en el flujo de trabajo de alto nivel. El inter´es de esta invesigaci´on se centra en el m´etodo de evaluaci´on estad´ıstica [21]. Esta inclusi´on viene impuesta por el deseo de poder evaluar ciertas caracter´ısticas como la robustez y calidad de los resultados. La mayor´ıa de estos m´etodos requieren que se establezca un n´umero de r´eplicas estad´ısticas rque se pasar´a a un generador de muestras (probablemente junto a otros par´ametros), el cual generar´a un n´umero igual de alineamientos derivados, cada uno de los cuales constituir´a un problema independiente de la misma magnitud al original, que ser´a resuelto independientemente y agrupado posteriormente con el resto para poder obtener el resultado final. En este punto se encuentra el primer nivel de concurrencia originado por la independencia de los datos. Adem´as, se puede identificar en la selecci´on de modelos [42] una tarea que la precede, y que adem´as es apta para un tratamiento determinista. Mientras los par´ametros a usar con el algoritmo de construcci´on del ´arbol son proporcionados por el usuario, es conveniente, de manera general, utilizar procesos de selecci´on que eval´uen un amplio rango de modelos My que seleccionen aqu´el que se ajuste mejor al alineamiento de los datos. De nuevo, encontramos un flujo de trabajo muy sencillo compuesto por tantas tareas como modelos se tengan en consideraci´on m; siendo cada uno independiente del resto. Al finalizar la evaluaci´on de todos los modelos, se recopilan los resultados y se escoge aquel que haya obtenido mejor valoraci´on. Esta tarea tiene lugar despu´es del alineamiento y antes de que ninguna de los procesos asignados a cada muestra estad´ıstica comience su ejecuci´on. Esto en lo referente a concurrencia por independencia de los problemas. Sin embargo, se puede encontrar dentro de los propios datos una fuente de tareas independientes realizable por medios automatizados. Para empezar, los datos biol´ogicos presentan estructura, y, de hecho, empiezan a aparecer 11 CAP´ ITULO 3. FLUJOS DE TRABAJO CON SELECCI´ ON DE MODELOS: UNA APROXIMACI´ ON MULTILOCUS. . . estudios multilocus cada vez m´as completos que tienen como objetivo explotar los beneficios derivados del uso de informaci´on extra´ıda a priori [30, 29]. Esto quiere decir que es beneficioso hacer uso de la informaci´on que se puede extraer de los datos a priori. En este punto de vista, la preclasificaci´on de acuerdo a hechos establecidos o hip´otesis (a pesar de la prueba de las mismas por cualquier medio necesario) permite establecer particionados del conjunto de datos con dos principales beneficios: la generaci´on de tareas independientes y el reducido tama˜no de las mismas. A continuaci´on se van a examinar los efectos de esta propuesta en las dos dimensiones fundamentales. 1. Considerese la naturaleza de los distintos caracteres clad´ısticos, representados por l. A pesar de la naturaleza homog´enea de los alineamientos de las secuencias, el genoma est´a compuesto en realidad por zonas codificantes (genes) y no codificantes, a menudo afectadas de distinta forma e intensidad por la presi´on evolutiva. Por tanto, un alineamiento deber´ıa estar dividido en subconjuntos de columnas correspondientes a cada unidad gen´etica, habiendo gen total. As´ı pues, basta con que el alineamiento tenga una secuencia con los gumbrales definidos para poder realizar la divisi´on de forma adecuada. Cada subalineamien- to generado puede ser procesado como se ha descrito anteriormente (creando las muestras estad´ısticas y aplicando la selecci´on de modelos), para combinar en una ´ultima etapa los resultados obtenidos con los de sus compa˜neros mediante alg´un m´etodo de coalescencia [13]. 2. Mientras que las secuencias representan organismos individuales (sen t´erminos de la complejidad del problema) y es tarea del an´alisis filogen´etico determinar su relaci´on, ser´ıa posible identificar grupos de secuencias relacionados (evolutivamente hablando) de forma anticipada (haplogrupos en nuestro caso), y manejar estos grupos como unidades indivisibles de cara a la construcci´on del ´arbol. Un esquema de clasificaci´on jer´arquica ser´ıa apropiado para cumplir este objetivo, clasificando la pertenencia de cada secuencia a uno de los hgrupos de sub´arboles, los cuales pueden ser construidos de forma independiente y congregados en la ´ultima etapa mediante alg´un algoritmo de generaci´on de super´arboles [7]. Se ha mostrado anteriormente c´omo esta t´ecnica ofrece grandes mejoras cuando se manejan grandes conjuntos de datos [9]. Obviamente, ambas estrategias pueden ser combinadas para conseguir un mayor efecto. N´otese que el problema del alineamiento m´ultiple, el cual precede a todos los dem´as, puede beneficiarse de la aplicaci´on de los mismos principios expuestos en los dos ´ultimos p´arrafos. De hecho, ambas fuentes de informaci´on de la partici´on –una secuencia anotada para ly un clasificador de secuencias para s— son aplicadas habitualmente a las secuencias de forma individual, alineadas por pares con la secuencia de referencia. 12 3.3. DISE˜ NO DEL FLUJO DE TRABAJO 3.3 Dise˜no del flujo de trabajo El n´umero y variedad de tareas independientes implicadas es realmente grande y ofrece muchas posibilidades de cara a una ejecuci´on concurrente. La naturaleza de estas divisiones es a la vez simple y homog´enea, formada por pasos de generaci´on, de ejecuci´on de m´ultiples instancias y de combinaci´on, en los que la zona intermedia comprende la mayor´ıa de la carga computacional. Los diferentes tipos de tareas est´an dispuestos de forma anidada: las operaciones de reducci´on generan grupos de datos relacionados, aunque independientes, hasta el momento en que los problemas m´as simples son resueltos y el proceso de clasificaci´on se revierte mediante la combinaci´on adecuada de algoritmos de agrupaci´on (ver Figura 3.1). Partiendo de estas consideraciones, se propone un flujo de trabajo modular basado en la definici´on de “cajas negras” reutilizables para cada una de las capas significativas de trabajo concurrente. La construcci´on presentada tiene un formato com´un, aunque ´unico, y las variaciones que sea necesario incluir se pueden aplicar de forma sencilla. Cada una de las capas puede ser substituida por clasificadores triviales si se desea. El problema b´asico de la computaci´on de filogenias est´a formado por una etapa central y tres capas paralelas, como se explica a continuaci´on. Algoritmos de ´arboles. La etapa b´asica de la computaci´on puede tener, de hecho, el mismo estilo que un flujo de trabajo, debido a la combinaci´on de varios programas (los m´etodos de distancias son ejemplos muy comunes) o a las implementaciones de paralelismo a bajo nivel. Sea cual sea el caso, se puede suponer que esta “caja” toma el alineamiento Ay un conjunto de par´ametros y genera el ´arbol T. Muestreo estad´ıstico. La capa paralela fundamental comprende el muestreo estad´ıstico y la soluci´on concurrente a un n´umero de problemas base, seguido de la aplicaci´on de un algoritmo de consenso. La interface es semejante a la de la caja b´asica, exceptuando que es necesario indicar el n´umero de r´eplicas rque determinan la magnitud total de la carga de trabajo asociada. ´ Arboles de genes. El preprocesamiento en lopera sobre alineamientos de secuencias sencillos y genera una serie de subalineamientos para ser posteriormente muestreados de acuerdo a las instrucciones de divisi´on representadas por S, las cuales son requeridas junto con los par´ametros de sus tareas subordinadas. Es conveniente darse cuenta de que se puede suministrar un conjunto de modelos M, en vez de un ´unico modelo µ, para elegir los par´ametros m´as adecuados para cada gen, tal y como se explicar´a a continuaci´on. Como en las anteriores cajas, su prop´osito es generar un ´unico ´arbol Tque d´e explicaci´on al alineamiento A. 13 CAP´ ITULO 4. DETECCI´ ON AUTOM´ ATICA DE ERRORES DE SECUENCIACI´ ON A PARTIR DE INFORMACI´ ON FILOGEN´ ETICA las bases de datos de secuencias p´ublicas ha sido exponencial en los ´ultimos 30 a˜nos, teniendo como referencia la duplicaci´on actual del n´umero de registros en GenBank cada 35 meses aproximadamente (encontrando alg´un caso de hasta 18 meses) [6]. En cualquier caso, no se espera una disminuci´on en la velocidad de crecimiento del n´umero de animales y genoma humano. Todo esto permite recoger conjuntos de datos gen´eticos razonablemente completos y representativos para grandes grupos de especies, teniendo disponible una cantidad elevada de informaci´on gen´etica humana. De todas formas, la situaci´on no es tan negativa como parece en un principio, dado que la abundancia de material gen´etico podr´ıa, de hecho, convertirse en una ventaja. La compatibilidad parcial en estudios multilocus se pueden consolidar a trav´es de super´arboles y m´etodos de coalescencia [7, 13], y se pueden formar filogenias robustas, al menos considerando un nivel medio-alto, incluso aunque necesiten ser sometidas a procesos de revaluaci´on y refinamiento. Por lo tanto, es a los elementos que conforman las filogenias, es decir, las secuencias, a las que debemos dirigir nuestra atenci´on. Muchos de los contenidos de las bases de datos p´ublicas no pueden someterse a curaciones independientes, en parte, debido al crecimiento que est´an sufriendo. Esto se traduce en que la mayor´ıa de sus metadatos no est´an ni estandarizados ni son homog´eneos. Por este motivo se puede afirmar, con cierta confianza, que la calidad individual de las secuencias de ADN, es decir, la precisi´on de la copia respecto a la original, es desconocida a priori. Los errores de secuenciaci´on pueden ocurrir debido a la contaminaci´on de la muestra –como fue el caso en una de las secuencias de referencia m´as conocidas, la del ADNmt [4]– como tambi´en a causa de las t´ecnicas modernas de secuenciaci´on de alta productividad, las cuales replican segmentos muy peque˜nos del ADN (al contrario de lo que hac´ıan otros algoritmos m´as antiguos y menos eficientes) y despu´es reconstruyen la secuencia mediante alineamientos locales de estos fragmentos, en los que los trozos m´as cortos son m´as susceptibles a generar falsos positivos. Los ratios de error de las tecnolog´ıas actuales son desconocidos, y la contaminaci´on es claramente un factor no medible de forma aislada. Por el contrario, la disponibilidad de un gran n´umero de secuencias muy estructuradas y mayormente correctas permite detectar y descartar los errores m´as claros, as´ı como revelar caracter´ısticas inusuales para ser confirmadas como novedades representativas o descartadas como elementos espurios. El objetivo planteando en este cap´ıtulo es sistematizar la evaluaci´on de nuevas secuencias de forma que ´esta sea a la vez significativa e informativa. Para ello se van a utilizar colecciones de secuencias homog´eneas como base para poder determinar la calidad de las nuevas secuencias, tomando como indicador el ajuste que ´estas tengan a la estructura actual del conjunto. Este ajuste proporcionar´a tambi´en informaci´on suficiente como para poder se˜nalar los errores potenciales para que sean revisados. Inicialmente, estos conjuntos pueden representar cualquier cosa, desde genes individuales hasta 20 4.1. INTRODUCCI´ ON genomas completos a trav´es de un rango variable de organismos. Existen algunas comprobaciones b´asicas que pueden ser corroboradas directamente en los alineamientos de las secuencias, aunque ser´ıa muy beneficioso partir de una filogenia que ofrezca una clasificaci´on natural y una jerarqu´ıa de las variaciones gen´eticas conocidas, por lo que se puede conseguir una verificaci´on mucho m´as completa dentro de un contexto evolutivo. En esencia, una filogenia provee un medio para agrupar cada nueva secuencia con sus parientes putativos y permite medir la similitud entre ´esta y sus ancestros. Las mutaciones representativas de cada grupo tienen un alto ´ındice de conservaci´on, por lo que las excepciones a dichas mutaciones son ocasionadas en la mayor´ıa de los casos por errores en el proceso de secuenciaci´on (una vez se ha alcanzado una cierta profundidad en el ´arbol base). Es posible descubrir nuevos subgrupos y variaciones inusuales, y deber´ıan ser tratados como excepciones leg´ıtimas, pero tan pronto como sean integradas en la filogenia (mejorando adem´as la aproximaci´on del ´arbol evolutivo) los descubrimientos extremadamente divergentes se ir´an volviendo cada vez m´as raros, hasta que desaparezcan. La interpretaci´on de la alarma depender´a del estado actual de la filogenia gu´ıa y de las secuencias en consideraci´on. Un hecho interesante es que el mismo proceso seguido para determinar lo buena que es una secuencia est´a conectado al proceso de construcci´on de filogenias, pues se analiza dicha secuencia seg´un el lugar m´as adecuado en el que se sit´ua dentro de la filogenia. Esta es una caracter´ıstica adicional muy interesante, pues se ha mostrado que aunque las actualizaciones peri´odicas de grandes filogenias puede ser viable a trav´es de m´etodos est´andar [10], todav´ıa resulta un proceso muy costoso cuya frecuencia para actualizar la filogenia debe ser limitada, lo que se deriva en una degradaci´on y empobrecimiento de la misma en ese lapso de tiempo. Por otro lado, un ´arbol bien fundamentado y con una estructura principal estable a largo plazo, facilita el proceso de localizaci´on de la nueva secuencia, adem´as de que su incorporaci´on al ´arbol, si no exacta, afectar´ıa, como mucho, de forma localizada. Sin duda, es dif´ıcil lograr un m´aximo beneficio en el compromiso entre estas dos perspectivas. En resumen, en el resto del cap´ıtulo se va a presentar un algoritmo inspirado en la filogen´etica que tiene como objetivo evaluar la precisi´on de las secuencias biol´ogicas a trav´es de la historia evolutiva incorporada en el ´arbol de referencia. Estableciendo el valor de los umbrales de forma adecuada, el algoritmo ser´a capaz de determinar cu´ando una secuencia es correcta (sus mutaciones caracter´ısticas no entran en contradicci´on con la estructura evolutiva manifestada en la filogenia), cu´ando puede ser potencialmente correcta aunque inusual y cu´ando ´esta se desv´ıa demasiado de la filogenia, por lo que se requerir´a de una revisi´on cuidadosa. 21 CAP´ ITULO 4. DETECCI´ ON AUTOM´ ATICA DE ERRORES DE SECUENCIACI´ ON A PARTIR DE INFORMACI´ ON FILOGEN´ ETICA 4.2 Detecci´on de errores de secuenciaci´on Para poder determinar si una mutaci´on es real o puede haber sido generada en el proceso de secuenciaci´on, el algoritmo requiere alg´un elemento con el que poder comparar la secuencia. Un ´arbol filogen´etico contiene una gran cantidad de informaci´on referente a la evoluci´on a lo largo del tiempo de un tipo de estructura biol´ogica determinada (ADN nuclear, ADNmt, proteinas, . . . ), por lo que es un candidato id´oneo con el que realizar la comparaci´on de la nueva secuencia. Como es obvio, es necesario que el ´arbol seleccionado haya sido construido bajo un modelo evolutivo conocido y aceptado, y con secuencias verificadas minuciosamente, de forma que no se hayan introducido errores en la fase de construcci´on de la filogenia. Esta ´ultima hip´otesis es muy importante debido a que, si una rama queda definida por ciertas mutaciones err´oneas, se podr´ıan validar secuencias con errores en esas mismas posiciones. Como se ha mencionado anteriormente, el proceso principal del algoritmo presentado se basa en la localizaci´on del lugar en el ´arbol filogen´etico donde la secuencia de entrada se ajusta mejor. Antes de exponer su funcionamiento de forma m´as detallada, se van a explicar dos operaciones b´asicas que ser´an necesarias a lo largo de la ejecuci´on del algoritmo. La primera es el Filtro Hamming. Esta operaci´on toma como entrada dos secuencias de la misma longitud y provee como salida su distancia Hamming, es decir, el n´umero total de posiciones en las que las dos secuencias no comparten el mismo valor (incluyendo en las posibilidades tanto los huecos como el car´acter N o “desconocido”). La segunda es el Filtro de Referencia. Antes de ejecutar el algoritmo, se debe seleccionar una secuencia de referencia (que tendr´a que estar incluida en el ´arbol filogen´etico). Como entrada, la operaci´on toma dos secuencias, de las cuales una ser´a la que se da como secuencia a analizar por el algoritmo y la otra pertenecer´a al ´arbol. Primero, la operaci´on obtiene una lista de las posiciones donde la secuencia de referencia difiere de la secuencia del ´arbol. Dado el hecho de que los huecos introducidos en la secuencia de referencia son debidos al alineamiento previo para la construcci´on de la filogenia, estas posiciones ser´an ignoradas. A continuaci´on, se comparan los valores de la secuencia del ´arbol y de la secuencia de entrada en las posiciones incluidas en la lista obtenida en el paso anterior. Como salida se obtiene el n´umero de diferencias resultantes de esta comparaci´on. El primer paso del algoritmo es alinear la secuencia de entrada con las secuencias del ´arbol. Esto significa que las secuencias del ´arbol filogen´etico deben estar previamente alineadas. Para este paso se ha utilizado MUSCLE [14], una herramienta para procesos de multialineamiento, a˜nadiendo la nueva secuencia al alineamiento del ´arbol. Despu´es, el algoritmo localiza el nodo m´as pr´oximo a la secuencia nueva. Para ello toma un nodo, al que denominaremos padre, y todos los nodos a 22 4.3. PRUEBAS Y RESULTADOS un nivel de distancia por debajo, sus nodos hijo. A continuaci´on aplica el Filtro de Referencia a los pares generados con la secuencia de entrada y cada una de las secuencias de los nodos seleccionados. Normalmente, uno de los nodos hijo es seleccionado como el m´as cercano de todos los pares evaluados, por lo que se establece este nodo como un nuevo nodo padre y se repite el proceso hasta que el algoritmo alcance a un nodo hoja. Como es obvio, el primer nodo seleccionado es la ra´ız del ´arbol. Existen otros casos que podr´ıan darse en vez del caso com´un presentado. En vez de un ´unico nodo, se pueden obtener dos o m´as como nodos m´as cercanos a la secuencia. Si todos los nodos son nodos hijo, el algoritmo explora cada una de las posibilidades independientemente, aplicando el proceso de forma individual. Las pruebas mostradas en la siguiente secci´on demuestran que esta situaci´on con m´ultiples caminos no se mantiene normalmente m´as all´a de dos o tres iteraciones. Si uno de los nodos cercanos es el nodo padre, ´este es descartado en la medida en que preferimos un resultado mas cercano a las hojas. Las pruebas han revelado algunas situaciones a las que se han denominado m´ınimos locales, donde el padre resulta ser el ´unico nodo m´as cercano, pero, como su nombre indica, se trata de una situaci´on local: existen otros nodos, m´as pr´oximos a las hojas, que son m´as cercanos a la secuencia que el nodo padre. Para poder eludir estos m´ınimos locales, el algoritmo aplica el Filtro Hamming a los mismos pares manejados en el Filtro de Referencia. Los diferentes resultados son procesados como en los casos expuestos anteriormente, exceptuando el caso en el que se obtenga de nuevo el padre como ´unico nodo cercano. En este caso, se ha alcanzado un m´ınimo global, por lo que ese nodo del ´arbol es el m´as cercano a la secuencia nueva. En un ´ultimo paso, el algoritmo aplica nuevamente el Filtro de Referencia al nodo seleccionado para obtener el n´umero total de diferencias con la nueva secuencia. Dos umbrales determinar´an si la secuencia es correcta (Right), si existen posibles errores (Alarm), o si es, con bastante seguridad, err´onea (Wrong). Es importante saber que estos umbrales no funcionar´an de forma adecuada si la secuencia de entrada corresponde a alguna especie desconocida para la filogenia, u otros casos similares, que se reflejan como “agujeros” en el ´arbol filogen´etico. Intuitivamente, debido a las situaciones de m´ultiples caminos, el algoritmo podr´ıa mostrar m´as de un nodo soluci´on. Mirando de forma detenida al ´arbol se ha observado que todas estas soluciones con parientes cercanos entre s´ı, es decir, nodos con el mismo nodo padre o sobrino. 4.3 Pruebas y resultados Como se ha destacado anteriormente, la detecci´on de variantes inusuales requiere de una filogenia relativamente estable y muy poblada. Por el momento hay pocas instancias donde la cobertura sea suficientemente elevada 23 CAP´ ITULO 4. DETECCI´ ON AUTOM´ ATICA DE ERRORES DE SECUENCIACI´ ON A PARTIR DE INFORMACI´ ON FILOGEN´ ETICA como para asegurar robustez. Posiblemente el m´as significativo de todos es el ADNmt humano, dado que es f´acil de secuenciar, no recombinante (por lo que puede ser utilizado en bloque en la reconstrucci´on filogen´etica) y altamente informativo [5]. La mayor´ıa de las ´areas de la filogenia mitocondrial humana est´an fielmente representadas y las mutaciones caracter´ısticas por las cuales grandes grupos de individuos son relacionados est´an organizadas en jerarqu´ıas extensivas de haplogrupos mitocondriales [44]. De hecho, el algoritmo toma la inspiraci´on de los procedimientos aplicados en la construcci´on incremental de la filogenia de MITOMAP [35]. Mientras se combinan estos procesos con reconstrucciones estrictas, es posible aplicar los resultados del algoritmo para construir copias de la filogenia mitocondrial que reflejen las ´ultimas incorporaciones al conjunto de secuencias de ADNmt humano publicadas. Esto es de gran utilidad debido a que la mayor parte de las secuencias pueden ser situadas con gran precisi´on en el ´arbol, lo que tiene el valor a˜nadido de tener filogenias permanentemente actualizadas y siempre disponibles. En referencia, las actualizaciones actuales a la filogenia mitocondrial humana ZARAMIT [8] requieren por encima de un a˜no de tiempo de CPU secuencial (tiempo que puede ser reducido a semanas mediante el uso apropiado de procesamiento paralelo). Es claramente inviable y muy ineficiente producir actualizaciones cada vez que unas pocas secuencias son publicadas, pues su valor incremental no justifica el gasto adicional de recursos computacionales. No obstante, el n´umero de adiciones entre reconstrucciones, por ejemplo cada pocos meses, puede ser extremadamente significativo. Por ejemplo, durante los primeros seis meses del 2011 aproximadamente 1000 secuencias nuevas fueron publicadas en GenBank, lo que implica un crecimiento del 14 % de la colecci´on de secuencias mitocondriales completas acumuladas principalmente a lo largo de la ´ultima d´ecada. Para los experimentos realizados se ha utilizado el ´ultimo ´arbol filogen´etico creado por el proyecto ZARAMIT, compuesto por 7390 secuencias de ADNmt humano obtenidas de GenBank. Como secuencia de referencia, se ha utilizado la secuencia de referencia revisada de Cambridge (rCRS). Los umbrales para el Filtro de Referencia se han establecido en 0 para discriminar entre los estados Right yAlarm, y 3 para la distinci´on entre Alarm yWrong. 4.3.1. Estudio de comportamiento Para poder estudiar el comportamiento del algoritmo, se han dividido los experimentos en tres grupos, cada uno centrado en obtener unos resultados espec´ıficos dentro de todos los casos posibles. 1. Localizaci´on correcta de las hojas: El primer experimento tiene como objetivo localizar correctamente algunas de las secuencias que for- 24 4.3. PRUEBAS Y RESULTADOS Tabla 4.1: Secuencias pertenecientes al ´arbol filogen´etico y su clasificaci´on por el algoritmo. Accession Ref. Clasificaci´on Localizaci´on Distancia DQ246811 [33] RIGHT DQ246811 0 DQ246826 [33] RIGHT DQ246826 0 DQ246828 [33] RIGHT DQ246828 0 DQ246830 [33] RIGHT Anc3521, Anc3534 271 AY738944 [1] RIGHT AY738944 0 AY738945 [1] RIGHT AY738945 0 AY738946 [1] RIGHT AY738946 0 AY738947 [1] RIGHT AY738947 0 AY738948 [1] RIGHT AY738948 0 AY738949 [1] RIGHT AY738949 0 AY738957 [1] RIGHT AY738957 0 AY738958 [1] ALARM(1) Anc4104, Anc3956 7 AY738980 [1] RIGHT AY738980 0 AY738981 [1] RIGHT Anc4051, . . . 2 AY738982 [1] RIGHT AY738982 0 AY738990 [1] RIGHT AY738990 0 AY738991 [1] RIGHT AY738991 0 AY738992 [1] RIGHT AY738992 0 AY738993 [1] RIGHT AY738993 0 AY738994 [1] RIGHT AY738994 0 man parte de las hojas del ´arbol filogen´etico. Espec´ıficamente, se han seleccionado 20 secuencias del conjunto de hojas del ´arbol. Los accession as´ı como los resultados obtenidos por el algoritmo se presentan en la Tabla 4.1. Como se puede ver en la Tabla 4.1, 17 secuencias se han localizado correctamente. A pesar de que la clasificaci´on de dos de las tres restantes ha sido Right, el algoritmo no ha sido capaz de localizarlas en el ´arbol. La primera, DQ246830, tiene una distancia enorme con los nodos m´as cercanos, hecho que no sucede con la segunda, AY738981. Si observamos la primera de estas secuencias, se puede ver que el primer fragmento de la regi´on de control no est´a, por lo que es normal haber 25 CAP´ ITULO 4. DETECCI´ ON AUTOM´ ATICA DE ERRORES DE SECUENCIACI´ ON A PARTIR DE INFORMACI´ ON FILOGEN´ ETICA Tabla 4.2: Secuencias de distintos animales y su clasificaci´on por el algoritmo. Accession Animal Clasificaci´on Distancia NC 001643 Chimpanc´e RIGHT 1966 NC 001941 Oveja RIGHT 4720 NC 007402 Serpiente de rayo de sol RIGHT 6324 NC 006160 Mosca blanca RIGHT 9303 NC 005313 At´un bala WRONG(7) 5779 NC 009885 Nematodo RIGHT 9077 NC 006281 Cangrejo azul RIGHT 9251 NC 005805 Rana de ´arbol moteada RIGHT 6191 NC 008159 Coral de setas WRONG(6) 8239 NC 009684 Pato silvestre WRONG(6) 5822 obtenido esa distancia y que el algoritmo no haya sido capaz de localizar la secuencia. En el segundo caso, el algoritmo no ha podido localizar la secuencia, pero en los resultados est´a incluido el nodo Anc4051 y otros nodos y hojas que son hijos del mismo nodo. El Anc4051 es el bisabuelo de AY738981, por lo que el resultado es muy pr´oximo a la soluci´on ideal. Esto puede suceder con algunas secuencias debido a la dependencia existente entre el algoritmo y el ´arbol. Algo parecido sucede con la secuencia AY738958, solo que en este caso una de las mutaciones no se ajusta exactamente al cl´uster donde ha sido situada, por lo que se ha clasificado como Alarm. 2. ADN mitocondrial no humano: En estos experimentos se han utilizado secuencias procedentes de otros animales, de forma que se pueda comprobar el comportamiento del algoritmo con secuencias que no se ajustan a un ´arbol filogen´etico de ADNmt. Los animales espec´ıficos y los accessions de las secuencias se muestran en la Tabla 4.2. El primer elemento que llama la atenci´on de los resultados es que la mayor´ıa de las secuencias se han clasificado como Right. Este hecho est´a relacionado con el prealineamiento realizado por el algoritmo con la secuencia rCRS, lo que puede alterar el valor de algunas posiciones para obtener el mejor resultado posible. Dada esta situaci´on, el resultado m´as importante de estos experimentos es el campo Distancia. No es sorprendente que la secuencia del chimpanc´e sea la m´as cercana de entre todas las secuencias examinadas. En la mayor´ıa de los casos, la distancia implica que m´as del 30 % de los nucle´otidos est´an mal, lo que es otra se˜nal inequ´ıvoca de que la secuencia no se ajusta al ´arbol 26 4.3. PRUEBAS Y RESULTADOS Tabla 4.3: Secuencias sint´eticas creadas a partir de la secuencia AY738958. Accession Mutaci´on SEQ00001 -3106A SEQ00002 -3106A, G8859A SEQ00003 -3106A, G8859A, G15325A SEQ00004 T6775C SEQ00005 T6775C, G1437A Tabla 4.4: Clasificaci´on de las secuencias sint´eticas por el algoritmo. Accession Clasificaci´on Localizaci´on Distancia AY738958 ALARM(1) Anc4104, Anc3956 7 SEQ00001 ALARM(2) Anc4104, Anc3956 8 SEQ00002 ALARM(3) Anc4104, Anc3956 9 SEQ00003 WRONG(4) Anc4104, Anc3956 10 SEQ00004 RIGHT Anc4076 7 SEQ00005 ALARM(1) Anc4076 8 y ser´ıa aconsejable comprobar si pertenece a un Homo sapiens. 3. Mutaciones sint´eticas: Con estos experimentos se pretende probar que el algoritmo es capaz de detectar cada mutaci´on relevante de forma individual. Se ha tomado la secuencia AY738958 como secuencia base a la que se han introducido de forma artificial algunas mutaciones para observar c´omo cambia el resultado del algoritmo. Estas mutaciones se muestran en la Tabla 4.3. Como suele ser el formato habitual en biolog´ıa, las mutaciones se codifican mostrando el valor anterior, la posici´on y el nuevo valor asignado a dicha posici´on. La Tabla 4.4 contiene los resultados de cada una de las secuencias analizadas, mostrando de nuevo el resultado de la secuencia AY738958 para poder ver c´omo las mutaciones introducidas han afectado a los resultados. Las primeras tres secuencias sint´eticas demuestran c´omo una sola mutaci´on, aplicada en el sitio adecuado, puede cambiar una clasificaci´on de Alarm aWrong. Las ´ultimas dos reflejan c´omo, obviamente, una mutaci´on puede cambiar el nodo m´as cercano. Normalmente, como en este caso, el resultado cambiar´a de un nodo a otro nodo hermano, por lo que no ser´a un cambio demasiado relevante. Pero si tres o cuatro 27 CAP´ ITULO 4. DETECCI´ ON AUTOM´ ATICA DE ERRORES DE SECUENCIACI´ ON A PARTIR DE INFORMACI´ ON FILOGEN´ ETICA mutaciones o errores se acumulan a lo largo de la secuencia en las posiciones adecuadas, se pueden obtener como nodos cercanos unos verdaderamente alejados de la localizaci´on real que tendr´ıa la secuencia en el ´arbol filogen´etico. 4.3.2. Estudio de rendimiento Todos los experimentos se han ejecutado en un ordenador con un procesador Core 2 Duo E6750 y 8 GB de RAM. El coste temporal que implica cargar el ´arbol filogen´etico y toda la informaci´on necesaria por la aplicaci´on supone, como m´aximo, 20 segundos, teniendo en cuenta que estos datos solo deber´an ser cargados la primera vez, cuando se inicie la aplicaci´on. El programa tarda 20 segundos de media en alinear y localizar la secuencia de entrada. Por tanto, el programa tiene un rendimiento excelente y el usuario puede obtener los resultados en “tiempo real”, generando adem´as una realimentaci´on hacia el usuario al mostrar las posiciones que se han marcado como malas (si las hay) con respecto a los nodos m´as cercanos. El peor caso de los experimentos realizados se ha dado al utilizar como entrada la secuencia de la serpiente de rayo de sol, donde al programa ha tardado 32 segundos alinearla y localizarla. 28 5 Complejidad param´etrica en bioinform´atica: filogenias casi perfectas En este cap´ıtulo se va a mostrar la investigaci´on actualmente en curso sobre complejidad param´etrica aplicada a la bioinform´atica. En concreto, se centra el caso de algoritmos te´oricos que tienen como objetivo construir filogenias perfectas y casi perfectas. 5.1 Introducci´on Los problemas clasificados como NP-duros son realmente dif´ıciles de resolver de forma eficiente y ´optima debido, principalmente, a su alta explosi´on combinatoria. Es com´un encontrar algoritmos con heur´ısticas o aproximaciones para este tipo de problemas, que obtienen la soluci´on con un coste computacional aceptable, aunque no aseguran hallar la soluci´on ´optima. Un subconjunto de estos problemas NP-duros tiene una caracter´ıstica especial: su complejidad depende de la estructura de la entrada. Este tipo de casos permiten una descomposici´on del problema en un conjunto de par´ametros, los cuales afectan de forma distinta a la complejidad del problema. El Tratamiento por Fijaci´on de Par´ametros, en ingl´es, Fixed-Parameter Tractablility (FPT), es una t´ecnica aplicada a este tipo de casos que permite estudiar subconjuntos del problema fijando aquellos par´ametros que afecten de forma m´as severa a su complejidad. De esta forma, se pueden llegar a desarrollar algoritmos eficientes para subproblemas, a´un interesantes, del problema inicial, que adem´as obtengan la soluci´on ´optima. Pese a ser una herramienta puramente te´orica, su uso est´a muy extendido. Como ejemplo, el FPT ha sido utilizado en problemas de teor´ıa de grafos, como el cubrimiento de v´ertices, entre otros [22]. Como ya se ha comentado en los cap´ıtulos anteriores, en bioinform´atica existen muchos problemas que a´un no se han podido resolver debido a su coste computacional. El enfoque FPT constituye, en estos casos, una herramienta muy ´util para poder obtener un algoritmo que resuelva aquellos subproblemas “sencillos” pero igualmente interesantes. Dentro de la filogen´etica, se ha utilizado tanto para problemas relacionados con su construcci´on, como 29 CAP´ ITULO 5. COMPLEJIDAD PARAM´ ETRICA EN BIOINFORM´ ATICA: FILOGENIAS CASI PERFECTAS 5.4.3. Lemas A continuaci´on se van a mostrar los lemas que han tenido que ser adaptados al nuevo criterio y funci´on de penalizaci´on. Lema 4: Sea Tuna filogenia Cp-perfecta tal que penalty (T)≤q. Entonces Ttiene, como m´aximo, q log(P 1−P)+mraristas malas. Demostraci´on: Sea C0=C−Cp. Para cada cC0, sea lcel n´umero de aristas(u, v) en Ttal que c(u)6=c(v). Por el Lema 3 [17], para cada arista mala (u, v) tiene que existir un cC0tal que c(u)6=c(v). Adem´as, 1 ≤C0< m, dado que Tes una filogenia casi perfecta. Por lo tanto, el n´umero de aristas malas es, como m´aximo, PcC0lc. Por otro lado, X cC0 log P(lc+1−rc)(1 −P)(rc−1−lc)≤q, log P 1−PX cC0 lc+|C0|(1 −r) log P 1−P≤q, X cC0 lc≤q log P 1−P+|C0|(r−1) ≤    q log P 1−P+mr    Por lo tanto, el n´umero de aristas malas est´a limitado por q log(P 1−P)+mr. Lema 5: Supongamos que Qes una subfamilia de caracteres la cual tiene una filogenia Cp-perfecta T, siendo xla ra´ız de Ty que cumple que ∀cCp,c(x) = c(Sv (Q)) si c(Sv (Q)) 6=∗. Sea αla asingaci´on de penalizaci´on de x. Entonces, (Q, α) tiene una subfilogenia con un valor de verosimilitud, como m´aximo, igual al de T. La demostraci´on sigue siendo v´alida para el enunciado presentado. Lema 6: Supongamos que el par (Q, α) tiene una subfilogenia. Entonces, (Q, α) tiene una subfilogenia de m´axima verosimilitud en forma normal. La demostraci´on sigue siendo v´alida para el enunciado presentado. 36 6 Conclusiones Se ha desarrollado un sistema donde se ha aplicado una metodolog´ıa divide y vencer´as para el dise˜no de flujos de trabajo (empleando el principio de transparencia de caja negra) que integra selecci´on de modelos y reconstrucci´on filogen´etica, y que puede lidiar con an´alisis de filogenias extensivas de forma eficiente. El sistema ha sido probado con un conjunto grande de datos de mtDNA, aunque acepta cualquier tipo de dato biol´ogico como entrada. Adem´as, el criterio para la divisi´on de los datos de entrada en subconjuntos puede ser personalizado para reflejar correctamente la naturaleza los mismos; por supuesto, el n´umero de muestras estad´ısticas tambi´en puede ser modificado. El sistema obtiene speedups mayores a 50 compar´andolo con su equivalente secuencial en estudios filogen´eticos grandes; se han obtenido grandes mejoras en la fase de selecci´on de modelos compar´andola de forma independiente con herramientas de objetivo espec´ıfico como jModelTest. Se ha presentado un nuevo algoritmo para evaluar los errores cometidos por el proceso de secuenciaci´on, proporcionando como salida el nivel de veracidad de la secuencia dada una filogenia. Actualmente, esta comprobaci´on de las secuencias se realiza de forma manual, lo que implica una gran inversi´on de tiempo, sin tener en cuenta los posibles errores humanos. La soluci´on propuesta ofrece un detector autom´atico de posibles errores cometidos en el proceso de secuenciaci´on, con un rendimiento muy bueno, que tambi´en muestra como resultado los nodos m´as cercanos de la filogenia a la nueva secuencia. Si ´esta se considera lo suficientemente buena, el algoritmo autom´aticamente la agrega a la filogenia, por la que la informaci´on est´a siempre actualizada. Finalmente, como trabajo futuro del sistema se buscar´a obtener mejoras en el speedup conseguido, que parece degradarse a medida que el tama˜no del problema crece. Tambi´en se contemplar´an formas de integrar, como tareas previas al sistema actual, la obtenci´on autom´atica de datos a partir de la entrada y el proceso de alineamiento de las secuencias. En este contexto, se espera mejorar los costes computacionales investigando nuevos criterios que permitan explotar m´as a´un la naturaleza inherentemente paralela de este tipo de procesos. En lo concerniente al algoritmo de verificaci´on de secuencias, como fu- 37 CAP´ ITULO 6. CONCLUSIONES turas mejoras se plantea desarrollar un nuevo proceso de verificaci´on de las mutaciones detectadas como posibles errores, a˜nadiendo un nuevo nivel de viabilidad biol´ogica. Esta verificaci´on consistir´a en tener en cuenta las reversiones, por lo que una mutaci´on que ya haya aparecido anteriormente en el ´arbol implicar´a que no se trata de una mala mutaci´on. Adem´as, el ratio de conservaci´on entre distintas especies proporcionar´a un criterio extra para evaluar la viabilidad biol´ogica de la mutaci´on. Por otra parte, el proceso de alineamiento siempre es una tarea que implica un gran consumo de tiempo, as´ı que cualquier mejora en esta fase mejorar´a el rendimiento global del algoritmo. Respecto al trabajo ya planteado sobre complejidad param´etrica, resta finalizar la demostraci´on de que la generalizaci´on al caso de MV propuesta tiene unas garant´ıas de prestaciones y eficiencias razonables. Al tratarse de un trabajo de ´ındole te´orico pueden surgir dificultades adicionales en estas demostraciones, pero a su vez se puede contar con una evaluaci´on experimental del algoritmo que se plantea. Por otro lado, se propone tambi´en extender el algoritmo a casos generales de m´axima verosimilitud y m´axima verosimilitud ancestral, adem´as de abordar pruebas de completitud para los problemas cercanos para los que no se encuentren soluciones param´etricas eficientes. 38 Bibliograf´ıa [1] Achilli, A., Rengo, C., Magri, C., Battaglia, V., Olivieri, A., Scozzari, R., Cruciani, F., Zeviani, M., Briem, E., Carelli, V., Moral, P., Dugoujon, J., Roostalu, U., Loogv¨ oli, E., Kivisild, T., Bandelt, H., Richards, M., Villems, R., Santachiara- Benerecetti, A., Semino, O., and Torroni, A. The molecular dissection of mtDNA haplogroup H confirms that the Franco-Cantabrian glacial refuge was a major source for the European gene pool. American Journal of Human Genetics 75 (Nov 2004), 910–918. [2] Agarwala, R., and Fern´ andez-Baca, D. A polynomial-time algorithm for the perfect phylogeny problem when the number of character states is fixed. SIAM Journal on Computing 23, 6 (1994), 1216–1224. [3] ´ Alvarez-Jarreta, J. An´alisis te´orico-pr´actico de m´etodos de inferencia filogen´etica basados en selecci´on de modelos y m´etodos de super´arboles. Master’s thesis, Centro Polit´ecnico Superior, Universidad de Zaragoza, 2010. [4] Andrews, R., Kubacka, I., Chinnery, P., Lightowlers, R., Turnbull, D., and Howell, N. Reanalysis and revision of the Cambridge reference sequence for human mitochondrial DNA. Nature Genetics 23 (Oct. 1999), 147. [5] Bandelt, H., Macaulay, V., and Richards, M., Eds. Human mitochondrial DNA and the evolution of Homo sapiens. Springer, Berlin, Germany, 2006. [6] Benson, D., Karsch-Mizrachi, I., Lipman, D., Ostell, J., and Sayers, E. GenBank. Nucleic Acids Research 38 (Jan. 2010), D46– D51. [7] Bininda-Emonds, O., Gittleman, J., and Steel, M. The (super)tree of life: procedures, problems and prospects. Annual Review of Ecology and Systematics 33 (Dec. 2002), 265–289. [8] Blanco, R., and Mayordomo, E. ZARAMIT: a system for the evolutionary study of human mitochondrial DNA. In IWANN 2009, Part II (2009), vol. 5518 of Lecture Notes in Computer Science, pp. 1139–1142. [9] Blanco, R., Mayordomo, E., Montes, E., Mayo, R., and Alberto, A. Advances in Bioinformatics. Springer, Heidelberg, 2010, ch. Scalable phylogenetics through input preprocessing, pp. 123–130. [10] Blanco, R., Mayordomo, E., Montoya, J., and Ruiz-Pesini, E. Rebooting the human mitochondrial phylogeny: an automated and scalable methodology with expert knowledge. BMC Bioinformatics 12 (May 2011), 174. 39 BIBLIOGRAF´ IA [11] Bowers, S., McPhillips, T., Riddle, S., Anand, M., and Lud¨ ascher, B. Provenance and Annotation of Data and Processes. Springer, Heidelberg, 2008, ch. Kepler/pPOD: scientific workflow and provenance support for assembling the tree of life, pp. 70–77. [12] Couvares, P., Kosar, T., Roy, A., Weber, J., and Wenger, K. Workflows for e-Science. Springer, 2006, ch. Workflow management in Condor, pp. 357–375. [13] Degnan, J., and Rosenberg, N. Gene tree discordance, phylogenetic inference and the multispecies coalescent. Trends in Ecology & Evolution 24 (Jun. 2009), 332–340. [14] Edgar, R. MUSCLE: multiple sequence alignment with high accuracy and high throughput. Nucleic Acids Research 32 (Mar. 2004), 1792– 1797. [15] Felsenstein, J. Inferring Phylogenies. Sinauer Associates, 2004. [16] Fern´ andez-Baca, D. The perfect phylogeny problem. In Steiner Trees in Industries (2000), D. Du and X. Cheng, Eds., Kluwer Academic Publishers. [17] Fern´ andez-Baca, D., and Lagergren, J. A polynomial-time algorithm for near-perfect phylogeny. SIAM Journal on Computing 32 (2003), 1115–1127. [18] Georgakopoulos, D., Hornick, M., and Sheth, A. An overview of workflow management: from process modeling to workflow automation infrastructure. Distributed and Parallel Databases 3, 2 (1995), 119–153. [19] Gramm, J., Nickelsen, A., and Tantau, T. Fixed-parameter algorithms in phylogenetics. The Computer Journal (2007). [20] Herrnstadt, C., Elson, J., Fahy, E., Preston, G., Turnbull, D., Anderson, C., Ghosh, S., Olefsky, J., Beal, M., Davis, R., and Howell, N. Reduced-median-network analysis of complete mitochondrial DNA coding-region sequences for the major African, Asian, and European haplogroups. American Journal of Human Genetics 70, 5 (2002), 1152–1171. [21] Holder, M., and Lewis, P. Phylogeny estimation: traditional and bayesian approaches. Nature Reviews Genetics 4 (2003), 275–284. [22] H¨ uffner, F., Niedermeier, R., and Wernicke, S. Techniques for practical fixed-parameter algorithms. The Computer Journal 51, 1 (2008). 40 BIBLIOGRAF´ IA [23] Kim, S., Tang, H., and Mardis, E., Eds. Genome sequencing technology and algorithms. Artech House, Norwood, MA, 2007. [24] Luscombe, N., Greenbaum, D., and Gerstein, M. What is bioinformatics? an introduction and overview. In Yearbook of Medical Informatics 2001 (2001). [25] Oinn, T., Addis, M., Ferris, J., Marvin, D., Senger, M., Greenwood, M., Carver, T., Glover, K., Pocock, M., Wipat, A., and Li, P. Taverna: a tool for the composition and enactment of bioinformatics workflows. Bioinformatics 20, 17 (2004), 3045–3054. [26] Olsen, G., Matsuda, H., Hagstrom, R., and Overbeek, R. fastdnaml: a tool for construction of phylogenetic trees of dna sequences using maximum likelihood. Computer Applications in the Biosciences 10 (1994), 41–48. [27] Piontkivska, H. Efficiencies of maximum likelihood methods of phylogenetic inferences when different substitution models are used. Molecular Phylogenetics and Evolution 31 (2004), 865–873. [28] Polanski, A., and Kimmel, M. Bioinformatics. Springer, 2007. [29] Porter, M., P´ erez-Losada, M., and Crandall, K. Model-based multi-locus estimation of decapod phylogeny and divergence times. Molecular Phylogenetics and Evolution 37 (2005), 355–369. [30] Porter, M., P´ erez-Losada, M., and Crandall, K. Advantages of multilocus sequence analysis for taxonimic studies: a case study using 10 housekeeping genes in the genus ensifer (including former sinorhizobium). International Journal of Systematic and Evolutionary Microbiology 58 (2008), 200–214. [31] Posada, D. jModelTest: phylogenetic model averaging. Molecular Biology and Evolution 25, 7 (2008), 1253–1256. [32] Posada, D., and Buckley, T. Model selection and model averaging in phylogenetics: Advantages of akaike information criterion and bayesian approaches over likelihood ratio tests. Systematic Biology 53 (2004), 793–808. [33] Rajkumar, R., Banerjee, J., Gunturi, H., Trivedi, R., and Kashyap, V. K. Phylogeny and antiquity of M macrohaplogroup inferred from complete mt DNA sequence of Indian specific lineages. BMC Evolutionary Biology 5 (Apr. 2005), 26. [34] Richards, M., Macaulay, V., Bandelt, H., and Sykes, B. Phylogeography of mitochondrial DNA in western Europe. Annals of Human Genetics 62, 3 (1998), 241–260. 41 BIBLIOGRAF´ IA [35] Ruiz-Pesini, E., Lott, M., Procaccio, V., Poole, J., Brandon, M., Mishmar, D., Yi, C., Kreuziger, J., Baldi, P., and Wallace, D. An enhanced mitomap with a global mtdna mutational phylogeny. Nucleic Acids Research 35 (2007), D823–D828. [36] Ruiz-Pesini, E., Mishmar, D., Brandon, M., Procaccio, V., and Wallace, D. Effects of purifying and adaptive selection on regional variation in human mtdna. Science 303 (2004), 223–226. [37] Salas, A., Lareu, V., Calafell, F., Bertranpetit, J., and Carracedo, A. mtdna hypervariable region ii (hvii) sequences in human evolution studies. European Journal of Human Genetics 8 (2000), 964– 974. [38] Sanderson, M., Purvis, A., and Henze, C. Phylogenetic supertrees: assembling the trees of life. Trends in Ecology and Evolution 13, 3 (1998), 105–109. [39] Stamatakis, A., Ludwig, T., and Meier, H. RAxML-III: a fast program for maximum likelihood-based inference of large phylogenetic trees. Bioinformatics 21 (2005), 456–463. [40] Steel, M. The maximum likelihood point for a phylogenetic tree is not unique. Systematic Biology 43 (1994), 560–564. [41] Strimmer, K. Maximum Likelihood Methods in Molecular Phylogenetics. PhD thesis, Ludwig-Maximilians-Universit¨at M¨unchen, 1997. [42] Sullivan, J., and Joyce, P. Model selection in phylogenetics. Annual Review of Ecology, Evolution, and Systematics 36 (2005), 445–466. [43] Torroni, A., Achilli, A., Macaulay, V., Richards, M., and Bandelt, H. Harvesting the fruit of the human mtDNA tree. Trends in Genetics 22, 6 (2006), 339–345. [44] van Oven, M., and Kayser, M. Updated comprehensive phylogenetic tree of global human mitochondrial DNA variation. Human Mutation 29 (Feb. 2008), E386–E394. 42