Evaluación de 19 genomas recién secuenciados de la sección Flavi

En este estudio, presentamos las secuencias del genoma completo de 19 especies de Aspergillus de la sección Flavi (Fig. 1b). Dos de ellas (A. nomius y A. arachidicola18,19) también fueron publicadas por otros grupos en paralelo a este trabajo. Comparamos estas 19 con especies de la sección Flavi previamente secuenciadas (A. oryzae, A. flavus, A. sojae y A. luteovirescens3,12,13,14) así como con ocho especies de referencia: seis del resto del género Aspergillus más Neurospora crassa y Penicillium digitatum como grupos externos (Fig. 1a, b).

Como primera prueba de base, se comparó la calidad de los ensamblajes del genoma basándose en el tamaño del genoma, el contenido de GC y el número de proteínas predichas (Fig. 1c). Esto mostró una calidad razonable del borrador del genoma con 13 de los 18 genomas ensamblados en menos de 500 andamios (Fig. 1c, columna 5). Una causa de alarma fue A. coremiiformis con 2728 andamios, lo que nos hizo preocuparnos por la calidad del contenido génico. Sin embargo, el genoma cubre el 99,78% de los Benchmarking Universal Single-Copy Orthologs (BUSCO20), y el 96% de los grupos de etiquetas de secuencias expresadas (EST) pueden ser mapeados en el genoma. Por lo tanto, concluimos que la anotación del genoma es de una calidad lo suficientemente alta para las comparaciones del contenido de los genes a pesar del gran número de andamios.

Las especies de la sección Flavi generalmente tienen genomas expandidos

Los tamaños del genoma de Aspergillus sección Flavi son generalmente grandes en comparación con otros Aspergilli representativos (promedio de 37,96 Mbp frente a 31,7 Mbp (Fig. 1c)), como se informó anteriormente para A. oryzae21. Una excepción importante es A. coremiiformis, que tiene tanto menos genes como un genoma notablemente más pequeño, haciéndolo único en la sección.

La filogenia multigénica muestra la compleja herencia de A. oryzae

A continuación examinamos las relaciones evolutivas en la sección Flavi basándonos en una filogenia derivada de 200 genes (Fig. 1a). El soporte de la ramificación dentro del árbol es alto (100 de 100 bootstraps en la mayoría de las ramas). El árbol confirma que la sección Flavi es un grupo monofilético. Los clados de la Fig. 1a corresponden a un árbol filogenético previamente reportado basado en el gen de la beta-tubulina10,11,22 y las distancias entre secciones corresponden a trabajos anteriores23.

Un posible error en el árbol es que A. sojae se encuentra más cerca de A. flavus, ya que A. sojae se percibe como una versión domesticada de A. parasiticus. De hecho, esta ramificación también tiene el valor bootstrap más bajo del árbol. La explicación más probable es que, dado que las predicciones de los genes de A. sojae se basan en las anotaciones del genoma de A. flavus y A. oryzae24,25, se crea un sesgo en los genes predichos y este sesgo se refleja probablemente en el árbol. Como prueba, hemos generado árboles filogenéticos utilizando métodos alternativos que no dependen de la anotación genética (CVTree26,27). Estos muestran claramente que A. sojae es el más cercano a A. parasiticus, tanto cuando se utilizan secuencias del genoma completo como del proteoma (Fig. 1 y Fig. 2 suplementarias). Por lo tanto, pensamos que A. sojae debería colocarse junto a A. parasiticus en el árbol filogenético como la flecha indicada en la Fig. 1a.

Además, A. oryzae, percibido como una versión domesticada de A. flavus10,28,29,30, no está directamente junto a él en el árbol. Sin embargo, se ha sugerido previamente que A. oryzae desciende de un ancestro que fue el antecesor de A. minisclerotigenes o A. aflatoxiformans31. La filogenia (Fig. 1a, zoom) apoya esta sugerencia, mostrando que A. minisclerotigenes y A. aflatoxiformans son parientes más cercanos de A. oryzae que de A. flavus.

El análisis de las proteínas compartidas confirma la alta diversidad genética

Con el fin de examinar las características centrales compartidas por todas las especies de la sección Flavi, los clados, así como las características de las especies individuales, hicimos un análisis de los genes homólogos compartidos dentro y a través de las especies16, y los clasificamos en familias de proteínas homólogas (Fig. 2). Esto permitió identificar (1) Las familias de proteínas del genoma central con al menos un miembro en todas las especies comparadas. Se espera que esto cubra las proteínas esenciales. (2) Genes específicos de sección y de clado: genes que tienen homólogos en todos los miembros de un clado/sección, pero no en ninguna otra especie. (3) Genes específicos de especie-genes que no tienen homólogos en ninguna otra especie en la comparación.

El núcleo del genoma de las 31 especies en este conjunto de datos es de 2082 familias de proteínas. Para las 29 especies de Aspergillus este número es de 3853, y para la sección de especies de Flavi solo constituye 4903 familias de proteínas. Así, más de la mitad del genoma de las especies de la sección Flavi varía a través de las especies.

Examinando las familias de proteínas específicas del clado, sólo se encuentran muy pocas (27-54) (Fig. 2a), lo que es bajo comparado con la sección Nigri examinada previamente16. Como las secciones Nigri y Flavi son aproximadamente igual de ricas en especies, esto podría indicar que las especies de la sección Flavi son más distintas. Esto se ve apoyado por el hecho de que el número de genes específicos de las especies es muy alto (166-2181), donde vemos que 166 (A. sojae) es un número artificialmente bajo, debido a que la llamada de genes en este genoma se basa en los genomas de A. flavus y A. oryzae.

Los genes específicos de cada especie suelen codificar la regulación y los P450

Queríamos ver si los genes específicos de cada especie podían estar relacionados con funciones conocidas de Flavi, como la fermentación de alimentos y la patogenicidad en plantas y humanos. Para ello, examinamos las funciones previstas de los genes específicos de la especie utilizando las anotaciones de InterPro, GO y KOG32,33,34,35. La porción con una anotación funcional fue baja; 20, 12 y 9% para InterPro, GO y KOG, respectivamente; en total el 21% tenía una anotación (Figs. Suplementarias 3-5). Este es un porcentaje muy alto-pero no inusual-de funciones no identificables.

Nos centraremos en InterPro ya que cubre más genes: las funciones más comunes de InterPro incluyen factores de transcripción, proteínas quinasas, transportadores y P450s (Fig. Suplementaria 3), que también están significativamente sobrerrepresentados. Aunque estos rasgos no pueden vincularse directamente con la fermentación de alimentos y la patogenicidad, la regulación está implicada en la adaptación y los P450s desempeñan funciones tanto en la degradación de sustratos como en la producción de compuestos bioactivos, ambos relevantes para la patogenicidad de los hongos.

Los genes de las especies están sobrerrepresentados en las regiones subteloméricas

Se ha demostrado que las secuencias subteloméricas son regiones ampliamente reordenadas en A. nidulans, A. oryzae y A. fumigatus21. Esto también se observa en mamíferos, nematodos y levaduras36. Estudios anteriores37,38 mostraron que las regiones sub-teloméricas tienen un sesgo de genes únicos, divergentes o ausentes. Otro estudio ha mostrado que los grupos de genes de metabolitos secundarios (SMGCs) están enriquecidos en regiones sub-teloméricas en A. nidulans y A. fumigatus21.

Por lo tanto, examinamos la densidad genética y la localización de los genes específicos de la especie, los grupos de metabolitos secundarios y el genoma central, utilizando el genoma de A. oryzae como referencia para evaluar la posible sobrerrepresentación de estos genes en las regiones subteloméricas (Fig. 3).

Tanto la inspección visual como la prueba exacta de Fisher confirmaron que tanto la especie específica (valor p = 7.266e-07) como los SMGC (valor p < 2,2e-16) están enriquecidos hacia las regiones subteloméricas (100 kbp desde los extremos cromosómicos), mientras que los genes centrales se encuentran con menos frecuencia en las regiones subteloméricas. El hecho de que los genes específicos de una especie no estén distribuidos de forma aleatoria, desaconseja que se trate simplemente de errores de anotación o de modelización de genes, indicando por tanto que se trata, efectivamente, de genes legítimos. La distribución de los genes específicos de la especie sugiere que los nuevos genes se incorporan con más frecuencia a las regiones subteloméricas que a otros lugares. Si esto es el resultado de una selección para la región subtelomérica, o una contra-selección contra otras regiones, o ambas cosas, los datos no lo revelan.

El análisis de sintenia revela islas de contenido génico muy variable

Las regiones sintéticas y no sintéticas son otro factor a considerar cuando se analiza la localización del genoma. Se ha demostrado que el genoma de A. oryzae tiene un patrón de mosaico de regiones sinténicas y no sinténicas en relación con los Aspergilli lejanamente relacionados1,2. Examinamos la sintenia a través de la sección Flavi y en A. nidulans y A. fumigatus utilizando A. oryzae RIB40 como referencia (Tabla 1). Este análisis apoya nuestro hallazgo anterior de que A. oryzae está estrechamente relacionado con A. aflatoxiformans que con A. flavus.

Una visión general de los genes sintéticos compartidos se ilustra en la Fig. 6 suplementaria. En general, hay menos regiones de sintenia hacia los extremos teloméricos como se vio previamente1,2 en una comparación de A. nidulans, A. fumigatus y A. oryzae. Además, observamos que los cromosomas 1 y 2 tienen un grado muy alto de conservación de la sintenia, mientras que los cromosomas 6 y 8 tienen una conservación mucho menor de la sintenia.

Encontramos densas islas de genes no sintéticos en regiones no subteloméricas en los cromosomas 4, 6 y 8. Estas islas podrían ser causadas por la transferencia horizontal de genes (HGT), por el barajado de genes o por la formación de genes de novo. Hemos investigado la presencia de HGT utilizando BLASTp para examinar las mejores coincidencias en la base de datos no redundante del NCBI. Se espera que los HGT recientes tengan una alta identidad de secuencia con otro grupo de especies de las que se habría transferido, y que no se encuentren en las especies estrechamente relacionadas39. Ninguna de estas islas mostró signos de HGTs recientes. Además, sólo 23 de los 80 genes de los bloques no sintéticos eran específicos de A. oryzae. Por lo tanto, parece probable que estas islas no sintéticas sean causadas por una mezcla de reordenamientos significativos, eventos de duplicación y la aparición de genes específicos de A. oryzae.

Tomado en conjunto, el hecho de que observemos algunos cromosomas muy conservados y algunos bloques no sintéticos altamente reordenados podría indicar una presión evolutiva para la estabilidad en algunas regiones mientras que otras regiones están frecuentemente sujetas a barajados y reordenamientos de genes, es decir, puntos calientes de reordenación.

La sección Flavi es una rica fuente de enzimas activas en carbohidratos

Las enZimas activas en carbohidratos (CAZymes) son esenciales para las fuentes de carbono que una especie puede degradar y utilizar. Dentro de la sección Flavi la utilización de CAZimas/carbono se describe principalmente para A. oryzae1,2,40 y en menor medida para A. flavus41,42,43,44,45 y A. sojae46,47, mientras que sólo se han realizado estudios incidentales con otras especies de este grupo48,49,50,51,52,53,54, a menudo describiendo la producción o caracterización de una determinada actividad o proteína CAZyme, respectivamente.

Utilizamos la base de datos CAZy para predecir el contenido de CAZyme en los genomas de la sección (Fig. 4). Se predijo un total de 13.759 CAZimas para las 23 especies de Flavi (una media de 598/especie). Esto es bastante rico comparado con la referencia incluida Aspergilli (508/especie).

De este análisis se desprende que hay una clara diferencia entre los clados de la sección Flavi (Fig. 4b), mostrando de nuevo una variación en el contenido de genes en la sección.

El contenido variable de CAZyme no refleja la capacidad de degradar la biomasa vegetal

Para evaluar la capacidad real de utilización del carbono en toda la sección Flavi, realizamos un perfil de crecimiento de 31 especies (29 Aspergilli, incluyendo 23 especies de la sección Flavi) en 35 sustratos relacionados con la biomasa vegetal (Fig. 5, Datos suplementarios 1) y lo comparamos con la predicción del contenido del gen CAZyme relacionado con la degradación de la biomasa vegetal (Datos suplementarios 2). En un estudio anterior, la variación en el crecimiento entre Aspergilli relacionados a distancia podría estar relacionada con las diferencias en el contenido del gen CAZyme55, pero este no fue el caso de las especies relacionadas más estrechamente de la sección Nigri de Aspergillus16.

La glucosa dio lugar al mejor crecimiento de todos los monosacáridos para todas las especies y, por tanto, se utilizó como referencia interna para el crecimiento (Fig. suplementaria 7). El crecimiento en otras fuentes de carbono se comparó con el crecimiento en d-glucosa y esta diferencia relativa se comparó entre las especies. El crecimiento sobre monosacáridos fue en gran medida similar entre las especies de la sección Flavi (Fig. 5, Fig. Suplementaria 7 y Datos Suplementarios 1).

Los conjuntos de CAZyme relacionados con la degradación de la biomasa vegetal son en general muy similares para la sección Flavi (Fig. 5), con la excepción de A. coremiiformis, que tiene un conjunto de genes muy reducido. Esto se debe principalmente a la reducción de las familias de glucósidos hidrolasas, pero también a un número de familias relacionadas con la degradación de la pectina, el xilano y el xiloglucano. Sorprendentemente, esta especie mostró un mejor crecimiento relativo sobre el xilano que la mayoría de las otras especies, mientras que el crecimiento sobre otros polisacáridos fue principalmente similar al de la sección Flavi. Por lo tanto, el conjunto reducido de genes no ha reducido su capacidad para degradar la biomasa vegetal. Esto podría ser similar al caso de T. reesei, que también tiene un conjunto reducido de genes CAZyme, pero produce las enzimas correspondientes a niveles muy altos56. Sin embargo, el origen de este enfoque es probablemente muy diferente, ya que su contenido de CAZyme fue formado por la pérdida y luego la ganancia masiva de HGT de las enzimas que degradan la pared celular de las plantas57, mientras que no hay indicios de esto para A. coremiiformis.

Las diferencias hidrolíticas son específicas del clado dentro de la sección Flavi (Datos Suplementarios 2). El clado A. togoensis tiene un conjunto reducido de genes xilanolíticos y xiloglucanolíticos, pero esto no se refleja en el crecimiento. Por el contrario, los genes GH115 (alfa-glucuronidasa) están expandidos en los clados A. flavus, A. tamarii y A. nomius (se ha informado de enzimas o actividad xilanolítica en varias especies de estos clados49,50,51,53,58,59,60,61,62), GH62 (arabinoxilano arabinofuranohidrolasa) se expandió en el clado A. leporis, y los clados A. leporis y A. avenaceus fueron los únicos clados con CE15 (glucuronoil esterasas), que también se encontraron en especies de Aspergillus fuera de la sección Flavi.

La capacidad de degradación del galactomanano se conservó casi por completo en la sección Flavi, pero curiosamente el crecimiento en goma guar que consiste principalmente en galactomanano fue variable entre las especies. Del mismo modo, la capacidad amilolítica reducida de los clados A. togoensis y A. avenaceus no dio lugar a un crecimiento reducido sobre el almidón o la maltosa.

Se observaron variaciones en el número de genes pectinolíticos. Las diferencias más pronunciadas fueron la ausencia de genes PL11 (ramnogalacturonano liasa) en la mayoría de las especies de la sección Flavi, y la expansión de GH78 (alfa-ramnosidasa) en los clados A. flavus y A. tamarii. Sin embargo, estas diferencias y las más pequeñas en otras familias no dieron lugar a una gran variación en el crecimiento sobre la pectina.

Diferencias más obvias se presentaron durante el crecimiento sobre la celobiosa, la lactosa y la lignina. La mayoría de las especies crecieron poco en celobiosa a pesar de que el número de genes que codifican la beta-glucosidasa es similar en la mayoría de las especies (Datos suplementarios 2). Del mismo modo, sólo A. arachidicola y, en menor medida, A. albertensis crecieron bien en lactosa, mientras que el número de beta-galactosidasas en estas especies es similar al de las otras especies. Lo más interesante fue el hallazgo de que A. albertensis creció tan bien en lignina como en d-glucosa, sugiriendo aplicaciones potenciales en la producción de biocombustible.

En resumen, el potencial CAZyme en la sección Flavi se conserva en gran medida (con la excepción de A. coremiiformis) con algunas variaciones en el número de copias, pero el potencial genómico y las variaciones no se reflejan necesariamente en el crecimiento. Por lo tanto, es probable que, como se sugirió anteriormente55, las diferencias observadas sean en gran medida a nivel de regulación.

La familia de CAZimas GH28 está inflada en el clado A. flavus

Estábamos particularmente interesados en las CAZimas GH28, ya que son importantes para la fermentación de alimentos y la calidad del producto final fermentado63. Se creó un árbol filogenético de todos los miembros de GH28 de la sección Flavi (Fig. 8 suplementaria). El árbol consta de 429 proteínas, una media de 18,7 por especie.

Dentro del árbol hay diferentes agrupaciones. Cinco grupos tienen miembros de las 23 especies, nueve grupos carecen de una a cuatro especies (normalmente A. coremiiformis y A. caelatus), y dos grupos son específicos de los clados A. flavus, A. tamarii y A. nomius. Por último, hay ocho grupos que contienen de 2 a 13 especies, que no siguen la filogenia, lo que sugiere que son fuentes de variación de GH28.

En general, las especies del clado A. flavus tienen un alto número de miembros de GH28. Se sabe que A. sojae tiene un alto número de GH28, que también se ve aquí con 24 miembros; sin embargo, A. sergii tiene un número aún mayor con 25 miembros. Podría ser interesante investigar si esto podría explotarse utilizando A. sergii como una nueva especie en la fermentación de alimentos y/o como fuente de nuevas enzimas.

Análisis del metabolismo secundario

Se sabe que el género Aspergillus produce un gran número de SMs y el número de SMGCs predichos es aún mayor. La mayoría de los SMGCs predichos no están caracterizados y por lo tanto tienen el potencial de producir una diversidad de compuestos bioactivos novedosos. Examinamos la diversidad y el potencial de producción de SM en la sección Flavi, tanto cuantitativamente en términos de número de clusters, como cualitativamente en términos de los compuestos que estos clusters podrían producir potencialmente.

El metabolismo secundario en la sección Flavi es diverso y prolífico

Para evaluar cuantitativamente el potencial de producción de SM, se predijeron los SMGCs utilizando una herramienta de predicción similar a SMURF64 para todas las especies excepto N. crassa y A. sojae, ya que éstas fueron secuenciadas por otros métodos y con métodos de llamada de genes diferentes (Fig. 6c). Dentro de las 28 especies de Aspergillus, hay un total de 1972 SMGCs predichos y para los genomas de la sección Flavi, el total es de 1606 SMGCs (73/especie). Esto supone más de 15 extra por especie en comparación con el muy prolífico género Penicillium65.

Queríamos examinar cómo de únicas son las SMGCs, y por ello construimos familias de SMGCs (Datos Suplementarios 3). Para todo el conjunto de datos, pudimos dividirlo en 477 familias de SMGC, y para la sección Flavi 308 familias de SMGC. De éstas, 150 clusters de SMGCs sólo se encuentran en una especie de la sección Flavi (Fig. 6a), mostrando un gran número de clusters únicos en cada especie (6,8 SMGCs únicos/especie). En comparación con Aspergillus sección Nigri, el número de clusters por especie en este estudio es ligeramente inferior, pero el número de miembros en cada familia de SMGC también es menor, lo que demuestra una mayor diversidad en el metabolismo secundario en la sección Flavi en comparación con la sección Nigri.

La desreplicación del metabolismo secundario predice productores de toxinas

Para evaluar cualitativamente el potencial de producción de SM, utilizamos una línea de «desreplicación genética» en la que los clusters predichos se asocian con clusters caracterizados verificados (de la base de datos MIBiG66) en un método de culpabilidad por asociación67. En base a esto, 20 familias de clusters fueron acopladas a una familia compuesta (Fig. 6b). Algunas familias de clústeres se encontraron en todos o casi todos los genomas de Flavi, por ejemplo, los similares a los clústeres de naftopirona68, nidulanina A69, azanigerona70, 4,4′-piperazina-2,5-diildimetil-bis-fenol y aflavarina71/endocrocina72,73. La mayoría de las familias suelen seguir los grupos filogenéticos, lo que sugiere un patrón de distribución basado en la pérdida, pero algunas, como las familias de SMGC similares a los grupos de asperfuranona74, pseurotina A75 o fumagilina76 no siguieron la filogenia. Además, se identificaron productores potenciales de toxinas conocidas como la aflatoxina y la aspiroclorina (Fig. 6b).

La combinación de datos y análisis vincula un compuesto a un clúster

A partir de los clústeres de SMGC conocidos, nos interesaba vincular compuestos y clústeres basándonos en el patrón de presencia/ausencia de compuestos producidos y clústeres predichos. Por lo tanto, creamos un mapa de calor de todas las familias de clusters encontradas en al menos cinco especies, añadimos las familias de compuestos predichas a partir de la desreplicación del MIBiG, además de las familias de compuestos curadas manualmente a partir de un estudio de la literatura (Fig. 9 suplementaria). Además de esto, medimos la producción de SM de las especies de Flavi (Datos Suplementarios 4).

De particular interés fueron las miyakamidas. Se aislaron originalmente de un aislado de A. flavus y se demostró que tenían propiedades antibióticas77, pero no se conoce el grupo de genes biosintéticos. Nuestro análisis químico mostró la producción en A. sojae, A. nomius, A. parasiticus, A. novoparasiticus y A. transmontanensis.

Realizamos la retro-biosíntesis a partir de la estructura química y predijimos que el clúster de genes biosintéticos debería contener una sintetasa peptídica no ribosomal (NRPS) con 2-3 dominios de adenilación (ya que dos de los tres aminoácidos son similares), una N-metiltransferasa, una acetiltransferasa y, potencialmente, una descarboxilasa/deshidrogenasa (Fig. suplementaria 10A). Buscando familias de grupos con miembros en todas las especies productoras de miyakamida que tuvieran columnas vertebrales de NRPS con 2-3 dominios de adenilación y un dominio de metiltransferasa, sólo una familia de grupos cumplía los requisitos. La familia de grupos tiene una columna vertebral NRPS con un dominio de metiltransferasa, tres dominios A en la mayoría de las especies, y dos en A. novoparasiticus. La predicción de sólo dos dominios A se debe probablemente a un error de anotación, ya que la similitud de la secuencia se conserva antes del inicio del gen (Fig. Suplementaria 10B). El tamaño del cluster predicho es de 1-9 genes, la diferencia es probablemente causada por errores de predicción del SMGC (Gráfico de sintenia en la Fig. Suplementaria 10B). El gráfico de sintenia muestra que el NRPS y dos pequeños genes con función desconocida están ampliamente conservados. Por lo tanto, proponemos que el NRPS identificado junto con los dos genes conservados de función desconocida son probables candidatos para la biosíntesis de la miyakamida.

El clúster de genes biosintéticos de la aflatoxina está muy conservado

Tal vez el metabolito secundario más conocido de la sección Flavi es la aflatoxina, altamente cancerígena. Se sabe que las aflatoxinas son producidas por muchas especies de la sección Flavi (A. arachidicola, A. luteovirescens, A. flavus, A. minisclerotigenes, A. nomius, A. aflatoxiformans, A. pseudocaelatus, A. pseudonomius, A. pseudotamarii y algunos aislados de A. oryzae)4,10.

El análisis de desreplicación (Fig. 6b) identificó una familia de SMGC que se predice que está implicada en la producción de esterigmatocistina y aflatoxina, que son todas las especies de los clados de A. flavus, A. nomius y A. tamarii, excepto A. tamarii. Un gráfico de sintenia de la familia SMGC (Fig. 11 suplementaria) muestra que el clado está extremadamente bien conservado, sin reordenamientos y con una alta identidad de alineación para los genes de aflatoxina. Sólo A. caelatus tiene una forma truncada con sólo los genes aflB, aflC y aflD y A. tamarii parece tener una pérdida completa del cluster. Curiosamente, la mayoría de los clusters predichos no incluían los genes aflP y aflQ que son responsables del último paso de la biosíntesis de las aflatoxinas. Buscamos en los genomas la aflP (Fig. 12 suplementaria), y la encontramos en todos los genomas, pero con diferentes sitios de inicio y una secuencia extra en el centro de las proteínas. Los datos de RNA-seq apoyan estos modelos (Supplementary Fig. 13) y sugieren errores en los modelos de genes de A. flavus. Del mismo modo, el gen aflQ se encuentra en todas las demás especies, pero a 5-10 genes de distancia de los grupos predichos. Así, el análisis detallado muestra que todas estas especies tienen los genes necesarios para la biosíntesis de aflatoxinas.