Construction and Analysis of a Long Non-Coding RNA (lncRNA)-Associated ceRNA Network in β-Thalassemia and Hereditary Persistence of Fetal Hemoglobin.

BACKGROUND Higher fetal hemoglobin (HbF) levels can ameliorate the clinical severity of ß-thalassemia. The use of integrative strategies to combine results from gene microarray expression profiling, experimental evidence, and bioinformatics helps reveal functional long noncoding RNAs (lncRNAs) in ß-thalassemia and HbF induction. MATERIAL AND METHODS In a previous study, a microarray profiling was performed of 7 individuals with high HbF levels and 7 normal individuals. Thirteen paired samples were used for validation. lncRNA NR_001589 and uc002fcj.1 were chosen for further research. The quantitative reverse transcription-PCR was used to detect the expression levels of 2 lncRNAs. The Spearman correlation test was employed. The nuclear and cytoplasmic distribution experiment in K562 cells was used to verify the subcellular localization of 2 lncRNAs. Potential relationships among lncRNAs, predicted microRNAs (miRNAs), and target gene HBG1/2 were based on competitive endogenous RNA theory and bioinformatics analysis. RESULTS Average expression levels of NR_001589 and uc002fcj.1 were significantly higher in the high-HbF group than in the control group. A positive correlation existed between NR_001589, uc002fcj.1, and HbF. The expression of NR_001589 was in both the cytoplasm and the nucleus, mostly (77%) in the cytoplasm. The expression of uc002fcj.1 was in both the cytoplasm and the nucleus; the cytoplasmic proportion was 43% of the total amount. A triple lncRNA-miRNA-mRNA network was established. CONCLUSIONS Novel candidate genetic factors associated with the HBG1/2 expression were identified. Further functional investigation of NR_001589 and uc002fcj.1 can help deepen the understanding of molecular mechanisms in ß-thalassemia.

Background

Genome-wide sequencing shows that about 93% of the DNA sequence in the human genome is transcribed into RNA, but only about 2% of the DNA sequence eventually encodes a protein [1]. Some nonprotein-coding RNAs, formerly known as “transcriptional noise,” now are known to serve as significant regulators of target gene expression [2]. Among them, long noncoding RNAs (lncRNAs) are longer than 200 nucleotides and are involved in a variety of biological processes, such as cell proliferation, differentiation, and chromosomal variation [3]. MicroRNAs (miRNAs), 22–25 nucleotides in length, negatively regulate target genes at the post-transcriptional level and participate in hematopoietic landscape shaping [4]. lncRNAs can function in various diseases, including hematopoiesis and other blood diseases, by interacting with miRNAs [5,6]. lncRNAs act directly against miRNAs to antagonize the expression and function of miRNAs, or are degraded by miRNAs, thus affecting target lncRNAs in the pathophysiological process [7,8]. lncRNAs compete with miRNAs for direct binding to mRNAs [9]. Some lncRNAs can also cleave miRNAs from sequences of intronic or exonic region during maturation [10,11]. However, few studies have reported the important regulatory roles of lncRNAs in β-thalassemia and fetal hemoglobin (HbF) induction.

β-thalassemia is a genetic and hemolytic disease caused by the dysfunction of globin synthesis [12]. Severe β-thalassemia probably accounts for more than 50 000 deaths per year of all deaths of children in tropical and subtropical areas [13]. In Guangxi province of southern China, the mutation gene frequency of β-thalassemias is up to 6.43% [14]. HbF, composed of 2 α chains and 2 γ chains, is the major hemoglobin type during fetal life and is replaced by adult hemoglobin after birth [15]. Accumulating evidence has shown that increased HbF levels effectively ameliorate the clinical symptoms and improve the prognosis of β-thalassemia. Genetic regulation of HbF levels has been of particular therapeutic interest in recent years [16,17]. Focusing on individuals with high levels of HbF in the geographic regions where β-thalassemias are prevalent with specific molecular pathology and racial/ethnic characteristics may provide valuable insights into the mechanisms underlying the expression of HBG1/2 genes. So far, detailed studies on mRNAs and miRNAs have helped guide the diagnosis and therapy of β-thalassemia. However, few studies have been conducted on the function of lncRNAs in β-thalassemia.

In a previous study, a microarray profiling of individuals with high HbF levels and normal individuals was performed, but the lncRNA function was poorly clarified. In the present study, lncRNAs NR_001589 [18] and uc002fcj.1 were selected to explore their regulatory mechanisms based on the competitive endogenous RNA (ceRNA) theory [19]. NR_001589 was of interest because it is located upstream of the β-globin locus. The previous study suggested that NR_001589 might activate HBE1 and regulate HbF expression. Uc002fcj.1, located on chromosome 16, was the most upregulated lncRNA in the high-HbF group compared with the normal group. The quantitative reverse transcription-polymerase chain reaction (qRT-PCR) was performed for validating these 2 differentially expressed lncRNAs. Their subcellular localization was confirmed using K562 cell lines. Putative miRNA-HBG1/2 sites in lncRNAs were predicted by bioinformatics analysis. A triple lncRNA-miRNA-mRNA network was established. The findings of this study offer new insights into the role of lncRNAs in HbF induction in patients with β-hemoglobinopathies, although deeper explorations are needed on this novel regulatory mechanism.

Material and Methods

Study participants and microarray analysis

The details are available in Reference [18]. Thirteen paired samples (13 subjects in the high-HbF group and 13 subjects in the control group) were used for validation. This study approved by the First Affiliated Hospital of Guangxi Medical University (2013-KY-007).

RNA extraction from nucleated erythrocytes and reticulocytes

Isolation of nucleated red blood cells an reticulocytes was shown in our previous study [18]. Total RNA was extracted from reticulocytes using TRIzol (Invitrogen Life Technologies, USA) in accordance with the manufacturer’s protocol. The quantity and quality of the total RNA were assessed using a NanoDrop ND1000 spectrophotometer (NanoDrop, USA).

qRT-PCR validation of differentially expressed lncRNAs: NR_001589 and uc002fcj.1

Based on a previous study [18], qRT-PCR was performed to further confirm whether lncRNA NR_001598 and uc002fcj.c had differential expression. RNA was reverse-transcribed into cDNA using SuperScript III Reverse Transcriptase (Invitrogen Life Technologies, CA, USA) according to the manufacturer’s protocols. Table 1 presents the sequences of qRT-PCR primers used. β-actin was used as the control gene.

Table 1

Sequences of the qRT-PCR primers.

Primer name	Sequence (5′–3′)		Annealing temperature (°C)	Product length (bp)
uc002fcj.1	Forward	GTCTGGTCCTTCTTACTGATGTG	60	159
	Reverse	CTTTCTTTTCGGTGTTCCTG
NR_001589	Forward	TTCAGCGAGGATTTTACCC	60	90
	Reverse	CTCAGTGGTCTTGTGGGCTA
β-actin (H)	Forward	GTGGCCGAGGACTTTGATTG	60	73
	Reverse	CCTGTAACAACGCATCTCATATT
U48	Forward	GATGATGACCCCAGGTAACTCT	60	50
	Reverse	TGCGGTGATGGCATCAGCGACAC

Subcellular localization of lncRNAs: NR_001589 and uc002fcj.1

Cell culture

K562 cell is a widely used human erythroid-like cell line capable of undergoing erythroid differentiation [20]. Numerous studies have used K562 cells to elucidate the regulatory mechanism of HbF expression in β-thalassemia in vitro [21]. K562 cells were purchased from the Stem Cell Bank of the Chinese Academy of Sciences (Shanghai, China) and cultured in Roswell Park Memorial Institute (RPMI) 1640 medium (Sigma-Aldrich, MO, USA) with 10% fetal bovine serum (Gibco, South America), 100 IU/mL penicillin, and 100 μg/mL streptomycin (Solarbio, China) in a 5% CO2 humidified atmosphere.

Nuclear and cytoplasmic separation experiments

Frozen K562 cells at −80°C were slowly thawed on ice and centrifuged at 500 g for 5 min to collect cells. The cells were washed by adding 500 μL of 1×PBS (phosphate-buffered saline) and collected by centrifugation at 500 g for 5 min. Then, 20 volumes of cell lysis buffer were added to the cell pellet, mixed well, and placed on ice for 5 min. After centrifugation at 1500 g for 5 min, the supernatant was collected as a cytoplasmic crude extract. Attempts were made to remove the supernatant. An equal volume of cell lysis buffer was added, mixed well, and placed on ice for 10 min. After centrifugation at 1500 g for 5 min, the precipitate comprised the separated nucleus. The cytoplasmic crude extract was centrifuged at 16 000 g for 5 min, and the supernatant was finally isolated as a cytoplasmic fraction.

RNA extraction from the nucleus and the cytoplasm

Nuclear and cytoplasmic RNAs of K562 cells were extracted separately using TRIzol (Invitrogen Life Technologies) according to the manufacturer’s protocol. The quantity and quality of the extracted RNA were tested on a NanoDrop ND-1000 spectrophotometer (NanoDrop, NY, USA). Denaturing agarose gel electrophoresis was used to assess the integrity of the RNA. Total RNA was reverse-transcribed into cDNA using SuperScript III Reverse Transcriptase (Invitrogen, NY, USA) in accordance with the manufacturer’s instructions. The amount of input RNA used was 500 ng, and the final volume for all reactions was adjusted to 20 μL with ddH2O. cDNA was stored at −20°C overnight and then used for qRT-PCR.

qRT-PCR validation of subcellular localization for NR_001589 and uc002fcj.1

A qRT-PCR was performed using the ViiA 7 Real-Time PCR System (ABI, NY, USA). A reaction volume of 10 μL was mixed, consisting of 5 μL of 2×Master Mix (ArrayStar, MD, USA), 0.5 μL of PCR Forward Primer, 0.5 μL of PCR Reverse Primer, 2 μL of template cDNA, and 2 μL of double-distilled water. The following cycling conditions were applied: 95°C for 10 min followed by 40 cycles of 95°C (10 s) and 60°C (60 s). The lncRNA PCR results were quantified using the 2−ΔΔct method, with normalization using β-actin and U6.

Statistical analysis

All statistical data were analyzed with SPSS 20.0 software (SPSS, Inc., IL, USA). Data are shown as the mean ± standard deviation. The t test was used to analyze the statistical significance of the microarray and qRT-PCR results. The Spearman correlation coefficient analysis was performed to assess correlations between the levels of lncRNAs verified by qRT-PCR and HbF levels. Statistical differences were considered significant at P<0.05.

Bioinformatics analysis

The HBG1/2 mRNA sequences and lncRNA sequences were obtained from the UCSC (). The predicted miRNA target sites of NR_001589 were collected from miRcode () and RegRNA2.0 (). The predicted miRNA targets sites of uc002fcj.1 were collected from miRcode () and TargetScan (). The interactions of miRNA–HBG1/2 in the present study were acquired from 5 highly reliable miRNA-mRNA reference databases: miRcode, DIANA Tools (), miRDB (), microRNA.org, and TargetScan.

Results

Validation of dysregulated lncRNAs

The results of qRT-PCR showed that lncRNAs NR_001589 and uc002fcj.1 were all upregulated in the high-HbF group compared with the control group (Figure 1).

Figure 1

Validation of lncRNA using qRT-PCR. lncRNA levels were measured in peripheral reticulocytes derived from 13 high-HbF patients and 13 age- and sex-matched controls. The relative expression level of each RNA was normalized. Data displayed in histograms are expressed as means ± standard deviation. ** P<0.01; *** P<0.001, compared with the high-HbF and normal-HbF groups. The qRT-PCR data of NR_001589 were derived from a previous study [18].

Correlation analysis between lncRNAs and HbF levels

Significant positive correlations were observed between NR_001589 and HbF (r=0.501, P=0.009) and uc002fcj.1 and HbF (r=0.779, P=0.000) (Figure 2).

Figure 2

Associations between lncRNA and HbF levels. Spearman correlation coefficient analysis was used to evaluate associations between (A) NR_001589 and HbF and (B) uc002fcj.1 and HbF. The correlation data of NR_001589 and HbF were derived from a previous study [18].

Subcellular localization of lncRNAs

The subcellular distribution of lncRNA determines its possible ways of functioning. Subcellular localization in K562 cells is necessary for subsequent mechanistic studies. The expression of NR_001589 was found in both cytoplasm and nucleus, mostly (77%) in the cytoplasm. The expression of uc002fcj.1 was seen in both the cytoplasm and the nucleus; the cytoplasmic proportion was 43% of the total amount (Figure 3).

Figure 3

Relative expression of lncRNA in the cytoplasm and the nucleus of K562 cells. The expression of uc002fcj.1 was found in both the cytoplasm and the nucleus; the cytoplasmic proportion is 43%. The expression of NR_001589 was seen mostly (77%) in the cytoplasm. β-actin, a cytoplasmic marker; U48, a nuclear marker.

Establishment of the lncRNA-miRNA-mRNA network

A lncRNA-associated ceRNA network was constructed by combining lncRNA-miRNA interactions and miRNA-HBG1/2 interactions. The network was visualized, and was composed of 2 lncRNA nodes, 2 mRNA nodes, and 14 miRNA nodes (Figure 4). Table 2 presents the specific miRNAs and source online platforms.

Figure 4

lncRNA-miRNA-mRNA ceRNA network. Differentially expressed lncRNA-mediated ceRNA network. Nodes in red represent miRNA. Rectangles in yellow represent HBG1/2. V shape in blue represents lncRNA. All lines denote interactions. Inhibition is denoted.

Table 2

Prediction and online platforms.

NR_001589	miRDB	miR-3132	HBG1	miRDB
	miRDB	miR-4291		miRDB, TargetScan
	Regrna2.0	miR-3619-5p		TargetScan, miRcode
	miRcode	miR-214/761		miRcode
	miRcode	miR-148ab-3p/152		miRcode
	miRcode	miR-18ab/4735-3p		miRcode
	miRcode	miR-137		microRNA.org, miRcode, DIANA Tools
	miRcode	miR-27abc/27a-3p		miRcode
	miRcode	miR-203		miRcode
	miRcode	miR-27abc/27a-3p	HBG2	miRcode
	miRDB	miR-3132		miRDB
	miRcode	miR-137		microRNA.org, miRcode
	Regrna2.0	miR-3619-5p		miRcode
	miRcode	miR-18ab/4735-3p		miRcode
	miRcode	miR-148ab-3p/152		miRcode
	miRDB	miR-4291		miRDB
	miRcode	miR-214/762		miRcode
	miRcode	miR-203		miRcode
uc002fcj.1	miRcode	miR-193/193b/193a-3p	HBG1	miRcode
	miRcode	miR-3619-5p		TargetScan, miRcode
	miRcode	miR-214/761		miRcode
	miRcode	miR-338/338-3p		miRcode
	miRcode	miR-34ac/bc-5p/449abc/449c-5p		miRcode
	miRcode	miR-27abc/27a-3p		miRcode
	miRcode	miR-128/128ab		miRcode
	miRcode	miR-148ab-3p/152		miRcode
	miRcode	miR-203		miRcode
	miRcode	miR-128/128ab	HBG2	miRcode
	miRcode	miR-148ab-3p/152		miRcode
	miRcode	miR-27abc/27a-3p		miRcode
	miRcode	miR-34ac/bc-5p/449abc/449c-5p		miRcode
	miRcode	miR-338/338-3p		miRcode
	miRcode	miR-214/761		miRcode
	miRcode	miR-3619-5p		miRcode
	miRcode	miR-193/193b/193a-3p		miRcode
	miRcode	miR-203		miRcode

Discussion

Great efforts have been made to elucidate the molecular mechanism underlying β-thalassemia. Previous studies focused mainly on mRNA and miRNAs. HBS1L-MYB, BCL11A, and KLF1 regulate γ-globin gene (HBG1/2) expression and influence HbF levels [22-24]. Additional, several miRNAs have been identified as critical factors regulating HbF expression, such as miR-15a, miR-16-1, miR-96, miR-210, miR-221, miR-222, miR-486-3p, and the let-7 family [25-29]. Accumulating evidence suggest the roles of lncRNAs in a variety of biological processes. Dysregulation of lncRNA has been found in genetic diseases, including hematopoiesis and the pathogenesis of blood diseases [30,31]. Studying the relationship of lncRNAs with miRNAs and/or mRNAs, whose functions have been annotated, might help infer the potential functions of lncRNAs. Reportedly, lncRNA has a natural “sponge” role as a ceRNA, thus affecting the inhibitory effects of miRNAs on target genes [32]. miRNAs regulate lncRNAs through similar interactions with the highly conserved region of lncRNAs and vice versa [33]. Therefore, it is crucial to learn the regulatory role of lncRNAs and their functional relationship with miRNAs as ceRNA in β-thalassemia and HbF induction.

This novel study confirmed that NR_001589 and uc002fcj.1 were significantly upregulated in the high-HbF group. The interplay data from databases and a previous study were combined to generate a triple network based on the ceRNA theory. Based on the results (Figures 1, 2, and 4), it was hypothesized that NR_001589 and uc002fcj.1 could interact with miRNAs and alter the expression of γ-globin gene. The miRcode, RegRNA2.0, and TargetScan were used to obtain NR_001589- and uc002fcj.1-targeting miRNAs, so as to find more relevant miRNAs and their potential regulatory role in HbF induction. The results showed that these miRNAs also interacted with HBG1/2. The miRNAs related to NR_001589, uc002fcj.1, and HBG1/2 (miR-3619-5p and miR-137) gained attention.

miR-3619-5p has been proved to be a cancer suppressor in prostate cancer and non-small cell lung cancer. It is associated with proliferation, invasion, and autophagy [34-36]. Bioinformatics databases (TargetScan and miRcode; Table 2) found binding sites for HBG1 and miR-3619-5p. Currently, whether miR-3619-5p is involved in HbF regulation is unknown. The association between miR-3619-5p and HBG1 warrants further investigation. The biological roles of miR-137 in cell proliferation, migration, invasion, and apoptosis have been reported. miR-137 is also involved in human cord blood–derived CD34+ cell erythropoiesis [37]. Complementary sequences of HBG1 and miR-137 were detected by the bioinformatics software (microRNA.org, miRcode, and DIANA Tools; Table 2). The regulatory effect of miR-137 in HbF induction needs further investigation.

This novel study detected the subcellular distributions of NR_001589 and uc002fcj.1 in K562 cells (Figure 3) and the advantage of lncRNA analysis. Nuclear and cytoplasmic lncRNAs can regulate gene expression in different ways [38]. Intranuclear lncRNAs bind the transcription factors and recruit related proteins. Histone trimethylation is induced, and the expression of nearby gene mRNAs is regulated. In addition, lncRNAs can directly bind with the promoter to regulate gene expression [39]. The ceRNA theory indicates that all types of RNA transcripts can crosstalk with each other through miRNA-binding sites. A recent study showed that cytoplasmic lncRNAs can serve as ceRNAs, function as precursors of miRNAs, and participate in mRNA and protein modifications [40]. Based on these results, it was speculated that lncRNA NR_001589, distributed mainly in the cytoplasm, might function as a ceRNA by sponging some microRNAs (including miR137), affecting the expression of HBG1/2. However, lncRNA uc002fcj.1, distributed in both the nucleus and the cytoplasm, affects transcription of the HBG1/2 gene and also influences post-transcriptional modification, thereby affecting HbF levels. All these topics need further exploration.

The present study has some limitations. First, the available microarray data on β-thalassemia and hereditary persistence of fetal hemoglobin (HPFH) was lacking. Second, the lncRNA microarray research is still in its infancy compared with mRNA and miRNA microarray testing. Finally, further experimental studies should be conducted to analyze the complex regulatory patterns underlying β-thalassemia and HPFH.

Conclusions

This study shows that NR_001589 and uc002fcj.1 can act as a ceRNA to promote the expression of HBG1/2 by sponging miRNA during β-thalassemia and HbF induction. These results might help in designing a series of in vivo and in vitro experiments to explore the functions of NR_001589 and uc002fcj.1 through the ceRNA language. After establishing multiple lncRNA-miRNA-mRNA relationships, it can be presumed that these genetic factors are involved in β-thalassemia, thus laying the theoretical foundation for subsequent investigations.