Next Article in Journal
Characterization of ETFDH and PHGDH Mutations in a Patient with Mild Glutaric Aciduria Type II and Serine Deficiency
Next Article in Special Issue
miR-21-5p Regulates the Proliferation and Differentiation of Skeletal Muscle Satellite Cells by Targeting KLF3 in Chicken
Previous Article in Journal
Pediatric and Adolescent Oncofertility in Male Patients—From Alpha to Omega
Previous Article in Special Issue
Characterization of Four Orphan Receptors (GPR3, GPR6, GPR12 and GPR12L) in Chickens and Ducks and Regulation of GPR12 Expression in Ovarian Granulosa Cells by Progesterone
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Identification and Functional Annotation of Genes Related to Bone Stability in Laying Hens Using Random Forests

1
Institute of Farm Animal Genetics, Friedrich-Loeffler-Institut, 31535 Neustadt, Germany
2
Breeding Informatics Group, Department of Animal Sciences, University of Göttingen, 37075 Göttingen, Germany
3
Center for Integrated Breeding Research (CiBreed), University of Göttingen, 37075 Göttingen, Germany
4
Livestock Center of the Faculty of Veterinary Medicine, Ludwig-Maximilians-University Munich, 85764 Oberschleissheim, Germany
5
Animal Breeding and Genetics Group, Department of Animal Sciences, University of Göttingen, 37075 Göttingen, Germany
*
Author to whom correspondence should be addressed.
Submission received: 19 April 2021 / Revised: 5 May 2021 / Accepted: 6 May 2021 / Published: 8 May 2021
(This article belongs to the Special Issue Poultry Genetics, Breeding and Biotechnology)

Abstract

:
Skeletal disorders, including fractures and osteoporosis, in laying hens cause major welfare and economic problems. Although genetics have been shown to play a key role in bone integrity, little is yet known about the underlying genetic architecture of the traits. This study aimed to identify genes associated with bone breaking strength and bone mineral density of the tibiotarsus and the humerus in laying hens. Potentially informative single nucleotide polymorphisms (SNP) were identified using Random Forests classification. We then searched for genes known to be related to bone stability in close proximity to the SNPs and identified 16 potential candidates. Some of them had human orthologues. Based on our findings, we can support the assumption that multiple genes determine bone strength, with each of them having a rather small effect, as illustrated by our SNP effect estimates. Furthermore, the enrichment analysis showed that some of these candidates are involved in metabolic pathways critical for bone integrity. In conclusion, the identified candidates represent genes that may play a role in the bone integrity of chickens. Although further studies are needed to determine causality, the genes reported here are promising in terms of alleviating bone disorders in laying hens.

1. Introduction

The very high incidence of skeletal disorders in laying hens, including brittle and fractured bones, is undoubtedly one of the most serious problems facing the egg production industry [1,2]. Bone demineralisation associated with eggshell calcification favours the loss of structural bone tissue and ultimately predisposes the birds to osteoporosis in the course of the laying period [3,4]. Besides dramatic effects on animal welfare [5,6,7], bone weakness also has an economic impact [4,8]. According to a widespread assumption, the reduction in bone stability is primarily the result of selection for high laying performance [9,10,11]. However, the role of genetic selection on egg production is now seen in a more differentiated view, with recent studies pointing to factors other than egg number alone [12,13,14].
In the urgently needed improvement of the skeletal health of laying hens, genetics play an important role alongside husbandry and feeding of the birds [3,15,16]. To date, a number of quantitative trait loci (QTL) have been mapped to skeletal traits in chickens [17,18,19,20,21,22]. Dunn et al. [23] discovered a QTL on chromosome 1 that was recently fine-mapped leading to the identification of a promising region around the cystathionine beta synthase gene associated with osteoporosis [24]. The discovery of candidate positions for bone integrity is inevitably linked to technical advances in genotyping and bioinformatics. Today, testing hundreds of thousands of single nucleotide polymorphisms (SNP) by means of genome-wide association studies (GWAS) has become common practice [15,25,26].
Despite its widespread use, GWAS has some potential pitfalls. In addition to population stratification, these include the identification of gene loci with small effect sizes, which rarely reach the statistical significance level due to their low strength of association [27]. At this point, machine learning algorithms represent a promising advance. Several studies have demonstrated their potential in identifying genes with small effect sizes [28,29]. The Random Forests (RF) models in particular seem to have a great potential for analysing a large number of loci simultaneously and identifying corresponding associations [29,30,31]. Recently, this approach has been used to identify genes associated with eggshell strength [27].
The aim of the current study was to identify genomic positions associated with bone stability traits, i.e., breaking strength and mineral density of the tibiotarsus and the humerus, in laying hens. The animal model used comprised four layer lines that differed in their phylogenetic origin (brown-egg vs. white-egg layers) and their egg production level (high vs. moderately performing lines) [32]. Jansen et al. [14] have recently reported promising heritability estimates for bone traits in this set of populations, supporting the assumption of an inherited component of hens’ susceptibility to osteoporosis. In the study reported here, we took a deeper look into the underlying genetic architecture of these hens. This includes the adoption of RF-based feature selection in order to find potentially important SNPs. Subsequently, we performed a series of functional analyses including gene set enrichment analysis. Furthermore, SNP effects were estimated to confirm candidate genes known from the literature to be associated with bone metabolism.

2. Materials and Methods

2.1. Population and Experimental Setup

The population consisted of four purebred chicken layer lines (Gallus gallus domesticus), which are phylogenetically distinct (brown- vs. white-egg lines). Within each of these phylogenetic groups, the two lines differed in terms of egg-laying rate (high- vs. moderate-performing lines) [32,33]. The set of populations was previously subjected to phenotypic analysis and the estimation of genetic parameters [14]. The data set only comprised hens whose total egg number was within the line specific threefold interquartile range and who laid at least one egg from 67 to 69 weeks of age [14]. For the statistical analyses done in this study, we combined the four chicken lines into one set as described below.
For the current research, we used the bone breaking strength (BBS) and bone mineral density (BMD) measurements previously reported by Jansen et al. [14]. A summary of these measurements, taken from reference [14], is presented in Table S1. Briefly, BBS and BMD of the tibiotarsus and humerus were determined by the three-point bending test and dual-energy X-ray absorptiometry, respectively, using dissected bones after the hens were sacrificed at 69 weeks of age.
The experimental setup is shown in Figure 1. We applied the machine learning-based approach of Random Forests to identify genomic positions potentially associated with the given phenotypes. Subsequent functional analyses included gene set enrichment analysis and retrospective SNP effects analysis.

2.2. Genotyping

Initially, deoxyribonucleic acid (DNA) samples from the hens and sires were extracted from blood samples. The hens were genotyped for 51,837 SNPs with a custom-made SNP array (Affymetrix Inc., Santa Clara, CA, USA). From the same chicken lines, in total 80 sires were genotyped for 580,961 SNP markers using the Affymetrix® Axiom® Genome-Wide Chicken Genotyping Array [34]. Quality control was applied to both data sets using the SNP & Variation Suite (SVS) v8.9 [35]. We only considered SNPs from autosomal chromosomes 1 to 28. The genotypes were filtered for a SNP call rate of ≥99% and an animal call rate of ≥95%. Furthermore, missing genotypes were imputed in a two-step procedure using BEAGLE 5.0 [36]. Missing markers within the sire data set were imputed using the default settings. After this, the female genotypes were imputed from 37,606 SNPs left after quality control to 497,041 SNPs. Here, the sire genotypes served as a reference population and the effective population size was set to ne = 5000. After imputation, 524 hens and 497,041 SNPs remained, of which 490,745 SNPs were finally annotated using the genome assembly GRCg6a (galGal6) [37], with duplicated SNPs and those with ambiguous chromosome annotation being removed.

2.3. Random Forests Classification

We applied the machine learning algorithm of Random Forests (RF) to identify SNPs associated with bone characteristics, i.e., BBS and BMD of the tibiotarsus and humerus. Briefly, the RF algorithm constructs a multitude of classifying decision trees assigning importance values to each SNP, thus determining those SNPs that explain variation in the response variable [29]. As shown by Ramzan et al. [27,38], we performed SNP selection by applying the Boruta algorithm, which works as a wrapper around the classification algorithm [39]. This algorithm is based on the idea that an unimportant attribute is not more useful for classification than a random one. Hence, if an attribute shows lower importance than a random attribute, it can be deemed irrelevant. The second idea is that importance measures get more accurate with less irrelevant attributes, such that iteratively removing unimportant attributes increases the accuracy of the importance measure. The procedure of the algorithm is as follows: The dataset is first expanded by adding shuffled copies of the original values of each SNP, called shadow attributes. RF classification is then applied iteratively, assigning a value to each SNP, which is considered as the importance of the SNP. At each iteration, SNPs whose importance is less than the best of their shadow attributes are removed.
We used the Python (v3.8.3) [40] implementation from Homola [41] that specifies the proportion of the shadow attributes by which a SNP has to be better in order to be selected as important. Embedded in the Boruta algorithm, the RF classification itself was carried out with the ‘RandomForestRegressor’ from the Scikit-learn package [42] using default settings. The parameter perc was set to 99, representing a threshold of 99%, as no SNP has been confirmed as important at the 100% level. RF classification was performed separately for each bone trait. The input file consisted of the SNP genotypes, coded as ‘0’ (AA), ‘1’ (AB), or ‘2’ (BB), and the phenotypic values of the respective bone trait. To account for possible confounding effects due to population stratification, residuals representing adjusted phenotypes were analysed instead of the raw values [29]. The following model was used to estimate the residuals:
γ i j k l = μ + G i + L L j + S k + ε i j k l
where γ i j k l is the observation for a bone trait, μ is the general mean, G i is the fixed effect of generation ( i = 1, 2), L L j is the fixed effect of layer line ( j = 1 to 4), S k is the random effect of sire ( k = 1 to 145), and ε i j k l is the residual error. The model was computed using JMP v14.0 (SAS Institute Inc. Cary, NC, USA, 2018). Normal distribution of the residuals was assumed (Figure S1).
The output of the RF classification was a list of confirmed SNPs, i.e., markers that are more than coincidentally associated with a given bone trait.

2.4. Functional Analyses

2.4.1. Gene Extraction

All steps of the functional analyses were carried out using R v4.0.3 [43]. Extraction of genes associated with SNPs identified by the RF classification from the Ensembl database v102 [44] was performed using BioMart [45]. All protein-coding genes that are located within 5 kb upstream and downstream of the respective SNPs were considered for the gene lists. Information on the biological functions of these genes was obtained from both the NCBI [46] and Ensembl databases, as well as from the literature. The gene lists were then screened for genes known to be associated with bone stability traits. In this way, we identified a number of genes that were henceforth regarded as candidate genes.

2.4.2. SNP Effects Analysis

The genotypic effect was analysed for those SNPs located in intragenic or in flanking genomic regions of candidate genes, which have previously been shown to be significantly associated with a bone trait (see Table 1). SNP effects for each locus were analysed as described by Wiedemann et al. [47]. For this purpose, the actual SNP genotypes were coded as ‘0’ (AA), ‘1’ (AB), or ‘2’ (BB), with the B allele representing the minor allele. The minor allele was considered the effect allele, whereas the major allele was termed ‘other allele’. All models were computed with the R package lme4 [48].
A linear regression model adjusted for fixed factors was applied to estimate the allele substitution effects by single marker regression (SMR):
γ i j k l m = μ + G i + L L j + b 1 S N P k + S l + ε i j k l m
where γ i j k l m is the observation for a bone trait, μ is the overall mean effect, G i is the fixed effect of generation ( i = 1, 2), L L j is the fixed effect of layer line ( j = 1 to 4),   b 1 is the regression coefficient of the SNP genotype ( S N P k ), S l is the random effect of sire ( l = 1 to 145), and ε i j k l m is the residual error. Standardised allele substitution effects were calculated according to model (2) after both the dependent variable and the SNP genotypes coded as ‘0’, ‘1’, or ‘2’ were standardised to have a mean of 0 and a standard deviation of 1.
To calculate the additive and dominance effects, a dominant-recessive model (DRM) was applied considering the SNP genotype as a fixed class variable. The statistical model was as follows:
γ i j k l m = μ + G i + L L j + S N P k + S l + ε i j k l m
where γ i j k l m is the observation for a bone trait, μ is the overall mean effect, G i is the fixed effect of generation ( i = 1, 2), L L j is the fixed effect of layer line ( j = 1 to 4), S N P k is the fixed effect of SNP genotype ( k = 1 to 3), S l is the random effect of sire ( l = 1 to 145), and ε i j k l m is the residual error. Least squares means (LSM) for the different genotypes were estimated with the emmeans package [49]. Significant differences between LSM were tested using a t-test and adjusted by the Bonferroni method. Additive and dominance effects were estimated by contrasting the respective genotypes according to the following formulas.
a = μ A A μ B B 2
d = μ A B μ A A + μ B B 2
where a is the additive effect, d is the dominance effect, μ A A and μ B B are the phenotypic mean values of the homozygous genotypes, and μ A B is the phenotypic mean value of the heterozygous genotype.

2.4.3. Gene Set Analysis

With the gene sets including all genes extracted, we performed gene set analysis (GSA) using g:Profiler2 [50]. This involved the Gene Ontology (GO) (Ensembl v102) and the Kyoto Encyclopedia of Genes and Genomes (KEGG) [51] (FTP release 2020-09-07) databases. The GSA was carried out considering all known genes obtained from Ensembl for the calculation of statistical significance and applying the default g: SCS algorithm [52] for computing the multiple testing correction. Only GO- and pathway terms with significant enrichment (p < 0.05) were considered for further analyses. Tree maps of the GO terms were generated using rrvgo [53].
Table 1. General information for 17 loci associated with the bone breaking strengths (BBS) or bone mineral densities (BMD) of the tibiotarsus (Tib) and humerus (Hum) selected for the SNP effects analysis.
Table 1. General information for 17 loci associated with the bone breaking strengths (BBS) or bone mineral densities (BMD) of the tibiotarsus (Tib) and humerus (Hum) selected for the SNP effects analysis.
SNPTraitLocationGGA 1Position 2GenotypesN IndividualsGenotype FrequenciesEA/OA 3EA FrequencyCandidate GeneReference 4
AX-75268181Tib_BMDintragenic1139,001,157CC/CT/TT392/96/360.75/0.18/0.07T/C0.16MCF2L[54]
AX-76044166Tib_BBSintragenic215,440,861AA/AG/GG421/63/400.80/0.12/0.08G/A0.14MPP7[55]
AX-80813610Tib_BMDdownstream223,056,581CC/CG/GG339/113/720.65/0.22/0.13G/C0.25CALCR[56]
AX-76099065Tib_BMDintragenic246,101,680GG/GA/AA392/77/550.75/0.15/0.10A/G0.18SFRP4[57]
AX-76601713Tib_BBSintragenic310,617,925AA/AG/GG265/102/1570.51/0.19/0.30G/A0.40ACTR2[15]
AX-77276717Tib_BBSintragenic319,498,104GG/GA/AA322/145/570.61/0.28/0.11A/G0.25TGFB2[58]
AX-76491534Tib_BBSintragenic349,027,160AA/AG/GG432/62/300.82/0.12/0.06G/A0.12CCDC170[59]
AX-76772658Tib_BBS/Hum_BBSintragenic511,438,677TT/TC/CC219/199/1090.41/0.38/0.21C/T0.40SOX6[60]
AX-77113061Tib_BMDupstream85,889,886GG/AG/AA202/156/1660.38/0.30/0.32A/G0.47TMCO1[61]
AX-77091655Hum_BBS/Hum_BMDupstream824,931,025CC/CA/AA286/139/990.54/0.27/0.19A/C0.32PODN[15]
AX-75597497Hum_BBSdownstream1019,108,829AA/AG/GG376/124/240.72/0.24/0.04G/A0.16SMAD6[62]
AX-75677174Tib_BMDintragenic1110,044,055CC/CT/TT377/107/400.72/020/0.08T/C0.18GPATCH1[55]
AX-75711229Tib_BBSintragenic123,804,145GG/AG/AA459/58/70.88/0.11/0.01A/G0.07ASPN[63]
AX-75913642Tib_BBSupstream188,793,585GG/AG/AA451/61/120.86/0.12/0.02A/G0.08SOX9[64]
AX-76351785Hum_BBSintragenic273,497,444CC/CT/TT316/138/700.61/0.26/0.13T/C0.26WNT9B[65]
AX-76351898Hum_BMDdownstream273,518,924GG/GA/AA483/31/100.92/0.06/0.02A/G0.05WNT3[55]
AX-76351899Hum_BMDdownstream273,519,091TT/TC/CC483/31/100.92/0.06/0.02C/T0.05WNT3[55]
1 GGA, Gallus gallus chromosome; 2 Physical position (bp) according to the GRCg6a (galGal6) genome assembly; 3 EA, effect allele (minor allele); OA, other allele (major allele); 4 References from the literature suggesting an association of the gene with bone stability traits.

3. Results

3.1. Identified Single Nucleotide Polymorphisms

Lists of confirmed SNPs were obtained from the RF classifier for each of the phenotypic traits. For the tibiotarsus, 358 (BBS) and 374 (BMD) SNPs were confirmed as important, whereas for the humerus 188 (BBS) and 178 (BMD) markers were identified, respectively. There were no confirmed SNPs on GGA (Gallus gallus chromosome) 16 for any of the four traits studied (Figure S2). In the case of the tibiotarsus, the majority of SNPs were located on GGA 1. In general, there were fewer markers for the humerus, with no markers found on GGA 28. Comparing the two bone types, more than twice as many SNPs were identified for the tibiotarsus.

3.2. Functional Analyses

3.2.1. Extracted Gene Sets

We identified 240 (BBS) and 220 (BMD) genes within an interval of 5 kb upstream and downstream of SNPs that were found to be significant for the tibiotarsus. In contrast, gene sets for the humerus included 115 (BBS) and 113 (BMD) genes.
A Venn diagram was drawn to find overlaps and differences between the genes identified for the BBS and BMD of the two bone types (Figure 2). The corresponding detailed gene list is given in Table S2. It was found that the overlaps of loci between the individual traits were rather small. It ranged from 1.7% (six genes) between BBS of tibiotarsus and humerus up to an overlap of 6.7% (31 genes) between BBS and BMD within the tibiotarsus. No gene was found in all bone and trait combinations. Rather, they were mainly unique genes.
Based on the information on their biological functions and from the literature review, we reduced the gene lists to genes that are known to be related to bone stability. We found 16 genes with an already described association (Table 1). These genes are located on GGA 1 (MCF2L), GGA 2 (MPP7, CALCR, and SFRP4), GGA 3 (ACTR2, TGFB2, and CCDC170), GGA 5 (SOX6), GGA 8 (TMCO1, PODN), GGA 10 (SMAD6), GGA 11 (GPATCH1), GGA 12 (ASPN), GGA 18 (SOX9), and GGA 27 (WNT9B, WNT3).

3.2.2. SNP Effects Analysis

To reveal the biological significance of the candidate genes, we analysed their associations with the corresponding phenotypic bone traits. To this end, we performed SNP effects analyses of all markers detected by the RF classifier and then assigned to genes (Table 1). Since the SNPs AX-77091655 (PODN) and AX-76772658 (SOX6) were associated with two traits each and, in addition, two further markers were assigned to the WNT3 gene, SNP effects were estimated for 19 SNP and bone trait combinations.
Results from the SMR model are shown in Table 2. Analysis of variance revealed significant effects of SNP genotypes on the respective bone traits. Only the SNPs AX-77276717 (TGFB2) and AX-75711229 (ASPN) had no significant effect. Locus AX-76099065 (SFRP4) had the greatest effect on tibiotarsus BMD, with the substitution of allele G for allele A leading to a reduction of 0.016 g/cm2. In contrast, increasing the number of the copies of the effect alleles at loci AX-76351898 and AX-76351899, both assigned to the WNT3 gene, would yield an increase in humerus BMD of 0.016 g/cm2. Of all loci significantly associated with BBS, the SNP AX-76491534 (CCDC170) showed the largest effect, whereby substitution of one copy of allele A with allele G would result in a 15.63 N decrease of tibiotarsus BBS. The counterpart is the SNP AX-76351785 (WNT9B) in which the T allele would presumably cause an increase of 11.51 N of humerus BBS. This is also the largest effect among all significant SNPs, with a change of 0.21 expressed in SD units.
Results obtained from the DRM are shown in Table 3. Comparison of the genotypic values (LSM) revealed significant differences among the genotypes. This applies to all loci studied, with exception of the SNPs AX-77276717 (TGFB2) and AX-80813610 (CALCR), where only a tendency towards a higher value for the homozygote genotype of the effect allele was observed. For the other loci, the effects indicated by the direction of the beta coefficients (SMR) were also reflected in the genotypic values. A significantly higher LSM was found for the homozygote genotype of the effect allele of the SNP AX-75711229 (ASPN). However, this estimate might be biased as the corresponding genotype had a frequency of only 0.01 (Table 1) and no significant allele substitution effect was detected for this locus (Table 2).
Significant additive effects of the respective other allele (major allele) were accounted for all loci with exception of the SNPs AX-77276717 (TGFB2) and AX-76772658 (SOX6) (Table 3). The estimates ranged from −0.02 to 0.017 g/cm2 for the BMD-related SNPs and from −16.70 to 15.70 N for the markers associated with the BBS. Effects of complete dominance were observed for the SNPs AX-76044166 (MPP7), AX-75711229 (ASPN), AX-75597497 (SMAD6) and AX-76099065 (SFRP4), with one copy of the major allele masking the recessive allele, thus leading to full trait expression. In contrast, complete dominance in favour of the effect allele was seen for the SNP AX-77113061 (TMCO1).

3.2.3. Gene Set Analysis

GSA was performed considering the total gene sets. We restricted the results presented to the GO biological process (BP) category, as we sought to determine overarching biological objectives to which the gene products of the extracted genes contribute. Furthermore, the genes were grouped according to their KEGG pathways. Full lists of significantly enriched GO terms, including those from the cellular component and molecular function categories, are given in Table S3.
A large number of genes were involved in common processes. The analysis reported 81 (BBS) and 51 (BMD) significantly enriched BPs for the tibiotarsus and 33 (BBS) and 42 (BMD) BPs for the humerus, respectively (Table S3). Of these, Figure 3 (tibiotarsus) and Figure 4 (humerus) show the top 15 significantly enriched GO BP terms with the highest −log10 p-values and all significantly enriched KEGG pathways obtained from the RF classifier. Although certain BPs overlapped between the bone and trait combinations, no relation to the skeletal system was evident in the enriched BPs. Visualizing the results using tree maps to investigate redundancy based on semantic similarity between different GO terms also did not yield any biologically relevant findings (Figures S3 and S4). However, the literature points to the involvement of Wnt- and MAPK signalling pathways in the pathogenesis of osteoporosis [25]. GSA revealed the Wnt signaling pathway (KEGG:04310) to be significantly enriched in both BMD gene sets (Figure 3B and Figure 4B). In addition, significant enrichment for the MAPK signaling pathway (KEGG:04010) was identified in the genes for BMD of the tibiotarsus.

4. Discussion

The objective of the present study was to identify genomic positions potentially associated with skeletal integrity in a laying hen population. There is solid evidence that osteoporosis is a polygenic disorder, i.e., determined by multiple functional genes acting conjointly rather than a few major genes [15,25]. For this reason, we applied RF classification, an approach known to be able to detect genes with modest effects [29,30]. To our knowledge, this is the first study applying a machine-learning approach to bone data in chickens. Using RF classification, we identified a high number of potentially informative SNPs. Although a large number of genes were adjacent to these SNPs, only 16 candidate genes related to skeletal disorders were identified; of these, many had human orthologues. However, for the vast majority of genes, no involvement in bone metabolism has been suspected so far, which is in line with previous reports [15,66]. From the 16 identified candidates discussed below, we first focus on genes that have previously been linked to BBS or BMD (n = 10), followed by genes for which an association with osteoarthritis is suggested (n = 3). Finally, genes are discussed that are functionally related to the Wnt signalling pathway (n = 3).
Ten of our candidate genes can be grouped as having previously been associated with BBS or BMD traits in the literature. Of these, the membrane palmitoylated protein 7 gene (MPP7) was associated with vertebral BMD in humans [55]. Its strong functional role in osteoblast biology was demonstrated by means of in vivo and in vitro studies [67]. Based on these reports, we consider MPP7 to be a good candidate for bone disorders in chickens. In our study, the calcitonin receptor gene (CALCR) was identified as a strong candidate for BMD. Calcitonin plays a role in calcium homeostasis and is primarily an inhibitor of bone resorption [68]. Our observations are in line with previous reports, as CALCR polymorphisms were associated with site-specific BMD in humans [56,69], and alpha-calcitonin gene-related peptide deficient mice were shown to have a lower bone mass [70]. One of the candidates for BBS located on GGA 3 is the actin related protein 2 gene (ACTR2), which was recently identified by Raymond et al. [15] as being associated with BBS in laying hens. ACTR2 is functionally linked to bone via its importance for cilia formation, as cilia are known to play an integral role in skeletal development [15,71]. Although no significant effect of the variant corresponding to the transforming growth factor beta 2 gene (TGFB2) was observed in our study, TGFB2 is considered a very promising candidate for skeletal integrity in the chicken. As a cytokine, the protein encoded by TGFB2 has important functions in many biological processes related to bone remodelling [19,58]. Analyses in different chicken populations including broilers and layers suggest TGFB2 to be associated with various bone characteristics [19,21,58]. In this context, the SMAD family member 6 gene (SMAD6) has to be mentioned, which we identified as a candidate for BBS. Its protein acts as a regulator of the TGF-beta family and inhibits bone morphogenetic protein pathways, which are integral parts of osteoblast and chondrocyte differentiation [72,73]. A study on mice revealed their essential role in bone formation, as SMAD3 knockout resulted in osteopenia [62]. The coiled-coil domain containing 170 gene (CCDC170) is our third candidate for BBS located on GGA 3. The region around this locus has been linked to BMD in humans [59,65]. However, since the function of the protein is yet unclear, it has been speculated whether associations attributed to CCDC170 do not rather belong to the adjacent estrogen receptor 1 gene [59]. In a follow-up study, CCDC170 polymorphisms were in turn associated with osteoporosis-relevant phenotypes [74]. Only one of our candidates was located on GGA 5. The corresponding variant is located in the intron of the SRY-box 6 gene (SOX6), which encodes a transcription factor known to affect developmental processes and skeletal formation in humans [60,65]. In addition, the gene was linked with BMD of the femoral neck [75], and skeletal abnormalities have previously been observed in SOX6 knockout mice, suggesting an integral role in cartilage formation [76]. We identified the transmembrane and coiled-coil domains 1 gene (TMCO1), located on GGA 8, as a candidate for BMD. TMCO1 plays an important role in bone formation-mediating calcium homeostasis within the endoplasmic reticulum [61]. Disruption of the endoplasmic reticulum of an osteoblast can lead to severe bone disorders [77]. Recently, Li et al. [61] demonstrated that TMCO1 deficiency leads to reduced bone formation and osteoblast differentiation in humans and mice. In addition to SOX6, the podocan gene (PODN) is another candidate that was associated with two traits, namely BBS and BMD of the humerus. PODN encodes a proteoglycan that was shown to bind type 1 collagen, suggesting a potential role in growth regulation [78]. At this point, the great influence of collagen on mechanical properties of bones should be mentioned, which is assumed to apply equally to humans [79] and chickens [80]. That PODN could be a promising candidate for bone integrity in laying hens is supported by findings of Raymond et al. [15]. Although the G-patch domain containing 1 gene (GPATCH1), identified as BMD candidate, is considered a candidate gene for osteoporosis in humans [55], functional information is limited and its role in skeletal pathophysiology is not yet clear.
For a group of three candidates, the literature suggests a functional relationship with osteoarthritis, a pathological condition of cartilage degradation [81]. Osteoarthritis and osteoporosis are closely related and characterised by subchondral bone loss and excessive bone resorption [20,81,82]. It is assumed that both diseases are partly determined by common genes [83]. One of the candidates found in our study is the MCF.2 cell line derived transforming sequence like gene (MCF2L), shown to be expressed in cartilage tissue, and linked to joint osteoarthritis in humans [54,84]. In addition, Mao et al. [85] recently pointed out the relevance of MCF2L for osteoporosis, which underlines the link between both disorders. The asporin gene (ASPN), also known as biglycan (BGN), is assumed to regulate chondrogenesis. While the results of Mishra et al. [63] point to a functional role of ASPN in osteoarthritis, other studies reported only a marginal relationship or contradict such an association [86,87]. Given these contradictory results and the fact that the association with ASPN was not significant in our study, we consider ASPN a suggestive candidate that requires further investigation. The SRY-box 9 gene (SOX9) is our third candidate linked to osteoarthritis [88]. SOX9 is considered a pivotal player in chondrogenesis, as its protein, the transcription factor SOX-9, was shown to stimulate chondrocyte differentiation [64,89]. In addition, SOX9 mediates the Wnt signalling pathway, abnormalities of which are correlated with cartilage degradation [64].
The remaining candidates, i.e., the SFRP4, WNT3, and WNT9B genes, are functionally linked to the Wnt signalling pathway, which plays a key role in various basic developmental processes [90]. The secreted frizzled related protein 4 gene (SFRP4) encodes a protein that primarily antagonizes Wnt polypeptides [90] and is one of the BMD candidates. A mutation in SFRP4 was shown to cause pathological reduction of cortical bone tissue in mice and humans [57]. The Wnt signalling pathway is crucial for bone metabolism and to date, several Wnt genes are known to be associated with traits such as bone mass and BMD [55,91]. This also includes the Wnt family member 3 gene (WNT3), which was identified in this study [55,88]. The Wnt family member 9B gene (WNT9B), located adjacent to WNT3 on GGA 27, was identified as a candidate for BBS. Although its role in skeletal biology is less explored than that of other Wnt genes, we consider WNT9B a susceptibility gene for bone strength due to its association with femur BMD [65]. The high importance of the Wnt signalling pathway for bone strength is supported by the significant enrichment that was shown in the GSA for this functional pathway. Furthermore, the mitogen-activated protein kinase (MAPK) signalling pathway was enriched, which is also very important for skeletal development and, in particular, for chondrogenesis [92]. These observations are in accordance with recent results from pathway analyses [25,93].
Taken together, we identified a number of genetic loci associated with the bone traits studied. Based on these findings, we can confirm the assumption that bone stability is determined by multiple genes, each of which has a rather small effect size. The genes presented here represent suggestive susceptibility genes of bone integrity in chickens, some of which are nonetheless very promising based on what is known so far. Follow-up studies will be required to determine causalities and further uncover the biological significance of these genes. Here, the use of an F2 mapping population for high-resolution mapping of loci is recommended [94]. Considering the animal model, a follow-up study should also investigate the influence of phylogenetic origin on bone phenotypic plasticity, which was not done here, as we focused on finding loci that are significant for laying hens across phylogenetically divergent layer lines.

5. Conclusions

In this study, RF classification was performed to identify loci related with bone integrity in laying hens. In the subsequent functional analyses, a set of 16 promising candidate genes was identified, although in some cases rather small SNP effect estimates were observed. Some of the genes were shown to be involved in pivotal pathways that regulate bone metabolism. Our results strongly support genetics as a crucial factor that contributes significantly to the regulation of bone strength and thus offers great opportunities to improve bone health in laying hens. Further functional analyses on the candidate genes identified at a suggestive level have to follow in order to confirm their biological significance.

Supplementary Materials

The following are available online at https://0-www-mdpi-com.brum.beds.ac.uk/article/10.3390/genes12050702/s1, Figure S1. Histograms of the residuals, Figure S2: Number of annotated SNPs per chromosome, Figure S3: Tree maps of significantly enriched Gene Ontology terms for genes associated with bone breaking strength, Figure S4: Tree maps of significantly enriched Gene Ontology terms for genes associated with bone mineral density, Table S1: Summary of the phenotypic measurements, Table S2: List of all extracted genes, Table S3: List of all enriched terms/pathways.

Author Contributions

Conceptualisation, S.J., U.B. and S.W; data curation, A.W.; formal analysis, S.J., U.B. and J.S.; funding acquisition, S.W.; investigation, S.J., C.H. and A.W.; resources, A.M.S.; supervision, A.R.S. and S.W.; visualisation, S.J.; writing—original draft preparation, S.J.; writing—review and editing, S.J., U.B., C.H., F.R., J.S., A.O.S., A.M.S., A.R.S., A.W. and S.W. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the AdaptHuhn Initiative of the Friedrich-Loeffler-Institut.

Institutional Review Board Statement

The study was conducted according to the guidelines of the Declaration of Helsinki and approved by the Lower Saxony State Office for Consumer Protection and Food Safety (LAVES) (33.19-42502-04-15/1988).

Informed Consent Statement

Not applicable.

Data Availability Statement

The data presented in this study are available on request from the corresponding author.

Acknowledgments

We like to thank Torsten Pook (Animal Breeding and Genetics Group, University of Göttingen) for performing the imputation of the chicken genotypes. Furthermore, we greatly appreciate the support in data collection and DNA preparation provided by the technical staff of the Institute of Farm Animal Genetics.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Sandilands, V. The laying hen and bone fractures. Vet. Rec. 2011, 169, 411. [Google Scholar] [CrossRef]
  2. Petrik, M.T.; Guerin, M.T.; Widowski, T.M. On-farm comparison of keel fracture prevalence and other welfare indicators in conventional cage and floor-housed laying hens in Ontario, Canada. Poult. Sci. 2015, 94, 579–585. [Google Scholar] [CrossRef]
  3. Fleming, R.H.; McCormack, H.A.; McTeir, L.; Whitehead, C.C. Relationships between genetic, environmental and nutritional factors influencing osteoporosis in laying hens. Br. Poult. Sci. 2006, 47, 742–755. [Google Scholar] [CrossRef]
  4. Cransberg, P.H.; Parkinson, G.B.; Wilson, S.; Thorp, B.H. Sequential studies of skeletal calcium reserves and structural bone volume in a commercial layer flock. Br. Poult. Sci. 2001, 42, 260–265. [Google Scholar] [CrossRef]
  5. McCoy, M.A.; Reilly, G.A.C.; Kilpatrick, D.J. Density and breaking strength of bones of mortalities among caged layers. Res. Vet. Sci. 1996, 60, 185–186. [Google Scholar] [CrossRef]
  6. Nasr, M.A.F.; Nicol, C.J.; Murrell, J.C. Do laying hens with keel bone fractures experience pain? PLoS ONE 2012, 7, e42420. [Google Scholar] [CrossRef]
  7. Webster, A.B. Welfare implications of avian osteoporosis. Poult. Sci. 2004, 83, 184–192. [Google Scholar] [CrossRef]
  8. Nasr, M.A.F.; Murrell, J.; Nicol, C.J. The effect of keel fractures on egg production, feed and water consumption in individual laying hens. Br. Poult. Sci. 2013, 54, 165–170. [Google Scholar] [CrossRef]
  9. Eusemann, B.K.; Patt, A.; Schrader, L.; Weigend, S.; Thöne-Reineke, C.; Petow, S. The role of egg production in the etiology of keel bone damage in laying hens. Front. Vet. Sci. 2020, 7, 81. [Google Scholar] [CrossRef] [Green Version]
  10. Habig, C.; Baulain, U.; Henning, M.; Scholz, A.; Sharifi, A.; Janisch, S.; Simianer, H.; Weigend, S. How bone stability in laying hens is affected by phylogenetic background and performance level. Eur. Poult. Sci. 2017, 81, 1–17. [Google Scholar]
  11. Hocking, P.M.; Bain, M.; Channing, C.E.; Fleming, R.; Wilson, S. Genetic variation for egg production, egg quality and bone strength in selected and traditional breeds of laying fowl. Br. Poult. Sci. 2003, 44, 365–373. [Google Scholar] [CrossRef]
  12. Gebhardt-Henrich, S.G.; Fröhlich, E.K. Early onset of laying and bumblefoot favor keel bone fractures. Animals 2015, 5, 1192–1206. [Google Scholar] [CrossRef] [Green Version]
  13. Dunn, I.C.; De Koning, D.-J.; McCormack, H.A.; Fleming, R.H.; Wilson, P.W.; Andersson, B.; Schmutz, M.; Benavides, C.; Dominguez-Gasca, N.; Sanchez-Rodriguez, E.; et al. No evidence that selection for egg production persistency causes loss of bone quality in laying hens. Genet. Sel. Evol. 2021, 53, 11. [Google Scholar] [CrossRef]
  14. Jansen, S.; Baulain, U.; Habig, C.; Weigend, A.; Halle, I.; Scholz, A.M.; Simianer, H.; Sharifi, A.R.; Weigend, S. Relationship between Bone Stability and Egg Production in Genetically Divergent Chicken Layer Lines. Animals 2020, 10, 850. [Google Scholar] [CrossRef]
  15. Raymond, B.; Johansson, A.M.; McCormack, H.A.; Fleming, R.H.; Schmutz, M.; Dunn, I.C.; De Koning, D.J. Genome-wide association study for bone strength in laying hens. J. Anim. Sci. 2018, 96, 2525–2535. [Google Scholar] [CrossRef]
  16. Bishop, S.C.; Fleming, R.H.; McCormack, H.A.; Flock, D.K.; Whitehead, C.C. Inheritance of bone characteristics affecting osteoporosis in laying hens. Br. Poult. Sci. 2000, 41, 33–40. [Google Scholar] [CrossRef]
  17. Rubin, C.J.; Brändström, H.; Wright, D.; Kerje, S.; Gunnarsson, U.; Schutz, K.; Fredriksson, R.; Jensen, P.; Andersson, L.; Ohlsson, C.; et al. Quantitative trait loci for BMD and bone strength in an intercross between domestic and wildtype chickens. J. Bone Miner. Res. 2007, 22, 375–384. [Google Scholar] [CrossRef] [Green Version]
  18. Zhou, H.; Deeb, N.; Evock-Clover, C.M.; Mitchell, A.D.; Ashwell, C.M.; Lamont, S.J. Genome-Wide Linkage Analysis to Identify Chromosomal Regions Affecting Phenotypic Traits in the Chicken. III. Skeletal Integrity. Poult. Sci. 2007, 86, 255–266. [Google Scholar] [CrossRef]
  19. Bennett, A.K.; Hester, P.Y.; Spurlock, D.M. Relationships of a transforming growth factor-beta2 single nucleotide polymorphism and messenger ribonucleic acid abundance with bone and production traits in chickens. Poult. Sci. 2007, 86, 829–834. [Google Scholar] [CrossRef]
  20. Johnsson, M.; Gustafson, I.; Rubin, C.-J.; Sahlqvist, A.-S.; Jonsson, K.B.; Kerje, S.; Ekwall, O.; Kämpe, O.; Andersson, L.; Jensen, P.; et al. A Sexual Ornament in Chickens Is Affected by Pleiotropic Alleles at HAO1 and BMP2, Selected during Domestication. PLoS Genetics 2012, 8, e1002914. [Google Scholar] [CrossRef] [Green Version]
  21. Schreiweis, M.A.; Hester, P.Y.; Moody, D.E. Identification of quantitative trait loci associated with bone traits and body weight in an F2 resource population of chickens. Genet. Sel. Evol. 2005, 37, 677. [Google Scholar] [CrossRef] [Green Version]
  22. Podisi, B.K.; Knott, S.A.; Dunn, I.C.; Burt, D.W.; Hocking, P.M. Bone mineral density QTL at sexual maturity and end of lay. Br. Poult. Sci. 2012, 53, 763–769. [Google Scholar] [CrossRef] [Green Version]
  23. Dunn, I.C.; Fleming, R.H.; McCormack, H.A.; Morrice, D.; Burt, D.W.; Preisinger, R.; Whitehead, C.C. A QTL for osteoporosis detected in an F2 population derived from White Leghorn chicken lines divergently selected for bone index. Anim. Genet. 2007, 38, 45–49. [Google Scholar] [CrossRef] [Green Version]
  24. De Koning, D.-J.; Dominguez-Gasca, N.; Fleming, R.H.; Gill, A.; Kurian, D.; Law, A.; McCormack, H.A.; Morrice, D.; Sanchez-Rodriguez, E.; Rodriguez-Navarro, A.B.; et al. An eQTL in the cystathionine beta synthase gene is linked to osteoporosis in laying hens. Genet. Sel. Evol. 2020, 52, 13. [Google Scholar] [CrossRef] [Green Version]
  25. Guo, J.; Sun, C.; Qu, L.; Shen, M.; Dou, T.; Ma, M.; Wang, K.; Yang, N. Genetic architecture of bone quality variation in layer chickens revealed by a genome-wide association study. Sci. Rep. 2017, 7, 45317. [Google Scholar] [CrossRef] [Green Version]
  26. Rothammer, S.; Bernau, M.; Kremer-Rücker, P.V.; Medugorac, I.; Scholz, A.M. Genome-wide QTL mapping results for regional DXA body composition and bone mineral density traits in pigs. Arch. Anim. Breed. 2017, 60, 51–59. [Google Scholar] [CrossRef]
  27. Ramzan, F.; Klees, S.; Schmitt, A.O.; Cavero, D.; Gültas, M. Identification of Age-Specific and Common Key Regulatory Mechanisms Governing Eggshell Strength in Chicken Using Random Forests. Genes 2020, 11, 464. [Google Scholar] [CrossRef]
  28. Romagnoni, A.; Jégou, S.; Van Steen, K.; Wainrib, G.; Hugot, J.-P.; Peyrin-Biroulet, L.; Chamaillard, M.; Colombel, J.-F.; Cottone, M.; D’Amato, M.; et al. Comparative performances of machine learning methods for classifying Crohn Disease patients using genome-wide genotyping data. Sci. Rep. 2019, 9, 10351. [Google Scholar] [CrossRef]
  29. Brieuc, M.S.O.; Waters, C.D.; Drinan, D.P.; Naish, K.A. A practical introduction to Random Forest for genetic association studies in ecology and evolution. Mol. Ecol. Resour. 2018, 18, 755–766. [Google Scholar] [CrossRef]
  30. Brieuc, M.S.O.; Ono, K.; Drinan, D.P.; Naish, K.A. Integration of Random Forest with population-based outlier analyses provides insight on the genomic basis and evolution of run timing in Chinook salmon (Oncorhynchus tshawytscha). Mol. Ecol. 2015, 24, 2729–2746. [Google Scholar] [CrossRef]
  31. Nguyen, T.-T.; Huang, J.Z.; Wu, Q.; Nguyen, T.T.; Li, M.J. Genome-wide association data classification and SNPs selection using two-stage quality-based Random Forests. BMC Genom. 2015, 16, S5. [Google Scholar] [CrossRef] [Green Version]
  32. Lieboldt, M.-A.; Halle, I.; Frahm, J.; Schrader, L.; Baulain, U.; Henning, M.; Preisinger, R.; Dänicke, S.; Weigend, S. Phylogenic versus selection effects on growth development, egg laying and egg quality in purebred laying hens. Eur. Poult. Sci. 2015, 79, 1–16. [Google Scholar]
  33. Dudde, A.; Weigend, S.; Krause, E.T.; Jansen, S.; Habig, C.; Schrader, L. Chickens in motion: Effects of egg production level and pen size on the motor abilities and bone stability of laying hens (Gallus gallus forma domestica). Appl. Anim. Behav. Sci. 2020, 227, 104998. [Google Scholar] [CrossRef]
  34. Kranis, A.; Gheyas, A.A.; Boschiero, C.; Turner, F.; Yu, L.; Smith, S.; Talbot, R.; Pirani, A.; Brew, F.; Kaiser, P.; et al. Development of a high density 600K SNP genotyping array for chicken. BMC Genom. 2013, 14, 59. [Google Scholar] [CrossRef] [Green Version]
  35. SNP & Variation Suite ™ (Version 8.9); Golden, Helix, Inc.: Bozeman, MT, USA; Available online: https://www.goldenhelix.com (accessed on 6 May 2021).
  36. Browning, B.L.; Zhou, Y.; Browning, S.R. A One-Penny Imputed Genome from Next-Generation Reference Panels. Am. J. Hum. Genet. 2018, 103, 338–348. [Google Scholar] [CrossRef] [Green Version]
  37. National Center for Biotechnology Information (NCBI). Genome Assembly GRCg6a. Available online: https://0-www-ncbi-nlm-nih-gov.brum.beds.ac.uk/assembly/GCF_000002315.6 (accessed on 17 June 2020).
  38. Ramzan, F.; Gültas, M.; Bertram, H.; Cavero, D.; Schmitt, A.O. Combining Random Forests and a Signal Detection Method Leads to the Robust Detection of Genotype-Phenotype Associations. Genes 2020, 11, 892. [Google Scholar] [CrossRef]
  39. Kursa, M.B.; Rudnicki, W.R. Feature selection with the Boruta package. J. Stat. Softw. 2010, 36, 1–13. [Google Scholar] [CrossRef] [Green Version]
  40. Python Software Foundation. Python Language Reference, Version 3.8.9. Available online: https://docs.python.org/3.8/reference/index.html (accessed on 6 May 2021).
  41. Homola, D. Python Implementations of the Boruta All-Relevant Feature Selection Method. Available online: https://github.com/scikit-learn-contrib/boruta_py (accessed on 14 October 2020).
  42. Pedregosa, F.; Varoquaux, G.; Gramfort, A.; Michel, V.; Thirion, B.; Grisel, O.; Blondel, M.; Prettenhofer, P.; Weiss, R.; Dubourg, V. Scikit-learn: Machine learning in Python. J. Mach. Learn. Res. 2011, 12, 2825–2830. [Google Scholar]
  43. R Core Team. R: A Language and Environment for Statistical Computing; R Foundation for Statistical Computing: Vienna, Austria, 2020. [Google Scholar]
  44. Yates, A.D.; Achuthan, P.; Akanni, W.; Allen, J.; Allen, J.; Alvarez-Jarreta, J.; Amode, M.R.; Armean, I.M.; Azov, A.G.; Bennett, R.; et al. Ensembl 2020. Nucleic Acids Res. 2020, 48, 682–688. [Google Scholar] [CrossRef] [PubMed]
  45. Kinsella, R.J.; Kähäri, A.; Haider, S.; Zamora, J.; Proctor, G.; Spudich, G.; Almeida-King, J.; Staines, D.; Derwent, P.; Kerhornou, A.; et al. Ensembl BioMarts: A hub for data retrieval across taxonomic space. Database 2011, 2011, bar030. [Google Scholar] [CrossRef] [PubMed]
  46. Brown, G.R.; Hem, V.; Katz, K.S.; Ovetsky, M.; Wallin, C.; Ermolaeva, O.; Tolstoy, I.; Tatusova, T.; Pruitt, K.D.; Maglott, D.R.; et al. Gene: A gene-centered information resource at NCBI. Nucleic Acids Res. 2015, 43, 36–42. [Google Scholar] [CrossRef] [PubMed]
  47. Wiedemann, I.; Maehlmeyer, A.; Jansen, S.; Sharifi, A.R.; Knorr, C. SNP g.1007A>G within the porcine DNAL4 gene affects sperm motility traits and percentage of midpiece abnormalities. Reprod. Domest. Anim. 2018, 53, 401–413. [Google Scholar] [CrossRef] [PubMed]
  48. Bates, D.; Mächler, M.; Bolker, B.; Walker, S. Fitting Linear Mixed-Effects Models Using lme4. J. Stat. Softw. 2015, 67, 1–48. [Google Scholar] [CrossRef]
  49. Van Lenth, R. Emmeans: Estimated Marginal Means, aka Least-Squares Means. R Package Version 1.5.3. Available online: https://CRAN.R-project.org/package=emmeans (accessed on 8 January 2021).
  50. Raudvere, U.; Kolberg, L.; Kuzmin, I.; Arak, T.; Adler, P.; Peterson, H.; Vilo, J. g:Profiler: A web server for functional enrichment analysis and conversions of gene lists (2019 update). Nucleic Acids Res. 2019, 47, 191–198. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  51. Kanehisa, M.; Furumichi, M.; Sato, Y.; Ishiguro-Watanabe, M.; Tanabe, M. KEGG: Integrating viruses and cellular organisms. Nucleic Acids Res. 2021, 49, D545–D551. [Google Scholar] [CrossRef]
  52. Reimand, J.; Kull, M.; Peterson, H.; Hansen, J.; Vilo, J. g:Profiler--a web-based toolset for functional profiling of gene lists from large-scale experiments. Nucleic Acids Res. 2007, 35, 193–200. [Google Scholar] [CrossRef]
  53. Sayols, S. Rrvgo: A Bioconductor Package to Reduce and Visualize Gene Ontology Terms. Available online: https://ssayols.github.io/rrvgo (accessed on 7 February 2021).
  54. Day-Williams, A.G.; Southam, L.; Panoutsopoulou, K.; Rayner, N.W.; Esko, T.; Estrada, K.; Helgadottir, H.T.; Hofman, A.; Ingvarsson, T.; Jonsson, H.; et al. A variant in MCF2L is associated with osteoarthritis. Am. J. Hum. Genet. 2011, 89, 446–450. [Google Scholar] [CrossRef] [PubMed]
  55. Estrada, K.; Styrkarsdottir, U.; Evangelou, E.; Hsu, Y.-H.; Duncan, E.L.; Ntzani, E.E.; Oei, L.; Albagha, O.M.E.; Amin, N.; Kemp, J.P.; et al. Genome-wide meta-analysis identifies 56 bone mineral density loci and reveals 14 loci associated with risk of fracture. Nat. Genet. 2012, 44, 491–501. [Google Scholar] [CrossRef] [Green Version]
  56. Zhang, H.; Tao, X.; Wu, J. Association of calcitonin receptor gene polymorphism with bone mineral density in postmenopausal Chinese women: A meta-analysis. Arch. Gynecol. Obstet. 2015, 291, 165–172. [Google Scholar] [CrossRef]
  57. Kiper, P.O.S.; Saito, H.; Gori, F.; Unger, S.; Hesse, E.; Yamana, K.; Kiviranta, R.; Solban, N.; Liu, J.; Brommage, R.; et al. Cortical-Bone Fragility--Insights from sFRP4 Deficiency in Pyle’s Disease. N. Engl. J. Med. 2016, 374, 2553–2562. [Google Scholar] [CrossRef]
  58. Li, H.; Deeb, N.; Zhou, H.; Mitchell, A.D.; Ashwell, C.M.; Lamont, S.J. Chicken quantitative trait loci for growth and body composition associated with transforming growth factor-beta genes. Poult. Sci. 2003, 82, 347–356. [Google Scholar] [CrossRef] [PubMed]
  59. Mullin, B.H.; Walsh, J.P.; Zheng, H.-F.; Brown, S.J.; Surdulescu, G.L.; Curtis, C.; Breen, G.; Dudbridge, F.; Richards, J.B.; Spector, T.D.; et al. Genome-wide association study using family-based cohorts identifies the WLS and CCDC170/ESR1 loci as associated with bone mineral density. BMC Genom. 2016, 17, 136. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  60. Yang, T.L.; Guo, Y.; Liu, Y.J.; Shen, H.; Liu, Y.Z.; Lei, S.F.; Li, J.; Tian, Q.; Deng, H.W. Genetic variants in the SOX6 gene are associated with bone mineral density in both Caucasian and Chinese populations. Osteoporos. Int. 2012, 23, 781–787. [Google Scholar] [CrossRef]
  61. Li, J.; Liu, C.; Li, Y.; Zheng, Q.; Xu, Y.; Liu, B.; Sun, W.; Li, Y.; Ji, S.; Liu, M.; et al. TMCO1-mediated Ca(2+) leak underlies osteoblast functions via CaMKII signaling. Nat. Commun. 2019, 10, 1589. [Google Scholar] [CrossRef] [Green Version]
  62. Borton, A.J.; Frederick, J.P.; Datto, M.B.; Wang, X.F.; Weinstein, R.S. The loss of Smad3 results in a lower rate of bone formation and osteopenia through dysregulation of osteoblast differentiation and apoptosis. J. Bone Miner. Res. 2001, 16, 1754–1764. [Google Scholar] [CrossRef]
  63. Mishra, A.; Awasthi, S.; Raj, S.; Mishra, P.; Srivastava, R.N. Identifying the role of ASPN and COMP genes in knee osteoarthritis development. J. Orthop. Surg. Res. 2019, 14, 337. [Google Scholar] [CrossRef]
  64. Szwedowski, D.; Szczepanek, J.; Paczesny, Ł.; Pękała, P.; Zabrzyński, J.; Kruczyński, J. Genetics in Cartilage Lesions: Basic Science and Therapy Approaches. Int. J. Mol. Sci. 2020, 21, 5430. [Google Scholar] [CrossRef]
  65. Karasik, D.; Rivadeneira, F.; Johnson, M.L. The genetics of bone mass and susceptibility to bone diseases. Nat. Rev. Rheumatol. 2016, 12, 323–334. [Google Scholar] [CrossRef]
  66. Johnsson, M.; Jonsson, K.B.; Andersson, L.; Jensen, P.; Wright, D. Genetic Regulation of Bone Metabolism in the Chicken: Similarities and Differences to Mammalian Systems. PLoS Genetics 2015, 11, e1005250. [Google Scholar] [CrossRef] [Green Version]
  67. Xiao, S.-M.; Kung, A.W.C.; Gao, Y.; Lau, K.-S.; Ma, A.; Zhang, Z.-L.; Liu, J.-M.; Xia, W.; He, J.-W.; Zhao, L.; et al. Post-genome wide association studies and functional analyses identify association of MPP7 gene variants with site-specific bone mineral density. Hum. Mol. Genet. 2012, 21, 1648–1657. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  68. Xie, J.; Guo, J.; Kanwal, Z.; Wu, M.; Lv, X.; Ibrahim, N.A.; Li, P.; Buabeid, M.A.; Arafa, E.-S.A.; Sun, Q. Calcitonin and Bone Physiology: In Vitro, In Vivo, and Clinical Investigations. J. Endocrinol. 2020, 2020, 3236828. [Google Scholar] [CrossRef]
  69. Zofková, I.; Zajícková, K.; Hill, M.; Krepelová, A. Does polymorphism C1377T of the calcitonin receptor gene determine bone mineral density in postmenopausal women? Exp. Clin. Endocrinol. Diabetes 2003, 111, 447–449. [Google Scholar] [CrossRef] [PubMed]
  70. Schinke, T.; Liese, S.; Priemel, M.; Haberland, M.; Schilling, A.F.; Catala-Lehnen, P.; Blicharski, D.; Rueger, J.M.; Gagel, R.F.; Emeson, R.B.; et al. Decreased Bone Formation and Osteopenia in Mice Lacking α-Calcitonin Gene-Related Peptide. J. Bone Miner. Res. 2004, 19, 2049–2056. [Google Scholar] [CrossRef] [PubMed]
  71. Goetz, S.C.; Anderson, K.V. The primary cilium: A signalling centre during vertebrate development. Nat. Rev. Genet. 2010, 11, 331–344. [Google Scholar] [CrossRef]
  72. Ali, I.H.A.; Brazil, D.P. Bone morphogenetic proteins and their antagonists: Current and emerging clinical uses. Br. J. Pharmacol. 2014, 171, 3620–3632. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  73. Hata, A.; Lagna, G.; Massagué, J.; Hemmati-Brivanlou, A. Smad6 inhibits BMP/Smad1 signaling by specifically competing with the Smad4 tumor suppressor. Genes Dev. 1998, 12, 186–197. [Google Scholar] [CrossRef] [Green Version]
  74. Hidalgo-Bravo, A.; Parra-Torres, A.Y.; Casas-Avila, L.; Jimenez-Ortega, R.F.; Ramírez-Salazar, E.G.; Patiño, N.; Rivera-Paredez, B.; Salmerón, J.; Valdés-Flores, M.; Velázquez-Cruz, R. Association of RMND1/CCDC170–ESR1 single nucleotide polymorphisms with hip fracture and osteoporosis in postmenopausal women. Climacteric 2019, 22, 97–104. [Google Scholar] [CrossRef]
  75. Rivadeneira, F.; Styrkársdottir, U.; Estrada, K.; Halldórsson, B.V.; Hsu, Y.-H.; Richards, J.B.; Zillikens, M.C.; Kavvoura, F.K.; Amin, N.; Aulchenko, Y.S.; et al. Twenty bone-mineral-density loci identified by large-scale meta-analysis of genome-wide association studies. Nat. Genetic. 2009, 41, 1199–1206. [Google Scholar]
  76. Smits, P.; Li, P.; Mandel, J.; Zhang, Z.; Deng, J.M.; Behringer, R.R.; de Crombrugghe, B.; Lefebvre, V. The transcription factors L-Sox5 and Sox6 are essential for cartilage formation. Dev. Cell 2001, 1, 277–290. [Google Scholar] [CrossRef] [Green Version]
  77. Sammels, E.; Parys, J.B.; Missiaen, L.; De Smedt, H.; Bultynck, G. Intracellular Ca2+ storage in health and disease: A dynamic equilibrium. Cell Calcium 2010, 47, 297–314. [Google Scholar] [CrossRef]
  78. Shimizu-Hirota, R.; Sasamura, H.; Kuroda, M.; Kobayashi, E.; Saruta, T. Functional characterization of podocan, a member of a new class in the small leucine-rich repeat protein family. FEBS Lett. 2004, 563, 69–74. [Google Scholar] [CrossRef] [Green Version]
  79. Viguet-Carrin, S.; Garnero, P.; Delmas, P.D. The role of collagen in bone strength. Osteoporos. Int. 2006, 17, 319–336. [Google Scholar] [CrossRef] [PubMed]
  80. Sparke, A.J.; Sims, T.J.; Avery, N.C.; Bailey, A.J.; Fleming, R.H.; Whitehead, C.C. Differences in composition of avian bone collagen following genetic selection for resistance to osteoporosis. Br. Poult. Sci. 2002, 43, 127–134. [Google Scholar] [CrossRef] [PubMed]
  81. Bultink, I.E.; Lems, W.F. Osteoarthritis and osteoporosis: What is the overlap? Curr. Rheumatol. Rep. 2013, 15, 328. [Google Scholar] [CrossRef] [PubMed]
  82. Oliveira, M.C.; Vullings, J.; van de Loo, F.A.J. Osteoporosis and osteoarthritis are two sides of the same coin paid for obesity. Nutrition 2020, 70, 110486. [Google Scholar] [CrossRef] [PubMed]
  83. Cong, Y.; Ru, J.Y.; Bao, N.R.; Guo, T.; Zhao, J.N. A single nucleotide polymorphism in the TGF-β1 gene (rs1982073 C>T) may contribute to increased risks of bone fracture, osteoporosis, and osteoarthritis: A meta-analysis. Clin. Rheumatol. 2016, 35, 973–985. [Google Scholar] [CrossRef] [PubMed]
  84. Shepherd, C.; Skelton, A.J.; Rushton, M.D.; Reynard, L.N.; Loughlin, J. Expression analysis of the osteoarthritis genetic susceptibility locus mapping to an intron of the MCF2L gene and marked by the polymorphism rs11842874. BMC Med. Genet. 2015, 16, 108. [Google Scholar] [CrossRef] [Green Version]
  85. Mao, J.H.; Sui, Y.X.; Ao, S.; Wang, Y.; Liu, Y.; Leng, H. miR-140-3p exhibits repressive functions on preosteoblast viability and differentiation by downregulating MCF2L in osteoporosis. In Vitro Cell. Dev. Biol. Anim. 2020, 56, 49–58. [Google Scholar] [CrossRef] [PubMed]
  86. González-Huerta, N.C.; Borgonio-Cuadra, V.M.; Zenteno, J.C.; Cortés-González, S.; Duarte-Salazar, C.; Miranda-Duarte, A. D14 repeat polymorphism of the asporin gene is associated with primary osteoarthritis of the knee in a Mexican Mestizo population. Int. J. Rheum. Dis. 2017, 20, 1935–1941. [Google Scholar] [CrossRef]
  87. Sobhan, M.R.; Mehdinejad, M.; Jamaladini, M.H.; Mazaheri, M.; Zare-Shehneh, M.; Neamatzadeh, H. Association between aspartic acid repeat polymorphism of the asporin gene and risk of knee osteoarthritis: A systematic review and meta-analysis. Acta Orthop. Traumatol. Turc. 2017, 51, 409–415. [Google Scholar] [CrossRef]
  88. Hsu, Y.H.; Kiel, D.P. Clinical review: Genome-wide association studies of skeletal phenotypes: What we have learned and where we are headed. J. Clin. Endocrinol. Metab. 2012, 97, 1958–1977. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  89. Juhász, T.; Matta, C.; Somogyi, C.; Katona, É.; Takács, R.; Soha, R.F.; Szabó, I.A.; Cserháti, C.; Sződy, R.; Karácsonyi, Z.; et al. Mechanical loading stimulates chondrogenesis via the PKA/CREB-Sox9 and PP2A pathways in chicken micromass cultures. Cell. Signal. 2014, 26, 468–482. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  90. Kawano, Y.; Kypta, R. Secreted antagonists of the Wnt signalling pathway. J. Cell Sci. 2003, 116, 2627–2634. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  91. Johnson, M.L. LRP5 and bone mass regulation: Where are we now? BoneKEy Rep. 2012, 1, 1. [Google Scholar] [CrossRef] [Green Version]
  92. Zhang, Y.; Pizzute, T.; Pei, M. A review of crosstalk between MAPK and Wnt signals and its impact on cartilage regeneration. Cell Tissue Res. 2014, 358, 633–649. [Google Scholar] [CrossRef] [Green Version]
  93. Gu, H.; Huang, Z.; Chen, G.; Zhou, K.; Zhang, Y.; Chen, J.; Xu, J.; Yin, X. Network and pathway-based analyses of genes associated with osteoporosis. Medicine 2020, 99, e19120. [Google Scholar] [CrossRef] [PubMed]
  94. Yi, G.; Shen, M.; Yuan, J.; Sun, C.; Duan, Z.; Qu, L.; Dou, T.; Ma, M.; Lu, J.; Guo, J.; et al. Genome-wide association study dissects genetic architecture underlying longitudinal egg weights in chickens. BMC Genom. 2015, 16, 746. [Google Scholar] [CrossRef] [Green Version]
Figure 1. Schematic illustration of the study design. The data set included 524 laying hens phenotyped for bone stability traits. The corresponding genotypes included 490,745 SNP markers. Association analysis was performed applying Random Forests (RF) classification. Genes harbouring significant SNPs were extracted and screened for links to bone stability. Gene set analyses were performed considering all genes obtained from the RF classification. Retrospectively, SNP effects were estimated for a subset of candidate genes identified in gene sets obtained from the RF classifier.
Figure 1. Schematic illustration of the study design. The data set included 524 laying hens phenotyped for bone stability traits. The corresponding genotypes included 490,745 SNP markers. Association analysis was performed applying Random Forests (RF) classification. Genes harbouring significant SNPs were extracted and screened for links to bone stability. Gene set analyses were performed considering all genes obtained from the RF classification. Retrospectively, SNP effects were estimated for a subset of candidate genes identified in gene sets obtained from the RF classifier.
Genes 12 00702 g001
Figure 2. Venn diagram showing the overlap of genes for the bone breaking strengths (BBS) and bone mineral densities (BMD) of the tibiotarsus and humerus.
Figure 2. Venn diagram showing the overlap of genes for the bone breaking strengths (BBS) and bone mineral densities (BMD) of the tibiotarsus and humerus.
Genes 12 00702 g002
Figure 3. Significantly enriched Gene Ontology terms of the category biological processes (GO:BP; top 15 with the highest -log10 p-values) and KEGG pathways for the bone breaking strength (A) and bone mineral density (B) of the tibiotarsus. The dot size represents the absolute number of genes enriched in the term. The proportion of enriched genes in all queried genes is represented on the x-axis. The colour represents the −log10 transformed adjusted p-values.
Figure 3. Significantly enriched Gene Ontology terms of the category biological processes (GO:BP; top 15 with the highest -log10 p-values) and KEGG pathways for the bone breaking strength (A) and bone mineral density (B) of the tibiotarsus. The dot size represents the absolute number of genes enriched in the term. The proportion of enriched genes in all queried genes is represented on the x-axis. The colour represents the −log10 transformed adjusted p-values.
Genes 12 00702 g003
Figure 4. Significantly enriched Gene Ontology terms of the category biological processes (GO:BP; top 15 with the highest -log10 p-values) and KEGG pathways for the bone breaking strength (A) and bone mineral density (B) of the humerus. The dot size represents the absolute number of genes enriched in the term. The proportion of enriched genes in all queried genes is represented on the x-axis. The colour represents the −log10 transformed adjusted p-values.
Figure 4. Significantly enriched Gene Ontology terms of the category biological processes (GO:BP; top 15 with the highest -log10 p-values) and KEGG pathways for the bone breaking strength (A) and bone mineral density (B) of the humerus. The dot size represents the absolute number of genes enriched in the term. The proportion of enriched genes in all queried genes is represented on the x-axis. The colour represents the −log10 transformed adjusted p-values.
Genes 12 00702 g004
Table 2. SNP effects analysis—Part 1: Analysis of variance table and allele substitution effect obtained from the single marker regression model.
Table 2. SNP effects analysis—Part 1: Analysis of variance table and allele substitution effect obtained from the single marker regression model.
SNPTrait 1Candidate GeneGenerationLayer LineSNP GenotypeAllele Substitution Effect 2
F-Statisticsp-ValueF-Statisticsp-ValueF-Statisticsp-ValueBeta (SE 3)Standardised Beta 4 (SE)t-Valuep-Value
AX-76044166Tib_BBSMPP780.92<0.000146.34<0.00014.050.04488.22 (4.09)0.10 (0.05)2.010.0448
AX-76601713Tib_BBSACTR286.02<0.0001106.86<0.000113.330.0003−10.19 (2.79)−0.18 (0.05)−3.650.0003
AX-77276717Tib_BBSTGFB281.07<0.0001102.16<0.00013.320.06964.67 (2.57)0.06 (0.04)1.820.0696
AX-76491534Tib_BBSCCDC17091.49<0.000184.86<0.000112.580.0004−15.63 (4.41)−0.17 (0.05)−3.550.0004
AX-76772658Tib_BBSSOX681.50<0.0001117.84<0.000110.710.00127.63 (2.33)0.12 (0.04)3.270.0012
AX-75711229Tib_BBSASPN79.24<0.000184.23<0.00012.080.15036.66 (4.62)0.05 (0.04)1.440.1503
AX-75913642Tib_BBSSOX983.08<0.0001111.94<0.00019.670.0019−12.87 (4.14)−0.11 (0.04)−3.110.0019
AX-76772658Hum_BBSSOX636.26<0.000152.59<0.00015.670.0177−5.32 (2.23)−0.10 (0.04)−2.380.0177
AX-77091655Hum_BBSPODN39.91<0.000141.64<0.00018.350.00416.69 (2.31)0.13 (0.04)2.890.0041
AX-75597497Hum_BBSSMAD636.38<0.000153.40<0.00014.620.0321−7.13 (3.32)−0.10 (0.05)−2.150.0321
AX-76351785Hum_BBSWNT9B37.27<0.000167.22<0.000121.57<0.000111.51 (2.48)0.21 (0.04)4.64<0.0001
AX-75268181Tib_BMDMCF2L4.300.0401106.46<0.000113.530.0003−0.015 (0.004)−0.15 (0.05)−3.670.0003
AX-80813610Tib_BMDCALCR4.240.041556.10<0.00014.860.02980.008 (0.004)0.10 (0.05)2.210.028
AX-76099065Tib_BMDSFRP44.310.040065.23<0.00018.550.0036−0.016 (0.006)−0.18 (0.06)−2.920.0036
AX-77113061Tib_BMDTMCO14.450.036999.26<0.00015.270.02210.008 (0.003)0.11 (0.05)2.300.0221
AX-75677174Tib_BMDGPATCH14.270.040661.13<0.000110.840.00110.013 (0.004)0.13(0.04)3.290.0011
AX-77091655Hum_BMDPODN20.70<0.000151.56<0.000111.530.00080.007 (0.002)0.14 (0.04)3.390.0008
AX-76351898Hum_BMDWNT319.82<0.000177.58<0.000113.810.00020.016 (0.004)0.15 (0.04)3.720.0002
AX-76351899Hum_BMDWNT319.82<0.000177.58<0.000113.810.00020.016 (0.004)0.15 (0.04)3.720.0002
1 BBS, bone breaking strength; BMD, bone mineral density; Tib, tibiotarsus; Hum, humerus; 2 Allele substitution effect per copy of the effect allele (minor allele); 3 SE, standard error; 4 Standardised regression coefficients expressed in SD unit.
Table 3. SNP effects analysis—Part 2: Genotypic values (least squares means) and additive and dominance effects obtained from the dominant-recessive model.
Table 3. SNP effects analysis—Part 2: Genotypic values (least squares means) and additive and dominance effects obtained from the dominant-recessive model.
SNPTrait 1Candidate GeneGenotypic ValuesHomozygous Additive Allele Effect 5Dominance Effect 5
AA 2,3 (SE 4)AB 2,3 (SE)BB 2,3 (SE)Estimate (SE)t-Valuep-ValueEstimate (SE)t-Valuep-Value
AX-76044166Tib_BBSMPP7155.33 (2.26) ab145.80 (5.85) b172.76 (7.25) a−8.71 (4.05)−2.150.0320−18.20 (5.45)−3.350.0009
AX-76601713Tib_BBSACTR2162.77 (3.08) a156.79 (3.81) a143.10 (3.79) b9.83 (2.82)3.490.00053.86 (4.03)0.960.3392
AX-77276717Tib_BBSTGFB2153.42 (2.25) a157.05 (3.06) a163.72 (5.10) a−5.15 (2.83)−1.820.0694−1.52 (3.73)−0.410.6843
AX-76491534Tib_BBSCCDC170159.13 (2.19) a144.04 (6.28) ab127.83 (8.09) b15.70 (4.42)3.540.00040.56 (5.88)0.0960.9239
AX-76772658Tib_BBSSOX6149.06 (2.65) b158.53 (2.58) a163.29 (3.90) a−7.11 (2.43)−2.930.00352.36 (3.13)0.750.4520
AX-75711229Tib_BBSASPN155.14 (1.94) b154.78 (5.29) b188.53 (13.11) a−16.70 (6.62)−2.520.0120−17.10 (8.02)−2.130.0340
AX-75913642Tib_BBSSOX9157.50 (1.93) a148.13 (4.83) ab124.13 (10.37) b16.70 (5.31)3.140.00187.40 (6.44)1.150.2506
AX-76772658Hum_BBSSOX6127.04 (2.51) a116.24 (2.46) b119.38 (3.71) ab3.83 (2.31)1.660.0984−6.96 (3.02)−2.310.0215
AX-77091655Hum_BBSPODN118.01 (2.31) b120.73 (3.04) b132.21 (3.73) a−7.10 (2.31)−3.070.0023−4.38 (3.44)−1.270.2043
AX-75597497Hum_BBSSMAD6122.16 (2.08) a123.48 (3.64) a98.16 (6.97) b12.0 (3.71)3.230.001313.30 (4.61)2.880.0040
AX-76351785Hum_BBSWNT9B115.73 (2.19) c124.86 (3.05) b139.61 (4.34) a−11.90 (2.54)−4.70<0.0001−2.81 (3.49)−0.800.4215
AX-75268181Tib_BMDMCF2L0.263 (0.003) a0.253 (0.005) a0.228 (0.008) b0.017 (0.004)3.920.00010.008 (0.006)1.350.1768
AX-80813610Tib_BMDCALCR0.256 (0.003) a0.258 (0.005) a0.273 (0.006) a−0.009 (0.004)−2.240.0257−0.007 (0.005)−1.270.2051
AX-76099065Tib_BMDSFRP40.261 (0.003) ab0.265 (0.008) a0.235 (0.009) b0.013 (0.006)2.320.02060.018 (0.006)2.710.0071
AX-77113061Tib_BMDTMCO10.246 (0.005) a0.267 (0.004) a0.266 (0.004) a−0.01 (0.004)−2.820.00500.011 (0.004)2.510.0125
AX-75677174Tib_BMDGPATCH10.254 (0.003) b0.269 (0.005) a0.278 (0.007) a−0.012 (0.004)−3.050.00240.004 (0.005)0.560.5739
AX-77091655Hum_BMDPODN0.164 (0.002) b0.167 (0.003) b0.178 (0.003) a−0.007 (0.002)−3.530.0005−0.004 (0.003)−1.250.2117
AX-76351898Hum_BMDWNT30.166 (0.002) b0.176 (0.006) b0.206 (0.010) a−0.02 (0.005)−3.840.0001−0.009 (0.007)−1.290.1991
AX-76351899Hum_BMDWNT30.166 (0.002) b0.176 (0.006) b0.206 (0.010) a−0.02 (0.005)−3.840.0001−0.009 (0.007)−1.290.1991
1 BBS, bone breaking strength; BMD, bone mineral density; Tib, tibiotarsus; Hum, humerus; 2 AA or BB represents the homozygote of the other allele or effect allele, respectively. AB denotes the heterozygote (see Table 1 for the actual genotypes); 3 Means with different letters within a column differ significantly at p < 0.05; 4 SE, standard error; 5 Effect of the other allele (major allele).
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Jansen, S.; Baulain, U.; Habig, C.; Ramzan, F.; Schauer, J.; Schmitt, A.O.; Scholz, A.M.; Sharifi, A.R.; Weigend, A.; Weigend, S. Identification and Functional Annotation of Genes Related to Bone Stability in Laying Hens Using Random Forests. Genes 2021, 12, 702. https://0-doi-org.brum.beds.ac.uk/10.3390/genes12050702

AMA Style

Jansen S, Baulain U, Habig C, Ramzan F, Schauer J, Schmitt AO, Scholz AM, Sharifi AR, Weigend A, Weigend S. Identification and Functional Annotation of Genes Related to Bone Stability in Laying Hens Using Random Forests. Genes. 2021; 12(5):702. https://0-doi-org.brum.beds.ac.uk/10.3390/genes12050702

Chicago/Turabian Style

Jansen, Simon, Ulrich Baulain, Christin Habig, Faisal Ramzan, Jens Schauer, Armin Otto Schmitt, Armin Manfred Scholz, Ahmad Reza Sharifi, Annett Weigend, and Steffen Weigend. 2021. "Identification and Functional Annotation of Genes Related to Bone Stability in Laying Hens Using Random Forests" Genes 12, no. 5: 702. https://0-doi-org.brum.beds.ac.uk/10.3390/genes12050702

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop