CEP samples were obtained from three patients that underwent discectomy operations at the Xinqiao Hospital, Third Military Medical University (Chongqing, China; Table I). All procedures in the present study were approved by the Ethics Committee of Xinqiao Hospital, Third Military Medical University and were conducted in accordance with the Declaration of Helsinki; written informed consent was obtained from each patient.

CEP samples were cut into pieces and digested in Dulbecco's modified Eagle's medium (DMEM)/Ham's F12 medium (Hyclone; GE Healthcare Life Sciences, Logan, UT, USA) containing 1% fetal calf serum (FCS; Gibco; Thermo Fisher Scientific, Inc., Waltham, MA, USA) and 0.2% collagenase II (Sigma-Aldrich; Merck KGaA, Darmstadt, Germany) for 12 h at 37°C. The suspension was passed through a cell filter (70 µm) and centrifuged at 110 × g for 5 min at room temperature. Following removal of the supernatant, the pellet was resuspended in 3 ml DMEM/F12 containing 10% FCS and 1% penicillin-streptomycin (Hyclone; GE Healthcare Life Sciences). The cells were then cultured in a 25-cm² cell-culture flask at 37°C and 5% CO₂ for 6 weeks. Cells were subjected to agarose selection following the first passage of expansion, as previously described (9). Briefly, Costar culture dishes (Corning Inc., Corning, NY, USA) were pre-coated with 1% low melting-point agarose, and a mixture containing DMEM/F12 (0.5 ml), 2% low melting-point agarose (0.5 ml) and 5×10⁴ suspended CEP cells in culture medium (1 ml) was transferred to the culture dishes and incubated in humidified atmosphere in 5% CO₂ and 37°C. The culture medium was changed twice per week and, 6 weeks later, cell aggregates were carefully aspirated with a sterile Pasteur pipette and cultured in a 25 cm² cell culture flask in 5% CO₂ and 37°C for 2 weeks to reach 80–100% confluence. CESCs which were sub-cultured to the third passage, after ~6 weeks of culture, were used.

For hypoxic culture conditions, CESCs at a density of 2.5×10⁴/ml were cultured in a 1% O₂ microenvironment at 37°C. For normoxic culture conditions, CESCs were cultured in a 21% O₂ microenvironment at 37°C. The culture medium was changed twice per week over a 21-day period.

Total RNA was extracted from cells (~5×10⁴) using TRIzol (Invitrogen; Thermo Fisher Scientific, Inc.) according to the manufacturer's protocol. RNA was purified with an RNeasy kit (Qiagen GmbH, Hilden, Germany) according to the manufacturer's protocol. The RNA purity was determined by the ratio of optical density (OD) value in 260 nm/OD value in 280 nm using a NanoDrop spectrophotometer (NanoDrop Technologies; Thermo Fisher Scientific, Inc., Wilmington, DE, USA). RNA quality was assessed by 1% agarose gel electrophoresis. The samples with a bright band of ribosomal 28S/18S RNA at a ratio >1.5:1 were used for the microarray analysis. For gene expression profiling, RNA was hybridized using the GeneChip HTA 2.0 system (Affymetrix; Thermo Fisher Scientific, Inc.). A total of 100 ng RNA was used to generate biotinylated and amplified sense-strand cDNA from the whole expressed genome using GeneChip WT PLUS Reagent kit (Applied Biosystems; Thermo Fisher Scientific, Inc.) according to the manufacturer's protocol. GeneChip HTA 2.0 cartridges were hybridized in a 45°C incubator (GeneChip Hybridization Oven 640; Thermo Fisher Scientific, Inc.) for 16 h under rotation at 60 rpm. After the hybridization, the microarrays were stained using the GeneChip Fluidics Station 450 (Thermo Fisher Scientific, Inc.) and screened using the GeneChip Scanner 3000 7G (Thermo Fisher Scientific, Inc.) according to the manufacturer's protocol. Gene probes of the HTA 2.0 target exons and junctions and are able to detect gene expression and AS events simultaneously. Microarray screening was performed by CapitalBio Corporation (Beijing, China) according to the manufacturer's protocol. Image signals of the probes were obtained and saved as DAT files. The GeneChip Command Console (AGCC) software version 3.2 (Affymetrix; Thermo Fisher Scientific, Inc.) was used to convert the image signals (DAT files) into digital signals (CEL files), which were converted to CHP files for probeset signal integration, quantile normalization and background correction through the robust multichip analysis (RMA) algorithm using the Affymetrix Expression Console software version 1.3 (Affymetrix; Thermo Fisher Scientific, Inc.). The CHP files were then analyzed by Affymetrix Transcriptome Analysis Console software version 3.0 (Affymetrix; Thermo Fisher Scientific, Inc.) to determine the differentially expressed genes (DEGs) and alternatively spliced genes (ASGs). GO terms and related signaling pathway enrichments were determined using web tools, such as DAVID (https://david.ncifcrf.gov), Molecule Annotation System (http://bioinfo.capitalbio.com/mas3) and KEGG (http://www.genome.jp/kegg. The workflow of the present study procedures was as previously described (20). The microarray data have been submitted to the Gene Expression Omnibus (accession no. GSE84484).

Data from the cell samples cultured under normoxia were used as a control to calculate the fold change in gene expression. The default filter criteria were set as fold change ≤-2 or ≥2 (linear) for significant DEGs. For ASG identification, the splicing index (SI) calculation mode was used to analyze the level of exon inclusion/exclusion. SI indicates the ratio of the signal intensity of an exon normalized to that of the target gene between two groups. SI values were calculated using the following formulas:

Where NI(i,j)_A represents the normalized intensity (NI) of the i^th exon signal normalized to that of the j^th gene signal under condition A, N represents the normoxic culture condition and H represents the hypoxic culture condition. The default filter criteria were defined as SI (linear) ≤-2 or ≥2.

Total RNA was extracted from cells (~5×10⁴ cells/ml) using TRIzol reagent (Invitrogen; Thermo Fisher Scientific, Inc.) according to the manufacturer's protocol. RNA was purified with an RNeasy kit (Qiagen GmbH, Hilden, Germany) according to the manufacturer's protocol. The RNA purity was determined by the ratio of optical density (OD) value in 260 nm/OD value in 280 nm using a NanoDrop spectrophotometer (NanoDrop Technologies; Thermo Fisher Scientific, Inc.). RNA quality was assessed by 1% agarose gel electrophoresis. The samples with a bright band of ribosomal 28S/18S RNA at a ratio >1.5:1 were used. RNA was reverse transcribed into cDNA with the PrimeScript RT Master Mix kit (cat no. RR047A; Takara Bio, Inc., Otsu Japan). qPCR was performed with the ABI 7900H Real-Time PCR system and the SYBR Premix Ex Taq II kit (cat no. RR820A; Takara Bio, Inc.) with the following conditions: At 95°C for 30 sec, then 40 cycles at 95°C for 5 sec and at 60°C for 34 sec; a dissociation curve analysis was performed in a temperature range of 60–95°C. The mRNA expression levels of each gene were normalized to the level of β-actin expression and analyzed using the 2^−∆∆Cq method (21). Data were analyzed using StepOne software version 2.3 (Applied Biosystems; Thermo Fisher Scientific, Inc.). The primers used for qPCR are listed in Table II-A.

Total RNA extraction and reverse transcription of cDNA were carried out as aforementioned. For semi-quantitative PCR, cDNA (2 µl) was combined with Premix Taq (cat no. RR901A; Takara Bio, Inc.) and primers that were designed flanking the constitutively expressed exons; primers are listed in Table II-B. PCR thermocycling conditions were as follows: At 94°C for 30 sec, then 30 cycles at 94°C for 30 sec, at 55°C for 30 sec and at 72°C for 60 sec. PCR products were resolved on a 1.5% agarose gel, using a DNA ladder (Thermo Fisher Scientific, Inc.) and gold view nucleic acid gel stain regent (Beijing Solarbio Science & Technology Co., Ltd., Beijing, China). Gels were imaged using a fluorescence imaging analysis system (Bio-Rad Laboratories, Inc., Hercules, CA, USA). The mRNA expression levels of each ASG were normalized to β-actin. Differences of ASG between two groups were reflected in the ratio of exon-inclusion fragment expression level vs. β-actin expression level/exon-exclusion fragment expression level vs. β-actin expression level. The ASGs used for validation were in accordance with the following selection criteria: i) The whole exon exclusion/inclusion was included; ii) relatively higher SI absolute value was included; iii) the first or the last AS exon was excluded.

The data of 3 independent experiments were expressed as the mean ± standard deviation. Independent-samples t-test was used to determine the significance between groups. P<0.05 was considered to indicate a statistically significant difference. Statistical analysis was performed using SPSS software version 19.0 (IBM Corp., Armonk, NY, USA).

HTA 2.0 data analysis identified a total of 448 DEGs, of which 330 (73%) were upregulated and 118 (27%) were downregulated. Based on these results, 10 DEGs were chosen for validation by RT-qPCR, and 8 of these were validated successfully. Since the GeneChip signal increases linearly as target gene expression increases, but the qPCR signal increases exponentially as target gene expression increases, the ratio calculated by GeneChip and qPCR was different. Thus, a similar tendency for increase or decrease in the results from both methods was considered to indicate successful validation of the DEGs. In the present study, the expression of 5 genes [vascular cell adhesion molecule 1, cartilage oligomeric matrix protein, fibroblast growth factor 9, interleukin 1 receptor-associated kinase 3 (IRAK3) and complement factor D (CFD)] was upregulated in the hypoxia group compared with the normoxia group. Conversely, the expression of 3 genes (transferrin receptor, phosphoglycerate mutase 2 and calpain 6) was downregulated in the hypoxia group compared with the normoxia group (Fig. 1).

Enriched biological processes, molecular functions and cellular components in CESCs were identified by GO analysis of the DEGs under normoxia and hypoxia. GO terms were highly enriched in CESCs under hypoxia, such as multicellular organismal process, cargo receptor activity and extracellular region part; the top 10 GO terms from each category are indicated in Fig. 2A-C.

The KEGG mapping tool was used to identify the enriched functional signaling pathways for the identified DEGs under hypoxia in CESCs. The results revealed that several signaling pathways were significantly enriched, such as drug metabolism-cytochrome P450, hematopoietic cell lineage and rheumatoid arthritis; the top 10 KEGG pathways are in Fig. 2D.

Genome-wide AS analysis between the normoxic group and hypoxic group identified 9,475 alternatively spliced exons, which belonged to 1,808 ASGs. In addition, 3,677 (39%) of the alternatively spliced exons with an SI value ≥2 were defined as ‘exon inclusion’ events, whereas the remaining 5,798 (61%) exons were considered as ‘exon exclusion’ events (with an SI value ≤-2). Under hypoxic conditions, each ASG had an average of 5.2 (9,475/1,808) alternatively spliced exons, which confirmed the occurrence of multiple alternative splicing events in the same gene. For example, nine alternatively spliced exons were identified for nuclear receptor subfamily 4 group A member 3 (NR4A3), which suggested that AS regulation was complex. In addition, 114 of the 1,808 ASGs were also identified as DEGs, indicating an underlying link between AS and differential gene expression. A total of 10 ASGs were chosen for validation by semi-quantitative PCR analysis, of which 7 were successfully validated. Successful validation was confirmed when the HTA 2.0 array results and semi-quantitative PCR demonstrated that an exon was more likely to be spliced or preserved under hypoxic conditions compared with under normoxic conditions. In the present study, AS events were observed and validated in GDNF family receptor α2, cell adhesion molecule 1, NR4A3, relaxin/insulin-like family peptide receptor 1, NADH:ubiquinione oxidoreductase subunit C1, transmembrane protein 106B and fibronectin 1 under hypoxic conditions compared with under normoxic conditions (Fig. 3).

GO term enrichment analysis of ASGs in CESCs under normoxic and hypoxic conditions was conducted to determine the enriched biological processes, molecular functions and cellular components. The results revealed enrichment in GO terms, such as antigen processing and presentation of peptide antigen via major histocompatibility complex (MHC) class I, protein binding and cytoplasm; the top 10 GO terms for ASGs under hypoxic conditions in CESCs are presented in Table III.

The 1,808 identified ASGs were further examined by KEGG pathway analysis to identify the enriched functional pathways. These results revealed that several vital pathways were significantly enriched, such as focal adhesion, cell adhesion molecules and axon guidance; a list of the top 10 KEGG pathways identified for ASGs under hypoxic conditions in CESCs is presented in Table IV.