LTA4H extensively associates with mRNAs and lncRNAs indicative of its novel regulatory targets

Cloning and plasmid construction

Hot fusion primer pairs were created using CE Design V1.04. Each primer contained a 17–30 bp sequence from the pIRES-hrGFP-1a vector and a gene-specific sequence fragment.

F-primer: agcccgggcggatccgaattcATGCCCGAGATAGTGGATACCTG

R-primer: gtcatccttgtagtcctcgagATCCACTTTTAAGTCTTTCCCCAC.

At 37 °C for 2 to 3 h, we digested the pIRES-hrGFP-1a vector with EcoRI and XhoI (NEB). The enzyme-digested vector was purified on a Qiagen column kit using 1.0% agarose gel. HeLa cells’ total RNA was obtained using Trizol. Oligo dT primers were used to transcribe the purified RNA for cDNA. Following that, PCR amplification was used to synthesize the inserted fragment. PCR insert and linearized vector digested with EcoRI and XhoI (NEB) were combined in a PCR microtube and ligated with Vazyme’s ClonExpress® II One Step Cloning Kit (Vazyme, Nanjing, Jiangsu, China). Chemical transformation was used to introduce plasmids into Escherichia coli strains. We incubated cells overnight at 37 °C on LB agar plates containing 1uL/ml ampicillin. Finally, 28 cycles of colony PCR were performed on the backbone vectors using universal primers to screen colonies.

Cell culture and transfections

The China Center for Type Culture Collection (CCTCC), Wuhan, Hubei, China provided human cervical carcinoma (CC) cell lines, HeLa (CCTCC@GDC0009). In Dulbecco’s modified Eagle’s medium (DMEM), which contains 10% fetal bovine serum (FBS), 100 ug/mL streptomycin and 100 U/mL penicillin, we cultivated HeLa cells at 37 °C and 5% CO2. Following the manufacturer’s instructions, HeLa cells were transfected using Lipofectamine 2000 (Invitrogen, Carlsbad, CA, USA). Transfected cells were collected after 48 h for RT-qPCR and western blot analyses.

Assessment of gene overexpression

To evaluated the effect of LTA4H overexpression, we used GAPDH (glyceraldehyde 3-phosphate dehydrogenase) as a control gene. The synthesis of cDNA was carried out according to standard procedures, and we performed RT-qPCR using Bestar SYBR Green RT-PCR Master Mix (DBI Bioscience, Shanghai, China) on a Bio-Rad S1000. Additional file contains primer information. After normalizing to GAPDH mRNA concentration level, each transcript was quantified using 2-ΔΔCT method (Livak & Schmittgen, 2001). The GraphPad Prism software (San Diego, CA, USA) was then used for comparison with the paired Student’s t-test.

Immunoprecipitation

HeLa cells were lysed on ice for 5 min with ice-cold lysis buffer (1 ×PBS, 0.5% sodium deoxycholate, 0.1% SDS, 0.5% NP40) containing RNase inhibitors (Takara, 2313) and protease inhibitors (329-98-6; Solarbio). In order to remove cell debris, the mixture was forcefully shaken and centrifuged for 20 min at 13,000 × g at 4 °C. The centrifuged supernatant was incubated overnight at 4 °C with DynaBeads protein A/G bound to normal IgG or anti-Flag LTA4H antibody. The beads were washed twice with low salt washing buffer, high salt washing buffer and 1X PNK buffer solution respectively, and the samples were suspended in the Elution Buffer to extract RNA from the LTA4H-RNA complex.

Western blot

Resuspend sample with 40 ul Elution Buffer 50 mM Tris-Cl (PH = 8.0), 10mM EDTA (PH = 8.0), 1%SDS; incubate it at 70 °C, 1,400 rpm for 20 min. The supernatant was put in a fresh EP tube. The complex was separated on a 10% SDS-PAGE gel after being boiled in boiling water with 1X SDS sample buffer for 10 min. with TBST buffer (20 mM Tris-buffered saline and 0.1% Tween-20) contained 5% non-fat milk power, we diluted the primary antibody: flag antibody (1:2,000, F7425; Sigma), actin (1:2,000, 66CUSABIO). The membranes were soaked in the primary antibody incubation solution and incubated at room temperature for 1 h. The membranes were then soaked in the HRP-conjugated secondary antibody incubation solution and incubated at room temperature for 1 h. The enhanced chemiluminescence (ECL) reagent (170506; Bio-Rad, Hercules, CA, USA) was used to detect the binding secondary antibody (anti-mouse or anti-rabbit 1:10,000) (Abcam).

iRIP-seq library preparation and sequencing

TRIzol (Invitrogen) was used to isolate the RBP-bound RNAs from the immunoprecipitation of anti-Flag. In accordance with the manufacturer’s instructions, complementary DNA (cDNA) libraries were prepared using KAPA RNA Hyper RNA binding protein connects the future Prep Kit (KK8541; KAPA). On the Illumina HiSeq X Ten platform, the cDNA libraries were sequenced for 150 bp paired-ends.

Data analysis

Only uniquely mapped reads were used for the subsequent analysis after reads were matched onto the genome using TopHat 2 (Kim et al., 2013) .“ABLIRC” strategy was utilized to determine the genomic locations where LTA4H binds (Xia et al., 2017). Peaks were formed from reads that had at least one base pair of overlap. Using computational simulation, reads with the same number and lengths as reads in peaks were generated randomly for each gene. For the purpose of generating random max peak height from overlapping reads, the outputting reads were further mapped to the same genes. The whole procedure was done 500 times. All observed peaks with heights greater than those of random maximum peaks (p-value 0.05) were chosen. The simulation independently analyzed the IP and input samples, removing the IP peaks that overlapped the input peaks. The peaks were used for motifs analysis with the Hypergeometric Optimization of Motif Enrichment (HOMER) software (Heinz et al., 2010).

Functional enrichment analysis

GO term and KEGG path enrichment analysis was performed using KOBAS 2.0 server (Xie et al., 2011). According to the annotation information of peak associated gene, the GO Term of each gene was counted, and significance of each Term was analyzed by Benjamini–Hochberg FDR (BH) and hypergeometric test to determine the degree of enrichment.

Reverse transcription qPCR validation

RT-qPCR was performed using total RNA from the iRIP-seq library preparation. Using the M-MLV Reverse Transcriptase (Vazyme), RNA was reverse transcribed into cDNA. Real-time PCR was carried out with the StepOne RealTime PCR System using the HieffTM qPCR SYBR® Green Master Mix (Low Rox Plus; Yeasen, Pudong, China). Denaturation at 95 °C for 5 min was followed by 40 cycles of denaturation at 95 °C for 15 s, annealing and extension at 60 °C for 30 s under PCR cycling conditions. PCR amplifications were carried out in triplicate for each sample.

Statistical analysis

The statistical software SPSS 16.0 (Chicago, IL, USA) was used to manipulate the experimental data, which were all presented as mean standard deviation (SD). All experiments were run at least three times independently, and P < 0.05 was considered significant.

Deregulated expression of LTA4H in various cancers

Previous researches have shown that LTA4H is significantly expressed in several malignancies and affects the initiation and growth of tumors (Chen et al., 2004; Guo et al., 2011; Jeong et al., 2009; Sun et al., 2006). In order to explore the relationship between LTA4H and laryngeal squamous cell carcinoma (LSCC), or more broadly head and neck squamous cell carcinoma (HNSCC), we first studied the expression level of LTA4H in LSCC tissues and normal tissues through The Cancer Genome Atlas (TCGA) database (Fig. 1A). Box plot and scatter diagram were used to display the expression levels of LTA4H (Transcripts Per Million (TPM)). According to the finding, LSCC tissues had lower levels of LTA4H than normal control tissues. The results obtained by Gao et al. (2019) with TCGA database are also LTA4H downregulation in LSCC. However, the author then detected LTA4H expression in LSCC and normal tissue by qPCR, and have proved that TA4H in LSCC tissues was significantly up-regulated (Gao et al., 2019). In addition, we also studied the association between LTA4H expression level and prognosis in HNSCC and normal tissues (Fig. 1B). A blue curve represented the low expression group, while a red curve represented the high expression group. A significant finding was that HNSCC patients with high expression of LTA4H had a poor prognosis. Therefore, the potential function of LTA4H in laryngeal squamous cell carcinoma deserves further study.

Characterization of the LTA4H-RNA interaction map by iRIP-seq analysis

To explore the potential function of LTA4H in LSCC, we obtained a LTA4H-bound RNA profile in modal HeLa cells by applying theiRIP-seq approach. The iRIP-seq is an advanced technique for studying RBPs, which achieves the precision of CLIP-seq to obtain both direct and indirect binding sites of protein and RNA accurately, whereas maintains the simplicity of RIP-seq. Labelled antibody and control antibody were used for immunoprecipitation, and two separately replicate experiments were performed. For immunoprecipitation, two separate iRIP repetitions were carried out using flag-tagged LTA4H. The western blots of both IP samples showed the presence of the protein Flag- LTA4H, but the IgG control did not (Fig. 2A). Then, we performed paired-end sequencing for the cDNA libraries using the Illumina HiSeq X Ten platform, and obtained the high-quality clean reads. After removing the adapter sequences and low-quality reads, we were left with 30,619,638 and 49,122,624 reads for IP group and input control of replicate 1, and 28,885,522 and 40,002,548 reads from those of the replicate 2 (Table S1). Next, using TopHat 2, we mapped the sequencing reads to reference genomes GRCh38 (Kim et al., 2013). About 77.23–78.94% of them were aligned and about 39.85–83.65% were matched uniquely. The uniquely mapped reads are overwhelmingly from mature mRNAs. The percentage of splice reads among the uniquely aligned reads was substantially greater in the IP sample compared with the control sample, indicating that LTA4H might have an involvement in the splicing role (Table S2).

Correlation analysis of IP and input samples by comparing Reads per kilo base of a gene per million reads (RPKM) (Mortazavi et al., 2008) of the same gene revealed transcripts were obviously enriched in IP samples than input control, which indicated that the specificity of the LTA4H-bound RNA was good (Fig. 2B). We also made correlation analysis between the two IP replicates, and the results showed that R = 0.917, which indicated that the two groups of IP samples had good repeatability (Fig. 2B). The results of these two groups of samples are almost the same, which indicated that the iRIP-seq experiment is reliable. The reads distribution across reference genomic regions showed LTA4H binding reads tend to concentrate in the CDS, the intron regions than input control reads, as well as in 3′ UTR, 5′ UTR and noncoding exons regions (Fig. 2C).

In order to eliminate the interference caused by gene expression quantity for predicting LTA4H specific binding sites, we adopted the ABLIRC method (Chi et al., 2009) to identify LTA4H-bound peaks precisely. There were 29,242 overlapping peaks in the two replication groups in Hela cells, indicating the overlap of peaks from the two sets of experiments is relatively high (Fig. 2D). Interesting, after sorting according to the number of reads on overlapping peaks, the top peaks were mainly located in mRNAs and lncRNAs (Fig. 2E). The results demonstrated that LTA4H has an extensive capability for RNA binding and may function as a regulator by interacting to mRNAs and lncRNAs. Peak distribution across reference genomic regions revealed that the LTA4H binding peaks located in the intron region accounted for a large proportion (66.93% and 63.20%), followed by CDS region (Fig. 2F).

Then, HOMER was employed to obtain sequence motif enriched within LTA4H peaks. The results showed UG-rich motif and GA-rich motif were presented as the first two motifs of LTA4H peaks of two replicates, which may be the key sites of LTA4H binding to its target (Fig. 2G). We found a high frequency of motif AAGG in both repeats in LTA4H binding peaks. It has recently been reported that TRA2B interacts with motif AAGG to promote cancer cell growth by disrupting gene expression processes associated with aging (Kajita et al., 2016). Our results suggest that LTA4H may interact with TRA2B for binding of the motif AAGG to regulate gene expression in cancerous cells.

In conclusion, the obtained LTA4H-binding RNA map will help our understanding of the overall regulatory mechanism of LTA4H-RNA association during gene expression in Hela cells.

Analysis of pre-mRNA and mRNAs associated by LTA4H

Further, the LTA4H overlapped peak associated genes were compared to the Gene Ontology database for enrichment biological process analysis. We found that LTA4H-bound genes were involved in gene expression, mitotic cell cycle, viral replication and DNA repair (Fig. 3A). Next, the DEseq package was used to identify the LTA4H-bound differentially enriched genes (DEGs) (Anders & Huber, 2010). Among the 14,170 DEGs, there were 2,776 enriched genes and 11,394 non-enriched genes related to LTA4H. We constructed a volcanic map to show the significantly enriched genes associated with LTA4H, all of which are associated with oncogenesis, including lncRNAs NEAT1 and LINC00657, and mRNAs ROR2, LTBP3, HSP90B1 and EGFR (Fig. 3B). To explore the potential biological role of these enriched genes, we continued to analyze the enriched genes using the GO database, and it showed they were mainly involved in negative regulation of transcription, mitotic cell cycle, gene expression, and viral replication (Fig. 3C).

Next, we performed an overlap analysis of LTA4H bound genes from ABLIRC algorithm and DEG from DEseq. Running DEseq identified fewer enriched genes, which were well overlapped by the LTA4H-bound genes by ABLIRC and resulted in 2425 overlapped genes (Fig. 3D). The results demonstrated a significant association between LTA4H-bound and enriched gene expression. GO analysis showed the 2425 overlapped genes were mainly clustered at gene expression, mitotic cell cycle, viral replication and DNA repair (Fig. 3E).

To further verify the presence of Flag-LTA4H protein binding on target genes, we next showed the distribution of reads binding location and coverage depth compared to Peak associated genes. The results across LTBP3 show the two replicates were consistent, and the IP groups were obviously biased towards the intron and exon regions than input control, which was the potential binding region of LTA4H on LTBP3 (Fig. 3F, left panel). And further, we used this gene to verify that it directly interacts with mRNA by RT-qPCR (Fig. 3F, right panel). In comparison to the control group, LTBP3 was considerably higher in the IP group. Similarly, mRNAs EGFR, ROR2 and HSP90B1 were distinctly enriched in IP samples compared to the input samples and the results of subsequent validation were as expected (Fig. S1). Taken together, our results suggest that LTA4H and mRNAs are closely interacted in Hela cells.

Analysis of LTA4H-bound lncRNAs

We also conducted the reads density landscape for lncRNAs that LTA4H highly enriched. The results across the cancer-promoting gene LINC00657 show there are many LTA4H-bound peaks in IP compared with input (Fig. 4A, left panel). In the validation experiment, LINC00657 was found to be significantly enriched in IP1, and there was no obvious bias in IP2 compared with the control group, which was considered to be related to experimental error (Fig. 4A, right panel). It may also be that the peak site is not enriched in IP2, that is, LTA4H may not bind LINC00657 specifically, they just bind randomly. The exact mechanism of LTA4H binding to LINC00657 has not been fully clarified and needs further study. NEAT1 had an obvious bias in IP groups, and we know that dysregulation of NEAT1 plays a key carcinogenic role (Chen et al., 2015; Wang et al., 2016) (Fig. 4B, left panel). Besides that, NEAT1 has been verified to be significantly enriched in the IP group (Fig. 4B, right panel). In summary, we hypothesized that LTA4H is preferentially bound to genes related to tumor formation and progression.