Identification of lncRNA and mRNA expression profiles in rat spinal cords at various time‑points following cardiac ischemia/reperfusion.

The identification of the expression patterns of long non‑coding RNAs (lncRNAs) and mRNAs in the spinal cord under normal and cardiac ischemia/reperfusion (I/R) conditions is essential for understanding the genetic mechanisms underlying the pathogenesis of cardiac I/R injury. The present study used high‑throughput RNA sequencing to investigate differential gene and lncRNA expression patterns in the spinal cords of rats during I/R‑induced cardiac injury. Male Sprague Dawley rats were assigned to the following groups: i) Control; ii) 2 h (2 h post‑reperfusion); and iii)v0.5 h (0.5 h post‑reperfusion). Further mRNA/lncRNA microarray analysis revealed that the expression profiles of lncRNA and mRNA in the spinal cords differed markedly between the control and 2 h groups, and in total 7,980 differentially expressed (>2‑fold) lncRNAs (234 upregulated, 7,746 downregulated) and 3,428 mRNAs (767 upregulated, 2,661 downregulated) were identified. Reverse transcription‑quantitative polymerase chain reaction analysis was performed to determine the expression patterns of several lncRNAs. The results indicated that the expression levels of lncRNA NONRATT025386 were significantly upregulated in the 2 and 0.5 h groups when compared with those in the control group, whereas the expression levels of NONRATT016113, NONRATT018298 and NONRATT018300 were elevated in the 2 h group compared with those in the control group; however, there was no statistically significant difference between the 0.5 h and control groups. Furthermore, the expression of lncRNA NONRATT002188 was significantly downregulated in the 0.5 and 2 h groups when compared with the control group. The present study determined the expression pattern of lncRNAs and mRNAs in rat spinal cords during cardiac I/R. It was suggested that lncRNAs and mRNAs from spinal cords may be novel therapeutic targets for the treatment of I/R‑induced cardiac injury.

Introduction

Cardiac ischemia/reperfusion (I/R) injury is associated with various etiological factors (1-3), including primary heart conditions, and neuropathic, vascular or systemic disorders. However, the pathological origin of I/R-induced cardiomyocyte death remains poorly understood. The neuropathic accumulation associated with myocardial I/R is thought to directly or indirectly activate sensory or sympathetic fibers that innervate the heart. However, as the underlying heart-specific neuronal pathway and mechanism are unknown, neurogenic therapeutic interventions often only have limited success.

It is well established that the autonomic nervous system serves a crucial physiological role in regulating cardiac function (4-9). Previous studies have demonstrated that some myocardial ischemic events are triggered by the autonomic nervous system, and a sympathetic-parasympathetic imbalance may lead to the pathophysiological development of myocardial ischemia (10,11). Our previous study revealed that an injection of retrograde tracer pseudorabies virus (PRV)-614 into the left ventricle wall in mice resulted in the retrograde infection of neurons in the intermediolateral nucleus of the spinal cord and the rostral ventrolateral medulla of the brainstem via the sympathetic pathway (12). The spinal cord has been implicated in the pathogenesis of cardiac injury caused by I/R (4,13-16). Furthermore, there is convincing evidence that the heart receives dense innervation from sensory afferent fibers, which peripherally release a variety of vasodilator neuropeptides, such as calcitonin gene-related peptide (CGRP) (17) and substance P (SP), in response to local stimuli (18). Despite extensive research in this area, the mechanisms underlying cardiac I/R injury are largely unknown. Therefore, there is an urgent requirement for more information and a genomic approach may prove helpful.

High-throughput RNA sequencing (RNA-seq) is a powerful tool that has been used to identify novel protein-coding and non-coding RNA transcripts involved in the regulation of gene expression (19-24). Recent research has focused on the cardiac long non-coding RNAs (lncRNAs) implicated in cardiac I/R injury (25-29); however, few studies have explored the important role of the spinal cord during focal cardiac I/R. Therefore, the present study was designed to identify the expression patterns of differentially expressed lncRNAs and mRNAs in the spinal cord under normal and cardiac I/R conditions, with the aim to gain a better understanding of the genetic mechanisms underlying the pathogenesis of cardiac I/R injury. The present study also determined the expression levels of various genes in the spinal cord at different time-points during cardiac I/R injury.

Materials and methods

Animals

A total of 24 male Sprague Dawley (SD) rats aged 8-10 weeks (200-240 g; specific pathogen-free grade; no. 42000600010250) were supplied by the Experimental Animal Research Center of Hubei Province (Hubei, China). The present experiment protocol was approved by the Institutional Ethical Committee of Tongji Hospital, Tongji Medical College, Huazhong University of Science and Technology (Hubei, China; no. TJ-A20150804). All animals were humanely treated according to the National Institutes of Health Guide for the Care and Use of Laboratory Animals (revised 2011) and the Guide for the Care and Use of Laboratory Animals (National Academic Press, USA; revised 2011) (30). Animals were housed in compliance with the Specific Pathogen-Free Animal Criteria, and maintained at a standard temperature of 21-23°C and 65±5% humidity under a 12-h light/dark cycle condition, with ad libitum access to artificial feed (food and water).

Myocardial I/R injury model

Surgical procedures to establish the myocardial I/R injury model were performed according to previously described methods (31-33).

Experimental protocol

Experiment A

A total of 6 SD rats were randomly allocated into 2 groups. The model group (n=3) received 2 h reperfusion following 30-min occlusion of the left anterior descending coronary artery by pulling the reversible trap (I/R group). The control group (n=3) received the same surgical procedure, without any occlusion of the coronary artery and reperfusion (sham group). T1-T4 spinal cord samples were collected for RNA-seq and reverse transcription-quantitative polymerase chain reaction (RT-qPCR) analysis.

Experiment B

A total of 18 SD rats were randomly assigned into 3 groups: i) A control group (n=6); ii) a 0.5 h group (0.5 h reperfusion following 30 min ischemia; n=6); and iii) a 2 h group (2 h reperfusion following 30 min ischemia; n=6). T1-T4 spinal cord segments were obtained for RT-qPCR analysis.

Tissue preparation and microarray gene expression analysis

The rats were sacrificed following the completion of the aforementioned experiments. Briefly, upon completion of the experiments, the rats were anesthetized by intraperitoneal injection with 100 mg/kg body weight ketamine and 10 mg/kg body weight xylazine (34). Then the rats were quickly decapitated to limit animal suffering and minimize the effects on the experimental results, and the T1-T4 spinal cord segments were immediately removed and frozen with liquid nitrogen for 1 min, then stored at −80°C until required. Total RNA from each animal was quantified using a mirVana miRNA Isolation kit (Ambion; Thermo Fisher Scientific, Inc., Waltham, MA, USA), and RNA integrity was assessed according to the manufacturer's protocol, which included standard denaturing agarose gel electrophoresis (35-37). For RNA-seq, the micro-array work was performed by CapitalBio Technology Co. Ltd. (Beijing, China), whereby 6 tissue samples (3 model group samples and 3 control group samples) were used for mRNA and lncRNA microarray analysis (38).

The present study used an Agilent Array platform for microarray analysis. Tissue preparation and microarray hybridization were performed based on the manufacturer's standard protocols (Agilent Technologies, Inc., Santa Clara, CA, USA) with minor modifications. Briefly, the mRNA was purified from the total RNA once the ribosomal RNA was removed using an mRNA-ONLY Eukaryotic mRNA Isolation kit (Epicentre; Illumina, Inc., San Diego, CA, USA). Each sample was then amplified and transcribed into fluorescent complementary RNA (cRNA) along the entire length of the transcripts using a random priming method (39). The labeled cRNAs were hybridized onto the mouse lncRNA Array v2.0 (8×60 K; Arraystar, Inc., Rockville, MD, USA). Then the arrays were scanned using a G2505C Scanner (Agilent Technologies, Inc.).

Bioinformatics analysis

Gene Ontology (GO) annotations were employed to investigate the differentially expressed mRNAs and lncRNAs in the T1-T4 spinal cord segments, according to the GO database (www.geneontology.org/). For the functions of genes and their products, the GO database describes 3 biological functional groups: Biological process, Cellular compartment, and Molecular function. The present study conducted GO functional enrichment analysis on the differentially expressed mRNAs involved in protein-protein interaction (PPI) networks. In addition, the key regulatory pathways in the spinal cord that respond to I/R-induced cardiac injury were also investigated using Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis (www.genome.jp/kegg) and STRING database (string-db.org).

RT-qPCR analysis

The present study extracted total RNA from the upper thoracic spinal cord segments (T1-T4) (40) using TRIzol® reagent (Invitrogen; Thermo Fisher Scientific, Inc.) according to our previous research (35-37). The primers for RT-qPCR were designed based on the lncRNA sequences (Table I), and were synthesized and purified at Invitrogen (Thermo Fisher Scientific, Inc.). The RT reactions were performed using a iScript™ cDNA Synthesis kit (Bio-Rad Laboratories, Inc., Hercules, CA, USA). RT-qPCR was performed using a ABI StepOnePlusä Multicolor system with the SsoAdvanced™ Universal SYBR®-Green Supermix (Bio-Rad Laboratories, Hercules, Inc.). The PCR thermocy-cling conditions were as follows: Initial denaturation at 95°C for 1 min followed by 40 cycles of 95°C for 15 sec, 60°C for 15 sec and 72°C for 45 sec. Compared with the averages for the housekeeping gene (GAPDH), the data were quantified using the 2−∆∆Cq method for relative fold-change, as described previously (41).

Table I

Primer sequences for reverse transcriptase-quantitative polymerase chain reaction.

A, Upregulated
Gene	Length (bp)	Forward (5′-3′)	Reverse (5′-3′)
NONRATT025386	86	GGGTCTGGGGTGGGCTAA	GGAGGTTTCTGAGTGGGATGTG
NONRATT016113	96	CCACAAGCGTCTCGGGATT	AGCGAAAACAGTCATTTTAACCAA
NONRATT018298	166	GACAGTCAACGGAACCAAACTAA	CGTGAACAAAAGCAAGCAAAC
NONRATT018300	178	GCCAACAACCAGTAAGAACCAC	CCATACCTTTGCTACTTTGGAGA
NONRATT020994	150	GAACGCCACCCCACCAT	CCTTGAAGTCTGAGGCAGGAA
B, Downregulated
Gene	Length (bp)	Forward (5′-3′)	(5′-3′)
XR_590210.1	128	TTTCAGCCCATCAATGGTTTC	TCCTCAGGAGTGCCCTTTCT
NONRATT002188	105	TTCCTACATACTGAGCAACGACC	CCTACCTGTAGCTGCCACTCC
XR_589980.1	115	GGATGCCCACTCAAGGGTC	GATGATAAATGCTTGCCCACC
XR_598701.1	145	AACAATGGGGACGGTAGTGC	GAAATGAACCTGGGAGAAACG
XR_590197.1	92	ACTTCCCTGGATTCTGCTCTG	GGGTCCCCTAACACTATTGCTT
GAPDH	68	CGCTAACATCAAATGGGGTG	TTGCTGACAATCTTGAGGGAG

Statistical analysis

Data are expressed as the mean ± standard error of the mean. The data and statistical graphs were analyzed using the GraphPad Prism v6.02 package (GraphPad Software, Inc., La Jolla, CA, USA). Between-group counts were compared using a Student's t-test (Mann-Whitney U), and the data from three groups were analyzed by one-way analysis of variance followed by Dunn's post hoc test. P<0.05 was considered to indicate a statistically significant difference.

Results

Characteristics of ischemic myocardial tissues

With regard to the animals in the model group, there were evident cyanotic changes in the myocardium of the occluded area 30 min post-ischemia, and a reactive hyperemic response after refilling of the left anterior descending coronary artery. With regard to the samples in the model group, examination under a dissecting microscope revealed discoloration of the occluded distal myocardium at the early stage of reperfusion in the infarct region (data not shown).

Expression profiling of lncRNAs in the spinal cord 2 h post-cardiac I/R

To select possible targets of lncRNAs in the model and control groups, the present study detected up to 16,987 coding transcripts in the T1-T4 spinal cords 2 h post-reperfusion. A total of 3,621 upregulated and 13,366 downregulated lncRNAs were identified in the spinal cords. On average, 234 lncRNAs were upregulated in the spinal cords of the model group, compared with those in the control group, whereas an average of 7,746 lncRNAs were downregulated (>2.0-fold-change; P<0.05). The distributions of the log2 ratios of the lncRNAs in the model and control samples were nearly identical; Fig. 1 presents the heatmaps of the expression ratios (log2 scale) of the lncRNAs in the spinal cords. The top 20 up- and downregulated mRNAs are listed in Tables II and III.

Figure 1

Differential expression of lncRNA in T1-T4 spinal cord sections in the cardiac I/R injury model. (A) Scatter plot comparing global lncRNA expression profiles in the spinal cords of cardiac I/R and control rats. Red indicates upregulated and green indicates downregulated genes. Black indicates the absence of a statistically significant difference, meaning it failed to pass the cutoff values of 1 and -1 in the log2 scale and P<0.05. (B) Heatmap presenting hierarchical clustering of lncRNAs with expression changes of >2-fold. In the clustering analysis, up- and downregulated genes are represented by red and green, respectively. I/R, ischemia/reperfusion; lncRNA, long non-coding RNA.

Table II

Top 20 upregulated lncRNAs in the spinal cord at 2 h post-reperfusion.

LncRNAs (sequence name)	Fold-change (R/N)	RNA length	Chromosome log2
gi|672017878|ref|XR_345533.2|	83.01321	1,727	Chr3
NONRATT025386	23.99574	563	Chr6
NONRATT024318	23.93667	473	Chr6
gi|672024701|ref|XR_599241.1|	13.07902	728	Chr10
NONRATT025509	10.08656	553	Chr7
gi|672017768|ref|XR_600487.1|	9.916187	634	Chr3
NONRATT025839	9.860808	508	Chr7
NONRATT023339	9.858076	706	Chr5
NONRATT000120	9.325455	550	Chr1
NONRATT002260	8.604115	692	Chr1
gi|672055933|ref|XR_592974.1|	8.306979	8,102	Chr6
uc.126	7.209348	271	-
gi|672086728|ref|XR_597427.1|	6.502431	4,903	Chr20
NONRATT008414	6.458688	518	Chr13
NONRATT026470	6.452916	655	Chr7
NONRATT015818	6.317123	255	Chr2
gi|672080453|ref|XR_596511.1|	6.286473	1,683	Chr16
gi|672027556|ref|XR_340041.2|	6.236808	633	Chr13
gi|672030740|ref|XR_598338.1|	6.169331	1,360	Chr18
NONRATT016237	6.087859	817	Chr2

Values are presented as fold-changes in the reperfusion groups (reperfusion 2 h) over the control group (N >2.0-fold; P<0.05). lncRNA, long non-coding; Chr, chromosome; R/N, reperfusion/normal.

Table III

Top 20 downregulated lncRNAs in the spinal cord at 2 h post-reperfusion.

LncRNAs (sequence name)	Fold-change (R/N)	RNA length	Chromosome log2
NR_130708.1	−27.7049	1,379	Chr3
NONRATT028627	−25.2873	410	Chr8
NONRATT021959	−22.8492	709	Chr4
gi|672034655|ref|XR_590005.1|	−21.8020	1,290	Chr1
NONRATT023191	−19.8345	1,977	Chr5
NONRATT025830	−19.4635	525	Chr7
NONRATT027253	−19.4441	1,233	Chr7
NONRATT023189	−19.2427	2,086	Chr5
NONRATT014248	−18.4873	1,049	Chr19
NONRATT008322	−17.5179	451	Chr13
NONRATT014862	−16.4839	619	Chr2
NONRATT016808	−16.4007	518	Chr2
NONRATT017256	−16.2487	786	Chr20
NONRATT011649	−15.8013	566	Chr16
gi|672021532|ref|XR_347699.2|	−15.494	537	Chr7
NONRATT024978	−14.3807	1,486	Chr6
NONRATT012913	−13.9501	338	Chr17
NONRATT004220	−13.7216	2,171	Chr10
NONRATT018550	−13.4908	365	Chr3
NONRATT008489	−13.4256	1,076	Chr13

Values are present as fold-changes in the reperfusion groups (reperfusion 2 h) over the control group (N >2.0-fold; P<0.05). lncRNA, long non-coding; Chr, chromosome; R/N, reperfusion/normal.

Expression profiling of mRNAs in the spinal cord 2 h post-cardiac I/R

To explore the potential role of mRNAs in the T1-T4 spinal cords 2 h post-cardiac I/R, the present study determined the expression profiles of mRNAs by high-throughput RNA-seq. Of the 26,466 mRNAs that were quantified using reads per kilobase per million mapped reads (RPKM) values, 3,428 mRNAs were deregulated by 2-fold following I/R-induced cardiac injury, of which 767 mRNAs were upregulated and 2,661 mRNAs were downregulated; Fig. 2 presents the heatmaps of the expression ratios (log2 scale) of the mRNAs in the spinal cords. The top 20 up- and downregulated mRNAs are listed in Tables IV and V.

Figure 2

Differential expression of mRNA genes in the T1-T4 spinal cord sections in the cardiac I/R injury model. (A) Scatter plot comparing global mRNA gene expression profiles in the spinal cords of cardiac I/R and control rats. Red indicates upregulated and green indicates downregulated genes. Black indicates the absence of a statistically significant difference, meaning it failed to pass the cutoff values of 1 and -1 in the log2 scale and false discovery rate (corrected P-value) <0.05. (B) Heatmap presenting hierarchical clustering of mRNAs with expression changes of >2-fold. In the clustering analysis, up- and downregulated genes are represented by red and green, respectively. I/R, ischemia/reperfusion.

Table IV

Top 20 upregulated mRNAs in the spinal cord at 2 h post-reperfusion.

mRNAs (sequence name)	Fold-change (R/N)	GENE_SYMBOL	GENE_NAME
A_64_P004112	25.97118	-	-
A_64_P151353	25.46518	Acsm5	Acyl-CoA synthetase medium-chain family member 5
A_64_P181171	21.23528	-	-
A_64_P149280	21.07919	Vegfb	Vascular endothelial growth factor B
A_64_P273771	15.23514	-	-
A_64_P260129	14.28584	Gzmc	Granzyme C
A_44_P122912	13.21478	Ces1c	Carboxylesterase 1C
A_44_P553341	13.01391	Lrfn5	Leucine rich repeat and fibronectin type III domain containing 5
A_64_P094055	10.5432	-	-
A_44_P401110	10.10015	Prss40	'Protease, serine, 40'
A_64_P147769	8.814655	Gucy1b2	'Guanylate cyclase 1, soluble, β2'
A_64_P082082	8.515041	-	-
A_64_P045902	8.093776	-	-
A_64_P118367	8.069224	Lrrd1	Leucine-rich repeats and death domain containing 1
A_64_P156605	7.85615	-	-
A_64_P135295	7.737585	Cd300c	CD300c molecule
A_64_P042127	7.403613	Akr1c3	'Aldo-keto reductase family 1, member C3'
A_44_P445575	7.245778	Irx2	Iroquois homeobox 2
A_64_P162476	7.220482	Art1	ADP-ribosyltransferase 1
A_64_P186630	7.156018	Vom2r60	'Vomeronasal 2 receptor, 60'

Values are fold-changes in the reperfusion groups (reperfusion 2 h) over the control group (N >2.0-fold; P<0.05). R/N, reperfusion/normal.

Table V

Top 20 downregulated mRNAs in the spinal cord at 2 h post-reperfusion.

mRNAs (sequence name)	Fold-change (R/N)	GENE_SYMBOL	GENE_NAME
A_64_P141763	−38.1959	Muc16	'Mucin 16, cell surface associated'
A_42_P473398	−33.0016	Cxcl1	Chemokine (C-X-C motif) ligand 1
A_44_P550907	−29.4102	RGD1306750	LOC362451
A_44_P702019	−23.3124	Rfx4	'Regulatory factor X, 4 (influences HLA class II expression)'
A_64_P021855	−22.6481	LOC499643	Similar to hypothetical protein FLJ25371
A_44_P367541	−21.9947	Olr671	Olfactory receptor 671
A_64_P036090	−20.4897	-	-
A_44_P371339	−20.0377	IL6	Interleukin 6
A_64_P141762	−19.8627	-	-
A_64_P018554	−19.1504	-	-
A_64_P101228	−16.1995	Tas2r120	'Taste receptor, type 2, member 120'
A_44_P563447	−15.8884	Wnt10a	'Wingless-type MMTV integration site family, member 10A'
A_44_P698466	−15.406	Lrrtm2	Leucine rich repeat transmembrane neuronal 2
A_44_P461456	−15.3963	Prl4a1	'Prolactin family 4, subfamily a, member 1'
A_44_P378749	−14.7924	LOC100912608	Homeobox protein Hox-A10-like
A_44_P547892	−14.7229	Olr1345	Olfactory receptor 1345
A_64_P001947	−14.5408	-	-
A_64_P009237	−14.3488	-	-
A_64_P029912	−13.3859	Cldnd2	Claudin domain containing 2
A_64_P074460	−13.1395	-	-

Values are fold-changes in the reperfusion groups (reperfusion 2 h) over the control group (N >2.0-fold; P<0.05 by analysis of variance). R/N, reperfusion/normal.

Hierarchical clustering analysis of the differentially expressed mRNAs and lncRNAs

Following the determination of the expression values of the differentially expressed genes (DEGs), the present study carried out hierarchical clustering analysis on the DEGs. As shown in Figs. 1 and 2, the differentially expressed mRNAs and lncRNAs clearly distinguished the cardiac I/R tissues from the control samples. In the cardiac I/R tissues, there were more downregulated genes than upregulated genes (Figs. 1 and 2).

Functional and pathway enrichment analyses

The significantly enriched GO terms (top 30) were comprised of 16 biological processes, 6 cellular compartments, and 8 molecular functions (Fig. 3). It was revealed that the differentially expressed biological processes in the spinal cords were primarily involved in the serotonin receptor signaling pathway, regulation of protein kinase B signaling, regulation of keratinocyte migration, and skeletal muscle satellite cell differentiation. The significantly enriched pathway terms (top 30) primarily involved KEGG pathways including 'olfactory transduction', 'arachidonic acid metabolism', the 'phosphoinositide 3-kinase-protein kinase B signaling pathway', 'extracellular matrix-receptor interaction', 'cytokine-cytokine receptor interaction' and 'neuroactive ligand-receptor interaction' (Fig. 4). The DEGs were analyzed with GO background significant enrichment, which demonstrated the number of genes associated with biological processes, cellular compartments and molecular functions (Fig. 5).

Figure 3

Top 30 significantly enriched GO terms as determined by GO annotation analysis. The blue bars represent biological process classification; the yellow bars represent cellular components; and the red bars represent molecular function. GO, Gene Ontology.

Figure 4

Top 30 significantly enriched pathway terms determined by Gene Ontology annotation analysis. Blue bars represent the KEGG pathway and orange bars represent PANTHER. KEGG, Kyoto Encyclopedia of Genes and Genomes.

Figure 5

Differentially expressed genes analyzed by GO background significant enrichment. Solid bars represent the target gene set and slashed bars represent the gene set. GO, Gene Ontology.

The results of the biological process analysis revealed that the DEGs involved in PPI networks were mainly enriched in 'neurological system processes' (P=1.92×10−11), 'sensory perception' (P=1.92×10−11), the 'detection of chemical stimuli' (P=5.10×10−10) and 'cell surface receptor signaling pathways' (P=4.23×10−10; Fig. 6). A total of 7 cellular components from GO terms were significantly enriched for the DEGs involved in PPI networks: 'Extracellular regions' (P=0.02), 'membrane parts' (P=1.54×10−5), 'extracellular spaces' (P=1.52×10−6), 'intrinsic components of the membrane' (P=1.29×10−8), the 'cell periphery' (P=1.76×10−7), 'integral components of the membrane' (P=1.61×10−8) and the 'plasma membrane' (P=3.22×10−7; Fig. 7). A total of 10 molecular functions from GO terms were significantly enriched for the DEGs involved in PPI networks: 'Molecular transducer activity' (P=6.13×10−14), 'signal transducer activity' (P=9.42×10−14), 'receptor activity' (P= 4.77×10−14), 'signal receptor activity' (P=4.77×10−15), 'transmembrane signal receptor activity' (P=1.31×10−15), 'serine-type peptidase activity' (P=9.57×10−5), 'endopeptidase activity' (P=0.02), 'olfactory receptor activity' (P=5.07×10−9), 'G-protein-coupled receptor activity' (P=3.03×10−14) and 'serine-type endopeptidase activity' (P=1.08×10−5) (Fig. 8).

Figure 6

Enrichment GO hierarchy of Biological Process. The GO terms enriched for the differentially expressed genes involved in the protein-protein interaction network constructed from the STRING database. GO, Gene Ontology.

Figure 7

Enrichment GO hierarchy of Cellular Component. A total of 7 Cellular Components from the GO terms were significantly enriched for the differentially expressed genes involved in the protein-protein interaction network. GO, Gene Ontology.

Figure 8

Enrichment GO hierarchy of Molecular Function. A total of 10 molecular functions from GO terms were significantly enriched for the differentially expressed genes involved in the protein-protein interaction network. GO, Gene Ontology.

RT-qPCR validation of lncRNA expression in the spinal cords 2 h post-cardiac I/R injury

To validate the reliability of the RNA sequencing results in the rats, the present study analyzed the differentially expressed lncRNAs, including 5 upregulated lncRNAs and 4 downregulated lncRNAs, by RT-qPCR analysis. T1-T4 spinal cord tissues were collected from the control and I/R groups 2 h post-reperfusion. The expression levels of 4 upregulated lncRNAs (NONRATT025386, NONRATT016113, NONRATT018298 and NONRATT018300) increased significantly in the I/R group when compared with those in the control group, whereas the expression level of one downregulated lncRNA (NONRATT002188) decreased significantly (Fig. 9). The RT-qPCR results for 3 lncRNAs (XR_589980.1, XR_598701.1 and XR_590197.1) were not consistent with the data from the RNA sequencing (XR_589980.1 decreased, and XR_598701.1 and XR_590197.1 increased post-reperfusion when compared with the control; Fig. 10).

Figure 9

Expression levels of 9 lncRNAs in the spinal cord (T1-T4) at different time-points (0.5 and 2 h) following cardiac ischemia/reperfusion injury. The expression levels of NONRATT025386, were upregulated in the 0.5 and 2 h groups. (A) The expression levels of NONRATT016113, NONRATT018298 and NONRATT018300 were only significantly upregulated in the 2 h group. (B) The expression levels of lncRNA NONRATT002188 were significantly downregulated in the 0.5 and 2 h groups. *P<0.05, **P<0.01 and ***P<0.001 vs. control. lncRNA, long non-coding RNA.

Figure 10

Validation of the differential expression of lncRNAs in spinal cord sections at different time-points (0.5 and 2 h) following cardiac isch-emia/reperfusion injury. A total of 7 upregulated lncRNAs and 2 downregulated lncRNAs were confirmed by RT-qPCR. The levels of lncRNAs were normalized to GAPDH and expressed as fold-changes compared with the sham group. RT-qPCR, reverse transcription-quantitative polymerase chain reaction; lncRNA, long non-coding RNA.

Expression levels of 4 lncRNAs in the spinal cord 0.5 h post-cardiac I/R injury

The present results indicated that the expression levels of the lncRNA NONRATT025386 were significantly upregulated in the 0.5 h group when compared with the control group, whereas the expression levels of the lncRNAs NONRATT002188 and XR_590197.1 were signifi-cantly downregulated in the 0.5 h group compared with the control group. Furthermore, the lncRNA NONRATT016113 showed no significant difference between the two groups (Figs. 10 and 11).

Figure 11

Reverse transcription-quantitative polymerase chain reaction confirmation of 3 deregulated lncRNAs in the rat spinal cord sections in the ischemia/reperfusion injury model. The expression levels of the lncRNA (A) NONRATT025386 was significantly upregulated in the 0.5 h group and the (B) lncRNA NONRATT016113 showed no significance between the two groups. The expression levels of the lncRNAs (C) NONRATT002188 and (D) gi|672 036840|ref|XR_590197.1| were significantly downregulated. Mann-Whitney test. *P<0.05 and **P<0.01 vs. control. lncRNA, long non-coding RNA.

Expression levels of 9 lncRNAs in the spinal cord at different time-points (0.5 and 2 h) following cardiac I/R injury

The present study collected samples from spinal cord tissues at 0.5 and 2 h post-cardiac I/R injury for RT-qPCR validation. The results indicated that the expression levels of the lncRNA NONRATT025386 were significantly upregulated in the 0.5 and 2 h groups compared with the control group. In addition, the expression levels of NONRATT016113, NONRATT018298 and NONRATT018300 were also significantly increased in the 2 h group, but not at 0.5 h, compared with the control group (Figs. 9–11). By contrast, the expression levels of the lncRNA NONRATT002188 were significantly downregulated in the 0.5 and 2 h groups when compared with the control group (Figs. 9–11).

Discussion

With regard to ischemic cardiac tissues, previous research has focused on multiple signaling pathways that regulate the critical balance between cell death and survival during cardiac I/R injury (42-47). The present study, to the best of our knowledge, for the first time provides evidence that suggests that many DEGs, pathways and biological processes of the T1-T4 spinal cord are implicated in myocardial I/R. Based on high-throughput RNA seq, 16,987 lncRNAs in the T1-T4 spinal cord tissues were identified, of which 3,621 were upregulated and 13,366 were downregulated (>2.0-fold-change; P<0.05). Among the 26,466 mRNAs that were quantified using RPKM values, 3,428 were deregulated by 2-fold following I/R-induced cardiac injury, of which 767 were upregulated and 2,661 were downregulated. According to these results, some differentially expressed lncRNAs were verified by RT-qPCR analyses.

Previous studies have shown that the spinal cord serves an important role in the pathogenesis of cardiac disease (48-51). The present study used a virally mediated trans-synaptic tracing method, by injecting the PRV virus into the rat heart and kidney, and these viruses were subsequently found in the lateral medial column of the spinal cord in the corresponding segment (12,52). This revealed the characteristics of the transcriptome in the T1-T4 spinal cord following cardiac I/R injury and the specific spinal segment that innervates the heart serves a significant role in cardiac I/R injury. Cheng et al (53) reported that melatonin regulation of the transcriptome was associated with the reversal of morphine tolerance by transcriptomic analysis of the L5-S3 segmental spinal cord. Mohrman et al (54) revealed the spinal cord transcriptomic and metabolic patterns in a excitotoxic injection injury model of syringomyelia. Niu et al (55) demonstrated the upregulation of tumor necrosis factor (TNF)-α in spinal cord neurons during coronary artery occlusion in rats, suggesting that TNF-α in the spinal cord may be associated with the nociception initiated by acute myocardial ischemia/infarction. Schulz et al (56) reported that connexin 43 in the spinal cord serves an important role in providing protection from cardiac I/R injuries. It is well known that myocardial ischemia creates an autonomic nervous system imbalance, and can trigger cardiac arrhythmias (57,58). Howard-Quijano et al (59) indicated that neuromodulation by spinal cord stimulation (SCS), in which a 4-pole lead was placed percutaneously in the T1-T4 epidural space, attenuated local cardiac sympathoexcitation from ischemia-induced increases in afferent signaling, reduced ventricular arrhythmias and improved myocardial function during acute ischemia. Liao et al (60) reported that chronic thoracic SCS at the T1-T3 level induced significant remodeling of cardiac sympathetic innervation over the peri-infarct and infarct regions, and was associated with improved left ventricular function and reduced myocardial norepinephrine spillover. The present study revealed that the differential expression of certain mRNAs and lncRNAs in the spinal cord affects the myocardial ischemic regions, suggesting that in the spinal T1-T4 segment these genes are involved in the response to cardiac injury. Although the functions of mRNAs and lncRNAs in the spinal cord are unclear, the present findings provide a novel paradigm for cardioprotection against I/R-induced myocardial injury.

The present results also revealed that some mRNAs in the spinal cord, including chemokine C-X-C motif ligand 1 (CXCL1), regulatory factor X4 (RFX4), WNT10a and interleukin (IL)-6, were differentially expressed following cardiac I/R. Haider et al (61) reported that the angiogenic potential of the mononuclear cell (MNC) secretome is regulated by CXCL-1 upregulation in spinal cord tissue, and factors in the MNC secretome may mitigate the pathophysiological processes of secondary damage following spinal cord injury, and may also improve functional outcomes in rats. Ashique et al (62) reported that the spinal cords of RFX4 mutants were correlated with defects in patterning and cilia formation, suggesting that RFX4 is a regionally specific transcriptional regulator of Sonic hedgehog signaling during the development of the central nervous system. Zhao et al (63) demonstrated that hyperbaric oxygen (HBO) reverses Wnt-10a upregulation induced by chronic constriction (CCI) injury in the dorsal root ganglion (DRG), spinal cord and hippocampus, suggesting that HBO attenuates CCI-induced rat neuropathic pain and inflammatory responses, potentially through regulation of the Kindlin-1/Wnt-10a signaling pathway. Although the detailed functions of many mRNAs from spinal cords are not fully understood, the present results provide novel insight into the molecular mechanisms underlying cardiac I/R injury.

It is well known that regulations of gene expression are varied over the time course (64-66). Our previous research associated with Itchy E3 ubiquitin protein ligase demonstrated that gene expression was significantly different in the C5-C8 spinal cord at 0.5 and 2 h following compound 48/80 injection when compared with the control group (36). Similar to the above method, we also screened key genes in myocardial tissues under 30 min cardiac ischemia following 2 h reperfusion compared with the sham group (67). Li et al (68) also indicated that impairment of sensory nerves with significant reductions in CGRP and SP in the DRG, ventricular myocardium and serum may be associated with an increase in myocardial vulnerability in acute cardiac I/R injury in diabetic rats. It was revealed that the injury was relatively evident at 2 h post-myocardial reperfusion, which may be considered as an acute cardiac I/R injury. Notably, this time point could be equivalent to the patients who received early percutaneous coronary intervention following myocardial infarction (69,70). Therefore, the present study chose 0.5 and 2 h post-reperfusion as the different time points of cardiac reperfusion injury. The present study revealed that there are significant differences in lncRNA and mRNA expression patterns at different time points of myocardial I/R, suggesting that the neural modulation of cardiac I/R injury may be temporal- and spatial-dependent.

In particular, our previous study also demonstrated that proton magnetic resonance spectroscopy was able to simultaneously detect and quantify the absolute concentrations of multiple metabolites within the spinal cord underlying α-Me-5-HT-evoked pruritus (71). Using the above method, we can also detect the changes of various metabolites in the spinal cord following cardiac I/R injury, which may deepen our understandings of the pathophysiology and pharmacological therapies for acute myocardial infarction.

The present study screened several differentially expressed mRNAs under cardiac I/R injury. In the future, whether the proteins encoded by the mRNAs are consistent with these mRNAs will be verified by immunoblotting. If so, the effects on cardiac I/R injury may be observed by activating or silencing the expression of one specific mRNA. Li et al (68) revealed that the regulation of the spinal cord served a significant role under cardiac I/R injury in diabetic neuropathic rats. Thus, it can be hypothesized that intervention on the spinal cord may have an important influence on cardioprotection in the future.

The present study used high-throughput RNA seq, coupled with RT-qPCR analysis, to demonstrate that the expression profiles of lncRNAs and mRNAs in spinal cords differed markedly between the control and 2 h groups, and ultimately identified 7,980 differentially expressed (>2-fold) lncRNAs (234 upregulated and 7,746 downregulated) and 3,428 mRNAs (767 upregulated and 2,661 downregulated). The expression patterns of several lncRNAs were confirmed by RT-qPCR. The results also indicated that the expression levels of the lncRNA NONRATT025386 were significantly upregulated in the 0.5 and 2 h groups compared with the control group, whereas the expression levels of NONRATT016113, NONRATT018298 and NONRATT018300 were significantly increased in the 2 h group when compared with the control group, although there was no statistically significant difference between the expression levels in the 0.5 h and control groups. Furthermore, the expression levels of the lncRNA NONRATT002188 were significantly downregulated in the 0.5 and 2 h groups compared with the control group.

In conclusion, this study revealed that high-throughput RNA seq can facilitate the systematic exploration of gene expression on a genome-wide scale, and can be used to investigate DEGs and lncRNA expression patterns in the spinal cords of rats during I/R-induced cardiac injury. In the search for better treatments for cardiac I/R injury, expanded sets of differentially expressed mRNAs and lncRNAs may prove very useful for identifying novel therapeutic targets.