Bioinformatics tools for the analysis of plant-associated bacterial genomes

Memoria presentada por Pedro Manuel Martínez García para optar al grado de Doctor por la Universidad de Málaga Bioinformatics tools for the analysis of plant-associated bacterial genomes Directores: Dr. Cayo J. Ramos Rodríguez Dr. Pablo Rodríguez Palenzuela Catedrático. Área de Genética Departamento de Biología Celular, Genética y Fisiología Universidad de Málaga Instituto de Hortofruticultura Subtropical y Mediterránea “La Mayora” (IHSM) Catedrático. Área de Bioquímica y Biología Molecular Departamento de Biotecnología Universidad Politécnica de Madrid Centro de Biotecnología y Genómica de Plantas (CBGP) Universidad de Málaga Málaga, 2015 AUTOR: Pedro Manuel Martínez García EDITA: Publicaciones y Divulgación Científica. Universidad de Málaga Esta obra está sujeta a una licencia Creative Commons: Reconocimiento - No comercial - SinObraDerivada (cc-by-nc-nd): Http://creativecommons.org/licences/by-nc-nd/3.0/es Cualquier parte de esta obra se puede reproducir sin autorización pero con el reconocimiento y atribución de los autores. No se puede hacer uso comercial de la obra y no se puede alterar, transformar o hacer obras derivadas. Esta Tesis Doctoral está depositada en el Repositorio Institucional de la Universidad de Málaga (RIUMA): riuma.uma.es COMITÉ EVALUADOR Presidente Dr. José Manuel Palacios Alberti Departamento de Biotecnología Centro de Biotecnología y Genómica de Plantas (CBGP) Universidad Politécnica de Madrid Secretario Dr. Javier Ruíz Albert Departamento de Biología Celular, Genética y Fisiología Instituto de Hortofruticultura Subtropical y Mediterránea (IHSM) Universidad de Málaga Vocales Dr. Miguel Redondo Nieto Departamento de Biología Universidad Autónoma de Madrid Dr. Antonio Jesús Pérez Pulido Departamento de Biología Molecular e Ingeniería Bioquímica Universidad Pablo de Olavide Dr. José María Vinardell González Departamento de Microbiología Universidad de Sevilla Suplentes Dra. Carmen Beuzón López Departamento de Biología Celular, Genética y Fisiología Instituto de Hortofruticultura Subtropical y Mediterránea (IHSM) Universidad de Málaga Dr. Francisco Javier López Baena Departamento de Microbiología Universidad de Sevilla Área de Genética. Departamento de Biología Celular, Genética y Fisiología. Instituto de Hortofruticultura Subtropical y Mediterránea (IHSM) Universidad de Málaga-Consejo Superior de Investigaciones Científicas Dr. CAYO J. RAMOS RODRÍGUEZ, Catedrático del Área de Genética del Departamento de Biología Celular, Genética y Fisiología de la Universidad de Málaga, y Dr. PABLO RODRÍGUEZ PALENZUELA, Catedrático del Departamento de Biotecnología de la Universidad Politécnica de Madrid y Centro de Biotecnología y Genómica de Plantas (CBGP), INFORMAN: Que, PEDRO MANUEL MARTÍNEZ GARCÍA ha realizado en este Departamento y bajo su dirección el trabajo titulado “Bioinformatics tools for the analysis of plant-associated bacterial genomes”, que constituye su memoria de Tesis Doctoral para aspirar al grado de Doctor por la Universidad de Málaga. Y para que así conste, y tenga los efectos que correspondan, en cumplimiento de la legislación vigente, extienden el presente informe. En Málaga, a 10 de Abril de 2015. Fdo. Cayo J. Ramos Rodríguez Fdo. Pablo Rodríguez Palenzuela AGRADECIMIENTOS Durante estos dos años y medio son muchas las personas que han intervenido directa o indirectamente en el presente trabajo, con lo que no puedo más que dedicarles a ellos esta Tesis Doctoral. En primer lugar, quiero agradecer a mis dos directores, Cayo y Pablo, por confiar en mi para desarrollar este proyecto. Gracias por todo lo que me habéis enseñado y por todas y cada una de las sugerencias que han hecho que esta Tesis llegue a buen puerto. Gracias a los dos por valorar mi trabajo y considerar siempre mis opiniones. Cayo, gracias por elegir el papelito sin foto que era mi currículum cuando los azares de la administración andaluza lo hicieron llegar a tu despacho. Te la jugaste conmigo, y no sabes cuánto te lo agradezco. Pablo, gracias por tu respaldo constante en todo lo que hago y por las conversaciones y elucubraciones que hemos compartido en este tiempo. Quiero también dar las gracias a Emilia, mi jefa number three. Gracias por contar conmigo en tus proyectos, por tus consejos científicos y por dedicarme tiempo siempre que lo he necesitado. No puedo olvidar al resto de compis del 285. A nuestras alumnas Bea y Silvia, por refrescar el labo con esa juventud que tanto envidiamos. A Mariela, por venir todos los días con un sonrisón latino que quita el sentío. A Chechu, por tantas conversaciones, trascendentales o cotidianas, pero siempre interesantes. A Saray, la burgalesa menos siesa que haya conocido un andaluz. Ha sido un gustazo compartir contigo estos meses, me llevo una compi de las buenas. A mi Isabella, que hasta mi llegada fue la más guapísima del CBGP. Gracias por tu compañía estos años, no sabes cómo te echamos de menos. Y sobre todo a mi Piluca de mi corazón! Mi super secre, mi confidente y mi reina mora...gracias por todo lo que has hecho por mi y por ser taaaaaaaaan buenísima gente. A todos mis compañeros del CBGP, en especial a Bea, Marco, Adri y los Alex, los “otros” bioinformáticos. Gracias por formar parte del rato de ocio diario que han sido nuestros almuerzos. Gracias a Bea, mi gurú del aprendizaje automático. A Alex Junior, por tu complicidad salmantina y tantas frikadas que he aprendido contigo. Y gracias también a Alex Senior, por tus consejos técnicos, procedimentales, y cómo no, por alguna que otra noche loca que nos hemos pegao por Lega y Madrid. A los compañeros y jefazos de las áreas de Genética y Microbiología de la Universidad de Málaga, en especial a Antonio de Vicente. Gracias por tu eficiencia en la gestión de todo lo referente al genoma de la UMAF0158 y por tus expertas indicaciones necrosis-apicalianas. Y por supuesto gracias a Conchita, mi única compi de despacho en estos años. Por efímera que fuera tu estancia, me encantó tenerte conmigo esos meses. Los cromosomas circulares van por ti! No querría pasar sin dar las gracias a Jesús Mercado por confiar en nosotros para la anotación y análisis del genoma de la PICF7. Ha sido un placer colaborar contigo. Gracias a todos los amig@s que he ido dejando en Tomares, Sevilla, Salamanca, Florencia, Madrid, Copenhage y Barcelona, sin excepción. Gracias a Luisito, Mateo, Juan, Jorge, Enrique y Pedro por el pasado y el presente juntos. Gracias a Arbe, Julito, Chico, Rafa, Javi, María, Mariwa, Pedro Ángel, Esther, Gloria, Juanma, Grego, Magneto, Antuan, Quinito, Carlete, Carmela, Recu, Gema y Joselito, mis tomareños del alma. Por todos los momentos juntos. Los fines de año locos en la sierra, las infinitas ferias y carnavales. Por los Al Rumbos habidos y por haber. Siempre todos. Siempre juntos. Todos sois parte de mi y como tal todos y cada uno sois parte de este trabajo. Sin vosotros no sería el piltrafa en el que me he convertido. Gracias a Benavente, Luisma, Moro y Elvi, mis europeos perpetuos. A mis primos, en especial a Chio y Alberto, sois más que familia. A Rodri, Agus, Morta y Antonio, mi familia salmantina. A Fernan y Laura, mi familia fiorentina. A Peter, Thomas y Steven, mi familia danesa. A Ro y Arlette, mi familia catalana. A Mari, Patas y Aza, mi familia madrileña (gracias por la portada Patorras!!). A Fidel, Tito, Sulas y Madaleno, por amenizar los inamenizables años de facultad. A mis padres y hermana. Por vuestro eterno apoyo y cariño, tan incondicionales que sólo la genética puede explicarlos. A Julia. Por todo. Por ser el mejor público para mis chistes, la mejor comensal cuando cocino y la mejor receptora involuntaria de cuantas parrafadas y comeduras de tarro me han surgido en el desarrollo de esta Tesis. Porque tu amor incondicional no lo explica la genética. Te quiero. Gracias a todos los políticos de esta nuestra madre patria por vuestra inconmensurable labor, tan necesaria en estos momentos difíciles. No hay dinero en el mundo que recompense el esfuerzo y dedicación que prestáis al servicio público. Vuestro constante compromiso para con la sociedad os eleva a la condición de mártires del bien común, ÍNDICE Índice Prefacio 1 Resumen 5 Introducción General 9 1. Bioinformática y análisis genómico . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 1.1. Análisis de secuencias biológicas . . . . . . . . . . . . . . . . . . . . . . . . . . 12 1.2. Bases de datos biológicas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13 1.3. Herramientas de anotación . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14 1.4. Minería de datos y aprendizaje automático . . . . . . . . . . . . . . . . . . . . 15 2. Interacciones planta-bacteria . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17 2.1. Enfermedades vegetales producidas por bacterias fitopatógenas . . . . . . 20 Mancha bacteriana . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20 Marchitez bacteriana . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 Tumores y agallas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 Costra bacteriana . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23 Podredumbre blanda . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23 Cancro bacteriano . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 Necrosis apical . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 Pseudomonas syringae pv. syringae UMAF0158, modelo de estudio de la necrosis apical del mango . . . . . . . 25 2.2. Bacterias beneficiosas y control biológico de enfermedades vegetales . . . 26 Rizobacterias, biocontrol y promoción del crecimiento vegetal . . . . 27 Pseudomonas fluorescens PICF7, agente de biocontrol de la verticilosis del olivo . . . . . . . . . . . . . . . . . . . 28 3. Bioinformática para el análisis genómico de bacterias asociadas a plantas . . . . . 30 Objetivos 33 I Índice Chapter I. Bioinformatics analysis of the complete genome sequence of the mango tree pathogen Pseudomonas syringae pv. syringae UMAF0158 reveals traits relevant to virulence and epiphytic lifestyle 37 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62 Material and Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63 Supplementary Material . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66 Chapter II. Complete genome sequence of Pseudomonas fluorescens strain PICF7, an indigenous root endophyte from olive (Olea europaea L.) and effective biocontrol agent against Verticillium dahliae 67 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79 Material and Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80 Chapter III. T346Hunter: A novel web-based tool for the prediction of type III, type IV and type VI secretion systems in bacterial genomes 83 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 90 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96 Material and Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97 Supplementary Material . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99 Chapter IV. A supervised machine-learning approach for the prediction of bacterial associations with plants 101 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 105 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 107 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 120 Material and Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 121 Supplementary Material . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 125 II Índice Conclusiones 127 Bibliografía 133 Anexo: Otras publicaciones 163 III PREFACIO Prefacio Esta Tesis Doctoral se ha dirigido al estudio bioinformático de los genomas de dos cepas bacterianas asociadas a cultivos de mango y olivo, ambos de relevancia en la agricultura andaluza, así como al desarrollo de herramientas computacionales para analizar genomas bacterianos en el contexto de las interacciones planta-bacteria. Durante el desarrollo de la misma se contribuyó en varios trabajos como consecuencia de la colaboración con el Grupo de Bacterias Fitopatógenas del Centro de Biotecnología y Genómica de Plantas de la Universidad Politécnica de Madrid (CBGP-UPM-INIA), coordinado por los Dres. Emilia López Solanilla y Pablo Rodríguez Palenzuela, este último director de la presente Tesis Doctoral junto el Dr. Cayo Ramos de la Universidad de Málaga. En uno de estos trabajos se estudió el papel de la percepción de la luz en la virulencia de Pseudomonas syringae pv. tomato DC3000, agente causal de la mancha bacteriana del tomate (Río-Álvarez et al., 2013). Asimismo, se colaboró en el estudio funcional de los quimiorreceptores implicados en el proceso de entrada en la planta a través de heridas de Dickeya dadantii 3937, agente causal de la podredumbre blanda de la patata (Río-Álvarez et al., 2014). Los resultados obtenidos en ambos trabajos no se han incorporado en el cuerpo de esta Tesis Doctoral, si bien las publicaciones científicas derivadas de los mismos se han incluido en el Anexo. Un tercer trabajo se centró en analizar la resistencia de P. syringae pv. tomato DC3000 frente a compuestos antimicrobianos de tomate mediada por bombas de extrusión multidroga, cuyos resultados se incluyeron en un manuscrito que actualmente se encuentra en proceso de revisión. Por otro lado, y en consecuencia de la estrecha colaboración con el resto de integrantes de la línea de investigación “Biología y control de enfermedades de plantas” del Instituto de Hortofruticultura Subtropical y Mediterránea “La Mayora” – Universidad de Málaga (IHSM-UMA-CSIC), se analizaron los genomas de dos cepas de Bacillus amyloliquefaciens, potenciales agentes de biocontrol frente al oídio de cucurbitáceas (Romero et al., 2007). El manuscrito que describe los resultados de este estudio ha sido recientemente aceptado para su publicación en una revista de relevancia en el área y se encuentra en proceso de producción. En lo que respecta al cuerpo de esta Tesis Doctoral, éste se divide en cuatro capítulos que cubren sus dos objetivos principales. Los capítulos 1 y 2 abordan el análisis bioinformático de los genomas de las cepas P. syringae pv. syringae UMAF0158 y Pseudomonas fluorescens PICF7, respectivamente. P. syringae pv. syringae UMAF0158 es el agente causal de la necrosis apical del mango, enfermedad que ha sido objeto de intenso estudio por parte de nuestro grupo de investigación desde su descripción en 1998 (Cazorla et al., 1998). Los resultados obtenidos tras el análisis de su genoma están siendo revisados en una revista de impacto para su publicación. Por su parte, P. fluorescens PICF7 es una cepa probada 3 Prefacio experimentalmente como agente efectivo contra la verticilosis del olivo, y fue aislada en 2003 por el grupo de Etiología y Control de Enfermedades de Cultivos del Instituto de Agricultura Sostenible (IAS-CSIC) (Mercado-Blanco et al., 2004), dirigido por el Dr. Jesús Mercado-Blanco, con quien nuestro grupo ha establecido diversas colaboraciones. El estudio bioinformático del genoma de esta cepa se publicó en la revista Standards in Genomic Sciences (SIGS) con título Complete genome sequence of Pseudomonas fluorescens strain PICF7, an indigenous root endophyte from olive (Olea europaea L.) and effective biocontrol agent against Verticillium dahliae (Martínez-García et al., 2015a). Por su parte, los capítulos 3 y 4 se focalizan en herramientas bioinformáticas de anotación y análisis de genomas bacterianos, con especial énfasis en las interacciones planta-bacteria. El capítulo 3 describe una herramienta online para la identificación genómica de sistemas de secreción bacterianos, publicada en la revista PLOS ONE con título T346Hunter: A novel web-based tool for the prediction of type III, type IV and type VI secretion systems in bacterial genomes (Martínez-García et al., 2015b). En el capítulo 4 se aborda el problema de la predicción de estilos de vida bacterianos usando información genómica. Para ello se implementó un pipeline que combina búsquedas de secuencias por homología con técnicas de aprendizaje automático supervisado para clasificar genomas bacterianos asociados a plantas. El correspondiente manuscrito que describe dicha herramienta se encuentra actualmente en proceso de revisión para su publicación en una revista relevante en patogénesis. 4 Introducción General 1. Bioinformática y análisis genómico En la actualidad, la cantidad ingente de datos generados por las nuevas tecnologías es sencillamente inabordable sin la ayuda de herramientas computacionales, que permiten automatizar la extracción de información relevante para generar nuevo conocimiento. Cuando la naturaleza de dichos datos es de índole biológica, el conjunto de herramientas que los tratan se enmarcan en el ámbito de la bioinformática. Según Luscombe et al. (2001), la bioinformática conceptualiza la biología en términos de moléculas (en un sentido físico-químico) y aplica técnicas informáticas para entender y organizar a gran escala la información asociada a dichas moléculas. Es, pues, en el contexto de la biología molecular en el que la bioinformática toma sentido y se presenta como una disciplina indispensable. En la era de la biotecnología (Ganguly et al., 2014), en que las técnicas de secuenciación masiva producen a diario inmensas cantidades de información genómica, la necesidad de implementar herramientas bioinformáticas para procesar y analizar dicha información es, cuanto menos, crucial. En genómica, la secuenciación consiste en la aplicación de métodos bioquímicos con objeto de determinar el orden de los nucleótidos (A, C, G y T) que componen una secuencia de ADN. Desde que Sanger y Gilbert (Sanger et al., 1977; Maxam, Gilbert, 1977), desarrollaron las primeras aproximaciones en los años 70, las tecnologías de secuenciación han experimentado una evolución notable. Como consecuencia, en 1995 se publicó por vez primera el genoma secuenciado de un organismo de vida libre, la bacteria Haemophilus influenzae (Fleischmann et al., 1995). Seis años después vio la luz el primer borrador del genoma humano (Lander et al., 2001), gracias a un proyecto coordinado de investigación sin precedentes. Este hito, aunque excepcional, supuso un gasto de 3.000 millones de dólares (un dólar por nucleótido leído), lo que implicaba un coste prohibitivo para los laboratorios. La necesidad de abaratar la secuenciación llevó al desarrollo de lo que hoy conocemos como secuenciación de nueva generación (next generation sequencing), que permite secuenciar un genoma humano por 1.800 dólares y un genoma bacteriano por 500. No es de extrañar que sólo en la base de datos del NCBI (http://www.ncbi.nlm.nih.gov) haya depositados a día de hoy (Marzo de 2015) cerca de 200 millones de secuencias, de las cuales alrededor de 10.000 corresponden a genomas de cepas bacterianas. Esta desmesurada producción de datos genómicos traslada el problema de costes de secuenciación al de ensamblar, procesar y manejar dichos datos de manera que se pueda extraer información útil de ellos. Es así como los métodos computacionales adquieren una importancia sustancial en la biología moderna, adaptando sus paradigmas a problemas específicos a nivel molecular. La etapa actual de progreso en la generación de datos biológicos coincide, no por 11 Introducción General casualidad, con el momento de más auge en lo que a tecnologías de la información se refiere. Irremediablemente, el efecto de la sociedad de la información también ha sido notable en el desarrollo de la electrónica, la ingeniería del software y las telecomunicaciones, lo que ha supuesto un estímulo en el avance de las técnicas de procesamiento y análisis inteligente de datos (Peek, Swift, 2012; Pea, Vityaev, 2010). Así pues, si la generación masiva de datos biológicos deriva en el problema evidente de procesar dicha información, ésto, por fortuna, sucede en un momento en que la informática está preparada para dar soluciones efectivas. Por su relevancia en esta Tesis Doctoral, a continuación se revisan algunas de las técnicas computacionales más relevantes para el tratamiento y análisis de datos biológicos. 1. 1. Análisis de secuencias biológicas El análisis de secuencias biológicas consiste en la aplicación de métodos analíticos a secuencias de ADN, ARN o aminoácidos con objeto de inferir propiedades biológicas de las mismas, como función, estructura o evolución (Durbin et al., 1998). En este contexto, los métodos de comparación de secuencias son de especial importancia. Por un mero proceso evolutivo, secuencias nuevas son necesariamente resultado de adaptaciones de otras ya existentes. Si una secuencia tiene un alto grado de similitud con otra conocida, es razonable inferir que su función será a su vez semejante. Cuando esto ocurre, se dice que ambas secuencias son homólogas, y se asume que tienen un origen evolutivo común. Aquí toma sentido el concepto de alineamiento, que permite comparar dos o más secuencias resaltando zonas de alta similitud, indicando posibles relaciones funcionales entre los genes o proteínas comparados. Desde el punto de vista de la programación, una secuencia (por biológica que sea) no es más que una sucesión de caracteres, o formalmente, cadena de caracteres. Típicamente, los algoritmos computacionales de alineamiento implementan métodos de optimización que buscan el mejor alineamiento entre dos cadenas en función de un sistema de puntuación determinado. Las cadenas pueden alinearse local o globalmente. En el primer caso, se buscan subsegmentos de una secuencia A en otra secuencia B, lo que resulta especialmente útil para búsquedas en bases de datos de gran tamaño. Un ejemplo de programa basado en este tipo de alineamiento es BLAST (Altschul, 2005). Los alineamientos locales permiten identificar secuencias con cierto grado de similitud, pudiendo ser o no homólogas. Por su parte, los algoritmos de alineamiento global comparan las secuencias enteras, forzando al alineamiento a ocupar la longitud total de las mismas. Este tipo de alineamiento es computacionalmente más caro, siendo útil cuando se comparan un número no demasiado grande de secuencias similares, como ocurre con los alineamientos 12 Introducción General múltiples. Una aplicación directa de estos alineamientos es el estudio filogenético de un conjunto de secuencias. Dado que las secuencias que forman el conjunto a alinear tienen una supuesta relación evolutiva, mediante el alineamiento de secuencias sucesivamente menos emparentadas se infiere el árbol filogenético del conjunto. Así es como operan los llamados métodos jerárquicos o de árbol, como ClustalW (Larkin et al., 2007) o T-Coffee (Notredame et al., 2000). Otro uso común de los alineamientos múltiples es el de, dado el alineamiento de un conjunto de secuencias pertenecientes a una familia determinada, comprobar si una secuencia nueva pertenece a dicha familia en función de lo bien que se alinee con el resto. Si queremos afirmar que esa secuencia pertenece a la familia en cuestión, sería deseable que al menos los patrones más conservados en el alineamiento múltiple estén presentes. Para capturar esta información se suele hacer uso de modelos probabilísticos, en especial los llamados modelos ocultos de Markov o HMM (hidden Markov models) (Zhang, Wood, 2003; Eddy, 2004). Así, los distintos algoritmos de alineamiento a aplicar dependen en gran medida del conocimiento biológico que se quiera extraer, haciendo de la comparación de secuencias una parte esencial de la bioinformática. De hecho, una vez obtenida la secuencia del genoma de un organismo concreto, el siguiente paso no es otro que el de compararla con secuencias ya existentes, para lo que el investigador tiene a su disposición una suerte de repositorios públicos con información sobre genes, proteínas, y otros componentes genéticos conocidos. 1.2. Bases de datos biológicas Las bases de datos biológicas son a día de hoy herramientas indispensables para los científicos a la hora de abordar cualquier tipo de fenómeno biológico, ya sea la estructura de las biomoléculas y sus interacciones o la propia evolución de los organismos. El conocimiento biológico moderno está en gran medida almacenado en infinidad de repositorios generalistas o especializados que crecen en proporción a la capacidad de generación de información, haciendo especialmente laboriosa la tarea de asegurar la coherencia de los datos. Por ello, iniciativas como la de la revista Nucleic Acids Research (NAR), que publica anualmente una exhaustiva clasificación con las bases de datos biológicas y bioinformáticas más relevantes (http://www.oxfordjournals.org/our_journals/nar/database/c/), son particularmente útiles. El DDBJ (DNA Data Bank of Japan), el EMBL-EBI (European Molecular Biology Laboratory) y el GenBank del NCBI (The National Center for Biotechnology Information) forman el consorcio INSDC (International Nucleotide Sequence Database Collaboration), e intercambian información a diario para formar el actual mayor banco de secuencias de ADN existente (Nakamura et al., 2013). GenBank, en concreto, almacena una colección 13 Introducción General anotada de todas las secuencias de ADN disponibles públicamente (Benson et al., 2014), hasta un total de 181.336.445 en Febrero de 2015. En lo que respecta a repositorios de proteínas, UniProtKB (Universal Protein Knowledgebase, Boutet et al. 2007) surge con la misión de proveer a la comunidad científica de una base de datos pública, exhaustiva y de alta calidad de secuencias de proteínas e información funcional de las mismas. A día de hoy, UniProtKB es sin duda la base de datos de proteínas más prominente y usada, y consta a su vez de dos repositorios: Swiss-Prot (547.599 secuencias anotadas manualmente) y TrEMBL (90.860.905 secuencias anotadas de forma automática). Otras bases de datos tienen como objeto almacenar familias de proteínas. Tal es el caso de PFAM (Finn et al., 2014), que consiste en una extensa colección de alineamientos múltiples y modelos ocultos de Markov, y que resulta especialmente útil en la identificación de regiones funcionales (dominios) en secuencias de proteínas. Estas y otras bases de datos biológicas son rastreadas a diario por los investigadores, y algunas de ellas constituyen el núcleo de información en que se basan los programas de anotación automática. 1.3. Herramientas de anotación Una vez obtenida la secuencia de un genoma pasamos a la fase de anotación del mismo. En genómica se distinguen dos tipos de anotaciones, la estructural y la funcional. La primera se encarga de identificar qué genes están presentes, dónde se ubican en la secuencia y qué proteínas codifican. Ésto puede inferirse ab initio directamente de las propiedades específicas de la secuencia de ADN, como la existencia de codones de inicio y terminación donde empiezan y acaban los marcos de lectura. Este tipo de información se integra en detectores de contenido, en base a los cuales se identifican regiones codificantes e intergénicas y se determina la posición de los genes. Un ejemplo de programas de predicción de genes ab initio es Glimmer (Delcher et al., 2007), que se basa en los llamados modelos de Markov interpolados para identificar marcos de lectura en bacterias, arqueas y virus. Otro forma de predecir genes es en base a posibles homologías con otros genes conocidos. Haciendo uso de herramientas de comparación de secuencias como las descritas en el apartado 1.1, se puede analizar el genoma en busca de genes razonablemente parecidos a otros ya estudiados. En la práctica, los pipelines de anotación automática de genomas suelen usar modelos híbridos, de manera que en primera instancia predicen las posiciones de los genes y las proteínas que codifican, y a continuación éstos se comparan con otros conocidos mediante búsquedas en bases de datos biológicas. Una de las herramientas de anotación de genomas bacterianos más usadas es PGAP (Prokaryotic Genome Annotation Pipeline), 14 Introducción General diseñada por el NCBI (Angiuoli et al., 2008). En base a una secuencia de ADN bacteriano, PGAP combina predicción de genes con búsquedas por homología para producir una anotación general que almacena en archivos ASN.1 (Abstract Syntax Notation One), avalados por ISO. NCBI ToolKit ofrece una serie de programas para, a partir de estos archivos, extraer la información sustancial del genoma anotado, generando a su vez ficheros en formato FASTA y tablas fácilmente procesables. Tanto los ASN.1 como las tablas y archivos FASTA son los formatos más extendidos entre la comunidad científica y suelen ser compatibles con la mayoría de programas de análisis de secuencias, lo que supone una ventaja para su procesamiento. Es por ello que han sido los formatos usados en el desarrollo de esta Tesis Doctoral para anotar y tratar genomas bacterianos. Otra de las herramientas más usadas de anotación de genomas bacterianos es RAST (Rapid Annotation using Subsystem Technology; Overbeek et al. 2014). RAST ofrece una serie de ventajas respecto a otros sistemas de anotación, lo que le ha hecho adquirir una popularidad notable en los últimos años. Entre estas ventajas, destaca la rapidez con que el usuario obtiene la anotación del genoma. En cuestión de minutos, RAST produce una anotación pormenorizada tomando como entrada la secuencia de un genoma bacteriano junto con información básica sobre la cepa. RAST se basa en un repositorio propio que integra anotaciones de una gran variedad de fuentes distintas, curadas a mano por expertos microbiólogos. A raíz de esa base se generan anotaciones automáticas de bastante precisión, generando a su vez anotaciones específicas de colecciones de proteínas funcionalmente relacionadas, o subsistemas. Este tipo de anotación aporta un componente de especificidad que no ofrecen otros sistemas como PGAP. Pero esta búsqueda de subsistemas es dependiente del repositorio de genomas sobre el que RAST trabaja, con lo que a efectos de cepas de reciente descripción, la anotación que este servidor ofrece no dista demasiado de la que produce PGAP. Es por ello que en la actualidad, y especialmente cuando se trabaja con bacterias no descritas previamente en la literatura científica, resulta más práctico complementar la anotación general con anotaciones específicas, que producen outputs con información exhaustiva de sistemas moleculares concretos. 1.4. Minería de datos y aprendizaje automático La minería de datos utiliza métodos estadísticos y de inteligencia artificial, entre otros, para interpretar grandes volúmenes de datos y estructurarlos de manera que la información subyacente sea comprensible. Su fin último es explotar y analizar datos para ayudar a la toma de decisiones, así como generar sistemas inteligentes capaces de entenderlos (Hastie et al., 2005; Gorunescu, 2011). En bioinformática, una de las aplicaciones de minería de 15 Introducción General datos más usadas es la generación de modelos predictivos, en especial aquellos enfocados a problemas de clasificación (Larrañaga et al., 2006). Para ello se utilizan técnicas de aprendizaje automático, que inducen modelos probabilísticos en base a un conjunto de datos llamados de entrenamiento. Estos datos de entrenamiento constituyen la información previa sobre la que se infiere el modelo de representación, que posteriormente puede aplicarse a nuevos datos y observar si éstos se comportan como los datos de entrenamiento. El aprendizaje de los datos puede ser supervisado o no supervisado (Love, 2002). En el primer caso, el entrenamiento se lleva a cabo con conocimiento de qué subconjuntos de datos son representativos de aquello que se quiere clasificar. En el segundo, este conocimiento es desconocido, por lo que la tarea del programa es encontrar patrones que ayuden a definir y clasificar los datos en distintos grupos, en función de su contenido. Para ilustrar ambos tipos de aprendizaje supongamos que disponemos de una base de datos de secuencias de ADN humano y queremos generar un modelo para discernir si esas secuencias corresponden a promotores. Los datos de entrenamiento contendrían secuencias ubicadas tanto en promotores como fuera de ellos. Para aplicar aprendizaje supervisado, tendríamos que etiquetar las secuencias que corresponden a promotores y las que no, de manera que el programa tenga información previa sobre los tipos de datos a clasificar. Si no dispusiéramos de esa información, al aplicar aprendizaje no supervisado esperaríamos que el modelo infiriera los tipos de secuencias, por lo que en este caso toma especial relevancia la interpretación de la información obtenida. Este tipo de aprendizaje no podría usarse para decidir si una secuencia es un promotor o no, pero sí podría, por ejemplo, clasificar las secuencias entre aquellas que tienen un alto contenido en GC, AT, y el resto. De esta información se pueden inferir conocimiento, si por ejemplo tenemos en cuenta que el 70% de los promotores humanos tienen un alto contenido en islas CpG (Saxonov et al., 2006). Como se ha señalado, aquí la interpretación del investigador es fundamental. 16 Introducción General 2. Interacciones planta-bacteria La extraordinaria capacidad de adaptación de las bacterias a diferentes condiciones ambientales les ha permitido colonizar prácticamente cada rincón de la Tierra. Como en cualquier organismo vivo, los mecanismos que dictan el estilo de vida de las bacterias se basan en su material genético y las modificaciones que éste sufre como resultado de la interacción de las mismas con el medio. Así, entender esta capacidad de adaptación y los componentes moleculares subyacentes a la misma pasa en gran medida por el estudio de sus genomas. El genoma de una bacteria está esencialmente formado por un cromosoma, generalmente circular y cerrado por enlace covalente. Muchas bacterias disponen también de una o varias moléculas de ADN extracomosómico, normalmente también cerrado y circular, llamado plásmido. Aunque estos últimos suelen codificar información para funciones no esenciales para la vida de la bacteria, sí albergan genes que les pueden proporcionan propiedades fenotípicas útiles en un contexto de adaptación al crecimiento en determinados medios. Por lo general, el tamaño de un plásmido en una bacteria es significativamente menor que el del cromosoma. El tamaño de los genomas bacterianos conocidos varían entre 139 kilobases (kb) y 14 megabases (Mb) (López-Madrigal et al., 2011; Chang et al., 2011), aunque la mayoría de ellos oscila entre 2 y 6 Mb (Koonin, Wolf, 2008). Como hemos dicho, la ubiquidad de las bacterias se explica por su capacidad de desarrollar estrategias adaptativas que les permiten sobrevivir en prácticamente cualquier ambiente posible. Parte esencial de esas estrategias son las distintas interacciones que las bacterias establecen con otros seres vivos, ya sean más complejos, como plantas, animales y hongos, como con otros microorganismos. A menudo, las bacterias son capaces de penetrar en otros organismos y proliferar en su interior. Esta capacidad de invasión reside en que cuentan con los componentes genéticos necesarios para colonizar los tejidos del huésped, invadirlos y en ocasiones sintetizar factores de virulencia que provocan enfermedad. En general, las bacterias que habitan en organismos eucariotas se benefician de los nutrientes y del hábitat protegido que éstos les proporcionan. Desde el punto de vista del huésped, esta relación puede ser de mutualismo, comensalismo o parasitismo, en función de si es beneficiosa, neutral o perjudicial, respectivamente. Se cree que la mayoría de las interacciones que los seres humanos tenemos con los microorganismos que habitan en nuestro cuerpo son mutualistas o comensalistas (Fukada, 2014) y, desde un punto de vista evolutivo, dicha afirmación tiene pleno sentido. Dado que las bacterias han colonizado el mundo billones de años antes de la existencia del ser humano, en el proceso de evolución de los organismos complejos, éstas han debido desempeñar un papel fundamental en el desarrollo de sus funciones vitales. No es de extrañar que el número estimado de células 17 Introducción General bacterianas presentes en el cuerpo humano sea diez veces superior al de las propias células humanas (Pappas, 2009). De hecho, sólo el microbioma intestinal alberga entre 500 y 1000 especies distintas de bacterias, cuyo número total de genes es cien veces mayor que el que contiene todo el genoma humano (Gill et al., 2006). De entre la batería de genes presentes en el genoma de Lactobacillus casei, por ejemplo, aquellos implicados en la utilización de nutrientes como el azúcar, cuya concentración es escasa en el lumen intestinal, parecen jugar un papel fundamental en la colonización y el establecimiento de la bacteria en este tejido (Licandro-Seraut et al., 2014). Asimismo, el hospedador obtiene una serie de ventajas inducidas por el la actividad de las bacterias que lo habitan. Un estudio in vivo demostró que Acetobacter pomorum dispone de genes que modulan la ruta de señalización de insulina en Drosophila melanogaster, cuya disrupción afecta significativamente al crecimiento de dicho organismo (Shin et al., 2011). Pero el mutualismo entre bacterias y organismos complejos no se reduce, ni mucho menos, al reino animal. De hecho, uno de los ejemplos más paradigmáticos de relación mutualista es la asociación que establecen las bacterias del género Rhizobium con las plantas leguminosas, que proporcionan a las bacterias un hábitat seguro y nutritivo mientras se benefician de la capacidad de estas para fijar nitrógeno, y así crecer en suelos con bajas concentraciones de este elemento (Gourion et al., 2015). Las bacterias rizobiales se establecen endosimbióticamente dentro de nódulos que se forman en las raíces de las plantas leguminosas como consecuencia de esta interacción. Estas plantas producen flavonoides, que son reconocidos por receptores bacterianos y que inducen la expresión de genes de nodulación, implicados en la formación de los nódulos radiculares (Haag et al., 2013). Por su especial relevancia en sanidad, la forma de interacción bacteriana tradicionalmente más estudiada es la que establecen las bacterias patógenas con sus distintos huéspedes. No en vano, las enfermedades infecciosas siguen siendo una de las principales causas de mortalidad en los países en vías de desarrollo, más de la mitad de las cuales están provocadas por bacterias (Molicotti et al., 2014). Tuberculosis, peste, neumonía o diarrea son causa directa de la infección de diversos patógenos bacterianos, que cuentan con un arsenal de genes implicados en la síntesis de distintos factores de virulencia, necesarios para establecerse en nuestro organismo y provocar enfermedad. Bacterias del género Streptococcus, como Streptococcus pneumoniae, agente causal de la neumonía, albergan genes relacionados con el metabolismo de carbohidratos complejos, que les sirven para la adquisición de nutrientes esenciales, adherirse a los tejidos e interferir con las funciones del sistema inmune (Shelburne et al., 2008). Ciertas especies de Streptococcus portan además elementos genéticos móviles, como los transposones del tipo Tn916/Tn1545, que contienen genes responsables 18 Introducción General de la resistencia a antibióticos (Santoro et al., 2014). La relativa facilidad con que este tipo de determinantes génicos se transmiten a otras bacterias por transferencia horizontal contribuye en gran medida a la propagación de resistencia a antibióticos de unas bacterias a otras, ya sean o no del mismo género, con el correspondiente impacto en salud pública que ello conlleva. Otros patógenos como Bacillus anthracis, agente causal del ántrax, basan su capacidad infecciosa principalmente en la secreción de compuestos tóxicos. Típicamente, las cepas de B. anthracis portan, además del cromosoma, dos megaplásmidos (pXO1 y pXO2) que albergan genes encargados de la síntesis y secreción de toxinas, como las llamadas toxina causante de edema (edema factor) y toxina letal (lethal factor) (Keim et al., 2009; Brossier, Mock, 2001). Por otro lado, estudios recientes han mostrado que enteropatógenos comúnmente asociados a animales y humanos, como Salmonella, pueden tambier infectar plantas, dando lugar a intoxicaciones alimentarias provocadas por la ingesta de frutas y vegetales (Schikora et al., 2012; Wiedemann et al., 2015). Este género bacteriano es el principal causante de infecciones en mamíferos como la gastroenteritis o la fiebre tifoidea, y se estima que, sólo en los Estados Unidos, provoca alrededor de un millón y medio de casos de enfermedad al año (Mead et al., 1999). Se ha demostrado que algunos genes de Salmonella typhimurium, como los encargados de la síntesis de celulosa, juegan un papel fundamental en el proceso de adhesion a los tejidos de las plantas (Lapidot et al., 2006). Asimismo, la translocación de proteínas efectoras (o efectores) al interior del huésped resulta esencial en el proceso de infección. El cromosoma de Salmonella suele contener al menos dos operones que codifican las subunidades estructurales de los llamados sistemas de secreción de tipo III o T3SS (type III secretion system), cuya implicación en la asociación bacteria-huésped se ha descrito ampliamente en la bibliografía (Gophna et al., 2003; Izoré et al., 2011). Sin embargo, los mecanismos específicos del papel del T3SS en la interacción de Salmonella con plantas, como la activación de la expresión del mismo, o la función de las proteínas exportadas, están aún por dilucidarse (Wiedemann et al., 2015). La determinación de dichos mecanismos ayudará a comprender el tipo de asociación que establecen las enterobacterias con las plantas, en especial si atendemos a la importancia que tienen los T3SS como factores de virulencia en ciertas bacterias fitopatógenas. Pseudomonas syringae, por ejemplo, transloca efectores a través del T3SS que alteran el correcto desarrollo y suprimen el sistema inmune en plantas de tabaco y tomate, entre otras (Block et al., 2008). Tanto el T3SS como sus efectores asociados y otras herramientas moleculares definen el tipo de asociación bacteriana con las distintas plantas que infectan, y son objeto de intenso estudio actualmente. 19 Introducción General Dado su peso específico en este trabajo, a continuación se desarrollan en mayor detalle algunas de las interacciones establecidas entre bacterias fitopatógenas y sus plantas huésped. 2.1. Enfermedades vegetales producidas por bacterias fitopatógenas Como parte de su mecanismo de adaptación, algunas bacterias son capaces de invadir los tejidos de ciertas plantas y, en ocasiones, producir enfermedades. La pérdida o contaminación de cultivos por enfermedades agrícolas no sólo tiene un impacto agroalimentario evidente, sino también económico. La agricultura es en muchos países el principal motor de la economía, con lo que la propagación de plagas vegetales pueden tener consecuencias devastadoras. Desde que en 1878 T. J. Burril demostrara que el fuego bacteriano del peral era causa directa de la acción de la bacteria Erwinia amylovora (Burrill, 1878), se han descrito incontables infecciones causadas por diversas bacterias fitopatógenas, así como muchos de los mecanismos moleculares implicados en los procesos de infección de plantas (Vidaver, Lambrecht, 2004; Mansfield et al., 2012). En los siguientes apartados se revisan algunas enfermedades relevantes de plantas y sus correspondientes agentes causales. Mancha bacteriana Es el tipo más frecuente de enfermedad vegetal, caracterizándose por la aparición de manchas de distintos tamaños ya sea en hojas, flores, tallo o frutos. Cuando las manchas se propagan con cierta celeridad la infección da lugar a la llamada quema bacteriana, que en ocasiones puede destruir prácticamente la totalidad de la superficie de la planta, marchitándola y provocando la muerte de muchos tejidos (Agrios, 2005). La mayoría de las manchas y quemas bacterianas tienen su origen en la infección de cepas pertenecientes a los géneros Pseudomonas yXanthomonas, especialmente las causadas por patovares (pvs.) de las especies P. syringae yXanthomonas campestris. P. syringae pv. tomato es el agente causal de la mancha bacteriana del tomate, una enfermedad propagada a nivel mundial que genera significativas pérdidas económicas desde la década de los 70, y que parece favorecerse en climas húmedos y frescos (Agrios, 2005). Su sintomatología se caracteriza por manchas pequeñas y oscuras en las hojas, denominadas “manchas pardas”, rodeadas por un halo clorótico amarillo, que con frecuencia producen defoliación. Las manchas aparecen también en flores, tallo y frutos. La patogenicidad de P. syringae pv. tomato se basa en el T3SS y sus efectores, así como, entre otros factores, en la síntesis y secreción de distintas moléculas fitotóxicas. Tal es el caso de la coronatina, una fitotoxina producida por varios patovares de P. syringae que induce la formación de los halos cloróticos en la planta (Brooks et al., 2005). Se cree que la coronatina mimetiza al ácido 20 Introducción General que combinan formas diferentes de control, como control biológico, físico, químico y cultural, para paliar enfermedades vegetales de forma viable económicamente y minimizando el impacto ecológico (Lewis, Papavizas, 1991). En lo que respecta al control biológico, éste introduce un elemento adicional en el triángulo de interacción planta-patógeno-ambiente, los organismos antagonistas. De entre estos, cabe destacar determinados microorganismos, ya sean naturales o modificados, que permiten reducir los efectos indeseables de los patógenos vegetales, al tiempo que favorecen el desarrollo de otros microorganismos beneficiosos para las plantas, como los enemigos naturales. Dado su peso específico en este trabajo, en lo que sigue hacemos énfasis en el uso de bacterias como agentes de biocontrol y promotores del crecimiento vegetal. Rizobacterias, biocontrol y promoción del crecimiento vegetal La rizosfera comprende una zona del suelo en que las interacciones entre las raíces de las plantas y los microorganismos existentes son de un dinamismo y especificidad excepcionales. La acumulación de exudados vegetales genera un entorno rico en nutrientes, que sirven de reclamo para los microorganismos que habitan el suelo, resultando en un aumento en la concentración de biomasa y actividad microbiana (Philippot et al., 2013; McNear, 2012). Se cree que entre el 1 y el 2% de las bacterias que habitan la rizosfera promueven directa o indirectamente el crecimiento vegetal (Bonfante, Anca, 2009), destacando los géneros Bacillus yPseudomonas como los más prominentes (Kloepper, 1981). Una forma de contribución directa a este crecimiento viene dada por la capacidad de algunas bacterias de incrementar la disponibilidad de nutrientes que la planta necesita. Un ejemplo de ello es la fijación biológica de nitrógeno que, como hemos visto, llevan a cabo determinadas bacterias, que incorporan este elemento a la rizosfera y proveen a las plantas del sustento necesario para crecer (Rodrigues et al., 2008; Hatayama et al., 2005). En ocasiones la contribución bacteriana viene dada por la síntesis de biomoléculas, típicamente hormonas, implicadas directamente en el crecimiento de la planta. Así como la producción de hormonas por parte de patógenos vegetales puede alterar perniciosamente el crecimiento de la planta, en bacterias promotoras del crecimiento vegetal (PGPR, del inglés, plant growth promoting rhizobacteria) la síntesis de estas moléculas ayuda al correcto desarrollo de la misma. Algunas rizobacterias sintetizan citoquininas y giberelinas que intervienen en la madurez de los brotes (Van Loon, 2007), o auxinas, que promueven la formación de raíces laterales (Yang et al., 2009). Por otro lado, los efectos indirectos de las PGPR se centran en su capacidad de paliar la severidad de determinadas infecciones, así como de presentar actividad antagonista frente a determinados patógenos. Para ello, ciertas bacterias producen compuestos como 27 Introducción General antibióticos, toxinas o enzimas líticas que generan un efecto deletéreo en microorganismos fitopatógenos. Un estudio reciente mostró que la síntesis de los antibióticos bacilomicina y macrolactina por parte de la cepa Bacillus amyloliquefaciens NJN-6 presenta efectos antagonistas contra Fusarium oxysporum yR. solanacearum, respectivamente (Yuan et al., 2012). Por otro lado, en el proceso de competencia por los nutrientes en la rizosfera, la síntesis de sideróforos puede suponer una ventaja respecto a potenciales competidores. Una capacidad superior de captar hierro, cuya disponibilidad es limitada en ambientes aeróbicos, puede tener un efecto de antibiosis frente a otros microorganismos, cuyo acceso a este nutriente se verá reducido significativamente, así como su capacidad de proliferar. Se ha demostrado que la producción de pioverdinas, un sideróforo fluorescente sintetizado por algunas Pseudomonas (Schalk, Guillon, 2013), ayuda a controlar los efectos de los patógenos F. oxysporum en patata y Gaeumannomyces graminis en tabaco (Voisard et al., 1989; Schippers et al., 1990). Asimismo, algunas rizobacterias, a medida que colonizan la raíz, inducen en las plantas un mecanismo de resistencia llamado ISR (Induced Systemic Resistance), mediante el que la parte aérea de las mismas se protege contra un amplio rango de patógenos. Esta propiedad se describió inicialmente al observarse que la cepa Pseudomonas fluorescens WCS417r protegía sistemáticamente plantas de clavel contra F. oxysporum f. sp. dianthi (Peer van et al., 1991), así como una selección de rizobacterias inducían resistencia en pepino contra Colletotrichum orbiculare (Wei, 1991). Se han descrito diversos determinantes bacterianos capaces de desencadenar ISR, como la flagelina (Ausubel, 2005), el antígeno O de los lipopolisacáridos (Leeman et al., 1995), algunos sideróforos (Ran et al., 2005; Loon van et al., 2008) y ciertos lipopéptidos. Tal es el caso de los lipopéptidos surfactina y fengicina, sintetizados por Bacillus subtilis, que activan este mecanismo en plantas de judía contra el patógeno Botrytis cinerea (Ongena et al., 2007), agente causal de la podredumbre gris. Dada su relevancia en esta Tesis Doctoral, el siguiente apartado se centra en una cepa bacteriana probada como agente de control biológico en plantas de olivo, P. fluorescens PICF7. Pseudomonas fluorescens PICF7, agente de biocontrol de la verticilosis del olivo P. fluorescens PICF7 es una cepa endofita natural de raíces de olivo y eficaz agente de biocontrol de la verticilosis (Verticillium dahliae Kleb.), una de las enfermedades más devastadoras que afecta al cultivo de dicha leñosa en la cuenca mediterránea. El hongo penetra en el olivo por heridas o a través de las raíces y llega al xilema, desde el que se expande al resto de la planta. Cuando la enfermedad está muy avanzada, el hongo crece 28 Introducción General fuera de los tejidos vasculares, dando lugar a una flacidez diurna y una defoliación intensa de la planta, que acaba por marchitarse de forma permanente hasta llegar en muchos casos a la muerte (López-Escudero, Mercado-Blanco, 2011). P. fluorescens PICF7 se aisló del suelo de olivares junto con otras siete Pseudomonas en un estudio llevado a cabo por el Dr. Jesús Mercado-Blanco y sus colaboradores (MercadoBlanco et al., 2004). En dicho trabajo, todos los aislados presentaron actividad supresora contra V. dahliae, siendo la cepa PICF7 la que exhibió mayor capacidad de biocontrol. Estudios posteriores han confirmado que PICF7 es capaz de controlar eficientemente infecciones provocadas por el patotipo defoliante (D) de V. dahliae en plantones de la variedad Picual, que presenta una alta susceptibilidad a la enfermedad. Asimismo, mediante microscopía confocal, disección tridimensional de tejidos y marcaje autofluorescente, se demostró la capacidad de esta cepa de colonizar endofíticamente raíces de olivos de la variedad Arbequina (Prieto, Mercado-Blanco, 2008). Sin embargo, los determinantes bacterianos de PICF7 relacionados tanto con su estilo de vida endofítico como su actividad de biocontrol están aún por dilucidarse. Una primera aproximación se ha centrado en estudiar, mediante mutagénesis dirigida, la implicación de fenotipos como la motilidad tipo “swimming” o la producción del sideróforo pioverdina en las propiedades naturales de PICF7. Mutantes en ambos fenotipos no mostraron diferencias significativas en la colonización de las raíces ni el antagonismo de PICF7 contra V. dahliae (MaldonadoGonzález et al., 2013). Por otro lado, la planta parece responder a la colonización radicular de PICF7, lo que podría explicar su capacidad de biocontrol. Mediante la metodología Suppression Subtractive Hybridization (SSH), se pudo confirmar que PICF7 induce ISR en los tejidos aéreos (Gómez-Lama et al., 2014), así como un amplio espectro de respuestas defensivas en las propias raíces (Schilirò et al., 2012). Con objeto de comprobar el potencial de PICF7 como agente de biocontrol contra otros patógenos de olivo, un estudio reciente examinó la interacción entre esta cepa y P. savastanoi NCPPB 3335, cepa modelo de estudio de infección bacteriana en plantas leñosas (Ramos et al., 2012). El estudio mostró que, aunque PICF7 no es capaz de suprimir el desarrollo de tumores, su co-inoculación con el patógeno produce una disminución de la población de este último in planta, dando lugar a síntomas necróticos reducidos y a una alteración de la colonización del tejido hiperplásico (Maldonado-González et al., 2013). La secuencia del genoma de PICF7 se ha caracterizado recientemente, y el estudio de la misma ha sido abordado en esta Tesis Doctoral. Como consecuencia, se han identificado una serie de componentes genéticos potencialmente implicados en las propiedades fenotípicas de este agente de biocontrol. 29 Introducción General 3. Bioinformática para el análisis genómico de bacterias asociadas a plantas Como se ha señalado, la evolución de las metodologías de secuenciación genómica ha sido notable en los últimos años (Zhao, Grant, 2011). Los avances en estas técnicas han supuesto tal reducción en tiempo y costes que han cambiado el paradigma de la investigación biológica (Delanty, Goldstein, 2013; Salto-Tellez, Gonzalez De Castro, 2014). La facilidad con que se puede acceder al genoma de prácticamente cualquier organismo vivo ha hecho que especialistas de todas las áreas secuencien, casi rutinariamente, los genomas de aquellos organismos que les son de interés (Rothberg, Leamon, 2008). La fitobacteriología no es una excepción, y en la actualidad existe un repertorio extenso de genomas de bacterias asociadas a plantas accesibles en bases de datos públicas (Hamilton et al., 2011; Winsor et al., 2011). Los mecanismos que dictan el estilo de vida de estas bacterias responden en gran medida a su información genética, por lo que profundizar en el conocimiento de los mismos pasa necesariamente por el estudio de sus genomas. Aunque el desarrollo de la secuenciación masiva es relativamente reciente (Margulies et al., 2005), los bacteriólogos llevan décadas caracterizando genes involucrados en las interacciones que las bacterias establecen con sus plantas huésped. Los procedimientos de secuenciación usados, rudimentarios en comparación con los actuales, han permitido caracterizar un amplio repertorio de componentes genéticos bacterianos, que han ayudado en gran medida a describir mecanismos específicos de patogénesis, de control biológico, de promoción del crecimiento vegetal, etc. En estos procesos de asociación intervienen multitud de factores, muchos de ellos esenciales para que la interacción surta efecto. Factores como adhesinas o exopolisacáridos son fundamentales en la fase de adherencia a los tejidos de la planta (Pizarro-Cerdá, Cossart, 2006; Mhedbi-Hajri et al., 2011). Asimismo, para subsistir en dichos tejidos muchas bacterias disponen de sideróforos, antibióticos o bombas de extrusión multidroga, necesarios para competir por los nutrientes y protegerse de compuestos tóxicos antimicrobianos (Ahmed, Holmström, 2014; Wang, Raaijmakers, 2004; Martínez et al., 2009). Por otro lado, factores como fitotoxinas, fitohormonas o enzimas degradadoras de la pared celular vegetal juegan un papel esencial en la colonización del huésped (Dudler, 2014; Costacurta, Vanderleyden, 1995; Barras et al., 1994). La determinación de las secuencias genéticas encargadas de la síntesis de estos factores ha contribuido de forma significativa al estudio de las interacciones planta-bacteria, y con los avances en las técnicas de secuenciación, la identificación de nuevos factores se ha facilitado notablemente. 30 Introducción General Algunos de los factores anteriormente mencionados son transportados al interior del huésped, para lo que la bacteria dispone de herramientas específicas encargadas de la secreción de moléculas al exterior. Tal es el caso de los sistemas de secreción, que permiten a las bacterias transportar proteínas al medio extracelular, así como al interior de otras células eucariotas o procariotas (Tseng et al., 2009). En el contexto de las interacciones planta-bacteria, son de especial relevancia los sistemas de secreción de tipos III, IV y VI (T3SS, T4SS y T6SS, respectivamente). Algunos ejemplos del papel que juegan tanto el T3SS como el T4SS en las asociaciones que establecen ciertas bacterias con sus plantas huésped se mencionaron en el apartado 2.1. Por otro lado, el T6SS es un sistema de reciente descripción, cuyos componentes se asemejan estructuralmente al bacteriófago T4 (Silverman et al., 2012). Como el T3SS y el T4SS, el T6SS transloca proteínas al exterior, y aunque está presente en una proporción considerable de bacterias Gram-negativas (Bingle et al., 2008; Boyer et al., 2009), existe todavía poca evidencia experimental que describa el modo de ensamblaje de este sistema, así como las funciones específicas de las proteínas secretadas. Estudios funcionales de sus efectores han mostrado que el T6SS no sólo juega un papel importante en la virulencia de ciertas bacterias, sino también en las relaciones mutualistas y comensalistas que éstas establecen con determinados organismos eucariotas (Jani, Cotter, 2010). Asimismo, trabajos recientes han indicado que efectores del T6SS tienen función antibiótica frente a otros microorganismos, ayudando en el proceso de competencia de algunas bacterias (Dong et al., 2013). Dada la relevancia que tienen los sistemas de secreción en el estilo de vida de las bacterias, se han caracterizado también muchos de los genes implicados en la síntesis tanto de sus componentes estructurales como de sus efectores asociados. En suma, y aunque algunos de los mecanismos involucrados en las interacciones planta-bacteria están aún por dilucisarse, el número de genes bacterianos conocidos que intervienen de una forma u otra no es nada desdeñable. Si consideamos además el conjunto creciente de genomas bacterianos disponibles, la cantidad de información resultante es inmensa. Desde un punto de vista computacional, dicha información resulta muy útil, en tanto que cubre un amplio espectro del conocimiento existente de las interacciones planta-bacteria, y puede usarse para generar nuevas hipótesis en este campo. Una aplicación directa de estos datos es el desarrollo de herramientas específicas de anotación, capaces de escanear secuencias de entrada en busca de genes homólogos a los ya conocidos y descritos en bacterias asociadas a plantas. La anotación generada podría usarse para direccionar posteriores análisis experimentales en base a los genes identificados. El desarrollo de este tipo de herramientas ha sido abordado en esta Tesis Doctoral. 31 Introducción General Otro aspecto interesante a la hora de tratar la ingente cantidad de información disponible es que ésta ofrece la oportunidad de abordar la genómica de las interacciones planta-bacteria desde un punto de vista global. Si el modo de vida de una bacteria está de alguna forma codificado en su material genético, sería interesante comprobar hasta qué punto la información genómica disponible es suficiente para identificar las distintas asociaciones que establecen las bacterias con las plantas. En este contexto, los métodos de aprendizaje automático son particularmente útiles, ya que, como se mencionó en el apartado 1.4, permiten usar la información existente para inducir modelos de representación. De este modo, podríamos entrenar modelos en base a los datos genómicos conocidos para generar clasificadores que nos permitan agrupar genomas de bacterias con estilos de vida similares, como patógenas o comensalistas. Este tipo de aproximación ya ha sido aplicada a genomas de bacterias asociadas a humanos con resultados prometedores, aunque el límite parece estar en discernir bacterias patógenas y no patógenas (Iraola et al., 2012; Andreatta et al., 2010; Barbosa et al., 2014). En lo que respecta a las interacciones planta-bacteria, una herramienta de predicción de estilos de vida basada en información genómica sería de gran utilidad, no sólo en el ámbito de la fitopatología, sino también en el de la seguridad alimentaria. En la última década, se han dado diversas epidemias relacionadas con el consumo de productos vegetales, afectando a miles de personas en todo el mundo. Sorprendentemente, los patógenos responsables han sido cepas de bacterias tradicionalmente asociadas al consumo de alimentos de origen animal, como Escherichia coli oSalmonella (Deering et al., 2012; Lim et al., 2014). Teniendo en cuenta este escenario, en el presente trabajo se ha abordado el diseño de un método que, en base a la información genómica disponible, identifica interacciones potenciales entre bacterias y plantas. 32 OBJETIVOS Objetivos Esta Tesis Doctoral tiene un doble cometido: por un lado, analizar las secuencias genómicas de dos cepas bacterianas asociadas a cultivos de relevancia en la agricultura andaluza, haciendo hincapié en los factores genéticos implicados en la asociación de las mismas con sus plantas huésped; por otro, desarrollar nuevas herramientas computacionales basadas en la información existente de las interacciones planta-bacteria que permitan extraer nuevo conocimiento en este campo. Para ello se establecieron los siguientes objetivos específicos: 1. Analizar el genoma secuenciado de la cepa Pseudomonas syringae pv. syringae UMAF0158, agente causal de la necrosis apical del mango, poniendo especial énfasis en genes potencialmente implicados en la síntesis de factores de virulencia. 2. Analizar el genoma secuenciado de la cepa Pseudomonas fluorescens PICF7, endofita de olivo y agente de biocontrol contra la verticilosis, destacando genes potencialmente implicados en biocontrol y endofitismo. 3. Desarrollar herramientas bioinformáticas de anotación de genes bacterianos implicados en las interacciones planta-bacteria. 4. Combinar dichas herremientas con técnicas de aprendizaje automático para implementar un sistema automatizado de clasificación de genomas bacterianos que permita identificar potenciales asociaciones planta-bacteria. 35 Chapter I but they grouped into different pv. syringae phylotypes (Gutiérrez-Barranquero et al., 2013a) and clades of phylogroup 2 (Berge et al., 2014). Our analysis reveals the existence of genetic differences between these two strains, which may confer UMAF0158 its ability to infect mango trees. In addition to the presence of the mangotoxin biosynthetic operon mbo, the UMAF0158 genome differs from that of B728a in the codification of a cellulose synthase operon and in harboring two additional secretion systems, i.e., a T3SS and a T6SS. Moreover, UMAF0158 displays a different repertoire of T3Es, which may be a determinant of its association with mango trees. Our data provide the basis for further functional studies of the virulence mechanisms and host specificity determinants of the P. syringae pv. syringae UMAF1058, a representative strain of phylotype 1 of this pathovar, which is the fourth P. syringae strain whose complete genome sequence has been made available. 43 Chapter I Results and Discussion General features The P. syringae pv. syringae UMAF0158 genome is composed of one circular chromosome of 5,787,986 bp (Table 1; Figure 1) and one plasmid. The pPSS158 plasmid of 63,004 bp has an average GC content of 54.6% and 71 predicted coding sequences (CDSs) (Table 1). Among the latter, repA, a T4SS conjugative system and the rulAB genes appear as the most relevant features (Table S1). Conjugative plasmids harboring rulAB genes have been shown to contribute to UV and solar radiations tolerance, as well as to epiphytic fitness (Cazorla et al., 2008). In total, 5017 CDSs were identified within the UMAF0158 chromosome, which has an average GC content of 59.3% (Table 1). Among the predicted chromosomal CDSs, a putative function was assigned to 4030 (80%), while the remaining 987 CDSs were designated as hypothetical proteins. A total of 18 genes were predicted to be pseudogenes. The classification of the UMAF0158 CDSs into functional categories according to the COG (Clusters of Orthologous Groups) database is summarized in Table 2 in comparison with P. syringae pv. syringae B728a, P. syringae pv. tomato DC3000 and P. syringae pv. phaseolicola 1448A. With the exception of category L (replication, recombination and repair), which includes a reduced number of UMAF0158 CDSs (137) in comparison with those of the other three genomes (187, 272 and 370 CDSs for B728a, 1448A and DC3000, respectively), no significant differences were found regarding the remaining functional categories, further supporting the relatedness of these four strains. In addition to CDSs, a total of 63 tRNAs and 5 rRNA operons were found on the UMAF0158 chromosome. Phylogeny In order to establish the phylogenetic relationship between UMAF0158 and other related P. syringae strains, we selected 25 genome sequenced strains (Table S2) and compared a set of five protein-coding house-keeping genes, namely gapA,gltA,recA,rpoA and rpoB. We created an alignment of the proteins and reconstructed the phylogenetic tree shown in Figure 2, using neighbor-joining methods. The strain Pseudomonas fluorescens Pf-5 was used as an out-group. The resulting phylogeny clustered UMAF0158 with P. syringae Cit 7, a strain originally isolated from a healthy orange tree (Lindow, 1985; Baltrus et al., 2011), and more separate from B728a, the model strain for pv. syringae. This result is in agreement with previous phylogenetic analyses, which clustered both UMAF0158 and Cit 7 in phylotype 1 (Carrión et al., 2013; Gutiérrez-Barranquero et al., 2013a). 44 Chapter I Such a phylotype of the pathovar syringae is mainly associated with the mango host and characterized by mangotoxin production. Additionally, UMAF0158 and other strains of phylotype 1 are pathogenic on mango, lilac, tomato or pear (Cazorla et al., 1998), but not on bean; in contrast, B728a is pathogenic on bean but not on mango (Feil et al., 2005; Gutiérrez-Barranquero et al., 2013a). Given that UMAF0158 is the only strain belonging to this group whose complete genome is available, it could be taken as a representative of phylotype 1 for pv. syringae. Figure 1. Features of the P. syringae pv. syringae UMAF0158 chromosome. From the outside in, the outermost circle (black) shows the scale line; circles 2 and 3 represent predicted coding regions on the plus and minus strand, respectively, which are color coded based on COG categories; circles 4 and 5 show tRNA (blue) and rRNA (red), respectively; circle 6 depicts ORFs associated with virulence (see Virulence factors section and Table S6). 45 Table 1. General features of the P. syringae pv. syringae UMAF0158 genome and comparison with P. syringae pv. syringae B728a, P. syringae pv. phaseolicola 1448A and P. syringae pv. tomato DC3000. UMAF0158 B728a 1448A DC3000 Molecule Chromosome pPSS158 Chromosome Chromosome p1448A-A p1448A-B Chromosome pDC3000A pDC3000B Size (bp) 5,787,986 63,004 6,093,698 5,928,787 131,950 51,711 6,397,126 73,661 67,473 G+C content (%) 59.3 54.6 59.2 58 54.1 56 58.4 55.1 56.1 No. of predicted CDSs 5,017 71 5,089* 4,985* 127* 60* 5,481* 68* 70* No. of rRNAs 16 - 16 16 - - 15 - - No. of tRNAs 63 - 64 64 - - 63 - - Reference This study Feil et al. (2005) Joardar et al. (2005) Buell et al. (2003) * The number of predicted CDSs corresponds to those indicated at NCBI for the corresponding genome sequences (March 1st, 2015). Chapter I Table 2. Number of CDSs associated with COG functional categories in the P. syringae pv. syringae UMAF0158 genome and comparison with P. syringae pv. syringae B728a, P. syringae pv. phaseolicola 1448A and P. syringae pv. tomato DC3000. Functional category UMAF0158 B728a 1448A DC3000 A RNA processing and modification 1 1 1 1 B Chromatin structure 1 1 1 1 C Energy production and conversion 223 231 210 234 D Cell cycle control, cell division, chromosome partitioning 38 46 41 44 E Amino acid transport and metabolism 459 467 458 461 F Nucleotide transport and metabolism 89 87 88 81 G Carbohydrate transport and metabolism 276 268 263 264 H Coenzyme transport and metabolism 177 178 179 173 I Lipid transport and metabolism 163 162 166 176 J Translation, ribosomal structure and biogenesis 198 205 204 201 K Transcription 360 360 343 367 L Replication, recombination and repair 137 187 272 370 M Cell wall/membrane/envelope biogenesis 270 288 270 263 N Cell motility 163 166 168 160 O Posttranslational modification, protein turnover, chaperones 156 157 152 157 P Inorganic ion transport and metabolism 272 277 287 278 Q Secondary metabolites biosynthesis, transport, and catabolism 119 128 112 118 R General function prediction only 521 522 520 544 S Function unknown 385 393 356 406 T Signal transduction mechanisms 343 349 339 358 U Intracellular trafficking, secretion, and vesicular transport 150 148 156 136 V Defense mechanisms 48 53 48 58 Z Cytoskeleton 1 1 1 0 Total 4,550 4,675 4,635 4,851 47 Chapter I Regarding the other P. syringae strains in the phylogeny with complete genome sequences, B728a was the closest to UMAF0158. Accordingly, this strain shares the highest number of CDSs predicted in UMAF0158 (see next section). These data, together with the fact that both B728a and UMAF0158 belong to P. syringae pv. syringae, prompted us to pay special attention to the genomic differences between these two strains. Figure 2. Phylogenetic analysis of P. syringae pv. syringae UMAF0158 and 25 selected strains of the P. syringae complex. Multilocus sequence analysis were performed using a concatenated dataset for gapA,gltA,recA,rpoA and rpoB. The evolutionary history was inferred using the Maximum Likelihood method based on the JTT matrix-based model. The percentage of trees in which the associated taxa clustered in the bootstrap test (1000 replicates) is shown next to the branches. Some strains are included in phylotypes 1 and 3 of pv. syringae according to Carrión et al. (2013) and Gutiérrez-Barranquero et al. (2013a). 48 Chapter I Comparative genomics The sequence of the UMAF0158 chromosome was compared to that of selected P. syringae strains (Figure 3). Of the 5017 CDSs predicted in UMAF0158, 4912 (98%) have orthologs (BLASTP E-value ≤1e−10) in other P. syringae and 3570 (71%) are present in all strains. All of the 105 genes found to be unique to UMAF0158 are heavily enriched in hypothetical proteins (103). The other two genes include a membrane protein (PSYRMG_17725) and a flavodoxin (PSYRMG_09680). Among the selected strains, Cit 7, BRIP39023 and 642 have the closest number of CDSs compared to UMAF0158. These three strains, which belong to pv. syringae or are closed to it and whose complete genome sequences are not yet available, share 93.6, 93.4 and 91.2% of the CDSs predicted in UMAF0158, respectively (Figure 3). Regarding P. syringae with complete genome sequences, B728a shares 90.7% of the CDSs followed by 1448A and DC3000, which share 88.4 and 88.1%, respectively (Figure 3). Figure 3 shows some sequence features associated with mechanisms of horizontal transfer, including regions with differential distributions of trinucleotides and GC-content, predicted prophages and putative horizontally transferred genes. In most cases, these features match with non-conserved regions of the UMAF0158 chromosome (white-colored in the six most outer rings). Comparison of CDSs between UMAF0158 and B728a In order to compare the genomes of UMAF0158 and B728a, we proceeded to identify regions enriched in coding genes in either strain that are not present in the other. The search was performed so that only regions spanning at least 4000 kb were retained (the whole set of differential protein coding genes are listed in Tables S3 and S4). These regions are summarized in Tables 3 and 4. Thirteen regions were identified in UMAF0158 with sizes ranging from 6051 to 20822 kb. Eight of such regions are highly enriched in hypothetical proteins (at least 80% of their CDSs). Two of the remaining five regions contain a combination of mobile genetic elements and hypothetical proteins. The remaining regions correspond to three potential operons: an additional T3SS, a cellulose production operon (PSYRMG_20805-20845), and the well-described mangotoxin biosynthetic operon mbo (PSYRMG_10110-10135) (Carrión et al., 2012). This operon is present in only a limited number of strains belonging to genomospecies 1, and it has been acquired once during evolution by horizontal transfer (Carrión et al., 2013). In addition, ten regions were identified in B728a with sizes ranging from 4968 to 43402 kb. Most of these regions contain mobile genetic elements. It is worth noting a region containing the streptomycin 49 Chapter I resistance transposon Tn5393 (Feil et al., 2005). Two other regions are enriched in secretion components with one of them corresponding to a T4SS, which is addressed in the next section. Figure 3. Conservation analysis of the P. syringae pv. syringae UMAF0158 chromosome. From the outside in, the outermost circle (black) shows the scale line. Circles 2 to 4 display CDSs homology (E-value ≤1e−10) among UMAF0158 and the three P. syringae with complete genome sequences: DC3000 (grey), 1448A (orange) and B728a (red). Circles 5 to 7 display CDSs homology (E-value ≤1e−10) among UMAF0158 and the draft genomes of the three phylogenetically closest P. syringae strains among the 25 selected in this study: 642 (purple), BRIP39023 (green) and Cit 7 (blue). Circles 8 and 9 display putative horizontally transferred regions (red) and prophages (purple), respectively; circle 10 shows G+C in relation to the mean G+C in 2 kb windows (red); circle 11 shows trinucleotide composition (black). 50 Table 3. Regions of P. syringae syringae pv. syringae UMAF0158 chromosome with low homology to P. syringae pv. syringae B728a. Location (bp) Length No. of CDSs No. of hypothetical No. of CDSs not present in B728a Relevant features 248304-257525 9221 kb 16 9 (56%) 15 (94%) mobile genetic elements 258593-272300 13707 kb 21 11 (52%) 19 (90%) mobile genetic elements 519318-537403 18085 kb 20 8 (40%) 15 (75%) T3SS components 996632-1011913 15281 kb 17 17 (100%) 17 (100%) hypothetical proteins 1133844-1143487 9643 kb 20 19 (95%) 18 (90%) mobile genetic elements 1357822-1378644 20822 kb 5 4 (80%) 5 (100%) hemolysin secretion/activation 2201422-2221045 19623 kb 22 18 (82%) 19 (86%) mobile genetic elements 2233554-2244999 11445 kb 33 28 (85%) 31 (94%) thiamin biosynthesis, peptidase 2324538-2330589 6051 kb 6 2 (33%) 6 (100%) mangotoxin biosynthetic operon 2710838-2721328 10490 kb 10 10 (100%) 10 (100%) hypothetical proteins 3032080-3050678 18598 kb 14 13 (93%) 14 (100%) hypothetical proteins 4684145-4696189 12044 kb 7 0 6 (86%) cellulose synthase 5668482-5685647 17165 kb 21 21 (100%) 21 (100%) hypothetical proteins Table 4. Regions of P. syringae syringae pv. syringae B728a chromosome with low homology to P. syringae pv. syringae UMAF0158. Location (bp) Length No. of CDSs No. of hypothetical No. of CDSs not present in UMAF0158 Relevant features* 248304-257525 9221 kb 16 9 (56%) 15 (94%) HP/MGE 258593-272300 13707 kb 21 11 (52%) 19 (90%) HP/MGE T3 effector 519318-537403 18085 kb 20 8 (40%) 15 (75%) T4SS components 996632-1011913 15281 kb 17 17 (100%) 17 (100%) HP membrane transport 1133844-1143487 9643 kb 20 19 (95%) 18 (90%) MGE secretion pilus proteins 1357822-1378644 20822 kb 5 4 (80%) 5 (100%) HP/Virulence protein 2201422-2221045 19623 kb 22 18 (82%) 19 (86%) streptomycin resistance 2233554-2244999 11445 kb 33 28 (85%) 31 (94%) HP/phage-related proteins 2324538-2330589 6051 kb 6 2 (33%) 6 (100%) HP/MGE/T3 effector 2710838-2721328 10490 kb 10 10 (100%) 10 (100%) HP membrane transport plasmid-related proteins phage-related proteins T3 effector * HP and MGE refer to hypotetical proteins and mobile genetic elements, respectively. Chapter I Arrebola et al., 2003). UMAF0158 contains orthologs of genes participating in the synthesis of syringopeptin and syringomycin. These two toxins induce necrosis in plant tissues and have been shown to be the major virulence determinants of P. syringae pv. syringae (Bender et al., 1999; Scholz-Schroeder et al., 2001). The two clusters encoding these toxins form a larger cluster (PSYRMG_03860 – 03910), which is consistent with previously reported data (Feil et al., 2005; Scholz-Schroeder et al., 2001). The production of syringomycin by UMAF0158 has been experimentally validated in previous studies by our group (Arrebola et al., 2003). Syringolins are another family of phytotoxins synthesized by a number of P. syringae pv. syringae (Krahn et al., 2011). UMAF0158 has a gene cluster resembling that of the production of syringolin A (PSYRMG_24250-24275), a toxin that has been shown to counteract stomatal innate immunity in beans and Arabidopsis (Schellenberg et al., 2010). Phaseolotoxin and coronatine are two chlorosis-inducing toxins that also represent major virulence factors for some P. syringae isolates (Aguilera et al., 2007; Zheng et al., 2012). UMAF0158 lacks orthologs for most of the genes involved in the production of coronatine; however, analysis of the 23 genes required for the synthesis of phaseolotoxin (Aguilera et al., 2007) showed that orthologs of 17 of these genes are included in its genome. Given that no inhibition halos were observed in the bioassay for toxin detection when ornithine was added (Arrebola et al., 2003) and no chlorosis was detected among the symptoms of UMAF0158 infection, it is likely that the lack of the other six genes prevent the synthesis of phaseolotoxin by this strain. The production of mangotoxin by UMAF0158 and its contribution to the virulence of this strain has been widely described (Arrebola et al., 2007; Carrión et al., 2012). Two operons are involved in the synthesis of this toxin, namely, mgo (PSYRMG_1582015835) and mbo (PSYRMG_10110-10135) (Carrión et al., 2012; Arrebola et al., 2012). The latter has been shown to be absent in B728a (Carrión et al., 2013). Unlike other toxins, mangotoxin is thought to be associated with host specificity, as it has been found to be synthesized by strains of phylotypes 1 and 2 of pv. syringae, which were mainly isolated from mango trees and other woody crops (Gutiérrez-Barranquero et al., 2013a). It is worth noting that, even though the metabolic cost of toxin production generally prohibits bacteria from producing more than one, the synthesis of at least two of them has been experimentally validated in UMAF0158 (i.e., syringomicin and mangotoxin; Arrebola et al., 2003). Whether this strain is capable of producing the rest of the phytotoxins mentioned above is a question that requires further investigation. 59 Chapter I Phytohormones Bacterial-produced phytohormones are typically transported to the plant cell to regulate plant biological processes, providing a beneficial context for the pathogen (Costacurta, Vanderleyden, 1995). Such is the case of auxin, which is predominantly represented by indole-3-acetic acid (IAA), a key plant growth regulator that is also involved in plant-bacteria interactions. The downregulation of this hormone in plants has been shown to restrict P. syringae growth in Arabidopsis (Navarro et al., 2006), suggesting that bacteria may have evolved the production of auxin to overcome this plant response. Accordingly, auxin production has been demonstrated to promote susceptibility to P. syringae (Mutka et al., 2013). The genome of UMAF0158 contains orthologs of two of the genes involved in the biosynthesis of IAA, namely, iaaH and iaaM. Further analyses will elucidate whether this strain is able to produce IAA. Detoxifying compounds As a response to bacterial infections, plants synthesize reactive oxygen species (ROS), such as hydrogen peroxide (H2O2), superoxide (O2-) and hydroxyl radical (OH). These molecules have a toxic effect on invading bacteria (Cabiscol et al., 2000), which have evolved mechanisms to counterattack by detoxification. DC3000 typically makes use of the catalases KatB and KatE together with the catalase-peroxidase KatG to detoxify plantproduced H2O2. Interestingly, orthologs of these three gene products are present in the UMAF1058 genome. Moreover, an ortholog of the gene that encodes Dps, a ferritinlike protein that has been shown to protect plant-associated bacteria against oxidative stress (Colburn-Clifford et al., 2010), has been also identified (PSYRMG_23380). Other genes possibly implicated in detoxification found in UMAF0158 correspond to a cluster presumably encoding a cbb(3)-type cytochrome C oxidase (PSYRMG_07490-07505) and a proline iminopeptidase (PSYRMG_18265). The latter has been reported to be required for pathogenicity of Xanthomonas campestris (Zhang et al., 2007) and to have dealanylating activity toward ascomycin, an antibiotic produced by Streptomyces that inhibits protein synthesis (Sudo et al., 1996). Genes involved in copper resistance have also been identified, including a cluster containing copA and copB (PSYRMG_23630 and PSYRMG_23625, respectively), and a locus with high similarity to the cueAR system in Pseudomonas putida (Adaikkalam, Swarup, 2005). This system consists of a copper-transporting ATPase transmembrane protein and its transcriptional regulator (PSYRMG_19635 and PSYRMG_19630, respectively). Additionally, the copABCD operon described in other P. syringae (Cazorla et al., 2002) is absent in the UMAF0158 chromosome, and also any 60 Chapter I copper resistant genes are present in the UMAF0158 plasmid, in agreement with the copper sensitivity of this strain (Gutiérrez-Barranquero et al., 2013b). PCWDEs Some phytopathogenic bacteria need to overcome the plant cell wall in the process of accessing the host cytoplasm. Therefore, many plant pathogens harbor a collection of genes encoding PCWDEs, which are considered important virulence determinants (Barras et al., 1994). The UMAF1058 genome contains genes predicting several PCWDEs, such as a cellulase (PSYRMG_06950), a lipoyl synthase (PSYRMG_12655), a xylanase (PSYRMG_13355) and a pectin lyase (PSYRMG_10750). 61 Chapter I Conclusions Summarizing, bioinformatics analysis of the complete genome of P. syringae pv. syringae UMAF0158, a pathogen of mango trees, revealed a high degree of conservation with other pseudomonads belonging to the P. syringae complex, including the model strain P. syringae pv. syringae B728a. However, the resulted phylogeny clustered UMAF0158 with P. syringae Cit 7 and more separately from B728a. Indeed, our data revealed a number of genetic factors that could be involved in the differential pathogenic and epiphytic lifestyle of UMAF0158, in comparison with the model strain B728a. The mangotoxin biosynthetic operon mbo is included among these factors, whose role in the pathogenicity of UMAF0158 has been previously reported (Carrión et al., 2012, 2013). Moreover, UMAF0158 harbors an operon presumably involved in cellulose production and two clusters predicting additional secretion systems of types III and IV, as well as displays a particular T3Es repertoire. Additionally, the conjugative plasmid pPSS158 contains rulAB genes, involved in UV resistance and epiphytic fitness (Cazorla et al., 2008). This work provides the basis for further analyses on the specific mechanisms that enable this strain to infect mango trees as well as for functional analyses of the factors governing host specificity in pv. syringae strains from different phylotypes. 62 Chapter I Material and Methods Bacterial growth and DNA methods The bacterial strain UMAF0158 (CECT 7752) belonging to Pseudomonas syringae pv. syringae was routinely grown in KB medium with 48 h of incubation at 28ºC. UMAF0158 was inoculated on 100 ml of LB medium and grown for 15 h at 28ºC with shaking (150 r.p.m.). After this period, the OD600nm of the culture was 1.8. Serial dilutions of this culture and plating on LB plates yielded 1.3×109 CFU/ml of pure bacterial culture. The rest of the culture was divided into 54 aliquots of 1.5 ml, and DNA was extracted from all of the cultures using the Jet-Flex genomic DNA purification kit (Genomed GmbH, Germany). DNA samples were collected together and further purified by extraction with 1:1 phenol:chloroform and precipitated with 4 M NaCl and 13% PEG. DNA was suspended in 500 µl MilliQ H2O, and NanoDrop measurements indicated 3.6 µg/µl (in total 1800 µg of DNA) with an A260/A280 of 1.85. The extracted DNA was visualized in 1% agarose after digestion with the restriction enzymes EcoRI and PstI. Whole Genome Sequencing The finished UMAF1058 genome was generated at the Beijing Genomics Institute (BGI-HK) using an Illumina HiSeq 2500 system. Briefly, the isolated DNA was used to generate three libraries of 500 bp, 2000 bp and 6000 bp, producing 1336, 1312 and 1352 Mb of raw data, respectively. These were passed through a filtering pipeline that removed known sequencing and library preparation artifacts. After data treatment, the SOAPdenovo 1.05 software package (Li et al., 2008, 2010) was used for sequence assembly and quality assessment. Assembly results were then combined and mapped to the genome of P. syringae B728a, yielding to two scaffolds corresponding to one chromosome and one plasmid. Finally, a PCR gap closure and three circle PCR verification were performed to obtain the final complete sequences of the chromosome and the pPSS158 plasmid. Genomic data and annotation The assembled genome of UMAF0158 was submitted to the NCBI Prokaryotic Genome Annotation Pipeline for automatic annotation and manually reviewed. Gene locations and protein products were generated from above annotation (ASN.1 file) using the script “asn2all” belonging to the NCBI ToolKit (http://www.ncbi.nlm.nih.gov/toolkit). Genome sequences (DNA, proteins and predicted genes locations) of DC3000, B728a and 1488A were downloaded from the NCBI respository of completely sequenced bacterial 63 Chapter I genomes (ftp://ftp.ncbi.nlm.nih.gov/genomes/Bacteria), while the corresponding sequences of P. savastanoi pv. savastanoi strain NCPPB3335 were obtained from our own sequencing project. The rest of the P. syringae genomes were downloaded from the NCBI draft bacterial genome repository (ftp://ftp.ncbi.nlm.nih.gov/genomes/Bacteria_DRAFT/). In this case, proteins and gene locations were generated from DNA sequences (fna files) using Glimmer v3.02 (Delcher et al., 2007). All genomes were downloaded on July 15, 2014. Accession numbers and references for all genome sequences used in this work are summarized in Table S2. The UMAF0158 annotation of COGs was performed by aligning the set of predicted protein sequences against the COG PSSM of the CDD (Marchler-Bauer et al., 2014) using rps-BLAST. Hits with an E-value ≤0.001 were first retained. Then, only the best hit was selected for each protein. The same procedure was used to assign COG categories to the repertoire of predicted proteins of the B728a, DC3000 and 1448A strains. Predictions of horizontally transferred regions and prophages were computed using Alien Hunter 1.7 (Vernikos, Parkhill, 2006) and Prophage Finder (Bose, Barber, 2006), respectively. T346Hunter (Martínez-García et al., 2015b) was used to identify secretion systems clusters. Virulence factors were predicted using a customized annotation tool developed by our group (not published). Further details on such a tool are provided in Chapter IV. The finished genomic sequences of UMAF0158 were deposited in GenBank under accession numbers CP005970 (chromosome) and CP005971 (plasmid pPSS158) Trinucleotide composition The distribution of all 64 trinucleotides was determined for the whole chromosome and 2 kb sub-windows. Then, the χ2statistic of the difference between the trinucleotide composition of each window and that of the whole chromosome was computed. Large values for χ2in a given window denote different trinucleotide compositions from the rest of the chromosome. Probability values were computed assuming uniform distribution of the DNA composition along the genome. Because this assumption may have been incorrect, large χ2values were interpreted as indicators of unusual regions on the chromosome that require further investigation. Phylogenetic tree Phylogeny was determined by multilocus sequence analysis using a concatenated dataset for the housekeeping genes gapA,gltA,recA,rpoA and rpoB. Then, the Maximum Likelihood method based on the JTT (Jones-Taylor-Thornton) matrix-based model (Jones et al., 1992) was applied. The percentage of trees in which associated taxa clustered in the 64 Chapter I bootstrap test (1000 replicates) is shown next to the branches in Figure 2 (Felsenstein, 1985). Multiple alignments and evolutionary analyses were conducted using MEGA5 software (Tamura et al., 2011). Comparative genomics Each of the predicted proteins in UMAF0158 was compared to those of the other P. syringae strains using BLASTP (E-value ≤1e−10). The same procedure was used to compare predicted proteins in B728a with those in UMAF0158. Distribution of T3Es We performed BLASTP searches of the T3Es in http://Pseudomonas-syringae.org/ against the 26 P. syringae predicted proteomes. First, only hits with an E-value ≤1e−10 were retained. If no hits were found for a given T3E, it was considered absent. For a given strain, when a gene product was found to match with several T3Es, the one with the best E-value was selected. If there were more than one T3E with the best E-value, the alignment with the greatest number of identities was retained. Then, lengths of the query T3Es were compared to those of the alignments. We labelled a potential T3E as incomplete when the alignment was ≥25% smaller than the length of the original T3E. Otherwise, the T3E was labelled as complete. Based on the presence of complete T3Es, a binary matrix was created and used to generate a dendrogram by means of the R package APE (Paradis et al., 2004). Further information on this analysis can be found in Table S5 and File S1. Circular genome visualization Circular layouts were generated using Circos (Krzywinski et al., 2009). 65 Chapter I Supplementary Material Table S1. Predicted ORF in PssUMAF0158 plasmid (pPSS158). ORFs were annotated by the NCBI Prokaryotic Genome Annotation Pipeline and manually curated. Available: http://bacterial-virulence-factors.cbgp.upm.es/supplementary/chapter1/Table_S1.docx Table S2. Accession numbers and references for the genome sequences of 26 Pseudomonas syringae strains used in this work. Available: http://bacterial-virulence-factors.cbgp.upm.es/supplementary/chapter1/Table_S2.docx Table S3. Regions of UMAF0158 chromosome with low homology to B728a. Available: http://bacterial-virulence-factors.cbgp.upm.es/supplementary/chapter1/Table_S3.docx Table S4. Regions of B728a chromosome with low homology to UMAF0158. Available: http://bacterial-virulence-factors.cbgp.upm.es/supplementary/chapter1/Table_S4.docx Table S5. T3Es repertoires found in 26 Pseudomonas syringae strains. Columns provide information on BLASTp alignments between T3Es from http://pseudomonas-syringae.org/ and strains gene products, such as E-value, fraction and number of identical positions, alignment length, coordinates (start-end) for the query effector and subject gene product in the alignment and number of gaps. Other relevant features are provided, such as lengths of both the effector and the gene product, the rate between effector length and alignment length and whether the searched effectors are considered complete. Available: http://bacterial-virulence-factors.cbgp.upm.es/supplementary/chapter1/Table_S5.xlsx Table S6. Relevant virulence factors found in UMAF0158 chromosome. Available: http://bacterial-virulence-factors.cbgp.upm.es/supplementary/chapter1/Table_S6.xlsx File S1. BLASTp alignments (stockholm format) generated by the T3Es searches. Available: http://bacterial-virulence-factors.cbgp.upm.es/supplementary/chapter1/File_S1.zip 66 CHAPTER II Complete genome sequence of Pseudomonas fluorescens strain PICF7, an indigenous root endophyte from olive (Olea europaea L.) and effective biocontrol agent against Verticillium dahliae Pedro Manuel Martínez-García, David Ruano-Rosa, Elisabetta Schilirò, Pilar Prieto, Cayo Ramos, Pablo Rodríguez-Palenzuela and Jesús Mercado-Blanco. Complete genome sequence of Pseudomonas fluorescens strain PICF7, an indigenous root endophyte from olive (Olea europaea L.) and effective biocontrol agent against Verticillium dahliae.Standards in Genomic Sciences. 2015. doi: 10.1186/1944-3277-10-10. Chapter II Table 2. Project information. MIGS ID Property Term MIGS-31 Finishing quality Finished MIGS-28 Libraries used Three libraries of 500 bp, 2,000bp and 6,000bp, respectively MIGS-29 Sequencing platforms Solexa MIGS-30 Assemblers SOAPdenovo 1.05 MIGS-32 Gene calling method NCBI Prokaryotic Genome Annotation Pipeline Locus Tag PFLUOLIPICF7 Genbank ID CP005975 Date of Release May 31, 2017 GOLD ID Gi0079402 BIOPROJECT PRJNA203247 NCBI taxon ID 1334632 Project relevance Plant-bacteria interaction Model for endophytic lifestyle Agricultural, Environmental Table 3. Genome statistics. Attribute Genome (total) Value % of total Genome size (bp) 6,136,735 100 DNA coding region (bp) 5,439,499 88.6 DNA G+C content (bp) 3,706,588 60.4 DNA scaffolds 1 - Total genes 5,655 100 Protein-coding genes 5,567 98.4 RNA genes 68 1.6 Pseudo genes 30 0.8 Protein-coding genes with function prediction 4,573 82.1 Protein-coding genes assigned to COGs 4,581 82.3 Proteins with signal peptides 644 11.6 Proteins with transmembrane helices 1,319 23.7 75 Chapter II Table 4. Number of genes associated with general COG functional categories. Functional category Value % of total A RNA processing and modification 1 0.02 B Chromatin structure 5 0.09 C Energy production and conversion 280 5.03 D Cell cycle control, cell division, chromosome partitioning 41 0.74 E Amino acid transport and metabolism 554 9.95 F Nucleotide transport and metabolism 96 1.72 G Carbohydrate transport and metabolism 307 5.51 H Coenzyme transport and metabolism 196 3.52 I Lipid transport and metabolism 219 3.93 J Translation, ribosomal structure and biogenesis 200 3.59 K Transcription 501 9 L Replication, recombination and repair 156 2.8 M Cell wall/membrane/envelope biogenesis 267 4.8 N Cell motility 162 2.9 O Posttranslational modification, protein turnover, chaperones 177 3.18 P Inorganic ion transport and metabolism 301 5.41 Q Secondary metabolites biosynthesis, transport, and catabolism 151 2.71 R General function prediction only 592 10.63 S Function unknown 446 8.01 T Signal transduction mechanisms 366 6.57 U Intracellular trafficking, secretion, and vesicular transport 153 2.75 V Defense mechanisms 67 1.2 Z Cytoskeleton - - W Extracellular estructures - Not in COGs 986 17.7 76 Chapter II Figure 3. Graphical map of the chromosome. From outside to the centre: genes on forward strand (coloured by COG categories), genes on reverse strand (coloured by COG categories), RNA genes: tRNAs - blue, rRNAs – pink, G+C in relation to the mean G+C in 2kb windows and trinucleotide distribution in 2kb windows. The latter was defined as the χ2statistic on the difference between the trinucleotide composition of 2kb windows and that of the whole chromosome. Insights from the genome sequence The genome contains a complete canonical type III secretion system and two known effector proteins, namely, AvrE1 and HopB1. In addition, two complete type VI secretion system (T6SS) clusters were identified. T6SS has been described to promote antibacterial activity against a wide range of competitor bacteria (Hood et al., 2010). PICF7 genome also encodes gene clusters for the synthesis of the siderophores pyochelin and pyoverdine and the hemophore HasAp. A repertoire of cell adhesion proteins has been also identified, 77 Chapter II including two filamentous hemagglutinin proteins and several fimbrial proteins clustered together with a number of pilus assembly proteins. Notably, two genes have been found to show high similarity with attC and attG genes from Agrobacterium, whose mutation leads to lack of attachment on tomato, carrot, and Bryophyllum daigremontiana (Matthysse et al., 2008). It is worth mentioning the presence of genome components presumably involved in the synthesis of detoxifying compounds. Such is the case of two clusters containing genes for copper resistance and for production of a cbb(3)-type cytochrome C oxidase, respectively. An ortholog of the gene that codes for Dps, a ferritin-like protein reported to protect plantassociated bacteria against oxidative stress (Colburn-Clifford et al., 2010), has also been found. Additional identified traits involved in detoxification are orthologs of catalase KatB and hydroperoxidase KatG, which detoxify plant-produced H2O2(Guo et al., 2012), and a gene coding for a proline iminopeptidase, which has been shown to have dealanylating activity toward the antibiotic ascamycin (Sudo et al., 1996). A gene predicting a salycilic hydroxylase has been also identified in PICF7 genome. This gene could be involved in the degradation of the plant defence hormone salicylic acid, thus disrupting the systemic response against colonizing bacteria. In addition, all genes required for biosynthesis of the exopolysaccharide alginate (Vázquez Peñaloza et al., 1997) are present in a gene cluster. Genes predicting volatile components are present in PICF7 genome as well. Volatile components have been shown to act as antibiotics and to induce plant growth (Ryu et al., 2003; Ren et al., 2010). An example is hydrogen cyanide (HCN), an inorganic compound with antagonistic effects against soil microbes (Ahmad et al., 2008). Orthologs of genes required for the biosynthesis of other volatile components such as 2,3-butanediol and acetoin were also found. Further genome analysis revealed other factors presumably involved in the endophytic fitness of PICF7. Such is the case of enzymes like a cellulase and a phytase, as well as the gene coding for aminocyclopropane-1-carboxylate deaminase suggested to be key in the modulation of ethylene levels in plants by bacteria (Hardoim et al., 2008). 78 Chapter II Conclusions In this report we describe the complete genome sequence of Pseudomonas fluorescens strain PICF7, a “Pseudomonadales” in the order Gammaproteobacteria that was originally isolated from the roots of healthy nursery-produced olive plants cv. Picual in Córdoba province, Spain. This strain was selected for sequencing based on its ability to exert biocontrol against Verticillium wilt of olive and to develop an endophytic lifestyle within olive root tissues. Such properties likely have origins in a repertoire of genes including a putative T3SS, two putative T6SS, and several genes presumably implicated in siderophore production. It also has a collection of genes predicting adhesion proteins, detoxifying compounds, volatile components and enzymes such as a cellulase, aphytase and a deaminase. Further functional studies and comparative genomics with related isolates will provide insights into biocontrol and endophytism. 79 Chapter II Material and Methods Growth conditions and DNA isolation P. fluorescens strain PICF7 was grown in 50 ml of LB medium and incubated for 16 h at 28ºC. After this period of time, the OD600 of the culture was 1.2. Serial dilutions from this culture and plating on LB plates yielded 2.8 x 108 CFU/mL of a pure bacterial culture (colonies showed uniform morphology and kanamycin resistance). The culture was divided into two 25-ml aliquots and total genomic DNA was extracted using the ’Jet-Flex genomic DNA purification’ kit (Genomed GmbH, Löhne, Germany), according to the manufacturer’s indications. DNA samples were further purified by extraction with phenol:chloroform and precipitation with ethanol. DNA quality and quantity were checked by agarose gel electrophoresis, spectrophotometry using a ND1000 spectrophotometer (NanoDrop Technologies, Wilmington, DE), and digestion with different restriction enzymes. Two DNA aliquots (0.6 µg/µL, ∼20 µg each) were sent in a dry ice container to the sequencing service. Genome sequencing and assembly The genome of PICF7 was sequenced at the Beijing Genomics Institute (BGI) using Solexa paired-end sequencing. Draft assemblies were based on 3,482,351 reads with a length of 500 bp resulting in 1,200 Mb, 2,456,221 reads with a length of 2,000 bp resulting in 1,209 Mb and 1,924,515 reads with a length of 6,000 bp resulting in 1,309 Mb. The SOAPdenovo 1.05 software package (Li et al., 2008, 2010) developed by BGI was used for sequence assembly and quality assessment. Genome annotation Automatic annotation was performed using the NCBI Prokaryotic Genome Annotation Pipeline (Angiuoli et al., 2008). Identification of known type III effectors effectors was conducted by BLASTP searches of the effectors described in http://pseudomonassyringae.org/against the predicted protein sequences of PICF7. Functional annotation was performed by aligning the latter CDSs against the COG PSSM of the CDD (MarchlerBauer et al., 2014) using RPS-BLAST. Hits with an E-value <= 0.001 were first retained. Then, only the best hit was selected for each protein. Signal peptides and transmembrane helices were predicted using SignalP (Emanuelsson et al., 2007) and TMHMM (Krogh et al., 2001), respectively. T346Hunter (Martínez-García et al., 2015b) was used to identify secretion systems clusters. Genetic factors potentially implicated in bacteria-plant associations 80 Chapter II were predicted using a customized annotation tool developed by our group (not published). See Chapter IV for further details on such a tool. 81 CHAPTER III T346Hunter: A novel web-based tool for the prediction of type III, type IV and type VI secretion systems in bacterial genomes Pedro Manuel Martínez-García, Cayo Ramos and Pablo Rodríguez-Palenzuela. T346Hunter: A novel web-based tool for the prediction of type III, type IV and type VI secretion systems in bacterial genomes. PLoS ONE. 2015. doi: 10.1371/journal.pone.0119317. Chapter III mechanisms has recently emerged. Consequently, hundreds of genome sequences are nowadays available for B. pseudomallei (Nandi et al., 2014), and the role of T3SS and T6SS in the virulence of this species has been previously reported (D’Cruze et al., 2011; Burtnick et al., 2011). Figure 1B shows the whole-genome graphical overview generated by T346Hunter displaying the predicted secretion systems clusters for B. pseudomallei 668 chromosome 2 (RefSeq NC_009075). The system localises three clusters of NF-T3SS, one cluster of flagellar T3SS and five clusters of T6SS. These predictions are consistent with previously reported in silico analyses (Boyer et al., 2009; Abby, Rocha, 2012). The bsa NFT3SS of B. pseudomallei has been shown to be an important part of the virulence armoury of this strain (Stevens et al., 2004). Figure 1C shows the output generated by T346Hunter containing detailed information of such cluster, including a tabulated output and a gene map graphic representing its genomic context. Core components The sets of core components used for T3SS, T4SS and T6SS were as described by Abby and Rocha (2012), Bi et al. (2013) and Shrivastava and Mande (2008), respectively. However, there is no consensus on the definition of core in terms of secretion systems components. As long as we understand, “core” is the minimum set of components experimentally proven to be necessary for a secretion system to be functional. That appears to be the meaning used by Abby and Rocha (2012) and Shrivastava and Mande (2008) when they suggest a set of T3SS core components and T6SS components of major requirement, respectively. On the other hand, Bi et al. (2013) do not explicitly describe minimum required sets of components, but suggest a list of core components for each of the 18 T4SS they collect in their database. It is not clear though whether such proteins are indispensable for these T4SS to be functional. For instance, the trb T4SS encoded by A. tumefaciens C58 (Li et al., 1998) lacks TrbN, which belongs to the above core list. Given this controversy, T346Hunter makes no discrimination regarding the different uses of the core set, leaving to the user the role of interpreting the results. We chose 4 as the minimum threshold of core components after trying different values and manually screening the predicted secretion systems, since it offers a trade-off between false positives and false negatives. Again, it is the user who has to sift through the predictions. Further information about core components can be found in Table S1. 91 Chapter III Please leave your email, you will be notified once the job is done. DNA sequence A [email protected] T346Hunter: A novel web-based tool for the prediction of T3SS, T4SS and T6SS Upload your sequence files for secretion systems prediction E-mail NC_009075.fna E-value (HMMER) <= 0.0005 Sequence shape Circular Secretion systems to predict E-val T3SS T4SS T6SS Home About Methods Contact us (b) Secretion Systems Loci in NC_009075 General View General View T3SS T6SS flagella NF-T3SS T4SS T6SS Complete = YES: 9 of 9 core components9 of 9 core components 18 9 9 (Abby & Rocha, 2012) T3SS_cluster_type Num_compon Core_comp_in_cluster Total_core_NF-T3SS gene start end st len product SS_typ comp E_val BLAST_DNA BLAST_aa Non-flagellar surface presentation of antigens protein SpaS 3e-104 100 % Perc_core NCBI_BLASTN NCBI_BLASTP sctuT3SS 406-121029832101763 BURPS668_A2164 type III secretion protein SpaR/YscT/HrcT 6e-60 NCBI_BLASTN scttT3SS 247 -121037302102987 NCBI_BLASTP HrpO family type III secretion protein 8.4e-31 NCBI_BLASTN NCBI_BLASTP sctsT3SS 84 -12104030 2103776 surface presentation of antigens protein SpaP 2.4e-85 NCBI_BLASTN sctrT3SS 226-12104746 2104066 NCBI_BLASTP T3SS apparatus protein YscQ/HrcQ 9.5e-24 NCBI_BLASTN NCBI_BLASTP sctq T3SS 323-12105707 2104736 BsaU protein 2e-71 NCBI_BLASTN bsauT3SS 434 -121070082105704 NCBI_BLASTP 6e-77 NCBI_BLASTN NCBI_BLASTP bsatT3SS 154 -121074502106986 ATP synthase SpaL 3.4e-179 sctn T3SS -12108754 2107447 NCBI_BLASTP NCBI_BLASTN 435 type III secretion system protein 8e-99 NCBI_BLASTN NCBI_BLASTP bsar T3SS 135 -121091582108751 HrcV family type III secretion protein 7.5e-260 NCBI_BLASTN sctv T3SS 690 -121112422109170 NCBI_BLASTP T3S regulator YopN/LcrE/InvE/MxiC 5.5e-45 NCBI_BLASTN NCBI_BLASTP sctw T3SS 373 -121123992111278 YscC/HrcC family T3S outer membrane protein 8.7e-121 NCBI_BLASTN sctc T3SS 617 -121142492112396 NCBI_BLASTP AraC family transcriptional regulator NCBI_BLASTN NCBI_BLASTP none 221-121149282114263 hypothetical protein NCBI_BLASTN none 431 21151362115005 NCBI_BLASTP type III secretion system protein PrgH/EprH 4.3e-07 NCBI_BLASTN NCBI_BLASTP sctd T3SS 4281 21167012115415 type III secretion system needle protein 6.2e-30 sctf T3SS 121169672116698 NCBI_BLASTPNCBI_BLASTN putative type III secretion system protein 2.5e-15 NCBI_BLASTN sctiT3SS 1 21173292117027 NCBI_BLASTP YscJ/HrcJ family T3S apparatus lipoprotein 2.1e-60 NCBI_BLASTN NCBI_BLASTP sctj T3SS 315121182812117334 T3S apparatus protein OrgA/MxiK 8.8e-15 sctkT3SS 1 21188652118263 NCBI_BLASTPNCBI_BLASTN 89 200 HrpE/YscL family T3S apparatus protein 1.8e-20 sctlT3SS 12119562 2118834 NCBI_BLASTP NCBI_BLASTN 242 BsaT protein none none none none 100 Genomic context BURPS668_A2164 BURPS668_A2165 BURPS668_A2166 BURPS668_A2167 BURPS668_A2168 BURPS668_A2169 BURPS668_A2170 BURPS668_A2171 BURPS668_A2172 BURPS668_A2173 BURPS668_A2174 BURPS668_A2175 BURPS668_A2179 BURPS668_A2180 BURPS668_A2181 BURPS668_A2182 BURPS668_A2183 BURPS668_A2184 none NCBI_BLASTN NCBI_BLASTP nonenone hypothetical protein BURPS668_A2165 BURPS668_A2166 BURPS668_A2167 BURPS668_A2168 BURPS668_A2169 BURPS668_A2170 BURPS668_A2171 BURPS668_A2172 BURPS668_A2173 BURPS668_A2174 BURPS668_A2175 BURPS668_A2176 BURPS668_A2177 BURPS668_A2178 BURPS668_A2179 BURPS668_A2180 BURPS668_A2181 BURPS668_A2182 BURPS668_A2183 BURPS668_A2184 2115215 2115418 1 67 B C All hits Predicted clusters Download profiles Figure 1. Example of execution of T346Hunter using the sequence of B. pseudomallei 668 chromosome 2 as input. A. Web interface of T346Hunter. B. Genome-wide graphical view showing the predicted secretion systems of B. pseudomallei 668 chromosome 2. C. Genomic representation of one of the three NF-T3SS clusters identified, including a graphical gene map and a tabulated gene list with detailed information of each component. Hyperlinks to NCBI for direct execution of BLASTn and BLASTp against the non-redundant nucleotide and protein databases are provided for each gene within the loci. Some other relevant information is also included, such as the percentage of core components found in the cluster and PubMed hyperlinks to the studies we have based our methods on to build the component profiles found in such a cluster. 92 Chapter III Identification of T3SS, T4SS and T6SS gene clusters in sequenced bacterial genomes Complete bacterial genomic sequences of 2,997 chromosomes and 2,164 plasmids sequenced available as of 14 February 2014 were downloaded from the NCBI RefSeq project Pruitt et al. (2012). T346Hunter was executed on these sequences and localised clusters enriched in either NF-T3SS, T4SS or T6SS components. In total, 2,814 clusters were identified (512 NF-T3SS, 1,466 T4SS and 836 T6SS) across 1,121 organisms. Predicted clusters are summarised in Table S2 and can be queried at the Predicted Clusters interface from the T346Hunter website. Sequences with negative predictions are listed in Table S3. Comparison with currently available tools In order to validate the performance of T346Hunter, systematic comparisons with other available applications for secretion systems prediction would certainly be the best choice. By crossing our predictions with loci identified by other servers one could have a measure of the relative accuracy of our tool. But, in practise, this is not straightforward to carry out. On the one hand, predicted data are not always available in a format that are ready to be systematically analysed. Servers usually provide their data in a way that either they have to be queried using some kind of criteria (e.g. strain name) or they are just embedded in the webpage. On the other hand, few servers are available to predict secretion systems clusters as such. To our knowledge, only SecReT4 (Bi et al., 2013) provides a specific tool for genomic localisation of T4SS clusters. Despite these difficulties, and given the need for assessing the accuracy of our predictions relative to others’ work, we attempted to accomplish comparisons either by systematic processing, when possible, or by manual inspection, when not. First, we aimed to compare our predictions of T3SS with those of T3SSscan-FLAGscan (Abby, Rocha, 2012). Such server does not localise T3SS clusters in user-submitted sequences, but does keep a repository of predicted NF-T3SS loci that can be accessed using different features. Therefore, we randomly queried clusters for 100 genomic sequences and manually compared them with our predictions. Since negative predictions are not provided in the server, this selection was restricted to positive predictions. We found that T346Hunter predicted any NF-T3SS cluster in the 100 sequences examined. More precisely, T346Hunter and T3SSscan-FLAGscan identified the same number of clusters in 95 sequences (Table S4). For each of the resting five sequences, the number of predicted clusters just differs in one. This difference is probably explained by the different searching criteria used by both methods, particularly in defining contiguous genes within a cluster and setting the minimum required number of core components. In total, 125 clusters were 93 Chapter III identified by T3SSscan-FLAGscan and 128 by T346Hunter. To test how T346Hunter performs in predicting T4SS, we compared it with SecReT4 (Bi et al., 2013). This server provides a summary of T4SS predictions on a number of sequences by means of a table in the webpage. We could, then, proceed to systematically compare the two methods. Again, such a comparison was necessarily restricted to positive predictions. In total, 387 sequences were examined and all of which were found to contain at least one T4SS cluster. Of 387 such replicons, 324 (84%) were predicted by both methods to contain the same number of T4SS (Table S5). In this case, differences in searching criteria may have had a stronger impact in the predictions. SecReT4, for instance, does not restrict T4SS genes to be located in a specific cluster. In contrast, T346Hunter identifies a cluster whenever orthologues of 4 core components are found in a window of up to 70 kb. Despite such divergent parametrisation, 53 out of the 63 sequences with discordant predictions (84%) only differ in one cluster. Accounting for the 387 sequences analysed, SecReT4 and T346Hunter identified 522 and 595 T4SS, respectively. We went ahead to inspect the accuracy of our T6SS predictions. As a specific server for the prediction of T6SS is not available, tool-by-tool comparison was precluded in this case. However, systematic localisation of T6SS loci has been previously performed (Shrivastava, Mande, 2008; Boyer et al., 2009; Barret et al., 2011), reporting data we could use to compare our predictions with. We focused on Boyer et al. (2009), which completed the widest analysis. Since no readily processable summary of the predictions was provided, comparisons had to be manually performed. Among the 100 sequences reported with positive predictions, T346Hunter identified at least one T6SS in 98 of them. Furthermore, both approaches predicted the same number of T6SS clusters in 93 sequences (Table S6). Such a subtle difference is explained by the requirement of 4 core components in a cluster imposed by T346Hunter, which was not applied by Boyer et al. (2009). Nonetheless, predictions on the resting 7 sequences of both approaches differed in only one cluster. Summing up all identified clusters in the 100 sequences examined, T346Hunter and Boyer et al. (2009) predicted 170 and 175, respectively. Regarding general annotation engines, some of the most widely used tools are the NCBI Prokaryotic Genome Annotation Pipeline (Angiuoli et al., 2008) and the RAST server system (Overbeek et al., 2014). In the last few years, RAST has become particularly popular and is now frequently used to rapidly annotate bacterial genomes against its comprehensively curated subsystem database. Due to its constant growth, RAST automated annotations are nowadays of a great quality, having reached a high degree of specificity. Indeed, T3SS, T4SS and T6SS are included among the subsystems 94 Chapter III collection of RAST, and thorough reports of related genes are provided within general annotations. Such reports include visual and tabular information of the corresponding genomic clusters, thus offering an exhaustive output. However, bacterial strains that has not been incorporated into RAST database are reported with no subsystems, and users need to manually inspect individual features to infer the existence of secretion systems clusters. This makes RAST subsystems search dependent on its database of bacterial isolates, and makes it particularly not suitable for the analysis of newly characterised bacteria. Furthermore, when subsystems are reported, the number of secretion systems clusters are not directly shown in the output and it rather needs to be derived from reported tables. Besides, no information regarding core components is provided, and some conjugal T4SS are not categorised as such. Therefore, even though RAST performs quite well in detecting genes encoding secretion systems when compared to other general annotation tools, its annotations lack some relevant information on T3SS, T4SS and T6SS, and do not directly offer the whole picture of the underlying genomic clusters. 95 Chapter III Conclusions The development of web-based tools for the prediction of virulence factors is crucial for allowing researchers to identify the bacterial pathogenic arsenal. Here, we present T346Hunter, an online tool for annotation and localisation of secretion systems clusters in sequenced bacterial genomes. Because they are distinctive features of pathogenesis, T346Hunter searches for T3SS, T4SS and T6SS, whose identification is of particular interest in the development of strategies against bacterial-mediated diseases. The server will be continuously updated as new experimental and bioinformatics information on secretion systems becomes available. We believe T346Hunter will help researchers uncover the mechanisms of bacterial secretion as a virulence trait. 96 Chapter III Material and Methods Protein sequences and profiles Sequence profiles of secretion systems components were generated by selecting orthologues of each component in order to capture the diversity of the T3SS, T4SS and T6SS. Following the approach described by Abby and Rocha (2012), we selected protein sequences corresponding to the components of the flagellar and non-flagellar T3SS (NFT3SS). Meanwhile, to build profiles that represent the variety of T4SS, protein sequences of the components of 18 archetypal T4SS (Bi et al., 2013) were also selected. In both cases, sequences of each component were extracted from a set of model organisms representative of the diversity of all these types of systems. We based the construction of T6SS components profile on a previously reported list of coding sequences belonging to several bacterial genomes; sequences that were found to be orthologues of the components of the first described T6SS, namely, V. cholerae,Pseudomonas aeruginosa and Burkholderia mallei (Shrivastava, Mande, 2008). The sequences of these orthologues were also included in our set of protein families. All these sequences were downloaded from Uniprot (Boutet et al., 2007) and the NCBI website (http://www.ncbi.nlm.nih.gov/), and selected based on their corresponding genome annotations. Then, sequences corresponding to each component were aligned with Muscle (Edgar, 2004) and manually adjusted with Seaview (Gouy et al., 2010). Finally, protein profiles were built with HMMER3 (Eddy, 2011). When fewer than five representative model organisms were found to code for a given component, no profile was built; instead, files were generated in multi-FASTA format. We extended the above set by collecting protein sequences from AtlasT4SS (Souza et al., 2012) and proceeding in the same way as above to generate profiles of orthologue clusters described in that database. Secretion systems loci identified in this study were manually screened, and additional profiles were incorporated based on the RefSeq genome annotation. Consequently, our final dataset comprises sequence information for a total of 364 components of the T3SS, T4SS and T6SS (65, 449 and 20, respectively). Further information about each of the component profiles can be found in Table S1. Identification of secretion systems clusters T346Hunter performs BLASTp (Altschul, 2005) and HMMER3 (Eddy, 2011) searches of the protein sequences and profiles described above against user-supplied genomic sequences. Regions containing homologous genes (E-value ≤0.0005 by default) of at least 4 different core components of T3SS, T4SS or T6SS and spanning up to 70 kb are retained 97 Chapter III and included in the output report. We consider core components of T3SS, T4SS and T6SS as described by Abby and Rocha (2012), Bi et al. (2013) and Shrivastava and Mande (2008), respectively (see Table S1 for details). Implementation T346Hunter runs on a Linux platform with an Apache web server. The web interface was implemented using HTML and CSS, and data pipelines were developed using PHP, Perl, R and shell scripts. DNA and protein sequences are processed by means of Bioconductor (Gentleman et al., 2004) and BioPerl (Stajich et al., 2002). Circular genome images are generated using circos (Krzywinski et al., 2009), and gene maps are produced using the R package genoPlotR (Guy et al., 2010). Open reading frame predictions are generated with Glimmer v3.02 (Delcher et al., 2007). 98 Chapter III Supplementary Material Figure S1. Overview of T346Hunter prediction workflow. Number of sequences Concatenate sequences >1 1 GLIMMER Search secretion systems components (HMMER3/BLASTp) T3SS/T4SS/T6SS DATABASE Please leave your email, you will be notified once the job is done. DNA sequence [email protected] T346Hunter: A novel web-based tool for the prediction of T3SS, T4SS and T6SS Upload your sequence files for secretion systems prediction E-mail NC_009075.fna Home About Methods Predicted clusters Contact us ...or upload NCBI sequence files for a faster execution [email protected] NC_009075.fna E-value (HMMER) <= 0.0005 Sequence shape E-value (BLASTp) <= Circular Secretion systems to predict T3SS 0.0005 DNA sequence NC_009075.fna Protein sequences NC_009075.ptt Genes Location E-mail T4SS T6SS Identify clusters (b) Secretion Systems Loci in NC_009075 General View General View T3SS T6SS flagella NF-T3SS T4SS T6SS ... Gene Hit SS Eval ... id0320 sctu T3SS 3e-104 id0321 sctt T3SS 6e-60 id0322 scts T3SS 8.4e-31 id2144 vasf T6SS 1.2e-10 ... ... ... id2147 vca0107 T6SS 3.3e-60 ... … … ... … … ... … … ... Tab-delimited output (all hits) HTML-formatted output (clusters) All hits User manual Download profiles 99 Chapter III Table S1. Protein profiles and sequences used in this study. Available: http://journals.plos.org/plosone/article/asset?unique&id=info:doi/10.1371/journal.pone.0119317.s002 Table S2. Summary of predicted clusters. Available: http://journals.plos.org/plosone/article/asset?unique&id=info:doi/10.1371/journal.pone.0119317.s003 Table S3. List of complete bacterial genomic sequences and plasmids not found to encode NF-T3SS, T4SS or T6SS. Available: http://journals.plos.org/plosone/article/asset?unique&id=info:doi/10.1371/journal.pone.0119317.s004 Table S4. Summary of NF-T3SS clusters predicted by T3SSscan-FLAGscan and T346Hunter on 100 sequences. Available: http://journals.plos.org/plosone/article/asset?unique&id=info:doi/10.1371/journal.pone.0119317.s005 Table S5. Summary of T4SS clusters predicted by SecReT4 and T346Hunter on 387 sequences. Available: http://journals.plos.org/plosone/article/asset?unique&id=info:doi/10.1371/journal.pone.0119317.s006 Table S6. Summary of T6SS clusters identified in Boyer et al. (2009) and by T346Hunter on 100 sequences. Available: http://journals.plos.org/plosone/article/asset?unique&id=info:doi/10.1371/journal.pone.0119317.s007 Table S7. Predictions of secretion systems clusters not fulfilling the restriction of 4 core components. Available: http://journals.plos.org/plosone/article/asset?unique&id=info:doi/10.1371/journal.pone.0119317.s008 Data S1. Hidden Markov Model profiles and sequences used by T346Hunter. Available: http://journals.plos.org/plosone/article/asset?unique&id=info:doi/10.1371/journal.pone.0119317.s009 100 Chapter IV Results and Discussion PIFAR, an open-access web-based resource for plant-bacteria interaction factors Although our aim for collecting bacterial factors was to test their predictive potential for the classification of plant-associated bacteria, we have also deposited them in a web server, called PIFAR (Plant-bacteria Interaction FActors Resource), so that they can be accessed by external users (Figure 1). PIFAR is freely available at http://bacterial-virulencefactors.cbgp.upm.es/PIFAR. The website has two main purposes: to hold comprehensive information on gene products described as implicated in bacterial interactions with plants and to help researchers to identify such products in input genome sequences. PIFAR holds a curated set of products that cover a broad spectrum of the current genomic information on plant-associated bacteria, including pathogenic (i.e. Pseudomonas syringae) and non pathogenic species (i.e. Rhizobium leguminosarum). Entries in the database are structured in terms of the molecular processes they are involved. To account for adhesion and attachment to plant tissues, genes implicated in the synthesis of factors such as adhesins or exopolysaccharides (EPSs) were included (Pizarro-Cerdá, Cossart, 2006; Mhedbi-Hajri et al., 2011). Regarding bacterial survival, genes involved in the production of factors like siderophores, antibiotics, lipopolysaccharides (LPSs) or multidrug efflux pumps (MDRs) were considered (Erbs, Newman, 2003; Wang, Raaijmakers, 2004; Martínez et al., 2009; Ahmed, Holmström, 2014). Meanwhile, to account for plant-host colonization, genes implicated in the synthesis of factors such as phytotoxins, phytohormones or plant cell wall degrading enzymes (PCWDEs) were also incorporated (Dudler, 2014; Costacurta, Vanderleyden, 1995; Barras et al., 1994). The final collection comprised 170 factors involving 603 gene products and 16 PFAM domains (Finn et al., 2014). Factors were selected by manual inspection of the available scientific literature. The curation based on selecting relevant genes reported to affect the outcome of bacteria-plant interactions. We are aware that the literature is somehow biased towards virulence factors, since pathogenesis has been more extensively studied than any other bacterial association to plants. Similarly, information on some bacterial species and/or genera is also overrepresented. However, we do not conceive our repertoire of factors as finished. On the one hand, it will be updated as new experimental data becomes available. On the other hand, we expect it to grow by means of the contribution of other experts in plant bacteriology. To this end, PIFAR allows data submissions of missing factors or newly charaterized ones. Further information on how to use the web server can be found in File S1. 107 Chapter IV Figure 1. Front page of PIFAR. Several interfaces are available for the user to navigate throught the website (left; see File S1). The pie chart represents the distribution of bacterial factors included in the database. By clicking on the different segments, users can access the database by factor type. A supervised machine learning strategy to classify plant-associated bacterial lifestyles Three classifiers were trained to discriminate between different bacterial lifestyles (Figure 2). To that aim, we selected 420 sequenced bacterial strains (Figure 2A; Table S1) classified into three classes: plant-pathogenic (PP, 109 strains), plant-associatednon-pathogenic (PANP, 93 strains) and non-plant-associated (NPA, 218 strains). This classification was based on literature searches and NCBI annotations (see Material and Methods). Then, we performed a systematic screening of our database against such strains (Figure 2B). For each genome analyzed, we obtained its repertoire of plant-associated bacterial factors. Based on presence/absence of such factors, a vector was generated for each genome so that each position contained the count of the identified factors by type (Figure 2C). Therefore, a given feature in a vector represents the total number of factors of a certain type (i.e. phytotoxins) found in a given bacterial strain. These features were extended using T346Hunter (Martínez-García et al., 2015b). T346Hunter is an annotation tool, developed by our group, which identifies genomic clusters potentially involved in the synthesis of secretion systems of types III, IV and VI (T3SS, T4SS and T6SS, respectively). Thus, three more features were incorporated to the vectors consisting of the numbers of T3SS, T4SS and T6SS genomic clusters that a given genome harbors (T4SS involved in conjugation and 108 Chapter IV DNA uptake/release were excluded). The result was a 420x19 matrix (Figure 2D) with no repeated rows (a unique vector resulted from each strain) that can be checked in Table S2. Finally, we used random forests (Breiman, 2001) to train three supervised machine-learning classifiers (Figure 2E and Figure 3). Sequenced bacterial genomes Manual selection PLANT-ASSOCIATED NON-PATHOGENIC (PANP) 93 strains PLANT-PATHOGENIC (PP) 109 strains NON PLANT-ASSOCIATED (NPA) 218 strains PIFAR + T346Hunter Identification of factors involved in bacterial interactions with plants Generation of vectors based on counts of identified factors Strain Toxins PCWDEs … T3SS T4SS T6SS P. syringae DC3000 1 7 ... 1 0 2 D. dadantii 3937 0 14 ... 1 1 1 A. tumefaciens C58 0 4 ... 1 2 1 ... ... ... ... Strain Toxins PCWDEs … T3SS T4SS T6SS P. fluorescens PICF7 0 2 ... 1 0 2 A. radiobacter K84 0 4 ... 1 2 0 E. billingiae Eb661 0 2 ... 0 0 2 ... ... ... ... Strain Toxins PCWDEs … T3SS T4SS T6SS P. marinus MIT 9215 0 1 ... 0 0 0 S. arenicola CNS-205 0 2 ... 0 0 0 P. staleyi DSM 6068 0 1 ... 0 0 1 ... ... ... ... Toxins PCWDEs … T3SS T4SS T6SS CLASS 1 7 ... 1 0 2 PP 0 14 ... 1 1 1 PP 0 2 ... 1 0 2 PANP 0 4 ... 1 2 0 PANP 0 1 ... 0 0 0 NPA 0 2 ... 0 0 0 NPA - - - - - - - - Create labeled matrix Train supervised models Random forests classifiers PP PANP A C B E D NPA - - - - Feat. 1 Feat. 2 … Feat. 19 Label Figure 2. Schematic representation of the bioinformatics pipeline used for the generation of the training set. A. The genomes of a selection of 420 bacterial strains were retrieved from NCBI. Strains were selected according to three categories: plant-pathogenic (PP), plant-associated non-pathogenic (PANP) and non-plant-associated (NPA). B. PIFAR annotation tool was used to screen these strains for factors in the database. T346Hunter (Martínez-García et al., 2015b) was used to identify genomic clusters encoding T3SS, T4SS and T6SS. C. The results of above searches were used to create vectors containing the repertoire of factors of each genome. D. Each vector was labeled as PP, PANP or NPA and a matrix was created consisting of all 420 vectors. E. The resulting matrix was used as the training set of three random forests classifiers. 109 Chapter IV First, we generated a 3-class model, which classified bacterial genomes into PP, PANP and NPA with a precision rate of 90% (Figure 3A; Table 1). The classifier seems to perform quite better in predicting NPA strains (212 out of 218; 97.3%) than PP (94 out of 109; 86.2%) or PANP (72 out of 93; 77.4%). The latter presents the highest error rate in this study (22.6%), which may be explained by the bias towards virulence factors existing in our database. This feature could be showing that, to date, more data on pathogenicity than on other kind of plant-bacteria interactions is available in the literature. Anyway, given its error rate, we consider this PP-PANP-NPA classifier as limited. By combining PANP and NPA classes into a new class consisting of non-plantpathogenic (NPP) strains, we generated a second model that classified bacterial genomes into PP and NPP with a precision rate of 94% (Figure 3B, Table 2). PP classes were correctly classified in 85.3% of the cases (93 out of 109), while NPP classification reached an accuracy of 97.7% (304 out of 311). Such a significant difference in error rates could barely be explained by the lightly unbalanced dataset (109 PP versus 311 NPP). More likely, it may be regarded as a consequence of the database content itself. It is possible that we are missing unknown factors that could help to better capture the pathogenic lifestyle of several phytobacteria, particularly those misclassified by our models. Such factors, yet to be described, would complement the repertoire of well-known virulence factors in these organisms, and their characterization would help improving the performance of this classifier. In any case, the available information on bacterial genetic factors seems to be limited to accomplish PP-NPP classification. Finally, a third model was built by merging PP and PANP classes into a new class that harbors all plant-associated (PA) bacterial strains included in the dataset. Such model classified bacterial genomes into PA and NPA with a precision rate of 94.3% (Figure 3C, Table 3). In this case, the general error rate is the most similar to that of false positives and false negatives of the three classifiers. PA classes were correctly assigned in 188 out of 202 strains (93.1%), while NPA classes were correctly classified in 208 out of 218 strains (95.4%). This slight difference may be due to the fact that factors included in our database are by definition involved in associations to plants, and thus capture plant-associated bacterial traits. Since this is the positive class to be predicted, the unavoidable lack of information may be prejudicing PA classification. We expect this tendency will decrease as new experimental data is incorporated into the database. Nevertheless, the classifier shows a good performace, with a consistent number of correctly/incorrectly classified genomes (396/24) and an acceptable balance between false positives and false negatives (6.9%-4.6%). In this case, our repertoire of genetic components has enough predictive strength to separate 110 Chapter IV PA and NPA bacteria with high accuracy (94%). Therefore, we focused on this classifier for further analyses. PP-PANP-NPA classifier Toxins PCWDEs … T3SS T4SS T6SS CLASS 1 7 ... 1 0 2 PP 0 14 ... 1 1 1 PP 0 2 ... 1 0 2 PANP 0 4 ... 1 2 0 PANP 0 1 ... 0 0 0 NPA 0 2 ... 0 0 0 NPA - - - - - - - - - - - - Toxins PCWDEs … T3SS T4SS T6SS CLASS 1 7 ... 1 0 2 PP 0 14 ... 1 1 1 PP 0 2 ... 1 0 2 PANP 0 4 ... 1 2 0 PANP 0 1 ... 0 0 0 NPA 0 2 ... 0 0 0 NPA - - - - - - - - - - - - Toxins PCWDEs … T3SS T4SS T6SS CLASS 1 7 ... 1 0 2 PP 0 14 ... 1 1 1 PP 0 2 ... 1 0 2 PANP 0 4 ... 1 2 0 PANP 0 1 ... 0 0 0 NPA 0 2 ... 0 0 0 NPA - - - - - - - - - - - - Precision = 0.87 Accuracy = 0.93 Sensitivity = 0.87 PP-NPP classifier Precision = 0.94 Accuracy = 0.95 Sensitivity = 0.92 PA-NPA classifier Precision = 0.94 Accuracy = 0.94 Sensitivity = 0.94 109 strains 93 strains 218 strains 109 strains 311 strains 202 strains 218 strains B A C Figure 3. Training data assembly to train the three classifiers. A. PP VS PANP VS NPA classifier. Assembly of training data as generated in Figure 2. B. PP VS NPP classifier. PANP and NPA are combined to form the NPP class. C. PA VS NPA classifier. PP and PANP are combined to form the PA class. 111 Chapter IV Table 1. Confusion matrix showing classification performance of the three-classes classifier. PP, PANP and NPA refer to plant-pathogenic, plant-associated non-pathogenic and nonplant-associated, respectively. Total Predicted as PP Predicted as PANP Predicted as NPA Plant-pathogenic 109 94 (86.2%) 8 (7.3%) 7 (6.4%) Plant-associated (NP) 93 5 (5.4%) 72 (77.4%) 16 (17.2%) Non-plant-associated 218 2 (0.9%) 4 (1.8%) 212 (97.3%) Table 2. Confusion matrix showing classification performance of the PP-NPP classifier. PP and NPP refer to plant-pathogenic and non-plant-pathogenic, respectively. Total Predicted as PP Predicted as NPP Plant-pathogenic 109 93 (85.3%) 16 (14.7%) Non-plant-pathogenic 311 7 (2.3%) 304 (97.7%) Table 3. Confusion matrix showing classification performance of the PA-NPA classifier. PA and NPA refer to plant-associated and non-plant-associated, respectively. Total Predicted as PA Predicted as NPA Plant-associated 202 188 (93.1%) 14 (6.9%) Non-plant-associated 218 10 (4.6%) 208 (95.4%) Adhesion, plant cell wall degradation and detoxification as the most predictive features for identifying plant-associated bacteria An attractive aspect of random forests (Breiman, 2001) is that it assesses an importance measure to each input variable in the final predictions. This importance is computed by measuring the decrease of classification performance when the values of a given variable in a given tree of the forest are randomly permuted (Chen, Ishwaran, 2012). The use of random forests importances to rank variables in our prediction model can give us hints on which repertoire of bacterial factors are the most relevant in the classification, and therefore, are of biological interest. Figure 4 shows a ranking with the importance measure assigned to each of the features in the PA-NPA classifier. Factor types having the largest importance are those related to bacterial adhesion. Adhesion mechanisms are widely conserved in both pathogenic and non-pathogenic plant associated bacteria, and are used by multiple 112 Chapter IV species as a step previous to colonization or biofilm aggregation (Pizarro-Cerdá, Cossart, 2006; Rigano et al., 2007). The second feature in the ranking of importances correspond to PCWDEs. The secretion of PCWDEs is a widely distributed machanism of plantassociated bacteria, and are used to penetrate the cell wall prior to plant infection. PCWDEs have been extensively studied in soft-rot enterobacteria, such as Pectobacterium and Dickeya (Barras et al., 1994; Toth et al., 2003). The third type of factors in the ranking correspond to those involving bacterial detoxification. In the last two decades, thorough research has been dedicated to the mechanisms that enable plant-associated bacteria to resist the action of plant products such as antimicrobial agents (López-Solanilla et al., 1998; GarcíaOlmedo et al., 2001) or active oxygen species (Hassouni et al., 1999; Miguel et al., 2000; Colburn-Clifford et al., 2010). These findings show that defense against plant-derived toxic compounds plays a key role for bacterial survival in host cells. Toxins Pigments LPSs Antibiot. T4SS Biofilm Volatiles Hormones MAMPs Sideroph. Metabol. T3SS Proteases MDRs EPSs T6SS Detoxif. PCWDEs Adhesion % of random forests importance 0 5 10 15 Figure 4. Ranking of predictive features based on random forests importances. Random forests importances (Gini index) were normalized dividing individual values by the sum of all importances. Horizontal axis correspond to the percentage of such normalized values. Regarding features with low importance measures, it is worth mentioning that phytotoxins are assigned the smallest value. Phytotoxins are products of the host-pathogen interaction that directly injure plant cells, but they are not required for pathogenicity . 113 Chapter IV Moreover, synthesis of phytotoxins is not a common trait, and only some bacteria such as P. syringae are well-known to produce them. Another significant factor that is assigned a low importance is the T4SS. Note that we have only considered in the analysis T4SS involved in effector translocation. Although the T4SS is an important machinery in plantmicrobe interactions (e.g. Agrobacterium), it is not a common feature of plant-associated bacteria. Nevertheless, we must be cautious when interpreting importance measures, since it is likely that multiple sets of lowly predictive factors are in the end jointly predictive (Chen, Ishwaran, 2012). The application of the PA-NPA classifier to sequenced bacterial genomes reveals potential associations between human-pathogenic bacteria and plants The application of our PA-NPA classifier offers the possibility of systematically determining whether a given sequenced bacterial strain could potentially associate with plants. With the increasing number of foodborne poisoning outbreaks that have been linked to fresh produce (Holden et al., 2009; Lim et al., 2014), prediction of potential interactions between human pathogens and vegetables can be very useful for clinical and industrial purposes. Therefore, we went ahead by screening with our model the whole set of complete and draft bacterial genomes available as of 22 December 2014 in the NCBI Genome database (http://www.ncbi.nlm.nih.gov/genome). A total of 9,446 strains were analyzed (see Material and Methods), obtaining their corresponding outcomes (plantassociated/non-plant-associated) together with a probability (0-1; Table S3). All assigned probabilities separated by genus were deposited in File S2. Figure 5 shows the distribution of probabilities assigned to all strains. We observe that most of them are given quite low plant-associated probabilities, 50% of them being lower than 0.15 (median) and 75% lower than 0.47 (3rd quartile). Random forests’ default threshold is 0.5 when labeling a strain as PA; however, for purposes of discussion, and aiming to maintain a conservative approach, we will not discriminate basing on random forests’ labels and will just focus on assigned probabilities. All in all, paradigmatic pathogenic and non-pathogenic plant-associated bateria are classified quite accurately, with associated probabilities ranging from 0.9 to 1. That is the case of species such as Agrobacterium fabrum,Agrobacterium radiobacter,P. syringae, Pseudomonas fluorescens,Burkholderia cepacia,Xanthomonas campestris,R. leguminosarum, Dickeya dadantii,Xylella fastidiosa,Bacillus amyloliquefaciens and Bacillus subtilis, among others. Since strains belonging to these species have been used to train the classifier, it seems reasonable that plant-associated bacterial variations of them are identified as such. It 114 Chapter IV is worth mentioning the resulted outcomes relative to the genus Pseudomonas.Pseudomonas are well-known for their ability to occupy a wide spectrum of environments as freeliving soil organisms or as prominent opportunistic pathogen (Auling, 2001). Among the 221 analyzed Pseudomonas strains (Figure 6), 200 (91%) are assigned plant-associated probabilities above 0.85, and only 10 (5%) below 0.7. As a case example, the 48 P. aeruginosa isolates are given values ranging from 0.93 to 1. This bacterium is a well documented opportunistic human pathogen, but is also capable of associate and even cause infections to plants (He et al., 2004; Walker et al., 2004).The lowest probability (0.06) among the Pseudomonas genus is assigned to the strain Pseudomonas thermotolerans DSM 14292, which was isolated from a industrial cooking water and is able to growth at 55º C (Manaia, Moore, 2002). Plant−associated probability Number of strains 0.0 0.2 0.4 0.6 0.8 1.0 0 500 1000 1500 Figure 5. Distribution of probabilities assigned by the PA-NPA classifier to all sequenced bacterial strains in NCBI. It is also worth to note the outcomes corresponding to the genus Bacillus.Bacillus comprises a number of species and exhibits a wide range of physiologic abilities, which make them being present in any natural environment (Todar, 2005). This feature is somehow observed when analyzing Bacillus’ plant-associated probabilities (Figure 6). The genus is divided into two: one cluster containing strains being assigned considerable high PA-probabilities (ranging from 0.78 to 1), suchas B. subtilis,B. amyloliquefaciens and Bacillus 115 Chapter IV licheniformis, and other cluster harbouring strains with relatively low plant-associated probabilities (0.2-0.5), such as Bacillus thuringiensis,Bacillus cereus and Bacillus anthracis. Human pathogens deserve special attention. We observe that several well-known bacteria such as Staphylococcus,Streptococcus,Helicobacter or Legionella are sharply assigned low PA-probabilities, which generally range from 0 to 0.25 (Figure 6; Table S3). In contrast, other human pathogens present significantly higher values, specially strains belonging to Enterobacteriaceae. Such is the case of E. coli strain O157:H7, a fatal enterohemorrhagic bacterium whose infection has been traditionally associated to bovine food (Ayscue et al., 2009). Back in 1996, this pathogen caused an important outbreak of foodborne disease, with thousands of cases reported, apparently because of the ingestion of contaminated radish sprout salad (Watanabe et al., 1999). Interestingly, all the 32 E. coli O157:H7 strains analyzed were assigned significant plant-associated probabilities, ranging from 0.75 to 0.82. Similarly, another serotype of E. coli (strain O104:H4) that recently provoked a large epidemic of the hemolytic-uremic syndrome in Germany (Stockman, 2013) seems to present relatively large values. The total of 49 sequenced isolates of E. coli O104:H4 show an average probability of 0.69. In general, the 756 bacteria of the Escherichia genus show a spread distribution of their assigned probabilities (Figure 6), with two prominent clusters with values ranging from 0.5 to 0.6 and from 0.7 to 0.9, respectively. This is in agreement with an increasing opinion that plants are secondary reservoirs for commensal and pathogenic E. coli strains (Méric et al., 2013). Other well-known pathogen that has been suggested not to be exclusively devoted to human hosts is Salmonella.Salmonella is the causal agent of several worldwide diseases such as gastroenteritis or typhoid fever, provoking around 1.4 million human illnesses annually only in the United States (Mead et al., 1999). A number of studies have shown that this bacterium is able to associate with plants, colonize the phyllosphere and even infect and cause death to plant organs (Kirzinger et al., 2011; Schikora et al., 2012). Like E. coli, many cases have been recently reported linking food poisoning with the ingestion of Salmonellacontaminated raw fruits and vegetables. It is not surprising, then, that around half of the 325 analyzed Salmonella strains show plant-associated probabilities above 0.7 (Figure 6). Although strains of the Vibrio genus are mostly assigned low probabilities (around 0.3), members of the species Vibrio vulnificus,Vibrio harveyi,Vibrio parahaemolyticus or Vibrio nigripulchritudo generally show quite high values. Interestingly, V. parahaemolyticus has been described as associated with plant roots as well as to fix nitrogen, suggesting this organism could use the rhizosphere as a refugium until an eventual transfer to the overlying waters, potentially contributing to infectious outbreaks (Criminger et al., 2007). 116 Chapter IV learning classification can be found in (Sokolova, Lapalme, 2009). Feature selection We tested whether there are features in the vectors that are uninformative and can be removed. We obtained random forests (Breiman, 2001) importance measures and ranked all features in the vectors. We kept the whole set of features and run 100 random forest classifiers. The classifier with the best precision rate was recorded. This process was repeated by decreasing the number of features based on their rank. In the first step, the feature with lowest importance measure was removed. In the second step, the two features with lowest importance measures were removed, and so on. After the whole process, all classifiers were compared and the one with the best precision rate was selected. This was repeated for the three types of classifiers. In the three cases, the best classifiers were those based on the whole set of features, and consequently, no feature was removed. Machine Learning: algorithms, configuration and assessment Cross-validation analyses were used to assess the performance of the trained classifiers. We carried out several methods in order to cover different supervised machine learning algorithms, such as decision trees, bayesian or lazy learners. Six classification algorithms were evaluated, namely, bayesian networks, simple logistic regression, J48, naive bayes, IBk and random forests. Bayesian network and simple logistic regression were used as base classifiers, while J48, naive bayes and IBk were performed as ensemble classifiers. Random forests is a combination of bagging and a random selection of features (Breiman, 2001). These algorithms are reviewed in detail in Larrañaga et al. (2006). Each supervised classifier was run 100 times using 10-fold cross-validation, each iteration with a randomly generated seed. In the specific case of random forests, default parameters were used: the number of variables randomly sampled as candidates at each split was set as the square root of the number of features (√19 ∼4), the number of total trees in each iteration as 500 and the percentage of OOB (out-of-bag) data used in each tree as 37% of the training data. Performance statistics for the best classifiers are listed in Table 4. ROC curves of binary classifiers are shown in Figure 7. Basing on the precision rate, random forests was chosen as the best strategy. Training and evaluation of the classifiers were performed using the R packages RWeka (Hall et al., 2009; Hornik et al., 2009) and randomForest (Liaw, Wiener, 2002). ROC curves were produced using the R package ROCR (Sing et al., 2005). For the ranking of random forests’ importances, the Gini index was used (Chen, Ishwaran, 2012). 123 Chapter IV Table 4. Classification performance of each learned classifier. Classifiers Precision Accuracy Error rate Sensitivity F-score PP-PANP-NPA Bayesian Network 0.78 0.86 0.14 0.77 0.77 Simple logistic regression 0.82 0.89 0.11 0.81 0.81 Bagging (J48) 0.81 0.88 0.12 0.78 0.80 Bagging (Naive bayes) 0.78 0.87 0.13 0.78 0.78 Bagging (IBk) 0.81 0.88 0.12 0.77 0.79 Random forests 0.90 0.93 0.07 0.87 0.88 PP-NPP Bayesian Network 0.79 0.84 0.16 0.84 0.82 Simple logistic regression 0.87 0.90 0.10 0.86 0.87 Bagging (J48) 0.88 0.91 0.09 0.87 0.88 Bagging (Naive bayes) 0.83 0.88 0.12 0.86 0.85 Bagging (IBk) 0.88 0.89 0.11 0.82 0.85 Random forests 0.94 0.95 0.05 0.92 0.93 PP-PANP-NPA Bayesian Network 0.89 0.89 0.11 0.89 0.89 Simple logistic regression 0.89 0.89 0.11 0.89 0.89 Bagging (J48) 0.91 0.91 0.10 0.91 0.91 Bagging (Naive bayes) 0.90 0.90 0.10 0.90 0.90 Bagging (IBk) 0.90 0.89 0.11 0.89 0.89 Random forests 0.94 0.94 0.06 0.94 0.94 Figure 7. Comparison of ROC curves for the machine-learning schemes used in this study. 124 Chapter IV Supplementary Material Table S1. Bacterial strains used to train the classifiers. Available: http://bacterial-virulence-factors.cbgp.upm.es/supplementary/chapter4/Table_S1.xlsx Table S2. Training data (420x19 matrix). Available: http://bacterial-virulence-factors.cbgp.upm.es/supplementary/chapter4/Table_S2.xlsx Table S3. Plant-associated probabilities assigned by the PA-NPA classifier to the whole set of sequenced bacterial genomes from NCBI. Available: http://bacterial-virulence-factors.cbgp.upm.es/supplementary/chapter4/Table_S3.xlsx File S1. PIFAR User Manual. Available: http://bacterial-virulence-factors.cbgp.upm.es/supplementary/chapter4/userManual.pdf File S2. Tables containing plant-associated probabilities assigned by the PA-NPA classifier separated by genus (compressed folder). Available: http://bacterial-virulence-factors.cbgp.upm.es/supplementary/chapter4/File S2.zip 125 CONCLUSIONES Conclusiones 1. El conjunto de secuencias codificantes predichas en el genoma de Pseudomonas syringae pv. syringae (Pss) UMAF0158 presenta una homología del 91% respecto al del genoma de la cepa filogenétcamente más cercana de la que se conoce su secuencia completa, Pss B728a, patógena de judía. 2. Pss UMAF0158 alberga en su genoma determinantes bacterianos diferenciales respecto a Pss B728a que podrían explicar su capacidad para infectar mango, como el operón mbo de síntesis de mangotoxina, un cluster potencialmente implicado en la producción de celulosa, dos sistemas de secreción diferentes de tipos III y VI, así como un repertorio particular de efectores del sistema de secreción de tipo III. 3. El genoma de Pss B728a contiene regiones no presentes en el genoma de Pss UMAF0158, la mayoría de ellas correspondientes a elementos genéticos móviles. 4. La anotación del genoma de Pseudomonas fluorescens PICF7 reveló genes potenciales que podrían explicar su asociación con olivo, como los codificantes de sideróforos, enzimas detoxificadoras, compuestos volátiles, así como de los sistemas de secreción de tipos III y VI. 5. T346Hunter es una herramienta web, pública y de fácil manejo, que permite identificar regiones genómicas relacionadas con la síntesis de los sistemas de secreción de tipos III, IV y VI. 6. PIFAR (Plant-bacteria Interaction FActors Resource) es una aplicación web de acceso público que permite consultar y descargar una base de datos curada de factores de interacción planta-bacteria, así como la anotación de los mismos en genomas bacterianos. 7. La combinación de T346Hunter, PIFAR y el método de aprendizaje automático random forests ha permitido generar un clasificador de genomas bacterianos asociados a plantas. 8. La aplicación del clasificador descrito en la conclusión 7 sobre aproximadamente 9.500 genomas bacterianos ha revelado posibles asociaciones con plantas de patógenos humanos, particularmente enterobacterias. 129 BIBLIOGRAFÍA