Variance components estimation of complex traits including microbiota information

Abstract

The influence of the microbiome on relevant complex traits for dairy cattle, such as feed efficiency or methane emissions has been well established. Further, recent studies have released evidences on the control of the genetic background of the animal over the microbiota composition. However, until now most analyses have focused on single microorganism approaches instead of the joint microbiome as a whole, including underlying relationships. The joint analysis of the genetic background of the host and its microbiota requires accounting for the distance (or dissimilarity) between communities of microorganisms in different hosts. Therefore, it is necessary to incorporate the whole microbiome into the statistical models to assess its association with complex traits. Microbiome relationship matrix (MRM) allow considering the microbiota as a whole. Several methods have been proposed to ordinate these matrices; those differ on the metric used to account for the distance (or dissimilarities) between microbial communities (e.g. Euclidean, Bray-Curtis, χ2). These distances account for alpha and beta diversity in different ways. Consensus on what method is the most appropriate hasn’t been reached yet, and might depend on data singularities and the purpose of the study. The aim of this study was to compare several microbiota relationship matrices, within a variance component estimation framework. Five ordination methods to build the MRM were tested: metric multidimensional scaling (MDS), detrended correspondence analysis (DCA), non-metric multidimensional scaling (NMDS), redundancy analysis (RDA) and constrained correspondence analysis (CCA). The log transformed and standardized relative abundances matrix described in Ross et al. (2013) was used as a benchmark matrix. Simulated (n=1000) data were used to estimate variance components including phenotypes, genotypes and rumen microbiota information. Data were analysed considering two possible models. First, the genomic effect and the microbiota effect were included independently. Second, an interaction effect between the genomic and microbiota effects was added. All models were implemented within a Bayesian framework using the BGLR package in R. A total of 100 replicates were generated. Real data were analysed using the same models. Similar or slightly better estimation of simulated h2 (0.30) and m2 (0.50) in the independent effects models resulted from ordination methods of MDS (0.307 and 0.493), RDA (0.307 and 0.501) and CCA (0.305 and 0.500) compared to the benchmark MRM 19 (0.304 and 0.480), while poor performance of the DCA (0.249 and 0.349) and NMDS (0.217 and 0.266) methods were obtained at estimating those parameters. The correlation coefficients between genomic estimated breeding values (GEBV) and true breeding values (TBV), from higher to lower, were: the obtained with the benchmark matrix (ρ = 0.633), CCA (ρ = 0.631), RDA (ρ = 0.624), DCA (ρ = 0.598), MDS (ρ =0.592) and NMDS (0.557). Likewise, correlations for predicted microbiota effect in the same order were: the benchmark matrix (ρ = 0.975), CCA (ρ = 0.966), RDA (ρ = 0.949), MDS (ρ = 0.845), DCA (ρ = 0.807) and NMDS (ρ = 0.517). Similar results, in terms of matrices performance, were obtained for the interaction effects model. A real data set (n=70) was also analysed under the same frameworks. Low heritability estimates for feed efficiency (from 0.077 to 0.083) and microbiability (from 0.073 to 0.103) were observed; however, consistent values for the microbiability were obtained with the MRM that performed better in the simulations (from 0.073 to 0.077). Besides, high correlations (ρ > 0.85) between the genetic effect of the host and the phenotypes were obtained for all methods, as well as high correlations between the microbiota effect and the phenotypes for the RDA (ρ = 0.91) and CCA (ρ = 0.91) matrices. Both models were compared using the deviance information criteria (DIC), the effective number of parameters (pD), and the posterior mean of the log likelihood (PostMeanLogLik), resulting in slightly lower values for the independent effects model (DIC: 183.9 to 189.3) than the interaction effects model (DIC: 187.5 to 191.7), those results indicate that it might be a relationship linking genotype-microbiome-phenotype which could be used in prediction of complex traits. The analyses performed in this thesis suggest that canonical ordination methods of RDA and CCA to create MRM are preferred when whole microbiota information is included in the statistical models to analyse complex traits.

La influencia del microbioma sobre los rasgos complejos relevantes para el ganado lechero, tales como la eficiencia alimenticia o las emisiones de metano, ha sido bien establecida. Además, estudios recientes han publicado evidencias sobre el control de los antecedentes genéticos del animal sobre la composición de la microbiota. Sin embargo, hasta ahora la mayoría de los análisis se han centrado en enfoques de microorganismos únicos en lugar del conjunto de microbiomas que incluyen las relaciones subyacentes. El análisis conjunto de los antecedentes genéticos del huésped y su microbiota requiere tener en cuenta la distancia (o disimilitud) entre las comunidades de microorganismos en diferentes hospedadores. Por lo tanto, es necesario incorporar el microbioma completo en los modelos estadísticos para evaluar su asociación con rasgos complejos. La matriz de relación del microbioma (MRM) permite considerar la microbiota como un todo. Se han propuesto varios métodos para ordenar estas matrices; los cuales difieren en la métrica utilizada en la distancia (o disimilitud) entre las comunidades microbianas (por ejemplo, Euclidiana, Bray-Curtis, χ2). Estas distancias representan la diversidad alfa y beta de diferentes maneras. Aún no se ha llegado a un consenso sobre qué método es el más apropiado y podría depender de las singularidades de los datos y del propósito del estudio. El objetivo de este estudio fue comparar varias matrices de relación de microbiota, dentro de un marco de estimación de componentes de varianza. Se probaron cinco métodos de ordenación para construir la MRM: escalado multidimensional métrico (MDS), análisis de correspondencias sin tendencia (DCA), escalado multidimensional no métrico (NMDS), análisis de redundancia (RDA) y análisis de correspondencia restringida (CCA). La matriz de abundancias relativas log-transformada y estandarizada descrita en Ross et al. (2013) se utilizó como matriz de referencia. Se utilizaron datos simulados (n = 1000) para estimar los componentes de la varianza, incluidos los fenotipos, los genotipos y la información de la microbiota del rumen. Los datos fueron analizados considerando dos posibles modelos. Primero, el efecto genómico y el efecto de la microbiota se incluyeron de forma independiente. Segundo, se agregó un efecto de interacción entre los efectos genómico y de microbiota. Todos los modelos se implementaron dentro de un marco Bayesiano utilizando el paquete BGLR en R. Se generaron un total de 100 repeticiones. Los datos reales se analizaron utilizando los mismos modelos. Una estimación similar o ligeramente mejor en la simulación de h2 (0,30) y m2 (0,50) para los modelos de efectos independientes resultó de los métodos de ordenación MDS (0,307 y 0,493), RDA (0,307 y 0,501) y CCA (0,305 y 0,500) en comparación con la MRM de referencia (0,304 y 0,480), mientras que un pobre desempeño de los métodos DCA (0.249 y 0.349) y NMDS (0.217 y 0.266) se obtuvieron al estimar esos parámetros. Los coeficientes de correlación entre los valores genómicos de cría estimados (GEBV) y los valores verdaderos de cría (TBV), de mayor a menor, fueron: los obtenidos con la matriz de referencia (ρ = 0.633), CCA (ρ = 0.631), RDA (ρ = 0.624), DCA (ρ = 0.598), MDS (ρ = 0.592) y NMDS (0.557). Del mismo modo, las correlaciones para el efecto predicho de la microbiota en el mismo orden fueron: la matriz de referencia (ρ = 0.975), CCA (ρ = 0.966), RDA (ρ = 0.949), MDS (ρ = 0.845), DCA (ρ = 0.807) y NMDS (ρ = 0.517). Resultados similares, en términos de rendimiento de matrices, se obtuvieron para el modelo que incluyó interacción entre efectos. Un conjunto de datos reales (n = 70) también se analizó bajo los mismos estándares. Se observaron bajas estimas de heredabilidad para la eficiencia alimenticia (de 0.077 a 0.083) y la microbiabilidad (de 0.073 a 0.103); sin embargo, se obtuvieron valores consistentes para la microbiabilidad con las MRM que se desempeñaron mejor en las simulaciones (de 0,073 a 0,077). Además, se obtuvieron altas correlaciones (ρ > 0.85) entre el efecto genético del huésped y los fenotipos para todos los métodos, así como altas correlaciones entre el efecto de la microbiota y los fenotipos para las matrices RDA (ρ = 0.91) y CCA (ρ = 0.91). Ambos modelos se compararon utilizando los criterios de información de desviación (DIC), número efectivo de parámetros (pD) y la media posterior del logaritmo de la verosimilitud (PostMeanLogLik), resultando en valores ligeramente inferiores para el modelo de efectos independientes (DIC: 183.9 a 189.3) en comparación con el modelo de efectos de interacción (DIC: 187.5 a 191.7), esos resultados indican que podría existir una relación que vincula genotipo-microbioma-fenotipo que podría usarse en la predicción de rasgos complejos. Los análisis realizados en esta tesis sugieren que los métodos de ordenación canónica de RDA y CCA para crear MRM son preferidos cuando la información completa de la microbiota se incluye en los modelos estadísticos para el análisis de rasgos complejos.