Skip to main content
Advertisement
Browse Subject Areas
?

Click through the PLOS taxonomy to find articles in your field.

For more information about PLOS Subject Areas, click here.

  • Loading metrics

Genomic Regions Associated with Feed Efficiency Indicator Traits in an Experimental Nellore Cattle Population

  • Bianca Ferreira Olivieri ,

    bianca.olivieri@hotmail.com

    Affiliation Universidade Estadual Paulista, Faculdade de Ciências Agrárias e Veterinárias, Departamento de Zootecnia, Via de acesso Prof. Paulo Donato Castellane, s/no, CEP 14884-900 Jaboticabal, SP, Brazil

  • Maria Eugênia Zerlotti Mercadante,

    Affiliation Instituto de Zootecnia, Centro Avançado de Pesquisa Tecnológica do Agronegócio de Bovinos de Corte, Rodovia Carlos Tonanni, km 94, CEP 14.174-000, Sertãozinho, SP, Brazil

  • Joslaine Noely dos Santos Gonçalves Cyrillo,

    Affiliation Instituto de Zootecnia, Centro Avançado de Pesquisa Tecnológica do Agronegócio de Bovinos de Corte, Rodovia Carlos Tonanni, km 94, CEP 14.174-000, Sertãozinho, SP, Brazil

  • Renata Helena Branco,

    Affiliation Instituto de Zootecnia, Centro Avançado de Pesquisa Tecnológica do Agronegócio de Bovinos de Corte, Rodovia Carlos Tonanni, km 94, CEP 14.174-000, Sertãozinho, SP, Brazil

  • Sarah Figueiredo Martins Bonilha,

    Affiliation Instituto de Zootecnia, Centro Avançado de Pesquisa Tecnológica do Agronegócio de Bovinos de Corte, Rodovia Carlos Tonanni, km 94, CEP 14.174-000, Sertãozinho, SP, Brazil

  • Lucia Galvão de Albuquerque,

    Affiliation Universidade Estadual Paulista, Faculdade de Ciências Agrárias e Veterinárias, Departamento de Zootecnia, Via de acesso Prof. Paulo Donato Castellane, s/no, CEP 14884-900 Jaboticabal, SP, Brazil

  • Rafael Medeiros de Oliveira Silva,

    Affiliation Universidade Estadual Paulista, Faculdade de Ciências Agrárias e Veterinárias, Departamento de Zootecnia, Via de acesso Prof. Paulo Donato Castellane, s/no, CEP 14884-900 Jaboticabal, SP, Brazil

  • Fernando Baldi

    Affiliation Universidade Estadual Paulista, Faculdade de Ciências Agrárias e Veterinárias, Departamento de Zootecnia, Via de acesso Prof. Paulo Donato Castellane, s/no, CEP 14884-900 Jaboticabal, SP, Brazil

Correction

6 Feb 2017: Olivieri BF, Mercadante MEZ, Cyrillo JNdSG, Branco RH, Bonilha SFM, et al. (2017) Correction: Genomic Regions Associated with Feed Efficiency Indicator Traits in an Experimental Nellore Cattle Population. PLOS ONE 12(2): e0171845. https://doi.org/10.1371/journal.pone.0171845 View correction

Abstract

The objective of this study was to identify genomic regions and metabolic pathways associated with dry matter intake, average daily gain, feed efficiency and residual feed intake in an experimental Nellore cattle population. The high-density SNP chip (Illumina High-Density Bovine BeadChip, 777k) was used to genotype the animals. The SNP markers effects and their variances were estimated using the single-step genome wide association method. The (co)variance components were estimated by Bayesian inference. The chromosome segments that are responsible for more than 1.0% of additive genetic variance were selected to explore and determine possible quantitative trait loci. The bovine genome Map Viewer was used to identify genes. In total, 51 genomic regions were identified for all analyzed traits. The heritability estimated for feed efficiency was low magnitude (0.13±0.06). For average daily gain, dry matter intake and residual feed intake, heritability was moderate to high (0.43±0.05; 0.47±0.05, 0.18±0.05, respectively). A total of 8, 17, 14 and 12 windows that are responsible for more than 1% of the additive genetic variance for dry matter intake, average daily gain, feed efficiency and residual feed intake, respectively, were identified. Candidate genes GOLIM4, RFX6, CACNG7, CACNG6, CAPN8, CAPN2, AKT2, GPRC6A, and GPR45 were associated with feed efficiency traits. It was expected that the response to selection would be higher for residual feed intake than for feed efficiency. Genomic regions harboring possible QTL for feed efficiency indicator traits were identified. Candidate genes identified are involved in energy use, metabolism protein, ion transport, transmembrane transport, the olfactory system, the immune system, secretion and cellular activity. The identification of these regions and their respective candidate genes should contribute to the formation of a genetic basis in Nellore cattle for feed efficiency indicator traits, and these results would support the selection for these traits.

Introduction

Costs associated with animal feeding can be up to 50% of the total cost in beef cattle production systems, and there is growing interest concern in adopting strategies to reduce these costs [1]. So improving feed efficiency would increase both sustainability and profitability in the beef cattle industry. Selection for feed efficiency indicator traits in beef cattle can reduce production cost, decrease use of natural resources, and reduce impacts on the environment, optimizing production efficiency [2].

Feed efficiency (Kg gain/ Kg feed intake) showed a moderate positive correlation with weight gain and mature weight [3]. Thus, selection to improve feed efficiency could increase mature weight and increase energy maintenance requirements [4]. Also, residual feed intake (RFI), as proposed by Koch et al. [5], defined as the difference between actual feed intake and the feed intake required to meet maintenance requirements and growth, offers advantages over G:F since there is no correlation with performance traits [2, 6, 7, 8, 9].

Heritability estimates for RFI have been reported as moderate by many authors, ranging from 0.30 to 0.45 [3, 10, 11, 12, 13], suggesting that this trait can be improved by selection. However, RFI is difficult and expensive to measure, which often limits its implementation a selection criteria in beef cattle breeding programs. Recently, genome wide association studies (GWAS), using a high-density genotyping array, have been applied aiming to discover genomic regions associated with feed efficiency traits [11, 14, 15, 16]. It is important to highlight that most of these studies have been implemented with taurine breeds (Bos taurus), and there are few studies for indicine breeds (Bos indicus) under tropical conditions.

Therefore, in order to identify genomic regions associated with feed efficiency indicator traits, as well as to elucidate the genetic basis of them, it is important to encourage genomic studies with zebu animals, since zebu breeds are prevalent in herds under tropical and subtropical conditions. The objective of this study was to identify genomic regions and metabolic pathways associated with dry matter intake (DMI), average daily gain (ADG), feed efficiency (G:F) and residual feed intake (RFI) in an experimental Nellore cattle population.

Material and Methods

Data

This study was approved by ethics committee of the Faculty of Agrarian and Veterinary Sciences, Sao Paulo State University (UNESP).

The data set used in this study was provided by the APTA Beef Cattle Center—Institute of Animal Science (IZ), Sertãozinho, São Paulo, Brazil. Phenotypic information is animals born from 2004 to 2012. These animals belong to three experimental lines of Nellore cattle, which have been selected since 1978 for yearling weight: selection line (NeS) is a closed herd selected for higher yearling weight; traditional line (NeT) is submitted to the same selection criterion as NeS but, eventually, receives animals from NeS; and a control line (NeC) selected for average of yearling weight [17]. The analyzed data was obtained by feed efficiency tests performed from 2005 to 2013, consisting of 541 males and 355 females. Some of these animals (n = 683) were restricted to individual troughs, which offered daily feed and refusal was controlled, while the others (n = 213) were held in two collective pens equipped with the GrowSafe® feeding system.

After weaning, animals were kept in the test for during 83±15 days, preceded by 28 days of adaptation, for evaluation the feed intake and average daily gain (ADG). Animals were weighed every 14 days after 12 hours of fasting. The diet was formulated with 67% of total digestible nutrients (TDN) and 13% of crude protein (CP), allowing ADG of 1.0 kg/day.

Traits

In order to ensure ad libitum feed intake the food supply was adjusted daily, allowing refusals varying from 5–10% of offer. The following feed intake data was not considered in the analyses: for days when animals were handled outside of the facilities for a number hours, during failure (GrowSafe) and when no refusals were found. Diet dry matter percentage was determined from weekly samples of offer and refusals. The ADG in each test was considered as the linear regression coefficient of body weight (BW) on days in test (DIT): where, yi = BW in ith observation; α = intercept of regression equation corresponding to the initial BW; β = linear regression coefficient corresponding to ADG; DITi = days in test for ith observation; and εi = random error associated with each observation. Metabolic weight (BW0.75) was calculated as: , with α and β assuming the values obtained by the equation described above.

Feed efficiency (G:F) was calculated as the ratio of ADG to DMI. The residual feed intake (RFI) was considered as error of linear regression equation of dry matter intake on average daily gain and metabolic weight within each contemporary group (CG: sex, year of birth, and pen), as shown below: where, βT, βTA, and βTB are regression coefficients of classificatory variable test group (TG) and of interactions between TG and covariates ADG and BW0.75, respectively; and ε is RFI. The descriptive statistics for DMI, ADG, G:F and RFI are presented in Table 1.

thumbnail
Table 1. Descriptive statistics for dry matter intake (DMI), average daily gain (ADG), feed efficiency (G:F) and residual feed intake (RFI).

https://doi.org/10.1371/journal.pone.0164390.t001

DNA extraction

The extraction of DNA from blood samples was performed using the DNeasy Blood & Tissue Kit (Qiagen). The DNA purification was performed using a column containing silica fragments (column purification). Firstly, DNA binds to the membrane of the extraction column and then it was washed until it has high purity. At the end of the protocol, the DNA was eluted with buffer AE (blood). The amount of scanning and quality of the material obtained was taken with the use of a spectrophotometer apparatus (NanoDrop 1000, Thermo Scientific, USA, 2008). Quality was measured by absorption ratio A260 / A280. A ratio of less than 1.8 suggests contamination from protein, and more than 2.0 suggests RNA contamination.

Genotyped Animals

Animals were genotyped using the high-density SNP chip (Illumina High-Density BovineBeadChip 777,000 SNP). Markers with minor allele frequency (MAF) and call rate higher than 5% and 95%, respectively, were considered, as well as samples with a call rate higher than 93%. After quality control of markers, 438,874 SNPs for 689 animals were available.

Quantitative genetic analyses

The contemporary groups (CGs) included animals born on the same farm in the same year, and from the same management group as yearlings. The CGs with fewer than 3 records were eliminated from the analyses to maintain variability in the CGs. Records exceeding 3 standard deviations above or below the mean of each CG were excluded, to avoid the inclusion of possible measurement error or outliers. The model included the random additive animal effect, the fixed effects of CG, calving month, age of animal at beginning of the test (linear effect), and the dam age as co-variable (linear and quadratic effect). The (co)variance components were estimated using the single step genomic BLUP (ssGBLUP), under Bayesian inference [18]. The ssGBLUP is a modified version of the animal model (BLUP) with additive relationship matrix A-1 replaced by H-1 [19]: where A22 is a numerator relationship matrix for genotyped animals and G is the genomic relationship matrix created as described by VanRaden et al. (2009) [20]: where Z is the gene matrix containing allele frequency adjustment; D is the matrix that have the SNP weight (initially D = I); and, q is a weighting / standardization factor. According to Vitezica et al. (2011) [21], such factors can be obtained by ensuring that the G average diagonal is next to A22 The model can be represented by the following matrix equation: where y is the observations vector; β is the vector of fixed effects; a is the additive direct vector; X is known as incidence matrix; Z is the incidence genetic random effects additive direct matrix (the β vector associated with the y vector); e is the residual effect vector. The priori distributions of vectors y, a and e were given by: where H is the relationship coefficients matrix among animals obtained from the single-step analyzes (single-step); R is the residual variance matrix; I is the identity matrix; G is the genetic additive variance matrix and ⊗ is the Kronecker product. An inverted qui-square distribution was used for the prior values of the direct and residual genetic variances. A uniform distribution was used for the priori for the fixed effects a uniform distribution. A total of 3,342 animals were considered in the pedigree file. Analyses were performed using GIBBS2F90 [19, 22]. The a posteriori estimates were obtained using the application POSTGIBBSF90 [22].

The analyses were originated from chain lengths of 1,000,000 interactions, where the first 200,000 interactions were discarded. For parameter estimation, the samples were stored at each 100 cycles, building samples with 80,000 information. The data convergence was verified with the interactions versus evaluation graphic of sampled values and using the criteria proposed by Geweke; Heidelberger, Welch; and Raftery, Lewis [2325] through analysis package Bayesian Output Analysis (BOA) in the software R 2.9.0 (The R Development Core Team, 2009).

Genome-wide association analysis

Genome-wide association analysis for each trait was performed using the single-step GWAS (ssGWAS) methodology [26]. The same linear animal model used to estimate the (co)variance components was applied. The effects were decomposed in genotyped (ag) and ungenotyped (an) animals, as describe by Wang et al. [26], considering the effect of genotyped animals as: where Z is a matrix that relates genotypes of each locus and u is a vector of marker effects, and the variance of animal effects was assumed as: where D is a diagonal matrix of weights for variances of markers (D = I for GBLUP), σ2u is the genetic additive variance captured by each SNP marker when no weights are present, and G* is the weighted genomic relationship matrix.

The ratio of covariance of genetic effects (ag) and SNPs (u) is: sequentially: where λ is a variance ratio or a normalizing constant. According to VanRaden et al. [20], where M is the number of SNP and pi is the allele frequency of the second allele in the ith SNP. According to Stranden e Garrick [27], the markers effects can be described by:

The estimated SNP effects can be used to estimate the variance of each individual SNP effect [28] and apply a different weighting for each marker, such as:

The following iterative process described by Wang et al. [26] was used considering D to estimate the SNP effects:

  1. D = I,
  2. To calculate the matrix G = ZDZ’q
  3. To calculate GEBVs for all animals in data set using ssGBLUP;
  4. To calculate the SNP effect: ;
  5. To calculate the variance of each SNP:, where I is the ith marker;
  6. To normalize the values of SNPs to keep constant the additive genetic variance;
  7. Exit, or loop to step 2.

The effects of markers were obtained by 2 iterations from step 2 to 7. The percentage of genetic variance that is responsible for ith region was calculated as described by Wang et al. [26]: where ai is the genetic value of the ith region that consists of contiguous 10 consecutive SNPs, σ2a is the total genetic variance, Zj is the vector of gene content of the jth SNP for all individual, and ûj is the marker effect of the jth within the ith region. The results were presented by the proportion of variance expressed by each window of 10 SNPs. In addition, the genes located at ±500 Kb of each window were considered.

Search for genes

The chromosome segments that are responsible for more than 1.0% of additive genetic variance were selected to explore and determine possible quantitative trait loci. The bovine genome Map Viewer was used for identification of genes, available at "National Center for Biotechnology Information" (NCBI - http://www.ncbi.nlm.nih.gov) [29] in UMD3.1 version bovine genome and Ensembl Genome Browser (http://www.ensemble.org/index.html) [30]. The classification of genes for biological function, identification of metabolic pathways and enrichment of genes was performed on the website “The Database for Annotation, Visualization and Integrated Discovery (DAVID) v. 6.7” (http://david.abcc.ncifcrf.gov/) [31] and GeneCards (http://www.genecards.org/) [32].

Results and Discussion

Genetic parameter estimates

The (co)variance components and heritability estimates for DMI, ADG, G:F and RFI are showed in Table 2. The criteria used to diagnose the chain convergence indicated convergence of all estimated parameters. Thus, the burn-in period considered was sufficient to reach the convergence in all parameter estimates. The posterior marginal distributions of heritability estimates for feed efficiency indicator traits showed accurate values, tending to normal distribution. The symmetric distributions of central tendency statistics of analyzed traits indicated that the analyses are reliable.

thumbnail
Table 2. Heritability (h2) and (co)variance parameters estimate for dry matter intake (DMI), average daily gain (ADG), feed efficiency (G:F) and residual feed intake (RFI).

https://doi.org/10.1371/journal.pone.0164390.t002

The heritability estimated for G:F showed low magnitude. Grion et al. [12] and Ceacero et al. [33] found low heritability for the same trait (0.17 ± 0.07, 0.14 ± 0.06, respectively). For DMI and ADG, the estimated heritability was moderate to high. Baldi et al. [34] in a study with Nellore animals and modeling the weights with random regression models, found lower heritability estimates for average daily gain during the performance test (0.21). Similar heritability estimates (0.47) to those obtained in the present study for DMI was reported by Bolormaa et al. [11] with taurine and zebu animals. The heritability obtained for RFI showed moderate magnitude. Recently, Grion et al. [12], working with Nellore animals obtained a higher RFI heritability estimate (0.33± 0.10). Bolormaa et al. [11] also reported higher heritability estimates for RFI (0.36), in a study with nine herds of Bos Taurus, Bos indicus, and crossbreeds. Also, Silva et al. [35] reported similar heritability estimates for ADG (0.39±0.08) and DMI (0.43±0.08) in a study using the same experimental population. The results of this study pointed out that there is genetic variability in selecting for feed efficiency indicator traits in Nellore cattle. Thus, it is important to known whether there are more genes involved to better understand the genetic architecture of these traits.

The known genes found in the regions that accounted for more than 1.0% of additive genetic variance are presented in tables according to the studied trait. The results indicated a total of 8, 17, 14 and 12 different windows with known genes that are responsible for more than 1.0% of the genetic variance for DMI, ADG, G:F, and RFI, respectively.

Genomic regions

For DMI, eight genomic regions that are responsible for more than 1.0% of the additive genetic variance were found (Table 3; Fig 1). The window that is responsible for the most part of additive genetic variance for DMI was located in chromosome BTA4 where one candidate gene, called NXPH1, was found associated with the DMI. The RFX6, GPRC6A, FAM162B, KPNA5 and ZUFSP genes were identified in the window located in chromosome BTA9 at the 34 Mb position that is responsible for 3.40% of the additive genetic variance. The RFX6 gene is related to the regulation of transcription. Results of a study in rats have suggested that this gene acts on the differentiation of cells in insulin production [36].

thumbnail
Table 3. Genomics regions associated with dry matter intake (DMI) in Nellore cattle, percentage of additive genetic variance and candidate genes.

https://doi.org/10.1371/journal.pone.0164390.t003

thumbnail
Fig 1. Manhattan plot of additive genetic variance expressed by windows of 10 adjacent SNPs for DMI.

https://doi.org/10.1371/journal.pone.0164390.g001

In the associated window located in BTA15 at the 46 Mb position the OR2D2, OR2D3 and OR10A4 genes associated with DMI were found, which have been reported to play roles of olfactory receptors, coupled to G proteins, to recognize and mediate the olfactory signals in humans and rats [37, 38]. The ZNF214, ZNF215 and NLRP14 genes were found in the same chromosome. The NLRP14 gene encodes a protein with activity related to the immune system of the animal [39]. Three genomic regions located on chromosome BTA18 were found associated with DMI. Among the genes found in those regions, the CACNG7 and CACNG6 genes are associated with calcium channels [40].

The results of gene enrichment analysis and functional clusters showed that genes associated to DMI are involved in functions related to ion transport (Table A in S1 Tables). The metabolic pathways (Table 4) showed that those significant genes (p-value < 0.05) encode for a protein domain known as a zinc finger that acts as a structural element in proteins. Zinc is essential for several biochemical and cellular signaling pathways, participating in DNA and RNA synthesis and division and cell activation. Also, Zinc is essential for immune response [41].

thumbnail
Table 4. Enriched GO terms and KEGG pathways from DAVID software for DMI.

https://doi.org/10.1371/journal.pone.0164390.t004

For ADG, seventeen genomic regions located in 14 different chromosomes are responsible for more than 1.0% of the additive genetic variance (Table 5; Fig 2). In the window located in BTA1 at position 75 Mb, two candidate genes associated with ADG (MB21D2 and FGF12) were found. The OR52J3 and OR51A7 genes found in a window in BTA15 encode olfactory receptors. The CAPN8, CAPN2 and TP53BP2 genes are located in BTA16, where the first two encode subunits of the calpain enzyme. Calpain is an enzyme related to the tenderness of beef after slaughter [42]. Ten genes were located in the window, in BTA18 at position 49 Mb, that are responsible for the greatest proportion of additive variance (6.07%) of ADG. Among these genes, three of them (MAP3K10, CNTD2 and AKT2) encode kinases protein, which belongs to the largest family of proteins in eukaryotes. These genes also play an important part in intracellular communication, regulation and signal transduction, and catalyzed the phosphorylation of proteins by ATP transfer [43].

thumbnail
Table 5. Genomics regions associated with average daily gain (ADG) in Nellore cattle, percentage of additive genetic variance and candidate genes.

https://doi.org/10.1371/journal.pone.0164390.t005

thumbnail
Fig 2. Manhattan plot of additive genetic variance expressed by windows of 10 adjacent SNPs for ADG.

https://doi.org/10.1371/journal.pone.0164390.g002

The enrichment analyses for ADG revealed functional clusters related to the catabolism of macromolecules, transcription, protein kinase and binding nucleotides (Table B in S1 Tables). The gene enrichment and metabolic pathways showed that those significant genes (p <0.05) are related to endocytosis (extracellular transport) and merger of myoblasts, which are precursor cells of muscle fibers, and formation of syncytia (multinucleated cells) (Table 6).

thumbnail
Table 6. Enriched GO terms and KEGG pathways from DAVID software for ADG.

https://doi.org/10.1371/journal.pone.0164390.t006

Several genes were found in the 14 genomic regions that are responsible for more than 1% of the additive genetic variance of G:F (Table 7, Fig 3). The GOLIM4 and SERPINI1 genes were found in the window located in BTA1 at position 100 Mb. The GOLIM4 gene encodes an integral membrane protein localized in the Golgi apparatus, which is the main organelle in the secretory pathway in eukaryotic cells [44]. In chromosome BTA4 at position 26 Mb, the candidate genes PRPS1L1 and HDAC9 were found. The HDAC9 gene is a histone deacetylase enzyme, which is related the transcription factor. The position 78 Mb (BTA7) made an important contribution to the additive genetic variance for G:F (7.74%), but no candidate gene was found in this genomic region. The RASEF and FRMD3 genes were found in the window located in BTA8 at position 77 Mb. The RASEF gene encodes the GTP binding protein and binding calcium ions involved in the regulation of membrane traffic.

thumbnail
Table 7. Genomic regions associated with feed efficiency (G:F) in Nellore cattle, percentage of additive genetic variance and candidate genes.

https://doi.org/10.1371/journal.pone.0164390.t007

thumbnail
Fig 3. Manhattan plot of additive genetic variance expressed by windows of 10 adjacent SNPs for G:F.

https://doi.org/10.1371/journal.pone.0164390.g003

The associated candidate genes NREP, YTHDC2, and KCNN2 were found at the window at 2 Mb (BTA10). The KCNN2 gene is related to the activity in calcium/potassium channels. The Ca/K channels are present in the cytosol modulate tissue concentrations regulating the potential membrane and when present in the liver it is believed that they play a role in the response to metabolic stress [45]. Nine genes were found in the window located at the position 1.7 Mb in BTA11 (ACOXL, BUB1, TPC3, NPHP1, MALL, MRPS5, ZNF514, ZNF2 and PROM2). At the position 9 Mb in the same chromosome (BTA11), seven other candidate genes were listed: MRPS9, GPR45, TGFBRAP1, FHL2, TACR1, POLE4 and HK2. The GPR45 gene belongs to a family of receptors present in G proteins, and some studies in humans and mice have demonstrated that the GPR gene family is responsible for and is functional in the brain [46].

The SF3B3, COG4, FUK, ST3GAL2, and AARS genes were identified in BTA18, close to the 1700 kb region. The COG4 gene is responsible for in the Golgi apparatus, performing functions related to secretion [44]. In the same chromosome at position 60 Mb two genes were identified: ZNF677 and ZNF729. Both of them encode type Zinc Fingers binding proteins, which is characterized by the coordination and stabilization of zinc ions in several processes of ion exchange [16].

Seven candidate genes for G:F were found in BTA20 at the position 71 Mb: CEP72, SLC9A3, EXOC3, PDCD6, SDHA, CCDC127, and LRRC14B. In humans the SLC9A3 gene has the function of pH regulation, eliminating the acids produced by the metabolism and has proton antiporter activity, and solute carrier family 9 is involved in the exchange of sodium ions and protons, because comprises Na+/H+ exchanger proteins [47]. The window located at 5.6 Mb (BTA21) presented five candidate genes GABRG3, VIMP, CHSY1, LRRK1, and ALDH1A3. BTA22 presented six associated candidate genes FRMD4B, LMOD3, ARL6IP5, UBA3, TMF1, and EOGT. In BTA27 the associated genomic region did not show any gene because it is an intergenic region.

The results of gene enrichment and functional analysis reported clusters related to nucleotides, Golgi apparatus, protein transport and acetylation (Table C in S1 Tables). The metabolic pathways linked to carbohydrate metabolism, such as six-carbon polysaccharides fructose and mannose. These polysaccharides are intermediates in the glucose degradation process and glycolytic pathway, which are the main energy source for all cell types from mammals, being responsible for the ATP supply in aerobic and anaerobic conditions [48] (Table 8).

thumbnail
Table 8. Enriched GO terms and KEGG pathways from DAVID software for G:F.

https://doi.org/10.1371/journal.pone.0164390.t008

A total of 12 SNP windows distributed in BTA1, BTA4, BTA7, BTA8, BTA10, BTA18, BTA21, and BTA24 are responsible for more than 1% of the additive variance for RFI (Fig 4). Mujibi et al. [49], in a study of crossbreed beef cattle, reported 11 SNP windows associated with RFI. Santana et al. [16] found two SNP windows (located in BTA8 and BTA21) associated with RFI in a Nellore cattle population using two different densities of SNP markers.

thumbnail
Fig 4. Manhattan plot of additive genetic variance expressed by windows of 10 adjacent SNPs for RFI.

https://doi.org/10.1371/journal.pone.0164390.g004

Many candidate genes identified in regions associated with RFI are coding proteins involved in the ion transportation system. This system can consume more than 10% of the total energy used by ruminants. Thus, animals with a reduced energy expenditure on this of this system can redirect energy to spend on other processes, which may influence feed efficiency [16]. In the present study, the window located at 121 Mb (BTA1) is responsible for 5.67% of the additive genetic variance for RFI. This window harbored the ZIC1 e ZIC4 genes, which are related to ion transport, since codified protein subunits which depend on zinc ions for stability. Other genes are also related to ion transportation system, such as GPR98, KCNV2 (acting on potassium channels) and ZDHHC7.

A total of eight candidate genes was found in BTA4, of which three of them are linked to the perception of food: TAS2R3 and TAS2R4, which are related to taste through taste receptors, and OR9A, which encodes odor receptors. In another region of BTA4 (118 Mb) six candidates genes were found. The LMBR1 gene, located at 118 Mb (BTA4), has the function to encode membrane proteins. In addition, other candidate genes encodes proteins with biological functions associated with many cellular processes such as transcription (POLR3G and TBPL2), cellular secretion (KTN1 and GOLIM4), and transport across membranes (TMEM178B and AGK). No candidate gene was found at the position 18 Mb (BTA21) that is responsible for 2.18% of additive genetic variance for RFI (Table 9).

thumbnail
Table 9. Genomic regions associated with residual feed intake (RFI) in Nellore cattle, percentage of additive genetic variance and candidate genes.

https://doi.org/10.1371/journal.pone.0164390.t009

The results of enrichment pathway analysis indicated that genes associated with RFI (p <0.05) are related to sensory receptors that operate in food taste perception and receptors coupled to G proteins (Table 10). The enrichment by functional cluster for RFI showed that those genes are linked to cell membranes that are related to the thermodynamic equilibrium of cells. (Table D in S1 Tables). The most common mechanism to maintain the thermodynamic equilibrium of cells is ion exchange activity, maintaining cell differentiation potential through the sodium-potassium pump. For adequate pump operation, in order to maintain thermodynamic equilibrium, energy is needed (ATP to keep the sodium-potassium pump equilibrium). This process requires approximately 25% of the basal energy expenditure of an individual [50].

thumbnail
Table 10. Enriched GO terms and KEGG pathways from DAVID software for RFI.

https://doi.org/10.1371/journal.pone.0164390.t010

The GOLIM4 and SERPINI1 genes, located in BTA1, have been linked to more than one trait, such as G:F and RFI. It was observed that few SNP windows located nearby in the same chromosome were associated with more than one trait like in BTA8, where two SNP windows were associated with RFI and G:F. Seven genes related the zinc finger protein domain located in BTA18 were found in several nearby SNP windows associated with G:F and DMI. These results could be due to pleitropic effects, which means that the expression of different traits could be influenced by the same set of genes which acts in a coordinated manner to contribute to feed efficiency.

It is important to highlight that the results obtained in this study are also supported by previous studies. Rolf et al. [14] working with an Angus population, reported genomic regions associated with DMI, RFI and ADG close to those obtained in this study for the same traits. Bolormaa et al. [11] reported seven genomic regions near (at a maximum distance of 3 Mb) those identified in this study for RFI. Recently, Karisa et al. [51] reported a candidate gene (CYP2B) associated with RFI in BTA18 (49 Mb), at the same position that the present study identified candidate genes for ADG. Recently, in a study with a Nellore cattle population, De Oliveira et al. [52] reported differently located candidate genes for feed efficiency traits than found in this study (HRH4, ALDH7A1, APOA2, LIN7C, CXADR, ADAM12 and MAP7). However, the genes described in that study have similar gene ontology (immune system, energy and ion metabolism) to the genes reported in the present study.

The large number of genomic regions associated with feed efficiency traits obtained in this study should support a better understanding the genetic and physiological mechanisms that determine growth, feed intake and feed efficiency in zebu animals. The results demonstrate the probability that these traits have their expression controlled by many QTL with small individual effects, confirming their polygenic nature. The identification of relevant genes might be a difficult task, since the additive genetic variance contribution from each region or SNP window for many traits was lower than expected. Thus, strategies such as genomic selection that take into consideration the variability among all markers might be a more adequate alternative to improve these traits.

In recent years there has been growing concern about the contribution of the beef industry to climate change. Livestock in particular has been identified as a major contributor to global warming according to the FAO report [53]. Climate change and variability impacts on livestock productivity, and especially on the economic and political behavior of international markets, constituting the main current threat to beef exports in many countries. The ability to become a supplier country of reliable and safe food to the world should be combined with environmental sustainability. The results obtained in this study show that it is possible to improve beef cattle feed efficiency, trough selection using genomic information, and reach additional benefits to the environment by reducing greenhouse gas emissions. Finally, the data set used in this study belongs to a beef cattle research farm, which provides selected sires to commercial herds in many regions of Brazil. Therefore, the information found in this study will contribute to the selection for animals with better feed efficiency and increased environmental and social sustainability.

Conclusion

The results of this study pointed out that selection for feed efficiency indicator traits is feasible in Nellore cattle under tropical conditions. It is expected that the response to selection would be higher for RFI than for G:F. Several genomic regions harboring possible QTL for feed efficiency indicator traits were identified. The candidate genes identified are involved in energy and protein metabolism, ion transport, transmembrane transport, the olfactory system, the immune system, secretion (Golgi apparatus) and cellular activity (cell multiplication). The identification of these regions and their respective candidate genes should contribute to the formation of a genetic basis for Nellore feed efficiency indicator traits, and these results would support the selection for these traits.

Supporting Information

S1 Tables.

Table A—Functional cluster enrichment analysis via DAVID database for DMI.

Table B—Functional cluster enrichment analysis via DAVID database for ADG.

Table C—Functional cluster enrichment analysis via DAVID database for G:F.

Table D—Functional cluster enrichment analysis via DAVID database for RFI.

https://doi.org/10.1371/journal.pone.0164390.s001

(DOCX)

Author Contributions

  1. Conceptualization: BFO MEZM JNSGC RHB SFMB LGA RMOS FB.
  2. Formal analysis: BFO FB.
  3. Funding acquisition: LGA.
  4. Methodology: BFO MEZM FB.
  5. Project administration: BFO MEZM FB.
  6. Resources: MEZM JNSGC RHB SFMB LGA FB.
  7. Supervision: MEZM JNSGC RHB SFMB LGA RMOS FB.
  8. Visualization: BFO MEZM FB.
  9. Writing – original draft: BFO FB.
  10. Writing – review & editing: BFO RMOS FB.

References

  1. 1. Anderson RV, Rasby RJ, Klopfenstein TJ, and Clark RT. An evaluation of production and economic efficiency of two beef systems from calving to slaughter. J Anim Sci. 2005;83: 694–704. pmid:15705767
  2. 2. Arthur JPF and Herd RM. Residual feed intake in beef cattle. R Bras Zootec. 2008;37: 269–279.
  3. 3. Arthur PF, Archer JA, Johnston DJ, Herd RM, Richardson EC, and Parnell PF. Genetic and phenotypic variance and covariance componentes for feed intake, feed efficiency, and other postweaning traits in Angus cattle. J Anim Sci. 2001;79: 2805–2811. pmid:11768108
  4. 4. Tedeschi LO, Fox DG, Baker MJ, and Kirschten DP. Identifying differences in feed efficiency among group-fed cattle. J. Anim. Sci. 2006;84: 767–776. pmid:16478969
  5. 5. Koch RM, Swiger LA, Doyle C, and Gregory KE. Efficiency of feed use in beef cattle. J Anim Sci. 1963;22: 486–494.
  6. 6. Archer JA, Richardson EC, Herd RM, and Arthur PF. Potential for selection to improve efficiency of feed use in beef cattle: a review. Aust J Agric Res. 1999;50: 147–161.
  7. 7. Herd RM and Bishop SC. Genetic variation in residual feed intake and its association with other production traits in British Hereford cattle. Livest Sci. 2000;63: 111–119.
  8. 8. Arthur JPF, Renand G, and Krauss D. Genetic and phenotypic relationships among different measures of growth and feed efficiency in young Charolais bulls. Livest Prod Sci. 2001;68: 131–139.
  9. 9. Castilhos AM, Branco RH, Razook AG, Bonilha SFM, Mercadante MEZ, and Figueiredo LA. test post-weaning duration for performance, feed intake and feed efficiency in Nellore cattle. R Bras Zootec. 2011;40: 301–307.
  10. 10. Berry DP and Crowley JJ. Residual intake and body weight gain: A new measure of efficiency in growing cattle. J Anim Sci. 2011;90: 109–115. pmid:21890504
  11. 11. Bolormaa S, Pryce JE, Kemper KE, Hayes BJ, Zhang Y, Tier B, et al. Detection of quantitative trait loci in Bos indicus and Bos taurus cattle using genome-wide association studies. Genet Sel. Evol. 2013; 45:43. pmid:24168700
  12. 12. Grion AL, Mercadante MEZ, Cyrillo JNSG, Bonilha SFM, Magnani E, and Branco RH. Selection for feed efficiency traits and correlated genetic responses in feed intake and weight gain of Nellore cattle. J Anim Sci. 2014;92: 955–965. pmid:24492579
  13. 13. Santana MHA, Oliveira Junior GA, Gomes RC, Silva SL, Leme PR, Stella TR, et al. Genetic parameter estimates for feed efficiency and dry matter intake and their association with growth and carcass traits in Nellore cattle. Livest Sci. 2014;167: 80–85.
  14. 14. Rolf MM, Taylor JF, Schnabel RD, Mckay SD, Mcclure MC, Northcutt SL, et al. Genome-wide association analysis for feed efficiency in Angus cattle. Anim Genet. 2012;43: 367–374. pmcid: PMC3437496. pmid:22497295
  15. 15. Serão NVL, Beever JE, Faulkner DB, and Rogriguez-Zas SL. Integration of polygenic and individual SNP effects in genome-wide association analyses. Bioinformatics and Biomedicine Workshops (BIBMW), 2011 IEEE International Conference on. Atlanta, USA. 2011; 985–987.
  16. 16. Santana MHA, Utsunomiya YT, Neves HHR, Gomes RC, Garcia JF, Fukumasu H, et al. Genome-wide association analysis of feed intake and residual feed intake in Nellore cattle. BMC Genet. 2014;15: 21. pmid:24517472
  17. 17. Mercadante MEZ, Packer IU, Razook AG, Cyrillo JNSG, and Figueiredo LA. Direct responses to selection for yearling weight on reproductive performance of Nelore cows. J. Anim. Sci. 2003;81: 376–384. pmid:12643480
  18. 18. Gianola D and Fernando RL. Bayesian methods in animal breeding theory. J Anim Sci. 1986;63: 217–244.
  19. 19. Aguilar I, Misztal I, Johnson DL, Legarra A, Tsuruta S, and Lawlor TJ. Hot Topic: A unified approach to utilize phenotypic, full pedigree, and genomic information for genetic evaluation of Holstein final score. J Dairy Sci. 2010;93: 743–752. pmid:20105546
  20. 20. VanRaden PM, Van Tassell CP, Wiggans GR, Sonstegard TS, Schnabel RD, Taylor JF, et al. Invited review: Reliability of genomic predictions for North American Holstein bulls. J Dairy Sci. 2009;92: 16–24. pmid:19109259
  21. 21. Vitezica ZG, Aguilar I, Misztal I, and Legarra A. Bias in genomic predict for populations under selection. Genet Res. (Camb). 2011:93(5): 357–66. pmid:21767459.
  22. 22. Misztal I, Tsuruta S, Strabel T, Auvray B, Druet T, and Lee DH. BLUPF90 and related programs (BGF90). In The 7th World Congress Genetics Application Livestock Production. Montpellier, France. 2002: 28.
  23. 23. Geweke J. Evaluating the accuracy of sampling-based approaches to calculation of posterior moments. In: Bernardo JM, Berger JO, David AP, Smith AFM. (Ed.) Bayesian statistics. New York: Oxford University. 1992: pp 625–631.
  24. 24. Heidelberg P and Welch P. Simulation run lenght control in the presence of an initial transient. Oper Res. 1983;31: 1109–1114.
  25. 25. Raftery AL and Lewis S. Comment: one long run with diagnostics: implementation strategies for Markov chain Monte Carlo. Statist Sci. 1992;7: 493–497.
  26. 26. Wang H, Misztal I, Aguilar I, Legarra A, and Muir W. Genome-wide association mapping including phenotypes from relative without genotypes. Genet Res. 2012;94: 73–83.
  27. 27. Stranden I and Garrick DJ. Technical note: Derivation of equivalent computing algorithms for genomic predictions and reliabilities of animal merit. J Dairy Sci. 2009;92: 2971–2975. pmid:19448030
  28. 28. Zhang Z, Liu J, Ding X, Bijma P, de Koning DJ, and Zhang Q. Best linear unbiased prediction of genomic breeding values using a trait-specific marker derived relationship matrix. PLos One. 2010;9: 5.
  29. 29. NCBI—National Center for Biotechnology Information at http://www.ncbi.nlm.nih.gov.
  30. 30. Ensembl Genome Browser at http://www.ensembl.org/index.html.
  31. 31. DAVID—The Database for Annotation, Visualization and Integrated Discovery; v. 6.7” at http://david.abcc.ncifcrf.gov/.
  32. 32. GeneCards—The Human Gene Database at http://www.genecards.org/.
  33. 33. Ceacero TM, Mercadante MEZ, Cyrillo JNdSG, Canesin RC, Bonilha SFM, de Albuquerque LG. Phenotypic and Genetic Correlations of Feed Efficiency Traits with Growth and Carcass Traits in Nellore Cattle Selected for Postweaning Weight. PLoS ONE. 2016;11(8): e0161366. pmid:27537268
  34. 34. Baldi F, Albuquerque LG, Cyrillo JNSG, Branco RH, Oliveira BCJ, and Mercadante MEZ. Genetic parameter estimates for live weight and daily live weight gain obtained for Nellore bulls in a test station using different models. Livest Sci. 2012;144: 148–156.
  35. 35. Silva RMO, Fragomeni BO, Lourenco DAL, Magalhães AFB, Irano N, Carvalheiro R, et al. Accuracies of genomic prediction of feed efficiency traits using different prediction and validation methods in an experimental Nelore cattle population. J. Anim. Sci. 2016;94: 1.
  36. 36. Chandra V, Albagli-Curiel O, Hastoy B, Piccand J, Randriamampita C, Vaillant E, et al. RFX6 regulates insulin secretion by modulating Ca2+ homeostasis in human β cells. Cell Rep. 2014;9: 2206–2218. pmid:25497100
  37. 37. Malnic B, Godfrey PA, and Buck LB. The human olfactory receptor gene family. Proc Natl Acad Sci USA. 2004;101: 2584–2589. pmid:14983052 pmcid: PMC356993.
  38. 38. Godfrey PA, Malnic B, and Buck LB. The mouse olfactory receptor gene family. Proc Natl Acad Sci USA. 2004;101(7): 2156–2161. pmid:14769939
  39. 39. Tschopp J, Martinon F, and Burns K. NALPs: a novel protein family involved in inflammation. Nat Rev Mol Cell Biol. 2003;4: 95–104. pmid:12563287.
  40. 40. Burgess DL, Gefrides LA, Foreman PJ, and Noebels JL. A cluster of three novel Ca2+ channel gamma subunit genes on chromosome 19q13.4: evolution and expression profile of the gamma subunit gene family. Genom. 2001;71: 339–50. pmid:11170751.
  41. 41. Santos C and Fonseca J. ZINCO: fisiopatologia, clínica e nutrição. Associação Portuguesa de Nutrição Entérica e Parentérica. 2012;6: 1.
  42. 42. Goll DE, Thompson VF, Taylor RG, and Christiansen JA. Role of the calpain system in muscle growth. Biochim. 1992;74: 225–237.
  43. 43. Silva BV, Horta BAC, Alencastro RB, and Pinto AC. Proteínas quinases: características estruturais e inibidores químicos. Quím Nova. 2009;32: 453–462.
  44. 44. Starr T, Williams KF, and Storrie B. Both Post-Golgi and Intra-Golgi Cycling Affect the Distribution of the Golgi Phosphoprotein GPP130. Traffic. 2007;8: 1265–1279. pmid:17605763
  45. 45. Shakkottai VG, Regaya I, Wulff H, Fajloun Z, Tomita H, Fathallah M, et al. Design and characterization of a highly selective peptide inhibitor of the small conductance Calcium-activated K+ Channel, SkCa2. J Biol Chem. 2001;276: 43145–43151. pmid:11527975
  46. 46. Kawasawa Y, Kume K, Nakade S, Haga H, Izumi T, and Shimizu T. Brain-specific expression of novel G-protein-coupled receptors, with homologies to Xenopus PSP24 and human GPR45. Biochem Biophys Res Commun. 2000;276(3): 952–956. pmid:11027574.
  47. 47. Orlowski J and Grinstein S. Diversity of the mammalian sodium/proton exchanger SLC9 gene family. Pflügers Archiv. 2004;447: 549–565. pmid:12845533.
  48. 48. Machado UF. Transportadores de Glicose. Arquivo Brasileiro Endocrinologia e Metabologia. 1998;42(6): 413–421.
  49. 49. Mujibi FD, Nkrumah JD, Durunna OH, Stothard P, Mah J, Wang Z, et al. Accuracy of genomic breeding values for residual feed intake in crossbred beef cattle. J. Anim. Sci. 2011;89(11): 3353–3361. pmid:21642493
  50. 50. Leal D. P. Bioquímica Médica. Ed. Limusa:México. 2004.
  51. 51. Karisa B, Moore S, and Plastow G. Analysis of biological networks and biological pathways associated with residual feed intake in beef cattle. Anim Sci J. 2014;85: 374–387. pmid:24373146
  52. 52. De Oliveira PSN, Cesar ASM, Do Nascimento ML, Chaves AS, Tizioto PC, et al. Identification of genomic regions associated with feed efficiency in Nelore cattle. BMC Genet. 2014;15: 100. pmid:25257854
  53. 53. Steinfeld H, Gerber P, Wassenaar T, Castel V, Rosales M, and de Haan C. Livestock’s Long Shadow: Environmental Issues and Options. Rome: Food and Agriculture Organization of the United Nations; 2006.