BACKGROUND
brazilian livestock began with taurine cattle introduced from europe in the colonial period, around 450 years ago.
several years later, starting in  <dig> but mainly in  <dig> and  <dig>  a bos indicus cattle breed named nelore was introduced in brazil from india and the herd expanded to the majority of the territory due to its good adaptability to tropical climate.
it is estimated that  <dig> nelore cattle were introduced from india, the majority were descendants from  <dig> main sires named karvadi, taj mahal, kurupathy, golias, godhavari and rasta  <cit> .
the nelore, a breed with white or gray pelage and short horns, is nowadays the most representative breed in brazil and dominates the brazilian beef cattle production with more than 75 % of the total herd which equates to more than  <dig> million heads of nelore cattle.
artificial insemination is widely used in brazil  <cit>  and semen for only a few sires has been intensively used nationwide.
with the advent of the illumina bovinehd beadchip , a high density snp panel comprising over  <dig>  markers, we were able to measure the genomic diversity of the brazilian nelore from the main commercial sire families and characterize their genomic relationship and genetic profiles across a series of phenotypes.
traditional phenotype-based selection has been focused primarily on improving growth and maternal traits  <cit> .
however, marker-assisted or genomic selection is particularly valuable for traits that have late, onerous and expensive measurement such as meat quality, yield and feed efficiency  <cit> .
despite production advantages of bos indicus composite breeds in tropical and subtropical environments, the inferiority of beef quality, especially tenderness, is a worrying issue for producers  <cit> .
it is known that nelore breed show moderate heritability for meat quality  <cit>  as well as for growth-related traits  <cit> , which make the application of genomic tools feasible.
genome-wide association studies, or gwas, are today a common practice in order to find associations of genome regions to phenotypes of interest such as meat quality traits  <cit> .
however, gwas are sometimes elusive regarding the low power of association of the markers to the trait and also to the association of regions of the genome that lack apparent relation to the trait being studied  <cit> .
some meta-analysis approaches are appealing as they attempt to enrich gwas analyses with information from other sources.
one of these methods is the system’s biology awm/pcit methodology, which involves the construction of snp-based gene network interactions, integrating results from several gwas by means of the association weight matrices  and partial correlation and information theory  .
we took advantage of this methodology to explore growth and meat quality traits from  <dig> samples of brazilian nelore cattle, systematically sourced to represent  <dig> highly influential sires and their  <dig> half-sibs from across brazil.
further, we inferred a gene network focused on growth and meat quality traits and found several genes implicated in the regulation of these phenotypes.
RESULTS
genetic profile of brazilian nelore beef cattle
after quality control filters we retained genotypes of  <dig>  snps for linkage disequilibrium and haplotype block analysis, and  <dig>  tagsnps from  <dig> nelore animals, including  <dig> prominent sires and their  <dig> progeny .
these genotypes were coupled with genotypes of the same snp for a sample of  <dig> brahman, hereford and  <dig> angus cattle sourced from the bovine hapmap project  <cit> , and used to calculate a genomic relationship matrix  that was used in a principal component analysis   <cit> .
this analysis  shows that the first component, explaining 86 % of the variation in genotype profiles, separates the taurine  from the indicine  breeds, indicating that the genotyping data is reliable.
a grm considering only the  <dig> nelore animals and the  <dig>  tagsnps was also constructed and used for a pca analysis exclusive for the brazilian nelore sample.
the results for the pca  shows that although some half-sib families are well distinguished , the bulk of the samples were not well separated by the first component, indicating a lack of structuring.
the first and second components only explain  <dig>  % and  <dig>  % of the variation respectively.
it is interesting to note in fig.
1a that the distance from a given sire  to its corresponding cluster of half-sibs is visually symmetrical to the distance of the half-sib cluster to the center of the figure.
we suggest that this effect is due to the randomness diversity of the genetic effect of the dams that would bring more related individuals to the center of the figure, indicating that these sires are less related to its half-sib families and even less related to the bulk of more related individuals in the center.
the average genomic inbreeding coefficient from the grm  reports values of 2 %  for sires and 0 %  for progeny.
also, inbreeding coefficient derived from runs of homozygosity  reports values of 9 %  for sires, and 6 %  for progeny, showing that the  <dig> sires are more inbred than their progeny, and that the random genetic pool from the dams brought more variability into this sample.
the distribution of these inbreeding coefficients  can be seen in fig.
<dig>  additional file 1: figure s1c depicts the pedigree with the information of the dams showing that some dams participate in more than one family.table  <dig> genetic background of the sires
fig.
<dig> genetic profiling using the genomic relationship matrix.
a pca analysis of all  <dig> nelore .
inset shows a pca with other breeds from hapmap.
b heatmap and hierarchical clustering of the  <dig> nelore.
lateral palette colors represent the families ; upper color palette differentiates sires  from sibs ; shades of grey from the heatmap represent relationship similarities .
c pedigree view of the families showing the sires , sibs  and the lineage ancestral from father side
fig.
<dig> genomic inbreeding coefficients.
runs of homozygosity  estimations with a minimum length of  <dig> snps, and inbreeding coefficient derived from a genomic relationship matrix  for sires and its progeny
in addition to the pca, we performed a hierarchical cluster analysis of the grm and plotted the resulting heatmap .
this analysis shows that the  <dig> families are well distinguished in the heatmap  and the sire relatedness to its sibs can be distinguished in some cases by a white cross inside the central squares  or a white border line for some squares that follows the red lines in the upper palette representing the sires.
in overall the heatmap is greyish than what would be expected taking into account that the sires were selected to be more unrelated with each other, suggesting that the half-sib families are genetically related.
wright’s f statistic  was applied to evaluate the degree of genetic differentiation among each sire with its progeny as one population  and each lineage with its descendants as one population  yielding values of fst =  <dig>  and fst =  <dig> , respectively, which can be considered a small level of genetic differentiation in this nelore sample  <cit> .
the clustering analysis also put together the majority of the animals evidencing the relatedness of the families.
some darker grey shades can be seen within the progeny of individual ne <dig>  which is the individual that is more distant in the pca analysis, showing that this sire and its half-sibs are more unrelated to all other individuals.
although this individual sire is from the karvadi lineage and supposed to be related to the ne <dig>  ne <dig>  ne <dig> and ne <dig> , in the pca it is more distinguished from others.
overall, the lineage information obtained from pedigree was not useful to evidence any kind of relatedness between individuals with similar ancestry.
for example, nor the pca or the hierarchical clustering show any improvement in distinguishing related and unrelated animals regarding to its lineages .
the hierarchical clustering also does not cluster individuals according to its ancestral lineages .
one explanation for the failure of these analyses to differentiate the half-sib families is the fact that the families do not have the same number of animals and the same goes for the lineage ancestry of the families .
these distributions could affect the effectiveness of the grm to provide enough information to generate a good differentiation in the pca or hierarchical cluster analyses across lineages.
we obtained information of mitochondrial dna  of  <dig> of the  <dig> sires from a previous study  <cit> .
results show that  <dig>  of these have bos taurus maternal ancestry , which agrees with the hypothesis that most of the brazilian nelore herd was obtained by backcrossing with females of taurine origin  <cit> .
although two of the most divergent sires, ne <dig> and ne <dig> were of indicine maternal ancestry, in overall the maternal ancestry seems not to have influenced the pca results.
for a better characterization of the genomic data under study we analyzed the haplotype structure of the population and estimated the extent of the linkage disequilibrium , decay by physical genomic distance on autosomes, as the squared correlation coefficient  between snp pairs.
additional file 2: figure s <dig> shows that ld varies from r2 =  <dig>  between snps distance up to 1 kb for bta <dig> to r2 =  <dig>  between snps distance of 500 kb to 1 mb for bta <dig>  the average distance between snps for this ld study was  <dig>  kb, with large gaps of  <dig>  mb on bta <dig>  and of  <dig>  mb and  <dig>  mb on bta <dig> and on bta <dig>  respectively.
at the distance of 10 kb to 25 kb, the overall average r <dig> was ≈  <dig> , which is considered strong ld and can be used for qtl mapping  <cit> , while at distance of 25 kb to 50 kb, the overall average r <dig> was ≈  <dig> , which is sufficient to provide an accuracy of  <dig>  for estimation of genomic breeding value  <cit> .
a summary of the haplotype block description is provided on additional file 3: table s <dig>  from the snps considered for the ld and haplotype block study,  <dig>   were clustered into haplotype blocks constituted of  <dig> or more snps, spanning  <dig>  mb  of the total autosomal genome size.
there were  <dig>  haplotype blocks with overall average length of  <dig>  kb, bta <dig>  followed by bta <dig> , and bta <dig>  presented the longest haplotype blocks, the average number of snps into haplotype blocks was  <dig>  with a maximum of  <dig> snps .
the haplotype block size and distribution along autosomes is provided in additional file 4: figure s <dig>
the effective population size  was estimated to be around  <dig> animals in average per chromosome, with a confidence interval of  <dig> , one generation ago .
this number is higher than the ne obtained from pedigree records, varying from  <dig> to  <dig>  calculated from past periods of time   <cit>  and from another study also using pedigree records, that shows ne to be around  <dig> for a brazilian polled nelore breed  <cit> .
we estimated ne to be higher  <dig> generations ago , when the breed was introduced in brazil.
ne shows a decreasing pattern over the last generations especially for the last two .
gene networks for growth and meat quality traits
data inspection
we performed correlation analyses using the phenotype collection  with the vce software  <cit> .
although we obtained low heritability values  for lightness of meat , lightness of fat , ph , cooking loss  and dressing , higher estimates were observed for rib eye area , total carcass weight  and back fat thickness .
there was a negative correlation between bft and rea  and between bft and tcw  and a strong positive correlation between tcw and rea  which makes biological sense indicating that the data is consistent.
we calculated the influence  of inbreeding over the phenotypes, by estimating the inbreeding depression coefficients derived from froh and fgrm .
only bft, ph and clo appeared to be significantly affected by inbreeding although this was not a consensus among the two methods of estimation.
the results suggest that inbreeding did not strongly influentiate the traits.table  <dig> trait data and gwas results information
n
p ≤ 10−3
p ≤ 10−4
p ≤ 10−5
tcw, kg total carcass weight, dre, % dressing, rea, cm
<dig> rib eye area, bft, mm back fat thickness, lm lightness of meat, lf lightness of fat, ph ph, clo, % cooking loss.
n is the number of animals used, min.
is the minimum value, max is the maximum value, avg is the average value, sd is the standard deviation, h <dig> is the heritability , p ≤ 10− <dig> is the number of snps under this p-value, p ≤ 10− <dig> is the number of snps under this p-value, p ≤ 10− <dig> is the number of snps under this p-value.
inside parentheses is the fdr threshold  for each given significance  <cit>
*p <  <dig> ; **p <  <dig> ; aadjusted phenotype; sd standard deviation, se standard error, lm lightness of meat, lf lightness of fat, ph ph, clo cooking loss, dre dressing, rea rib eye area, tcw total carcass weight, bft back fat thickness
genome-wide association studies
the qxpak.
<dig> software was used  <cit>  to perform genome wide association studies  using the  <dig>  tagsnps for all eight phenotypes.
manhattan plots with the results are displayed in fig.
<dig>  when looking for markers with a significance threshold of p ≤ 10− <dig> we obtained several hundreds of significant tagsnps, varying from  <dig>  to  <dig> .
this high number of significant markers using this significance threshold indicates some inflation of the p-values .
if the threshold is pushed further to p ≤ 10− <dig> the number of significant tagsnps fall to a few hundreds and still if a more stringent threshold is used  we obtain just several dozens of significant tagsnps .
although we had some dozens of snps with very low p-values , the traits that had better heritabilities, like bft, rea and tcw showed a scattered distribution of significant snps over the chromosomes.
other traits showed more specific regions with significant snps like the ph trait, showing chromosome  <dig> with a visual peak of significant snps  in the region of  <dig> to 32 mb, pointing to the follicle stimulating hormone receptor  gene, but no evidence was found in the literature linking this gene to ph of meat.
other traits like lm have some diffuse peaks in chromosomes  <dig>   <dig>   <dig> and  <dig>  and lf with some less intense peaks in chromosomes  <dig> and  <dig>  but the lower heritability of these traits discouraged more profound analyses of these regions.fig.
<dig> gwas results for the eight traits.
a-h panels show the manhattan plots for: a total carcass weight ; b dressing % ; c rib eye area ; d back fat thickness ; e lightness of meat ; f lightness of fat ; g ph ; h cooking loss
we did a literature survey of the top  <dig> genes associated to the lowest p-value tagsnps for bft, rea and tcw.
for the later we searched for the involvement of the significant genes that could influence the total carcass weight phenotype.
we found the gene serotonin receptor  <dig>   associated to the marker bovinehd <dig>  which plays a role in appetite control  <cit>  and was found associated to birth weight in humans  <cit> .
we also found the leukemia inhibitory factor receptor  gene, associated to the marker bovinehd <dig>  which is involved in skeletal growth and bone formation and resorption in humans  <cit> .
for the bft phenotype we searched for genes involved in studies related to adipose tissue and fat.
we found the gene epha <dig> , associated to the marker bovinehd <dig>  found associated in a copy number variation  study to childhood obesity  <cit>  and the gene tmem <dig> , associated to the marker bovinehd <dig>  which is found up-regulated in brown adipose tissue during adipogenesis and myogenesis  <cit>  and also found associated in another gwas for bft  <cit> , discussed below.
we found no specific genes related to rea in the literature, that could be influencing muscle growth and myogenesis, associated to any tagsnps of the top  <dig> list.
in a previous work from our group  <cit> , gwas were performed using a bayesian approach, with a subset of the data used in this work .
we compared the set of snp, associated with significance below 10− <dig> in our work, against the ones associated in this previous work, to try to find genes simultaneously associated in both studies.
we found one common gene for lm,  <dig> common genes for lf,  <dig> common genes for bft and  <dig> genes for rea .
the accordance between both works for bft, rea and lf were of 26 %, 38 % and 38 % respectively, of all genes with significance below 10− <dig>  we could not estimate the accordance for the tcw trait as it was not used in this previous work.
to test for gene ontology  enrichment, we ordered all tagsnps according to the pleiotropy test  <cit>  used.
the ordered list of all corresponding genes was used as input to gorilla  <cit>  and results showed “cell adhesion” , “biological adhesion”  and “developmental growth”  as the three go terms from the biological process tree enriched with significant values after a fdr correction.
cell adhesion was also found significantly enriched for color of muscle and fat and for meat tenderness in our previous work  <cit> .
for the list of genes in each go, check additional file 8: spreadsheet s <dig>
gene networks
for a more holistic view of the genes influencing the eight phenotypes related to growth and meat quality we used the p-values and snp additive effects  obtained with the qxpak.
<dig> software to construct an awm .
tcw was chosen as key phenotype for the awm which was then used to calculate a pair-wise correlation matrix of all snp effects.
this correlation matrix was first filtered to maintain genes with correlation values equal or above  <dig> .
it was then used to run the pcit methodology, which generated a gene network that consisted of  <dig> genes  and  <dig> gene relationships  .
this is allegedly the network that contains genes more involved in growth traits because of the key phenotype , and also prone to be constituted of genes more involved with traits highly correlated to tcw, like rea and bft.
with less intensity, the network is supposed to contain genes involved with the other five phenotypes as they are less correlated to tcw.
additional file 10: figure s <dig> shows a heatmap and hierarchical cluster obtained with the standardized snp effects of the eight phenotypes that evidences that tcw, rea and bft are more correlated than ph, clo, lf, lm and dre.
for visualization we used cytoscape v <dig>  and applied several visual filters to identify highly connected nodes.
we suggest that the number of connections is a measure of the importance of the node to the function of the network as a whole.
the assumption is that the more connections a gene has the more likely is its role regulating and influencing the function of other genes in the network.
after applying the filters, the size of the node was changed to be proportional to its number of edges.
the smaller the number of connections of a node, the smaller was the size of the node and vice-versa.
this same rule is valid for the coloring: the more reddish, more connections and the more greenish, less connections.
the obtained network shows in its center, a cluster of nodes with larger size and stronger red color intensity showing that these nodes are more interconnected and are probably more involved in key biological mechanisms underlying the phenotypes involved in the network .
as an alternative type of analysis found in similar works , we decided to divide the network in a smaller sub-network to focus on a specific set of genes that have the highest degree of connections.
we used an algorithm to select the trio of transcription factors that mostly span the network with minimum redundancy  <cit> .
transcription factors  are known to regulate several cellular mechanisms in different degrees.
we arbitrarily decided to select three tf with most connections, which is a number that results in a sub-network with a reasonable sampling of the more important genes involved with the phenotypes.
we used cytoscape to create the sub-network by selecting the tf trio including the edges shared among its first degree nodes.
the resulting sub-network, with  <dig> nodes and  <dig> edges, is shown in fig.
<dig>  colors and node sizes represent the levels of connections inherited from the original network and triangular shaped nodes are representing the tf.
the trio of tf chosen by our algorithm is composed by the vdr , lhx <dig>  and zeb <dig>  genes.
the vdr, or vitamin d  receptor is a hormone receptor for vitamin d <dig> and is a transcription factor involved in a variety of metabolic pathways.
mutations in the vdr gene can lead to the rickets disease, characterized by growth retardation, bone deformity among other effects  <cit>  and this gene was also found related to bone density in mice  <cit> .
more interestingly, snps in the coding region of vdr were found to be associated with growth traits in bovine  <cit> .
the lhx <dig> gene is a transcription factor involved in developmental processes like neuronal and gonadal development  <cit>  and also to limb three-dimensional patterning and growth  <cit>  in mice.
taken together, the known functions of lhx <dig> are implicated in several important growth and developmental mechanisms which suggests it is an important player in the growth trait phenotype network.
the other tf from the most connected tf trio is the zeb <dig> , a zinc finger transcription factor which is involved in cancer and in regulation of the expression of several genes , this gene have been related to muscle differentiation and myogenesis  <cit> .fig.
<dig> gene network for growth and meat quality traits.
sub network for the more connected trio of transcription factors.
triangular shaped nodes show transcription factors.
greener nodes have lower connection levels than reddish nodes.
also, smaller nodes have lower connection levels than larger nodes.
labels represent gene symbols, exceptions are ensembl ids
from the set of genes with direct edges to the trio of tf, there are some genes with remarkable evidence in literature relating them to growth traits, specially bone and muscle growth and adipogenesis.
the asb <dig> , is a gene found to function in a negative regulation of muscle growth in salmon  <cit> .
the procr gene  was found to positively regulate other genes positively, in an osteoblastic cell line, functioning toward bone growth  <cit> .
the aldh9a <dig> , has a role on fetal growth in adult adipose tissue mass in bovine  <cit> .
cpeb <dig>  or cytoplasmic polyadenylation element binding protein  <dig>  is associated to human waist-hip ratio  <cit> .
the leucine-rich repeat lgi family, member  <dig>  lgi <dig>  is thought to have its function altered in obesity and to suppress adipogenesis  <cit> .
gpc <dig>  is implicated in bone growth  <cit> .
the fam3c gene, or family with sequence similarity  <dig>  member c was found associated to influence bone mineral density in humans  <cit> .
a gene, the iars , was found related to perinatal weak calf syndrome, which involves intrauterine growth retardation and low birth weights  <cit> .
there are several other genes related to cell growth and proliferation, like ergic <dig>   <cit> , khdrbs <dig>   <cit> , cep <dig>   <cit> , stim <dig>   <cit> , habp <dig>   <cit> .
from all  <dig> genes present in this sub-network,  <dig> were found also associated to rea  and bft  in the work from tizioto et al.,  <dig>  <cit> .
the rib eye area associated genes include the already discussed: vdr, asb <dig>  aldh9a <dig>  cpeb <dig>  gpc <dig>  fam3c, iars, khdrbs <dig>  cep <dig> and stim <dig>  only one gene was also found associated to bft, epb41l <dig> or erythrocyte membrane protein band  <dig> -like, which was found related to growth regulation  <cit> .
taken together, all the evidence found in literature of the relatedness of these genes to growth traits suggests that our methodology is reliable and the gene networks obtained shed a light on the biological mechanisms beyond growth and meat quality traits.
furthermore, as can be seen in table  <dig>  the snps that correspond to the best trio of tf have mostly p-values above 1e- <dig> meaning that this methodology captures snps that would not be taken as important by traditional gwas.table  <dig>
p-values of the snps assigned to the genes of the trio of tf
tcw, kg total carcass weight, dre, % dressing, rea, cm
<dig> rib eye area, bft, mm back fat thickness, lm lightness of meat, lf lightness of fat, ph ph, clo cooking loss
discussion
after the introduction of the nelore breed in brazil in the last century, there were some pedigree-based studies trying to characterize the founder lineages of the breed  <cit> , but no genomic approach has been made since then.
the introduction of the low and high density bovine genotyping panels, made possible the calculation and use of grm matrices  <cit>  in genomic similarity studies to create a genetic profile.
using only pedigree information from the  <dig> sires, the average inbreeding coefficient  was  <dig>  %  , which is similar to the fgrm reported for sires, but lower than the froh reported for sires.
comparing the inbreeding coefficients between sires and progeny, there is a clear decrease in inbreeding, which can only be attributed to the randomness of the genetic pool of the dams.
according to recent literature, fgrm was determined as the optimal method for estimating genomic inbreeding coefficients as it can distinguish between markers that are ibs  from markers that are ibd   <cit> .
the average froh for either sires or progeny estimated in this study are higher than the ones reported for a holstein population  which have shown effects of inbreeding depression for traits in dairy cattle  <cit> .
also, inbreeding depression appears not to have strongly affected the phenotypes .
this result agrees with the lower inbred estimation for the progeny compared to the sires.
we were not able to cross test this hypothesis by estimating inbreeding depression for the traits among sires, as we did not collect their phenotypes, altough we would expect a stronger influence of inbreeding in this case.
the average r <dig> values obtained in this study, considering snps with distance up to 150 kb apart, are higher than the values reported in the literature for nelore and indicine breeds.
beyond 150 kb distance, the average r <dig> values are similar to ones found in the literature  <cit> .
however, it is worth mentioning that it is difficult to make direct comparisons among ld studies, as ld levels vary due to many factors   <cit> .
in this study, both maternal and paternal haplotypes were considered in contrast to studies that only consider maternal haplotypes, as it is our understanding that the nelore population  presents higher influence from males than females and not considering paternal haplotypes would underestimate r <dig> values.
we believe that the backcrossing with taurine local breeds, evidenced by the mtdna origin, could have influenced the values of ld and ne estimated for this nelore sample.
there are some reports in the literature considering haplotype block structure in cattle, but none was found for nelore cattle.
these reports stated that with the increase in marker density, there will be more haplotype blocks of smaller sizes which will respond for higher genome coverage .
a recent study using the high density snp panel in crossbred cattle population reports 61 % of the genome covered by haplotype blocks, and a total of 78 % of snps clustered into haplotype blocks  <cit> .
the recent ne found for nelore in this study, is in agreement with reports of recent ne for holstein and swiss eringer breed varying from  <dig> to  <dig> individuals, and can go higher than  <dig> animals for past ne   <cit> .
our results show that the diversity among sire families and lineages is not high enough to correctly separate those in a pca, or hierarchical clustering analysis .
this is corroborated by the low fst values obtained for sire families and lineages subpopulations.
the uninformative values of ne and fst estimates were expected for the lineages as they do not represent isolated populations and should have been crosses among them.
some of the traits used in this work were already analyzed for genome wide associations elsewhere using part of the data and with a different approach  <cit> .
although the phenotype data is the same in both works, our sample size was larger and the software and algorithms used for association studies were different.
we have not found larger accordance than 38 %, between the set of associated genes with significance of 1e- <dig> to bft, rea and lm traits, studied in that work, although we could not verify this for tcw as this trait was not present in that work.
however, the quality and consistency of the genotyping data were found reliable since the pca against other breeds showed consistent results .
we also checked for the fluctuation of heterozygosity along the  <dig>  tagsnps and found reasonable smoothed values of heterozygosity fluctuating around  <dig>  to  <dig>  .
the traits correlations were also found reasonable .
furthermore some of the top tagsnps that had lowest p-values in the gwas, had the corresponding gene with functional biological correlation to the main growth traits, indicating that the gwas was able to extract some of the signal out of the noise from the data.
in addition, we also performed a gene enrichment test with gorilla, after ordering the tagsnps for the pleiotropic effect.
the pleiotropic test is a methodology to try to rank genes that are most influencing the traits in the awm.
the results showed go terms for biological and cell adhesion and developmental growth, which are again evidence that the terms are enriched for genes related mainly to growth which is the key phenotype of the awm.
cell adhesion is a biological process that is also found significantly enriched in the work of tizioto et al.,  <dig>  <cit> .
cell adhesion is a process involved in the binding of a cell to a substrate which sustains its growth and therefore is intrinsically related to the growth process.
with a similar objective of the pleiotropy test, the construction of the gene networks with the awm/pcit  is an alternative to analyze pure gwas results as this methodology gathers the  <dig> gwas results and constructs a gene network that influences the traits in a simultaneous way, but weighting accordingly to the correlations within the awm.
we decided to select tcw as the key trait for the awm because of its reasonable heritability estimate and its congruent correlations to rea and bft, together with the economic importance of growth traits for livestock.
as the resulting network is too large , we used a heuristic to divide the network into a sub-network in a way that it would be feasible to analyze and still gather most of the biological importance regarding these traits.
this heuristic is implemented in an algorithm that selects the trio of tf that most span the network and its first degree nodes.
as the coloring and shape of the nodes of the resulting sub-network are inherited from the original network, it can be noticed a majority of reddish and larger nodes  indicating that the objective of the heuristic was achieved .
also, a survey in the literature showed correlations between the trio of tf and several nodes from the sub-network to growth traits, suggesting that our methodology is correct.
it is remarkable to note that some genes like vdr  was already associated to growth in a bovine gwas  <cit>  and many others have relations to the biology of growth, like proliferation and growth of muscle and bones and adipogenesis .
although comparing the results of single snp association to the work of tizioto et al.
<dig> did not lead to more than 38 % of common associated genes, the genes from the sub-network are in 50 % of accordance with the genes associated to bft and rea from that work.
this suggests that the methodology that created this network was able to select genes that are more related to the traits involved with growth.
we compared the genes from the sub-network to genes found in other studies  that also used the awm/pcit methodology and similar datasets .
although there were no great correlations, some genes in these studies were found to have similarities or to encode the same protein domains to some of the genes from the sub-network: metalloprotease and ankyrin domains , the rab6b oncogene  <cit> , the mer proto-oncogene tyrosine kinase  and the slc solute carrier gene family   <cit> .
taken together, these results are suggesting that the gene networks obtained are related to growth and meat quality traits and their genes should be thoroughly inspected to try to discover the biological mechanism beneath growth and meat quality phenotypes.
CONCLUSIONS
we used high density snp panels to genotype  <dig> sires, and its half-sib families, totaling  <dig> animals.
the sires were previously selected in order to represent most of the variability of the nelore beef cattle in brazil.
we believe this is the first work that used a genomic approach in order to try to investigate the amount of diversity of the brazilian nelore breed by investigating its ancestral lineages.
results showed that the population studied was not structured enough in order to differentiate families using ancestral information of sires or lineages.
eight phenotypes related to growth and meat quality were used in whole genome association studies and the results were used to construct gene networks focused on growth, using the methodology of awm/pcit.
literature surveys showed that both the gwas and the gene networks constructed had significant snps associated to genes related to growth in former studies like the vdr, lhx <dig>  cpeb <dig>  asb <dig> and many others.
these results suggests that the methodology used to construct the gene networks can be used as an alternative approach to standard gwas, in order to reach for novel information and to try to understand the biological mechanisms and gene compositions that leads to complex phenotypes, like growth in beef cattle.

 brazilian livestock began with taurine cattle introduced from europe in the colonial period, around 450 years ago. several years later, starting in but mainly in and a bos indicus cattle breed named nelore was introduced in brazil from india and the herd expanded to the majority of the territory due to its good adaptability to tropical climate. it is estimated that nelore cattle were introduced from india, the majority were descendants from main sires named karvadi, taj mahal, kurupathy, golias, godhavari and rasta. the nelore, a breed with white or gray pelage and short horns, is nowadays the most representative breed in brazil and dominates the brazilian beef cattle production with more than 75 % of the total herd which equates to more than million heads of nelore cattle. artificial insemination is widely used in brazil and semen for only a few sires has been intensively used nationwide. with the advent of the illumina bovinehd beadchip, a high density snp panel comprising over markers, we were able to measure the genomic diversity of the brazilian nelore from the main commercial sire families and characterize their genomic relationship and genetic profiles across a series of phenotypes. traditional phenotype-based selection has been focused primarily on improving growth and maternal traits. however, marker-assisted or genomic selection is particularly valuable for traits that have late, onerous and expensive measurement such as meat quality, yield and feed efficiency. despite production advantages of bos indicus composite breeds in tropical and subtropical environments, the inferiority of beef quality, especially tenderness, is a worrying issue for producers. it is known that nelore breed show moderate heritability for meat quality as well as for growth-related traits, which make the application of genomic tools feasible. genome-wide association studies, or gwas, are today a common practice in order to find associations of genome regions to phenotypes of interest such as meat quality traits. however, gwas are sometimes elusive regarding the low power of association of the markers to the trait and also to the association of regions of the genome that lack apparent relation to the trait being studied. some meta-analysis approaches are appealing as they attempt to enrich gwas analyses with information from other sources. one of these methods is the system’s biology awm/pcit methodology, which involves the construction of snp-based gene network interactions, integrating results from several gwas by means of the association weight matrices and partial correlation and information theory. we took advantage of this methodology to explore growth and meat quality traits from samples of brazilian nelore cattle, systematically sourced to represent highly influential sires and their half-sibs from across brazil. further, we inferred a gene network focused on growth and meat quality traits and found several genes implicated in the regulation of these phenotypes.  genetic profile of brazilian nelore beef cattle after quality control filters we retained genotypes of snps for linkage disequilibrium and haplotype block analysis, and tagsnps from nelore animals, including prominent sires and their progeny. these genotypes were coupled with genotypes of the same snp for a sample of brahman, hereford and angus cattle sourced from the bovine hapmap project, and used to calculate a genomic relationship matrix that was used in a principal component analysis. this analysis shows that the first component, explaining 86 % of the variation in genotype profiles, separates the taurine from the indicine breeds, indicating that the genotyping data is reliable. a grm considering only the nelore animals and the tagsnps was also constructed and used for a pca analysis exclusive for the brazilian nelore sample. the results for the pca shows that although some half-sib families are well distinguished, the bulk of the samples were not well separated by the first component, indicating a lack of structuring. the first and second components only explain  % and  % of the variation respectively. it is interesting to note in fig. 1a that the distance from a given sire to its corresponding cluster of half-sibs is visually symmetrical to the distance of the half-sib cluster to the center of the figure. we suggest that this effect is due to the randomness diversity of the genetic effect of the dams that would bring more related individuals to the center of the figure, indicating that these sires are less related to its half-sib families and even less related to the bulk of more related individuals in the center. the average genomic inbreeding coefficient from the grm reports values of 2 % for sires and 0 % for progeny. also, inbreeding coefficient derived from runs of homozygosity reports values of 9 % for sires, and 6 % for progeny, showing that the sires are more inbred than their progeny, and that the random genetic pool from the dams brought more variability into this sample. the distribution of these inbreeding coefficients can be seen in fig.  additional file 1: figure s1c depicts the pedigree with the information of the dams showing that some dams participate in more than one family.table genetic background of the sires fig. genetic profiling using the genomic relationship matrix. a pca analysis of all nelore. inset shows a pca with other breeds from hapmap. b heatmap and hierarchical clustering of the nelore. lateral palette colors represent the families; upper color palette differentiates sires from sibs; shades of grey from the heatmap represent relationship similarities. c pedigree view of the families showing the sires, sibs and the lineage ancestral from father side fig. genomic inbreeding coefficients. runs of homozygosity estimations with a minimum length of snps, and inbreeding coefficient derived from a genomic relationship matrix for sires and its progeny in addition to the pca, we performed a hierarchical cluster analysis of the grm and plotted the resulting heatmap. this analysis shows that the families are well distinguished in the heatmap and the sire relatedness to its sibs can be distinguished in some cases by a white cross inside the central squares or a white border line for some squares that follows the red lines in the upper palette representing the sires. in overall the heatmap is greyish than what would be expected taking into account that the sires were selected to be more unrelated with each other, suggesting that the half-sib families are genetically related. wright’s f statistic was applied to evaluate the degree of genetic differentiation among each sire with its progeny as one population and each lineage with its descendants as one population yielding values of fst =  and fst = , respectively, which can be considered a small level of genetic differentiation in this nelore sample. the clustering analysis also put together the majority of the animals evidencing the relatedness of the families. some darker grey shades can be seen within the progeny of individual ne which is the individual that is more distant in the pca analysis, showing that this sire and its half-sibs are more unrelated to all other individuals. although this individual sire is from the karvadi lineage and supposed to be related to the ne ne ne and ne, in the pca it is more distinguished from others. overall, the lineage information obtained from pedigree was not useful to evidence any kind of relatedness between individuals with similar ancestry. for example, nor the pca or the hierarchical clustering show any improvement in distinguishing related and unrelated animals regarding to its lineages. the hierarchical clustering also does not cluster individuals according to its ancestral lineages. one explanation for the failure of these analyses to differentiate the half-sib families is the fact that the families do not have the same number of animals and the same goes for the lineage ancestry of the families. these distributions could affect the effectiveness of the grm to provide enough information to generate a good differentiation in the pca or hierarchical cluster analyses across lineages. we obtained information of mitochondrial dna of of the sires from a previous study. results show that of these have bos taurus maternal ancestry, which agrees with the hypothesis that most of the brazilian nelore herd was obtained by backcrossing with females of taurine origin. although two of the most divergent sires, ne and ne were of indicine maternal ancestry, in overall the maternal ancestry seems not to have influenced the pca results. for a better characterization of the genomic data under study we analyzed the haplotype structure of the population and estimated the extent of the linkage disequilibrium, decay by physical genomic distance on autosomes, as the squared correlation coefficient between snp pairs. additional file 2: figure s shows that ld varies from r2 =  between snps distance up to 1 kb for bta to r2 =  between snps distance of 500 kb to 1 mb for bta the average distance between snps for this ld study was  kb, with large gaps of  mb on bta and of  mb and  mb on bta and on bta respectively. at the distance of 10 kb to 25 kb, the overall average r was ≈ , which is considered strong ld and can be used for qtl mapping, while at distance of 25 kb to 50 kb, the overall average r was ≈ , which is sufficient to provide an accuracy of for estimation of genomic breeding value. a summary of the haplotype block description is provided on additional file 3: table s from the snps considered for the ld and haplotype block study, were clustered into haplotype blocks constituted of or more snps, spanning  mb of the total autosomal genome size. there were haplotype blocks with overall average length of  kb, bta followed by bta, and bta presented the longest haplotype blocks, the average number of snps into haplotype blocks was with a maximum of snps. the haplotype block size and distribution along autosomes is provided in additional file 4: figure s the effective population size was estimated to be around animals in average per chromosome, with a confidence interval of, one generation ago. this number is higher than the ne obtained from pedigree records, varying from to calculated from past periods of time and from another study also using pedigree records, that shows ne to be around for a brazilian polled nelore breed. we estimated ne to be higher generations ago, when the breed was introduced in brazil. ne shows a decreasing pattern over the last generations especially for the last two. gene networks for growth and meat quality traits data inspection we performed correlation analyses using the phenotype collection with the vce software. although we obtained low heritability values for lightness of meat, lightness of fat, ph, cooking loss and dressing, higher estimates were observed for rib eye area, total carcass weight and back fat thickness. there was a negative correlation between bft and rea and between bft and tcw and a strong positive correlation between tcw and rea which makes biological sense indicating that the data is consistent. we calculated the influence of inbreeding over the phenotypes, by estimating the inbreeding depression coefficients derived from froh and fgrm. only bft, ph and clo appeared to be significantly affected by inbreeding although this was not a consensus among the two methods of estimation. the results suggest that inbreeding did not strongly influentiate the traits.table trait data and gwas results information n p ≤ 10−3 p ≤ 10−4 p ≤ 10−5 tcw, kg total carcass weight, dre, % dressing, rea, cm rib eye area, bft, mm back fat thickness, lm lightness of meat, lf lightness of fat, ph ph, clo, % cooking loss. n is the number of animals used, min. is the minimum value, max is the maximum value, avg is the average value, sd is the standard deviation, h is the heritability, p ≤ 10− is the number of snps under this p-value, p ≤ 10− is the number of snps under this p-value, p ≤ 10− is the number of snps under this p-value. inside parentheses is the fdr threshold for each given significance *p < ; **p < ; aadjusted phenotype; sd standard deviation, se standard error, lm lightness of meat, lf lightness of fat, ph ph, clo cooking loss, dre dressing, rea rib eye area, tcw total carcass weight, bft back fat thickness genome-wide association studies the qxpak. software was used to perform genome wide association studies using the tagsnps for all eight phenotypes. manhattan plots with the results are displayed in fig.  when looking for markers with a significance threshold of p ≤ 10− we obtained several hundreds of significant tagsnps, varying from to. this high number of significant markers using this significance threshold indicates some inflation of the p-values. if the threshold is pushed further to p ≤ 10− the number of significant tagsnps fall to a few hundreds and still if a more stringent threshold is used we obtain just several dozens of significant tagsnps. although we had some dozens of snps with very low p-values, the traits that had better heritabilities, like bft, rea and tcw showed a scattered distribution of significant snps over the chromosomes. other traits showed more specific regions with significant snps like the ph trait, showing chromosome with a visual peak of significant snps in the region of to 32 mb, pointing to the follicle stimulating hormone receptor gene, but no evidence was found in the literature linking this gene to ph of meat. other traits like lm have some diffuse peaks in chromosomes and and lf with some less intense peaks in chromosomes and but the lower heritability of these traits discouraged more profound analyses of these regions.fig. gwas results for the eight traits. a-h panels show the manhattan plots for: a total carcass weight; b dressing %; c rib eye area; d back fat thickness; e lightness of meat; f lightness of fat; g ph; h cooking loss we did a literature survey of the top genes associated to the lowest p-value tagsnps for bft, rea and tcw. for the later we searched for the involvement of the significant genes that could influence the total carcass weight phenotype. we found the gene serotonin receptor associated to the marker bovinehd which plays a role in appetite control and was found associated to birth weight in humans. we also found the leukemia inhibitory factor receptor gene, associated to the marker bovinehd which is involved in skeletal growth and bone formation and resorption in humans. for the bft phenotype we searched for genes involved in studies related to adipose tissue and fat. we found the gene epha, associated to the marker bovinehd found associated in a copy number variation study to childhood obesity and the gene tmem, associated to the marker bovinehd which is found up-regulated in brown adipose tissue during adipogenesis and myogenesis and also found associated in another gwas for bft, discussed below. we found no specific genes related to rea in the literature, that could be influencing muscle growth and myogenesis, associated to any tagsnps of the top list. in a previous work from our group, gwas were performed using a bayesian approach, with a subset of the data used in this work. we compared the set of snp, associated with significance below 10− in our work, against the ones associated in this previous work, to try to find genes simultaneously associated in both studies. we found one common gene for lm, common genes for lf, common genes for bft and genes for rea. the accordance between both works for bft, rea and lf were of 26 %, 38 % and 38 % respectively, of all genes with significance below 10− we could not estimate the accordance for the tcw trait as it was not used in this previous work. to test for gene ontology enrichment, we ordered all tagsnps according to the pleiotropy test used. the ordered list of all corresponding genes was used as input to gorilla and results showed “cell adhesion”, “biological adhesion” and “developmental growth” as the three go terms from the biological process tree enriched with significant values after a fdr correction. cell adhesion was also found significantly enriched for color of muscle and fat and for meat tenderness in our previous work. for the list of genes in each go, check additional file 8: spreadsheet s gene networks for a more holistic view of the genes influencing the eight phenotypes related to growth and meat quality we used the p-values and snp additive effects obtained with the qxpak. software to construct an awm. tcw was chosen as key phenotype for the awm which was then used to calculate a pair-wise correlation matrix of all snp effects. this correlation matrix was first filtered to maintain genes with correlation values equal or above. it was then used to run the pcit methodology, which generated a gene network that consisted of genes and gene relationships. this is allegedly the network that contains genes more involved in growth traits because of the key phenotype, and also prone to be constituted of genes more involved with traits highly correlated to tcw, like rea and bft. with less intensity, the network is supposed to contain genes involved with the other five phenotypes as they are less correlated to tcw. additional file 10: figure s shows a heatmap and hierarchical cluster obtained with the standardized snp effects of the eight phenotypes that evidences that tcw, rea and bft are more correlated than ph, clo, lf, lm and dre. for visualization we used cytoscape v and applied several visual filters to identify highly connected nodes. we suggest that the number of connections is a measure of the importance of the node to the function of the network as a whole. the assumption is that the more connections a gene has the more likely is its role regulating and influencing the function of other genes in the network. after applying the filters, the size of the node was changed to be proportional to its number of edges. the smaller the number of connections of a node, the smaller was the size of the node and vice-versa. this same rule is valid for the coloring: the more reddish, more connections and the more greenish, less connections. the obtained network shows in its center, a cluster of nodes with larger size and stronger red color intensity showing that these nodes are more interconnected and are probably more involved in key biological mechanisms underlying the phenotypes involved in the network. as an alternative type of analysis found in similar works, we decided to divide the network in a smaller sub-network to focus on a specific set of genes that have the highest degree of connections. we used an algorithm to select the trio of transcription factors that mostly span the network with minimum redundancy. transcription factors are known to regulate several cellular mechanisms in different degrees. we arbitrarily decided to select three tf with most connections, which is a number that results in a sub-network with a reasonable sampling of the more important genes involved with the phenotypes. we used cytoscape to create the sub-network by selecting the tf trio including the edges shared among its first degree nodes. the resulting sub-network, with nodes and edges, is shown in fig.  colors and node sizes represent the levels of connections inherited from the original network and triangular shaped nodes are representing the tf. the trio of tf chosen by our algorithm is composed by the vdr, lhx and zeb genes. the vdr, or vitamin d receptor is a hormone receptor for vitamin d and is a transcription factor involved in a variety of metabolic pathways. mutations in the vdr gene can lead to the rickets disease, characterized by growth retardation, bone deformity among other effects and this gene was also found related to bone density in mice. more interestingly, snps in the coding region of vdr were found to be associated with growth traits in bovine. the lhx gene is a transcription factor involved in developmental processes like neuronal and gonadal development and also to limb three-dimensional patterning and growth in mice. taken together, the known functions of lhx are implicated in several important growth and developmental mechanisms which suggests it is an important player in the growth trait phenotype network. the other tf from the most connected tf trio is the zeb, a zinc finger transcription factor which is involved in cancer and in regulation of the expression of several genes, this gene have been related to muscle differentiation and myogenesis.fig. gene network for growth and meat quality traits. sub network for the more connected trio of transcription factors. triangular shaped nodes show transcription factors. greener nodes have lower connection levels than reddish nodes. also, smaller nodes have lower connection levels than larger nodes. labels represent gene symbols, exceptions are ensembl ids from the set of genes with direct edges to the trio of tf, there are some genes with remarkable evidence in literature relating them to growth traits, specially bone and muscle growth and adipogenesis. the asb, is a gene found to function in a negative regulation of muscle growth in salmon. the procr gene was found to positively regulate other genes positively, in an osteoblastic cell line, functioning toward bone growth. the aldh9a, has a role on fetal growth in adult adipose tissue mass in bovine. cpeb or cytoplasmic polyadenylation element binding protein is associated to human waist-hip ratio. the leucine-rich repeat lgi family, member lgi is thought to have its function altered in obesity and to suppress adipogenesis. gpc is implicated in bone growth. the fam3c gene, or family with sequence similarity member c was found associated to influence bone mineral density in humans. a gene, the iars, was found related to perinatal weak calf syndrome, which involves intrauterine growth retardation and low birth weights. there are several other genes related to cell growth and proliferation, like ergic, khdrbs, cep, stim, habp. from all genes present in this sub-network, were found also associated to rea and bft in the work from tizioto et al.,. the rib eye area associated genes include the already discussed: vdr, asb aldh9a cpeb gpc fam3c, iars, khdrbs cep and stim only one gene was also found associated to bft, epb41l or erythrocyte membrane protein band -like, which was found related to growth regulation. taken together, all the evidence found in literature of the relatedness of these genes to growth traits suggests that our methodology is reliable and the gene networks obtained shed a light on the biological mechanisms beyond growth and meat quality traits. furthermore, as can be seen in table  the snps that correspond to the best trio of tf have mostly p-values above 1e- meaning that this methodology captures snps that would not be taken as important by traditional gwas.table p-values of the snps assigned to the genes of the trio of tf tcw, kg total carcass weight, dre, % dressing, rea, cm rib eye area, bft, mm back fat thickness, lm lightness of meat, lf lightness of fat, ph ph, clo cooking loss discussion after the introduction of the nelore breed in brazil in the last century, there were some pedigree-based studies trying to characterize the founder lineages of the breed, but no genomic approach has been made since then. the introduction of the low and high density bovine genotyping panels, made possible the calculation and use of grm matrices in genomic similarity studies to create a genetic profile. using only pedigree information from the sires, the average inbreeding coefficient was  %, which is similar to the fgrm reported for sires, but lower than the froh reported for sires. comparing the inbreeding coefficients between sires and progeny, there is a clear decrease in inbreeding, which can only be attributed to the randomness of the genetic pool of the dams. according to recent literature, fgrm was determined as the optimal method for estimating genomic inbreeding coefficients as it can distinguish between markers that are ibs from markers that are ibd. the average froh for either sires or progeny estimated in this study are higher than the ones reported for a holstein population which have shown effects of inbreeding depression for traits in dairy cattle. also, inbreeding depression appears not to have strongly affected the phenotypes. this result agrees with the lower inbred estimation for the progeny compared to the sires. we were not able to cross test this hypothesis by estimating inbreeding depression for the traits among sires, as we did not collect their phenotypes, altough we would expect a stronger influence of inbreeding in this case. the average r values obtained in this study, considering snps with distance up to 150 kb apart, are higher than the values reported in the literature for nelore and indicine breeds. beyond 150 kb distance, the average r values are similar to ones found in the literature. however, it is worth mentioning that it is difficult to make direct comparisons among ld studies, as ld levels vary due to many factors. in this study, both maternal and paternal haplotypes were considered in contrast to studies that only consider maternal haplotypes, as it is our understanding that the nelore population presents higher influence from males than females and not considering paternal haplotypes would underestimate r values. we believe that the backcrossing with taurine local breeds, evidenced by the mtdna origin, could have influenced the values of ld and ne estimated for this nelore sample. there are some reports in the literature considering haplotype block structure in cattle, but none was found for nelore cattle. these reports stated that with the increase in marker density, there will be more haplotype blocks of smaller sizes which will respond for higher genome coverage. a recent study using the high density snp panel in crossbred cattle population reports 61 % of the genome covered by haplotype blocks, and a total of 78 % of snps clustered into haplotype blocks. the recent ne found for nelore in this study, is in agreement with reports of recent ne for holstein and swiss eringer breed varying from to individuals, and can go higher than animals for past ne. our results show that the diversity among sire families and lineages is not high enough to correctly separate those in a pca, or hierarchical clustering analysis. this is corroborated by the low fst values obtained for sire families and lineages subpopulations. the uninformative values of ne and fst estimates were expected for the lineages as they do not represent isolated populations and should have been crosses among them. some of the traits used in this work were already analyzed for genome wide associations elsewhere using part of the data and with a different approach. although the phenotype data is the same in both works, our sample size was larger and the software and algorithms used for association studies were different. we have not found larger accordance than 38 %, between the set of associated genes with significance of 1e- to bft, rea and lm traits, studied in that work, although we could not verify this for tcw as this trait was not present in that work. however, the quality and consistency of the genotyping data were found reliable since the pca against other breeds showed consistent results. we also checked for the fluctuation of heterozygosity along the tagsnps and found reasonable smoothed values of heterozygosity fluctuating around to. the traits correlations were also found reasonable. furthermore some of the top tagsnps that had lowest p-values in the gwas, had the corresponding gene with functional biological correlation to the main growth traits, indicating that the gwas was able to extract some of the signal out of the noise from the data. in addition, we also performed a gene enrichment test with gorilla, after ordering the tagsnps for the pleiotropic effect. the pleiotropic test is a methodology to try to rank genes that are most influencing the traits in the awm. the results showed go terms for biological and cell adhesion and developmental growth, which are again evidence that the terms are enriched for genes related mainly to growth which is the key phenotype of the awm. cell adhesion is a biological process that is also found significantly enriched in the work of tizioto et al.,. cell adhesion is a process involved in the binding of a cell to a substrate which sustains its growth and therefore is intrinsically related to the growth process. with a similar objective of the pleiotropy test, the construction of the gene networks with the awm/pcit is an alternative to analyze pure gwas results as this methodology gathers the gwas results and constructs a gene network that influences the traits in a simultaneous way, but weighting accordingly to the correlations within the awm. we decided to select tcw as the key trait for the awm because of its reasonable heritability estimate and its congruent correlations to rea and bft, together with the economic importance of growth traits for livestock. as the resulting network is too large, we used a heuristic to divide the network into a sub-network in a way that it would be feasible to analyze and still gather most of the biological importance regarding these traits. this heuristic is implemented in an algorithm that selects the trio of tf that most span the network and its first degree nodes. as the coloring and shape of the nodes of the resulting sub-network are inherited from the original network, it can be noticed a majority of reddish and larger nodes indicating that the objective of the heuristic was achieved. also, a survey in the literature showed correlations between the trio of tf and several nodes from the sub-network to growth traits, suggesting that our methodology is correct. it is remarkable to note that some genes like vdr was already associated to growth in a bovine gwas and many others have relations to the biology of growth, like proliferation and growth of muscle and bones and adipogenesis. although comparing the results of single snp association to the work of tizioto et al. did not lead to more than 38 % of common associated genes, the genes from the sub-network are in 50 % of accordance with the genes associated to bft and rea from that work. this suggests that the methodology that created this network was able to select genes that are more related to the traits involved with growth. we compared the genes from the sub-network to genes found in other studies that also used the awm/pcit methodology and similar datasets. although there were no great correlations, some genes in these studies were found to have similarities or to encode the same protein domains to some of the genes from the sub-network: metalloprotease and ankyrin domains, the rab6b oncogene, the mer proto-oncogene tyrosine kinase and the slc solute carrier gene family. taken together, these results are suggesting that the gene networks obtained are related to growth and meat quality traits and their genes should be thoroughly inspected to try to discover the biological mechanism beneath growth and meat quality phenotypes.  we used high density snp panels to genotype sires, and its half-sib families, totaling animals. the sires were previously selected in order to represent most of the variability of the nelore beef cattle in brazil. we believe this is the first work that used a genomic approach in order to try to investigate the amount of diversity of the brazilian nelore breed by investigating its ancestral lineages. results showed that the population studied was not structured enough in order to differentiate families using ancestral information of sires or lineages. eight phenotypes related to growth and meat quality were used in whole genome association studies and the results were used to construct gene networks focused on growth, using the methodology of awm/pcit. literature surveys showed that both the gwas and the gene networks constructed had significant snps associated to genes related to growth in former studies like the vdr, lhx cpeb asb and many others. these results suggests that the methodology used to construct the gene networks can be used as an alternative approach to standard gwas, in order to reach for novel information and to try to understand the biological mechanisms and gene compositions that leads to complex phenotypes, like growth in beef cattle. 
the gene networks constructed from the awm/pcit methodology were a useful alternative in characterizing genes and gene networks that were allegedly influential in growth and meat quality traits in nelore cattle. our results indicate a lack of structuring between the individuals studied since principal component analyses were not able to differentiate families by its sires or by its ancestral lineages. an alternative methodology to traditional gwas that involves the construction of gene network interactions, derived from results of several gwas is the awm /pcit. genome-wide association studies are often used to associate markers and genomic regions to growth and meat quality traits that can be used to assist selection programs. the awm/pcit methodology was employed to evaluate the genes that participate in a series of eight phenotypes related to growth and meat quality obtained from this nelore sample. the diversity of the nelore sample studied is not high enough to differentiate among families neither by sires nor by using the available ancestral lineage information. with the aim of evaluating the genetic architecture of brazilian nelore cattle, we used high-density snp genotyping data from nelore animals comprising half-sibling families derived from highly disseminated and unrelated sires from across brazil. the application of the awm/pcit methodology revealed a trio of transcription factors which in combination connected genes through edges and whose biological functions were inspected, some revealing to participate in biological growth processes in literature searches.
 nelore is the major beef cattle breed in brazil with more than million heads. genome-wide association studies are often used to associate markers and genomic regions to growth and meat quality traits that can be used to assist selection programs. an alternative methodology to traditional gwas that involves the construction of gene network interactions, derived from results of several gwas is the awm /pcit. with the aim of evaluating the genetic architecture of brazilian nelore cattle, we used high-density snp genotyping data from nelore animals comprising half-sibling families derived from highly disseminated and unrelated sires from across brazil. the awm/pcit methodology was employed to evaluate the genes that participate in a series of eight phenotypes related to growth and meat quality obtained from this nelore sample.  our results indicate a lack of structuring between the individuals studied since principal component analyses were not able to differentiate families by its sires or by its ancestral lineages. the application of the awm/pcit methodology revealed a trio of transcription factors which in combination connected genes through edges and whose biological functions were inspected, some revealing to participate in biological growth processes in literature searches.  the diversity of the nelore sample studied is not high enough to differentiate among families neither by sires nor by using the available ancestral lineage information. the gene networks constructed from the awm/pcit methodology were a useful alternative in characterizing genes and gene networks that were allegedly influential in growth and meat quality traits in nelore cattle. electronic supplementary material the online version of this article contains supplementary material, which is available to authorized users. keywords genotypingawmpcitgwasissue-copyright-statement© the author 2016 