Abstract

Objective: Heterotopic ossification (HO) arises from inflammation, microenvironmental disturbance, and abnormal osteogenic differentiation in soft tissues. MicroRNA-155 participates in inflammatory and bone metabolic processes, but its regulatory mechanism in HO remains unclear. This study adopted bioinformatics methods to screen functional target genes of miR-155 and explore their potential biological functions.

Methods: Differentially expressed genes (DEGs) related to HO were extracted from the Gene Expression Omnibus (GEO) database. Three databases, including Target Scan, miRDB, and miRTar Base, were jointly used to predict miR-155 target genes. Key genes were obtained via the intersection of predicted targets and HO-related DEGs. Functional enrichment, protein-protein interaction network, and miRNA binding prediction analyses were performed on candidate genes.

Results: A total of 435 HO-associated DEGs and 132 shared target genes of miR-155 were identified. Only two overlapping core genes, muscle blind-like splicing regulator 3 (MBNL3) and MYB proto-oncogene-like 1 (MYBL1), were screened out. MBNL3 mainly mediates post-transcriptional regulation via RNA splicing, while MYBL1 controls transcription and cell cycle processes. Both genes contain credible binding sites for hsa-miR-155-5p.

Conclusion: miR-155 could regulate HO through dual transcriptional and post-transcriptional pathways mediated by MYBL1 and MBNL3. This research provides candidate genes for follow-up experimental verification.

Keywords: miR-155; heterotopic ossification; bioinformatics analysis; regulation mechanisms

Abbreviations: HO: Heterotopic ossification; miR-155: microRNA-155; hsa-miR-155-5p: Homo sapiens microRNA-155-5p; DEGs: Differentially expressed genes; GEO: Gene Expression Omnibus; HGNC: HUGO Gene Nomenclature Committee; GO: Gene ontology; KEGG: Kyoto encyclopedia of genes and genomes; PPI: Protein–protein interaction; STRING: Search tool for the retrieval of interacting genes/proteins; 3'UTR: 3' untranslated region; ΔGhybrid: Hybridization free energy; BMP: Bone morphogenetic protein; TGF-β: Transforming growth factor beta; RUNX2: Runt-related transcription factor 2; ACVR1/ALK2: Activin A receptor type 1 / activin receptor-like kinase 2; Wnt/β-catenin: Wingless-related integration site / beta-catenin; PI3K/Akt: Phosphoinositide 3-kinase / protein kinase B; TNF-α: Tumor necrosis factor alpha; IL-1β/IL-6: Interleukin-1 beta/interleukin-6; SOCS1: Suppressor of cytokine signaling 1;SHIP1 – Src homology 2-containing inositol phosphatase 1; SMAD5 – Mothers against decapentaplegic homolog 5; MITF: Microphthalmia-associated transcription factor; qPCR: Quantitative polymerase chain reaction; MBNL3: Muscle blind-like splicing regulator 3; MBNL1 / MBNL2: Muscle blind-like splicing regulator 1/2; CHCR: Cys3His CCG1-required protein; CELF1/2/3/4: CUGBP Elav-like family member 1/2/3/4; DMPK: Dystrophia myotonica protein kinase; MYBL1 (A-MYB): MYB proto-oncogene-like 1; DREAM complex: DP, RB-like, E2F, and MuvB complex; MuvB: Multicellular vulval abnormal B complex; LIN9/37/52/54: Lin-9/37/52/54, DREAM complex subunit; RBBP4: Retinoblastoma binding protein 4; CREBBP: CREB binding protein; SASP: Senescence-associated secretory phenotype; RB: Retinoblastoma protein; CDK: Cyclin-dependent kinase; mTOR: Mechanistic target of rapamycin; FOXO: Forkhead box O transcription factor; ROS: Reactive oxygen species; CNKI: China National Knowledge Infrastructure; FDR: False discovery rate

Introduction

Heterotopic ossification (HO) refers to the formation of mature bone tissue within soft tissues where bone is not normally found, including muscles, fascia, tendons, ligaments, and peri-articular tissues. It is a process in which local cells undergo abnormal osteogenic differentiation in response to injury or other stimuli, ultimately forming osteoid-like tissue or even mature bone tissue [1]. Based on etiology and clinical origin, HO is generally classified into acquired HO and hereditary HO. Acquired HO is more common clinically with frequent triggers, including trauma, fractures, joint replacement surgery, burns, and central nervous system injuries. Its core mechanism involves an external injury provoking a local inflammatory response, which subjects local precursor cells, such as mesenchymal cells, to abnormal stimuli, driving them to differentiate into chondrogenic or osteogenic lineages. Hereditary HO is relatively rare, arising from genetic mutations that render osteogenesis-related signaling pathways hypersensitive, with a classic example named fibrodysplasia ossificans progressiva (FOP) [2].

Pathological Basis of Heterotopic Ossification

Following soft tissue injury, the repair process may deviate from its normal course, leading to ossification in abnormal locations. Inflammatory responses frequently occur in the injured area after trauma, surgery, or burns, while nervous system injuries can also trigger similar changes. Upon entering the injury site, inflammatory cells release various cytokines, altering the local microenvironment. Persistent stimulation disturbs tissue repair, redirecting mesenchymal precursors toward osteogenic differentiation and triggering ectopic bone formation [3].

Pathologically, the development of HO proceeds through several sequential stages. The early stage is primarily characterized by an inflammatory response and cell recruitment. Subsequently, signaling molecules such as BMP and TGF-β drive the proliferation of cellular components, including fibroblast-like and mesenchymal cells, and their gradual transformation into chondroid or osteoid tissues. Finally, these tissues mature further to form stable bone tissue [3]. While sharing similarities with normal endochondral ossification, it occurs outside the skeletal system and is thus classified as abnormal bone formation.

Molecular regulatory mechanisms of heterotopic ossification

Various molecular signals influence the progression of heterotopic ossification. Among the multiple pathways involved, the BMP pathway has garnered significant attention. Upon BMP binding to its receptor, Smad1/5/8 proteins are activated, leading to upregulated expression of RUNX2 and Osterix, which drives mesenchymal cells toward osteogenic differentiation [4]. BMP signaling abnormalities are particularly prominent in hereditary forms of heterotopic ossification. For instance, in fibrodysplasia ossificans progressiva (FOP), mutations in ACVR1/ALK2 result in sustained BMP signaling enhancement, a dysregulation linked to abnormal osteogenesis in soft tissues [5].

In addition, TGF-β, Wnt/β-catenin, and PI3K/Akt pathways are also implicated in heterotopic ossification. TGF-βis frequently involved in injury repair, while Wnt/β-catenin regulates the transcription of osteogenic genes and influences the differentiation fate of osteoprogenitor cells; meanwhile, the PI3K/Akt pathway modulates cell proliferation, survival, and differentiation [6].

The inflammatory response is another critical component of the early stages of heterotopic ossification. Following tissue injury, local expression levels of inflammatory cytokines such as TNF-α, IL-1β, and IL-6 increase, accompanied by the infiltration of inflammatory cells like macrophages. These inflammation-related changes alter the local microenvironment and affect the expression of osteogenesis-related genes [3]. Furthermore, recent studies suggest that microRNAs may also play a regulatory role in heterotopic ossification by targeting specific mRNAs regulating inflammation, cell differentiation, and osteogenesis-related processes at the post-transcriptional level [7]. Therefore, the formation of heterotopic ossification could not be simply attributed to an abnormality in a single signaling pathway; it is rather likely the result of the interplay among multiple factors, including inflammatory responses, osteogenic signaling pathways, and non-coding RNA regulation.

Biological functions of miR-155

MicroRNAs are a class of endogenous, small non-coding RNAs that primarily regulate gene expression at the post-transcriptional level by binding to target gene mRNAs, leading to the inhibition of protein translation and potential mRNA degradation [9]. miR-155 is a microRNA that has been extensively studied in recent years, playing a crucial regulatory role in various physiological and pathological processes such as inflammatory responses, immune regulation and cell differentiation [8]. Studies have identified SOCS1 and SHIP1 as key target genes of miR-155. By modulating these genes, miR-155 participates in inflammatory responses and immune cell activation, and is also associated with alterations in cellular function [8].

In the context of bone-related research, miR-155 regulates bone metabolism through effects on both osteogenesis and osteoclastogenesis. On one hand, miR-155 targets SMAD5, thereby inhibiting BMP-related signaling and reducing the expression of osteogenic molecules such as Runx2, which suppresses osteogenic differentiation [10]. On the other hand, miR-155 targets SOCS1 and MITF, influencing the process of osteoclast differentiation [11]. Given that heterotopic ossification involves both the above local inflammatory responses and aberrant osteogenic differentiation processes, miR-155 might play a potential regulatory function in heterotopic ossification and warrants further investigation.

Research on miR-155 in bone-related diseases includes osteoarthritis, rheumatoid arthritis, and osteogenic/osteoclastic differentiation. D'Adamo et al. found that in osteoarthritis, miR-155 is associated with impaired chondrocyte autophagy and cellular damage [12]. In rheumatoid arthritis, it is linked to an enhanced inflammatory response [13], and in studies regarding bone formation, miR-155 influences osteogenic differentiation by targeting molecules such as SMAD5[10]. These studies demonstrate that miR-155 is involved not only in the regulation of inflammation but also in bone metabolism and cellular differentiation.

Currently, there is a lack of direct research linking miR-155 to heterotopic ossification. Bioinformatic screening, including multi-database prediction and DEG intersection analysis, can narrow candidate genes and infer potential regulatory pathways.

Materials and Methods

A series of public bioinformatics databases and online platforms were adopted for data acquisition, miRNA target prediction, intersection screening, functional enrichment, protein interaction, and binding site prediction in this study.

Data Collection

High-throughput transcriptome datasets related to heterotopic ossification were retrieved and downloaded from the GEO database [14]. After screening, the dataset GSE94683 containing heterotopic ossification expression profiles was selected for subsequent analysis, with its sample metadata reviewed. Meanwhile, literature retrieval targeting miR-155, heterotopic ossification, inflammation, and osteogenic differentiation was conducted on CNKI and PubMed to clarify the regulatory association between miR-155 and bone differentiation/inflammatory response.

Prediction of miR-155 Target Genes

Three authoritative miRNA target prediction databases, Target Scan, miRDB, and miRTar Base, were jointly used to predict target genes of hsa-miR-155-5p. Target Scan predicts targets based on seed sequence conservation, miRDB provides machine learning scoring results, and miRTar Base collects experimentally validated miRNA-mRNA interactions. Raw prediction data downloaded from the three databases were standardized to unified gene symbols and deduplicated.

All processed gene lists were imported into an online Venn diagram analysis platform to extract intersecting genes shared by all three databases, which were regarded as reliable candidate targets of miR-155 and could reduce prediction bias caused by distinct database algorithms [15].

Intersection Analysis

Differentially expressed genes (DEGs) linked to heterotopic ossification were screened from GEO datasets and standardized to HGNC gene symbols. The candidate target gene list of miR-155 and the DEG list were uploaded to the online Venn tool to obtain overlapping hub genes, which were exported and preserved for follow-up functional analysis.

Functional Enrichment Analysis

The overlapping hub genes were subjected to GO and KEGG enrichment analysis with Homo sapiens as the reference species [16]. GO analysis annotates genes into biological processes, cellular components, and molecular functions, while KEGG analysis identifies enriched signaling pathways. Enrichment results focusing on heterotopic ossification, cell cycle, transcriptional regulation, and RNA processing were visualized as bar charts via the Micro-Bioinformatics online plotting platform.

Protein-Protein Interaction (PPI) Network and Pathway Analysis

Hub genes MBNL3 and MYBL1 were submitted to the STRING database to construct a human PPI network under medium confidence (interaction score ≥ 0.400). Basic network metrics, including node count, edge number, average node degree, and PPI enrichment p-value, were recorded to interpret gene interaction modules and biological functions.

Combined with Reactome pathway database annotation, the biological processes and signaling cascades participated by hub genes were further supplemented to support the GO/KEGG enrichment outcomes [17].

Target Binding Prediction

Transcript sequences and 3′ UTR regions of MBNL3 and MYBL1 were extracted from the NCBI database. The STarMir online tool was used to predict binding sites between hsa-miR-155-5p and target mRNA sequences under human species mode. Key predictive indicators, including seed match type, logistic binding probability, and ΔGhybrid, were recorded, and hybrid secondary structure diagrams were saved. Higher logistic probability and lower ΔGhybrid indicate more stable miRNA-mRNA binding.

Results

Prediction of miR-155 target genes

Target Scan, miRDB, and miRTar Base adopt distinct prediction algorithms, leading to inconsistent target counts. We extracted their 132 overlapping genes to reduce prediction bias (Figure 1).

Subsequent analysis did not rely directly on results from a single database; instead, it utilized these 132 common target genes as a foundation for further comparison with differentially expressed genes associated with heterotopic ossification.

The number of target genes identified by the three databases varied, a discrepancy attributable to the specific algorithms employed by each. Target Scan prioritizes the matching between the miRNA seed sequence and the target gene's 3' UTR [18]. MiRDB generates predictions primarily based on model scores [19], and miRTarBase catalogs miRNA-target gene interactions supported by existing experimental evidence [20]. Consequently, relying on any single database makes the results susceptible to the limitations of its specific algorithm or scope of coverage.

Screening and analysis of key overlapping genes

After compiling the list of candidates miR-155 target genes, the next step was to determine whether these genes were also associated with heterotopic ossification. Previously, a search of the GEO database for expression data related to heterotopic ossification, specifically the GSE94683 dataset, yielded 435 differentially expressed genes associated with the condition. Subsequently, these 435 differentially expressed genes and the 132 candidate miR-155 target genes were analyzed using an online Venn diagram platform to identify overlaps. This analysis yielded only two overlapping genes: MBNL3 and MYBL1 (Figure 2).

This step further narrowed the scope of the initial prediction; the intersection analysis retained MBNL3 and MYBL1, identifying them as both candidate miR-155 targets and differentially expressed genes in heterotopic ossification. Only two intersecting genes were identified under strict screening, so all subsequent analyses centered on MBNL3 and MYBL1. Functional and pathway analyses were conducted to further determine the biological processes in which these genes might be involved.

Bioinformatics analysis of MBNL3

Functional annotation of the MBNL3 gene

MBNL3 (muscle blind-like splicing regulator 3) is a protein-coding gene belonging to the muscle blind-like family. According to functional annotations from NCBI Gene and UniProt, MBNL3 is primarily associated with processes such as RNA binding and the regulation of alternative splicing. The protein it encodes contains a C3H1-type zinc finger domain-a structural motif typically involved in RNA recognition and binding [21]. Previous studies have identified the Muscle blind family as key regulators of alternative splicing during developmental stages. Structurally, MBNL family proteins typically contain pairs of zinc-finger RNA-binding domains. Research indicates that while the zinc-finger regions of MBNL1 and MBNL3 are primarily responsible for RNA binding, their splicing regulatory functions also involve other structural domains. This suggests that MBNL3’s role extends beyond mere RNA binding; alterations in its expression or activity can influence exon selection in certain pre-mRNAs and the composition of mature transcripts.

Clues regarding MBNL3's function extend beyond RNA splicing to include potential links with cellular differentiation states. Previous studies have associated the Muscle blind family with cell fate determination, noting roles in tissue-specific expression and the maintenance of differentiation states [21]. Consequently, it is hypothesized that MBNL3 activity may be influenced by the tissue environment and that its functional emphasis may shift across different developmental stages.

In earlier literature, MBNL3 was referred to as CHCR. Experiments have shown that sustained expression of CHCR/MBNL3 hinders the entry of C2C12 cells into the myogenic differentiation pathway [22]. In models of myogenesis, MBNL3 has been linked to the regulation of alternative splicing, which in turn affects cellular differentiation status [23]. While these findings do not directly explain heterotopic ossification, they establish a functional link between MBNL3 and the regulation of cellular states. This suggests that MBNL3 may influence RNA splicing and changes in cellular state at the post-transcriptional level.

Analysis of the MBNL3-related protein-protein interaction network

This study conducted a protein-protein interaction (PPI) analysis centered on MBNL3. The STRING database was used to construct the interaction network [24]. The interaction confidence level was set to "medium confidence," corresponding to a combined score of ≥ 0.400. As shown in Figure 3, the network comprises 21 nodes and 100 interaction edges, with an average node degree of 9.52; the database-expected number of edges was 22, and the PPI enrichment p-value was < 1.0×10-16. Proteins closely associated with MBNL3 in the network include MBNL1, MBNL2, CELF1, CELF2, CELF3, CELF4, QKI, ADAR, CNBP, and DMPK. An analysis of these nodes reveals that most are involved in RNA binding, RNA processing, or the regulation of splicing. This aligns with the earlier functional annotations and further supports the conclusion that MBNL3 primarily participates in RNA splicing and post-transcriptional regulatory processes.

GO/KEGG enrichment analysis of MBNL3

Based on the PPI network, a GO enrichment analysis was conducted on the gene set associated with MBNL3 [25]. The results regarding GO biological processes (Figure 4) show that the associated genes are primarily enriched in categories such as spliceosome-mediated mRNA splicing regulation, alternative mRNA splicing regulation, RNA splicing regulation, protein-RNA complex assembly, spliceosome complex assembly, mRNA splice site recognition, RNA processing, and the regulation of mRNA stability. These categories are consistent in direction, centering on RNA splicing and RNA processing.

Heterotopic ossification involves the recruitment and proliferation of precursor cells, as well as abnormal osteogenic differentiation. Alterations in RNA splicing regulation can modify the expression patterns of genes related to the maintenance of cellular states, inflammatory responses, or osteogenic differentiation.

The results for GO molecular functions are shown in Figure 5a; the enrichment patterns support the conclusions regarding RNA-related regulation. Enriched terms are concentrated on RNA-binding functions, including mRNA binding and pre-mRNA binding, with mRNA 3' UTR binding also being prominent. Terms related to miRNA binding and regulatory RNA binding also appear; collectively, these point toward RNA-binding activity. Results regarding Gene Ontology (GO) cellular components are shown in Figure 5b. The associated genes primarily localize to cytoplasmic stress granules and the nucleus. Cytoplasmic stress granules are frequently involved in mRNA storage, degradation, and translational regulation, while the nucleus is associated with RNA transcription, splicing, and processing. These results further support the involvement of MBNL3 in RNA processing and post-transcriptional regulation.

KEGG analysis yielded no pathways reaching statistical significance (Figure 6). This result does not imply that MBNL3 lacks biological significance; rather, it indicates that its associated network does not clearly cluster within any specific canonical signaling pathway. When considered alongside the GO results, MBNL3 appears to function primarily in RNA splicing, RNA processing, and post-transcriptional regulation. Consequently, subsequent discussions should avoid forcing an interpretation of MBNL3 as a core molecule of a canonical signaling pathway.

Bioinformatics analysis of MYBL1

Functional annotation of the MYBL1 gene

MYBL1, also known as A-MYB (MYB proto-oncogene-like 1), belongs to the MYB transcription factor family. Database annotations indicate that the protein encoded by MYBL1 is primarily localized in the nucleus and functions as a DNA-binding transcription factor. It is also involved in RNA polymerase II-dependent transcriptional activation [26]. As a transcription factor, MYBL1 operates at the level of gene transcription and likely alters downstream gene expression; it does not primarily serve structural maintenance functions, nor does it act as a signal-receiving membrane receptor.

Functional annotations and findings from the literature point in the same direction: MYBL1 functions primarily as a nuclear transcriptional regulator. Reactome analysis revealed terms related to the regulation of gene expression, with the DREAM complex also showing enrichment. These results suggest that MYBL1 may participate in cell cycle-dependent gene expression regulation; during the development of heterotopic ossification, it likely influences the proliferative capacity of progenitor cells and alters their differentiation status and osteogenic propensity. This inference also helps explain the subsequent pathway enrichment results.

Although both MYBL1 and MBNL3 are identified as intersection genes, their primary functional roles differ. Functional clues for MBNL3 point toward RNA processing and the regulation of alternative splicing, suggesting that MBNL3 is more closely associated with post-transcriptional regulation. In contrast, MYBL1 is more closely linked to the regulation of gene transcription. If MYBL1 is involved in heterotopic ossification, transcriptional regulation is likely the primary mechanism; it may affect genes related to osteogenic differentiation as well as those involved in cell proliferation. Cell fate determination and local pathological responses represent potential areas for future investigation.

Currently, there is a lack of direct experimental research linking MYBL1 to heterotopic ossification. Therefore, the assessment of MYBL1 in this section is primarily based on database annotations, PPI network analysis, GO/KEGG enrichment, and Reactome results.

Analysis of the MYBL1-related protein-protein interaction network

This study conducted a PPI analysis centered on MYBL1. The STRING database was used to construct the relevant protein-protein interaction network, with the interaction confidence set to "medium." The network diagram utilizes the visualization results generated directly by STRING. As shown in Figure 7, the network comprises 11 nodes and 22 interaction edges; the average node degree is 4, the expected number of edges is 11, and the PPI enrichment p-value is 0.00028. The actual number of edges exceeds the database's expected value. The average node degree reflects the connectivity between proteins in the network, with higher degrees indicating closer functional associations. This suggests a degree of concentration in network connectivity, implying that the associated proteins may belong to the same functional module.

In terms of node composition, MYBL1 is closely associated with proteins such as LIN9, LIN37, LIN52, LIN54, RBBP4, CREBBP, NCL, BCL7A, LMO2, and NFIB. Notably, LIN9, LIN37, LIN52, LIN54, and RBBP4 are closely linked to modules associated with the DREAM complex/MuvB [27].

GO enrichment analysis of MYBL1

The GO results align with the transcriptional regulatory characteristics of MYBL1 (Figure 8). Genes associated with MYBL1 largely point to nuclear regulation, with categories related to transcriptional regulation being particularly prominent. These include RNA polymerase II-mediated transcriptional regulation, as well as the regulation and negative regulation of DNA-templated transcription; positive regulation of RNA biosynthesis further indicates a role in expression control. Categories such as the negative regulation of developmental processes, regulation of stem cell population maintenance, and chromatin remodeling suggest potential involvement in changes to cell fate.

These GO terms can be understood at three levels. The first is transcriptional regulation, indicating that the MYBL1-associated network is linked to the switching of gene expression. The second involves developmental differentiation and stem cell population maintenance, suggesting an influence on cell state. The third is chromatin remodeling, implying that the associated effects extend beyond the expression of individual genes to potentially influence the broader transcriptional environment. In heterotopic ossification, the differentiation trajectory of progenitor cells is closely tied to the activation of osteogenesis-related genes; thus, the findings regarding transcriptional regulation offer clues as to how MYBL1 influences cell differentiation pathways. Consequently, MYBL1 is best analyzed as a node within transcriptional regulation.

Pathway enrichment analysis of MYBL1

The KEGG pathway enrichment results are shown in Figure 9; MYBL1-associated genes are implicated in processes such as cellular senescence, the Notch signaling pathway, adherens junctions, and chromatin remodeling. While these pathways do not correspond in a simple one-to-one manner to heterotopic ossification, they are all linked to changes in cellular state, cell fate determination, or the regulation of the local microenvironment. Notably, the Notch signaling pathway is closely associated with cell differentiation [28], whereas the cellular senescence pathway influences cell cycle arrest and SASP-related microenvironmental changes [29].

In addition to KEGG analysis, supplementary analyses were conducted using STRING local network clustering and Reactome pathway analysis. The results (Figure 10) indicate that the MYBL1-associated network primarily involves the Myb complex regulatory module. Events associated with the G1 phase and Polo-like kinase mediation also appear in the results. Transcriptional regulation of E2F target genes mediated by the DREAM complex plays a critical role; processes involving G1/S-specific transcription and the G1/S transition are also included. The DREAM complex regulates cell cycle-dependent gene expression, targeting E2F target genes, among others, and influences transcriptional programs during the G1/S and G2/M phases [27]. Collectively, these findings suggest that MYBL1 is likely involved in cell cycle-related transcriptional regulation.

When these results are contextualized within heterotopic ossification, a clearer hypothesis emerges. Heterotopic ossification entails the proliferation of progenitor cells, alterations in cellular state, and aberrant osteogenic differentiation. The MYBL1-associated network encompasses transcriptional regulation, the cell cycle, stem cell maintenance, chromatin remodeling, and the Notch signaling pathway. To further evaluate the relationship between these findings and established regulatory mechanisms of heterotopic ossification, the MYBL1 pathway analysis results were compared with known pathways associated with the condition, pathways that have been the focus of previous research. A key distinction to be made is whether these pathways represent direct osteogenic signaling or mechanisms regulating cellular state.

The results regarding MYBL1 do not directly point to major osteogenic pathways like BMP/Smad; instead, functional clues suggest a role in regulating cellular states. Both the Notch signaling pathway and cellular senescence appeared in the results. The involvement of the DREAM/E2F transcriptional program and G1/S phase regulation points toward cell cycle-related processes. This suggests that MYBL1 likely does not participate directly in heterotopic bone formation as a typical osteogenic signaling molecule; rather, it may indirectly contribute to the pathological processes of heterotopic ossification by regulating nuclear transcriptional programs, thereby influencing cell cycle progression, cell fate decisions, and changes in the local microenvironment.

Visualization of MYBL1-related pathway analysis

To interpret the pathway results associated with MYBL1, the cellular senescence pathway was selected for visualization in this section. The cellular senescence pathway stood out prominently in the KEGG results (Figure 9). This pathway encompasses more than just cellular aging; it also involves cell cycle arrest, alterations in cellular states, and the regulation of the microenvironment. This explains the cell cycle regulation clues observed in the MYBL1 results. Changes in transcriptional programs also align with this pathway.

The RB-E2F module is more closely associated with the regulation of transcriptional programs, whereas the TGFβ/SMAD module links to the tissue microenvironment. SASP-the senescence-associated secretory phenotype-is another key phenomenon in cellular senescence. SASP can influence the release of inflammatory factors and alter the levels of cytokines and growth factors, thereby affecting the local microenvironment [30]. The MYBL1-associated network points to the DREAM/E2F-mediated cell cycle transcriptional program; this program may alter cell states and potentially affect the local microenvironment. Thus, it is hypothesized that MYBL1 may indirectly contribute to the formation of an abnormal osteogenic environment by influencing cell cycle progression and cell states.

Prediction of targeted binding between miR-155 and key target genes

Prediction of targeted binding between miR-155 and MBNL3

To determine whether miR-155 might directly target MBNL3, this study utilized the STarMir platform to analyze the interaction between hsa-miR-155-5p and the MBNL3 3'UTR [31]. MBNL3 transcript information was obtained from NCBI, and the 3'UTR sequence was compiled based on this data. Subsequently, this sequence and hsa-miR-155-5p were input into STarMir, and the appropriate model was selected to predict target sites. STarMir identified two candidate binding regions within the MBNL3 3'UTR. The first candidate site, spanning regions 157–201, exhibited a logistic probability of 0.681, an 8-mer seed type, and a ΔGhybrid of -15.5 kcal/mol. The second candidate site, spanning regions 891–915, showed a logistic probability of 0.584, an offset-6mer seed type, and a ΔGhybrid of -17.1 kcal/mol.

The advantage of the 891–915 region lies in its lower free energy; the ΔGhybrid of -17.1 kcal/mol suggests a relatively stable hybridization structure. However, this site also has drawbacks: a lower logistic probability and an offset-6mer seed type, indicating a less typical seed region matching pattern. The 157–201 region warrants greater attention; with a logistic probability of 0.681, higher than that of the 891–915 region, and an 8-mer seed type, both the prediction probability and the seed region match support this choice. Conformational prediction plots indicate that miR-155 can form local complementary base pairs with the MBNL3 3'UTR within the 157–201 region. These results support the identification of MBNL3 as a potential target gene of miR-155.

Previous functional analysis indicated that MBNL3 is primarily involved in RNA splicing, RNA processing, and post-transcriptional regulation. Taken together, the results of this study suggest that the miR-155–MBNL3 axis may constitute a post-transcriptional regulatory pathway that helps explain the potential mechanism by which miR-155 contributes to heterotopic ossification; this pathway warrants further experimental validation.

Prediction of miR-155 targeting of MYBL1

The same methodology was applied to predict the targeting of MYBL1. The analysis involved extracting the MYBL1 3'UTR sequence and inputting both hsa-miR-155-5p and the target sequence into STarMir. The results revealed two potential binding sites within the MYBL1 3'UTR region: a high-confidence site located at positions 2093–2115 (logistic probability: 0.797; seed type: 8-mer; ΔGhybrid: -15.2 kcal/mol) and another site located at positions 577–594 (logistic probability: 0.497; seed type: offset-6mer; ΔGhybrid: -11.6 kcal/mol).

Two candidate binding regions were identified for MYBL1, with distinct differences in their parameters. The region at 2093–2115 exhibited a higher logistic probability, indicating greater prediction confidence; its 8-mer seed type represents a canonical miRNA seed-region match. While ΔGhybrid alone cannot determine a targeting relationship, the relatively low free energy of the 2093–2115 region suggests a potentially stable hybridization structure; thus, the three metrics must be evaluated collectively. Given that the prediction probability, seed-region match, and hybridization free energy all favor the 2093–2115 region, this site is the most suitable candidate for further investigation.

The binding conformation model shows that miR-155 can form local complementary base pairs with the MYBL1 3'UTR in the 2093–2115 region; however, this result remains a sequence-based prediction. It does not directly prove that miR-155 regulates MYBL1 intracellularly, but it suggests a potential targeting relationship between the two. The preceding GO, KEGG, and Reactome analyses point in a similar direction; genes associated with MYBL1 are primarily involved in transcriptional regulation, with cell cycle programs and the maintenance of cellular states also appearing in the results. If miR-155 modulates MYBL1 expression, the downstream effects would likely manifest in transcriptional regulation, cell cycle dynamics, and the maintenance of cellular states.

Should this targeting relationship be experimentally validated, MYBL1 could emerge as a key candidate molecule linking miR-155 to cell cycle transcriptional programs. Given that heterotopic ossification involves progenitor cell proliferation, differentiation, and alterations in the microenvironment, the miR-155–MYBL1 axis warrants further investigation.

Summary

In this section, STarMir was used to predict potential binding interactions between miR-155 and the 3'UTRs of MBNL3 and MYBL1. Predicted binding regions for miR-155 were identified within the 3'UTRs of both MBNL3 and MYBL1. Specifically, the 157–201 region of the MBNL3 3'UTR stands out as a key site, while the 2093–2115 region of the MYBL1 3'UTR also merits attention in future studies. Previous functional analyses suggest that MBNL3 and MYBL1 operate at distinct regulatory levels: MBNL3 is primarily associated with RNA processing, alternative splicing, and post-transcriptional regulation, potentially influencing the stability of gene expression and the state of progenitor cells during heterotopic ossification. In contrast, MYBL1 is implicated in the regulation of gene transcription and cell cycle dynamics, suggesting it may influence progenitor cell proliferation and the readiness for differentiation.

Discussion

Prior work has linked miR-155 to inflammation and osteogenesis, yet its HO-specific downstream cascade remains uncharacterized. Integrating multi-database miRNA prediction and HO transcriptomic data, we identified MBNL3 and MYBL1 as two novel miR-155 targets, mediating ectopic ossification via post-transcriptional splicing and transcriptional cell cycle regulation, respectively.

MBNL3 carries a conserved miR-155 binding site and governs RNA splicing rather than canonical BMP/TGF-β osteogenic signaling. Existing evidence confirms that MBNL3 stabilizes mesenchymal cell identity. We propose that injury-upregulated miR-155 inhibits MBNL3, disrupts splicing homeostasis, and removes barriers to aberrant osteogenic trans differentiation, offering a splicing-centered mechanism supplementing classic signaling models.

As a DREAM complex transcription factor, MYBL1 restrains progenitor proliferation via E2F-mediated cell cycle control. miR-155-mediated MYBL1 suppression loosens this checkpoint, expanding osteogenic progenitor populations, while crosstalk with senescence-related inflammation further accelerates HO progression. Unlike direct osteogenic regulators, MYBL1 shapes HO through transcriptional fate control.

Classical miR-155 targets SMAD5 and SOCS1 were not recovered here, demonstrating disease-specific miR-155 regulatory networks in HO. STarMir prediction validated direct miR-155 binding to the 3′UTR of both genes, supporting our dual regulatory axis model.

This study is limited to single-dataset bioinformatic analysis without experimental validation. Future dual-luciferase assays and cell functional tests will verify the miRNA-target interaction and their osteogenic regulatory roles, and multi-clinical cohorts will validate gene expression changes in HO patients.

In short, this study identifies two-layered miR-155 mediators controlling splicing and cell cycle, providing new molecular clues for HO mechanistic and therapeutic research.

Conclusion

Although MBNL3 and MYBL1 are both intersection genes, their functional directions differ. miR-155 may be involved in RNA splicing, RNA processing, and post-transcriptional regulation via MBNL3, a direction that can be used to explain the pathological changes in ectopic ossification. Therefore, MBNL3 appears to be more closely related to RNA-level regulation. The MYBL1 results point to the Notch signaling pathway, cellular senescence, and the DREAM/E2F cell cycle transcription program. These processes cannot be directly equated with classic osteogenic pathways such as BMP/Smad, but they can explain precursor cell differentiation, cell cycle changes, and alterations in the local microenvironment during ectopic ossification. This difference suggests that miR-155 may have multi-layered regulatory clues in ectopic ossification. MBNL3 corresponds to post-transcriptional regulation, while MYBL1 points to gene transcription and cell state regulation.

The STarMir targeting prediction results corroborate the above judgments. miR-155 has predicted binding regions in the 3'UTR of both MBNL3 and MYBL1, and both have appeared in the ectopic ossification intersection genes. The targeting prediction results further suggest that they may be direct regulatory targets of miR-155. Combining these three types of results, the regulatory clues of miR-155 can be divided into two levels. MBNL3 corresponds to post-transcriptional regulation, focusing on RNA processing and splicing. MYBL1 corresponds to transcriptional regulation and cell cycle changes. Together, they constitute a candidate regulatory network for miR-155's involvement in heterotopic ossification.

Conflict of Interest

The authors declare that there is no conflict of interest regarding the publication of this paper.

References

1. Meyers C, Lisiecki J, Miller S, Levin A, Fayad L, et al. (2019) Heterotopic Ossification: A Comprehensive Review. JBMR Plus 3(4): e10172.
2. Cao G, Zhang S, Wang Y, Quan S, Yue C, et al. (2023) Pathogenesis of acquired heterotopic ossification: risk factors, cellular mechanisms, and therapeutic implications. Bone 168: 116655.
3. Huang Y, Wang X, Zhou D, Zhou W, Dai F, Lin H (2021) Macrophages in heterotopic ossification: from mechanisms to therapy. NPJ Regen Med 6(1): 70.
4. Wu M, Wu S, Chen W, Li YP. (2024) The roles and regulatory mechanisms of TGF-β and BMP signaling in bone and cartilage development, homeostasis, and disease. Cell Res 34(2):101-123.
5. Shen Q, Little SC, Xu M, Haupt J, Ast C, et al. (2009) The fibrodysplasia ossificans progressiva R206H ACVR1 mutation activates BMP-independent chondrogenesis and zebrafish embryo ventralization. J Clin Invest 119(11):3462-3472.
6. Kan C, Chen L, Hu Y, Ding N, Lu H, et al. (2018) Conserved signaling pathways underlying heterotopic ossification. Bone 109: 43-48.
7. Pulik Ł, Mierzejewski B, Sibilska A, Grabowska I, Ciemerych MA, et al. (2022) The role of miRNA and lncRNA in heterotopic ossification pathogenesis. Stem Cell Res Ther 15 13(1): 523.
8. O'Connell RM, Rao DS, Baltimore D (2012) microRNA regulation of inflammatory responses. Annu Rev Immunol 30: 295-312.
9. Bartel DP. (2009) MicroRNAs: target recognition and regulatory functions. Cell 136(2): 215-233.
10. Gu Y, Ma L, Song L, Li X, Chen D, Bai X (2017) miR-155 Inhibits Mouse Osteoblast Differentiation by Suppressing SMAD5 Expression. Biomed Res Int 2017: 1893520.
11. Zhang J, Zhao H, Chen J, Xia B, Jin Y, et al. (2012) Interferon-β-induced miR-155 inhibits osteoclast differentiation by targeting SOCS1 and MITF. FEBS Lett 586(19): 3255-3262.
12. D'Adamo S, Alvarez-Garcia O, Muramatsu Y, Flamigni F, Lotz MK. (2016) MicroRNA-155 suppresses autophagy in chondrocytes by modulating expression of autophagy proteins. Osteoarthritis Cartilage 24(6): 1082-1091.
13. Blüml S, Bonelli M, Niederreiter B, Puchner A, Mayr G, et al. (2011) Essential role of microRNA-155 in the pathogenesis of autoimmune arthritis in mice. Arthritis Rheum 63(5): 1281-1288.
14. Clough E, Barrett T, Wilhite SE, Ledoux P, Evangelista C, et al. (2024) NCBI GEO: archive for gene expression and epigenomics data sets: 23-year update. Nucleic Acids Research 52(D1): D138-D144.
15. Tang D, Chen M, Huang X, Zhang G, Zeng L, et al. (2023) SRplot: A free online platform for data visualization and graphing[J]. PLoS ONE, 18(11): e0294236.
16. Kanehisa M, Furumichi M, Sato Y, Kawashima M, Ishiguro-Watanabe M. (2023) KEGG for taxonomy-based analysis of pathways and genomes. Nucleic Acids Research 51(D1): D587-D592.
17. Milacic M, Beavers D, Conley P, Gong C, Gillespie M, et al. (2024) The Reactome Pathway Knowledgebase. Nucleic Acids Research 52(D1): D672-D678.
18. Agarwal V, Bell GW, Nam JW, Bartel DP. (2015) Predicting effective microRNA target sites in mammalian mRNAs. Elife 4: e05005.
19. Chen Y, Wang X. (2020) miRDB: an online database for prediction of functional microRNA targets. Nucleic Acids Res 48(D1): D127-D131.
20. Huang HY, Lin YC, Li J, Huang KY, Shrestha S, et al. (2020) miRTar Base 2020: updates to the experimentally validated microRNA-target interaction database. Nucleic Acids Res 48(D1): D148-D154.
21. Pascual M, Vicente M, Monferrer L, Artero R. (2006) The Muscle blind family of proteins: an emerging class of regulators of developmentally programmed alternative splicing. Differentiation 74(2-3): 65-80.
22. Squillace RM, Chenault DM, Wang EH. (2002) Inhibition of muscle differentiation by the novel muscle blind-related protein CHCR. Dev Biol 250(1): 218-230.
23. Lee KS, Cao Y, Witwicka HE, Tom S, Tapscott SJ, Wang EH. (2010) RNA-binding protein Muscle blind-like 3 (MBNL3) disrupts myocyte enhancer factor 2 (Mef2) {beta}-exon splicing. J Biol Chem 285(44): 33779-33787.
24. Szklarczyk D, Kirsch R, Koutrouli M, Nastou K, Mehryary F, et al. (2023) The STRING database in 2023: protein-protein association networks and functional enrichment analyses for any sequenced genome of interest. Nucleic Acids Res 51(D1): D638-D646.
25. Aleksander S A, Balhoff J, Carbon S, Cherry J M, Drabkin H J, et al. (2023) The Gene Ontology knowledgebase in 2023. Genetics 224(1): iyad031.
26. Bolcun-Filas E, Bannister LA, Barash A, Schimenti KJ, Hartford SA, et al. (2011) A-MYB (MYBL1) transcription factor is a master regulator of male meiosis. Development 138(15): 3319-3330.
27. Sadasivam S, DeCaprio JA. (2013) The DREAM complex: master coordinator of cell cycle-dependent gene expression. Nat Rev Cancer 13(8): 585-595.
28. Thomas S, Jaganathan BG. (2022) Signaling network regulating osteogenesis in mesenchymal stem cells. J Cell Commun Signal 16(1): 47-61.
29. Zhang Q, Zhou D, Liang Y. (2022) Single-cell analyses of heterotopic ossification: characteristics of injury-related senescent fibroblasts. Journal of Inflammation Research 15: 5579-5593.
30. Matsuda S, Revandkar A, Dubash T D, Ravi A, Wittner B S, et al. (2023) TGF-β in the microenvironment induces a physiologically occurring immune-suppressive senescent state. Cell Reports 42(3): 112129.
31. Rennie W, Liu C, Carmack C S, Wolenc A, Kanoria S, et al. (2014) STarMir: a web server for prediction of microRNA binding sites. Nucleic Acids Research 42(W1): W114-W118.