Skip to main content

Heritability and recursive influence of host genetics on the rumen microbiota drive body weight variance in male Hu sheep lambs



Heritable rumen microbiota is an important modulator of ruminant growth performance. However, no information exists to date on host genetics-rumen microbiota interactions and their association with phenotype in sheep. To solve this, we curated and analyzed whole-genome resequencing genotypes, 16S rumen-microbiota data, and longitudinal body weight (BW) phenotypes from 1150 sheep.


A variance component model indicated significant heritability of rumen microbial community diversity. Genome-wide association studies (GWAS) using microbial features as traits identified 411 loci-taxon significant associations (P < 10−8). We found a heritability of 39% for 180-day-old BW, while also the rumen microbiota likely played a significant role, explaining that 20% of the phenotypic variation. Microbiota-wide association studies (MWAS) and GWAS identified four marker genera (Bonferroni corrected P < 0.05) and five novel genetic variants (P < 10−8) that were significantly associated with BW. Integrative analysis identified the mediating role of marker genera in genotype influencing phenotype and unravelled that the same genetic markers have direct and indirect effects on sheep weight.


This study reveals a reciprocal interplay among host genetic variations, the rumen microbiota and the body weight traits of sheep. The information obtained provide insights into the diverse microbiota characteristics of rumen and may help in designing precision microbiota management strategies for controlling and manipulating sheep rumen microbiota to increase productivity.

Video Abstract


Sheep weight is the most important growth indicator used in production and is directly related to meat production [1], fat deposition [2], and reproductive performances [3]. Global meat production and consumption continue to grow due to the increased demand driven by population growth, individual economic growth, and urbanization [4]. It is estimated that global meat production will increase 76% by 2050 to meet the increased demand, with half of the global demand for ruminants coming from developing countries, particularly China, the largest producer and importer of sheep meat [4,5,6]. This thus poses new challenges to increase the sheep meat production. At the same time, the interest of the consumer for low-fat and healthy meat is increasing due to the negative impact of meat fat on human health [7]; however, in sheep, higher weight is often associated with excessive fat deposition [2]. In light of this, to balance yield and healthiness in sheep meat production, there is an urgent need to comprehensively elucidate the mechanisms that control the body weight (BW) and to develop new breeding intervention strategies.

Although classical genetics studies suggested heritability of approximately 30 to 60% for sheep BW [8,9,10,11,12], current genetic breeding and selection for BW have been only partially successful because the phenotype represents a complex multifactorial trait. In addition, new research in genome-wide association analysis (GWAS) of sheep weight was limited so far due to the small sample sizes, the relative dispersion of SNP chips, and investigations conducted in complex non-laboratory settings such as grazing [13,14,15]. Therefore, it is necessary to conduct larger sample size GWAS and whole-genome sequencing studies. Furthermore, in addition to breeding and selection, body weight may also be controlled through the regulation of the gut microbiota. In recent years, several studies have supported the role of the commensal microbiota in the gastrointestinal tract in the regulation of the sheep metabolism [16] and immune response [17], and therefore, the hypothesis to influence the sheep BW through the manipulation of the gut microbiota has grown. In this regard, the rumen, which is primarily responsible for digestive and nutrient absorptive functions in sheep, is colonized by a highly complex and anaerobic microbial ecosystem that is able to convert low-nutrient plant material through fermentation into essential and readily metabolites, accounting for up to 70% of host energy requirements [18, 19].

Since the research studies were conducted by Pomp’s team [20], there is a growing evidence that in mammals [21, 22], including humans [23,24,25], and in agriculturally relevant avian hosts [26, 27], the host genetic variation may affect the composition and structure of their gut microbiota. Of particular note, empirical studies on ruminants were focused only on cattle [28,29,30,31,32,33,34], while an association between rumen microbiota and the sheep genome has not been demonstrated yet, further leading to a lack of empirical evidence on whether host genetic influences on the rumen microbiota can also affect sheep phenotypes, such as the traits related to BW.

We therefore hypothesized that host genetics may influence the rumen microbial community in sheep, and that the same host single-nucleotide polymorphism (SNP) marker could similarly influence sheep weight through a direct mechanism and an indirect one mediated by marker microbiota. To test these hypotheses, we performed whole-genome re-sequencing and rumen microbiota 16S ribosomal RNA (16S rRNA) amplicon sequencing from a cohort of 1150 sheep of the same sex, age, and breed, living on the same farm, and raised on the same diet. We characterized the composition of rumen microbiota, and we estimated the heritability. Additionally, the GWAS for rumen microbial features (microbiota GWAS, mbGWAS) was conducted to explore the relationship between host genetics and the rumen microbiota in sheep. The host additive genetics and rumen microbial effects were investigated in relation to the sheep BW by estimating heritability (h2) and microbiability (m2, the proportion of the contribution of gut microbiota on host phenotypes), identifying the genetic SNP markers and rumen marker microbiotas significantly associated with BW by GWAS and microbiota-wide association studies (MWAS). Finally, the recursive influence of host genomic variations on weight was identified by combining mbGWAS, GWAS, and MWAS results.

Material and methods

Animals, phenotypic data and sample collection

This study included 1150 male Hu lambs reared on the same farm and divided into four batches over a 24-month period to investigate animal performance and collect the samples (Fig. 1 and Table S1). Specifically, just after weaning, all lambs were transferred to the Minqin experimental farm of Lanzhou University (N38°43′41′′, E103°013′), and they were housed in individual feeding pens (0.8 m × 1.0 m, l × w) and fed on the same diet (Supplementary Table 1). For each individual, BW traits were measured with a calibrated livestock scales (at 06:00–08:00 am, before morning feeding) every 20 days, from the age of 100 days to 180 days. At the age of 180 days and after fasting for 12 h, whole blood samples were obtained via the jugular vein and stored at −20 °C until use; each animal was slaughtered according to standard commercial procedures following the requirements of the China Council on Animal Care, and whole rumen content samples were collected immediately and stored at −80 °C. The concentration of VFAs was measured using TRACE-1300 series GC ultra-gas chromatograph (Thermo Scientific, Milan, Italy) by following standard procedures.

Fig. 1
figure 1

Study design and workflow. This schematic representation highlights, for each step, the research question that we sought to answer, the analysis workflow, the data used, and the generalized result

Extraction of host and gut microbial DNA

Host genomic DNA was extracted from blood sample using the whole blood genomic DNA rapid extraction kit (EasyPure Blood Genomic DNA Kit; TIANGEN Bio Company, Beijing, China). All rumen content samples were separately thawed and homogenized again on ice, and then total microbial DNA was extracted from ~200 mg of each sample using EasyPure Stool Genomic DNA Kit (TransGen Biotech, EE301-01, Beijing, China) according to the manufacturers’ instructions. The quality of DNA was assessed on 1% agarose gel electrophoresis. All samples were included in the study, resulting in a total of 1150 host DNA and 1150 rumen microbial DNA samples.

16S rRNA gene sequencing and analysis

The V3–V4 regions of bacterial 16S rRNA gene were amplified using specific barcoded primers (341F: CCTAYGGGRBGCASCAG and 806R: GGACTACNNGGGTATCTAAT). Amplicons were sequenced NovaSeq PE250 platform of Illumina (Novogene Biotech Co., Ltd, Beijing, China). Raw sequences were assigned to samples based on their unique barcodes and then trimmed to remove the barcode and primer sequences. The pair-end reads of each sample were assembled using FLASH [35]. The clean sequences underwent quality control analysis using FastQC software (FastQC, and chimeric sequence removal using UCHIME ( The filtered data were further processed via DADA2 method [36] in QIIME2 ( to produce tables of amplicon sequence variants (ASVs), while the taxonomic assignment was performed using QIIME2 classify-sklearn algorithm by a pre-trained Naive Bayes classifier on 16S Silva database (version 138). To avoid the interference of contingent opportunistic factors and low abundance feature sequence/ASV, the table were filtered using the QIIME2 feature-table filter-features commands (--p-min-samples 3 and --p-min-frequency 5), ensuring that each feature sequence was present in at least 3 samples, and that the total sequencing depth for each feature sequence was greater than 5 reads. Prior to further downstream processing, the data were reduced to the minimum library size (QIIME2 feature-table rarefy commands: -p-sampling-depth 34,736) to obtain the final ASV count data for six taxonomic levels (from phylum to species). Finally, we obtained a preliminary ASV tables of rumen microbes containing 1150 individuals and 11,976 ASVs. These data generated 813 genera. The count data were further normalized to relative abundance by the total-sum scaling method. To further minimize the effect of background noise on downstream analysis, we removed all unclassified genera and then conducted quality control on genera based on the following criteria: (1) mean relative abundance > 0.0001% and (2) occurrence in more than 3 samples. The final genus-level dataset comprises 430 unique bacterial genera.

Alpha diversity metrics, including the Richness, Chao1, ACE, and Shannon index, were calculated based on the ASV table using the “diversity” function in the Vegan R package ( The principal coordinates analysis (PCoA) utilizing the Bray-Curtis dissimilarity matrix was calculated using the “pcoa’” function in the ape R package (

Whole-genome sequencing and data processing

The 1150 host qualified DNA samples were subjected to whole genome resequencing at the Novogene facility (Novogene Co., Ltd., China). Each DNA sample was randomly fragmented into 350-bp fragments using a Covaris crusher, followed by library preparation and repairing the ends of the DNA fragments, adding of polyA tails and sequencing adaptors, and PCR amplification according to the manufacturer’s instructions for the TruSeq Nano DNA HT Sample preparation kit (Illumina USA). The PCR amplification products were then purified using the AMPure XP system, initially quantified using Qubit3.0, and the library was diluted to 1 ng/μl. The insert size and effective concentration of the library were then measured using the Agilent 2100 Bioanalyzer and real-time fluorescence PCR, respectively. The selected libraries were sequenced using the Illumina HiSeq X Ten platform (PE150). After resequencing, low-quality reads were removed to obtain high-quality clean data using Trimmomatic (v0.36). We applied the following filtering criteria to eliminate adapters and low-quality bases: reads containing more than 10% unknown nucleotides (N), reads containing more than 50% low-quality bases (Q-value < 5), and reads containing more than 10 nucleotides aligned to the adaptor sequence with up to two mismatches. The clean reads were mapped against the sheep reference genome (Oar_v1.0) using Burrows-Wheeler-Alignment Tool (BWA) [37] with the command bwa mem -M. Duplicate reads were then marked and removed using SAMBAMBA ( and indexed in SAMtools ( Variant detection was performed using the Genome Analysis Toolkit (GATK, Specifically, first gVCF files for each sample using HAplotypeCaller were generated, and then genotypes were called using GenotypeGVCFs. Finally, the SNPs identified were subjected to rigorous quality control using the VariantFiltration module using the following criteria: (1) QD > 10.0, (2) MQ > 40.0, (3) FS < 60.0, (4) MQRankSum > −12.5, (5) ReadPosRankSum > −8.0. Subsequently, the SNP datasets generated above (71,403,155 unfiltered SNP loci) were controlled for quality using vcftools with the following parameters: --remove-indels, --minDP 3, --min-alleles 2, --max-alleles 2, --max-missing 0.3, and --maf 0.05. Following these steps, a total of 23,409,311 SNPs distributed over 27 chromosomes and 1150 sheep were obtained for subsequent analysis (n_autosomal SNPs = 23,112,008,Table S10).

Core and keystone microbiota

The rumen microbial genera with a prevalence [prevalence = (the number of sheep samples in which a specific genus was detected/total number of sheep samples) × 100%] equal to 100% were defined as the core microbiota of sheep, which are the specific rumen genus present in all individuals (n = 1150) when they are considered core microbiota. To unravel the patterns of interaction between rumen microbiotas at the genus level, a microbial co-occurrence network was constructed to identify keystone taxa by calculating for each node the within-module connectivity (Zi) and among-module connectivity (Pi) [38, 39]. In details, we excluded rumen genera with a summed relative abundance of less than 0.01% and a prevalence of less than 1.5%, and then we calculated the Spearman correlation factors between microbial genera based on centered log ratio (CLR)-transformed data using the R package Hmisc ( Correlation coefficients (in absolute value) greater than 0.6 and correlation coefficient matrices with P_adj (Benjamini-Hochberg) less than 0.05 were retained for network building, and networks were constructed in the igraph package ( followed by visualization as co-occurrence networks in Gephi software. Subsequently, we calculated the network modularity and module division using the igraph R-package, and we identified the network nodes according to their within-module connectivity (Zi) and among-module connectivity (Pi). We defined each node as one of four types: (1) peripheral nodes (Zi < 2.5, Pi < 0.62), (2) connectors (Zi < 2.5, Pi > 0.62), (3) module hubs (Zi > 2.5, Pi < 0.62), and (4) network hubs (Zi > 2.5, Pi > 0.62). Connectors, module hubs, and network hubs are generally considered to be the keystone genera.

Investigation of the association between host genetics and the rumen microbiota

Mantel test

Based on the filtered SNPs dataset, we used the method proposed by Yang et al. [40] to construct a genetic relationship matrix (GRM) of 1150 individuals using GCTA v1.94.1 software [41]. Later, a microbial relationship matrix (MRM) of these 1150 individuals was constructed based on a normalized (z-value) dataset of rumen microbial ASVs using the formula described by Wen et al. [42] and a customized R script by Tang et al. [43]. Based on this, we used the Mantel test to search for the correlation between GRM and MRM using Pearson correlation (9999 permutations).

Estimation of rumen microbial heritability

Subsequently, we conducted an estimation of rumen microbial heritability (h2) at the genus level. We first excluded rumen genera that were present in less than 1.5% of the sheep samples, resulting in the retention of 290 genera. Subsequently, we classified microbial traits using a threshold of 60% prevalence, as described by Wen et al. [26]. Genera present in less than 60% of the sheep samples were considered binary genera traits (present or absent), while genera present in more than 60% of the samples were considered quantitative genera traits (relative abundance). In addition, we produced a CLR transformation for quantitative traits data to avoid spurious results. We estimated h2 of 300 microbial features (5 alpha indexes, top 5 PCoA scores, 209 binary traits and 81 quantitative traits) based on GRM using restricted maximum likelihood (REML) analysis implemented in GCTA and using birthplace, rearing season, and the top five principal components (PCs) as covariates. The estimation model is as follows:


In the above equation, \(y\) is the vector of observations for rumen microbial traits;\(b\) is the vector of fixed effects;\(a\) is the vector of additive genetic effects following a distribution of \(\mathbf{N}(0,\mathbf{G}{\upsigma }_{{\varvec{a}}}^{2})\), where \(G\) is GRM and \({\sigma }_{a}^{2}\) is the additive genetic variance; and \(e\) is the vector of residual effects following a distribution of \(\mathbf{N}\left(0,\mathbf{I}{\upsigma }_{e}^{2}\right)\), where \(\mathbf{I}\) is an identity matrix and \({\sigma }_{e}^{2}\) is the residual variance. \(X\) and \(W\) are incidence matrices for \(b\) and \(a\), respectively. The h2 was estimated as \(\frac{{{\varvec{\sigma}}}_{{\varvec{a}}}^{2}}{{{\varvec{\sigma}}}_{{\varvec{p}}}^{2}}\), where \({\upsigma }_{p}^{2}\) is the phenotypic variance. A likelihood ratio test (LRT) was used to test whether the heritability of a given phenotype was significant (PLRT < 0.05).

Microbiota genome-wide association studies (mbGWAS)

For significantly heritable microbial features above (h2, PLRT < 0.05), a mbGWAS analysis was performed using the mixed linear model (MLM) implemented in the rMVP R package [44] ( To minimize potential sources of bias in the analyses, we excluded SNPs located on the sex chromosomes for mbGWAS. Furthermore, to avoid the risk of false positives, we excluded bacterial genera with a prevalence below 30% before conducting mbGWAS for binary traits. Finally, we carried out mbGWAS to detect SNP-microbial features associations for 80 heritable microbial features (h2, PLRT < 0.05,nine heritable microbial community traits, 64 heritable quantitative genera traits, and seven heritable binary genera with a prevalence ranging between 30 and 60%) using 23,112,008 autosomal SNPs and 1150 animals. For each mbGWAS, the model included birthplace, rearing season, and top three PCs as covariates. We first used a genome-wide suggestive significant P-value threshold of 1 × 10−6 to select marker genetic variants that showed association with microbial features. This threshold was selected to maximize the strength of genetic instruments. Subsequently, we set a genome-wide significant threshold of 1 × 10−8, along with a Bonferroni-adjusted study-wide significance level of 0.05/NSNPs (0.05/23,112,008 = 2.163378e-09) for significant associations. We used Ensemble Variant Effect Predictor (VEP; for variant annotation and functional annotation.

Evaluating effects of host genetics and the rumen microbiota on body weight

Heritability and microbiability

To gain a better understanding of the relative contribution of host genetics and rumen flora to the sheep BW variation, the h2 and microbiability (m2) were examined for all weight phenotypes using the GRM and MRM described above, respectively. Phenotypic heritability was first estimated using a same model to that used for the rumen microbial heritability, as described above. However, in this case, the trait of interest (\(y\)) in the model [31] represents the vector of observations for weight phenotypes. The m2 was the concept that corresponds to the heritability, and it was used to assess the proportion of the total phenotypic variance explained by the gut microbiota [45, 46]. The estimation model is as follows:


In the above equation, \(y\) is the vector of observations for weight traits; \(m\) represents the random effect of rumen microbiota following a multinomial distribution of \(\mathbf{N}(0,\mathbf{M}{\sigma }_{m}^{2})\), where \(M\) is the MRM and \({\sigma }_{m}^{2}\) is the rumen microbial variance; \(Z\) is incidence matric for \(m\). The m2 was estimated as \(\frac{{{\varvec{\sigma}}}_{{\varvec{m}}}^{2}}{{{\varvec{\sigma}}}_{{\varvec{p}}}^{2}}\). Other parameters were consistent with the parameters in model 1 described above. In this study, it referred specifically to the part of the weight phenotypic variation caused by rumen microbiota, as operationalized by using MRM instead of GRM in GCTA.

Here, we considered both sources, host genetics and rumen microbiome, as simultaneous components of the total phenotypic variation. Briefly, we expanded the model by incorporating both the GRM and MRM as input matrices. The extended estimation model is represented as follows:

$${\varvec{y}}=\mathbf{X}\mathbf{b}+\mathbf{W}\mathbf{a}+ \mathbf{Z}\mathbf{m}+\mathbf{e}$$

The parameters were consistent with the parameters described above.

GWAS and microbiota-wide association studies (MWAS)

Subsequently, we pinpointed the marker SNPs and marker microbiotas associated with body weight at 180 days of age in sheep. Phenotypic GWAS were conducted using the method described in mbGWAS, with a significant genome-wide P-value threshold of 1 × 10−8 and a suggestive significant genome-wide P-value of 1 × 10−6. As for marker rumen genera, the filtered genus-level dataset described above was used, and statistical analyses on the associations between these 290 genera and body weight at 180 days of age were performed in R using the two-part microbiota-wide association model described by Wen et al. [26] and Fu et al. [47]. The first part of the model in our study targeted binary traits, which were named as binary traits based on the presence or absence of the bacterial genus. Specifically, a relative abundance greater than zero was coded as 1 (present), while equal to zero was coded as 0 (absent), obtaining a sample prevalence of less than 60%. The second part of the model is focused on the quantitative traits, and it is commonly used for regression analyses between phenotype and abundance of bacterial genera with sample prevalence greater than or equal to 60%. The MWAS model is described as follows:

$$\left\{\begin{array}{c}y={\beta }_{1}b+e\\ y={\beta }_{2}q+e\end{array}\right.$$

where \(y\) is the 180-day-old weight value after adjustment for birthplace and rearing season. The first three PCs were also used as a covariate in the association analyses to adjust for population substructure. \({\beta }_{1}\) and \({\beta }_{2}\) represent the regression coefficients for the binary and quantitative models, respectively. \(b\) is the binary trait, \(q\) is the CLR-transformed bacterial genus abundance, and e is the residual effect. The final calculated P-values were subjected to multiple comparisons using a Bonferroni correction, setting the threshold P adj-value to 0.05.


Landscape of core and keystone rumen microbiota composition

To identify symbiotic interactions, if any, between the host and rumen microbiota, the microbial genera that could be critical for the organization and maintenance of the rumen ecosystem were investigated. After low-quality genera were filtered out, a total of 430 unique genera were identified, among which, 17 genera comprised the core rumen microbiota (genera present in 100% of individuals; animal cohort information, see Table S1) and the cumulative abundance accounted for more than 75% of the total microbial abundance in the rumen (including Christensenellaceae R-7 group, Clostridia UCG-014, Eubacterium coprostanoligenes group, Eubacterium ruminantium group, F082, Fibrobacter, Lachnospiraceae ND3007 group, Lachnospiraceae NK3A20 group, Muribaculaceae, NK4A214 group, Prevotella, Prevotellaceae UCG-001, Rikenellaceae RC9 gut group, Ruminococcus, Saccharofermentans, Treponema, and Erysipelatoclostridiaceae UCG-004; Fig. 2 and Table S2). A co-occurrence network of rumen microbial communities was analyzed using a Spearman rank correlation coefficient matrix based on significant (PBenjamini-Hochberg < 0.05) and strong (|r| > 0.60) values to identify co-occurrence patterns. These data were crucial to identify 15 connectors as keystone genera with high connectivity between microbial communities (Fig. 2 and Table S3). Altogether, the genera that were identified showed a pivotal role in the rumen microecosystem as potential drivers of microbiota structure and function, particularly Christensenellaceae R-7 group, Lachnospiraceae NK3A20 group, and NK4A214 group, who were identified in both the core microbiota and keystone taxa.

Fig. 2
figure 2

Compositional profile of rumen bacterial communities in sheep at the genus level

Heritability of rumen microbial features

To investigate whether host genetics may impact the rumen microbiota, the Mantel test was used to examine the possible associations between the genetic relationship matrix (GRM) and microbial relationship matrix (MRM). These analyses identified significant associations between the host genome and rumen microbiota (Mantel statistic r: 0.1122, P-value: 1e-04), which indicated that host genetics can modulate specific taxa of the rumen microbiota. Subsequently, h2 was estimated for 300 rumen microbiota features in 1150 sheep to explore the relative proportion of total variation in microbiota community regulated by host genetics. These included ten measures of microbial community trait [amplicon sequence variant (ASV) richness, Chao1, ACE, Shannon index, Firmicutes:Bacteroidetes (F:B) ratio, and the first five principal coordinates (PCoAs) of the Bray-Curtis dissimilarity metric; Table S4]. A total of 81 single-taxon traits (taxa as quantitative traits), representing the relative abundance of individual microbiota genera, were detected at high sample prevalence in over 60% of animals. The h2 value was also estimated for 209 binary traits that reflected the presence or absence of a bacterial genus in a sample (limited to genera detected in 1.5–60% of animals, treated as low-abundance rumen microbiota). All microbial genera included in the estimation of h2 accounted for over 91% of the total rumen microbial abundance, with single-taxon traits accounting for 89.88% (Fig. 3a).

Fig. 3
figure 3

a Cumulative relative abundance and number of bacterial genera in the rumen of sheep at different prevalence levels. b Percentage of heritable rumen microbial genera in sheep. The red dashed line indicates mean heritability. c Heritability estimates for microbial features

We found that 64% (52/81) of the quantitative genera and nine community traits were significantly heritable [likelihood ratio test (LRT); PLRT <0.05)]. The 76% (13/17) of the core genera and 80% (12/15) of keystone taxa were also found among heritable taxa, and a similar proportion was obtained using a centered log-ratio transformation (CLR) for the correction of potential constitutive artifacts (69%, the CLR transformed quantitative genera: 56/81). Ultimately, a total of 79% of the high-prevalence genera were heritable [64/81; 82% (14/17) core microbes and 87% (13/15) keystone taxa] (Fig. 3b and Table S5). Binary rumen genera only represented 1.25% of the total abundance but could also be connected to host genome, with 22% (47/209) identified as heritable but were largely limited by their prevalence. These heritable genera showed significant heritability estimates ranging from 0.09 to 0.68 (Fig. 3c); the h2 was the highest for Candidatus Liberibacter (0.68, PLRT = 2.42E-07) which had a very low relative abundance and prevalence, whereas heritability was lowest for Caldalkalibacillus (0.09, PLRT = 5.00E-02). Among the heritable core genera, Christensenellaceae R-7 group had the highest heritability (0.58, PLRT = 1.07E-08), while F082 had the lowest value (0.10, PLRT = 0.10). Four of these heritable core genera (F082, Muribaculaceae, Prevotella, and Rikenellaceae RC9 gut group) belonged to phylum Bacteroidetes, one (Fibrobacter) belonged to Fibrobacterota, and the rest were Firmicutes. Possibly the most intriguing observation was that rumen microbial community structure, diversity, richness, and composition overview were all moderate to high heritability traits, with heritability as high as 0.51 for both ASV richness and ASV Chao1 and the first four PCoAs (cumulatively explaining 54.72% of the overall variance in rumen microbiota composition) all having a heritability above 0.2. Members of Firmicutes in the rumen represent the most predominant heritable taxa, followed by Bacteroidota (Table S5). The cumulative total abundance of heritable genera was as high as 73%, with an average prevalence of 57%. Additionally, the F:B ratio, which can reflect obesity in humans and other mammals, was also heritable (0.15, PLRT = 0.02). Overall, these values indicate that rumen microbes are largely heritable, and that the possible effects of host genetics on the rumen core microbiota and keystone taxa appear to be nearly universal.

Large-scale GWAS identified host genetics that profoundly affect sheep rumen microbiota

To gain a deeper understanding of the impact of host additive genetics on rumen microbiota, we investigated the association between 23,112,008 autosomal genetic variants and 80 rumen microbiota features that were significantly heritable (h2, PLRT < 0.05) in a cohort of 1150 sheep (see “Material and methods”: “mbGWAS”). Finally, a total of 411 associations involving 405 loci that associated independently with 1 or more of the 80 rumen metabolites at significance (P < 1 × 10−8). Using a more conservative Bonferroni-corrected study-wide significant P-value of 2.163378e-09 (0.05/23,112,008), we were able to identify 171 associations, which involved a total of 171 genomic loci and 28 microbiota features (as shown in Fig. 4 and Table S7).

Fig. 4
figure 4

Genome-wide association of sheep genetics and heritable rumen microbial variations. Manhattan plot of host genomic associations with microbial features with at least one genome-wide significant association (P < 1 × 10−6). The y-axis shows the −log10 transformation of the association P-value observed at each tested variant. The x-axis shows the genomic position of variants. The thresholds of study-wide (0.05/NSNPs; P = 2.163378e-09) and genome-wide (P = 1 × 10−8) significance are shown with horizontal lines. The autosome variants were annotated by Ensembl Variant Effect Predictor

According to the genome-wide suggestive significant threshold of P < 1 × 10−6, we identified 6845 SNP loci, of 106 were associated with two or more microbiota features (Fig. 4 and Table S7). Of these, 2709 SNPs were located within the genes. These suggestive significant SNPs were widely and evenly distributed across all autosomal chromosomes (Table S7). The average number of genetic variants for each microbiota trait was 86. The most significant associations were Prevotellaceae Ga6A1 group with the SNP: rs405050318 (P = 3.13E-13). Both genetic variants, rs427324596 and rs417318619, exhibit associations with six distinct microbiota traits. For the four alpha-diversity indexes and F:B ratio index, a total of 217 SNP loci and 259 associations were observed. Among these, two SNPs (rs427324596 and rs417318619) exhibited significant associations with all four indexes simultaneously. Furthermore, 6423 SNPs were significantly associated with 71 rumen genera (6504 associations), with 51.66% (3318 no redundant loci) of these SNPs associated with Eubacterium nodatum group, Anaerovorax, Prevotellaceae UCG-003, Moryella, UCG-005, Ruminococcus gauvreauii group, Prevotellaceae Ga6A1 group, Pseudobutyrivibrio, RF39, and UCG-002. The selection of candidate genes for heritable rumen microbes further increases the feasibility of inducing a particular microbiota in the production changing the breeding techniques, including the regulation of host phenotypic and production performances based on the control that the heritable rumen taxa may have on animal performances.

Proportion of variation in sheep body weight explained by host genetics and rumen microbiota

To explore to which extent the host additive genetics and rumen microbiota contributed to weight traits in sheep, we calculated h2 and microbiability (m2) estimates for the six BW traits using the same GRM and MRM described above, respectively (Fig. 5a and Table S8). We found all BW traits to have upper moderate to high heritability, ranging from 36 for 100-day-old weight to 61% for birth weight, with a heritability of 39% for 180-day-old weight (Fig. 5a and Table S8). As with heritability, the m2 was between 0 and 1 and with increasing m2 values, a major contribution of rumen microbiota to the traits was demonstrated. After correction for host genetics, 180-day-old weight had a moderate microbiability estimate of 20%. The m2 for birth weight obtained using rumen microbiota data collected at 180-day-old was nonsignificant, and the mean m2 for BW at other ages was 20%, ranging from 17% for 100-day-old weight to 24% for 100-day-old weight.

Fig. 5
figure 5

a Heritability (h2) and microbiability (m2) for weight traits in sheep. b Manhattan plot shows the results of microbiota-wide association studies (MWAS). Four bacterial taxa were significantly associated with body weight at 180 days of age in sheep (Bonferroni-corrected P-value < 0.05). Red labels indicate that the bacterial genus is a core microbe, and green labels indicate that the bacterial genus is a keystone taxon. c and d A linear model was fitted to identify the association between rumen marker genera and BW. e Correlation patterns showing the rumen marker genera associated with rumen volatile fatty acids and rumen morphology

We further considered both host genetics and rumen microbiome as simultaneous components of the total phenotypic variation and calculated the h2 and m2 of body weight traits. The results from these models were largely consistent (Fig. 5a and Table S8). The combined contribution of host genetics and rumen microbiota to the 180-day-old weight phenotype was 52%. Specifically, the heritability was estimated to be 0.35, indicating the proportion of phenotypic variance attributed to host genetics. Additionally, the microbiota contribution was estimated to be 0.17, indicating the proportion of phenotypic variance attributed to the rumen microbiota.

Heritable rumen microbes affect sheep weight

Above, we demonstrated that host genetics can influence the rumen microbial community in sheep, and that together with rumen microbiota, may contribute to sheep BW. Later, a two-part model for the MWAS analysis between rumen microbiota and BW in sheep [26, 47] was used to identify which indicator of microbiota taxa may regulate this mechanism. From this analysis, four last genera were identified (Fig. 5b and Table S9), including Lachnospiraceae ND3007 group, Rikenellaceae RC9 gut group, Syntrophococcus, and Oribacterium, which were all heritable and significantly associated with 180-day-old weight (Bonferroni corrected P < 0.05). Among them, Lachnospiraceae ND3007 group and Rikenellaceae RC9 gut group emerged as moderately heritable core taxa, explaining 10.85% and 10.32% of variance on 180-day-old weight (Fig. 5c & d), respectively. By contrast, the linear fit correlation coefficients for the potential core and keystone genera Oribacterium (prevalence = 99.91%) and Syntrophococcus (prevalence = 99.57%) were both over 9.8%.

To investigate the potential relationships among BW-related microbial markers and rumen metabolic profiles, Spearman’s rank correlation analysis was performed between the selected bacterial genera and rumen volatile fatty acids (VFA) profile, rumen papillae length and width, and muscle thickness (Fig. 5e). This analysis showed that total VFAs (sum of all individual VFAs) shared a positive correlation with Syntrophococcus but were negatively correlated with Lachnospiraceae ND3007 group, Rikenellaceae RC9 gut group, and Oribacterium. Acetic acid was negatively correlated with Lachnospiraceae ND3007 group. Propionic acid levels were positively correlated with Syntrophococcus and Oribacterium but negatively correlated with Rikenellaceae RC9 gut group. Likewise, the abundance of candidate genera was significantly associated with rumen histomorphology. In particular, Lachnospiraceae ND3007 group and Oribacterium were positively associated with rumen papillae length, Rikenellaceae RC9 gut group was correlated with rumen papillae width, Syntrophococcus was negatively correlated with rumen papillae width, and Rikenellaceae RC9 gut group shared a negative relationship with rumen muscle thickness (Fig. 5e). Overall, these four weight-related indicator taxa and their association with metabolic profiles further emphasized the central role of heritable rumen microbiota, specifically core and keystone taxa, in the host phenotype and rumen microbiota metabolism in sheep.

Genome-wide association study reveals novel loci controlling sheep weight

We next conducted GWAS to identify the SNPs loci in the host genome linked to BW in sheep, using the phenotypic data from the same cohort of 1150 animals and the same autosomal SNPs dataset used for the mbGWAS. A total of five SNP loci were identified to be significantly associated with the body weight of sheep at 180 days old, at a significance level of P < 1 × 10−8, all located on OAR 9. With a more moderate genome-wide suggestive significant P-value of 1 × 10−6, we identified 94 SNPs, spanning OARs 1, 2, 4, 6, 9, 10, 13, 15, 16, and 23 (Fig. 6a and Table S10). The most significant association was observed with the intergenic SNP, rs413796993 (P = 6.99E-14), on OAR 9 at 38007144bp. The second most significant SNP (rs417240663; P = 3E-10) was located on OAR 9 at 38659928bp, in an intron of the XKR4 gene, which cumulatively harbored 3 significant and 16 suggestive significant SNPs. Notably, we observed that several regions in OAR 6 and OAR 9 together contained more than 84% (79/94) of the marker SNPs. Eight of these loci were located within an ~900-kb region on OAR 6 (33.97–34.88 Mb), including four located in introns of BMPRIB, one in introns of PDLIM5, and one in introns of UNC5C. Additionally, two SNPs were located in a 72.66-kb region between 42.21 and 42.28 Mb on OAR 6, both in the intron of LCORL. On OAR 9, a large region spanning ~4.8 Mb (37.42–42.23 Mb) was found to contain 60 SNPs, 33 of which were located within seven genes, including ATP6V1H, LYPLA1, RP1, XKR4, SDCBP, TOX, and CA8; a SNP located in KCNQ3 was also identified on OAR 9. In addition, four of the 15 SNPs on other autosomes were also located in introns, such as rs399594440 in UNC80 on OAR2, 10-74576409 (variant ID unknown) in GPC5 on OAR10, rs425587495 in NRP1 on OAR13, and rs417185097 in EXT2 on OAR15. Altogether, the analyses described, which were based on a large sample size and with a reduced environmental interference, led to the identification of genetic regions and novel candidate genes that significantly may influence BW in sheep, which provide new molecular markers for breeding efforts to regulate sheep growth and development and thus improve meat production.

Fig. 6
figure 6

a and c Genome-wide association analysis (GWAS) for body weight at 180 days of age and Lachnospiraceae ND3007 group (the marker genera most strongly associated with 180-day-old weight) in sheep. The horizontal red solid line indicates the genome-wide significance (P < 1 × 10−8) thresholds. b and d Magnified view of the region of interest. e Comparison of the body weight at 180 days old of different genotypes of the rs405307925 SNP of TOX. f Comparison of the Lachnospiraceae ND3007 group in rumen of different genotypes of the rs405307925 SNP of TOX

Recursive influence of host genetics on rumen microbiota drive sheep weight

Although host genetics and rumen microbiota jointly contribute to the sheep weight, our results showed that the host genetics and the rumen microbiota are not independent, and that their effects on weight traits are also not independent and unlinked, while there are some weight-associated microbial features that are controlled by the weight-associated genetic variants. In this case, the sheep genome may have both direct and indirect (marker microbiota-mediated) effects on weight trait, and this phenomenon was defined as a “recursive” scenario. To this regard, we investigated the host genetic aspects that may be responsible for this scenario and the SNPs loci that may have both direct and indirect effects on sheep weight. To address this question, the resultant data from GWAS, mbGWAS, and MWAS were combined and analyzed. Based on finding above showing that the four microbial features were significantly associated with 180-day weight (MWAS identified), we found that all had moderate to strong h2 (0.16–0.40, PLRT < 0.05), and 306 SNPs affecting these four marker genera were identified (mbGWAS identified; Table S11), indicating an “indirect’” scenario for the control of phenotypes by genetic variants, i.e., host genetics indirectly influencing sheep weight through the control of weight-associated microbiota. Subsequently, we compared the overlapping loci of marker microbiota-associated SNPs in “indirect” scenario (mbGWAS identified) and weight-associated SNPs (GWAS identified) to identify host genetics that influence weight through a “recursive” model. Intriguingly, we found that the rs405307925, located at 40,691,375 bp on OAR 9 (Fig. 6, Tables S7–S11), appeared in both the mbGWAS signals sets of the Lachnospiraceae ND3007 group (the marker genera most strongly associated with 180-day-old weight) and the 180-day-old weight GWAS signals sets. Therefore, this represents the first time that a “recursive” scenario of genomic markers influencing the phenotypes was identified using a locus overlap algorithm in sheep. To identify more potential “recursive” scenarios, we then used a Kruskal-Wallis test to examine the differences in BW within marker microbiota-associated SNPs (mGWAS identified) genotypes in the “indirect” scenario described above. Finally, we identified six genetic markers (PKruskal–Wallis test < 0.01; Table S12) that have the potential to recursive impact the body weight in sheep. These recursive loci consist of rs405307925, which is situated on the TOX gene, as well as rs425633529 and rs604806950, which are located on the GLRX3 and PRR14 genes, respectively.


Body weight represents the most economically important trait for sheep production, and it may be strongly influenced by host genetics and rumen microbial ecosystem. However, these studies were only focused on single effects, that is to say only the independent influence of host genome or the microbiota. Therefore, whether the host genetics may interact with the rumen microbiota to influence BW in sheep has remained largely unknown to date. Here, we described a relationship between host genome, rumen microbiota, and BW phenotype in 1150 sheep subjected to the same diet and management conditions. To gain better insights into this topic, we described the core and keystone structure of the rumen microbiota. We calculated the h2 of the rumen microbiota, and mbGWAS was performed to identify the influence of host genetics on individual rumen microbial genera. The h2 and m2 were also estimated for BW, identified host genetics loci, and rumen microbiota affecting BW. The different causal scenarios in which genetic markers may influence the BW phenotype were analyzed using the locus overlap algorithm and Kruskal-Wallis test. To our best knowledge, this is the first study that reports a large association between host genome and rumen microbiota and BW in sheep to date. Our study provided detailed guidelines on the combined use of microbiota and host genome information in order to predict the complex traits and the reliability of parameter inferences.

Mantel tests are a widely used method in microbiome analysis to examine the macroscopic correlations between host genetics and gut microecological variation. In this study, a significant correlation was observed between the host genome kinship matrix and the rumen microbiota relationship matrix, which is consistent with previous findings in human feces [48] and mice [49]. However, it is important to note that the Mantel test used to investigate the impact of host genetics on shaping the overall structure of the rumen microbiome has its limitations. The Mantel test is a matrix-based method that can only detect linear relationships within the matrix while disregarding nonlinear interactions [50]. Furthermore, the rumen microbiota community is highly complex and diverse, with intricate interactions between different species that can influence the magnitude of the correlation coefficient. Therefore, despite the weak but significant correlation observed, it still suggests that the host genetic background plays a role in shaping the structure of the rumen microbiota.

In this study for the first time, the heritability of the rumen microbiota has been investigated on sheep, and it was demonstrated that the rumen microbial community diversity was heritable, consistently with previously studies in cattle [29], demonstrating a higher h2 than those reported for cattle. This finding may be attributed to three main factors: firstly, differences between host species,secondly, the more control for environmental factors and covariates in the design and analysis of the current study [21, 51] and thirdly, higher detection of variants than SNP chips following the resequencing. All together, these observations also indicated that most of the variation in rumen microbial communities is due to host genetics in sheep. On the other hand, estimating bacterial heritability provided new perspectives in identifying taxa that were closely related to host but previously unidentified, such as Candidatus Liberibacter and Longimicrobiaceae, both of which showed a heritability above 50% but for which no host interactions were detected. However, we found that the heritability of the very low-abundance rumen taxa in sheep was limited by their prevalence in the population (most taxa were present in < 10% of samples). This may be due to which the effective sample size and efficacy of genetic analyses, and therefore very large sample sizes and more sensitive microbial community sequencing technologies (e.g., metagenomic and metatranscriptomic approaches), will be beneficial to address this issue in the future.

We next focused on the BW value, which was measured at 180 days of age, the typical age for slaughter in commercial meat lambs. Host genetics explained 39% of the total phenotypic variation, and indeed, the effect of the rumen microbial communities on sheep BW was 20% following the consideration of host genetics. These results highlighted the need for additional comprehensive analyses based on the genome and rumen microbiota. We therefore used MWAS to screen rumen microbes associated with BW, and we identified four heritable bacterial genera. Interestingly, these taxa derived from a core or keystone microbial subsets, suggesting that ecologically important core and keystone taxa play a role in the development of BW phenotypes in sheep. The increase in major VFAs is beneficial for the animal, and Syntrophococcus produced major VFAs that are absorbed through the rumen wall to sustain the energy requirements [52, 53]. Oribacterium has been found to be an obesogenic bacterium in human and rat studies [54, 55], and the increasing of Oribacterium in the rumen of cattle may promote the synthesis of a variety of fatty acids [56]. For example, the changes in propionate concentrations are correlated with the changes in the abundance of most marker genera. Propionate is the main precursor of hepatic gluconeogenesis in ruminants and is effectively involved in their energy balance [57].

By applying large-scale GWAS, two signal peaks on chromosomes 6 and 9 in sheep were observed. The OAR6 signalling peak region included an orthologous region that has been clearly associated with body size in a variety of mammals [58]. We replicated the relationship between the LCORL gene and BW in sheep found by Al-Mamun and colleagues [59]. Biologically relevant genes localized to the OAR9 signalling peak region were KCNQ3, ATP6V1H, LYPLA1, and XKR4. KCNQ3 modulates M-currents in NPY/AgRP neurons, affecting neuronal excitability and stimulating the physiological appetite, thereby it contributes to the energy balance [60]. The ATP6V1H gene is involved in mediating the regulatory bone formation [61], and SNPs in ATP6V1H were previously reported to affect feed efficiency in beef cattle [62]. LYPLA1 has been identified as a potential candidate gene associated with lipid-related biological processes [63]. XKR4 encodes an XK-related protein in the XK-Kyle blood group complex. XKR4 variants which were associated with economically important traits such as growth, and feed efficiency, have been widely detected in cattle, but not in sheep [64]. In addition, some interesting regulatory candidates were identified, with BMPR1B (also known as FECB and ALK-6), TOX, and SDCBP the ones of major interest. Several studies have demonstrated the key role of TOX and SDCBP on carcass weight, fertility, and feed efficiency traits in cattle [65]. An important gene involved in glucose metabolism is CA8 [66], while NRP1 is involved in the regulation of organ development and function [67]. Furthermore, PDLIM5 encodes a cytoskeleton-associated protein that plays an important role in cell proliferation and differentiation in a variety of tissues and cell types [68], while PTTG1IP was previously reported to be a selectively spliced gene associated with chicken muscle development [69]. The EXT2 gene encodes an essential component of the glycosyltransferase complex required for the biosynthesis of acetyl heparin, which in turn regulates the signalling involved in bone formation [70]. UNC80 is one of the subunits of the NALCN complex, which is associated with global dysplasia [71]. Although our analysis has identified several SNPs and genes previously not reported using the GWAS analysis, functional validation of these variants is required.

Following an integrated analysis of GWAS, mbGWAS, and MWAS results, our study described a recursive scenario of host genetics influencing phenotype, based on the hypothesis that the same SNPs are directly associated with animal phenotype able to influence the composition of the rumen flora and consequently also the host phenotype. We identified a recursive scheme in which rs405307925, located within the TOX gene, affects sheep BW, and its indirect effects are mediated through the Lactococcus ND3007 group. Notably, this is also the first report in which the SNPs discovered in the target phenotypic GWAS were overlapped with those in the GWAS of microbes associated with the target phenotype. Here, although the biological function of the TOX gene in sheep has not yet been reported, it is noteworthy that three previous studies described the key role of the TOX gene in regulating the immune T -cell function [72,73,74]. Although rumen microbes play an important role in host immunity and metabolism, multiple immunological mechanisms are put in place by the host to avoid microbial ecological dysregulation. Therefore, we hypothesized that the host genetics may regulate the immune system, controlling the complex rumen microbes. Future research will be needed to elucidate these complex mechanisms.

This study has some limitations. While it was conducted using the most dominant breed of sheep in China, Hu sheep, generalizing the findings to other breeds, may require additional considerations. As all study animals were males, the applicability to other genders of sheep may be limited. Additionally, the exact mechanism underlying the “recursive” scheme remains elusive. In future studies, a combination of metabolomics, metagenomics, and gene editing techniques will be employed to investigate the mechanisms involved in this recursive scenario.


In this study, we explored the association of host genetic variation and rumen microbiota and its impact on sheep BW. Using large-scale GWAS and MWAS, novel candidate genes and heritable rumen indicator taxa that potentially affect sheep BW were identified. Furthermore, our results supported the hypothesis that host genetics influences the composition of the rumen microbiota, with an intricate network involving host genetics, rumen microbiota, and sheep weight traits. Additionally, we identified six SNP loci that influence sheep BW through a “recursive” model. Of particular interest is rs405307925, located within the TOX gene, which can affect BW both directly and indirectly by influencing the relative abundance of the Lachnospiraceae ND3007 group, which in turn affect the BW through the production of VFAs. Our work is the first study that explored the host genetic-rumen microbiota-phenotype relationships on a large-scale sheep population, which was sufficiently attenuated from other factors such as environment and diet. These observations will provide a reliable and detailed guide that can be helpful for the combined use of rumen microbiota, and genomic information useful for predicting complex traits and parameter inferences, which may be useful for developing breeding strategies to improve sheep BW.

Availability of data and materials

Raw data were deposited in the National Center for Biotechnology Information database under BioProjectID PRJNA867677. Further data requests and informations can be obtained  by contacting the corresponding author Weimin Wang at 


  1. Li X, Yang J, Shen M, Xie XL, Liu GJ, Xu YX, Lv FH, Yang H, Yang YL, Liu CB, Zhou P, Wan PC, Zhang YS, Gao L, Yang JQ, Pi WH, Ren YL, Shen ZQ, Wang F, Deng J, Xu SS, Salehian-Dehkordi H, Hehua E, Esmailizadeh A, Dehghani-Qanatqestani M, Štěpánek O, Weimann C, Erhardt G, Amane A, Mwacharo JM, Han JL, Hanotte O, Lenstra JA, Kantanen J, Coltman DW, Kijas JW, Bruford MW, Periasamy K, Wang XH, Li MH. Whole-genome resequencing of wild and domestic sheep identifies genes associated with morphological and agronomic traits. Nat Commun. 2020;11:2815.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  2. Dos Santos ACS, Santos SA, Carvalho GGP, Mariz LDS, Tosto MSL, Valadares Filho SC, Azevedo JAG. A comparative study on the excretion of urinary metabolites in goats and sheep to evaluate spot sampling applied to protein nutrition trials. J Anim Sci. 2018;96:3381–97.

    Article  PubMed  PubMed Central  Google Scholar 

  3. McHugh N, Pabiou T, McDermott K, Wall E, Berry DP. A novel measure of ewe efficiency for breeding and benchmarking purposes. J Anim Sci. 2018;96:2051–9.

    Article  PubMed  PubMed Central  Google Scholar 

  4. Godfray, H.C.J., Aveyard, P., Garnett, T., Hall, J.W., Key, T.J., Lorimer, J., Pierrehumbert, R.T., Scarborough, P., Springmann, M., Jebb, S.A., 2018. Meat consumption, health, and the environment. Science 361

  5. Alexandratos N, Bruinsma J. World agriculture towards 2030/2050. ESA Working paper No. 12-03. Rome: Food and Agriculture Organization of the United Nations; 2012.

    Google Scholar 

  6. Du Y, Ge Y, Ren Y, Fan X, Pan K, Lin L, Wu X, Min Y, Meyerson LA, Heino M, Chang SX, Liu X, Mao F, Yang G, Peng C, Qu Z, Chang J, Didham RK. A global strategy to mitigate the environmental impact of China’s ruminant consumption boom. Nat Commun. 2018;9:4133.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  7. Koeth RA, Wang Z, Levison BS, Buffa JA, Org E, Sheehy BT, Britt EB, Fu X, Wu Y, Li L, Smith JD, DiDonato JA, Chen J, Li H, Wu GD, Lewis JD, Warrier M, Brown JM, Krauss RM, Tang WH, Bushman FD, Lusis AJ, Hazen SL. Intestinal microbiota metabolism of L-carnitine, a nutrient in red meat, promotes atherosclerosis. Nat Med. 2013;19:576–85.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  8. Ahmad SF, Khan NN, Ganai NA, Shanaz S, Rather MA, Alam S. Multivariate quantitative genetic analysis of body weight traits in Corriedale sheep. Trop Anim Health Prod. 2021;53:197.

    Article  PubMed  Google Scholar 

  9. Ghaderi-Zefrehei M, Safari A, Moridi M, Khanzadeh H, Dehsaraei AR. Bayesian estimate of genetic parameters for growth traits in Lori Bakhtiari sheep. Trop Anim Health Prod. 2021;53:457.

    Article  PubMed  Google Scholar 

  10. Mortimer SI, Hatcher S, Fogarty NM, van der Werf JHJ, Brown DJ, Swan AA, Greeff JC, Refshauge G, Edwards JEH, Gaunt GM. Genetic parameters for wool traits, live weight, and ultrasound carcass traits in Merino sheep. J Anim Sci. 2017;95:1879–91.

    Article  CAS  PubMed  Google Scholar 

  11. Tesema Z, Deribe B, Lakew M, Getachew T, Tilahun M, Belayneh N, Kefale A, Shibesh M, Zegeye A, Yizengaw L, Alebachew GW, Tiruneh S, Kiros S, Asfaw M, Bishaw M. Genetic and non-genetic parameter estimates for growth traits and Kleiber ratios in Dorper × indigenous sheep. Animal. 2022;16:100533.

    Article  CAS  PubMed  Google Scholar 

  12. Zhang, X., Li, G., Li, F., Zhang, D., Yuan, L., Zhao, Y., Zhang, Y., Li, X., Song, Q., Wang, W., 2021a. Effect of feed efficiency on growth performance, body composition, and fat deposition in growing Hu lambs. Anim Biotechnol, 1-16.

  13. Cao Y, Song X, Shan H, Jiang J, Xiong P, Wu J, Shi F, Jiang Y. Genome-wide association study of body weights in Hu sheep and population verification of related single-nucleotide polymorphisms. Front Genet. 2020;11:588.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  14. Ghasemi M, Zamani P, Vatankhah M, Abdoli R. Genome-wide association study of birth weight in sheep. Animal. 2019;13:1797–803.

    Article  CAS  PubMed  Google Scholar 

  15. Zhao B, Luo H, Huang X, Wei C, Di J, Tian Y, Fu X, Li B, Liu GE, Fang L, Zhang S, Tian K. Integration of a single-step genome-wide association study with a multi-tissue transcriptome analysis provides novel insights into the genetic basis of wool and weight traits in sheep. Genet Sel Evol. 2021;53:56.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  16. Zhang YK, Zhang XX, Li FD, Li C, Li GZ, Zhang DY, Song QZ, Li XL, Zhao Y, Wang WM. Characterization of the rumen microbiota and its relationship with residual feed intake in sheep. Animal. 2021;15:100161.

    Article  CAS  PubMed  Google Scholar 

  17. Yin X, Ji S, Duan C, Ju S, Zhang Y, Yan H, Liu Y. Rumen fluid transplantation affects growth performance of weaned lambs by altering gastrointestinal microbiota, immune function and feed digestibility. Animal. 2021;15:100076.

    Article  CAS  PubMed  Google Scholar 

  18. Bergman EN. Energy contributions of volatile fatty acids from the gastrointestinal tract in various species. Physiol Rev. 1990;70:567–90.

    Article  CAS  PubMed  Google Scholar 

  19. Morgavi DP, Forano E, Martin C, Newbold CJ. Microbial ecosystem and methanogenesis in ruminants. Animal. 2010;4:1024–36.

    Article  CAS  PubMed  Google Scholar 

  20. Benson AK, Kelly SA, Legge R, Ma F, Low SJ, Kim J, Zhang M, Oh PL, Nehrenberg D, Hua K, Kachman SD, Moriyama EN, Walter J, Peterson DA, Pomp D. Individuality in gut microbiota composition is a complex polygenic trait shaped by multiple environmental and host genetic factors. Proc Natl Acad Sci U S A. 2010;107:18933–8.

    Article  PubMed  PubMed Central  Google Scholar 

  21. Grieneisen L, Dasari M, Gould TJ, Björk JR, Grenier JC, Yotova V, Jansen D, Gottel N, Gordon JB, Learn NH, Gesquiere LR, Wango TL, Mututua RS, Warutere JK, Siodi L, Gilbert JA, Barreiro LB, Alberts SC, Tung J, Archie EA, Blekhman R. Gut microbiome heritability is nearly universal but environmentally contingent. Science. 2021;373:181–6.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  22. Yang H, Wu J, Huang X, Zhou Y, Zhang Y, Liu M, Liu Q, Ke S, He M, Fu H, Fang S, Xiong X, Jiang H, Chen Z, Wu Z, Gong H, Tong X, Huang Y, Ma J, Gao J, Charlier C, Coppieters W, Shagam L, Zhang Z, Ai H, Yang B, Georges M, Chen C, Huang L. ABO genotype alters the gut microbiota by regulating GalNAc levels in pigs. Nature. 2022;606:358–67.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  23. Lopera-Maya EA, Kurilshikov A, van der Graaf A, Hu S, Andreu-Sánchez S, Chen L, Vila AV, Gacesa R, Sinha T, Collij V, Klaassen MAY, Bolte LA, Gois MFB, Neerincx PBT, Swertz MA, Harmsen HJM, Wijmenga C, Fu J, Weersma RK, Zhernakova A, Sanna S. Effect of host genetics on the gut microbiome in 7,738 participants of the Dutch Microbiome Project. Nat Genet. 2022;54:143–51.

    Article  CAS  PubMed  Google Scholar 

  24. Qin Y, Havulinna AS, Liu Y, Jousilahti P, Ritchie SC, Tokolyi A, Sanders JG, Valsta L, Brożyńska M, Zhu Q, Tripathi A, Vázquez-Baeza Y, Loomba R, Cheng S, Jain M, Niiranen T, Lahti L, Knight R, Salomaa V, Inouye M, Méric G. Combined effects of host genetics and diet on human gut microbiota and incident disease in a single population cohort. Nat Genet. 2022;54:134–42.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  25. Sanna S, Kurilshikov A, van der Graaf A, Fu J, Zhernakova A. Challenges and future directions for studying effects of host genetics on the gut microbiome. Nat Genet. 2022;54:100–6.

    Article  CAS  PubMed  Google Scholar 

  26. Wen C, Yan W, Mai C, Duan Z, Zheng J, Sun C, Yang N. Joint contributions of the gut microbiota and host genetics to feed efficiency in chickens. Microbiome. 2021;9:126.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  27. Wen C, Yan W, Sun C, Ji C, Zhou Q, Zhang D, Zheng J, Yang N. The gut microbiota is largely independent of host genetics in regulating fat deposition in chickens. Isme j. 2019;13:1422–36.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  28. Difford GF, Plichta DR, Løvendahl P, Lassen J, Noel SJ, Højberg O, Wright AG, Zhu Z, Kristensen L, Nielsen HB, Guldbrandtsen B, Sahana G. Host genetics and the rumen microbiome jointly associate with methane emissions in dairy cows. PLoS Genet. 2018;14:e1007580.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  29. Li F, Li C, Chen Y, Liu J, Zhang C, Irving B, Fitzsimmons C, Plastow G, Guan LL. Host genetics influence the rumen microbiota and heritable rumen microbial features associate with feed efficiency in cattle. Microbiome. 2019;7:92.

    Article  PubMed  PubMed Central  Google Scholar 

  30. Wallace RJ, Sasson G, Garnsworthy PC, Tapio I, Gregson E, Bani P, Huhtanen P, Bayat AR, Strozzi F, Biscarini F, Snelling TJ, Saunders N, Potterton SL, Craigon J, Minuti A, Trevisi E, Callegari ML, Cappelli FP, Cabezas-Garcia EH, Vilkki J, Pinares-Patino C, Fliegerová KO, Mrázek J, Sechovcová H, Kopečný J, Bonin A, Boyer F, Taberlet P, Kokou F, Halperin E, Williams JL, Shingfield KJ, Mizrahi I. A heritable subset of the core rumen microbiome dictates dairy cow productivity and emissions. Sci Adv. 2019;5:eaav8391.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  31. Abbas W, Howard JT, Paz HA, Hales KE, Wells JE, Kuehn LA, Erickson GE, Spangler ML, Fernando SC. Influence of host genetics in shaping the rumen bacterial community in beef cattle. Sci Rep. 2020;10:15101.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  32. Zhang Q, Difford G, Sahana G, Løvendahl P, Lassen J, Lund MS, Guldbrandtsen B, Janss L. Bayesian modeling reveals host genetics associated with rumen microbiota jointly influence methane emission in dairy cows. Isme j. 2020;14:2019–33.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  33. Mahala S, Kala A, Kumar A. Host genetics associated with gut microbiota and methane emission in cattle. Mol Biol Rep. 2022;49:8153–61.

    Article  CAS  PubMed  Google Scholar 

  34. Martínez-Álvaro M, Auffret MD, Duthie CA, Dewhurst RJ, Cleveland MA, Watson M, Roehe R. Bovine host genome acts on rumen microbiome function linked to methane emissions. Commun Biol. 2022;5:350.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  35. Magoč T, Salzberg SL. FLASH: fast length adjustment of short reads to improve genome assemblies. Bioinformatics. 2011;27:2957–63.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  36. Callahan BJ, McMurdie PJ, Rosen MJ, Han AW, Johnson AJ, Holmes SP. DADA2: high-resolution sample inference from Illumina amplicon data. Nat Methods. 2016;13:581–3.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  37. Li H, Durbin R. Fast and accurate short read alignment with Burrows-Wheeler transform. Bioinformatics. 2009;25:1754–60.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  38. Olesen JM, Bascompte J, Dupont YL, Jordano P. The modularity of pollination networks. Proc Natl Acad Sci U S A. 2007;104:19891–6.

    Article  PubMed  PubMed Central  Google Scholar 

  39. Deng Y, Jiang YH, Yang Y, He Z, Luo F, Zhou J. Molecular ecological network analyses. BMC Bioinformatics. 2012;13:113.

    Article  PubMed  PubMed Central  Google Scholar 

  40. Yang J, Benyamin B, McEvoy BP, Gordon S, Henders AK, Nyholt DR, Madden PA, Heath AC, Martin NG, Montgomery GW, Goddard ME, Visscher PM. Common SNPs explain a large proportion of the heritability for human height. Nat Genet. 2010;42:565–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  41. Yang J, Lee SH, Goddard ME, Visscher PM. GCTA: a tool for genome-wide complex trait analysis. Am J Hum Genet. 2011;88:76–82.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  42. Wen C, Yan W, Sun C, Ji C, Zhou Q, Zhang D, Zheng J, Yang N. The gut microbiota is largely independent of host genetics in regulating fat deposition in chickens. Isme Journal. 2019;13(6):1422–36.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  43. Tang S, Xin Y, Ma Y, Xu X, Zhao S, Cao J. Screening of microbes associated with swine growth and fat deposition traits across the intestinal tract. Front Microbiol. 2020;11:586776.

    Article  PubMed  PubMed Central  Google Scholar 

  44. Yin L, Zhang H, Tang Z, Xu J, Yin D, Zhang Z, Yuan X, Zhu M, Zhao S, Li X, Liu X. rMVP: a memory-efficient, visualization-enhanced, and parallel-accelerated tool for genome-wide association study. Genomics Proteomics Bioinformatics. 2021;19:619–28.

    Article  PubMed  PubMed Central  Google Scholar 

  45. Difford, G., Lassen, J., Løvendahl, P., 2016. Genes and microbes, the next step in dairy cattle breeding. Book of Abstracts of the 67th Annual Meeting of the European Federation of Animal Science, 285-285.

  46. Difford GF, Plichta DR, Løvendahl P, Lassen J, Noel SJ, Højberg O, Wright A-DG, Zhu Z, Kristensen L, Nielsen HB. Host genetics and the rumen microbiome jointly associate with methane emissions in dairy cows. Plos Genetics. 2018;14:e1007580.

    Article  PubMed  PubMed Central  Google Scholar 

  47. Fu J, Bonder MJ, Cenit MC, Tigchelaar EF, Maatman A, Dekens JA, Brandsma E, Marczynska J, Imhann F, Weersma RK, Franke L, Poon TW, Xavier RJ, Gevers D, Hofker MH, Wijmenga C, Zhernakova A. The gut microbiome contributes to a substantial proportion of the variation in blood lipids. Circ Res. 2015;117:817–24.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  48. O’Connor A, Quizon PM, Albright JE, Lin FT, Bennett BJ. Responsiveness of cardiometabolic-related microbiota to diet is influenced by host genetics. Mamm Genome. 2014;25:583–99.

  49. Suzuki TA, Phifer-Rixey M, Mack KL, Sheehan MJ, Lin D, Bi K, Nachman MW. Host genetic determinants of the gut microbiota of wild mice. Mol Ecol. 2019;28:3197–207.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  50. Legendre P, Legendre L. Numerical ecology. Elsevier; 2012.

    Google Scholar 

  51. Speed D, Holmes J, Balding DJ. Evaluating and improving heritability models using summary statistics. Nat Genet. 2020;52:458–62.

    Article  CAS  PubMed  Google Scholar 

  52. Cui Z, Wu S, Liu S, Sun L, Feng Y, Cao Y, Chai S, Zhang G, Yao J. From maternal grazing to barn feeding during pre-weaning period: altered gastrointestinal microbiota contributes to change the development and function of the rumen and intestine of yak calves. Front Microbiol. 2020;11:485.

    Article  PubMed  PubMed Central  Google Scholar 

  53. Li M, Zhong H, Li M, Zheng N, Wang J, Zhao S. Contribution of ruminal bacteriome to the individual variation of nitrogen utilization efficiency of dairy cows. Front Microbiol. 2022;13:815225.

    Article  PubMed  PubMed Central  Google Scholar 

  54. Mayengbam S, Mickiewicz B, Trottier SK, Mu C, Wright DC, Reimer RA, Vogel HJ, Shearer J. Distinct gut microbiota and serum metabolites in response to weight loss induced by either dairy or exercise in a rodent model of obesity. J Proteome Res. 2019;18:3867–75.

    Article  CAS  PubMed  Google Scholar 

  55. Nardelli, C., Granata, I., D'Argenio, V., Tramontano, S., Compare, D., Guarracino, M.R., Nardone, G., Pilone, V., Sacchetti, L., 2020. Characterization of the duodenal mucosal microbiome in obese adult subjects by 16S rRNA sequencing. Microorganisms 8

  56. Sacks D, Baxter B, Campbell BCV, Carpenter JS, Cognard C, Dippel D, Eesa M, Fischer U, Hausegger K, Hirsch JA, Shazam Hussain M, Jansen O, Jayaraman MV, Khalessi AA, Kluck BW, Lavine S, Meyers PM, Ramee S, Rüfenacht DA, Schirmer CM, Vorwerk D. Multisociety Consensus Quality Improvement Revised Consensus Statement for Endovascular Therapy of Acute Ischemic Stroke. Int J Stroke. 2018;13:612–32.

    Article  PubMed  Google Scholar 

  57. Overton TR, Waldron MR. Nutritional management of transition dairy cows: strategies to optimize metabolic health. Journal of Dairy Science. 2004;87:E105–19.

    Article  Google Scholar 

  58. Xu L, Bickhart DM, Cole JB, Schroeder SG, Song J, Tassell CP, Sonstegard TS, Liu GE. Genomic signatures reveal new evidences for selection of important traits in domestic cattle. Mol Biol Evol. 2015;32:711–25.

    Article  CAS  PubMed  Google Scholar 

  59. Al-Mamun HA, Kwan P, Clark SA, Ferdosi MH, Tellam R, Gondro C. Genome-wide association study of body weight in Australian Merino sheep reveals an orthologous region on OAR6 to human and bovine genomic regions affecting height and weight. Genet Sel Evol. 2015;47:66.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  60. Stincic TL, Bosch MA, Hunker AC, Juarez B, Connors AM, Zweifel LS, Rønnekleiv OK, Kelly MJ. CRISPR knockdown of Kcnq3 attenuates the M-current and increases excitability of NPY/AgRP neurons to alter energy balance. Mol Metab. 2021;49:101218.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  61. Zhang Y, Huang H, Zhao G, Yokoyama T, Vega H, Huang Y, Sood R, Bishop K, Maduro V, Accardi J, Toro C, Boerkoel CF, Lyons K, Gahl WA, Duan X, Malicdan MC, Lin S. ATP6V1H deficiency impairs bone development through activation of MMP9 and MMP13. PLoS Genet. 2017;13:e1006481.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  62. de Las Heras-Saldana S, Clark SA, Duijvesteijn N, Gondro C, van der Werf JHJ, Chen Y. Combining information from genome-wide association and multi-tissue gene expression studies to elucidate factors underlying genetic variation for residual feed intake in Australian Angus cattle. BMC Genomics. 2019;20:939.

    Article  CAS  Google Scholar 

  63. Kalds P, Zhou S, Gao Y, Cai B, Huang S, Chen Y, Wang X. Genetics of the phenotypic evolution in sheep: a molecular look at diversity-driving genes. Genet Sel Evol. 2022;54:1–27.

    Article  Google Scholar 

  64. Ghoreishifar SM, Eriksson S, Johansson AM, Khansefid M, Moghaddaszadeh-Ahrabi S, Parna N, Davoudi P, Javanmard A. Signatures of selection reveal candidate genes involved in economic traits and cold acclimation in five Swedish cattle breeds. Genet Sel Evol. 2020;52:52.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  65. Naserkheil M, Mehrban H, Lee D, Park MN. Genome-wide association study for carcass primal cut yields using single-step Bayesian approach in Hanwoo cattle. Front Genet. 2021;12:752424.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  66. Kaya, N., Aldhalaan, H., Al-Younes, B., Colak, D., Shuaib, T., Al-Mohaileb, F., Al-Sugair, A., Nester, M., Al-Yamani, S., Al-Bakheet, A., Al-Hashmi, N., Al-Sayed, M., Meyer, B., Jungbluth, H., Al-Owain, M., 2011. Phenotypical spectrum of cerebellar ataxia associated with a novel mutation in the CA8 gene, encoding carbonic anhydrase (CA) VIII. Am J Med Genet B Neuropsychiatr Genet 156b, 826-834.

  67. Raimondi C, Brash JT, Fantin A, Ruhrberg C. NRP1 function and targeting in neurovascular development and eye disease. Prog Retin Eye Res. 2016;52:64–83.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  68. Huang X, Qu R, Ouyang J, Zhong S, Dai J. An overview of the cytoskeleton-associated role of PDLIM5. Front Physiol. 2020;11:975.

    Article  PubMed  PubMed Central  Google Scholar 

  69. Li Z, Xu Y, Lin Y. Transcriptome analyses reveal genes of alternative splicing associated with muscle development in chickens. Gene. 2018;676:146–55.

    Article  CAS  PubMed  Google Scholar 

  70. Wang X, Monteagudo S, Lories R. EXT1 and EXT2 regulate chondrogenesis by modulation of WNT signaling. Osteoarthritis & Cartilage. 2018;26:S97.

    Article  Google Scholar 

  71. Wie J, Bharthur A, Wolfgang M, Narayanan V, Ramsey K, Aranda K, Zhang Q, Zhou Y, Ren D. Intellectual disability-associated UNC80 mutations reveal inter-subunit interaction and dendritic function of the NALCN channel complex. Nat Commun. 2020;11:3351.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  72. Alfei F, Kanev K, Hofmann M, Wu M, Ghoneim HE, Roelli P, Utzschneider DT, von Hoesslin M, Cullen JG, Fan Y, Eisenberg V, Wohlleber D, Steiger K, Merkler D, Delorenzi M, Knolle PA, Cohen CJ, Thimme R, Youngblood B, Zehn D. TOX reinforces the phenotype and longevity of exhausted T cells in chronic viral infection. Nature. 2019;571:265–9.

    Article  CAS  PubMed  Google Scholar 

  73. Khan O, Giles JR, McDonald S, Manne S, Ngiow SF, Patel KP, Werner MT, Huang AC, Alexander KA, Wu JE, Attanasio J, Yan P, George SM, Bengsch B, Staupe RP, Donahue G, Xu W, Amaravadi RK, Xu X, Karakousis GC, Mitchell TC, Schuchter LM, Kaye J, Berger SL, Wherry EJ. TOX transcriptionally and epigenetically programs CD8(+) T cell exhaustion. Nature. 2019;571:211–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  74. Scott AC, Dündar F, Zumbo P, Chandran SS, Klebanoff CA, Shakiba M, Trivedi P, Menocal L, Appleby H, Camara S, Zamarin D, Walther T, Snyder A, Femia MR, Comen EA, Wen HY, Hellmann MD, Anandasabapathy N, Liu Y, Altorki NK, Lauer P, Levy O, Glickman MS, Kaye J, Betel D, Philip M, Schietinger A. TOX is a critical regulator of tumour-specific T cell differentiation. Nature. 2019;571:270–4.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

Download references


The author would like to thank Prof. Congjiao Sun, Prof. Xiaolei Liu, and Prof. Peng Hu for their great support for this study. Chaoliang Wen and Shi Tang provided invaluable help with the coding of the microbiability. We thank the Prof. Wen Yao team for developing open sourcing the shinyCircos software. The first author, Dr. Yukun Zhang, wants to thank his parents (Zhijun Zhang and Cunxiu Sun) and girlfriend Dr. Yuanyuan Kong for their continuous support. We also gratefully acknowledge the support of all members for their assistance in animal sampling, measurements, and data input.


This work was supported by the National Key R&D Program of China (2021YFD1300901), National Natural Science Foundation of China (31960653 and 32260818), West Light Foundation of the Chinese Academy of Sciences, and National Joint Research on Improved Breeds of Livestock and Poultry (19220038).

Author information

Authors and Affiliations



W.M.W., Y.K.Z. and F.D.L conceived and designed the study. W.M.W., Y.K.Z., X.X.Z., C.L., L.F.Y., X.P.Y., W.H.L., J.B.C., D.Y.Z., Y.Z., X.L.L., L.M.Z., D.X., J.H.W., C.C.L., and F.D.L. were responsible for sample collection. X.X.W. and J.H.W. was responsible for the volatile fatty acid determinations. L.F.Y. processed the genotyping-by-sequencing data and performed SNP calling. Y.K.Z. were responsible for genomic data and 16S rRNA sequencing data analysis. Y.K.Z. and W.M.W. were conducted data integration and bioinformatic analysis. Y.K.Z., F.D.L and W.M.W. wrote the initial draft of the manuscript. Y.K.Z., F.D.L, W.M.W., Z.H.J., X.Z.D., G.H.S. revised and edited the manuscript. All authors read and approved the final manuscript.

Corresponding authors

Correspondence to Weimin Wang or Fadi Li.

Ethics declarations

Ethics approval and consent to participate

All animal experiments and procedures were conducted under the approval and guidance of the Animal Ethics Committee of Lanzhou University (Nos.: 2020-01 and 2021-02).

Consent for publication

Not applicable.

Competing interests

The authors declare no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1:

Table S1. Diet and Animal cohort information. Table S2. Prevalence (detection rates) and mean relative abundance of rumen microbiota genera in a single large-scale homogeneous population-based cohort of 1,150 sheep. Table S3. Keystone taxa (rumen microbiota genera) in a single large-scale homogeneous population-based cohort of 1150 sheep. Table S4. The microbial alpha diversity indices and beta diversity. Table S5. Heritability of rumen microbiota. Table S6. Single nucleotide polymorphisms (SNPs) information. Table S7. Genome-wide association of sheep genetic and rumen microbial variations. Table S8. Heritability and Microbiability of Body Weight Traits in Sheep. Table S9. microbiome-wide association study (MWAS). Table S10. GWAS for sheep body weight at 180 days old. Table S11. Microbiota mediated indirect scenarios. Table S12. 'recursive' scenarios.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and Permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Wang, W., Zhang, Y., Zhang, X. et al. Heritability and recursive influence of host genetics on the rumen microbiota drive body weight variance in male Hu sheep lambs. Microbiome 11, 197 (2023).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI:


  • Sheep
  • Host genetics
  • Rumen microbiota
  • Body weight
  • Heritability
  • Microbiota GWAS
  • Microbiability