Screening for potential genes associated with bone overgrowth after mid-shaft femur fracture in a rat model

Background We investigated the underlying molecular mechanisms of bone overgrowth after femoral fracture by using high-throughput bioinformatics approaches. Methods The gene expression profile of GSE3298 (accession number) was obtained from the Gene Expression Omnibus database. Sixteen femoral growth plate samples, including nine samples without fracture and seven fracture samples for seven time points, were used for analysis. The Limma package was applied to identify differentially expressed genes (DEGs) between fractured and intact samples. The DAVID online tool was used for Gene ontology functional and pathway enrichment analysis. A protein-protein interaction (PPI) network established by String software was used to identify interactions between significant DEGs, and network modules were detected using plug-in MCODE. Additionally, a transcription regulatory network was constructed based on the ENCODE Project and PPI network. Results A total of 680 DEGs were screened in fractured femoral growth plate samples compared with controls, including 238 up- and 442 down-regulated genes. These DEGs were significantly involved in the calcium signaling pathway and cancer pathway. A PPI network was constructed with 167 nodes and 233 edges, and module analysis demonstrated that CCL2, CSF2, NOS2, and DLC1 may stimulate bone overgrowth after femoral fracture via anti-apoptosis-related functions. A transcription regulatory network was constructed with 387 interacting pairs, and overlapping nodes were significantly enriched in intracellular signaling cascade and regulation of cell proliferation, among others. Conclusions Bone overgrowth was associated with changes in the expression of identified DEGs such as CCL2, NOS2, CSF2, and DLC1 in the femoral head. They may be important in regulating bone overgrowth via the anti-apoptosis of osteoblasts. Electronic supplementary material The online version of this article (doi:10.1186/s13018-017-0510-6) contains supplementary material, which is available to authorized users.


Background
Femoral fracture, which is one of the most commonly occurring fractures during childhood, always results from casual falls, motor vehicle accidents, or sporting accidents [1]. Treatment of femoral fracture typically includes open-reduction, traction, and internal fixation. However, pediatric femoral fracture often results in the stimulation of bone overgrowth, particularly in children younger than 12 years [2,3]. Overgrowth is described as a universal phenomenon in patients with femoral shaft fractures and can elongate the lower limb by nearly 9 mm [4] or 11 mm [5]. It is crucial to explore the underlying molecular mechanism of bone overgrowth associated with femoral fracture.
In recent years, numerous studies have investigated the molecular mechanism of bone overgrowth after femoral fracture [6][7][8]. Bone homeostasis is thought to be maintained by a balance between bone formation by osteoblasts and bone resorption by osteoclasts in the growth plate. Various proteins such as β-catenin and triggering receptor expressed by myeloid cells-2 interact with each other by controlling the rate of osteoclastogenesis and further regulating bone homeostasis [9]. Several other factors were also found to be involved in this process. For instance, lipoprotein receptor-related protein 4 (LRP4) was found to be associated with the inhibitory function of sclerostin which is secreted by osteocytes and inhibits bone formation [6]. Additionally, the Wnt1/β-catenin signaling pathway is crucial for embryonic and bone homeostasis [10][11][12], and LRP4 may increase sclerostin secretion through Wnt1/β-catenin signaling [7]. Moreover, fibroblast growth factor receptor (FGFR) is also a critical gene in bone overgrowth and participates in FGFR3 signaling, further affecting chondrocyte proliferation [8]. Another report showed that osteocrin is highly expressed in osteoblasts and interacts with C-type natriuretic peptide receptors to modulate the action of the natriuretic system during bone elongation [13]. Therefore, expression changes of such related genes in cells after femoral fracture may provide insight into the physiological mechanisms of bone overgrowth.
In the past few years, DNA microarray technology has been increasingly utilized to comprehensively test for changes in the messenger RNA (mRNA) expression of genes and search for evidence of overgrowth after midfemoral fracture [3,14]. However, the potential molecular mechanism of bone overgrowth after mid-shaft femur fracture also remains unclear. The aim of this study was to explore potentially important genes associated with bone overgrowth after femoral fracture and clarify this phenomenon using high-throughput bioinformatics methods.

Data source
The gene expression profile of GSE3298 [3,14], which describes mRNA expression in the rat proximal femoral growth plate after mid-shaft fracture, was derived from the Gene Expression Omnibus (GEO, http://ncbi.nlm. nih.gov/geo/) database based on the GPL1355 Affy metrix Rat Genome 230 2.0 Array platform (Santa Clara, CA, USA). A total of 16 femoral growth plate samples were used for analysis, including nine samples without fracture and seven fracture samples for seven time points: 1 day, 3 days, 1 week, 2 weeks, 3 weeks, 4 weeks, and 6 weeks, after mid-shaft fracture.

Data preprocessing and differentially expressed genes (DEGs) screening
Each sample in the obtained dataset had a probe ID, which was converted into the corresponding gene name. Multiple probe IDs targeting the same gene were averaged as the gene expression value. After expression values were log 2 transformed, quantile normalization was carried out [15]. The Limma package (http://www.bioconductor.org/ packages/release/bioc/html/limma.html) [16] in R language was used to screen DEGs between intact and fractured samples. Gene expression differences were assessed using Student's t test, and expression changes were considered to be significant when by P < 0.05.

Functional annotation and pathway enrichment of DEGs
Functional annotation of genes was carried out using the Database for Annotation, Visualization, and Integrated Discovery (DAVID, http://david.abcc.ncifcrf.gov/) [17]. Gene ontology (GO, http://geneontology.org/) annotation and pathway enrichment analysis were performed to derive all associated functions with their enrichment scores and P values. Fisher's exact test was used to evaluate the differences between the intact and fractured femora. Only results showing enrichment scores of more than 2 and P values <0.05 were considered to be statistically significant.

Construction of protein-protein interaction (PPI) network
The PPI is thought to be important for understanding the potential functions of a certain protein. The up-and down-regulated genes identified as described above were respectively mapped to the Search Tool for the Retrieval of Interacting Genes (STRING, http://string-db.org/) software which is commonly used to predict PPI pairs [18]. The PPI network was constructed with interesting PPI pairs and visualized by Cytoscape 2.8 (http:// cytoscape.org/) [19].

Module detection
Molecular complex detection (MCODE) [20] is a clustering algorithm used to identify molecular complexes in the PPI network. Degree ≥2 and k-score ≥2 were selected as cutoff criteria. Next, Bingo [21] was applied to perform GO functional enrichment analysis with a threshold of adjusted P value <0.05, with multiple test adjustment conducted as described by Benjamini-Hochberg.

Transcription regulatory network
The ENCODE (ENCyclopedia of DNA Elements) Project is designed to identify all functional components in the human genome sequence, including proteincoding genes, non-protein-coding genes, sequences that mediate chromosome structure and dynamics, and transcriptional regulatory elements [22]. Based on transcription factor information determined using EN-CODE and the PPI network, a transcription regulatory network was constructed with interacting pairs using Cytoscape. Additionally, overlapping nodes in the network with node degrees of >2 were further analyzed for functional enrichment using DAVID. Only the results showing P values <0.05 were considered to be statistically significant.

Identification of DEGs
After data preprocessing and quantile normalization, the gene expression profile of GSE3298 was used to screen for DEGs in the proximal femoral growth plate between intact and fractured samples. A total of 680 DEGs were screened out with P < 0.05, including 238 up-and 442 down-regulated genes (Additional file 1: Table S1).

Functional analysis and pathway enrichment for DEGs
To evaluate DEG functions, GO and pathway analysis were performed for up-regulated and down-regulated genes, respectively. As shown in Table 1, up-regulated genes were mainly enriched in seven GO terms including cell fraction, response to organic substance, and response to wounding, among others, while downregulated genes mainly function in the extracellular region, vesicles, cytoplasmic vesicles, and membranebound vesicles, among others. Additionally, up-regulated genes were found to be mainly enriched in six pathways such as the calcium signaling pathway (P = 0.001) and neuroactive ligand-receptor interactions (P = 0.001), while the 442 down-regulated DEGs were significantly enriched in four pathways, including pathways in cancer (P = 0.008), calcium signaling pathway (P = 0.016), hedgehog signaling pathway (P = 0.025), and MAPK signaling pathway (P = 0.028) ( Table 2). These results suggest that perturbations in genes involved in these functions and/or pathways are associated with bone overgrowth following femoral fracture.

Construction of PPI network
To build the PPI network, significant protein interactions were predicted; the results are displayed in Additional file 2: Table S2. Based on these interaction pairs, a PPI network was constructed with 167 nodes and 233 edges (Fig. 1). Among all nodes, three DEGs showed relatively higher degrees, including chemokine (C-C motif ) ligand 2 (CCL2), nitric oxide synthase 2 (NOS2), and colony-stimulating factor 2 (CSF2). Moreover, GO analysis suggested that up-regulated CCL2 participated in the chemokine signaling pathway, and down-regulated NOS2 involved cancer and calcium signaling pathways. In addition, down-regulated CSF2 was mainly enriched in the hematopoietic cell lineage pathway.

Transcription regulatory network
Based on information of transcription factors from ENCODE and the PPI network, 387 interacting pairs (Additional file 3: Table S3) were screened out for construction of the transcription regulatory network (Fig. 3). Analysis of functional enrichment showed that the overlapping nodes were significantly enriched in GO terms such as intracellular signaling cascade, regulation of cell proliferation, and regulation of apoptosis, among others.

Discussion
Currently, the incidence of fractures of proximal femur have increased as industrial societies become older [23]. Femur fracture is thought to be associated with bone overgrowth, which is a common phenomenon, particularly in children. However, the underlying mechanism remains unclear. In this study, we identified potential genes involved in the molecular mechanism of bone overgrowth after femoral fracture in juvenile rats by using high-throughput bioinformatics. Based on gene expression profiles, a total of 680 DEGs were screened out, including 238 up-and 442 down-regulated DEGs. The up-regulated DEGs were found to be significantly enriched in six pathways, while down-regulated DEGs were strikingly enriched in the cancer pathway and calcium signaling pathway. PPI network construction accompanied by module detection revealed key genes such as CCL2, CSF2, NOS2, and DLC1 were identified to be potentially related with femoral overgrowth. CCL2, also known as monocyte chemoattractant protein-1 and small inducible cytokine A2, is a chemokine ligand and plays a crucial role in the recruitment and activation of macrophages/monocytes during inflammation after bone injury [24]. It is well known that activation of macrophages and monocytes can stimulate osteoclastic bone resorption or bone formation [25,26]. However, the recruitment of macrophages and monocytes to the inflamed bone by CCL2 is regulated by rhTNF to regulate bone formation and further improve fracture healing, which only occurs in the fractured environment [27][28][29]. Additionally, the production of CCL2 is stimulated by the receptor-activator of nuclear factor (NF)-κB ligand, which is regarded as an essential regulator of bone remodeling [30]. In the present study, CCL2 was found to be up-regulated at the proximal femoral growth plate of mid-shaft fracture samples compared with that of no fracture samples and was predicted to participate in cell functions; CCL2 was a key node in the PPI network constructed of significant interaction pairs of DEGs. Thus, CCL2 may be involved in bone overgrowth after femora fracture via recruiting macrophages and monocytes to fractured bone to accelerate bone formation.
CSF2, also known as granulocyte macrophage colonystimulating factor (GM-CSF), encodes a monomeric Fig. 1 Protein-protein interaction network for differentially expressed genes (DEGs). Red circles indicate protein products of significant DEGs; blue lines indicate interaction between different proteins glycoprotein which is regarded as a hemopoietic growth factor. CSF2 is released by osteoblast lineage cells [31]. CSF2 is associated with the production, differentiation, and function of granulocytes and macrophages in vitro [32]. The fusion of monocytes/macrophages can form osteoclasts, which were demonstrated to function in degrading bone [31,33]. Additionally, CSF2 was confirmed to be a target of NF-κB for inducing osteoclastogenesis and further promoting osteolytic bone metastasis [34]. CSF2 has also been reported to enhance osteoclast development which was mediated by tumor necrosis factor α [35]. In addition, CSF2 is an anti-apoptotic factor that can minimize the extent of cell death (such as osteoclasts) in tissues surrounding the injured areas [36]. In this study, CSF2 was identified to be down-regulated. As a result, reduced CSF2 may negatively regulate osteoclastogenesis, resulting in partial recovery of bone formation.
NOS2 is an isoenzyme of nitric oxide (NO) synthase and plays an important role in producing NO, a multifunctional signal molecule. Osteoclasts have been confirmed to express NOS2 and release NO in a regulated manner [37][38][39][40]. In addition, NO is involved in the mechanism of osteoclastic activity by releasing boneresorbing inflammatory cytokines [41]. Endogenously produced NO exerts potent biphasic actions that may significantly affect the proliferation, recruitment, differentiation, and/or survival of osteoblasts and osteoclasts [42][43][44]. Low levels of NO may support osteoblast bone formation and osteoclast-mediated bone remodeling as well as protect osteoblasts against apoptosis, while high NO levels inhibit osteoclastogenesis and prevent bone loss [45][46][47]. Additionally, NO is involved in the control of Ca 2+ dynamics and mediates Ca 2+ -inhibited bone  [39,48]. Furthermore, epidermal growth factor receptor/signal transducers and activators of transcription 3 can interact with the NOS2 promoter and activate NOS2 expression [49]. In this study, NOS2 was predicted to participate in anti-apoptosis-related functions and was found to be significantly enriched in the calcium pathway and cancer pathway. Therefore, NOS2 may participate in bone overgrowth after femora fracture via suppressing osteoblast apoptosis through the cancer pathway and calcium pathway. However, the regulation of NOS2 in bone overgrowth after fracture requires further investigation. DLC1 encodes a Rho GTPase-activating protein that regulates osteoclastogenesis via Rho protein signal transduction [50]. Moreover, Rho GTPases have been confirmed to play a critical role in regulating the actin cytoskeleton organization of osteoclasts [51]. In the present study, DLC1 was found to be down-regulated and enriched in the functions of regulation of Rho protein signal transduction in module C. Thus, this gene may also be important for reducing osteoclastogenesis and bone resorption after femur fracture.

Conclusions
The identified DEGs, particularly those in significant gene modules, including CCL2, CSF2, NOS2, and DLC1, may play a vital role in bone overgrowth after mid-shaft femur fracture. Experimental studies including samples of a larger size will be performed in the future. These data underscore the complexity of the regulation of bone overgrowth. Additionally, these findings form a basis for future studies focusing on the role of these key genes in the molecular mechanisms of bone growth disturbances with the longer-term goal of investigating proper treatment for children with fractured growing bones.

Acknowledgements
The study was supported by a special fund for medical service of the Jilin Finance Department, project number: SCZSY201507.

Funding
The study was supported by the special fund for medical service of Jilin Finance Department, project number: SCZSY201507.
Availability of data and materials All data generated or analyzed during this study are included in this published article and its supplementary information files.