This article is mentioned in:

Introduction

Low back pain (LBP) is a widely spread disease, which generally requires medical care and can lead to chronic disabilities (1). It is reported that ~84% of the general population suffer from LBP in their lifetime (2). Degenerative disc disease (DDD) is a commonly observed reason for LBP (3). However, it is very challenging to investigate the pathogenesis of DDD on the account of the vague DDD definitions and multiple interdependent factors involved, including changed mechanical loading (4), hampered nutrition supply (5), hereditary factors (6) and altered extracellular matrix (ECM) composition (7). The removal of protruding disc tissues is the prevalent therapeutic method for DDD, which relieves painful symptoms of patients; however, this overlooks the underlying biological changes of discs. In terms of drawbacks of current treatments, advanced novel therapies are urgently required, which directly deal with the underlying biochemical causes of DDD to both relieve symptoms and reverse disc degeneration. In recent years, cell-based therapies regenerating disc structure and function have aroused people's interests (8). In particular, mesenchymal stem cells (MSCs) are recognized as an eligible cell source for disc targeting tissue engineering. A great number of previous studies have confirmed the ability of MSCs to self-renew, expand and of multilineage differentiation (9–12). However, the microenvironment of degenerated discs does not suit exogenous MSCs due to the pathological changes occurring in degenerated discs. Hence, it may be a better choice to concentrate on stem cells in situ in degenerated discs. Numerous previous studies have reported evidence for existence of stem cells in degenerated intervertebral discs (IVDs) (11,13,14). Notably, our research group has previous isolated cartilage endplate-derived stem cells (CESCs) and confirmed their multilineage differentiation capacity (13). The cartilage endplate (CEP) refers to a thin layer of hyaline cartilage existing between the vertebral body and the disc, and it prevents nucleus pulposus (NP) from protruding out to the adjacent vertebrae. It is speculated that CEP degeneration may be involved in the initiation and development of DDD (15,16), since CEP is the predominant path through which nutrients get to IVDs and supply them. CEP degeneration has several manifestations, including CEP calcification (17), proteoglycan loss (18) and hampered ECM synthesis (19). It is reported that calcification or sclerosis of CEP inhibited nutrient diffusion into adjacent IVDs, resulting in DDD (20). Therefore, it may be helpful to modify the differentiation ability of CESCs to relieve CEP calcification and restore CEP structure. However, the mechanism of CESC differentiation remains to be fully understood.

Alternative splicing (AS) is a sophisticated regulatory process in which a single pre-mRNA generates different RNA isoforms, potentially leading to structurally and functionally diverse proteins (21,22). It is estimated that AS occurs to ~95% of multi-exonic genes in eukaryotic organisms (23,24). Generally speaking, AS is divided into the following types: Exon skip/inclusion, mutually exclusive exons, alternative 5′/3′ splice sites, intron retention, alternative promoters and polyadenylation sites (22). In addition, regulation of AS is elaborately achieved through cell type-, development- and extracellular signal-related pathways (25). Aberrated AS of genes is reported to serve roles in numerous human diseases, including autoimmune diseases, neurodegenerative diseases and cancer (26–28). Lately, the AS mechanism underlying stem cell differentiation has piqued people's interest. Kazantseva et al (29) reported that depletion of hTAF4-TAFH domain from TAF4 isoforms caused enhanced chondrogenic differentiation of human MSCs (29). In addition, a novel alternative trasncript-quantitative polymerase chain reaction method was created by McAlinden et al (30) to quantify isoforms of alternatively spliced Col2a1 gene, and the results indicated that the majority of ATDC5 cells were the chondroprogenitor cells induced by the standard chondrogenic differentiation method. Furthermore, Longo et al (31) found that in osteogenic differentiation of human MSCs, PTHrP isoforms became increasingly selective and were considered as novel molecular markers of stem cell state. Therefore, it is very meaningful and rewarding to elucidate AS mechanism during stem cell differentiation.

As mentioned above, it may be helpful to promote chondrogenic differentiation and inhibit osteogenic differentiation of CESCs to alleviate CEP calcification and rebuild nutrition provision, repairing and regenerating disc degeneration. The present study aimed to investigate the mechanism of AS underlying osteogenic differentiation of CESCs. The isolated CESCs were induced to undergo osteogenic differentiation and a genome-wide analysis was performed on both the undifferentiated and differentiated samples using Affymetrix Human Transcriptome Array (HTA) 2.0 system. Following data extraction and pre-treatment, a comparative analysis of alternative splicing events (ASE) was performed between the controlled CESCs (undifferentiated) and osteogenically differentiated CESCs. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analyses were performed to add functional annotation of genes of interest to exemplify AS mechanisms. The present study is the first, to the best of our knowledge, to elucidate the AS mechanisms in osteogenic differentiation of stem cells on a genome-wide scale.

Materials and methods

Ethics statement

The CEP tissues used in the present study was obtained from seven patients who underwent discectomy and fusion surgeries as a result of lumbar degenerative diseases at the Department of Orthopedics, Xinqiao Hospital, Third Military Medical University (Chongqing, China) (Table I). The present study was approved by the Ethics Committee of Xinqiao Hospital, Third Military Medical University. All procedures described below were in accordance with the Helsinki Declaration. Written informed consent was obtained from each patient and we extensive precautions were made to protect the privacy of each donor.

Table I Patient information.

Table I

Patient information.

Case no.	Gender	Age (years)	Diagnosis	Degenerated disc level	Surgery type
1	Female	56	Disc herniation	L4-L5	MEDa
2	Female	67	Disc herniation	L4-L5	MED
3	Male	68	Disc herniation	L5-S1	TLIFb
4	Female	50	Spondylolisthesis	L4-L5	Quadrant assisted TLIF
5	Male	52	Disc herniation	L5-S1	MED
6	Male	58	Disc herniation	L5-S1	TLIF
7	Male	60	Disc herniation	L4-L5	TLIF

a MED, microendoscopic discectomy;

b TLIF, transforaminal lumbar interbody fusion.

Agarose culture to select CESCs

The agarose culture method was used, according to a previous study (32). Briefly, a 60 mm-diameter culture dish was coated with 1% low-melting point agarose containing an equal volume of DMEM/F12 (37°C) and 2% low-melting point agarose. Then, 0.75 ml DMEM/F12, 0.75 ml 2% low-melting point agarose and 1.5 ml DMEM/F12 (20% FBS) containing 5×104 P1 CEP-derived cells were mixed and transferred to the coated culture dishes. The final concentration of FBS was 10%. The culture dishes were maintained at 4°C for 15 min until the gel solidified. The culture dishes were subsequently incubated in a humidified atmosphere containing 5% CO2 at 37°C. The culture medium was changed twice every week. After 6 weeks, cell clusters with a diameter >50 µm were isolated by sterile Pasteur pipette and were sub-cultured in 6-well plates (Corning Inc., Corning, NY, USA). Passage 3 CESCs were used in the present study.

Osteogenic differentiation assay

The CESCs were seeded at 3×104 cells/cm2 in 6-well culture plates pre-coated with gelatin. The complete osteogenic differentiation medium (cat. no. HUXMA-90021; Cyagen Biosciences, Inc., Guangzhou, China) consisted of 175 ml basal medium, 20 ml FBS, 2 ml penicillin-streptomycin, 2 ml glutamine, 400 µl ascorbate, 2 ml β-Glycerophosphate and 20 µl dexamethasone. After CESCs reached 50–70% confluence, the culture medium was replaced with 2 ml complete medium. The complete medium was replaced every 3 days and the cells were cultured for 21 days.

Affymetrix HTA 2.0

CESCs were induced to undergo osteogenic differentiation or left untreated in the undifferentiated state. Differentiated and undifferentiated samples were treated with TRIzol reagent and sent to Bioassay Laboratory of CapitalBio Corporation (Beijing, China). The alternative splicing events were analyzed using HTA 2.0, purchased from Affymetrix (Santa Clara, CA, USA). The Affymetrix HTA 2.0 contained ~339,000 probe sets (10 probes/exon and 4 probes/junction), covering ~67,000 transcript clusters and 573,000 probe selection regions (PSRs). Transcript clusters were referred to as genes in the present study for simplicity. The HTA 2.0 allowed probes to target exons and junctions within genes and provided AS information. The labeling, hybridization, scanning and data extraction of microarray were performed by Bioassay Laboratory of CapitalBio Corporation (Beijing, China), according to the recommended Affymetrix protocols. Briefly, the fluorescence signals of the microarray were scanned and saved as DAT image files. The Affymetrix GeneChip® Command Console software transformed DAT files into CEL files, to change image signals into digital signals, which recorded the fluorescence density of probes. Next, the Affymetrix Expression Console software to pre-treat CEL files through Robust Multichip Analysis algorithm (33), including background correction, probeset signal integration and quantile normalization. Following pre-treatment, the obtained chp files were analyzed by Affymetrix Transcriptome Analysis Console software to detect alternatively spliced genes (ASGs). The Expression Console and Transcriptome Analysis Console software were provided by Affymetrix. To identify significantly enriched GO terms and functional pathways, both the publicly available web-tools KEGG (http://www.genome.jp/kegg/) and DAVID (http://david.abcc.ncifcrf.gov/tools.jsp), and the commercial database Molecule Annotation System (MAS; CapitalBio Corporation) were used. The microarray data were submitted to NCBI's Gene Expression Omnibus (accession number, GSE63897). The overall workflow of the HTA data analysis is presented in Fig. 1.

Figure 1 Overall workflow of the HTA data analysis. Briefly, CESCs were isolated and induced into osteogenic differentiation. The total RNA was extracted, labeled and hybridized to HTA 2.0. The Affymetrix AGCC, EC and TAC software were used to scan and analyze the microarray data. HTA, human transcriptome array; CESCs, cartilage endplate-derived stem cells; AGCC, Affymetrix GeneChip® command console; EC, expression console; TAC, Transcriptome analysis console.

Criteria for detecting ASGs

The Spicing Index (SI) model (34,35) was used to identify ASGs. The SI represented the ratio of the exon signal intensities normalized against the gene signal intensities between two experimental conditions, and was used to detect the exon exclusion/inclusion level. The SI value was calculated in following ways:

The NI (i,j)A represented the signal intensity of i-th exon normalized against the j-th gene in condition A. The subscript U indicated undifferentiated condition and the subscript D indicated the differentiated condition. The default filter criteria was set as SI (linear) ≤−2 or ≥2.

ASG validation by semi-quantitative reverse transcription (RT)-PCR

RT-PCR was performed to identify the ASGs. The total RNA was extracted using TRIzol reagent and was used to generate cDNA using the the relevant kits (Takara Bio, Inc., Otsu, Japan), according to the manufacturer's protocol. The quality of the total RNA was examined using a spectrophotometer (Nanodrop 2000; Thermo Fisher Scientific, Inc., Waltham, MA, USA) at 260 and 280 nm. An oligo (dT) primer was used to reverse transcribe 1 µg total RNA into cDNA using the PrimeScrip RT reagent kit with gDNA Eraser (cat. no. RR047A; Takara Bio, Inc.) and 1 µl cDNA template was added for each reaction. The primers of genes of interest were designed in expressed constitutive exons flanking the target exon, using the Primer Premier 6.0 software (Premier Biosoft International, Palo Alto, CA, USA). GAPDH was used as the internal control. All primers are listed in Table II. The candidate genes for ASG validation were selected according to following criteria: i) Higher absolute value of SI was firstly considered; ii) whole exon gain/skip was privileged; iii) first and last alternative exons were excluded because of primer design difficulties.

Table II Primer sequences for alternatively spliced gene confirmation by RT-PCR.

Table II

Primer sequences for alternatively spliced gene confirmation by RT-PCR.

Gene symbol	Transcript ID	Primer sequences (5′-3′)
ADH1C	NM_000669	F: CGTTCAGATGAGCATGTGGTT
		R: AAGGTGCTGACGCCGAC
KDM5D	NM_004653	F: TTAAGGCCCGACATGGAACC
		R: CCGCTCCACATTGGGAATCT
PDE4B	NM_002600	F: GCAGGAGTGTGATGACGGTG
		R: GATCCAGTGGACTCCGACCT
USP9Y	NM_004654	F: ACTGTGCGTTCTTCTCCGTCA
		R: AAGACACAAGCATAAAGGTAGCAG
ADAMTSL3	NM_207517	F: TGTCCTGGACGTTGCATGGG
		R: GCAGCACCTTTGTTTGTAGCG
NTRK2	NM_006180	F: ACAATGCACGCAAGGACTTC
		R: AAATCTCCCACAACACGACCC
C7	NM_000587	F: CCTCAGGTTGGCATTTTGTCG
		R: GCAATGGCACAGACAATGGG
PKP2	NM_004572	F: TGTGTGGGGCCTTGAGAAAC
		R: CTCCGTCAGCGTAAGCAATG
RXFP1	NM_001253732	F: TGTAACGGTGTGGACGACTG
		R: ACCGATGGAACAGCTCGTAA
MLPH	NM_001042467	F: TGCTTGCCCCCATTATCCAG
		R: CTCGTTCAGATGGGCAGTGT
GAPDH	NM_002046.5	F: CTCTGCTCCTCCTGTTCGAC
		R: GCGCCCAATACGACCAAATC

[i] F, forward; R, reverse; RT-PCR, reverse transcription-polymerase chain reaction; ADH1C, alcohol dehyrogenase 1C; KDM5D, lysine-specific demethylase 5D; PDE4B, phosphodiesterase 4B; USP9Y, ubiquitin specific peptidase 9, Y-linked; ADAMTSL3, A disintegrin-like and me-talloprotease domain with thrombospondin type I motifs-like 3; NTRK2, neurotrophic tyrosine kinase, receptor, type 2; C7, complement component 7; PKP2, plakophilin-2; RXFP1, relaxin/insulin-like family peptide receptor 1; MLPH, melanophilin.

Statistical analysis

Student's t-test was used to determine the significance between groups. The data were expressed as the mean ± standard deviation. P<0.05 was considered to indicate a statistically significant difference.

Results

ASG detection and validation during osteogenic differentiation of CESCs

Analysis of HTA 2.0 data was performed using rigorous statistical methods to detect ASGs in the osteogenic differentiation of CESCs. According to the criteria and the SI algorithm mentioned above, this analysis of genome-wide AS identified 11,040 alternatively spliced exons, which belonged to 3,802 ASGs during osteogenic differentiation of CESCs. In addition, 6,149 (55%) alternatively spliced exons with ≥2 SI value were considered as 'general exon inclusion' events, while the remaining 4,891 (45%) exons were referred to as 'general exon exclusion' events. Furthermore, it was found that 52% (1,990/3,802) of the ASGs contained 83% (9,228/11,040) of the alternatively spliced exons, thus confirming that multiple AS events can occur to the same gene. During osteogenic differentiation of CESCs, each ASG had 4.6 (9,228/1,990) alternatively spliced exons on average. The ADH1C gene was a typical example, which had 9 alternatively spliced exons detected and suggested complicated splicing regulation. Based on these results, 10 ASGs were selected for RT-PCR validation (Table III). Fig. 2 showed that 6/10 the selected ASGs were validated successfully.

Figure 2 ASGs confirmed in osteogenic differentiation of CESCs by RT-PCR. A total of 6/10 ASGs were confirmed successfully by RT-PCR. The control sample refers to undifferentiated CESCs and GAPDH was used as internal control. ASG, Alternatively spliced gene; CESCs, cartilage endplate-derived stem cells; RT-PCR, reverse transcription-polymerase chain reaction; OD, osteogenically differentiated; ADH1C, alcohol dehyrogenase 1C; C7, complement component 7; MLPH, melanophilin; NTRK2, neurotrophic tyrosine kinase, receptor, type 2; PDE4B, phosphodiesterase 4B; RXFP1, relaxin/insulin-like family peptide receptor 1; GAPDH, glyceraldehyde 3-phosphate dehydrogenase; Incl, including; Excl, excluding; Ex, exon; bp, base pairs.

Table III Summary of alternatively spliced genes selected for validation by RT-PCR.

Table III

Summary of alternatively spliced genes selected for validation by RT-PCR.

Gene symbol	AS exon	PSR ID	Splicing index	Microarray results	RT-PCR results
ADH1C	5	PSR04019704.hg.1	8.1	Exon inclusion	Exon inclusion
C7	12	PSR05002182.hg.1	5.02	Exon inclusion	Exon inclusion
MLPH	5	PSR02023043.hg.1	5.33	Exon inclusion	Exon inclusion
NTRK2	22	PSR09004159.hg.1	4.94	Exon inclusion	Exon inclusion
PDE4B	4	PSR01011724.hg.1	4.47	Exon inclusion	Exon inclusion
RXFP1	7	PSR04011118.hg.1	−5.4	Exon exclusion	Exon exclusion
ADAMTSL3	29	PSR15007327.hg.1	−3.94	Exon exclusion	Exon inclusion
KDM5D	3	PSR0Y001860.hg.1	7.87	Exon inclusion	Exon exclusion
PKP2	7	PSR12017926.hg.1	−5.17	Exon exclusion	Exon inclusion
USP9Y	39	PSR0Y000610.hg.1	10.43	Exon inclusion	Exon exclusion

[i] AS, alternatively spliced; PSR, probe selection region; RT-PCR, reverse transcription-polymerase chain reaction; ADH1C, alcohol dehyrogenase 1C; C7, complement component 7; MLPH, melanophilin; NTRK2, neurotrophic tyrosine kinase, receptor, type 2; PDE4B, phosphodiesterase 4B; RXFP1, relaxin/insulin-like family peptide receptor 1; ADAMTSL3, A disintegrin-like and metalloprotease domain with thrombospondin type I motifs-like 3; KDM5D, lysine-specific demethylase 5D; PKP2, plakophilin-2; USP9Y, ubiquitin specific peptidase 9, Y-linked.

Molecular function analysis of ASGs during osteogenic differentiation of CESCs

GO enrichment analysis was performed on the ASGs during osteogenic differentiation of CESCs. The results suggested that numerous important GO terms were regulated by AS in osteogenic differentiation of CESCs, including regulation of transcription, metal ion binding, signal transduction and cell adhesion. Fig. 3 highlighted the top 10 GO functions regulated in the biological process, molecular function and cellular component categories during osteogenic differentiation of CESCs at the level of alternative splicing.

Figure 3 GO analysis at the AS level during osteogenic differentiation of CESCs. The top 10 GO terms regulated in (A) biological process, (B) molecular function and (C) cellular component during osteogenic differentiation of CESCs at the level of AS are shown. GO, gene ontology; CESCs, cartilage end-plate-derived stem cells; AS, alternative splicing.

The 3,802 ASGs were analyzed by KEGG pathway analysis in order to determine functional cellular pathways regulated during osteogenic differentiation of CESCs. Based on the results, several cellular pathways were affected, including the MAPK signaling pathway, insulin signaling pathway, cell adhesion molecules and calcium signaling pathway. Fig. 4 demonstrated the top 10 KEGG pathways regulated in osteogenic differentiation of CESCs at the level of AS.

Figure 4 KEGG pathway analysis of alternatively spliced genes in osteogenic differentiation of CESCs. The top 10 KEGG pathways regulated in osteogenic differentiation of CESCs at the level of alternative splicing are shown. KEGG, Kyoto Encyclopedia of Genes and Genomes; CESCs, cartilage endplate-derived stem cells.

Discussion

Tissue engineering is a bioengineering method to combine seed cells, biomaterials and biological factors to repair and regenerate damaged tissues and organs (36). The differentiation capacity is of great importance to seed cells, usually stem cells, to accomplish the repair and regeneration tasks. In recent years, the molecular mechanism of osteogenic differentiation of stem cells has been extensively studied (37–39). In particular, it is a very powerful approach to analyze gene transcription and translation on a genome-wide scale to fully elucidate the mechanisms of osteogenic differentiation of stem cells. This approach has been applied to several previous studies to investigate the global gene expression, and post-transcriptional and epigenetic changes in the differentiation process (40–43). However, little contribution has been made to obtain a comprehensive and coherent view of AS mechanisms of stem cell osteogenic differentiation in a genome-wide scale. The present study detected and verified ASEs in osteogenic differentiation of CESCs, and analyzed molecular functions and pathways using bioinformatics methods. To the best of our knowledge, the present study is the first to determine the ASEs in osteogenic differentiation of stem cells on the whole genome level.

The HTA 2.0 platform was used in the present study to cover the entire genome to both identify evidence-based sequences and discover novel ASEs. Previously, researchers tended to analyzed expression sequence tag (EST) data for the purpose of discovering novel ASEs (44). However, it should be noted that the number of sequenced ESTs is really small and the available EST data leave more ASEs undetectable. The present study detected 3,802 ASGs with 11,040 ASEs in osteogenic differentiation of CESCs. Additionally, 6 ASGs were validated successfully by RT-PCR. These novel validated isoforms of ASGs have not been recorded in the NCBI Reference Sequences database previously, hence further studies are required to be carried out to determine cellular and molecular functions of these isoforms. In terms of the molecular functions and pathways of these detected ASGs, GO analysis presented enrichment of numerous GO terms, including signal transduction, cell differentiation and cell cycle arrest. The enriched signal transduction process suggested that AS may modulate signaling pathways to exert its influence on cellular function networks. Besides, during osteogenic differentiation process, the state of CESCs shifted from proliferation to differentiation. Therefore, cell cycle arrest occurred and cell growth stopped slowly. Furthermore, the present study also performed KEGG pathway analysis and revealed several enriched signaling pathways, including focal adhesion, ECM-receptor interaction and the MAPK signaling pathway. The focal adhesion signaling pathway is critical for cell-matrix adhesion and serves important roles in numerous biological processes, including cell proliferation, cell differentiation and gene expression regulation (45–47). The significantly enriched focal adhesion pathway suggested an CESC and ECM interaction in osteogenic differentiation. Additionally, the MAPK signaling cascade is a highly conserved module, which is involved in diverse cellular processes, including cell differentiation, migration and proliferation (48–50). Therefore, the results of pathway analysis indicate that AS interacts with signaling pathways to regulate the osteogenic differentiation process.

The present stusdy used the HTA 2.0 platform to investigate AS events in osteogenic differentiation of CESCs. The results revealed various ASGs, and associated molecular functions and pathways. Further research is required to illuminate downstream mechanisms of AS modulation. The structure and function of novel isoforms may be potential targets in future studies.

Acknowledgments

The present study was financially supported by the National Natural Science Foundation of China (nos. 81401801 and 81301944).

References