Skip to main content
Download PDF
  • Research article
  • Open access
  • Published:

Multi-species transcriptome meta-analysis of the response to retinoic acid in vertebrates and comparative analysis of the effects of retinol and retinoic acid on gene expression in LMH cells

BMC Genomics volume 22, Article number: 146 (2021) Cite this article

Abstract

Background

Retinol (RO) and its active metabolite retinoic acid (RA) are major regulators of gene expression in vertebrates and influence various processes like organ development, cell differentiation, and immune response. To characterize a general transcriptomic response to RA-exposure in vertebrates, independent of species- and tissue-specific effects, four publicly available RNA-Seq datasets from Homo sapiens, Mus musculus, and Xenopus laevis were analyzed. To increase species and cell-type diversity we generated RNA-seq data with chicken hepatocellular carcinoma (LMH) cells. Additionally, we compared the response of LMH cells to RA and RO at different time points.

Results

By conducting a transcriptome meta-analysis, we identified three retinoic acid response core clusters (RARCCs) consisting of 27 interacting proteins, seven of which have not been associated with retinoids yet. Comparison of the transcriptional response of LMH cells to RO and RA exposure at different time points led to the identification of non-coding RNAs (ncRNAs) that are only differentially expressed (DE) during the early response.

Conclusions

We propose that these RARCCs stand on top of a common regulatory RA hierarchy among vertebrates. Based on the protein sets included in these clusters we were able to identify an RA-response cluster, a control center type cluster, and a cluster that directs cell proliferation. Concerning the comparison of the cellular response to RA and RO we conclude that ncRNAs play an underestimated role in retinoid-mediated gene regulation.

Background

RO and its derivative RA belong to the vitamin A group of compounds. Derivatives of RO, termed retinoids, are involved in cell proliferation, differentiation, cell adhesion, and apoptosis in different types of vertebrate tissues [1] and play an important role in immunity (reviewed in [2]), male and female reproduction, embryonic development, and barrier integrity (reviewed in [3]). Hence, an in-depth understanding of gene regulation by retinoids is essential to understand their involvement in processes that affect health and diseases. RA is thought to be the main mediator of these effects and is therefore the most studied fat-soluble vitamin [3]. RA binds to different nuclear receptors that regulate gene expression through the binding to certain canonical sequences termed retinoic acid response-elements (RAREs). RAREs are typically two direct repeats of the sequence motif PuG (G/T) TCA with a variable spacer of 0–8 bases length (DR0-DR8) or are inverted repeats with no spacer (IR0) [4,5,6]. In 2002 Balmer and Blomhoff compiled a list of over 500 genes that have been identified to be regulatory targets of RA in different species and categorized them in a hierarchical manner. They identified 27 direct targets and 105 genes that can be modulated by RA [7]. Since these results might be biased by individual assumptions, we intended to generate an unbiased set of core RA response genes independent of tissue or cell type, exposure time, and species. A direct comparison of the transcriptomic responses of different cell and tissue types from different species has not been conducted so far. Hence, we performed a meta-analysis of RNA-seq data sets from five different vertebrate tissues and cells from four different species treated with RA for different periods of time. This led to the discovery of 91 DE genes. We were able to identify three RA response core interaction clusters, comprising 27 proteins of which seven to our knowledge have not been linked to RA. We propose that these networks of proteins are species- and tissue-spanning and mark the starting point of tissue-dependent downstream gene regulation after RA-stimulation.

Little focus has been put on elucidating whether RA and RO differ in their effect on gene expression. The only study conducted so far that compared gene expression in response to RA and RO investigated the application of both compounds to human skin. By histological assessment and real-time quantitative PCR (qPCR) for 12 target genes, Kong et al. concluded that the response of skin to RA and RO is similar with RA being more potent in its effect on gene and protein expression [8]. We conducted an in-depth comparison of the transcriptomic responses to RA and RO in chicken LMH cells. We thereby confirmed that RA exerts a stronger effect on down-stream targets and found only a 76% overlap in differentially expressed (DE) genes between both treatments. Furthermore, we observed differences in the early response to RA, which indicates an involvement of ncRNAs in the RA response.

Results

Transcriptome and differential expression analyses

To gain insights into species and tissue-specific effects of RA, the transcriptomic responses of five different cell and tissue types from four different species were compared: (i) LMH cells exposed to 100 nM RA for 4 h (N = 3, this study), (ii) human neuroblastoma cells (SH-SY5Y) exposed to 1 μM RA for 24 h (N = 2, BioProject PRJEB6636) [9], (iii) murine embryonic stem cells (mESCs) exposed to 1 μM RA for 48 h (N = 3, BioProject PRJNA274740) [10], (iv) murine lymphoblasts (mLympho) exposed to 1 μM of RA for 2 h (N = 4, BioProject PRJNA282594) [11], and in vitro-generated pancreatic explants from Xenopus laevis (Xenopus) exposed to 5 μM RA for 1 h (N = 2, each sample contained ~ 50 pooled explants, BioProject PRJNA448780) [12]. Additionally, we performed a comparative analysis of the response of LMH cells to RA and RO after 1 h and 4 h treatments. Alignment metrics after mapping of RNAseq reads with TopHat are shown in Table 1 and detailed results per sample and dataset are summarized in Additional file 1.

Table 1 Summary statistics of transcriptome mappings of all datasets used in the study
Full size table

A meta-analysis of the effects of retinoic acid on gene expression in different vertebrate tissues

The results of all DE analyses are summarized in Additional file 2. DE analysis of the datasets by comparing untreated with RA treated cells or tissues led to the discovery of 139 DE genes in LMH cells (73.4% upregulated), 164 DE genes in SH-SY5Y cells (68.9% upregulated), 3967 DE genes in mESCs (56.8% upregulated), 679 DE genes in murine lymphocytes (57.4% upregulated), and 48 DE genes in Xenopus (97.9% upregulated; p-adj < 0.01, abs. LFC > 1). Concordance of DE genes between the five analyses is represented by a Venn diagram (Additional file 3) and summarized in Additional file 4. None of the discovered DE genes were common in all five systems and the majority of DE genes were limited to each respective cell/tissue type. An overlap in at least two systems could be observed for 262 out of all DE genes. Due to the little overlap between the five datasets, we conducted a meta-analysis with MetaVolcanoR. This led to the discovery of 91 DE genes with a p-value < 0.02 and abs. LFC > 1 (Fig. 1; complete results are summarized in Additional file 2), all of which were upregulated. The 20 highest ranked DE genes are shown in Table 2. Four transcription factors could be detected among DE genes with the PANTHER classification system [13]: HEYL (LFC = 1.130, p-value = 1.31 × 10− 2), HIC1 (LFC = 3.264, p-value = 1.49 × 10− 3), RARB (LFC = 3.539, p-value = 4.17 × 4− 3), and TWIST2 (LFC = 3.037, p-value = 1.99 × 10− 2).

Fig. 1
figure 1

Volcano plot of differentially expressed genes from a transcriptome meta-analysis that was conducted with MetaVolcanoR. The results of each respective differential expression analysis from chicken hepatocellular carcinoma (LMH) cells, human neuroblastoma cells (SHSY5Y), murine embryonic stem cells (mESCs), murine lymphoblasts (mLympho), and in vitro-generated pancreatic explants from Xenopus laevis (Xenopus) after exposure to retinoic acid were used as input data. Red dots represent transcripts with a p-value < 0.02 and a LFC > 1

Full size image
Table 2 Top 20 DE genes from a multi-species transcriptome meta-analysis. RNA-seq data from five different cell types from four different vertebrate species after retinoic acid exposure were subjected to differential expression analysis and used as input for a meta-analysis
Full size table

To identify potential functional protein clusters among DE gene from the meta-analysis we performed protein interaction network analysis with STRING. The analysis revealed significantly more interactions than expected (Fig. 2, number of edges: 36, expected number of edges: 13, PPI enrichment p-value: 2.2 × 10− 7). Three distinct interaction clusters were identified: Cluster (i) contains the proteins ADRA2C, CCDC80, CCL19, CNR1, GDNF, IL18, NTRK2, OXT, P2RX1, RET, SEMA3A, and TACR3, cluster (ii) consists of the proteins CYP26A1, CYP26B1, CYP26C1, DHRS3, HIC1, HOXA2, HOXB1, HOXB2, and RARB and cluster (iii) contains CLDN11, CLDN2, ERMN, GALNT5, IFNW1, and TSPAN10. To identify general functions of RA, which are common among the five analyzed datasets we performed a gene cluster analysis with clusterProfiler using DE genes with p-values < 0.05 and abs. LFC > 0.5 as input data. Results are shown in Fig. 3 (complete analysis output is summarized in Additional file 5). GO biological processes affected by DE genes from the meta-analysis (Fig. 3a) are mainly involved in morphogenesis, development, and extracellular organization as well as “axon guidance” and “neuron projection guidance”. In regard to GO cellular components (Fig. 3b), most of the terms involve synaptic and postsynaptic membranes. The term with the lowest p-value and highest GeneRatio is “collagen-containing extracellular matrix”. GO molecular functions, which are enriched in the meta-analysis (Fig. 3c) involve transcription activator activity, receptor activity, extracellular matrix structure, and binding of sulfur, heparin, and retinoic acid. In regard to KEGG pathways (Fig. 3d) only “Neuroactive ligand-receptor interaction” reached statistical significance (p-value < 0.05).

Fig. 2
figure 2

Protein interaction analysis of differentially expressed genes from a transcriptome meta-analysis that was conducted with differential expression data from chicken hepatocellular carcinoma cells, human neuroblastoma cells, murine embryonic stem cells, murine lymphoblasts, and in vitro-generated pancreatic explants from Xenopus laevis after exposure to retinoic acid . DE genes with p-values < 0.02 and LFC > 1 were used for the analysis

Full size image
Fig. 3
figure 3

Gene cluster analysis of differentially expressed genes from a transcriptome meta-analysis that was conducted with differential expression data from chicken hepatocellular carcinoma cells, human neuroblastoma cells, murine embryonic stem cells, murine lymphoblasts, and in vitro-generated pancreatic explants from Xenopus laevis after exposure to retinoic acid. DE genes with a p-value < 0.05 and an abs. LFC > 0.5 were used for the analysis. a GO biological processes, b GO cellular components, c GO molecular functions, and d KEGG pathways

Full size image

Comparison of early and late RA and RO response in LMH cells

To compare the response of hepatic cells to RA and RO we analyzed differential expression in LMH cells treated with RA and RO for time periods of 1 h and 4 h. This led to the discovery of 21 DE genes after 1 h of RA treatment, 139 DE genes after 4 h of RA treatment, 8 DE genes after 1 h of RO treatment, and 128 DE genes after 4 h of RO treatment (p-adj < 0.01, abs. LFC > 1). The majority of DE genes were upregulated (95% RA 1 h, 76% RA 4 h, 100% RO 1 h, 75% RO 4 h). Volcano plots of DE genes after RA and RO exposure for both time points are shown in Additional file 6 and the complete results of the DE analysis are summarized in Additional file 2. The numbers of common and discordant DE genes from all four treatments are summarized in a Venn diagram (Fig. 4, complete results Additional file 7). Only seven genes were commonly DE in all four treatments: AADACL4L5, BARL, CYP8B1, LEKR1, LOC107054076 (ncRNA), RBPMS, TBX21, and TNFRSF8. The genes ATF3, BAIAP2, LOC101749099 (ncRNA), and LOC101750589 (ncRNA) are exclusive for the early response to RA. RO specific genes are ADAMTS9, AFAP1L2, ARHGAP24, CDKL2, LOC112530664 (ncRNA), LOC112531076 (pseudo-gene), LOC112531755 (ncRNA), LOC112531791 (ncRNA), PALMD, RUNX1T1, and VSIG10L. A total number of 26 genes were DE in a RA-dependent manner whereas a major overlap of 101 DE genes between RA and RO treatment after 4 h of exposure was observed. Genes with differences in expression between RA and RO treatment (min. 1.2-fold difference in Fragments per kilobase of exon model per million reads mapped (FPKM) values) after 1 h or 4 h are depicted in a heatmap (Fig. 5). The majority of differences in FPKM values were found between the time points, which were not considered in the heatmap. The genes with the highest differences in FPKM values between RA and RO treatment are listed in Table 3. The most distinct genes between both treatments are AADACL4L3, CYP26B1, HIC1, and RARB, all of which differ most in the early response and show stronger upregulation after RA stimulation. The only genes with a stronger response to RO (FPKM fold-difference > 1.2) are ARHGAP8, CDKL2, HS3ST1, and SLC5A12 after 1 h of exposure as well as AFAP1L2, LOC101749099 (ncRNA), LOC112530664 (ncRNA), LOC112531076 (pseudo-gene), LOC112531791 (ncRNA), and VSIG10L after 4 h of exposure. Among the genes with the highest differences in FPKM values between RA and RO treatment are seven ncRNAs.

Fig. 4
figure 4

Venn diagram of differentially expressed genes in LMH cells after exposure to retinoic acid for 1 h (RA_1h), retinoic acid for 4 h (RA_4h), retinol for 1 h (RO_1h), and retinol for 4 h (RO_4h)

Full size image
Fig. 5
figure 5

Heatmap of DE genes that differ between retinoic acid and retinol treatment in LMH cells: Log(FPKM) values of genes with at least 1.2-fold difference in FPKM values between retinoic acid and retinol treatment after 1 h or 4 h hours are shown. Cells treated with retinoic acid for 1 h (RA_1h), were compared with cell treated with retinol for 1 h (RO_1h) and cells treated with retinoic acid for 4 h (RA_4h), were compared with cell treated with retinol for 4 h (RO_4h)

Full size image
Table 3 Genes that are differentially expressed in at least one treatment with a > 1.5 fold difference in FPKM values between retinoic acid and retinol treatment in LMH cells. Cells treated with retinoic acid for 1 h (RA_1h), were compared with cell treated with retinol for 1 h (RO_1h), and cells treated with retinoic acid for 4 h (RA_4h), were compared with cell treated with retinol for 4 h (RO_4h)
Full size table

To elucidate if the DE genes that we identified by exposing LMH cells to RA and RO might be RARE-regulated the chicken reference genome (GCF_000002315.5) was screened for RAREs (DR0-DR8 and IR0). The numbers of RAREs in the vicinity of DE genes (up to 10 kb upstream of transcript start and 10 kb downstream of transcript end) are summarized in Additional file 8. We detected RAREs in the vicinity of 103 out of 150 DE genes from all four treatments with an average of 2.07 RAREs per gene. The average occurrence of RAREs per gene in the genome is 0.77. Genes with ten or more RAREs close to the gene coding region are ARHGAP24, OBSCN, RARB, STARD13, and TOX.

To find out whether certain protein interaction networks are differentially affected by RA and RO treatment the products of DE genes after 4 h of exposure to RA and RO were subjected to protein interaction network analyses with STRING [14] (interaction graphs in Additional file 9). In both cases, the networks had significantly more interactions than expected (RA treatment: number of edges: 41, expected number of edges: 21. PPI enrichment p-value: 7.29 × 10− 5; RO treatment: number of edges: 28, expected number of edges: 17. PPI enrichment p-value: 0.0107). With a higher level of significance and a higher number of edges, we could observe a higher degree of protein interaction among RA-regulated genes. Among those genes is a cluster of HOX genes (HOXA1, HOXA3, HOXA5, HOXB3, and HOXB4) and a cluster of genes primarily involved in bone development (MSX2, RUNX2, THBS1, TNFRSF11B, TOR4A). The interaction cluster surrounding RARB is larger (15 proteins) in RA-treated cells compared to RO-treated cells (8 genes). One interaction cluster that both treatments have in common consists of four genes encoding proteins with G protein-coupled receptor activity: BDKRB2, GPR37L1, GRM8, and HTR2A. To investigate if short- and long-term RA and RO exposure have different effects on the cellular response we performed a cluster analysis of DE genes (p-adj < 0.01, abs. LCF > 0.5) with clusterProfiler (complete analysis output is summarized in Additional file 5). The analysis revealed that treatment with RA and RO leads to an increase in GO biological processes associated with embryo, organ and skeletal system development and morphogenesis. RA acts more potent on the GO terms “embryo organ morphogenesis”, “embryonic organ development”, “animal organ development”, and “embryo development ending in birth or egg hatching” (Fig. 6a). The impact of RA on GO molecular functions was significantly higher as compared to RO with the majority of GO terms related to transcription, DNA-binding, gene expression, and metal ion binding. Comparable p-values between cells treated with RA and RO were only found for the GO terms “DNA-binding transcription factor activity” and “transcription regulator activity” (Fig. 6b). Due to the limited amount of DE genes detected for the 1 h time point comparison of early and late response to RA and RO was only possible in the KEGG pathway analysis. KEGG pathways limited to the early response to RA and RO stimulation are “Cytokine-cytokine receptor interaction”, “Phosphatidylinositol signal system”, and “Primary bile acid biosynthesis”. “Apoptosis”, and “Glycosaminoglycan biosynthesis – heparin sulfate / heparin” were only affected after 1 h of RA stimulation and “Insulin signaling pathway” and “mTOR signaling pathway” after 1 h of RO exposure. An exposure of 4 h to RA and RO led to lower p-values in “Retinol metabolism” and “Adipocytokine signaling pathway”. Prominent effects of RA and RO limited to an exposure of 4 h include KEGG pathways related to lipid metabolism, “FoxO signaling pathway”, and “Wnt signaling pathway” (Fig. 6c).

Fig. 6
figure 6

Gene cluster comparison of differentially expressed genes in LMH cells with clusterProfiler after exposure to retinoic acid for 1 h (RA_1h), retinoic acid for 4 h (RA_4h), retinol for 1 h (RO_1h), and retinol for 4 h (RO_4h). Sufficient differentially expressed genes were found to analyze a GO biological processes, b GO molecular functions, and c KEGG pathways

Full size image

Discussion

A meta-analysis of the transcriptomic responses to retinoic acid from different species

To gain further insights into RA-dependent gene-regulation we acquired four RNA-seq datasets from the NCBI SRA and mapped them to the most recent genome assembly of each respective species (Homo sapiens, Mus musculus, and Xenopus laevis). To increase species and cell type variety we performed RNA-seq on chicken hepatocellular carcinoma (LMH) cells after RA exposure. We ended up with whole transcriptome DE data from five different systems: chicken LMH cells, human neuroblastoma cell line SH-SY5Y, murine embryonic stem cells, murine lymphoblasts, and in vitro-generated pancreatic explants from Xenopus laevis. Data quality regarding read length and coverage was mixed. Exon coverages around 50x were achieved with LMH cells and murine lymphoblasts. Coverages around 10x for the mESC and Xenopus mappings are acceptable whereas a 4x coverage and a multiple alignment frequency of 32.8% in SH-SY5Y cells might have introduced bias into the DE analysis of this dataset. The high frequency of multiple alignments is a result of the short read length and the absence of paired reads. Hence, accuracy of the results may be affected by the relatively low to medium quality of the SH-SY5Y, mESC and Xenopus data sets. The number of DE genes in response to RA-stimulation appears to stand in direct relation to the transcriptional activity of the respective cell- and tissue-types. mESCs are by far most susceptible to RA-stimulation with almost 4000 DE genes, followed by murine lymphoblasts with 679 DE genes. However, the overlap of DE genes between the five systems was not very prominent (Additional file 3). Hence, we conducted a transcriptome meta-analysis with MetaVolcanoR. By using the random effect model we circumvent the introduction of bias by differing p-value dimensions between the five datasets. It produces summary LFCs based on the variance, which are then used to estimate summary p-values. This is followed by perturbation ranking with the topconfects approach, which ranks results by confident effect sizes. By combining RNA-seq data of five different cell types from four different species taken at different time points the experimental conditions are highly adverse. This enabled us to compile a high confidence RA-response gene set in this meta-analysis. We discovered 91 DE genes (Fig. 1) of which 27 are part of three protein interaction clusters (Fig. 2), which we term retinoic acid response core clusters (RARCCs). Of those 27 RARCC genes, seven have not been previously associated with RA: P2RX1, TACR3, HIC1, ERMN, GALNT5, IFNW1, and TSPAN10 (Table 4).

Table 4 Summary of published data on genes that belong to the retinoic acid response protein interaction clusters
Full size table

Cluster (ii) is obviously on top of the RA response hierarchy since it contains proteins that directly metabolize RA, like CYP26A1, CYP26B1, CYP26C1, and DHRS3, or are well known RA receptors (RARB) and direct downstream targets (Hox genes). Cluster (i) appears to be a control center that exerts cell-type-dependent downstream responses. We come to that conclusion because it is an adverse mix of well-described RA responses, like neuronal development, pathogen response, and transcription initiation. Cluster (iii) is of special interest since it contains four genes that to our knowledge have not been associated with RA before. Based on their reported functions (UniProt) we assume that this RARCC is responsible for the stimulating effect that RA has on cell proliferation. Since the RARCCs are species and tissue spanning as well as time-point independent, we propose that these interaction networks are conserved among vertebrates and tissue/cell types. This knowledge might have implications for future research since these 27 genes can be considered on top of an RA-response hierarchy that gets more specialized downstream depending on cell and tissue type as well as time. The fact that 20 out of these 27 RARCC genes have been previously associated with the response to RA underscores the reliability of the meta-analysis approach. Furthermore, it makes the 7 newly discovered genes likely candidates to be RA-regulated, which needs to be further validated experimentally. Gene cluster analysis confirms previous findings concerning the cellular response to RA exposure (Fig. 3). These include neuron-specific GO terms and KEGG pathways as well as terms involved in cell proliferation, transcription, receptor activity, and binding of certain chemical compounds. RA is able to induce neuron-like phenotypes in stem cells and plays a major role in the switch between cell proliferation and neuronal differentiation (Reviewed in [29, 30]). That this function of RA is present on a higher level in the RA-response hierarchy is confirmed by our data. As Tang and Gudas outlined in their review article [31] RA is able to induce or inhibit cell proliferation in many cell types depending on the studied system. With our gene cluster analysis we can confirm that this is mostly a cell type-independent function of RA. This also holds true for the effect of RA on the activation of gene transcription in general. Further information on the effects of RA on gene transcription can be found in the overview article by Amann et al. [32]. In comparison with our findings for each respective dataset, these findings further indicate that we compiled a set of higher-level genes that orchestrate the RA-response further downstream depending on cell and tissue type. Only four DE transcription factors were detected in the meta-analysis of which three were previously described to be RA-responsive: HIC1 [33], RARB [26], and TWIST2 [34]. The transcriptional repressor HEYL, which is involved in cardiac gene expression [35], has not been linked to RA-mediated gene regulation and should be considered in further studies focusing on retinoid signaling in cardiac development (reviewed in [36]). The fact that virtually no genes were significantly downregulated in the meta-analysis of these five datasets confirms the notion that RA is mainly a transcriptional activator and not a repressor. Transcriptional repression seems to be regulated in a cell type-dependent manner. We saw downregulation of transcripts among the five individual datasets ranging from 43.2 to 2.1% of all DE genes. One gene that caught our attention is TOX, which is downregulated in three of the datasets: LMH cells, murine lymphoblasts, and mESCs. TOX is a transcriptional regulator that has not been linked to RA-biology. Our data indicate that TOX is a RARE regulated gene. We discovered 10 RAREs in the vicinity (up to 10 kb up- and down-stream) to the TOX gene in the chicken genome (Supporting Information 8). TOX is a diagnostic marker for cutaneous T-cell lymphomas (CTCL) and is overexpressed in affected CD4+ cells [37] Retinoids have been successfully used in the therapy of CTLC for over 30 years [38] and are still being used in T-cell lymphoma therapy [39]. The results from our analysis suggest that the success of retinoid cancer therapy is directly connected to TOX downregulation by RA. Identification of genetic polymorphisms related to TOX in patients that are susceptible to retinoid therapy and determination of the TOX expression state in CD4+ cells could aid in the tailored development of CTLC combination-therapy.

Comparison of the retinoic acid and retinol response in chicken LMH cells

Since most past studies focused on the effects of RA on gene expression we wanted to investigate if the response to its dietary form RO is different. In summary, we were able to confirm the conclusion by Kong et al. that RA is a more potent modulator of gene expression compared to RO [8]. After 1 h incubation time, twice as many genes were DE in the RA treatment (Additional file 6) although the concentration of RO was ten times higher. We explain that effect by a delay in the response to RO since it first has to be oxidized to RA via retinal, which takes place in the cytosol of hepatic cells [40] and in microsomes [41]. Furthermore, the RO treated cells might produce physiological amounts of RA, which are very likely lower than the applied RA concentration in the other experimental group. After 4 h a comparable number of genes were DE in RA and RO stimulated cells. With an overlap of 76% in DE genes, the transcriptional responses to the two chemicals are similar, but not identical. By real-time PCR Kong et al. identified six genes (COL1A1, COL3A1, CRABP2, FLG, TGM1, and TGM3) that respond to RA and RO stimulation in human skin, none of which were DE in LMH cells, which were derived from chicken liver carcinoma. Only four genes were limited to the early response to RA, two of which were ncRNAs. Among RO-specific genes three out of eleven were also ncRNAs. Furthermore, the majority of DE transcripts that showed the strongest difference between both treatments are ncRNAs (Table 3). We hypothesize that ncRNAs play an important role in early processes of RA mediated gene regulation and in RO metabolism, which has already been described in the context of RA-mediated TGM2 gene regulation [42]. BLAST search revealed that the ncRNA LOC112530664 has binding properties to TMEM168, a gene that promotes cell proliferation [43]. Furthermore, LOC112531791 shows sequence homology to ABCC8, which is a regulator of ATP-sensitive K+ channels and insulin release [44]. LOC112530664 shows stronger upregulation after RO treatment in comparison to RA and LOC112531791 is DE in an RO-dependent manner. Hence, these ncRNAs might be involved in RO metabolism by acting on those genes. To clarify the involvement of TMEM168 and ABBC8 in the cellular response to retinol requires further validation on the protein level.

Hox genes are crucial to vertebrate embryogenesis and are known to get activated by signaling cascades initiated by RA (reviewed in [45]). Our results indicate that RA has a stronger influence on Hox gene expression in comparison to RO. An interaction cluster consisting of five HOX proteins (HOXA1, HOXA3, HOXA5, HOXB3, and HOXB4) was detected with STRING (Additional file 9), whereas RO only led to the activation of three Hox genes. We assume that the conversion of RO to RA leads to lower RA levels in the cell compared to the direct administration of RA. Hence, the Hox gene response is less prominent after RO treatment. Another protein interaction cluster that we found in the RA response consists of MSX2, RUNX2, THBS1, TNFRSF11B, and TOR4A, all of which are involved in development. Dickson et al. demonstrated that embryonic chicken calvaria responds differently to RA and RO administration [46]. We can indirectly confirm these results since we did not find the mentioned bone development protein interaction cluster after RO administration. Other studies, which only focused on the RA response, found a stimulating effect on osteoclasts [47] and an inhibitory effect on osteoblasts [48], which suggests that RA leads to bone degradation rather than bone formation (discussed in this review [49]). The protein interaction cluster surrounding RARB is almost twice as large after RA treatment in comparison to RO treatment. Since the three direct interaction partners of RARB, namely NRIP1, ALDH1A3, and CYP26B1 are DE after both treatments, we assume that higher RA bioavailability in the cells leads to a higher downstream effect on gene expression of RARB targets. For instance, the RAR coregulatory NRIP1 [50], which is a transcription factor that regulates lipid and glucose metabolism in the liver [51], exhibits a stronger downstream effect after RA treatment. Furthermore, ALDH1A3 and CYP26B1, both of which are involved in retinoid metabolism [52, 53], may have more substrate to process in the presence of RA compared to RO. The interaction cluster containing the proteins BDKRB2, GPR37L1, GRM8, and HTR2A is present in the interaction maps of both treatments. Hence, G protein-coupled receptor activity seems to be equally responsive to RA and RO but has so far only been described for RA as a possible activator of Non-canonical Wnt Signaling [54].

Gene cluster analysis also revealed that RA is the more potent activator of gene expression in comparison to RO. Concerning GO biological processes (Fig. 6a) RA has a higher impact on terms that involve embryo and organ development a known function of RA as reviewed here [55]. The same holds true for terms belonging to GO molecular functions (Fig. 6b). Both compounds lead to an upregulation of terms affecting transcription, DNA-binding, gene expression, and metal ion binding with RA initiating a stronger response. KEGG pathway analysis (Fig. 6c) revealed some interesting insights into the early response of LMH cells to RA and RO. Pathways that are limited to 1 h of RA and RO treatment are “Cytokine-cytokine receptor interaction”, “Phosphatidylinositol signal system”, and “Primary bile acid biosynthesis”. An influence of RA on cytokines has been reported [56, 57] but not in a time-dependent manner. However, an immediate decrease in phosphatidylinositol turnover after RA exposure has been reported in neuroblastoma cells [58, 59]. Concerning the early response, we found differences in the response to RA and RO. Whereas early RA stimulation had an effect on “Apoptosis”, and “Glycosaminoglycan biosynthesis – heparin sulfate / heparin”, RO had an exclusive effect on “Insulin signaling pathway” and “mTOR signaling pathway”. With respect to insulin signaling, insulin was shown to regulate RA biosynthesis by upregulation of retinol dehydrogenase expression [60]. Our data suggest that vice versa, RO can upregulate genes that are involved in insulin metabolism in an immediate manner. Regarding mTOR signaling, synaptic RA receptors mediate hippocampal learning via mTOR dependent metaplasticity [61]. Since learning is an immediate and hippocampal consolidation a fast process [62], early activation of the mTOR signaling pathway after RO administration is conclusive.

Surprisingly, the response of LMH cells to RA and RO exposure with respect to KEGG “Retinol metabolism” genes was identical. After 1 h incubation time, only DHRS3 was significantly upregulated in both treatments with a higher LFC after RA (2.455) exposure compared to RO (1.874). DHRS3 reduces all-trans-retinal to all-trans-retinol or oxidizes all-trans-retinol to all-trans-retinal [25], most likely dependent on the stoichiometry between the two chemicals. Hence, lower expression in the presence of RO is conclusive, since it first has to be metabolized to RA via retinal as an intermediate product. After 4 h of RA and RO exposure, three additional genes were differentially expressed: CYP26B1, RDH10, and UGT1A1. CYP26B1 hydroxylates RA to 4-OH-RA, 4-oxo-RA, or 18-OH-RA [52]. RDH10, which catalyzes the conversion of all-trans-retinol to all-trans-retinal [63], was downregulated after both treatments. UGT1A1 activity leads to glucuronidation of RA [64], a detoxification process that takes place in the liver [65, 66]. Taken together upregulation of CYP26B1, and UGT1A1 and downregulation of RDH10 seem to be a detoxification mechanism to get rid of excess RA. If these genes play a role in retinoic acid syndrome [67] (reviewed in [68]) remains to be elucidated. The toxicity of RA in the treatment of acute promyelocytic leukemia has been described for the first time in a clinical case in 1992 [69]. We used a RA concentration of 100 nM and already observed a potentially toxic response. Hence, we conclude that a lower RA concentration or the application of RO in functional experiments might produce results that are closer to the natural response to retinoids. Furthermore, the studies discussed above utilized RA concentrations of 1 μM or higher, which may have introduced bias to the results by an overdose effect.

Conclusions

By conduction a meta-analysis of the transcriptomic responses to RA exposure of five different vertebrate systems we were able to identify a core RA response gene set. From our results, we conclude that on a higher hierarchical level RA is an activator of transcription and that RA mediates transcription repression in a cell type-dependent manner. Furthermore, we conclude that RA exerts its downstream functions via three distinct protein interaction clusters: The largest cluster comprises diverse downstream targets of RA and might function as a control hub, which acts cell type-dependent. One mid-sized cluster almost exclusively contains well described direct RA targets, which we consider on top of a general gene expression hierarchy in vertebrates. The smallest cluster seems to be cell-type independent since it mainly contains genes, which are involved in cell proliferation – a general function of RA. The comparative analysis of the influence of RA and RO on gene expression in LMH cells confirmed that RA is a more potent inducer of gene expression. However, a discordance of 24% in DE genes caught our attention. Among those are two RA- and three RO-specific ncRNAs from which we conclude that ncRNAs play a central role in the early response to retinoids.

Methods

Cell culture

LMH cells were kindly provided by Prof. Schwemmle’s lab and cultured in MEM supplemented with 10% FBS. Culture conditions were 37 °C, 5% CO2 and 95% RH. 2 × 106 Cells were seeded on a T25 flask and grown for 24 h. Cells were either treated with DMSO containing RA (final concentration 100 nM), Retinol (final concentration 1 μM) or DMSO only (control group) for 1 h or 4 h. The experiment was repeated three times independently.

RNA isolation

Cells were washed with PBS, detached from the surface by adding 1 ml of accutase for 5 min, resuspended in twice the amount of medium, transferred to a 15 ml tube, centrifuged for 5 min at 400 g, resuspended in PBS, centrifuged for 5 min at 400 g, and the resulting cell pellet was used for RNA isolation. RNA was isolated with the RNeasy Mini Kit and QIAshredder columns (QIAGEN) according to the manual.

NGS library preparation and sequencing

Procedures have already been described in a previous study [70]. Briefly: Total RNA input quality was evaluated on a TapeStation 4200 (Agilent, USA), and all samples showed a RIN score > 8. Samples were quantified with a fluorometric dye (Quant-IT, thermofisher, USA) and 500 ng per sample were used as input for the TruSeq stranded mRNA library kit (Illumina, USA) following the manufacturers manual. Resulting libraries showed a fragment size distribution of around 300 bp and were sequenced on a NovaSeq S2 Flowcell (Illumina, USA) with 50 bp paired-end reads.

Resource datasets

The following sample comparisons were carried out in the differential expression analysis: LMH cells exposed to RA or RO were compared to DMSO-treated cells. All three treatments were applied for 1 h and 4 h. SH-SY5Y exposed to 1 μM RA for 24 h were compared to DMSO treated cells (BioProject PRJEB6636; RA-treated samples: ERR550444, ERR550446; DMSO Samples: ERR550449, ERR550450) [9]. mESCs exposed to 1 μM RA for 48 h were compared to untreated control cells (BioProject PRJNA274740; RA-treated samples: SRR1792530, SRR1792529, SRR1792531; control samples: SRR1792526, SRR1792528, SRR1792527) [10]. Murine lymphoblasts exposed to 1 μM of RA for 2 h were compared to DMSO treated cells (BioProject PRJNA282594; RA-treated samples: SRR2001796, SRR2001794, SRR2001797, SRR2001795; control samples: SRR2001790, SRR2001793, SRR2001791, SRR2001792) [11]. In vitro-generated pancreatic explants from Xenopus laevis (Xenopus) exposed to 5 μM RA for 1 h were compared to DMSO treated cells where each sample contained ~ 50 pooled explants (BioProject PRJNA448780; RA-treated samples: SRR6941647, SRR6941648; control samples: SRR6941648, SRR6941644) [12].

Transcriptome analyses

DE analyses were carried out as previously described [70]. Briefly, quality control and trimming of raw sequencing reads was achieved with Trimmomatic version 0.36 (settings: PE -phred33 LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:36) [71]. Reads were aligned to the most recent reference genomes with TopHat version 2.1.0 (settings: --no-novel-juncs --min-isoform-fraction 0.0 --min-anchor-length 3 -r 192) [72]. Reference genomes (RefSeq assemblies) used for the alignment are Gallus gallus GCF_000002315.5, Mus musculus GCF_000001635.26, Homo sapiens GCF_000001405.39, and Xenopus laevis GCF_001663975.1. The R packages GenomicFeatures (Version 1.40.0) and summarizeOverlaps were used to count exon spanning reads [73]. DE analyses were conducted with DESeq2 (Version1.28.1) [74]. The R package EnhancedVolcano (Version 1.6.0) was used to generate volcano plots of DE results. Meta-analysis of DE datasets was carried out using the R package MetaVolcanoR (Version 1.2.0) using a random effect model. The following settings were applied: geneidcol = NULL, collaps = FALSE, cvar = FALSE, metathr = 0.01, ncores = 8. Gene Symbols of the input datasets were generalized by capitalizing all letters and removal of species specific pre- and suffixes.

Functional analyses

Gene cluster comparison and visualization was achived with the R package clusterProfiler [75]. Gene symbols were converted to ensemble IDs with the clusterProfiler Biological Id Translator (bitr). GO term analyses were performed with enrichGO (settings: pAdjustMethod = “fdr”, pvalueCutoff = 1, qvalueCutoff = 0.25, readable = TRUE, minGSSize = 10). KEGG pathtway analysis was done with enrichKEGG (settings: pvalueCutoff = 1, pAdjustMethod = “BH”,minGSSize = 10, maxGSSize = 500, qvalueCutoff = 0.25, use_internal_data = FALSE). Plots were created with the dotplot function. Protein interaction maps were done with STRING (Version 11.0) [14] using default settings.

Availability of data and materials

The datasets with the following BioProject IDs were aquired from the NCBI SRA: PRJEB6636 (SH-SY5Y cells) https://trace.ncbi.nlm.nih.gov/Traces/sra/sra.cgi?study=ERP006185, PRJNA274740 (mESCs) https://trace.ncbi.nlm.nih.gov/Traces/sra/sra.cgi?study=SRP053290, PRJNA282594 (murine lymphoblasts) https://trace.ncbi.nlm.nih.gov/Traces/sra/sra.cgi?study=SRP057791, PRJNA448780 (Xenopus laevis) https://trace.ncbi.nlm.nih.gov/Traces/sra/sra.cgi?study=SRP137258 . The raw sequencing data that was generated for this study is accessible under BioProject ID PRJNA667585 (reviewer link: https://dataview.ncbi.nlm.nih.gov/object/PRJNA667585?reviewer=koal33cfclsaj7d2c7nifjipl2) and will be made publicly available upon publication. For the mapping of RNA-seq reads from LMH cells chicken genome version GCF_000002315.5 was used (genome assembly file: https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/002/315/GCF_000002315.5_GRCg6a/GCF_000002315.5_GRCg6a_genomic.fna.gz, genomic features file: https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/002/315/GCF_000002315.5_GRCg6a/GCF_000002315.5_GRCg6a_genomic.gff.gz). For the mapping of RNA-seq reads from murine cells Mus musculus genome version GCF_000002315.5 was used (genome assembly file: https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/001/635/GCF_000001635.26_GRCm38.p6/GCF_000001635.26_GRCm38.p6_genomic.fna.gz, genomic features file: https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/001/635/GCF_000001635.26_GRCm38.p6/GCF_000001635.26_GRCm38.p6_genomic.gff.gz). For the mapping of RNA-seq reads from SH-SY5Y cells Homo sapiens genome version GCF_000001405.39 was used (genome assembly file: https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/001/405/GCF_000001405.39_GRCh38.p13/GCF_000001405.39_GRCh38.p13_genomic.fna.gz, genomic features file: https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/001/405/GCF_000001405.39_GRCh38.p13/GCF_000001405.39_GRCh38.p13_genomic.gff.gz). For the mapping of RNA-seq reads from Xenopus tissue Xenopus laevis genome version GCF_001663975.1 was used (genome assembly file: https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/001/663/975/GCF_001663975.1_Xenopus_laevis_v2/GCF_001663975.1_Xenopus_laevis_v2_genomic.fna.gz, genomic features file: https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/001/663/975/GCF_001663975.1_Xenopus_laevis_v2/GCF_001663975.1_Xenopus_laevis_v2_genomic.gff.gz). For the mapping of RNA-seq reads from SH-SY5Y cells Homo sapiens genome version GCF_000001405.39 was used (genome assembly file: https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/001/405/GCF_000001405.39_GRCh38.p13/GCF_000001405.39_GRCh38.p13_genomic.fna.gz, genomic features file: https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/001/405/GCF_000001405.39_GRCh38.p13/GCF_000001405.39_GRCh38.p13_genomic.gff.gz).

Abbreviations

CTCL:

Cutaneous T-cell lymphoma

DE:

Differentially expressed

DR:

Direct repeat

FDR:

False discovery rate

FPKM:

Fragments per kilobase of exon model per million reads mapped

LFC:

Log fold change

LMH cells:

Chicken hepatocellular carcinoma cells

ncRNAs:

Non-coding RNAs

p-adj:

Adjusted p-value

qPCR:

Real-time quantitative PCR

RA:

Retinoic acid

RARCC:

Retinoic acid response core cluster

RAREs:

Retinoic acid response-elements

RNA-seq:

RNA-sequencing

RO:

Retinol

SRA:

Sequence Read Archive

References

  1. Al Tanoury Z, Piskunov A, Andriamoratsiresy D, Gaouar S, Lutzing R, Ye T, et al. Genes involved in cell adhesion and signaling: a new repertoire of retinoic acid receptor target genes in mouse embryonic fibroblasts. J Cell Sci. 2014;127:521–33. https://doi.org/10.1242/jcs.131946.

    Article  CAS  PubMed  Google Scholar 

  2. Pino-Lagos K, Guo Y, Noelle RJ. Retinoic acid: a key player in immunity. Biofactors. 2010;36:430–6. https://doi.org/10.1002/biof.117.

    Article  CAS  PubMed  Google Scholar 

  3. Li Y, Wongsiriroj N, Blaner WS. The multifaceted nature of retinoid transport and metabolism. Hepatobiliary Surg Nutr. 2014;3:126–39. https://doi.org/10.3978/j.issn.2304-3881.2014.05.04.

    Article  PubMed  PubMed Central  Google Scholar 

  4. Leid M, Kastner P, Chambon P. Multiplicity generates diversity in the retinoic acid signalling pathways. Trends Biochem Sci. 1992;17:427–33. https://doi.org/10.1016/0968-0004(92)90014-Z.

    Article  CAS  PubMed  Google Scholar 

  5. Mangelsdorf DJ, Evans RM. The RXR heterodimers and orphan receptors. Cell. 1995;83:841–50. https://doi.org/10.1016/0092-8674(95)90200-7.

    Article  CAS  PubMed  Google Scholar 

  6. Moutier E, Ye T, Choukrallah M-A, Urban S, Osz J, Chatagnon A, et al. Retinoic acid receptors recognize the mouse genome through binding elements with diverse spacing and topology. J Biol Chem. 2012;287:26328–41. https://doi.org/10.1074/jbc.M112.361790.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  7. Balmer JE, Blomhoff R. Gene expression regulation by retinoic acid. J Lipid Res. 2002;43:1773–808. https://doi.org/10.1194/jlr.r100015-jlr200.

    Article  CAS  PubMed  Google Scholar 

  8. Kong R, Cui Y, Fisher GJ, Wang X, Chen Y, Schneider LM, Majmudar G. A comparative study of the effects of retinol and retinoic acid on histological, molecular, and clinical properties of human skin. J Cosmet Dermatol. 2016;15:49–57. https://doi.org/10.1111/jocd.12193.

    Article  PubMed  Google Scholar 

  9. Duffy DJ, Krstic A, Halasz M, Schwarzl T, Konietzny A, Iljin K, et al. Retinoic acid and TGF-β signalling cooperate to overcome MYCN-induced retinoid resistance. Genome Med. 2017;9:15. https://doi.org/10.1186/s13073-017-0407-3.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  10. Terranova C, Narla ST, Lee Y-W, Bard J, Parikh A, Stachowiak EK, et al. Global developmental gene programing involves a nuclear form of fibroblast growth factor Receptor-1 (FGFR1). PLoS One. 2015;10:e0123380. https://doi.org/10.1371/journal.pone.0123380.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  11. Churchman ML, Low J, Qu C, Paietta EM, Kasper LH, Chang Y, et al. Efficacy of Retinoids in IKZF1-mutated BCR-ABL1 acute lymphoblastic leukemia. Cancer Cell. 2015;28:343–56. https://doi.org/10.1016/j.ccell.2015.07.016.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  12. Gere-Becker MB, Pommerenke C, Lingner T, Pieler T. Retinoic acid-induced expression of Hnf1b and Fzd4 is required for pancreas development in Xenopus laevis. Development. 2018. https://doi.org/10.1242/dev.161372.

  13. Mi H, Muruganujan A, Ebert D, Huang X, Thomas PD. PANTHER version 14: more genomes, a new PANTHER GO-slim and improvements in enrichment analysis tools. Nucleic Acids Res. 2019;47:D419–26. https://doi.org/10.1093/nar/gky1038.

    Article  CAS  PubMed  Google Scholar 

  14. Szklarczyk D, Gable AL, Lyon D, Junge A, Wyder S, Huerta-Cepas J, et al. STRING v11: protein-protein association networks with increased coverage, supporting functional discovery in genome-wide experimental datasets. Nucleic Acids Res. 2019. https://doi.org/10.1093/nar/gky1131.

  15. Chen W-K, Chang N-CA, Chang Y-H, Chang K-L, Wu S-C, Yang T-S, et al. Characterization of the regulatory region of Adra2c, the gene encoding the murine α2C Adrenoceptor subtype. J Biomed Sci. 2004;11:886–901. https://doi.org/10.1159/000081836.

    Article  CAS  PubMed  Google Scholar 

  16. Chee LCY, Hendy J, Purton LE, McArthur GA. ATRA and the specific RARα agonist, NRX195183, have opposing effects on the clonogenicity of pre-leukemic murine AML1-ETO bone marrow cells. Leukemia. 2013;27:1369–80. https://doi.org/10.1038/leu.2012.362.

    Article  CAS  PubMed  Google Scholar 

  17. Mukhopadhyay B, Liu J, Osei-Hyiaman D, Godlewski G, Mukhopadhyay P, Wang L, et al. Transcriptional regulation of cannabinoid receptor-1 expression in the liver by retinoic acid acting via retinoic acid receptor-gamma. J Biol Chem. 2010;285:19002–11. https://doi.org/10.1074/jbc.M109.068460.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  18. Doxakis E, Davies AM. Retinoic acid negatively regulates GDNF and neurturin receptor expression and responsiveness in embryonic chicken sympathetic neurons. Mol Cell Neurosci. 2005;29:617–27. https://doi.org/10.1016/j.mcn.2005.04.011.

    Article  CAS  PubMed  Google Scholar 

  19. Iyer N, Grizotte-Lake M, Duncan K, Gordon SR, Palmer ACS, Calvin C, et al. Epithelium intrinsic vitamin a signaling co-ordinates pathogen clearance in the gut via IL-18. PLoS Pathog. 2020;16:e1008360. https://doi.org/10.1371/journal.ppat.1008360.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  20. Guo L, Lin W, Zhang Y, Wang J. A systematic analysis revealed the potential gene regulatory processes of ATRA-triggered neuroblastoma differentiation and identified a novel RA response sequence in the NTRK2 gene. Biomed Res Int. 2020;2020:6734048. https://doi.org/10.1155/2020/6734048.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  21. Richard S, Zingg HH. Identification of a retinoic acid response element in the human oxytocin promoter. J Biol Chem. 1991;266:21428–33.

    Article  CAS  PubMed  Google Scholar 

  22. Yamada S, Uchimura E, Ueda T, Nomura T, Fujita S, Matsumoto K, et al. Identification of twinfilin-2 as a factor involved in neurite outgrowth by RNAi-based screen. Biochem Biophys Res Commun. 2007;363:926–30. https://doi.org/10.1016/j.bbrc.2007.09.069.

    Article  CAS  PubMed  Google Scholar 

  23. Pezzini F, Bettinetti L, Di Leva F, Bianchi M, Zoratti E, Carrozzo R, et al. Transcriptomic profiling discloses molecular and cellular events related to neuronal differentiation in SH-SY5Y neuroblastoma cells. Cell Mol Neurobiol. 2017;37:665–82. https://doi.org/10.1007/s10571-016-0403-y.

    Article  CAS  PubMed  Google Scholar 

  24. Isoherranen N, Zhong G. Biochemical and physiological importance of the CYP26 retinoic acid hydroxylases. Pharmacol Ther. 2019;204:107400. https://doi.org/10.1016/j.pharmthera.2019.107400.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  25. Haeseleer F, Huang J, Lebioda L, Saari JC, Palczewski K. Molecular characterization of a novel short-chain dehydrogenase/reductase that reduces all-trans-retinal. J Biol Chem. 1998;273:21790–9. https://doi.org/10.1074/jbc.273.34.21790.

    Article  CAS  PubMed  Google Scholar 

  26. Mattei MG, de Thé H, Mattei JF, Marchio A, Tiollais P, Dejean A. Assignment of the human hap retinoic acid receptor RAR beta gene to the p24 band of chromosome 3. Hum Genet. 1988;80:189–90. https://doi.org/10.1007/BF00702867.

    Article  CAS  PubMed  Google Scholar 

  27. Katsushima K, Shinjo K, Natsume A, Ohka F, Fujii M, Osada H, et al. Contribution of microRNA-1275 to Claudin11 protein suppression via a polycomb-mediated silencing mechanism in human glioma stem-like cells. J Biol Chem. 2012;287:27396–406. https://doi.org/10.1074/jbc.M112.359109.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  28. Telgenhoff D, Ramsay S, Hilz S, Slusarewicz P, Shroot B. Claudin 2 mRNA and protein are present in human keratinocytes and may be regulated by all-trans-retinoic acid. Skin Pharmacol Physiol. 2008;21:211–7. https://doi.org/10.1159/000135637.

    Article  CAS  PubMed  Google Scholar 

  29. Maden M. Retinoic acid in the development, regeneration and maintenance of the nervous system. Nat Rev Neurosci. 2007;8:755–65. https://doi.org/10.1038/nrn2212.

    Article  CAS  PubMed  Google Scholar 

  30. Janesick A, Wu SC, Blumberg B. Retinoic acid signaling and neuronal differentiation. Cell Mol Life Sci. 2015;72:1559–76. https://doi.org/10.1007/s00018-014-1815-9.

    Article  CAS  PubMed  Google Scholar 

  31. Tang X-H, Gudas LJ. Retinoids, retinoic acid receptors, and cancer. Annu Rev Pathol. 2011;6:345–64. https://doi.org/10.1146/annurev-pathol-011110-130303.

    Article  CAS  PubMed  Google Scholar 

  32. Amann PM, Eichmüller SB, Schmidt J, Bazhin AV. Regulation of gene expression by retinoids. Curr Med Chem. 2011;18:1405–12. https://doi.org/10.2174/092986711795029618.

    Article  CAS  PubMed  Google Scholar 

  33. Burrows K, Antignano F, Chenery A, Bramhall M, Korinek V, Underhill TM, Zaph C. HIC1 links retinoic acid signalling to group 3 innate lymphoid cell-dependent regulation of intestinal immunity and homeostasis. PLoS Pathog. 2018;14:e1006869. https://doi.org/10.1371/journal.ppat.1006869.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  34. Yorgan TA, Heckt T, Rendenbach C, Helmis C, Seitz S, Streichert T, et al. Immediate effects of retinoic acid on gene expression in primary murine osteoblasts. J Bone Miner Metab. 2016. https://doi.org/10.1007/s00774-015-0666-2.

  35. Kathiriya IS, King IN, Murakami M, Nakagawa M, Astle JM, Gardner KA, et al. Hairy-related transcription factors inhibit GATA-dependent cardiac gene expression through a signal-responsive mechanism. J Biol Chem. 2004;279:54937–43. https://doi.org/10.1074/jbc.M409879200.

    Article  CAS  PubMed  Google Scholar 

  36. Perl E, Waxman JS. Retinoic Acid Signaling and Heart Development. In: Asson-Batres MA, Rochette-Egly C, editors. The biochemistry of retinoid signaling. III, Vitamin A and retinoic acid in embryonic development [electronic resource]. Cham: Springer; 2020. p. 119–49. https://doi.org/10.1007/978-3-030-42282-0_5.

    Chapter  Google Scholar 

  37. Dulmage BO, Geskin LJ. Lessons learned from gene expression profiling of cutaneous T-cell lymphoma. Br J Dermatol. 2013;169:1188–97. https://doi.org/10.1111/bjd.12578.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  38. Kessler JF, Meyskens FL, Levine N, Lynch PJ, Jones SE. Treatment of cutaneous T-cell lymphoma (mycosis fungoides) with 13-cis-retinoic acid. Lancet. 1983;1:1345–7. https://doi.org/10.1016/s0140-6736(83)92136-0.

    Article  CAS  PubMed  Google Scholar 

  39. Ma H, Abdul-Hay M. T-cell lymphomas, a challenging disease: types, treatments, and future. Int J Clin Oncol. 2016;22:18–51. https://doi.org/10.1007/s10147-016-1045-2.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  40. Rowland A, Miners JO, Mackenzie PI. The UDP-glucuronosyltransferases: their role in drug metabolism and detoxification. Int J Biochem Cell Biol. 2013;45:1121–32. https://doi.org/10.1016/j.biocel.2013.02.019.

    Article  CAS  PubMed  Google Scholar 

  41. Napoli J, Race K. Microsomes convert retinol and retinal into retinoic acid and interference in the conversions catalyzed by cytosol. Biochim Biophys Acta (BBA). 1990;1034:228–32. https://doi.org/10.1016/0304-4165(90)90081-7.

    Article  CAS  Google Scholar 

  42. Franzese O, Minotti L, Aguiari G, Corrà F, Cervellati C, Ferrari C, et al. Involvement of non-coding RNAs and transcription factors in the induction of transglutaminase isoforms by ATRA. Amino Acids. 2019;51:1273–88. https://doi.org/10.1007/s00726-019-02766-7.

    Article  CAS  PubMed  Google Scholar 

  43. Xu J, Su Z, Ding Q, Shen L, Nie X, Pan X, et al. Inhibition of proliferation by knockdown of Transmembrane (TMEM) 168 in Glioblastoma cells via suppression of Wnt/β-catenin pathway. Oncol Res. 2019;27:819–26. https://doi.org/10.3727/096504018X15478559215014.

    Article  PubMed  PubMed Central  Google Scholar 

  44. Harel S, Cohen ASA, Hussain K, Flanagan SE, Schlade-Bartusiak K, Patel M, et al. Alternating hypoglycemia and hyperglycemia in a toddler with a homozygous p.R1419H ABCC8 mutation: an unusual clinical picture. J Pediatr Endocrinol Metab. 2015;28:345–51. https://doi.org/10.1515/jpem-2014-0265.

    Article  CAS  PubMed  Google Scholar 

  45. Nolte C, de Kumar B, Krumlauf R. Hox genes: downstream "effectors" of retinoic acid signaling in vertebrate embryogenesis. Genesis. 2019;57:e23306. https://doi.org/10.1002/dvg.23306.

    Article  PubMed  Google Scholar 

  46. Dickson IR, Walls J, Webb S. Vitamin a and bone formation. Different responses to retinol and retinoic acid of chick bone cells in organ culture. Biochim Biophys Acta (BBA). 1989;1013:254–8. https://doi.org/10.1016/0167-4889(89)90143-2.

    Article  CAS  Google Scholar 

  47. Hu L, Lind T, Sundqvist A, Jacobson A, Melhus H. Retinoic acid increases proliferation of human osteoclast progenitors and inhibits RANKL-stimulated osteoclast differentiation by suppressing RANK. PLoS One. 2010;5:e13305. https://doi.org/10.1371/journal.pone.0013305.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  48. Roa LA, Bloemen M, Carels CEL, Wagener FADTG, von den Hoff JW. Retinoic acid disrupts osteogenesis in pre-osteoblasts by down-regulating WNT signaling. Int J Biochem Cell Biol. 2019;116:105597. https://doi.org/10.1016/j.biocel.2019.105597.

    Article  CAS  PubMed  Google Scholar 

  49. Henning P, Conaway HH, Lerner UH. Retinoid receptors in bone and their role in bone remodeling. Front Endocrinol (Lausanne). 2015;6:31. https://doi.org/10.3389/fendo.2015.00031.

    Article  Google Scholar 

  50. Wei L-N. Retinoids and receptor interacting protein 140 (RIP140) in gene regulation. Curr Med Chem. 2004;11:1527–32. https://doi.org/10.2174/0929867043365017.

    Article  CAS  PubMed  Google Scholar 

  51. Herzog B, Hallberg M, Seth A, Woods A, White R, Parker MG. The nuclear receptor cofactor RIP140 is required for the regulation of hepatic lipid and glucose metabolism by LXR. Mol Endocrinol. 2007;21:2687–97. https://doi.org/10.1210/me.2007-0213.

    Article  CAS  PubMed  Google Scholar 

  52. White JA, Ramshaw H, Taimi M, Stangle W, Zhang A, Everingham S, et al. Identification of the human cytochrome P450, P450RAI-2, which is predominantly expressed in the adult cerebellum and is responsible for all-trans-retinoic acid metabolism. Proc Natl Acad Sci U S A. 2000;97:6403–8. https://doi.org/10.1073/pnas.120161397.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  53. Moretti A, Li J, Donini S, Sobol RW, Rizzi M, Garavaglia S. Crystal structure of human aldehyde dehydrogenase 1A3 complexed with NAD+ and retinoic acid. Sci Rep. 2016;6:35710. https://doi.org/10.1038/srep35710.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  54. Harada Y, Yokota C, Habas R, Slusarski DC, He X. Retinoic acid-inducible G protein-coupled receptors bind to frizzled receptors and may activate non-canonical Wnt signaling. Biochem Biophys Res Commun. 2007;358:968–75. https://doi.org/10.1016/j.bbrc.2007.04.208.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  55. Kam RKT, Deng Y, Chen Y, Zhao H. Retinoic acid synthesis and functions in early embryonic development. Cell Biosci. 2012;2:11. https://doi.org/10.1186/2045-3701-2-11.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  56. Elias KM, Laurence A, Davidson TS, Stephens G, Kanno Y, Shevach EM, O'Shea JJ. Retinoic acid inhibits Th17 polarization and enhances FoxP3 expression through a Stat-3/Stat-5 independent signaling pathway. Blood. 2008;111:1013–20. https://doi.org/10.1182/blood-2007-06-096438.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  57. Babina M, Guhl S, Motakis E, Artuc M, Hazzan T, Worm M, et al. Retinoic acid potentiates inflammatory cytokines in human mast cells: identification of mast cells as prominent constituents of the skin retinoid network. Mol Cell Endocrinol. 2015;406:49–59. https://doi.org/10.1016/j.mce.2015.02.019.

    Article  CAS  PubMed  Google Scholar 

  58. Lanciotti M, Longone P, Cornaglia-Ferraris P, Ponzoni M. Retinoic acid inhibits phosphatidylinositol turnover only in RA-sensitive while not in RA-resistant human neuroblastoma cells. Biochem Biophys Res Commun. 1989;161:284–9. https://doi.org/10.1016/0006-291x(89)91593-3.

    Article  CAS  PubMed  Google Scholar 

  59. Ponzoni M, Lanciotti M. Retinoic acid rapidly decreases phosphatidylinositol turnover during neuroblastoma cell differentiation. J Neurochem. 1990;54:540–6. https://doi.org/10.1111/j.1471-4159.1990.tb01905.x.

    Article  CAS  PubMed  Google Scholar 

  60. Obrochta KM, Krois CR, Campos B, Napoli JL. Insulin regulates retinol dehydrogenase expression and all-trans-retinoic acid biosynthesis through FoxO1. J Biol Chem. 2015;290:7259–68. https://doi.org/10.1074/jbc.M114.609313.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  61. Hsu Y-T, Li J, Wu D, Südhof TC, Chen L. Synaptic retinoic acid receptor signaling mediates mTOR-dependent metaplasticity that controls hippocampal learning. Proc Natl Acad Sci U S A. 2019;116:7113–22. https://doi.org/10.1073/pnas.1820690116.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  62. Antony JW, Ferreira CS, Norman KA, Wimber M. Retrieval as a Fast Route to Memory Consolidation. Trends Cogn Sci (Regul Ed ). 2017;21:573–6. https://doi.org/10.1016/j.tics.2017.05.001.

    Article  Google Scholar 

  63. Wu BX, Chen Y, Chen Y, Fan J, Rohrer B, Crouch RK, Ma J-X. Cloning and characterization of a novel all-trans retinol short-chain dehydrogenase/reductase from the RPE. Invest Ophthalmol Vis Sci. 2002;43:3365–72.

    PubMed  Google Scholar 

  64. Rowbotham SE, Illingworth NA, Daly AK, Veal GJ, Boddy AV. Role of UDP-glucuronosyltransferase isoforms in 13-cis retinoic acid metabolism in humans. Drug Metab Dispos. 2010;38:1211–7. https://doi.org/10.1124/dmd.109.031625.

    Article  CAS  PubMed  Google Scholar 

  65. Perreault M, Białek A, Trottier J, Verreault M, Caron P, Milkiewicz P, Barbier O. Correction: Role of Glucuronidation for Hepatic Detoxification and Urinary Elimination of Toxic Bile Acids during Biliary Obstruction. PLoS One. 2014. https://doi.org/10.1371/annotation/ef13ed51-6848-419d-94d8-1bb62e7bcf52.

  66. Czernik PJ, Little JM, Barone GW, Raufman JP, Radominska-Pandya A. Glucuronidation of estrogens and retinoic acid and expression of UDP-glucuronosyltransferase 2B7 in human intestinal mucosa. Drug Metab Dispos. 2000;28:1210–6.

    CAS  PubMed  Google Scholar 

  67. Burney IA, Islam MU, Khursheed M. Retinoic acid syndrome: a potentially fatal side effect of retinoic acid therapy. J Pak Med Assoc. 1998;48:54–5.

    CAS  PubMed  Google Scholar 

  68. Tariq Z, Phinney RC, Mohamed I. A case of life-threatening retinoic acid syndrome and review of literature. Am J Ther. 2014;21:e28–30. https://doi.org/10.1097/MJT.0b013e31822aeece.

    Article  PubMed  Google Scholar 

  69. van der Hem KG, Ossenkoppele GJ. All-trans retinoic acid toxicity. Eur J Haematol. 1992;49:148–9. https://doi.org/10.1111/j.1600-0609.1992.tb00921.x.

    Article  PubMed  Google Scholar 

  70. Falker-Gieske C, Mott A, Preuß S, Franzenburg S, Bessei W, Bennewitz J, Tetens J. Analysis of the brain transcriptome in lines of laying hens divergently selected for feather pecking. BMC Genomics. 2020;21:595. https://doi.org/10.1186/s12864-020-07002-1.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  71. Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics. 2014;30:2114–20. https://doi.org/10.1093/bioinformatics/btu170.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  72. Trapnell C, Pachter L, Salzberg SL. TopHat: discovering splice junctions with RNA-Seq. Bioinformatics. 2009;25:1105–11. https://doi.org/10.1093/bioinformatics/btp120.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  73. Lawrence M, Huber W, Pagès H, Aboyoun P, Carlson M, Gentleman R, et al. Software for computing and annotating genomic ranges. PLoS Comput Biol. 2013;9:e1003118. https://doi.org/10.1371/journal.pcbi.1003118.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  74. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15:550. https://doi.org/10.1186/s13059-014-0550-8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  75. Yu G, Wang L-G, Han Y, He Q-Y. clusterProfiler: an R package for comparing biological themes among gene clusters. OMICS. 2012;16:284–7. https://doi.org/10.1089/omi.2011.0118.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgments

The DFG Competence Centre for Genome Analysis Kiel (CCGA) is greatly acknowledged for the cooperation in RNA-sequencing. We thank Dr. Sebastian Giese and Prof. Dr. Martin Schwemmle from the Institute of Virology Universtity Freiburg for the kind gesture of providing the LMH cell line. We further acknowledge support by the Open Access Publication Funds of the Göttingen University.

Funding

The publication fee was provided by the Open Access Publication Funds of the Göttingen University. Open Access funding enabled and organized by Projekt DEAL.

Author information

Authors and Affiliations

  1. Department of Animal Sciences, Georg-August-University, Burckhardtweg 2, 37077, Göttingen, Germany

    Clemens Falker-Gieske, Andrea Mott & Jens Tetens

  2. Institute of Clinical Molecular Biology, Christian-Albrechts-University of Kiel, Rosalind-Franklin-Straße 12, 24105, Kiel, Germany

    Sören Franzenburg

  3. Center for Integrated Breeding Research, Georg-August-University, Albrecht-Thaer-Weg 3, 37075, Göttingen, Germany

    Jens Tetens

Authors
  1. Clemens Falker-Gieske
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Andrea Mott
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Sören Franzenburg
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Jens Tetens
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

CFG perfomed cell culture experiments, all bioinformatic analyses and wrote the manuscript. AM performed cell culture experiments and RNA isolation. SF performed the sequencing. JT wrote the manuscript. The authors read and approved the final manuscript. The authors read and approved the final manuscript.

Corresponding author

Correspondence to Clemens Falker-Gieske.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1.

Complete alignment metrics after mapping with TopHat.

Additional file 2.

Differential expression analyses results of all datasets after exposure to retinoic acid.

Additional file 3.

Venn diagram of differentially expressed genes from all datasets after exposure to retinoic acid.

Additional file 4.

Common differentially expressed genes among all datasets after exposure to retinoic acid.

Additional file 5.

Results of clusterProfiler analyses of differentially expressed genes from the meta-analysis and LMH cells exposed to retinoic acid and retinol for 1 h and 4 h.

Additional file 6.

Volcano plots of differentially expressed genes in LMH cells after exposure to retinoic acid and retinol for 1 h and 4 h.

Additional file 7.

Common differentially expressed genes after exposure of LMH cells to retinoic acid and retinol for 1 h and 4 h.

Additional file 8.

The numbers of RAREs in the vicinity of DE genes (up to 10 kb upstream of transcript start and 10 kb downstream of transcript end) in LMH cells after exposure to retinoic acid and retinol for 1 h and 4 h. LFCs for each gene and treatment are listed and results with a p-adj < 0.01 are indicated by *.

Additional file 9.

Protein interaction network analysis results of genes that were differentially expressed in LMH cells after 4 h exposure to retinoic acid or retinol.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Falker-Gieske, C., Mott, A., Franzenburg, S. et al. Multi-species transcriptome meta-analysis of the response to retinoic acid in vertebrates and comparative analysis of the effects of retinol and retinoic acid on gene expression in LMH cells. BMC Genomics 22, 146 (2021). https://doi.org/10.1186/s12864-021-07451-2

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12864-021-07451-2

Keywords

BMC Genomics

ISSN: 1471-2164

Contact us