This article is mentioned in:

Introduction

Tendon injury is a frequent problem for recreational and competitive athletes. Individuals, who live sedentary lifestyles may also develop tendinopathy in the absence of any history of increased physical activity (1). An estimated 30–50% of all sports-associated injuries are caused by a disorder of the tendons (2). This injury may be caused by intrinsic or extrinsic factors, either alone or in combination. Extrinsic factors lead to the majority of acute tendon injuries, although in overuse syndromes, including tendinopathy, multifactorial combinations of intrinsic factors, such as age-associated cell activity changes and extrinsic factors, including overuse, repetitive strain injury and microtrauma may be the cause (3).

Previous studies have indicated that histopathological changes occur with tendinopathy and are associated with degeneration and disorganization of collagen fibers, increased cellularity and minimal inflammation (4). Macroscopic changes include thickening of the tendon, loss of mechanical properties and pain. Previous studies have demonstrated that several changes occur in response to overuse, including the production of matrix metalloproteinases, cytokines, tendon cell apoptosis, chondroid metaplasia of the tendon, collagen, glycosaminoglycan and expression of protective factors (5,6).

Currently, the non-surgical therapies available to patients who suffer from tendinopathies are exercise-based physical therapy, ultrasound, non-steroidal anti-inflammatory drugs and steroid or platelet-rich plasma injections (7). However, these therapies offer symptomatic relief, but do not result in definitive disease resolution. Through understanding the cellular and molecular mechanisms of causation and novel therapeutic targets, small molecules could potentially be identified for drug development. This may result in the development of more effective treatments, while minimizing side effects.

Microarray analysis supports the identification of drug-sensitive genes and the chemical substructures associated with specific genetic responses. It has become a powerful tool in drug development (8). In the present study, microarrays were utilized to identify differentially expressed genes (DEGs) between normal and degenerating tendon cells. The functions of DEGs were investigated by annotating to biological processes and pathways. Several target sites of the transcription factors and certain regulatory microRNAs were also screened. This information may assist in elucidating the molecular mechanism of tendon injuries. In addition, candidate small molecules were identified for the potential treatment of tendinopathy.

Materials and methods

Derivation of genetic data

The gene expression profile of GSE26051 (7) was downloaded from a public functional genomics data repository, the Gene Expression Omnibus (GEO; www.ncbi.nlm.nih.gov/geo/) database. A total of 46 specimens, including 23 normal samples and 23 tendinopathy specimens, were available based on the GPL570 platform. This information was approved by the ethics committee of the Hospital for Special Surgery (New York City, NY, USA).

DEG analysis

The derived genetic data was analyzed using the GEOquery (www.bioconductor.org/packages/release/bioc/html/GEOquery.html) and Limma (www.bioconductor.org/packages/release/bioc/html/limma.html) packages in the R programming language (v.2.13.0) (9). Geoquery can quickly access the expression profiling data on the GEO database, while Limma is the most popular method of statistical analysis to analyze the DEGs (10). The preprocessed microarray data were obtained by Geoquery package and then a log2 transformation was performed. The Limma package, a linear regression model, was applied to compare the normal samples and tendinopathy samples. Only the genes with P<0.05 were identified as DEGs.

Gene Ontology (GO) enrichment analysis

GO analysis has become a common approach for the functional annotation of large-scale genomic data (11). Gene ontology enrichment analysis software toolkit (GOEAST; omicslab.genetics.ac.cn/GOEAST/) is an easy-to-use web-based toolkit, which identifies statistically overrepresented GO terms within provided gene sets (12). GOEAST was utilized for GO enrichment analysis to identify the locations of DEGs within cellular compartments and molecular functions affected by DEGs, based on the hypergeometric distribution, with the false discovery rate (FDR) <0.001.

Biological pathway enrichment analysis

Biological pathways were investigated to examine the tendinopathy cell changes at the molecular level. All metabolic and non-metabolic pathways were downloaded from the open WikiPathways database (www.wikipathways.org.) (13,14) and WikiPathways cluster analysis was conducted (15,16) to the DEGs using the gene set analysis toolkit V2 platform. A count number >2 and P<0.05 were selected as the cut-off criteria.

Examining potential target sites of transcription factors and potential regulatory microRNAs

Well-annotated gene sets in the molecular signature database (MsigDB; www.broadinstitute.org/gsea/msigdb/index.jsp) were subject to gene set enrichment analysis (GSEA) (17). Subsequently, the GSEA results were statistically accounted for with the hypergeometric distribution. The consequences were adjusted for multiple testing using the Benjamini-Hochberg procedure. Finally, the target sites with FDR<0.01 were selected as the potential target sites that may regulate transcription factors. Similarly, the potential regulatory microRNAs were identified with an FDR<0.05.

Identification of candidate small molecules

The connectivity map (CMap) database contains data on 7,056 gene-expression profiles, involving 6,100 small molecule treatment-control pairs (18). The DEGs were divided into up- and downregulated groups. Subsequently, these genes were subjected to GSEA and compared with the DEGs in the CMap database. Finally, a correlation score for each perturbagen was calculated, ranging between −1 and +1 (19).

Results

DEG selection

In order to analyze differentially expressed genes between cells in tendinopathy and normal controls, a publicly available microarray dataset, GSE26051 was obtained and a classical t-test, corrected for multiple comparisons was performed. A total of 419 probes were considered to be differentially expressed in tendinopathy samples when compared with normal control tendons (P<0.001), which corresponded to 318 DEGs.

GO enrichment analysis of DEGs

To investigate the functional changes in the pathological process of tendinopathy, the DEGs were mapped to the GO database. The project provided three structured networks of defined terms to describe gene product attributes: biological process, molecular function and cellular compartment. Fig. 1 reveals the molecular function in which the majority of the DEGs were located, such as the cytoskeleton, actin cytoskeleton and sarcoplasm. In addition, Fig. 2 shows the biological processes of DEGs, for instance, protein complex binding and cytoskeletal protein binding. The majority of enriched GO biological processes of the DEGs between normal and pathological specimens were associated with a particular cellular compartment (Fig. 3), for example, multicellular organismal processes, developmental processes and single-multicellular organism processes.

Figure 1 Enriched GO terms of the cellular compartment of differentially expressed genes. White entries indicate nonsignificant aggregation (false discovery rate >0.05) and the grey entries indicate significant aggregation (false discovery rate <0.05). GO, gene ontology.

Figure 2 Enriched GO terms of the molecular function of the differentially expressed genes. The white entries indicate nonsignificant aggregation (false discovery rate >0.05) and the grey entries indicate significant aggregation (false discovery rate <0.05). GO, gene ontology.

Figure 3 Enriched GO terms of the molecular function of the differentially expressed genes. The white entries indicate nonsignificant aggregation (false discovery rate >0.05) and the grey entries indicate significant aggregation (false discovery rate <0.05). GO, gene ontology.

Pathway enrichment analysis

To gain further insights into the changes of biological pathways in cells of tendinopathy, the WikiPathways cluster analysis was used to identify the significant pathways associated with DEGs. P<0.05 and counts, which were >2 were selected as the cut-off criteria. A total of 10 pathways were identified and the main 8 pathways with a highly significant correlation are listed in Table I. The most significant pathway was focal adhesion with P=7.08E-5 and the genes enriched in focal adhesion were laminin, α4 (LAMA4), platelet-derived growth factor α (PDGFA), laminin γ1 (LAMC1) and Src homology 2 transforming protein 1 (SHC1).

Table I Enriched biological pathways (P<0.05). The main eight pathways are listed.

Table I

Enriched biological pathways (P<0.05). The main eight pathways are listed.

Pathway	Genes	P-value
Focal adhesion	LAMA4, PDGFA, LAMC1, SHC1	7.50 ×10−0.5
Integrin-mediated cell adhesion	ITGA5, SORBS1, ILK, SHC1	0.031
Myometrial relaxation and contraction pathways	YWHAZ, CRHR1, GRK5, PRKCG	0.074
Insulin signaling	PRKAA2, SORBS1, SGK2, CAP1	0.085
Osteoclast	TNFSF11, CTSK	0.0171
Serotonin receptor 2-> ELK-SRF/GATA4 signaling	HTR2A, ITPR1	0.0248
Delta-notch signaling pathway	YWHAZ, LAMC1, SHC1	0.0342
L-3 Signaling pathway	YWHAZ, CHEK1, SHC1,	0.0582

Examining potential target sites

As an important regulatory element, transcription factors can regulate gene expression. Taking upstream sequences of the DEGs as the analyzed object, the potential target sites of the transcription factor were examined. The main 20 target sites with a highly significantly correlation are listed in Table II. The most significant transcription factors were lymphoid enhancer-binding factor 1 (LEF1) and OCT1, in which LEF1 may regulate the pantothenate kinase 2 (PANK2) and G protein-coupled receptor kinase 5 (GRK5) by binding the target sequence CTTTGT.

Table II Enriched potential target sites of transcription factors.

Table II

Enriched potential target sites of transcription factors.

Target	Genes	P-value
hsa_CTTTGT_V$LEF1_Q2	PANK2, SEMA6D, GRK5, ATOH8	2.58 ×10−0.5
hsa_V$OCT1_06	CUGBP2, AGRP, LMO1, ANK3	2.58 ×10−0.5
hsa_V$PXR_Q2	PITX2, HOXD3, SLC6A6, LMO1	0.0001
hsa_AACTTT_UNKNOWN	PRKAA2, PFKFB1, SORBS1, TRDN	0.0001
hsa_CAGGTG_V$E12_Q6	EMX1, SORBS1, LMO1, CDKN1C	0.0001
hsa_TGACAGNY_V$MEISI_01	ELOVL5, DPF3, DNAJA4, ESRRG	0.0001
hsa_V$NKX25_02	PITX2, FOXP2, DPF3, LRP1B	0.0004
hsa_CATTGTYY_V$SOX9_B1	FOXP2, DPYSL5, HOXD3, GRK5	0.0004
hsa_V$MEISI_01	FOXP2, DPF3, SORBS1, EPHB2	0.0005
hsa_V$OCT1_02	FOXP2, PFKFB1, DNAJA4, TRDN	0.0005
hsa_V$POU6F1_01	PITX2, FOXP2, CUGBP2, EPHB2	0.0005
hsa_SMTTTTGT_UNKNOWN	DPYSL5, DPF3, LMO1, RCOR1	0.0006
hsa_V$FREAC3_01	FOXP2, MMP11, ATOH8, LEF1	0.0006
hsa_TTGTTT_VFOXO4_01	SORBS1, GRK5, ATOH8, MPZL1	0.0013
hsa_V$SMAD4_Q6	ELOVL5, MMP11, LMO1, CBFA2T3	0.0013
hsa_V$AP1_Q6_01	LRP1B, EPHB2, LAMC1, ATP6V1A	0.0013
hsa_V$WHN_B	FOXP2, LRP1B, LEF1, ARPC5	0.002
hsa_AAANWWTGC_UNKNOWN	FOXP2, DPYSL5, EPHB2, BNC2	0.0028
hsa_V$PITX2_Q2	PITX2, FOXP2, PDGFA, DPYSL5	0.0028
hsa_V$PTF1BETA_Q6	FOXP2, CRHR1, CALD1, EEF1A2	0.0038
hsa_YATGNWAAT_V$OCT_C	PITX2, PFKFB1, SORBS1, LLGL2	0.0038
hsa_TATAAA_V$TATA_01	PRKAA2, HOXD3, DNAJA4, SLC2A3	0.0038
hsa_V$EN1_01	FOXP2, NRG1, HOXD3, LMO1	0.0038
hsa_GCCNNNWTAAR_UNKNOWN	FOXP2, DPF3, ITPR1, ETV6	0.0038
hsa_V$CHX10_01	PITX2, FOXP2, LMO1, BACE2	0.0038
hsa_V$SMAD3_Q6	PFKFB1, CHDH, CBFA2T3, ESRRG	0.0038
hsa_V$SRF_C	FOXP2, DUSP2, CAP1, CALD1	0.0038
hsa_GCANCTGNY_V$MYOD_Q6	DPYSL5, HOXD3, DPF3, EMX1	0.0038
hsa_V$HP1SITEFACTOR_Q6	PITX2, FOXP2, HOXD3, ESRRG	0.0048
hsa_WGTTNNNNNAAA_UNKNOWN	CKC25A, DPF3, EMX1, ATOH8	0.0048
hsa_V$TBP_01	PRKAA2, DNAJA4, TRDN, ESRRG	0.0048
hsa_TGACATY_UNKNOWN	ELOVL5, MICAL2, DPF3, SORBS1	0.0048
hsa_V$ATF1_Q6	PDGFA, LEOVL5, ESRRG, CALD1	0.0048
hsa_RTAAACA_V$FREAC2_01	GRK5, CDKN1C, ESRRG, BNC2	0.0055
hsa_V$DR4_Q2	ARHGAP24, FOXP2, LAMA4	0.0055
hsa_TAATTA_V$CHX10_01	CUGBP2, LMO1, LAMC1, TRDN	0.0055
hsa_V$GATA6_01	PITX2, PFKFB1, SORBS1, COL4A3	0.0055
hsa_GGGAGGRR_V$MAZ_Q6	SORBS1, GRK5, MAPRE3, MTNR1B	0.0055
hsa_V$MF2_Q6_01	FOXP2, PRKAA2, DNAJA4, MYOG	0.0055
hsa_VTST1_01	PITX2, PDGFA, LAMA4, HOXD3	0.0065
hsa_GGGTGGRR_V$PAX4_03	DPYSL5, PFKFB1, ITPKB, ATOHB	0.0065
hsa_TGCCAAR_V$NF1_06	DUSP2, COL4A3, ESRRG, BNC2	0.0074
hsa_CAGCTG_V$AP4_Q5	PRKAA2, HOXD3, SORGS1, GRK5	0.0074
hsa_V$AR_Q6	PITX2, FOXP2, HOXD3, LMO1	0.0085
hsa_V$CEBP_Q2 _01	PITX2, CDKN1C, CALD1, BNC2	0.0093
hsa_V$ER_Q6	PRKAA2, ATP6V1A, ESRRG, CA4	0.0093

Examining the potential regulatory microRNA

MicroRNAs are involved in the regulation of numerous cellular processes by adjusting the stability of mRNA. The potential regulatory microRNAs were screened out based on the sequences of DEGs. The main 20 instances with a highly significant correlation were enumerated in Table III. The most significant microRNAs were in the miR-499 and miR-200 family, including miR-200B, miR-200C and miR-429. miR-499 may regulate the CUGBP2 and MYB genes by binding the target sequence AGTCTTA and the miR-200 family may regulate the LRP1B and SLC6A6 genes by binding CAGTATT.

Table III Enriched potential regulatory microRNAs.

Table III

Enriched potential regulatory microRNAs.

Target sequence	Potential microRNA	Genes	P-value
hsa_AGTCTTA	miR-499	CUGBP2, KLHDC5, FAM60A	0.0247
hsa_CAGTATT	miR-200B, miR-200C, miR-429	LRP1B, SLC6A6, LAMC1,	0.0247
hsa_GAGCCAG	miR-149	COL4A3, ACLY, RAP1B	0.0247
hsa_GTGCAAA	miR-507	SEMA6D, HECW1, LEF1	0.0296
hsa_ATACTGT	miR-144	CUGBP2, KPNA1, ESRRG	0.0296
hsa_GCAAGGA	miR-502	PANK2, HOXD3, LEF1	0.0317
hsa_TGCACTT	miR-519C, miR-519B, miR-519A	ARHGAP24, TNFSF11, WDR1	0.0317
hsa_TTGGGAG	miR-150	MMP19, NOTCH3, EPHB2	0.0444
hsa_CACTTTG	miR-520G, miR-520H	DPYSL5, TNFSF11, KPNA1	0.0444
hsa_ATAAGCT	miR-21	PITX2, ARHGAP24, CDC25A	0.049

Identification of candidate small molecules

In order to screen small molecule drugs, computational bioinformatics analysis of DEGs was performed using CMap. A total of 20 associated small molecules with a highly significant correlation are listed in Table IV, including 13 negatively-associated molecules and seven positively-associated small molecules. Among these molecules, Prestwick-1082 and Viomycin, with the highest negative correlation had the potential to treat the tendinopathy.

Table IV Enriched significant small molecules.

Table IV

Enriched significant small molecules.

Connectivity map name	Enrichment score	P-value
Propylthiouracil	0.91	0.00006
Sulfadimethoxine	−0.867	0.00008
Monensin	−0.815	0.00012
Viomycin	−0.876	0.00052
Nadolol	−0.872	0.00056
Cycloserine	−0.857	0.00056
Lisuride	−0.782	0.00088
Medrysone	0.728	0.00107
Luteolin	0.832	0.00121
Adiphenine	−0.752	0.00174
Diethylstilbestrol	−0.681	0.0028
Alpha-estradiol	0.433	0.00295
Podophyllotoxin	−0.802	0.003
Etiocholanolone	−0.678	0.003
Scopoletin	0.956	0.00344
Omeprazole	0.795	0.00346
Resveratrol	0.557	0.00348
Fuldrocortisone	−0.591	0.00357
Prestwick-1082	−0.878	0.00363
Prestwick-983	−0.874	0.00403

Discussion

Tendinopathy is a critical clinical problem as it is often asymptomatic at onset and during development, and is only recognized upon rupture of the tendon (20). Therefore, there is an urgent requirement to investigate the mechanism of tendinopathy and develop a mechanism to effectively prevent the condition or a treatment for it. In the present study, bioinformatics methods were used to investigate the molecular mechanism of tendinopathy and identify small molecule drugs, which have the potential to treat this condition. The results revealed that the expression of 318 genes were altered in the human samples of tendinopathy compared with normal tendons. These genes were mainly involved in pathways associated with adhesion. Furthermore, it was demonstrated that Prestwick-1082 and Viomycin may be effective for the treatment of tendinopathy.

The gene expression analysis, which focussed on identifying individual genes, which exhibited differences between two states, although useful, may be unable to detect biological processes, including metabolic pathways, transcriptional programs and stress responses, which are distributed across an entire network of genes and less detectable at the level of individual genes (17). Current approaches typically study entire pathways, whether through using singular enrichment analysis or by gene set enrichment analysis. In the present study, eight pathways were identified and focal adhesion was observed to be the most significant pathway in the development of tendinopathy. Focal adhesions lie at the convergence of integrin adhesion, signaling and the actin cytoskeleton (21). Genes in the integrin family, including LAMA4, PDGFA, LAMC1 and SHC1, are closely associated with focal adhesion. Among this family, the downregulation of LAMA4 may affect cell survival rate via lamin-integrin interaction (22) and LAMA4-deficient mice have previously been reported to develop a defect in endothelial cell viability, followed by cardiac hypertrophy and heart failure (23). For PDGFA, it may induce tyrosine phosphorylation of focal adhesion kinase, a member of the focal adhesion complex family. The PDGFA receptor acts as a high affinity binding site for several signaling molecules leading to activation of Ras, followed by activation of Raf, mitogen-activated protein kinase and extracellular signal-regulated kinase (24). This complex interacts with extracellular matrix proteins through integrin interactions, providing a direct sensor to the integrity and composition of the extracellular environment (25). Besides, LAMC1 belongs to the Lamins, a family of extracellular matrix glycoproteins, which are the major noncollagenous components of basement membranes. LAMC1 has been implicated in a wide variety of biological processes, including cell adhesion, differentiation, migration, signaling, neurite outgrowth and metastasis. SHCI has been reported to be involved in the aging process, a signaling pathway inducing elevation of extracellular oxidant levels, cytochrome c release and apoptosis, as well as the oxidative stress response (26). Consistent with the present findings, Riley (3) suggested that tendon matrix damage is the primary event, overwhelming the ability of the resident cell population to repair structural defects and degradation of the extracellular matrix may affect the structural properties of the tendon. Previous studies (27,28) have also reported that fibronectin is markedly increased following tendon injury when compared with the levels in the normal tendon and consequently has been implicated in cell adhesion, migration and differentiation at the site of injury. Therefore, the present study indicates that these integrin genes associated with focal adhesion have crucial roles in tendinopathy development.

Numerous studies have reported an abundance of transcription factors associated with human disease, thus making them targets for the investigation of the mechanisms of tendinopathy. In the present study, LEF1 was identified as one of the most significant transcription factors, which binds the target sequence: GTTTGT. A number of genes, including PANK2 and GRK5, which contain this sequence, can be identified by LEF1. PANK2 is a mitochondrial enzyme, which catalyzes the first regulatory step of coenzyme A synthesis and that is processed and active in the mitochondria (29). Mutations in PANK2 may lead to a variety of metabolic defects (30). Semaphorin 6D has been found to be involved in cardiac morphogenesis, cancer and immune responses (31).

MicroRNAs are small regulatory RNAs, which regulate the translation and degradation of target mRNAs and are extensively involved in human disease (32). The most significant microRNA in the present study was miR-499 and its targeting sequence was AGTCTTA. The genes, including CUGBP2 and MYB, which contained this sequence can be regulated by miR-499. CUGBP2 is an RNA-binding protein, which regulates mRNA translation and is abundant in the skeletal muscle (33). Ectopic overexpression of this protein may also induce apoptosis (34). Similar to MYB, it is an important regulator in the control of cell proliferation, apoptosis and differentiation, is highly expressed in immature, proliferating cells and is downregulated as cells become further differentiated. Tenocyte apoptosis has been observed to occur at an increased frequency in tendinopathy specimens (7).

There are several important implications of the present study. The identification of a group of small molecules with potential therapeutic efficacy for tendinopathy is an important observation. The data in Table IV show that the small molecules of Prestwick-1082 (enrichment score=-0.878) and viomycin (enrichment score=-0.876) were associated with significant negative scores, which suggest that these small molecules may be used as therapeutic drugs for tendinopathy.

Viomycin is an RNA-binding peptide antibiotic, which inhibits prokaryotic protein synthesis and group I intron self-splicing (35). It has a marked selectivity for RNAs, which form pseudoknots, a structure that may function as a ‘tag’ for recognition by this peptide and also induces interactions between RNA molecules (36). It has demonstrated promise in the search for drugs, which may be useful for treating tuberculosis (37). However, to the best of our knowledge, there are no previous studies investigating the use of these compounds as systemic therapies for tendinopathy. The present observations warrant further study and should generate hypotheses for laboratory, patient or population-based studies. The small molecule, Prestwick-1082 (enrichment score=-0.878) was associated with a significant negative score, which suggested that these small molecules are potential adjuvant drugs to improve the therapeutic effect in tendinopathy.

In conclusion, the present study has presented novel insights into the mechanism and treatment of tendinopathy. DEG profiles were analyzed using a computational bioinformatics approach. In addition, a group of small molecules were identified, which can be exploited as adjuvant drugs to improve treatment, including Prestwick-1082 and viomycin. Although it may be premature to suggest that these drugs may be ready for psychiatric clinical trials, it is clearly a direction that warrants additional consideration.

References