Recuperación, anotación y destilación del genoma
Navegamos en la base de datos GTDB [9], que contiene genomas bacterianos completos y en borrador con puntajes de completitud CheckM asociados, para filas bacterianas con al menos 100 genomas en cada una de las ventanas del 1 % que van del 70 al 100 % de compleción del genoma con <10 % de contaminación. Estos criterios fueron cumplidos por cuatro filos, a saber Actinobacteria, Bacteroidota, Firmicutes (sensu latoincluido bacilos y clostridios) y proteobacteria, que fueron considerados para el análisis. Los genomas dentro de cada ventana de completitud y phylum se seleccionaron aleatoriamente y sus secuencias se recuperaron de la base de datos NCBI. [14] usando sus códigos de acceso a la asamblea. Para el análisis del genoma completo submuestreado, se utilizó un estándar comunitario microbiano ZymoBIOMICS® secuenciado de Olm et al. [15] se descargó usando la descarga de kingfisher (https://github.com/wwood/kingfisher-download). Esta comunidad simulada contiene 8 especies bacterianas que tienen genomas de referencia completos disponibles. Las lecturas descargadas se submuestrearon aleatoriamente a diferentes profundidades (1 millón, 2 millones, 5 millones y TODAS) utilizando la muestra seqtk con una semilla de 1337 (https://github.com/lh3/seqtk), y cada uno procesado a través del ‘ Tubería Individual_Assembly_Binning.snakefile’ (https://github.com/earthhologenoma/EHI_bioinformatics). dRep [16] luego se usó en todos los MAG resultantes con la inclusión de las referencias del genoma bacteriano completo para determinar qué MAG era de qué genoma bacteriano. Un MAG representativo de cada especie (completitud de CheckM >90, contaminación <5, mínimo al menos 42 contigs) tuvo contigs submuestreados aleatoriamente diez veces a tres tasas diferentes (70, 80, 90% de contigs restantes) utilizando reformat.sh de BBMap (total de 240 submuestras MAG). Luego se ejecutó CheckM en estos MAG submuestreados para estimar la integridad.
Posteriormente, los genomas se anotaron y destilaron a los valores de plenitud de la ruta KEGG (definidos como la proporción de reacciones bioquímicas habilitadas por los genes presentes en un genoma para lograr una función metabólica) usando DRAM (1.2.4) [2]. El script resume_genomes.py se modificó para mostrar todos los módulos y números de módulos (https://github.com/EisenRa/DRAM_more_modules/tree/1.2.4_more_modules).
Estadísticas y visualización
Se realizaron análisis estadísticos sobre los datos de integridad del módulo KEGG. Solo los módulos KEGG generalizados presentes en al menos el 5% de los MAG se utilizaron para el modelado estadístico. Este filtrado dio como resultado 195 módulos KEGG observados en 11 842 MAG. Los análisis estadísticos se realizaron en tres pasos consecutivos.
En el primer paso, se utilizaron modelos lineales generalizados para estimar la relación entre la integridad de los módulos KEGG y la integridad de los genomas. Se utilizó una distribución binomial con la función de enlace logit, ya que la función plenitud representa la proporción de reacciones enzimáticas (o pasos) de un módulo presente en un genoma. El número total de pasos de cada módulo se utilizó como peso en los modelos. La completitud del genoma (variable numérica), el Phylum bacteriano (variable categórica con cuatro niveles) y su interacción se utilizaron como variables explicativas fijas, por lo que se estimó una pendiente específica de Phylum para cada módulo.
En el segundo paso, modelado lineal de efectos mixtos, tal como se implementa en el paquete R lme4 [17], se utilizó para explorar predictores de la fuerza de la relación plenitud-integridad entre módulos. Se utilizaron como variables de respuesta las pendientes específicas del Phylum estimadas en los modelos binomiales antes mencionados y el Phylum bacteriano (variable categórica con cuatro niveles), el dominio KEGG de cada módulo (variable categórica con diez niveles) y el número de pasos involucrados en un módulo ( variable numérica) se utilizaron como variables explicativas fijas. Dado que se estimaron cuatro pendientes para cada módulo (una para cada Phylum bacteriano) y se incluyeron en la variable de respuesta para este modelo, se incluyó un efecto aleatorio a nivel de módulo en el modelo como intercepto aleatorio (aleatorio = 1|módulo). Construimos intervalos de confianza bootstrap alrededor de los niveles de las variables categóricas a través de la función bootMER() con 999 simulaciones, y alrededor de la pendiente de la variable numérica usando la función confint.merMod(), ambas incluidas en el paquete lme4. Para realizar las predicciones marginales de cada variable categórica, las variables categóricas no focales (Phylum o dominio) se mantuvieron en su nivel de referencia y la variable numérica (pasos) en su valor medio. Los intervalos de confianza del 95% no superpuestos entre los niveles de los factores categóricos Phylum y el dominio se consideraron como evidencia contra la hipótesis nula de ausencia de diferencias entre los grupos. De igual manera, para la variable numérica número de pasos, un intervalo de confianza de la pendiente que no se superponga a cero se consideró como evidencia en contra de la hipótesis nula de no asociación.
En el tercer paso, primero usamos los perfiles funcionales de 8 genomas de una comunidad simulada y sus réplicas submuestreadas para realizar un Análisis de coordenadas principales (PCoA) para la evaluación visual del sesgo introducido por la reducción de la integridad del genoma a alrededor del 90%, 80 y 70%. Luego, usamos los modelos lineales generalizados binomiales específicos del módulo construidos en el paso uno (entrenados con los MAG recuperados de la base de datos GTDB) para corregir los perfiles funcionales de los genomas submuestreados. Para esto, generamos dos predicciones para cada genoma focal en función de su Phylum: la plenitud del módulo pronosticada dada su completitud observada y la plenitud del módulo pronosticada si el genoma estaba 100% completo. La diferencia entre ambas predicciones se sumó a la plenitud del módulo observado en cada genoma focal. Si la plenitud del módulo corregido era mayor que 1, se redondeaba a 1. Por último, se realizó una PCoA conjunta utilizando los perfiles funcionales de los genomas completos, los genomas submuestreados sin procesar y los genomas submuestreados corregidos, y se mostró en seis gráficos con fines de visualización. . .