Altered expression of differential gene and lncRNA in the lower thoracic spinal cord on different time courses of experimental obstructive jaundice model accompanied with altered peripheral nociception in rats

The spinal origin of jaundice-induced altered peripheral nociceptive response poorly understood. In the current study, we aimed to first validate rats with bile duct ligation (BDL) as a jaundice model accompanied by altered peripheral nociceptive response, and then to analyze differential gene and lncRNA expression patterns in the lower thoracic spinal cord on different time courses after BDL operation by using high-throughput RNA sequencing. The differentially expressed genes (DEGs) identified using reverse transcription-quantitative polymerase chain reaction (RT-qPCR) analysis, followed by clustering analysis, Gene Ontology analysis and pathway analysis. As a result, a total of 2033 lncRNAs were differentially expressed 28d after BDL, in which 1545 probe sets were up-regulated and 488 probe sets were down-regulated, whereas a total of 2800 mRNAs were differentially expressed, in which 1548 probe sets were up-regulated and 1252 probe sets were down-regulated. The RNAseq data of select mRNAs and lncRNAs was validated by RT-qPCR. 28d after BDL, the expressions of lncRNA NONRATT002335 and NONRATT018085 were significantly up-regulated whereas the expression of lncRNA NONRATT025415, NONRATT025388 and NONRATT025409 was significantly down-regulated. 14d after BDL, the expressions of lncRNA NONRATT002335 and NONRATT018085 were significantly up-regulated; the expression of lncRNA NONRATT025415, NONRATT025388 and NONRATT025409 was significantly down-regulated. In conclusion, the present study showed that jaundice accompanied with decreased peripheral nociception involved in the changes of gene and lncRNA expression profiles in spinal cord. These findings extend current understanding of spinal mechanism for obstructive jaundice accompanied by decreased peripheral nociception.

The introduction of high-throughput RNA sequencing (RNA-seq), where thousands of differential genes and lncRNA expression patterns can be studied in parallel [21][22][23][24], permits a broad assessment for altered expression of differential gene and lncRNA in the lower thoracic spinal cord on different time courses of obstructive jaundice model. In current study, we performed a comprehensive transcriptome analysis in obstructive jaundice model using RNA-seq, and identified lncRNAs with differential expression. To explore the function of lncRNAs, we predicted their potential targets with cis-regulatory effects, which were then put into gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) for further analysis. In addition, we compared the expression of some lncRNAs in thoracic spinal cord at different time points in obstructive jaundice model. group showed yellow ears (D), liver tissue(E) and rear claw(F) 28 days after surgery compared to control group (A, B, C). (G) The serum total bilirubin levels was significantly higher in the BDL group (112.8±12.88μml/L) than in the control group (0.45±0.06μml/L). (H) The mechanical PWT was tested to assess the jaundice-induced nociceptive threshold change prior to (day 0) and at days 3, 7, 14, 21, and 28 following BDL operation. The hindpaw of BDL rats displayed a dramatic decrease in mechanical PWT to von Frey filament stimulation from 4.28±0.59 g at baseline to 7.15 ±0.43 g 14d after BDL operation and further deteriorated to 7.02 ± 0.57 at day 21, and 7.76 ±0.58g at day 28. * P<0.05, ** P<0.01, *** P<0.001 compared with the each corresponding time point. www.impactjournals.com/oncotarget , whose expression changes were more than twofold. In clustering analysis, up-and down-regulated genes are colored in red and green, respectively. www.impactjournals.com/oncotarget

Evaluation of jaundice
We observed the general appearance of the animals used in this study (Figure 1), and found that rats in BDL group showed yellow ears (D), liver tissue(E) and rear claw(F) 28 days after surgery compared to control group (A, B, C). Otherwise, our results showed that the serum total bilirubin levels was significantly higher in the BDL group (112.8±12.88μml/L) than in the control group (0.45±0.06μml/L) ( Figure 1G).

The change of nociceptive threshold following BDL operation
The mechanical nociceptive stimulation was tested to assess the jaundice-induced nociceptive threshold change prior to (day 0) and at days 3, 7, 14, 21, and 28 following BDL operation. Our result showed that the hindpaw of BDL rats displayed a dramatic increase in mechanical nociceptive threshold to von Frey filament stimulation from 4.28±0.59 g at baseline to 7.15 ±0.43 g 14d after BDL operation and further deteriorated to 7.02 ± 0.57 at 21d and 7.76 ±0.58g at 28d after BDL operation ( Figure 1H). These data indicated that obstructive jaundice resulted in the decreased sensitivity in response to the mechanical nociceptive stimulation in the present study.

Differentially expressed genes of T6-T12 spinal cord 28d after BDL operation
To systematically identify jaundice-associated lncRNAs in spinal cord, the differential expression analysis was performed between BDL group and control group. The lncRNA and mRNA expressions of T6-T12 spinal cord of the animals were examined the HiSeq 2000 system (Illumina, Inc.) with a total of 1,198,903,526 raw reads from the nine libraries. Differentially expressed gene profiles in BDL group were compared to the corresponding data from control group. Figure 2 indicated the scatter plot comparing global mRNA (A) gene expression profiles in the spinal cord between BDL group (Model group) and sham group (Control group), and heat map showing hierarchical clustering of overall mRNA (B) expression pattern of reliably measured probe sets. In addition, heat map also showed hierarchical clustering of mRNA (C), whose expression changes were more than twofold.   We identified that a total of 2033 lncRNAs were differentially expressed between control group and BDL group, in which 1545 probe sets were up-regulated and 488 probe sets were down-regulated, whereas a total of 2800 mRNAs were differentially expressed between control group and BDL group, in which 1548 probe sets were upregulated and 1252 probe sets were down-regulated. The detail information of the top 20 up-regulated and 20 downregulated lncRNAs was in Table 1, and that of the top 20 up-regulated and 20 down-regulated mRNAs was in Table  2. The maximal and minimal fold change was 35.01 and 2.01, respectively.

Gene ontology annotation for differential expression genes
We had used the RNA-seq analyses to identify the differentially expressed genes which were annotated using the GO database (Gene Ontology, http://www. geneontology.org/). The present study showed three important results involved in biological functional groups including molecular function ( Figure 3A), biological process ( Figure 3B), and cellular component ( Figure  3C). Figure 3D indicated that the differential expression genes were analyzed with GO background significant enrichment, and the differentially expressed mRNAs in spinal cord were primarily involved in the biological processes GO functions (Figure 4).

Validation of the lncRNAs by RT-qPCR 28d after BDL operation
To validate the RNA-seq data, top 5 up-regulated lncRNAs (NONRATT025327, NONRATT000845, NONRATT001654, NONRATT002335, and NONRATT018085) and top 5 down-regulated lncRNAs (NONRATT025415, NONRATT025388, NONRATT025409, NONRATT02538 and NONRATT006517) that were differentially expressed between BDL group and control group were randomly selected (Table 3), and their relative expression levels were quantified by RT-qPCR. We found that the expressions of lncRNA NONRATT002335 and NONRATT018085 were significantly up-regulated underlying BDL group compared with control group, whereas the expression of lncRNA NONRATT025415, NONRATT025388 and NONRATT025409 was significantly down-regulated in model group ( Figure 5).

The expression of 10 lncRNAs in T6-T12 spinal cord at different time points (14d/28d) after BDL operation
Obviously, gene expressions after BDL operation injection are varied in different time points. We collected spinal tissue sample from T6-T12 spinal cord 14d vs 28d after BDL operation for RT-qPCR   and NONRATT018085 were significantly up-regulated underlying BDL group compared with control group; the expression of lncRNA NONRATT025415, NONRATT025388 and NONRATT025409 was significantly down-regulated underlying BDL group compared with control group whereas the expression of lncRNA NONRATT02538 and NONRATT006517 had not statistically different between control group and BDL group.
validation. 14d after BDL operation, the expressions of lncRNA NONRATT002335 and NONRATT018085 were significantly up-regulated underlying BDL group compared with control group whereas the expression of lncRNA NONRATT025327, NONRATT000845 and NONRATT001654 had not statistically different between control group and BDL group; The expression of lncRNA NONRATT025415, NONRATT025388 and NONRATT025409 was significantly down-regulated underlying BDL group compared with control group whereas the expression of lncRNA NONRATT02538 and NONRATT006517 had not statistically different between control group and BDL group (Figures 6 and 7).

DISCUSSION
In the present study, we confirmed that bile duct ligation resulted in obstructive cholestasis and jaundice and this model was consistent with previous reports describing marked elevation of the serum total bilirubin www.impactjournals.com/oncotarget in rats responding to bile duct ligation [5,25]. In addition, our result also indicated that obstructive jaundice induced the decreased sensitivity in response to the mechanical nociceptive stimulation 14d and 28d after BDL operation, which was in line with a previous experimental study showing that obstructive cholestasis by bile duct resection in rodents displayed decreased nociception [26][27][28].
We identified the spinal genes and lncRNAs that were differentially expressed between model group and control group by high throughput RNA-seq, providing an important view of spinal genetic heterogeneity for jaundice-induced decreased nociception in rats. The present study identified 2033 differentially expressed lncRNAs in the spinal cord, 488 of which were down- regulated >2-fold whereas 1545 of which were upregulated >2-fold. Similarly, a total of 2800 mRNAs were differentially expressed between control group and BDL group, 1252 of which were down-regulated whereas 1548 of which were up-regulated (FDR<0.05, |log2FC|>1). Among top 5 up-regulated lncRNAs and top 5 downregulated lncRNAs, 2 up-regulated (NONRATT002335 and NONRATT018085) and 3 down-regulated lncRNAs (NONRATT025415, NONRATT025388 and NONRATT025409) were identified. These results suggested that 50% differentially expressed lncRNAs were verified by RT-qPCR. Furthermore, we found that Serpina3n was an important up-regulated gene 28d after BDL operation, which were in agreement with a previous study that the serine protease inhibitor Serpina3n was upregulated in the dorsal root ganglia (DRG) after nerve injury, and attenuated neuropathic pain by inhibiting T cell-derived leukocyte elastase [29].
Previous studies have suggested that the different duration of cholestasis accompanies by some behavioral change in rodents [26,30]. Belghiti et al observed that chronic BDL rats displayed enhanced scratching behavior and thermal hyperalgesia already 48 h after BDL surgery, and this result was perfectly demonstrated by the upregulation and sensitization of the heat-sensitive TRPV1 channel [30]. Tian et al also demonstrated 5-HT-induced enhanced scratching and antinociception to mechanical and heat stimuli in BDL rats [31], suggesting that targeting 5-HT receptors may be an effective treatment for cholestatic itch. In this study, we found differential expression of some lncRNAs (Figure 7) in the spinal cord at different time points (14d/28d) after BDL operation, suggesting that the different lncRNAs in spinal cord segment may be involved in the neuronal response to obstructive jaundice. Although the functions of many lncRNAs in spinal cord are not fully known, our findings provide novel potential insights involving in the molecular mechanism of jaundice-induced altered peripheral nociception.
In conclusion, obstructive jaundice accompanied by altered peripheral nociception is always a crucial factor limiting the therapeutic efficacy of many drugs and quality of life for patients. Here we constructed the expression profiles of lncRNAs and potentially related mRNAs in rats with obstructive jaundice, and found some distinct lncRNA/mRNA expression profiles in spinal cord, suggesting that these unique noncoding transcripts may contribute to the acquisition of cholestasis-induced altered peripheral nociception. Although additional studies are needed to verify these lncRNAs/mRNAs mentioned above, our study provides important insights into novel indicators of treatment for patients with cholestasis.

Animal care
Male Sprague-Dawley rats were provided by the Center of Experimental Animal of Tongji Medical College (license number: 43004700019962). All experimental protocols, animal testing and surgeries were performed in adherence to the National Institute of Health Guide for the Care and Use of Laboratory Animals (NIH Publications No. 80-23) revised 1996. The experimental protocols were approved by the committee of experimental animals of Tongji Hospital, Tongji Medical College (IRBID:TJ-A0804). The animals were maintained in a climate controlled room on a 12-hlight/dark cycle (light on at 07:00 h). Rats were housed (2/cage), but they were

Obstructive jaundice model by bile duct ligation
Rats were placed in a temperature-controlled chamber before experimental operation. All surgical procedures were performed under sterile conditions. Common bile duct ligation (BDL) established the prolonged obstructive jaundice as described in previous studies [32][33][34]. After an intraperitoneal injection of ketamine hydrochloride (50 mg/kg) plus xylazine (5 mg/ kg), rats which underwent BDL were anesthetized. In those with BDL, the common bile duct was located and ligated using 4-0 silk at two points anterior to the pancreas and posterior to the hilum of the liver. The first ligation was made just above the duodenum and the second almost 2 mm above the first ligation. The bile duct was then transected at the midpoint between the two ligatures. Animals used in this study were sacrificed before signs of severe illness became apparent.

Experimental groups
Experiment A Rats were randomly assigned to two groups: (1) Control group (sham-ligation surgery, n = 9); (2) BDL group (common bile duct ligation, n = 9). The mechanical nociceptive thresholds were evaluated at different time points including day (d) 0 (baseline), 3, 7, 14, 21, and 28. 28d after operation, rats were sacrificed, and blood samples, spinal cord tissues were collected for further analysis. Serum total and direct bilirubin levels were measured. Thoracic segments of spinal cord (T6-T12) were prepared to analyze differential gene and lncRNA expression patterns by high-throughput RNA sequencing, and thoracic (T6-T12) segment of spinal cord for Real-Time quantitative PCR (RT-qPCR).
Experiment B Rats were randomly assigned to two groups: (1) Control group (n = 12); (2) BDL group (n = 12). 14 days after surgery the T6-T12 spinal cord were collected using a dissection microscope, rinsed with isotonic saline, dissected and fleshly frozen in liquid nitrogen for RT-qPCR.

Assessment of mechanical sensitivity
As in previous reports [35][36][37][38], mechanical paw withdraw threshold was examined using the blind method with a minimum of six animals/group. Calibrated von Frey filament apparatus (Stoelting, Wood Dale, IL, USA) was used to measure the paw withdrawal response for a mechanical stimulus. A range of filament comprises 1g, 1.4g, 2g, 4g, 6g, 8g, 10g, and 15g bending force. Rats were allowed to acclimate within a clear plastic chamber over an elevated mesh floor for 30 minutes before assessment at room temperature. Each filament was applied to the midplantar surface of the left hind paw until a withdrawal response occurred. To avoid tissue injury in refractory animals, stimulation was automatically terminated after 15 seconds. The interval between adjacent tests was >5 minutes. The results of 3 consecutive measurements were averaged for the mechanical paw withdraw threshold in BDL and control rats.

Analysis of jaundice index
28d following BDL operation, serum total bilirubin levels were determined in the Automated Blood Chemical Analyzer Vitro 350 (Orthoclinical Diagnostic Inc., Rochester, NY).

RNA extraction, library preparation and highthroughput RNA sequencing
RNA was isolated from the thoracic (T6-T12) spinal segment tissues in two groups using TRI reagent (Invitrogen, Carlsbad, CA, USA) according to the manufacturer's instructions. Tissue samples were determined using a NanoDrop 2000 (Thermo Fisher Scientific Inc, USA) for integrity, quality and purity.
Library preparation and high-throughput RNA sequencing were performed by the OE Biotechnology Corporation (Shanghai, China) [39][40][41]. Total RNA of all samples were sequenced by EMBL GeneCore (EMBL Genomics Core Facility, Heidelberg, Germany). The synthesis of polyadenylated transcriptome libraries for every sample was accompanied by deep sequencing in 3 lanes each generating >50 million reads pair end. RNAseq was performed in the two groups, each with three biological replicates.

Real-time quantitative-PCR
Thoracic (T6-T12) segment of spinal cord were dissected and isolated, and total RNA was extracted using TRIzol reagent (Beijing Tiangen Biotech CO.). 6 μg total RNA was used as a template for cDNA synthesis and amplification using the First-Stand Synthesis System (BioPhotometer, Eppendorf, Hamburg, Germany) according to the manufacturer's instructions [35,[42][43][44][45]. The primers were designed with the Primer Express 3.0 software (Applied Biosystems), and the specific forward (F) and reverse (R) primer sequences were in Table 4. Experiments were evaluated in triplicate and repeated at least three times. The threshold cycle (CT) was used to estimate the amount of target mRNA. The comparative CT method with the formula for relative fold-change = 2 −∆∆CT was used to quantify the amplified transcripts.

Data analysis
Results are expressed as the mean ± SEM. Statistical comparisons were performed with Mann-Whitney test. Statistical comparisons were performed using unpaired