Repositorio Institucional de Documentos

Proyecto Fin de Carrera Ingeniería Informática Curso 2013/2014 Diseño y Estudio de herramientas para el Análisis del Índice de Conservación del ADN mitocondrial Autor: Francisco Merino Casallo Bajo la dirección de: Jorge Álvarez Jarreta Elvira Mayordomo Cámara Departamento de Informática e Ingeniería de Sistemas Área de Lenguajes de Sistemas Informáticos Escuela de Ingeniería y Arquitectura Universidad de Zaragoza Junio de 2014 Diseño y Estudio de herramientas para el Análisis del Índice de Conservación del ADN mitocondrial RESUMEN El objetivo principal de este Proyecto Fin de Carrera es el estudio y desarrollo de herramientas software para calcular de manera automatizada el índice de conservación de conjuntos de secuencias biológicas. Este conjunto de herramientas, pese a poder utilizarse de forma independiente, constituyen un potente sistema cuando se utilizan de forma combinada. El índice de conservación es un estadístico muestral que tiene especial importancia en el estudio de la patogenicidad de mutaciones, el cual consiste en determinar si una mutación puede producir una enfermedad o por el contrario, es inocua al organismo en que se produce. Constituye una de las herramientas de las que disponen los biólogos para realizar los estudios evolutivos. En la actualidad, los esfuerzos realizados por la comunidad científica para automatizar este proceso han producido contadas herramientas web con algunas limitaciones. Sorprende especialmente el no haber encontrado ninguna que permita utilizar este valor estadístico sobre grandes volúmenes de datos como un valioso instrumento adicional en este tipo de estudios. Para resolver este problema, se han implementado una serie de algoritmos haciendo uso de técnicas de programación modular para desarrollar el núcleo del sistema; y pequeñas utilidades para automatizar tareas auxiliares como el tratamiento de los conjuntos de secuencias o los análisis estadísticos requeridos para validar los experimentos realizados. Entre las dificultades encontradas durante el desarrollo de este proyecto cabe destacar la necesidad de formación previa en temas biológicos, ya que no se habían tratado estos temas a lo largo de la carrera. Por otro lado, la necesidad de realizar un tratamiento sobre los conjuntos de secuencias biológicas también implicó importantes esfuerzos en distintas fases del trabajo realizado. Por último, la obtención de información adicional relativa a la división de dichas secuencias para facilitar su análisis supuso un importante obstáculo a resolver en el último tramo del proyecto. Los resultados del sistema diseñado para los conjuntos de secuencias biológicas han sido muy positivos, especialmente teniendo en cuenta el volumen de datos con el que se ha estado trabajando, muy superior al empleado en estudios anteriores. Gracias a la información extraída de los informes generados por el sistema desarrollado se ha podido constatar la presencia de mutaciones ya conocidas por la comunidad bióloga, lo que demuestra la validez de los resultados obtenidos. I Agradecimientos Me gustaría agradecérselo a muchas personas, pido perdón por no haber podido incluir a todos los que se lo merecen, ¡no me olvido de vosotros! En primer lugar quiero agradecer a mis dos directores, Jorge y Elvira, su inestimable ayuda al final de esta aventura que comenzó hace ya casi 7 años. Sin ellos no habría redescubierto mi pasión por la biología y comprendido que nuestra profesión es mucho más divertida y reconfortante cuando la pones al servicio de otras disciplinas, aunque muchos no comprendan por qué tu sueño no es trabajar en Google, Facebook o Twitter. Desde mis inicios en el antiguo CPS he podido disfrutar de muy valiosas amistades. Javier, que siempre ha querido estar a mi lado aunque dejásemos de compartir asignaturas, intentando enseñarme la importancia de los valores éticos pese a mi supuesto fanboyismo. Rafael, que además de aguantar mis lamentos en más de una práctica comparte conmigo la idea de “Mens sana in corpore sano”. Los compañeros de Púlsar e ISC, con los que he aprendido todo eso que no te enseñan dentro de las aulas y alguna otra cosa. Y a Andrés, César, Daniel, David, Jorge y Víctor, gracias por acompañarme de principio a fin, a pesar de todo. Estaré eternamente agradecido a mis padres, que siempre me animaron a estudiar aquello que me apasionaba, aunque no fuese medicina; y con sacrificio y paciencia han costeado mi formación y me han soportado en los momentos más complicados. A mi hermana Blanca, que me enseñó que en esta vida hay que ser valiente y luchar por lo que te importa, aunque sea en el extranjero. A mi hermana María, que me enseña todos los días que no hay que tirar la toalla, por muchos obstáculos que tenga el camino. A mi hermano Lucas, compañero en muchas penas madridistas y alguna alegría mundialista, que me demuestra a diario que hay otra forma de hacer las cosas. A mi abuelo Antonio, que me enseñó la importancia de cultivar la lectura y el conocimiento, y que si hay algo que te apasiona, no debe terminar el día de tu jubilación. A mi abuela Blanca, que se esforzó mucho por mantener esa “conexión especial” que teníamos cuando era pequeño, y que es la única que me acompaña las tardes de los domingos, por muy mala que sea la película que escoja. Y a mi tía Pilar, que a su manera, siempre intenta que no me aleje demasiado del camino. Y como dicen que los últimos serán los primeros... A Elena, que me recuerda todos los días lo maravillosa que puede llegar a ser la vida si encuentras a la persona adecuada y que no hay que tener miedo al cambio pues es la clave para continuar avanzando. V Índice general 1. Introducción 1 1.1. Contexto del proyecto . . . . . . . . . . . . . . . . . . . . . . 1 1.2. Objetivos ............................. 1 1.3. Metodología y herramientas . . . . . . . . . . . . . . . . . . . 2 1.4. Software.............................. 3 1.5. Entorno tecnológico . . . . . . . . . . . . . . . . . . . . . . . 3 1.6. Estructura de la memoria . . . . . . . . . . . . . . . . . . . . 4 2. Glosario biológico 5 3. Automatización del Estudio del Índice de Conservación 7 3.1. EstadodelArte.......................... 7 3.2. Diseño............................... 8 3.3. Implementación.......................... 14 3.4. Pruebas .............................. 18 3.5. Resultados............................. 21 4. Conclusiones 27 4.1. Trabajo realizado . . . . . . . . . . . . . . . . . . . . . . . . . 27 4.2. Con vistas al futuro . . . . . . . . . . . . . . . . . . . . . . . 27 4.3. De lo profesional a lo personal . . . . . . . . . . . . . . . . . . 28 A. Diagrama de Gantt 30 B. Fundamentos biológicos 32 B.1.Basebiológica........................... 32 B.1.1. Expresión Génica . . . . . . . . . . . . . . . . . . . . . 33 B.1.2. Mutaciones Génicas . . . . . . . . . . . . . . . . . . . 33 B.2. Introducción a la bioinformática . . . . . . . . . . . . . . . . 34 B.2.1. Alineamientos de secuencias . . . . . . . . . . . . . . . 35 B.2.2. Índice de Conservación . . . . . . . . . . . . . . . . . . 35 C. Selección del Lenguaje de Programación 36 D. Manual de usuario: Conjunto de herramientas desarrolladas para el análisis del IC 38 D.1. Manual de usuario: Herramienta de cálculo del IC ...... 38 D.2. Manual de usuario: Herramienta de traducción de nucleótidos a aminoácidos ........................... 39 D.3. Manual de usuario: Herramienta de combinación de informes 40 D.4.Consejosdeuso.......................... 40 E. Gráficas de resultados 42 VII CAPÍTULO 1. INTRODUCCIÓN 1.6 Estructura de la memoria La memoria se ha dividido en varias secciones y apéndices que se describen brevemente a continuación. El primer capítulo tras esta introducción se ha titulado Glosario biológico. Como su nombre sugiere, en esta sección se encontrarán todos aquellos conceptos y definiciones, sobre todo de índole biológica, que se han considerado básicos y necesarios para poder comprender los siguientes capítulos en su totalidad. Se recomienda encarecidamente al lector que dedique unos minutos a su lectura. En el siguiente capítulo se describe en detalle el conjunto de herramientas software desarrolladas en en este proyecto. Dicho capítulo comienza describiendo las motivaciones que han guiado su construcción y a continuación se explica el proceso que se ha seguido para la creación de dichas herramientas. Se termina detallando las pruebas realizadas y los resultados obtenidos. En el último capítulo antes de los apéndices se recogen todas las conclusiones relacionadas con el proyecto. El primero de los apéndices se centra en la descripción de la asignación de tiempo dedicada a cada parte del proyecto, incluyendo el diagrama de Gantt. Debido a las limitaciones de longitud del cuerpo de la memoria, se ha incluido en el segundo apéndice un desarrollo más extenso sobre los fundamentos biológicos en los que se sustenta el proyecto. El tercer apéndice se centra en describir de forma más detallada las motivaciones que guiaron la elección del lenguaje de programación empleado. En el penúltimo apéndice se ha incluido el manual de usuario de la implementación realizada de las distintas herramientas construidas. En el último apéndice se incluyen en un tamaño más legible la Figura 3.15 y la Figura 3.16. Finalmente se adjunta la bibliografía consultada durante el proyecto. 4 2 Glosario biológico En este capítulo se encuentra la definición de aquellos conceptos que se han considerado imprescindibles para la comprensión del proyecto. Para más detalles sobre aspectos biológicos el lector puede acudir al Apéndice B. ADN: siglas en castellano de ácido desoxirribonucleico. Es una macromolécula formada por una doble cadena de nucleótidos (por lo general siguiendo una estructura de doble hélice) que forma parte de todas las células, y es usada para su desarrollo y funcionamiento. En ella se encuentra toda la información genética y es, por tanto, el componente responsable de la transmisión hereditaria. ADN mitocondrial (ADNmt): en algunas células existen unos orgánulos denominados mitocondrias. En estos orgánulos se produce la oxidación de las moléculas de glucosa, obteniendo energía para la célula. Poseen ADN propio que gestiona sus funciones internas, siendo independiente del ADN nuclear. Este ADN posee unas características únicas que lo hacen idóneo para el estudio de la patogenicidad de las mutaciones, debido a su alta tasa de mutación y a su gran conservación entre organismos de la misma especie [14, 23, 29]. Gen: segmento de una secuencia de ADN o ARN que contiene toda la información necesaria para codificar un elemento funcional de la célula. Codifican proteínas o secuencias de ARN con objetivos muy específicos. El ADNmt animal contiene, salvo alguna excepción, 37 genes en su secuencia [6]. Proteína: macromolécula formada por cadenas de aminoácidos. Es fundamental para la vida, ya que puede desempeñar una gran cantidad de funciones básicas para el correcto funcionamiento del organismo (estructural, inmunológica, enzimática, etc.). 13 de los genes contenidos en el ADNmt animal codifican proteínas. Mutación génica: alteraciones producidas en la secuencia de nucleótidos de un gen. Estos cambios pueden provocar a su vez modificaciones en las cadenas de aminoácidos que tengan como resultado efectos patógenos sobre el individuo. Dichas alteraciones se clasifican en dos grupos: neutrales y no neutrales. El primer grupo está formado por aquellas que no afectan ni a la supervivencia del organismo ni a su reproducción mientras que las 5 CAPÍTULO 2. GLOSARIO BIOLÓGICO mutaciones no neutrales están formadas tanto por aquellas que tienen un efecto beneficioso como las que son perjudiciales. Índice de conservación: estadístico muestral que permite medir la frecuencia con la que cada nucleótido o aminoácido (dependiendo del tipo de secuencia manejada) aparece en una posición concreta para un conjunto de secuencias dado. Dicho valor estadístico tiene especial importancia en el estudio de la patogenicidad de mutaciones. Alineamiento: análisis y modificación de un conjunto de secuencias con el fin de conseguir la mayor cantidad posible de nucleótidos iguales en la misma posición. Para ello se realizan inserciones de huecos (habitualmente llamados gaps) en algunas de las secuencias. Secuencia de referencia: secuencia, normalmente de ADN, de un individuo específico de una especie determinada que ha sido ampliamente estudiada, y por lo tanto, es muy poco probable que contenga errores. Suele existir una por cada especie. Además contiene información referente al principio y final de las secciones de cada uno de sus genes. 6 3 Automatización del Estudio del Índice de Conservación En este capítulo se van a exponer de forma detallada tanto las motivaciones del proyecto como todo el trabajo realizado para el desarrollo de un sistema formado por distintas herramientas software que permitan el cálculo del índice de conservación (IC) de forma automática. 3.1 Estado del Arte La motivación principal para el desarrollo de un conjunto de herramientas que permitan la automatización del cálculo del IC era la de ofrecer un potente instrumento adicional para facilitar y dotar de mayor profundidad a los estudios de mutaciones a lo largo de la evolución. Este estadístico muestral se ha estado utilizando desde hace ya varios años para determinar si una mutación es neutral o no neutral [25, 26, 30]; e incluso se ha conseguido asociar algunas de estas mutaciones no neutrales con distintos tipos de cáncer, como el de mama o endometrio [7] o el gástrico [3]. También se han llevado a cabo otros estudios en los que se utiliza el IC junto con otros instrumentos para detectar una nueva mutación en un gen del ADNmt asociada a la pérdida auditiva hereditaria [20] o enfermedades cardiovasculares como la hipertensión [9] y la miocardiopatía no compactada [18, 28]. Aunque inicialmente el cálculo del IC se realizaba de forma manual, con el paso de los años se han ido desarrollando herramientas que facilitan el cálculo de este estadístico muestral de una forma cada vez más eficiente y automatizada. Ejemplos de estas herramientas son MITOMASTER [8] y MitoTool [12]. A pesar de ello, todavía existen limitaciones en dichas propuestas que sería beneficioso subsanar de cara a mejorar la precisión y eficiencia de los resultados. MITOMASTER se centra en el cálculo del IC interespecífico, es decir, entre individuos de distintas especies, obviando su estudio entre individuos de la misma especie; mientras que MitoTool realiza un cálculo del IC también únicamente entre individuos de distintas especies y además, fija esta comparación a 43 especies de primates, lo cual limita la 7 CAPÍTULO 3. AUTOMATIZACIÓN DEL ESTUDIO DEL ÍNDICE DE CONSERVACIÓN profundidad del estudio. El desarrollo de este tipo de sistemas es cada vez más importante debido al crecimiento exponencial que se está produciendo en cuanto al número y la variedad de las secuencias disponibles y que permite sacar conclusiones cada vez con una base más sólida. El objetivo de este proyecto es desarrollar un conjunto de herramientas que permitan realizar cálculos del IC no solo de forma interespecífica como MITOMASTER y MitoTool, sino también entre individuos de una misma especie. Además, el sistema que conforman este conjunto de herramientas permite realizar un análisis más profundo al ofrecer la posibilidad de trabajar con secuencias de aminoácidos. Los resultados de dicho estudio se presentan al personal investigador de una forma clara con el objetivo de facilitar su comprensión además de ofrecer la información necesaria para determinar cuáles deberían ser los objetivos de un estudio en mayor profundidad. 3.2 Diseño La construcción del sistema se ha llevado a cabo teniendo muy en cuenta el esquema de caja negra y la idea de modularidad. El objetivo no era otro que el facilitar la comprensión por parte del usuario de qué es lo que hace cada una de las herramientas desarrolladas, sin prestar atención a cómo lo hacen. Al desarrollar módulos independientes no solo se consigue agilizar la comprensión global del sistema, sino que además se obtiene mayor robustez y se facilita su mantenimiento. El sistema construido está compuesto por tres módulos principales: 1. Cálculo del IC y generación de informes. 2. Traducción de conjuntos de secuencias de nucleótidos a aminoácidos. 3. Combinación de informes para aquellos genes que codifican proteínas. Se considera importante mencionar el hecho de que gracias a la colaboración de investigadores expertos en estos temas como Eduardo Ruiz Pesini y su doctorando Antonio Martín Navarro, se han podido aplicar sus conocimientos en la materia para determinar tanto la información útil que debía incluirse en dichos informes como el formato de los mismos. A continuación se va a proceder a describir de forma detallada el diseño de cada uno de los módulos. 8 3.2. DISEÑO La Figura 3.1 muestra el diseño del primer módulo del sistema construido. Como ya se ha comentado anteriormente, cada módulo que se expone a lo largo de este capítulo en sus distintas fases puede operar de forma independiente. Durante la fase de diseño se ha visto la posibilidad de combinar todos los módulos para formar un sistema de análisis del IC en secuencias biológicas. Figura 3.1: Diseño del primer componente del sistema. En el apartado a) se indica el diseño global del módulo: tras el preprocesamiento del conjunto de datos de entrada, realiza un análisis de la variabilidad de dichos datos y genera un informe adaptado a las necesidades indicadas por el usuario. El apartado b) refleja la paralelización de distinta granularidad según los datos que se manejan en la fase de análisis. A continuación se va a proceder a detallar las operaciones internas que este primer módulo realiza. Como se puede apreciar en el apartado a) de la Figura 3.1, el primer paso consiste en realizar el preprocesamiento del conjunto de secuencias dado. Si se está trabajando con secuencias de nucleótidos, se ofrece la posibilidad de realizar una división por genes (en el caso de secuencias completas de ADNmt) o por secciones (tanto si se trata de secuencias completas de ADNmt como si son las asociadas a un gen concreto). Por otro lado, si se trabaja con secuencias de aminoácidos solo se permite la división por secciones ya que únicamente interesa estudiar las secuencias asociadas a los genes que codifican proteínas. En ambos casos, de forma 9 CAPÍTULO 3. AUTOMATIZACIÓN DEL ESTUDIO DEL ÍNDICE DE CONSERVACIÓN adicional se permite realizar un análisis global, es decir, sin realizar ningún tipo de división. La segunda fase del preprocesamiento consiste en realizar el alineamiento de las secuencias. Esta fase puede resultar muy compleja en función de varios aspectos (precisión, longitud de las secuencias, etc.) como se comenta en la Sección 3.4 y en la Sección 3.5. Posteriormente se realiza el análisis del conjunto de secuencias indicado. Puesto que se ofrece la posibilidad de trabajar con distintos tipos de secuencias (nucleótidos o aminoácidos) es necesaria la inclusión de dos bloques diferenciados, cada uno con operaciones específicas para realizar el análisis del IC sobre un alfabeto distinto como se puede apreciar en la Tabla 3.1 y en la Tabla 3.2. El algoritmo en el que se basa el análisis consiste en contabilizar el número de veces que cada nucleótido o aminoácido aparece en cada una de las posiciones del conjunto alineado de secuencias y dividir el resultado por el número total de secuencias para obtener su valor del IC. En aquellos casos en los que el nucleótido o aminoácido contemple varias posibilidades, como por ejemplo el símbolo R en nucleótidos que equivale a un nucleótido G o a un nucleótido A, se le asigna a cada una de estas equivalencias un peso de forma equitativa (0.5 en este caso). Con este comportamiento se pretende penalizar de algún modo el hecho de no conocer con certeza el nucleótido que aparece en dicha posición. Merece especial atención el caso del gap, que tiene un peso de 0 para conseguir un efecto más penalizante sobre el cálculo del IC que el caso anterior, ya que dicho elemento representa la falta de información. Como se puede ver en el apartado b) de la Figura 3.1, dicho análisis se puede realizar de forma paralela ya que cada columna o sección del conjunto de secuencias se puede analizar de forma independiente bajo este estadístico. Pese a ello hay que tener muy en cuenta que solo resulta interesante su paralelización en caso de que la longitud de las secuencias así lo recomiende, especialmente en aquellos alineamientos en los que no se haya realizado la división de las secuencias por genes o secciones. Una vez se ha llevado a cabo el análisis de la variabilidad de cada uno de los genes o secciones en los que se ha dividido el conjunto de secuencias, se generan los informes correspondientes. En caso de que no se haya realizado ninguna división previa se generará un único informe. Además, también se calcula la secuencia más frecuente (SMF), que incluye en cada posición el nucleótido (o aminoácido) que tiene un IC más elevado, es decir, aquel que aparece más veces en dicha posición. La Figura 3.2 muestra la estructura del segundo módulo del sistema. La motivación que guió la creación de este módulo consiste en realizar un análisis más profundo de los efectos que la mutación o mutaciones que ha sufrido el ADNmt provocan en los organismos. Dicho análisis permite prestar especial atención a las posiciones de las secuencias con un IC próximo al 100% 10 3.2. DISEÑO Símbolo Equivalencia Significado Peso G G Guanina 1 A A Adenina 1 T T Timina 1 C C Citosina 1 R G o A Purina 0.5 por nucleótido Y T o C Pirimidina 0.5 por nucleótido M A o C Amina 0.5 por nucleótido K G o T Cetona 0.5 por nucleótido S G o C Interacción fuerte (enlaces 3 H) 0.5 por nucleótido W A o T Interacción débil (enlaces 2 H) 0.5 por nucleótido H A o C o T no G 0.33 por nucleótido B G o T o C no A 0.33 por nucleótido V G o C o A no T 0.33 por nucleótido D G o A o T no C 0.33 por nucleótido N G o A o T o C Cualquier nucleótido 0.25 por nucleótido -gap Ningún nucleótido 0 Tabla 3.1: Elementos del alfabeto utilizado en secuencias de nucleótidos y el peso asignado. pero sin alcanzarlo. En estas posiciones existe una mínima probabilidad de mutación, lo que hace que sea más interesante analizar si dicha mutación afecta a la secuencia de aminoácidos, en cuyo caso la probabilidad de que esta sea mortal aumenta considerablemente. Como se ha comentado brevemente en el Capítulo 2 y se describe en mayor detalle en el Apéndice B, existen mutaciones en el ADNmt que no tienen efectos perjudiciales para el organismo, bien por el hecho de que no provocan cambios en la secuencia de aminoácidos correspondiente, o porque dichos cambios no afectan a las funcionalidades de las proteínas. En GenBank aunque el número de secuencias de ADNmt es muy elevado y sigue creciendo rápidamente, no se encuentran las traducciones a proteínas de todas estas secuencias. Por ello, se consideró necesario el desarrollo de una herramienta capaz de realizar la traducción de secuencias de nucleótidos y permitir un análisis más exhaustivo del IC. A la hora de diseñar esta herramienta se ha empleado como guía el proceso biológico que realiza la célula. Recordar que solo tiene sentido realizar esta traducción sobre el conjunto de secuencias asociadas a aquellos genes que codifican proteínas (13 de los 37 genes que se encuentran en el ADNmt animal). Al igual que en el primer módulo, se incluye una fase de preprocesamiento. Sin embargo, existen ciertas diferencias: la división por genes solo 11 CAPÍTULO 3. AUTOMATIZACIÓN DEL ESTUDIO DEL ÍNDICE DE CONSERVACIÓN Símbolo Equivalencia Significado Peso G G Glicina 1 A A Alanina 1 V V Valina 1 L L Leucina 1 I I Isoleucina 1 P P Prolina 1 F F Fenilalanina 1 Y Y Tirosina 1 C C Cisteína 1 M M Metionina 1 H H Histidina 1 K K Lisina 1 R R Arginina 1 W W Triptófano 1 S S Serina 1 T T Treonina 1 D D Ácido aspártico 1 E E Ácido glutámico 1 N N Asparagina 1 Q Q Glutamina 1 B D o N Ácido aspártico o Asparagina 0.5 por aminoácido Z E o Q Ácido glutámico o Glutamina 0.5 por aminoácido - - Terminador 1 X X Desconocido 1 =gap Ningún aminoácido 0 Tabla 3.2: Elementos del alfabeto utilizado en secuencias de aminoácidos y el peso asignado. selecciona aquellos que codifican proteínas, y no está permitida la división por secciones. Después se lleva a cabo la traducción a proteínas, que consiste en traducir cada triplete o codón (se denomina así a los conjuntos de tres nucleótidos) en su aminoácido correspondiente (los 64 tripletes posibles codifican solo 20 aminoácidos distintos). En la Tabla B.1 se indican los codones que codifican cada uno de los aminoácidos. Como se puede ver en el apartado b) de la Figura 3.2, dicha traducción se puede realizar de forma paralela ya que la traducción de cada secuencia del conjunto se puede llevar a cabo de forma independiente. A continuación se va a proceder a describir las operaciones internas que realiza el tercer módulo, cuya estructura se muestra en la Figura 3.3. El 12 3.2. DISEÑO Figura 3.2: Diseño del segundo componente del sistema. En el apartado a) se indica el diseño global del módulo: tras el preprocesamiento del conjunto de secuencias de ADNmt de entrada, realiza la traducción de dichas secuencias a proteínas. El apartado b) refleja la paralelización de distinta granularidad según los datos que se manejan en la fase de traducción. Figura 3.3: Diseño del tercer componente del sistema. Realiza la combinación en un único informe de los generados previamente tanto para las secuencias de nucleótidos como de aminoácidos. primer paso consiste en extraer la SMF de los informes tanto de nucleótidos como de aminoácidos. Para obtener dichos informes es necesaria la invo13 CAPÍTULO 3. AUTOMATIZACIÓN DEL ESTUDIO DEL ÍNDICE DE CONSERVACIÓN Figura 3.10: Tiempo de ejecución en segundos de los diferentes alineamientos de 22954 secuencias de ADNmt humano. Figura 3.11: Tiempo de ejecución en segundos de los diferentes alineamientos de 442 secuencias de ADNmt primate. El principal motivo por el que se han omitido los fragmentos en esta versión del sistema es porque su alineamiento supone un problema muy complejo y su tratamiento supera los objetivos de este proyecto. Cabe destacar que durante el proceso de construcción de las consultas, han sido detectados varios errores de notación en algunas de las secuencias disponibles en GenBank que han sido notificados para su futura corrección. La existencia de este tipo de secuencias inválidas en los primeros conjuntos de ADNmt primate descargados hizo que el alineamiento global no tuviese ninguna validez práctica y fuese necesaria una división por familias taxonó20 3.5. RESULTADOS micas como etapa intermedia. En estas primeras pruebas de alineamiento global se generaba la inserción de más de 300000 gaps, lo cual supone un incremento de 20 veces la longitud de la secuencia más larga en el conjunto de secuencias original. Tras resolver este problema, se han conseguido alinear las secuencias de ADNmt primate de forma global, sin que fuese necesaria la división por familias taxonómicas. Durante el proceso de obtención de conjuntos de secuencias de ADNmt primate se detectó que para algunas especies existen varias secuencias de referencia. Esto supone un gran problema de cara a la validez de los resultados ya que no existe un protocolo biológico para seleccionar una secuencia de referencia en el caso de que existan varias para la misma especie. Debido al coste que supone resolver este problema (formación y pruebas) ha quedado pendiente el desarrollo de mecanismos que permitan llevar a cabo dicho proceso. Por último, también se ha llevado a cabo la actualización de la herramienta de PhyloDAG encargada de descargar los conjuntos de secuencias de GenBank, para obtener información adicional referente a la localización de cada gen en las secuencias de ADNmt. Esto ha permitido simplificar el proceso de división por genes y alineamiento de las secuencias. Como se puede apreciar en la Figura 3.12 y la Figura 3.13, a pesar de que se sufre una penalización en cuanto al tiempo de descarga de las secuencias, el espacio en disco que ocupan los ficheros que las almacenaban y el consiguiente tiempo de lectura de dichos ficheros, el ahorro en costes posterior resulta muy superior y beneficioso, especialmente a la hora de resolver el problema que se ha comentado previamente sobre las múltiples secuencias de referencia asociadas a una misma especie. 3.5 Resultados Uno de los objetivos principales en la construcción de las herramientas desarrolladas a lo largo de este proyecto era el poder realizar análisis del IC con grandes cantidades de datos. El coste de realizar dicho análisis de forma manual resultaba inadmisible desde el punto de vista de los investigadores. Antes de que comenzase este proyecto, el alineamiento más grande realizado por el grupo de bioinformática de la Universidad de Zaragoza era el llevado a cabo en el proyecto ZARAMIT (centrado en la construcción de árboles filogenéticos), que contaba con 7390 secuencias de ADNmt humano. Además, dicho alineamiento se realizó mediante una reconstrucción bottomup con la ayuda de un árbol filogenético, lo que simplifica el problema pero exige que se disponga de un árbol filogenético que contenga todas las secuencias que se vayan a alinear. En este proyecto se ha estado trabajando con alineamientos de 22954 secuencias de ADNmt humano para las cuales no se ha construido todavía ningún árbol filogenético. Además en el tramo 21 CAPÍTULO 3. AUTOMATIZACIÓN DEL ESTUDIO DEL ÍNDICE DE CONSERVACIÓN Figura 3.12: Análisis de descarga y almacenamiento de 442 secuencias de ADNmt primate: a) Tiempos de descarga y lectura en segundos; b) Tamaño del fichero en MB. final se ha creado otro grupo de estudio en el que se incluyen 442 secuencias de ADNmt de distintas especies de primates, en las que se han excluido a los seres humanos. Los últimos esfuerzos se han enfocado en dotar al sistema de los mecanismos necesarios para realizar un tratamiento gen a gen no solo en humanos, sino también en el resto de las especies animales. Como se ha comentado en la Sección 3.4, la decisión de utilizar una herramienta de alineamiento u otra afecta, entre otras cosas, a la precisión del resultado obtenido, como queda reflejado en la Figura 3.14. En el caso de las secuencias de ADNmt humano, existen dos secuencias 22 3.5. RESULTADOS Figura 3.13: Análisis de descarga y almacenamiento de 1000 secuencias de ADNmt humano: a) Tiempos de descarga y lectura en segundos; b) Tamaño del fichero en MB. de referencia (rCRS y RSRS), lo que genera una disyuntiva biológica. La rCRS es una revisión de la primera secuencia de referencia publicada para este tipo de secuencias y la segunda es una secuencia sintética creada en el año 2012 que se sitúa en la raíz del árbol filogenético del ADNmt humano. Puesto que la disyuntiva todavía no se ha resuelto se han duplicado las pruebas realizadas con el objetivo de comprobar si existe algún tipo de diferencia entre los resultados obtenidos con cada una de ellas. A la vista de los resultados, se ha decidido descartar el obtener más resultados con la RSRS ya que no se han detectado diferencias que justifiquen el coste del 23 CAPÍTULO 3. AUTOMATIZACIÓN DEL ESTUDIO DEL ÍNDICE DE CONSERVACIÓN Figura 3.14: Longitud de los alineamientos de 442 secuencias de ADNmt primate. doble análisis, como se puede apreciar en la Figura 3.15 y la Figura 3.16. Como se ha comentado al comienzo de la Sección 3.3, se ha realizado un estudio estadístico del IC en cada gen del ADNmt humano con la ayuda de una de las herramientas incluidas en la librería PhyloDAG. En la Figura 3.15 y la Figura 3.16 se muestran los resultados obtenidos normalizados dividiendo los valores absolutos de las secciones asociadas a cada gen entre su longitud. Figura 3.15: Gráfica con los resultados del cálculo del IC del alineamiento normalizado de 22954 secuencias de ADNmt humano con división por genes utilizando Mafft - auto y la secuencia de referencia rCRS. Como se puede apreciar en la Figura 3.15 y la Figura 3.16 el gen de la Treonina (Thr) consta de un número de posiciones con una alta tasa de variabilidad muy superior al resto. Tras contactar con el personal investigador biólogo se confirmó que en dicho gen existen una serie de mutaciones génicas conocidas. Solo en la región de control (también conocida como región hipervariante) se pueden encontrar valores similares. La hipótesis más aceptada actualmente es que esta región no contiene información genética por lo que 24 3.5. RESULTADOS Figura 3.16: Gráfica con los resultados del cálculo del IC del alineamiento normalizado de 22954 secuencias de ADNmt humano con división por genes utilizando Mafft - auto y la secuencia de referencia RSRS. el encontrar posiciones altamente conservadas puede resultar interesante. Por otro lado, hay que tener en cuenta que la elección de un método de alineamiento concreto no solo afecta al tiempo de ejecución como se ha visto en la Figura 3.10 y Figura 3.11, sino que como se ha comentado anteriormente, la precisión del alineamiento también varía y por lo tanto, los IC son objeto de modificaciones debido a este problema. 25 4 Conclusiones 4.1 Trabajo realizado El trabajo realizado supone una nueva aportación al conjunto de herramientas software para estudios de patogenicidad de mutaciones, extendiendo su ámbito a las secuencias de aminoácidos y a especies distintas del ser humano con el objetivo de profundizar en el estudio de los efectos que las mutaciones tienen en los organismos. Se han cumplido satisfactoriamente no solo aquellos objetivos que se plantearon inicialmente, sino también aquellos que han ido surgiendo a lo largo del desarrollo de este proyecto. 4.2 Con vistas al futuro Aunque el sistema desarrollado se considera por si solo un potente instrumento adicional para el estudio de la patogenicidad de mutaciones, existen varios aspectos en los que se va a continuar trabajando. Se considera especialmente importante desarrollar los mecanismos necesarios para establecer una correlación entre las secuencias de las distintas especies y poder realizar estudios del IC interespecífico. Otro de los aspectos más interesantes es el de ampliar los estudios realizados a otras especies, no solo primates como se ha hecho en la segunda mitad de este proyecto, sino a otros grandes grupos como pueden ser los mamíferos. El hecho de manejar conjuntos de secuencias con mayor variabilidad entre ellas hace que las mutaciones en aquellas posiciones en las que el IC esté cercano al 100% sin llegar a alcanzarlo, tengan mayores probabilidades de ser mortales. Como se ha comentado con anterioridad, el hecho de utilizar GenBank como fuente de los conjuntos de secuencias biológicas provoca que la información que obtenemos no sea todo lo homogénea que nos gustaría, ya que se trata de un proyecto colaborativo. Una de las consecuencias de esto consiste en que no todas las secuencias cuentan con la sección hipervariante ya que entre la comunidad científica estuvo extendida durante mucho tiempo la creencia de que dicha sección no contenía información biológica relevante, por lo que no resultaba beneficioso su secuenciación tanto desde el punto de vista académico como desde el punto de vista aplicado. A este respecto, se considera interesante 27 CAPÍTULO 4. CONCLUSIONES añadir la capacidad de detectar aquellas secuencias que no cuenten con esta sección ya que la ausencia de información por este motivo no debería tener el mismo efecto penalizante en el cálculo del IC que aquellos generados en el proceso de alineamiento. En los estudios llevados a cabo hasta la fecha el número de secuencias sin región hipervariante ha sido ínfimo por lo que no se ha considerado necesario, pero con vistas a un futuro próximo en el que los estudios se amplíen a otras especies o se incluyan los fragmentos de secuencias disponibles, resultará conveniente. Con el objetivo de extender y facilitar el uso de dichas herramientas, también se ha pensado en ofrecer la capacidad de trabajar con conjuntos de secuencias almacenadas en formatos distintos a FASTA, que pese a ser el más extendido, no es el único manejado por la comunidad bióloga. Por último, pensando en adaptar el sistema para que pueda trabajar con secuencias de ADN nuclear, se considera también importante implementar aquellos mecanismos de paralelización que durante la etapa de diseño fueron valorados. Quedaron descartados en la implementación por resultar innecesarios al trabajar con secuencias de ADNmt, con una longitud en el caso peor inferior a los 22000 nucleótidos (tras la fase de preprocesamiento). Sin embargo, al plantear la aplicación del sistema a secuencias de ADN nuclear, que en el caso de los seres humanos tiene una longitud total aproximada de 3200 millones de nucleótidos, estos mecanismos resultan muy necesarios. Para que el lector se haga una idea más clara de la magnitud del problema que se está describiendo, se ofrece el siguiente ejemplo: en el peor caso antes mencionado, los tiempos de ejecución rondaban los 25 minutos para el análisis del IC. Puesto que el algoritmo tiene un coste temporal lineal, esto supone que en el caso de disponer del mismo número de secuencias de ADN nuclear humano, este mismo análisis tendría un tiempo de ejecución superior a 7 años. 4.3 De lo profesional a lo personal Desde el punto de vista profesional debo realizar una valoración muy positiva del trabajo que se ha llevado a cabo a lo largo de este proyecto ya que se ha conseguido desarrollar un sistema completo y robusto que pueda ser utilizado como herramienta de trabajo por los investigadores, a pesar de los distintos obstáculos que se han ido encontrando a lo largo del camino. Lo más valioso ha sido la posibilidad de demostrar la capacidad de trabajo, esfuerzo y aprendizaje, no solo en cuestiones relacionadas con la informática, sino aquellas que tienen un corte más biológico. Todo esto ha desembocado en la posibilidad de disfrutar de un contrato de 3 meses de duración en el mismo grupo donde se ha desarrollado este trabajo dentro del proyecto de investigación BASMATI y queda todavía pendiente la participación en la redacción de dos artículos de investigación de futura publicación, lo que ha hecho que se considere cursar un máster para continuar la formación en estos 28 4.3. DE LO PROFESIONAL A LO PERSONAL temas o incluso la realización de una tesis doctoral. Pese a que al empezar en la carrera de Ingeniería Informática mi interés por la biología quedó enterrado, debo agradecer a mis directores que me ofreciesen la posibilidad de realizar este proyecto, que me ha hecho recuperar mi pasión por la biología y descubrir el poder de la informática cuando la pones al servicio de otras especialidades. Por otro lado, desde el punto de vista personal, me he sentido muy cómodo, valorado y acogido en el grupo de bioinformática, permitiéndome colaborar de forma activa en otros proyectos de investigación. 29 C Selección del Lenguaje de Programación Tras realizar un primer análisis del problema al que se buscaba dar solución se decidió escoger como paradigma de programación la programación por procedimientos, que se deriva de la programación estructurada. La elección de dicho paradigma permitió explotar la técnica de programación modular, cuyas ventajas más destacables son la facilidad de depurar, actualizar y modificar el código. Entre los múltiples lenguajes de programación que se consideraron inicialmente se encontraban lenguajes de bajo (C) y alto nivel (C++, Java, Python). Pese a que lenguajes como C permiten una gestión de los recursos computacionales más eficiente, puesto que en este proyecto no se trabajaba con sistemas que tuviesen grandes limitaciones en este aspecto, se decidió elegir un lenguaje de alto nivel y aprovechar sus múltiples ventajas. Además, debido a la creciente complejidad de las arquitecturas de los microprocesadores modernos, los compiladores para lenguajes de alto nivel cada vez generan código más eficiente. Los lenguajes candidatos fueron C++, Java y Python por múltiples razones (experiencia previa, documentación, gestión de excepciones, robustez, etc.). Pese a que la experiencia previa utilizando C++ y Java era muy superior, varios factores fueron determinantes a la hora de seleccionar a Python. Dicha elección como lenguaje principal de este proyecto se basó en: legibilidad, librería estándar muy extensa y potente; existencia de la librería BioPython, que incluye herramientas para la bioinformática; y la utilización de dicho lenguaje en las herramientas desarrolladas por Jorge Álvarez, que permitieron simplificar la complejidad del problema a resolver. Por último, merece la pena comentar que pese a que Python se un lenguaje interpretado y esto penaliza el tiempo de ejecución de los algoritmos, debido a la alta interacción con herramientas externas como Mafft, dicha penalización se puede considerar despreciable. D Manual de usuario: Conjunto de herramientas desarrolladas para el análisis del IC D.1 Manual de usuario: Herramienta de cálculo del IC El programa de cálculo del IC dispone de los siguientes parámetros de entrada: fich_secs [obligatorio]: ruta de acceso al fichero que contiene el conjunto de secuencias de ADNmt o proteínas. El fichero tiene que estar en formato FASTA. dir_informe [opcional]: ruta del directorio donde se guardará el fichero con el informe generado. Por defecto, el directorio actual. De no existir el directorio se creará. umbral [obligatorio]: indica el límite inferior o superior del IC. En el informe se indican aquellas posiciones del conjunto de secuencias que tienen un IC superior o inferior al indicado por este umbral. tipo_umbral [opcional]: limite inferior o superior (por defecto el límite es superior). Permite realizar búsquedas de las posiciones con un IC tanto por encima como por debajo del umbral establecido. rango [opcional]: rango de la sección a analizar (por defecto la sección comprende la longitud total del de las secuencias, si se introducen alineadas, o de la longitud del alineamiento una vez efectuado). Permite realizar secciones no asociadas a los distintos genes del ADNmt. gen [opcional]: gen al que pertenece el conjunto de secuencias a analizar. Necesario para que en el informe generado se indique correctamente la posición absoluta de cada nucleótido. Por defecto el primer nucleótido de la secuencia comienza en 1, y en caso de indicarse un gen, D.2. MANUAL DE USUARIO: HERRAMIENTA DE TRADUCCIÓN DE NUCLEÓTIDOS A AMINOÁCIDOS la primera posición corresponde a la posición inicial del gen en la secuencia de ADNmt. En la versión actual el parámetro solo es aplicable en ADNmt humano. informe_detallado [opcional]: bandera (más conocida por el término anglosajón flag) que indica el deseo de obtener como salida un informe detallado del análisis del IC realizado (por defecto no se incluye indicando que se desea obtener un informe básico). En la Sección 3.5 se explican las principales diferencias entre estos dos tipos de informes. verbose [opcional]: bandera (o flag) que indica el deseo de obtener información por pantalla relativa al estado en el que se encuentra la ejecución del módulo (por defecto no se incluye ningún tipo de información del estado del proceso). En el caso de que las secuencias no estén alineadas, se solicitará por pantalla al usuario que indique la herramienta de alineamiento que se desea utilizar y su configuración. Dicho programa genera la siguiente salida: Fichero de texto que contiene el informe con los resultados del estudio del IC para cada posición del conjunto de secuencias indicado en el parámetro de entrada. D.2 Manual de usuario: Herramienta de traducción de nucleótidos a aminoácidos El programa de traducción de nucleótidos a aminoácidos dispone de los siguientes parámetros de entrada: fich_secs [obligatorio]: ruta de acceso al fichero que contiene el conjunto de secuencias de ADNmt. El fichero tiene que estar en formato FASTA. En el caso de que las secuencias no estén alineadas, se solicitará por pantalla al usuario que indique la herramienta de alineamiento que se desea utilizar y su configuración. Dicho programa genera la siguiente salida: Alineamiento del conjunto de proteínas resultante de aplicar el método desarrollado en este módulo que simula el proceso de traducción a proteínas que ocurre en las células. 39 APÉNDICE D. MANUAL DE USUARIO: CONJUNTO DE HERRAMIENTAS DESARROLLADAS PARA EL ANÁLISIS DEL IC D.3 Manual de usuario: Herramienta de combinación de informes El programa de combinación de informes dispone de los siguientes parámetros de entrada: ruta_nucleotidos [obligatorio]: ruta de acceso al fichero que contiene el informe del IC de las secuencias de nucleótidos. ruta_aminoacidos [obligatorio]: ruta de acceso al fichero que contiene el informe del IC de las secuencias de aminoácidos. Ambos informes deben pertenecer al mismo conjunto de secuencias de ADNmt. Dicho programa genera la siguiente salida: Fichero de texto que contiene el informe combinado en el que se muestra la información recopilada de los informes indicados en los parámetros de entrada. D.4 Consejos de uso Puesto que el número de ficheros que componen el sistema puede crecer de forma exponencial, tanto al añadir nuevos métodos de alineamiento como al ampliar el número de conjuntos de secuencias que se manejan (por ejemplo, porque se decide ampliar el estudio a otras especies), se propone establecer una estructura lo más sencilla posible en el sistema de ficheros. De este modo, se recomienda al usuario mantener en directorios separados los scripts que componen el núcleo del sistema, los ficheros que almacenan los distintos conjuntos de secuencias, y los informes que dichos scripts generan. Se muestra una propuesta de la estructura de ficheros en la Figura D.1. Otro de los detalles que se han tenido en cuenta es el de generar ficheros cuyos nombres sean lo más autoexplicativos posible, no solo en el caso de los scripts desarrollados, sino también en aquellos que almacenan los alineamientos y los informes generados. Para facilitar la comprensión del lector a este respecto, se incluyen a continuación tres ejemplos: hmtDNA.fasta: Conjunto de secuencias de ADNmt humano. hmtDNA_rCRS_aligned_mafft_auto_ATP6_peptide_[all]_0.5.txt: Informe básico del cálculo del IC con un umbral superior a 0.5, correspondiente al conjunto de secuencias de aminoácidos asociadas al gen ATP6 de ADNmt humano alineadas con Mafft - auto y la secuencia de referencia rCRS. 40 D.4. CONSEJOS DE USO Figura D.1: Propuesta de la estructura del sistema de ficheros. pmtDNA_aligned_mafft_auto_[all]_0.99_details.txt: Informe detallado del cálculo del IC con un umbral superior a 0.99, correspondiente al conjunto de secuencias de ADNmt primate alineadas con Mafft - auto. 41 E Gráficas de resultados Figura E.1: Figura 3.15 ampliada. Figura E.2: Figura 3.16 ampliada. 43 Bibliografía [1] J. Alvarez-Jarreta. Análisis teórico-práctico de métodos de inferencia filogenética basados en selección de modelos y métodos de superárboles. Master’s thesis, Centro Politécnico Superior, Universidad de Zaragoza, 2010. [2] Dennis A. Benson, Ilene Karsch-Mizrachi, David J. Lipman, James Ostell, and David L. Wheeler. Genbank. Nucleic Acids Research, 36(suppl 1):D25–D30, 2008. [3] Rui Bi, Wen-Liang Li, Ming-Qing Chen, Zhu Zhu, and Yong-Gang Yao. Rapid identification of mtdna somatic mutations in gastric cancer tissues based on the mtdna phylogeny. Mutation Research/Fundamental and Molecular Mechanisms of Mutagenesis, 709–710(0):15 – 20, 2011. [4] Roberto Blanco. Definición y prototipo de herramienta de análisis filogenético para el adn mitocondrial humano. Master’s thesis, Centro Politécnico Superior, Universidad de Zaragoza, 2008. [5] Roberto Blanco and Elvira Mayordomo. Zaramit: A system for the evolutionary study of human mitochondrial dna. In Sigeru Omatu, MiguelP. Rocha, José Bravo, Florentino Fernández, Emilio Corchado, Andrés Bustillo, and JuanM. Corchado, editors, Distributed Computing, Artificial Intelligence, Bioinformatics, Soft Computing, and Ambient Assisted Living, volume 5518 of Lecture Notes in Computer Science, pages 1139–1142. Springer Berlin Heidelberg, 2009. [6] Jeffrey L. Boore. Animal mitochondrial genomes. Nucleic Acids Research, 27(8):1767–1780, 1999. [7] M Brandon, P Baldi, and D C Wallace. Mitochondrial mutations in cancer. Oncogene, 25(34):4647–4662, 2006. [8] Marty C. Brandon, Eduardo Ruiz-Pesini, Dan Mishmar, Vincent Procaccio, Marie T. Lott, Kevin Cuong Nguyen, Syawal Spolim, Upen Patil, Pierre Baldi, and Douglas C. Wallace. Mitomaster: a bioinformatics tool for the analysis of mitochondrial dna sequences. Human Mutation, 30(1):1–6, 2009. [9] Hong Chen, Jing Zheng, Ling Xue, Yanzi Meng, Yan Wang, Bingjiao Zheng, Fang Fang, Suxue Shi, Quiaomeng Qiu, Pingping Jiang, Zhongqiu Lu, Jun Qin Mo, Jianxin Lu, and Min-Xin Guan. The 12s rrna a1555g mutation in the mitochondrial haplogroup d5a is responsible for maternally inherited hypertension and hearing loss in two chinese pedigrees. European Journal of Human Genetics, 20(6):607–612, 2012. 45