Skip to main content

Transcriptomics and miRNomics data integration in lymphoblastoid cells highlights the key role of immune-related functions in lithium treatment response in Bipolar disorder

Abstract

Background

Bipolar Disorder (BD) is a complex mental disease characterized by recurrent episodes of mania and depression. Lithium (Li) represents the mainstay of BD pharmacotherapy, despite the narrow therapeutic index and the high variability in treatment response. However, although several studies have been conducted, the molecular mechanisms underlying Li therapeutic effects remain unclear.

Methods

In order to identify molecular signatures and biological pathways associated with Li treatment response, we conducted transcriptome and miRNome microarray analyses on lymphoblastoid cell lines (LCLs) from 20 patients diagnosed with BD classified as Li responders (n = 11) or non-responders (n = 9).

Results

We found 335 mRNAs and 77 microRNAs (miRNAs) significantly modulated in BD responders versus non-responders. Interestingly, pathway and network analyses on these differentially expressed molecules suggested a modulatory effect of Li on several immune-related functions. Indeed, among the functional molecular nodes, we found NF-κB and TNF. Moreover, networks related to these molecules resulted overall inhibited in BD responder patients, suggesting anti-inflammatory properties of Li.

From the integrative analysis between transcriptomics and miRNomics data carried out using miRComb R package on the same samples from patients diagnosed with BD, we found 97 significantly and negatively correlated mRNA-miRNA pairs, mainly involved in inflammatory/immune response.

Conclusions

Our results highlight that Li exerts modulatory effects on immune-related functions and that epigenetic mechanisms, especially miRNAs, can influence the modulation of different genes and pathways involved in Li response. Moreover, our data suggest the potentiality to integrate data coming from different high-throughput approaches as a tool to prioritize genes and pathways.

Peer Review reports

Introduction

Bipolar disorder (BD) is a severe and disabling psychiatric condition characterized by intermitting states of mania and depression, affecting 1–3% of the population worldwide. The onset of the illness is around 20 years of age and it is associated with a reduced quality of life, substantial societal costs and the highest suicide rate among psychiatric disorders [1,2,3]. Management of patients diagnosed with BD requires both acute treatment of manic or hypomanic episodes above to a maintenance therapy to prevent relapses and further episodes.

Mood stabilizers are used as the first-line therapy in the treatment of BD, and among this class of compounds, lithium (Li), introduced by John Cade in 1949, still represents the gold standard treatment for stabilization, prophylaxis and suicide prevention [4, 5]. However, patients diagnosed with BD receiving Li need regular monitoring due to Li narrow therapeutic index and the risk to develop multiple side effects [5, 6]. Further complexity in the management of Li therapy is represented by the heterogeneity in Li response: approximately 30% of patients appear good responders, whereas 70% of them are classified as partial or non-responders [7,8,9].

A growing body of evidence suggests that the response to Li prophylaxis has a strong genetic background and familiar heritability [10,11,12,13,14,15,16]. Recent findings from the largest genome-wide association study (GWAS) conducted by the International Consortium on Lithium Genetics (ConLiGen) have indeed suggested the involvement of two long non-coding RNA (lncRNA) genes, AL157359.3 and AL157359.4, in Li response [7].

However, so far, findings from pharmacogenetic studies have been able to explain only a small proportion of the observed variability, suggesting that other factors could be involved, including for example epigenetic mechanisms, such as DNA methylation, histone modification and microRNAs (miRNAs). MiRNAs are single-stranded, non-coding RNA molecules, 18–25 nucleotides-long, with a key role in the regulation of messenger RNA (mRNA). Indeed, they act by inducing either degradation or translational silencing of their target mRNA, but in certain instances, miRNAs may activate translation or even act at the level of transcription by binding to specific gene promoters [17,18,19,20]. MiRNAs are involved in key processes such as neurogenesis [21], neural plasticity and higher brain functioning [22] and they have been also associated with neurodegenerative disorders [23, 24] and psychiatric illnesses [25]. Interestingly, although still limited, available data on miRNAs expression support the hypothesis that these small non-coding RNAs could influence the treatment response in patients diagnosed with BD [8, 26, 27].

Among human cellular models, lymphoblastoid cell lines (LCLs) represent a valid and useful experimental tool to study Li response in patients diagnosed with BD, including epigenetic mechanisms, as already demonstrated [28,29,30,31,32,33,34]. Despite some limitations associated with the Epstein-Barr Virus (EBV) transformation, LCLs show several advantages [35]. Indeed, the genomes of LCLs remain stable during subsequent cell divisions. This stability in part reflects the fact that the EBV genome is not incorporated into the germ-line genome but rather remains in the cell cytosol [35]. Interestingly, although some of the epigenomic signatures are lost during the EBV immortalization process and the subsequent in vitro propagation, many of them persist and may be useful to study disease-related epigenomic modifications [36,37,38,39,40,41,42].

Although Li treatment response in patients diagnosed with BD has been studying for decades, findings from already available studies have been elusive mainly due to: i) the complexity of the mechanism of action of Li [43, 44], ii) the involvement of multiple molecular processes [45, 46], and iii) the use of different animal models or in vitro human cell lines (for a review see [8, 34]). Moreover, the heterogeneity among all studies (in terms of analyzed samples, technologies and methods used, inclusion/exclusion criteria of patients’ selection, sample size), the lack of independent replication cohorts to validate initial findings and of standard operative procedures among different laboratories (i.e. methods of collecting, storing, processing and analyzing markers), as well as the high-throughput approaches across different studies have produced big data sets, contributing to make the picture even more complex to understand.

To overcome all these limitations and to reduce heterogeneity among the studies, new successful strategies have been introduced by combining modalities coming from different high-throughput methods, which can be useful to prioritize a stringent number of genes and pathways involved in Li treatment response. For example, Hunsberger and colleagues have conducted an integrative approach combining miRNAs and mRNAs expression profiling of LCLs derived from BD Li responder and non-responder patients before and after in vitro Li treatment [30]. Interestingly, the results of this study suggested that in vitro Li treatment down-regulated miRNA Let-7 family in both groups, but specifically in the group of BD responders [30]. However, these interesting findings have not been replicated by further studies [47]. Another example of high-throughput data integration has been proposed by Pisanu and colleagues who applied a convergent analysis of genome-wide genotyping and transcriptomic data from LCLs of patients diagnosed with BD identifying a zinc finger gene, ZNF493, as a potential Li-responsive target [48].

Within this scenario, in the present study, we aimed to investigate the effects of Li treatment response in LCLs obtained from patients diagnosed with BD characterized for Li response. Thus, we performed an integrative analysis on both mRNAs and miRNAs microarray data by combining differential expression analyses and correlation of expression levels with miRComb, an R package able to combine miRNA and mRNA expression data with hybridization information [49]. Our goal was to identify molecular signatures and biological pathways that could underlie the effects of Li therapeutic response. This could have a clinical utility helping clinicians to develop a personalized medicine approach. Moreover, we aimed to find biological signatures (mRNAs, miRNAs) involved in Li response not only to better understand the molecular mechanisms associated with Li therapeutic treatment, but also to identify novel targets for the development of new drugs.

Materials and methods

Patients

Twenty patients with a diagnosis of Bipolar Disorder type I were recruited at the Expert Center for bipolar disorder of Paris (France) as part of a research protocol (Clinical Trials Number NCT02627404) approved by the ethical committee (Comité de Protection des Personnes – La Pitié- Salpétrière Hospital – Paris – France) (reference: P111002-IDRCB2008-AO1465–50). All participants provided written informed consent prior to inclusion.

Diagnosis was made according to the Diagnostic and Statistical Manual of Mental Disorders (DSM-IV Edition), and the Li response was evaluated by using the Retrospective Criteria of Long-Term Treatment Response in Research Subjects with BD score, also known as ALDA scale.

Briefly, the ALDA scale is composed of two subscales: the A scale and the B scale. The A scale evaluates the improvement of symptoms during the treatment with Li, while the B scale describes five clinical factors that may potentially confound the pharmacological response.

The ALDA score assumes an integer value between 0 and 10 calculated by subtracting the B from the A score. Subjects with ALDA score ≥ 7 are classified as “responders”, whereas those with ALDA score < 7 represent “non-responders” [50]. Accordingly, in our study patients diagnosed with BD were defined as Li responders (R, n = 11) or non-responders (NR, n = 9).

Demographical information, details of the clinical characterization and assessment of Li response of all patients diagnosed with BD have been reported in Table 1.

Table 1 Demographical and clinical information of BD responder (R) and non-responder (NR) patients

Lymphoblastoid cell lines

Lymphoblastoid cell lines (LCLs) from patients diagnosed with BD were established from fresh blood by transforming lymphocytes with Epstein-Barr virus (EBV) following standard protocols [28, 51]. When reaching the confluence state, each cell line was stored in liquid nitrogen until use. For the present study, LCLs were thawed and cultured in RPMI-1640 medium supplemented with 10% foetal bovine serum, 1% penicillin/streptomycin, 1% L-Glutamine 200 mM (Thermo Fisher Scientific, USA) in a 5% CO2 humidified incubator at 37 °C. LCLs were seeded at 2 × 105 cells/ml and grown in T75 Flask. After 4 days, cells were harvested for RNA isolation.

All the cells were tested free of mycoplasma by using TransDetect® PCR Mycoplasma Detection Kit (Clinisciences, France).

RNA extraction

From each cell line, total RNA, including miRNAs, was isolated from 5 × 106 cell pellets with the miRNeasy mini kit according to the manufacturer’s protocol (QIAGEN, Hilden, Germany). RNA quantity and quality were determined by using a NanoDrop-1000 spectrophotometer (NanoDrop Technologies, Thermo Fisher Scientific, USA). RNA purity was considered adequate when the A260/280 ratio was in the range of 1.8–2.0 and when the A260/230 ratio was in the range of 2.0–2.2. RNA integrity number, assessed using the Agilent 2100 Bioanalyzer (Agilent, Santa Clara, CA, USA), was in the range of 7–10.

Microarrays

Transcriptome analysis

For the whole transcriptomic profiling, 250 ng of total RNA were processed with the WT PLUS Reagent Kit (Thermo Fisher Scientific, USA) and were subsequently hybridized onto the GeneChip Human Gene 2.1 ST Array Strips on a GeneAtlas platform (Thermo Fisher Scientific, USA). The comprehensive coverage of these array strips allows the simultaneously evaluation of > 30.000 coding transcripts, of > 11.000 long intergenic non-coding transcripts and alternative splicing events/transcript variants with probes designed to maximize exon coverage. Samples have been randomized and stratified in a way that each array strip included the same number of Li responder and Li non-responder patients diagnosed with BD. Washing and staining were conducted on the Fluidics station; scanning was performed on the Imaging station. All procedures were run following the manufacturer’s instructions.

miRNome analysis

250 ng of total RNA including miRNAs were processed with the FlashTag Biotin HSR RNA Labeling kit (Thermo Fisher Scientific, USA) and subsequently hybridized onto the GeneChip miRNA 4.1 Array Strip on a GeneAtlas platform (Thermo Fisher Scientific, USA). As reported before, samples have been randomized and stratified in a way that each array strip included the same number of Li responder and Li non-responder patients diagnosed with BD. Washing and staining were conducted on the Fluidics station, whereas scanning was performed on the Imaging station. All procedures were run following the manufacturer’s instructions.

Statistics

Microarray data analysis

Strip array data quality control was checked by GeneAtlas integrated analysis software (Thermo Fisher Scientific, USA) for both GeneChip Human Gene 2.1 ST Array Strips and GeneChip miRNA 4.1 Array Strips.

Briefly, raw microarray data (.CEL files) were imported and analyzed with the commercially available software Partek® Genomic Suite® (Partek, St. Louis, MO, USA). Probe set normalization and summarization were computed using Robust Multi-Array (RMA) algorithm. We checked for quality control and batch effects using the Principal Component Analyses (PCA), and we did not observe any outliers or batch effect.

Analysis of variance (ANOVA) was performed to assess, separately, the expression levels of both mRNAs and miRNAs in LCLs obtained from BD responders versus non-responders.

To identify the list of differentially expressed mRNAs and miRNAs, we filtered ANOVA results using both fold change |FC| ≥ 1.2 and p-value ≤0.05. Moreover, in the list of miRNAs we selected the probe sets only for the Human organism.

Pathway analysis of differentially expressed mRNAs

Differentially expressed genes with RefSeq annotations were analyzed using the Core Analysis functionality in Ingenuity Pathway Analysis (IPA) software (QIAGEN, Hilden, Germany) to describe key aspects of the molecular function, biological process and cellular component of gene products. This software is based on a proprietary database able to identify specific cellular pathways and to generate biological networks.

Pathway analysis of differentially expressed mature miRNAs

The identification of potentially altered molecular pathways as targeted by differentially modulated miRNAs was performed by using DIANA-miRPath (v.3) tool web-server [52]. Specifically, an in silico target prediction was performed on differentially expressed mature miRNAs by using DIANA-microT-CDS algorithm, which combines the analysis of 3′-untranslated regions (3′-UTRs) and coding sequence regions [53]. MiRNA-gene interactions were filtered using the default threshold equal to 0.8 and the False Discovering Rate (FDR) algorithm (FDR cut-off = 0.05) was applied.

An enrichment analysis of miRNA target genes was executed in Kyoto Encyclopedia of Genes and Genomes (KEGG) database and a list of significant pathways (p-value ≤0.05) showing all targeted genes and number of miRNAs was computed [52].

mRNAs-miRNAs combined analysis

In order to detect the relevant pairs of interacting mRNA-miRNA for Li response, an integrated analysis of mRNAs and miRNAs expression levels was performed. First, differentially expressed mRNAs and miRNAs between BD responders and non-responders were filtered to focus on gene expression regulation mechanisms associated with Li response. Among all identified mRNAs and miRNAs, a correlation analysis was performed to detect mRNA-miRNA pairs that were inversely correlated based on the biological assumption that a miRNA negatively regulates the expression of its mRNA targets. The analysis has been performed using miRComb, an R package able to combine miRNA and mRNA expression data with hybridization information [49].

We first used LIMMA package to perform differential analysis and to filter for differentially expressed mRNAs and miRNAs (|FC| ≥ 1.2 and nominal p-value ≤0.05). The intensity signals of microarray (after RMA normalization) among all the couples of identified differentially expressed mRNAs and miRNAs have been tested for association with Pearson correlation (that is to test for the presence of an inverse linear relationship between a given couple of mRNA-miRNA).

We also computed an S score = − 2 (logratiomiRNA x logratiomRNA) to have a weigh within each couple of mRNA and miRNA, representing the level of alteration among Li responders and non-responders (the higher the score the higher the deregulation). We considered a positive S score (> 0), to take into account only mRNA-miRNA pairs within which mRNA and miRNA showed an opposite direction of expression, in terms of FC [49].

We extracted the most relevant mRNA-miRNA pairs by filtering for Pearson correlation coefficient ρ ≤ − 0.38, showing an anti-correlation trend, a nominal p-value ≤0.05 and a positive S score.

Results

Transcriptomic signatures associated with lithium response

Our first aim was to identify changes at the transcriptome level that could be associated with Li response. From the comparison between BD responder and non-responder patients, we identified 335 differentially modulated genes (217 were up-regulated and 118 down-regulated) (Supplementary Table 1 for the entire list of significant genes) in LCLs of Li responders.

To strengthen our analysis, we compared our results to already available transcriptomic studies by performing a bibliographic search in PubMed, using the following keywords: “lithium; bipolar disorder; peripheral blood; gene expression” OR “lithium; bipolar disorder; lymphoblastoid; gene expression”. We found 15 studies conducted in peripheral blood samples and 21 performed in LCLs obtained from patients diagnosed with BD, responders or non-responders to Li, but also from controls. When we looked for common genes within our list and already available studies, we found an overlap of 42 Li-responsive genes modulated in the same direction upon treatment with Li (Supplementary Table 2) [29, 30, 32, 33, 54,55,56,57].

mRNAs pathway and network analyses

To identify possible deregulated biological processes underlying Li treatment response, we performed a pathway and network analysis on the list of 335 differentially expressed genes. We found 54 significantly modulated pathways represented in the pie chart (Fig. 1) and in Supplementary Table 3. As shown in the pie chart, 57% of the pathways are involved in cellular immune response, including Th1 and Th2 pathways, OX40 signalling, B cell development, T helper cell differentiation, Antigen Presentation, Calcium-induced T Lymphocyte Apoptosis, Inflammasome pathway and Nur77 signalling in T Lymphocytes. Other pathways relevant to BD pathophysiology and Li treatment response are related to neurotransmitters and other nervous system signalling (5%) (i.e. ErbB signalling), cell cycle regulation (4%) (i.e. Cdc42 and 14–3-3 mediated signalling) and metabolism (4%) (i.e. Type I Diabetes Mellitus Signalling).

Fig. 1
figure 1

Pie chart gathering 54 significant pathways, identified in patients diagnosed with BD who respond to Li, according to their biological functions

In addition, to predict how genes in our data set could interact with each other, we performed a network analysis on the same list of 335 differentially expressed genes. Interestingly, IPA revealed 21 relevant functional networks (Supplementary Table 4). As shown in Fig. 2, Nuclear Factor kappa B (NF-κB), Tumor Necrosis Factor (TNF) and Signal Transducer and Activator of Transcription 3 (STAT3) might represent important functional nodes associated with Li treatment response.

Fig. 2
figure 2

Networks implicated in response to Li treatment in BD responders. NF-kB, TNF and STAT3 are among the hub genes

Furthermore, to predict directional effects of the networks, we used Molecule Activity Predictor (MAP) and Path designer functionalities in IPA, which take into account the observed expression changes to compute functional effects on neighbouring molecules. Networks showing NF-κB, TNF and STAT3 among the hub genes resulted globally inhibited in patients diagnosed with BD who responded to Li treatment (Fig. 2).

Results of eQTL-mapping studies from blood and/or brain tissues

To better investigate the effects of Li response on gene expression levels in BD responders and non-responders, we performed an in-silico analysis based on the matching between expression Quantitative Trait Loci (eQTL) and the genetic signals available from BD GWAS [58].

Our aim was to identify to which extent the genetic regulation of gene expression was associated with the genetic signal of the analyzed phenotype, represented by BD symptomatology, hypothesizing that Li treatment in BD responders could counterbalance alterations in gene expression levels relevant for BD. For our purpose, we retrieved the list of differentially regulated genes in blood and brain tissues (according to GTEx eQTL data) considering BD GWAS summary statistics from the Psychiatric Genomics Consortium [59]. We took into account only those genes that reached significance (nominally significant p-value < 0.05). We then compared this list with the whole data set of differentially expressed genes between BD responder and non-responder patients found in our transcriptome microarray analysis (Table 2). Specifically, we were interested in those genes with an opposite directionality between case/control and Li response. The only gene that survived this analysis was fasciculation and elongation protein zeta (FEZ1). A significant association by TWAS suggests that a down-regulation of FEZ1 expression levels increases the risk for BD (TWAS Z = − 2.15; p-value = 0.0318). Interestingly, in our analysis, we found a significant increase of FEZ1 mRNA levels in BD responders versus non-responders (FC = 1.51; p-value = 0.0314), suggesting FEZ1 as a potential target of Li therapeutic effects.

Table 2 Common genes between GTEx data and differentially expressed genes in our analysis. TWAS: Transcriptome-Wide Association Studies; TWAS. Z is an estimation of the strength of association between the predicted expression of a gene and a complex trait. TWAS. P is the p-value associated with the TWAS test. On the left (a), genes that reached the significance in reference data set from GTEx; on the right (b), genes in common between GTEx eQTL data and our list of differentially expressed genes

miRNomic signatures associated with Li response

In order to evaluate which miRNAs could be differentially modulated by Li, we also performed a miRNome microarray analysis in LCLs. We found a list of 77 mature miRNAs that were significantly modulated in BD responders versus non-responders. Out of these, 46 were up-regulated and 31 down-regulated (Supplementary Table 5). To verify whether these differentially expressed mature miRNAs have been already associated with Li response, we conducted a search in PubMed using the following keywords “lithium; bipolar disorder; miRNAs” and we found 17 articles. We identified only two common miRNAs modulated in the same direction (up-regulated in LCLs of patients diagnosed with BD following chronic Li treatment) between our list and previous studies: hsa-miR-34a-5p (FC = 1.26, p-value = 0.0378) and hsa-miR-152-3p (FC = 1.54, p-value = 0.007), which were consistent with the literature [60].

MiRNAs pathway analysis

In order to investigate the biological systems regulated by differentially expressed miRNAs in BD responder patients, we conducted a pathway analysis by using Diana Tool software (miRPath v.3). Interestingly, we identified a list of 77 statistically significant pathways (Table 3 and Supplementary Table 6) mainly involved in: i) neurodevelopment, such as Hippo signalling pathway, Axon guidance, Glutamatergic synapse, Neurotrophin signalling pathway, FoxO signalling pathway, mTOR signalling pathway, GABAergic synapse, Cholinergic synapse, Wnt signalling pathway and Dopaminergic synapse; ii) intracellular signal transduction (i.e. ErbB signalling pathway, Ras signalling pathway, Rap1 signalling pathway, MAPK signalling pathway, cAMP signalling pathway, AMPK signalling pathway, PI3K-Akt signalling pathway and cGMP-PKG signalling pathway); iii) inflammatory/immune system response, including TGF-beta signalling pathway and T cell receptor signalling pathway.

Table 3 List of 77 statistically significant pathways targeted by differentially modulated miRNAs identified using Diana Tool software (miRPath v.3)

Predicted mRNA-miRNA interactions

The innovative approach of our study mainly consists in the integration analysis between transcriptomic and miRNomic data carried out on the same patients diagnosed with BD. Briefly, by using miRComb R package [49], we combined transcriptomic and miRNomic expression data with microarray hybridization information to compute all possible correlations between deregulated and differentially modulated mRNAs and miRNAs. Our goal was to highlight those mRNA-miRNA interactions that may play an important role in the outcome of Li treatment response.

We thus selected the previously mentioned 335 and 77 significantly deregulated mRNAs and miRNAs, respectively, and we computed all possible correlations. Among all mRNA-miRNA possible interactions, 7627 had a p-value ≤0.05. We considered only mature mRNA-miRNA pairs and we eliminated the isoforms of the same gene, obtaining 5450 significant mature mRNA-miRNA interactions with a Pearson correlation index ρ ≤ − 0.38 (Supplementary Table 7). In order to further reduce the number of interactions, we used a cut-off of p-value ≤0.05 and filtered mRNA-miRNA pairs by Pearson correlation index ρ ≤ − 0.5. We identified a list of 2027 significant mature mRNA-miRNA interactions (Supplementary Table 8). Moreover, as miRComb does not take into account competitivity among different miRNAs targeting the same gene, we also used the S score (> 0), an index representing the level of alteration among mRNAs and miRNAs, to further prioritize the mRNA-miRNA pairs in those cases where: i) mRNA and miRNA within the pair showed an opposite direction of expression and ii) different miRNAs shared the same mRNA target: higher score means that both mRNA and miRNA are highly deregulated in Li responders. We obtained 97 mRNA-miRNA interactions representing those that are more likely to occur in the context of good response to Li treatment (Table 4). As an example, we show some plots of the most significant pairs in that table, namely hsa-miR-574-3p and Ring Finger protein 125 (RNF125) (ρ = − 0.83), hsa-miR-3128 and ELAV Like RNA Binding Protein 3 (ELAVL3) (ρ = − 0.72), Creatine Kinase, Mitochondrial 1A (CKMT1A) (ρ = − 0.71), Keratin Associated Protein 9–1 (KRTAP9–1) (ρ = − 0.70), hsa-miR-3201 and Histone Cluster 3, H2a (HIST3H2A) (ρ = − 0.69) (Fig. 3).

Table 4 Top 97 mRNA-miRNA pairs with the most significant negative correlations predicted by miRComb, ranked by the Pearson correlation index and the S score. All interactions show a p-value ≤ 0.05, a Pearson correlation index ρ ≤ -0.5 and an S score > 0
Fig. 3
figure 3

Top 5 mRNA-miRNA miRComb interactions. Panel A shows the interaction between hsa-miR-574-3p and RNF125. Panel B shows interactions between hsa-miR-3128, ELAVL3, CKMT1A and KRTAP9–1, respectively. Panel C shows the correlation between hsa-miR-3201 and HIST3H2A. Li responders are depicted in blue circles, whereas non-responders are shown in red triangles. The graphs show the correlation between normalized values of log2 intensity signal. Regression line shows the negative correlation

Discussion

In this exploratory study, we investigated the effects of Li response in LCLs derived from patients diagnosed with BD, responders or non-responders to therapeutic treatment, to further elucidate the biological and molecular processes involved in Li therapeutic action. We used an integrative approach on both mRNAs and miRNAs microarray data by combining differential expression analyses with correlation of expression levels. The premise is that the expression changes of specific mRNAs and miRNAs will reflect alterations in specific target genes and biological pathways that underlie therapeutic action, thus, we sought to identify differentially regulated mRNAs and miRNAs in LCLs from BD responder and non-responder patients. This molecular signature may allow the identification of biomarkers able to shed light on both Li therapeutic effects and BD pathophysiology.

Through an exploratory transcriptome microarray analysis, we identified 335 differentially expressed genes between BD responder and non-responder patients and 42 of them were in common with previously published transcriptomic studies performed in LCLs or in peripheral blood samples of patients diagnosed with BD [29, 30, 32, 33, 54,55,56,57], supporting the reliability and reproducibility of our analysis. Among all common genes, the most studied and recurring ones were represented by Branched Chain Amino Acid Transaminase 1 (BCAT1), which was up-regulated in our data set and in three studies, all performed in LCLs [33, 54, 55] and by RAB11 Family Interacting Protein 1 (RAB11FIP1) whose up-regulation was reported by our study and by three different works, all performed in LCLs as well [29, 33, 55].

The pathway and network analyses performed on the list of 335 differentially expressed genes revealed that they are mainly involved in cellular immune response and, to a lesser extent, in neurotransmitters and other nervous system signalling, cell cycle regulation and metabolism. In support of these findings, our data suggest that the therapeutic effects of Li mainly affect the immune system and the inflammatory response: indeed, functional molecular nodes with a biological role in the treatment response to Li are represented by NF-κB, TNF and STAT3, known for their pro-inflammatory role in mental disorders [61, 62]. Interestingly, networks related to these molecules resulted globally inhibited in patients diagnosed with BD who responded to Li treatment, confirming the anti-inflammatory properties of Li [62].

Worth of note is that top canonical pathways identified in our study are also represented among the most significantly ones obtained from the cross-trait meta-GWAS and pathway analysis based on GWAS summary statistics and treatment response to Li coming from GWAS of schizophrenia and from ConLiGen, respectively [63]. Briefly, in this study, the authors tested whether a polygenic score for schizophrenia was associated with Li treatment response in patients affected by BD and explored the potential molecular underpinnings of this association [63]. According to their results, patients diagnosed with BD, who had a low polygenic load for schizophrenia, responded better to Li treatment. Moreover, the authors found that genetic variants in the Human Leukocyte Antigens (HLA) genes, the antigen presentation pathway, and inflammatory cytokines such as TNF, IL-4 and IFN-γ could represent central and functional hub genes, suggesting that pathways associated with immunity and inflammation could play a key role in Li treatment response [63].

In this regard, previous studies have reported modulatory anti-inflammatory effects of Li and highlighted the possibility that mechanisms involving pro-inflammatory cytokines might play a role in mediating the responsiveness of Li in patients suffering from BD (for a review see [62]). These mechanisms are not clearly understood, but it is commonly thought to be dependent upon the inhibition of Glycogen synthase kinase-3β (GSK-3β), the well-established pharmacological target of Li. GSK-3β, a serine/threonine kinase ubiquitously distributed in mammalian tissues, is able to enhance inflammation and immune responses [64] through the activation of NF-κB (for a review see [62, 65]). In a simplified model, Li might reduce the production of pro-inflammatory mediators by inhibiting GSK-3β, thus resulting in a decreased activity of NF-κB, which leads to an attenuated expression of inflammatory-associated molecules and enzymes [62].

According to this hypothesis, Guloksuz and colleagues found higher plasma levels of TNF-α in patients diagnosed with BD characterized for a poor response to Li compared to those showing a good response, suggesting that increased TNF-α levels may affect the clinical response to Li [66].

Furthermore, from the analysis of tissue specific eQTL (GTEx) and GWAS association data [59], we found that decreased expression levels of FEZ1, a well-recognized risk factor for schizophrenia [67, 68] also involved in neuronal development, neuropathology, and viral infections [69], increases the risk for BD. Interestingly, in our transcriptome analysis we observed a significant increase in the expression levels of FEZ1 in BD responders to Li treatment. Although no eQTL and GWAS data are already available for FEZ1 and Li treatment response, we can hypothesize that, if validated and replicated in larger cohorts of BD responder patients, this gene might be involved in Li therapeutic effects.

Using a similar approach, we also identified 77 miRNAs whose expression differentially changed between BD Li responders and non-responders. These findings indicate that epigenetic modifications may play a key role in modulating the effect of mood stabilizers, shed light on their molecular targets and possibly explain the high heterogeneity in the treatment response to Li. Among the most significant biological processes modulated by these miRNAs in BD responders, we identified several pathways involved in neurodevelopment, as the WNT signalling, highlighting the key role of this pathway in Li response, and in inflammatory/immune system response, corroborating transcriptomic findings.

Following the same study design previously used for the transcriptome analysis, we found that, in line with previous studies [60, 70, 71], hsa-miR-34a-5p and hsa-miR-152-3p were up-regulated in BD responders compared to non-responders. Interestingly, high levels of miR-34a have been found in postmortem cerebellar tissue from patients diagnosed with BD, as well as in BD patient-derived neuronal cultures generated by reprogramming human fibroblasts into induced pluripotent stem cells (iPSCs) subsequently differentiated into neurons [72]. Although contrasting results exist, several studies have suggested that Li has a modulatory effect on miR-34a expression levels [73, 74], which have been found consistently up-regulated in LCLs following Li treatment [60]. In non-neuronal cells, miR-34a has been shown to be a strong inhibitor of the WNT signalling and β-catenin-mediated transcription in response to p53 activation [75, 76]. Given the well-known effects of Li on the WNT signalling, through the inhibition of GSK-3β, a key molecule within this pathway [77, 78], we can speculate that miR-34a might contribute to the therapeutic effects of Li [60, 72,73,74, 79, 80].

In addition to miR-34a, we found that miR-152-3p expression levels were consistently up-regulated in LCLs from Li responders. To our knowledge, no study has investigated the effect of Li on the modulation of this miRNA. However, miR-152-3p targets a DNA methyltransferase (DNMT1) and could interfere with gene expression levels and with Li treatment response by acting through DNA methylation. Indeed, it is well established that Li monotherapy induces a global DNA hypomethylation in leukocytes and converging evidence from clinical and preclinical studies has suggested that different mood stabilizers may exert specific effects on the epigenome [81, 82]. According to all these findings, our hypothesis, which, however, needs to be further examined, is that Li response could be, at least in part, mediated by epigenetic mechanisms, which may interact each other in a miRNA-epigenetics regulatory circuit, influencing the expression of relevant genes possibly involved in Li response.

However, the core analysis of the present work is based on an innovative approach that enables correlation analysis between transcriptome and miRNome data, simultaneously carried out on the same patients, in order to highlight those mRNA-miRNA interactions that may play an important role in the outcome of Li treatment response. To date, only two studies have already integrated miRNAs and mRNAs expression profiling in LCLs from BD responder or non responder patients [30, 83]. However, their main findings have not been replicated by our study. Although the rationale for the experimental design is based on a similar hypothesis, these studies show differences with our study. Indeed, Hunsberger and collaborators applied GRANITE, an integrative genomic systems biology-based tool, to genome-wide mRNA and miRNA expression data obtained from LCLs of BD excellent responders or non-responders to Li treatment. They found that the Let-7 miRNA family was consistently downregulated by Li in the BD responder group. Conversely to our study, LCLs were cultured for 7 days with either a therapeutic dose of Li or vehicle in order to highlight genetic networks differentially influenced by the treatment in the responder and non-responder patients [30]. In the other study, Pisanu and colleagues sequenced small non-coding RNAs with next generation sequencing (NGS) in LCLs from patients diagnosed with BD and characterized for Li response. Then, they integrated these data with differentially expressed mRNAs identified by microarray-based techniques in the same subject. The results suggest that miR-320a, miR-155-3p and their target genes might represent relevant players in modulating clinical response to Li [83]. In our study, we used a microarray-based approach for both transcriptomic and miRNomic analyses, taking advantage of the same microarray platform (GeneAtlas platform, Thermo Fisher Scientific, USA). Therefore, the use of the same platform for both analyses (mRNAs and miRNAs) has allowed us to follow simpler and leaner bioinformatics and data integration analyses. In addition, while Pisanu and collaborators focused their attention only on the most promising genes, we decided to perform a more complex approach by investigating biological pathways modulated by differentially expressed mRNAs and miRNAs, for a better understanding of the biological systems involved in Li treatment response.

MiRNAs that most likely occur in our data set of mRNA-miRNA interactions are hsa-miR-574-3p, hsa-miR-3128 and hsa-miR-3201. Bibliographic references about these miRNAs are scarce, and mostly highlight a putative role in cancer or cardiovascular disorders [84, 85] as miRNAs could regulate a broad range of biological processes, like cell cycle, apoptosis, stress response, differentiation, proliferation and angiogenesis among the others [73, 86, 87]. Thus, it is possible that their involvement in such processes will in part explain the well-known neuroprotective and anti-apoptotic effects of Li [44, 77, 78, 88].

Regarding mRNAs potentially targeted by these miRNAs, we found genes previously associated with Li response. Among others, we found Inositol Monophosphatase 2 (IMPA2), which has been implicated in BD and Li response [44, 77, 89, 90], HLA genes that have been recently associated with Li response in BD [63] and Phosphatidylinositol-4,5-Bisphosphate 3-Kinase Catalytic Subunit Gamma (PIK3CG) whose expression has been already found to be modulated by Li treatment in LCLs of patients diagnosed with BD [33].

In order to identify mRNA-miRNA interactions that may play an important role in Li treatment response, we filtered the entire list of significant mature mRNA-miRNA correlations considering the strength of correlation and we obtained 97 mRNA-miRNA pairs. Interestingly, the best pair was represented by hsa-miR-574-3p and RNF125, because it showed the strongest negative correlation coefficient among all mRNA-miRNA interactions.

RNF125, also named T cell RING protein in activation (TRAC-1), is an E3 ubiquitin ligase that contains a RING finger domain in the N-terminus and three zinc-binding and one ubiquitin-interacting motif in the C-terminus. RNF125 protein may function as a positive regulator in the T-cell receptor signalling pathway [91,92,93]. Activation of T cells results in a pro-inflammatory response necessary to prevent the spread of infection. Limiting T cell signalling, however, is essential to prevent this protective response from causing injury to the host. Findings from our analysis suggest a down-regulation of RNF125, probably mediated by hsa-miR-574-3p that could result in a mitigation of immune response.

Therefore, all data seem to converge on the anti-inflammatory effects of Li in BD responders, supporting the hypothesis that pathways associated with cellular immune response could play a biological role in Li treatment response.

We are aware that our findings should be interpreted with caution in light of several limitations. In our study, we are not able to discriminate whether the observed effect is specifically associated with BD or is a general effect of Li treatment response. Indeed, in our dataset we have a few missing data and therefore we cannot provide information about the length of Li treatment. Moreover, we cannot exclude that the lack or the presence of Li treatment may have interfered with the results.

Conclusions

Taken together, our results strengthen the evidence that Li may act in modulating the inflammatory signaling and immune functions, and provide further evidence for the involvement of miRNAs in Li response. However, future studies should be conducted in larger cohorts of BD subjects by integrating data from different high-throughput approaches to offer a successful strategy for the identification of molecular signatures related both to mechanism of action and efficacy of Li. This would enable the personalization of treatment, define criteria for predicting whether patients diagnosed with BD will respond or not to Li, and improve their long-term management and prognosis, thus reducing the risk for suicidal behaviours.

Availability of data and materials

The datasets used and/or analysed during the current study are available from the corresponding author on reasonable request. We are willing to make the data public, but we are internally setting up the appropriate procedure, and in any case we are happy to share the data with researchers if they are interested, once the manuscript is published.

Abbreviations

BD:

Bipolar disorder

Li:

Lithium

GWAS:

Genome-wide association study

ConLiGen:

International Consortium on Lithium Genetics

lncRNA:

long non-coding RNA

miRNAs:

microRNAs

mRNA:

messenger RNA

LCLs:

Lymphoblastoid cell lines

EBV:

Epstein-Barr Virus

R:

Responders

NR:

Non-responders

RMA algorithm:

Robust Multi-Array

PCA:

Principal Component Analysis

ANOVA:

Analysis of variance

IPA:

Ingenuity Pathway Analysis

3′-UTRs:

3′-untranslated regions

FDR algorithm:

False Discovering Rate

KEGG database:

Kyoto Encyclopedia of Genes and Genomes

NF-κB:

Nuclear Factor kappa B

Akt/ERK:

AKT serine/threonine kinase 1/ extracellular regulated MAP kinase

TNF:

Tumor Necrosis Factor

MAP:

Molecule Activity Predictor

eQTL:

expression Quantitative Trait Loci

FEZ1:

Fasciculation and elongation protein zeta

RNF125:

Ring Finger protein 125

ELAVL3:

ELAV Like RNA Binding Protein 3

CKMT1A:

Creatine Kinase, Mitochondrial 1A

KRTAP9–1:

Keratin Associated Protein 9–1

HIST3H2A:

Histone Cluster 3, H2a

BCAT1:

Branched Chain Amino Acid Transaminase 1

RAB11FIP1:

RAB11 Family Interacting Protein 1

HLA:

Human Leukocyte Antigens

GSK-3β:

Glycogen synthase kinase-3β

iPSCs:

induced pluripotent stem cells

DNMT1:

DNA methyltransferase

IMPA2:

Inositol Monophosphatase 2

PIK3CG:

Phosphatidylinositol-4,5-Bisphosphate 3-Kinase Catalytic Subunit Gamma

TRAC-1:

T cell RING protein in activation

References

  1. Vigo D, Thornicroft G, Atun R. Estimating the true global burden of mental illness. Lancet Psychiatry. 2016;3(2):171–8. https://doi.org/10.1016/S2215-0366(15)00505-2.

    Article  PubMed  Google Scholar 

  2. Vieta E, Berk M, Schulze TG, Carvalho AF, Suppes T, Calabrese JR, et al. Bipolar disorders. Nature Reviews Disease Primers. 2018;4:18008. https://doi.org/10.1038/nrdp.2018.8.

  3. Grande I, Berk M, Birmaher B, Vieta E. Bipolar disorder. Lancet. 2016;387(10027):1561–72. https://doi.org/10.1016/S0140-6736(15)00241-X.

    Article  PubMed  Google Scholar 

  4. Won E, Kim YK. An Oldie but Goodie: Lithium in the treatment of Bipolar disorder through neuroprotective and neurotrophic mechanisms. Int J Mol Sci. 2017;18(12). https://doi.org/10.3390/ijms18122679.

  5. Geddes JR, Miklowitz DJ. Treatment of bipolar disorder. Lancet. 2013;381(9878):1672–82. https://doi.org/10.1016/S0140-6736(13)60857-0.

    Article  CAS  PubMed  Google Scholar 

  6. McKnight RF, Adida M, Budge K, Stockton S, Goodwin GM, Geddes JR. Lithium toxicity profile: a systematic review and meta-analysis. Lancet. 2012;379(9817):721–8. https://doi.org/10.1016/S0140-6736(11)61516-X.

    Article  CAS  PubMed  Google Scholar 

  7. Hou L, Heilbronner U, Degenhardt F, Adli M, Akiyama K, Akula N, et al. Genetic variants associated with response to lithium treatment in bipolar disorder: a genome-wide association study. Lancet. 2016;387(10023):1085–93. https://doi.org/10.1016/S0140-6736(16)00143-4.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  8. Bellivier F, Marie-Claire C. Molecular signatures of Lithium treatment: current knowledge. Pharmacopsychiatry. 2018;51(5):212–9. https://doi.org/10.1055/a-0650-4820.

    Article  CAS  PubMed  Google Scholar 

  9. Pisanu C, Katsila T, Patrinos GP, Squassina A. Recent trends on the role of epigenomics, metabolomics and noncoding RNAs in rationalizing mood stabilizing treatment. Pharmacogenomics. 2018;19(2):129–43. https://doi.org/10.2217/pgs-2017-0111.

    Article  CAS  PubMed  Google Scholar 

  10. Alda M, Grof P, Rouleau GA, Turecki G, Young LT. Investigating responders to lithium prophylaxis as a strategy for mapping susceptibility genes for bipolar disorder. Prog Neuro-Psychopharmacol Biol Psychiatry. 2005;29(6):1038–45. https://doi.org/10.1016/j.pnpbp.2005.03.021.

    Article  CAS  Google Scholar 

  11. Cruceanu C, Alda M, Turecki G. Lithium: a key to the genetics of bipolar disorder. Genome medicine. 2009;1(8):79. https://doi.org/10.1186/gm79.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  12. Cruceanu C, Ambalavanan A, Spiegelman D, Gauthier J, Lafreniere RG, Dion PA, et al. Family-based exome-sequencing approach identifies rare susceptibility variants for lithium-responsive bipolar disorder. Genome. 2013;56(10):634–40. https://doi.org/10.1139/gen-2013-0081.

    Article  CAS  PubMed  Google Scholar 

  13. Geoffroy PA, Bellivier F, Leboyer M, Etain B. Can the response to mood stabilizers be predicted in bipolar disorder? Front Biosci. 2014;6:120–38. https://doi.org/10.2741/e696.

    Article  Google Scholar 

  14. Pickard BS. Genomics of Lithium action and response. Neurotherapeutics. 2017;14(3):582–7. https://doi.org/10.1007/s13311-017-0554-7.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  15. Budde M, Degner D, Brockmoller J, Schulze TG. Pharmacogenomic aspects of bipolar disorder: an update. Eur Neuropsychopharmacol. 2017;27(6):599–609. https://doi.org/10.1016/j.euroneuro.2017.02.001.

    Article  CAS  PubMed  Google Scholar 

  16. Pisanu C, Heilbronner U, Squassina A. The role of pharmacogenomics in Bipolar disorder: moving towards precision medicine. Mol Diagn Ther. 2018;22(4):409–20. https://doi.org/10.1007/s40291-018-0335-y.

    Article  CAS  PubMed  Google Scholar 

  17. Hu Y, Ehli EA, Boomsma DI. MicroRNAs as biomarkers for psychiatric disorders with a focus on autism spectrum disorder: current progress in genetic association studies, expression profiling, and translational research. Autism Res. 2017;10(7):1184–203. https://doi.org/10.1002/aur.1789.

    Article  PubMed  Google Scholar 

  18. Younger ST, Corey DR. Transcriptional gene silencing in mammalian cells by miRNA mimics that target gene promoters. Nucleic Acids Res. 2011;39(13):5682–91. https://doi.org/10.1093/nar/gkr155.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  19. Dwivedi Y. Pathogenetic and therapeutic applications of microRNAs in major depressive disorder. Prog Neuro-Psychopharmacol Biol Psychiatry. 2016;64:341–8. https://doi.org/10.1016/j.pnpbp.2015.02.003.

    Article  CAS  Google Scholar 

  20. Mora C, Zonca V, Riva MA, Cattaneo A. Blood biomarkers and treatment response in major depression. Expert Rev Mol Diagn. 2018;18(6):513–29. https://doi.org/10.1080/14737159.2018.1470927.

    Article  CAS  PubMed  Google Scholar 

  21. Rajman M, Schratt G. MicroRNAs in neural development: from master regulators to fine-tuners. Development. 2017;144(13):2310–22. https://doi.org/10.1242/dev.144337.

    Article  CAS  PubMed  Google Scholar 

  22. Wei CW, Luo T, Zou SS, Wu AS. Research progress on the roles of microRNAs in governing synaptic plasticity, learning and memory. Life Sci. 2017;188:118–22. https://doi.org/10.1016/j.lfs.2017.08.033.

    Article  CAS  PubMed  Google Scholar 

  23. Qiu L, Tan EK, Zeng L. microRNAs and neurodegenerative diseases. Adv Exp Med Biol. 2015;888:85–105. https://doi.org/10.1007/978-3-319-22671-2_6.

    Article  CAS  PubMed  Google Scholar 

  24. Viswambharan V, Thanseem I, Vasu MM, Poovathinal SA, Anitha A. miRNAs as biomarkers of neurodegenerative disorders. Biomark Med. 2017;11(2):151–67. https://doi.org/10.2217/bmm-2016-0242.

    Article  CAS  PubMed  Google Scholar 

  25. Luoni A, Riva MA. MicroRNAs and psychiatric disorders: from aetiology to treatment. Pharmacol Ther. 2016;167:13–27. https://doi.org/10.1016/j.pharmthera.2016.07.006.

    Article  CAS  PubMed  Google Scholar 

  26. Pisanu C, Papadima EM, Del Zompo M, Squassina A. Understanding the molecular mechanisms underlying mood stabilizer treatments in bipolar disorder: potential involvement of epigenetics. Neurosci Lett. 2018;669:24–31. https://doi.org/10.1016/j.neulet.2016.06.045.

    Article  CAS  PubMed  Google Scholar 

  27. Fries GR, Carvalho AF, Quevedo J. The miRNome of bipolar disorder. J Affect Disord. 2018;233:110–6. https://doi.org/10.1016/j.jad.2017.09.025.

    Article  CAS  PubMed  Google Scholar 

  28. Sie L, Loong S, Tan EK. Utility of lymphoblastoid cell lines. J Neurosci Res. 2009;87(9):1953–9. https://doi.org/10.1002/jnr.22000.

    Article  CAS  PubMed  Google Scholar 

  29. Breen MS, White CH, Shekhtman T, Lin K, Looney D, Woelk CH, et al. Lithium-responsive genes and gene networks in bipolar disorder patient-derived lymphoblastoid cell lines. Pharmacogenomics J. 2016;16(5):446–53. https://doi.org/10.1038/tpj.2016.50.

    Article  CAS  PubMed  Google Scholar 

  30. Hunsberger JG, Chibane FL, Elkahloun AG, Henderson R, Singh R, Lawson J, et al. Novel integrative genomic tool for interrogating lithium response in bipolar disorder. Transl Psychiatry. 2015;5:e504. https://doi.org/10.1038/tp.2014.139.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  31. Geoffroy PA, Curis E, Courtin C, Moreira J, Morvillers T, Etain B, et al. Lithium response in bipolar disorders and core clock genes expression. World J Biol Psychiatry. 2018;19(8):619–32. https://doi.org/10.1080/15622975.2017.1282174.

    Article  PubMed  Google Scholar 

  32. Sugawara H, Iwamoto K, Bundo M, Ishiwata M, Ueda J, Kakiuchi C, et al. Effect of mood stabilizers on gene expression in lymphoblastoid cells. J Neural Transm. 2010;117(2):155–64. https://doi.org/10.1007/s00702-009-0340-8.

    Article  CAS  PubMed  Google Scholar 

  33. Fries GR, Colpo GD, Monroy-Jaramillo N, Zhao J, Zhao Z, Arnold JG, et al. Distinct lithium-induced gene expression effects in lymphoblastoid cell lines from patients with bipolar disorder. Eur Neuropsychopharmacol. 2017;27(11):1110–9. https://doi.org/10.1016/j.euroneuro.2017.09.003.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  34. Viswanath B, Jose SP, Squassina A, Thirthalli J, Purushottam M, Mukherjee O, et al. Cellular models to study bipolar disorder: a systematic review. J Affect Disord. 2015;184:36–50. https://doi.org/10.1016/j.jad.2015.05.037.

    Article  PubMed  Google Scholar 

  35. Gurwitz D. Human iPSC-derived neurons and lymphoblastoid cells for personalized medicine research in neuropsychiatric disorders. Dialogues Clin Neurosci. 2016;18(3):267–76.

    Article  PubMed  PubMed Central  Google Scholar 

  36. Wang L, Laing J, Yan B, Zhou H, Ke L, Wang C, et al. Epstein-Barr virus Episome physically interacts with active regions of the host genome in Lymphoblastoid cells. J Virol. 2020;94(24). https://doi.org/10.1128/JVI.01390-20.

  37. Harris-Arnold A, Arnold CP, Schaffert S, Hatton O, Krams SM, Esquivel CO, et al. Epstein-Barr virus modulates host cell microRNA-194 to promote IL-10 production and B lymphoma cell survival. Am J Transplant. 2015;15(11):2814–24. https://doi.org/10.1111/ajt.13375.

    Article  CAS  PubMed  Google Scholar 

  38. Forte E, Salinas RE, Chang C, Zhou T, Linnstaedt SD, Gottwein E, et al. The Epstein-Barr virus (EBV)-induced tumor suppressor microRNA MiR-34a is growth promoting in EBV-infected B cells. J Virol. 2012;86(12):6889–98. https://doi.org/10.1128/JVI.07056-11.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  39. Skalsky RL, Corcoran DL, Gottwein E, Frank CL, Kang D, Hafner M, et al. The viral and cellular microRNA targetome in lymphoblastoid cell lines. PLoS Pathog. 2012;8(1):e1002484. https://doi.org/10.1371/journal.ppat.1002484.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  40. Lee JE, Hong EJ, Nam HY, Kim JW, Han BG, Jeon JP. MicroRNA signatures associated with immortalization of EBV-transformed lymphoblastoid cell lines and their clinical traits. Cell Prolif. 2011;44(1):59–66. https://doi.org/10.1111/j.1365-2184.2010.00717.x.

    Article  CAS  PubMed  Google Scholar 

  41. Linnstaedt SD, Gottwein E, Skalsky RL, Luftig MA, Cullen BR. Virally induced cellular microRNA miR-155 plays a key role in B-cell immortalization by Epstein-Barr virus. J Virol. 2010;84(22):11670–8. https://doi.org/10.1128/JVI.01248-10.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  42. Leong MML, Lung ML. The impact of Epstein-Barr virus infection on epigenetic regulation of host cell gene expression in epithelial and lymphocytic malignancies. Front Oncol. 2021;11:629780. https://doi.org/10.3389/fonc.2021.629780.

    Article  PubMed  PubMed Central  Google Scholar 

  43. Malhi GS, Outhred T. Therapeutic mechanisms of Lithium in Bipolar disorder: recent advances and current understanding. CNS drugs. 2016;30(10):931–49. https://doi.org/10.1007/s40263-016-0380-1.

    Article  CAS  PubMed  Google Scholar 

  44. Alda M. Lithium in the treatment of bipolar disorder: pharmacology and pharmacogenetics. Mol Psychiatry. 2015;20(6):661–70. https://doi.org/10.1038/mp.2015.4.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  45. Robertson OD, Coronado NG, Sethi R, Berk M, Dodd S. Putative neuroprotective pharmacotherapies to target the staged progression of mental illness. Early Interv Psychiatry. 2019;13(5):1032–49. https://doi.org/10.1111/eip.12775.

    Article  PubMed  Google Scholar 

  46. Jakobsson E, Arguello-Miranda O, Chiu SW, Fazal Z, Kruczek J, Nunez-Corrales S, et al. Towards a unified understanding of Lithium action in basic biology and its significance for applied biology. J Membr Biol. 2017;250(6):587–604. https://doi.org/10.1007/s00232-017-9998-2.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  47. Milanesi E, Hadar A, Maffioletti E, Werner H, Shomron N, Gennarelli M, et al. Insulin-like growth factor 1 differentially affects Lithium sensitivity of Lymphoblastoid cell lines from Lithium responder and non-responder Bipolar disorder patients. J Mol Neurosci. 2015;56(3):681–7. https://doi.org/10.1007/s12031-015-0523-8.

    Article  CAS  PubMed  Google Scholar 

  48. Pisanu C, Congiu D, Costa M, Chillotti C, Ardau R, Severino G, et al. Convergent analysis of genome-wide genotyping and transcriptomic data suggests association of zinc finger genes with lithium response in bipolar disorder. Am J Med Genet B Neuropsychiatr Genet. 2018;177(7):658–64. https://doi.org/10.1002/ajmg.b.32663.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  49. Vila-Casadesus M, Gironella M, Lozano JJ. MiRComb: An R package to analyse miRNA-mRNA interactions Examples across Five Digestive Cancers. PloS one 2016;11(3):e0151127. doi:https://doi.org/10.1371/journal.pone.0151127.

  50. Grof P, Duffy A, Cavazzoni P, Grof E, Garnham J, MacDougall M, et al. Is response to prophylactic lithium a familial trait? The Journal of clinical psychiatry. 2002;63(10):942–7. https://doi.org/10.4088/jcp.v63n1013.

    Article  CAS  PubMed  Google Scholar 

  51. Neitzel H. A routine method for the establishment of permanent growing lymphoblastoid cell lines. Hum Genet. 1986;73(4):320–6. https://doi.org/10.1007/BF00279094.

    Article  CAS  PubMed  Google Scholar 

  52. Vlachos IS, Kostoulas N, Vergoulis T, Georgakilas G, Reczko M, Maragkakis M, et al. DIANA miRPath v.2.0: investigating the combinatorial effect of microRNAs in pathways. Nucleic acids Res. 2012;40(Web Server issue):W498–504. https://doi.org/10.1093/nar/gks494.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  53. Reczko M, Maragkakis M, Alexiou P, Grosse I, Hatzigeorgiou AG. Functional microRNA targets in protein coding sequences. Bioinformatics. 2012;28(6):771–6. https://doi.org/10.1093/bioinformatics/bts043.

    Article  CAS  PubMed  Google Scholar 

  54. Milanesi E, Voinsky I, Hadar A, Srouji A, Maj C, Shekhtman T, et al. RNA sequencing of bipolar disorder lymphoblastoid cell lines implicates the neurotrophic factor HRP-3 in lithium's clinical efficacy. World J Biol Psychiatry. 2019;20(6):449–61. https://doi.org/10.1080/15622975.2017.1372629.

    Article  PubMed  Google Scholar 

  55. Kittel-Schneider S, Hilscher M, Scholz CJ, Weber H, Grunewald L, Schwarz R, et al. Lithium-induced gene expression alterations in two peripheral cell models of bipolar disorder. World J Biol Psychiatry. 2019;20(6):462–75. https://doi.org/10.1080/15622975.2017.1396357.

    Article  PubMed  Google Scholar 

  56. Watanabe S, Iga J, Nishi A, Numata S, Kinoshita M, Kikuchi K, et al. Microarray analysis of global gene expression in leukocytes following lithium treatment. Human Psychopharmacol. 2014;29(2):190–8. https://doi.org/10.1002/hup.2381.

    Article  CAS  Google Scholar 

  57. Squassina A, Costa M, Congiu D, Manchia M, Angius A, Deiana V, et al. Insulin-like growth factor 1 (IGF-1) expression is up-regulated in lymphoblastoid cell lines of lithium responsive bipolar disorder patients. Pharmacol Res. 2013;73:1–7. https://doi.org/10.1016/j.phrs.2013.04.004.

    Article  CAS  PubMed  Google Scholar 

  58. Cloney R. Complex traits: integrating gene variation and expression to understand complex traits. Nat Rev Genet. 2016;17(4):194. https://doi.org/10.1038/nrg.2016.18.

    Article  CAS  PubMed  Google Scholar 

  59. Bipolar D. Schizophrenia working Group of the Psychiatric Genomics Consortium. Electronic address drve, Bipolar D, schizophrenia working Group of the Psychiatric Genomics C. genomic dissection of Bipolar disorder and schizophrenia, including 28 subphenotypes. Cell. 2018;173(7):1705–15 e16. https://doi.org/10.1016/j.cell.2018.05.046.

    Article  CAS  Google Scholar 

  60. Chen H, Wang N, Burmeister M, McInnis MG. MicroRNA expression changes in lymphoblastoid cell lines in response to lithium treatment. Int J Neuropsychopharmacol. 2009;12(7):975–81. https://doi.org/10.1017/S1461145709000029.

    Article  CAS  PubMed  Google Scholar 

  61. Garcia-Juarez M, Camacho-Morales A. Defining the role of anti- and pro-inflammatory outcomes of Interleukin-6 in mental health. Neuroscience. 2022;492:32–46. https://doi.org/10.1016/j.neuroscience.2022.03.020.

    Article  CAS  PubMed  Google Scholar 

  62. Nassar A, Azab AN. Effects of lithium on inflammation. ACS Chem Neurosci. 2014;5(6):451–8. https://doi.org/10.1021/cn500038f.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  63. International Consortium on Lithium G, Amare AT, Schubert KO, Hou L, Clark SR, Papiol S, et al. Association of Polygenic Score for schizophrenia and HLA antigen and inflammation genes with response to Lithium in Bipolar affective disorder: a genome-wide association study. JAMA Psychiatry. 2018;75(1):65–74. https://doi.org/10.1001/jamapsychiatry.2017.3433.

    Article  Google Scholar 

  64. Beurel E, Grieco SF, Jope RS. Glycogen synthase kinase-3 (GSK3): regulation, actions, and diseases. Pharmacol Ther. 2015;148:114–31. https://doi.org/10.1016/j.pharmthera.2014.11.016.

    Article  CAS  PubMed  Google Scholar 

  65. Brasier AR. The NF-kappaB regulatory network. Cardiovasc Toxicol. 2006;6(2):111–30. https://doi.org/10.1385/ct:6:2:111.

    Article  CAS  PubMed  Google Scholar 

  66. Guloksuz S, Altinbas K, Aktas Cetin E, Kenis G, Bilgic Gazioglu S, Deniz G, et al. Evidence for an association between tumor necrosis factor-alpha levels and lithium response. J Affect Disord. 2012;143(1–3):148–52. https://doi.org/10.1016/j.jad.2012.04.044.

    Article  CAS  PubMed  Google Scholar 

  67. Chen X, Ku L, Mei R, Liu G, Xu C, Wen Z, et al. Novel schizophrenia risk factor pathways regulate FEZ1 to advance oligodendroglia development. Transl Psychiatry. 2017;7(12):1293. https://doi.org/10.1038/s41398-017-0028-z.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  68. Vachev TI, Stoyanova VK, Ivanov HY, Minkov IN, Popov NT. Investigation of fasciculation and elongation protein zeta-1 (FEZ1) in peripheral blood reveals differences in gene expression in patients with schizophrenia. Balkan J Med Genet. 2015;18(1):31–8. https://doi.org/10.1515/bjmg-2015-0003.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  69. Maturana AD, Fujita T, Kuroda S. Functions of fasciculation and elongation protein zeta-1 (FEZ1) in the brain. ScientificWorldJournal. 2010;10:1646–54. https://doi.org/10.1100/tsw.2010.151.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  70. Banigan MG, Kao PF, Kozubek JA, Winslow AR, Medina J, Costa J, et al. Differential expression of exosomal microRNAs in prefrontal cortices of schizophrenia and bipolar disorder patients. PLoS One. 2013;8(1):e48814. https://doi.org/10.1371/journal.pone.0048814.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  71. Lenox RH, Wang L. Molecular basis of lithium action: integration of lithium-responsive signaling and gene expression networks. Mol Psychiatry. 2003;8(2):135–44. https://doi.org/10.1038/sj.mp.4001306.

    Article  CAS  PubMed  Google Scholar 

  72. Bavamian S, Mellios N, Lalonde J, Fass DM, Wang J, Sheridan SD, et al. Dysregulation of miR-34a links neuronal development to genetic risk factors for bipolar disorder. Mol Psychiatry. 2015;20(5):573–84. https://doi.org/10.1038/mp.2014.176.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  73. Alural B, Ozerdem A, Allmer J, Genc K, Genc S. Lithium protects against paraquat neurotoxicity by NRF2 activation and miR-34a inhibition in SH-SY5Y cells. Front Cell Neurosci. 2015;9:209. https://doi.org/10.3389/fncel.2015.00209.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  74. Tufekci KU, Alural B, Tarakcioglu E, San T, Genc S. Lithium inhibits oxidative stress-induced neuronal senescence through miR-34a. Mol Biol Rep. 2021;48(5):4171–80. https://doi.org/10.1007/s11033-021-06430-w.

    Article  CAS  PubMed  Google Scholar 

  75. Hashimi ST, Fulcher JA, Chang MH, Gov L, Wang S, Lee B. MicroRNA profiling identifies miR-34a and miR-21 and their target genes JAG1 and WNT1 in the coordinate regulation of dendritic cell differentiation. Blood. 2009;114(2):404–14. https://doi.org/10.1182/blood-2008-09-179150.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  76. Kim NH, Kim HS, Kim NG, Lee I, Choi HS, Li XY, et al. p53 and microRNA-34 are suppressors of canonical Wnt signaling. Sci Signal. 2011;4(197):ra71. https://doi.org/10.1126/scisignal.2001744.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  77. Machado-Vieira R, Manji HK, Zarate CA Jr. The role of lithium in the treatment of bipolar disorder: convergent evidence for neurotrophic effects as a unifying hypothesis. Bipolar Disord. 2009;11(Suppl 2):92–109. https://doi.org/10.1111/j.1399-5618.2009.00714.x.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  78. Malhi GS, Tanious M, Das P, Coulston CM, Berk M. Potential mechanisms of action of lithium in bipolar disorder Current understanding. CNS drugs. 2013;27(2):135–53. https://doi.org/10.1007/s40263-013-0039-0.

    Article  PubMed  Google Scholar 

  79. Zhou R, Yuan P, Wang Y, Hunsberger JG, Elkahloun A, Wei Y, et al. Evidence for selective microRNAs and their effectors as common long-term targets for the actions of mood stabilizers. Neuropsychopharmacology. 2009;34(6):1395–405. https://doi.org/10.1038/npp.2008.131.

    Article  CAS  PubMed  Google Scholar 

  80. Hunsberger JG, Fessler EB, Chibane FL, Leng Y, Maric D, Elkahloun AG, et al. Mood stabilizer-regulated miRNAs in neuropsychiatric and neurodegenerative diseases: identifying associations and functions. Am J Transl Res. 2013;5(4):450–64.

    CAS  PubMed  PubMed Central  Google Scholar 

  81. Backlund L, Wei YB, Martinsson L, Melas PA, Liu JJ, Mu N, et al. Mood stabilizers and the influence on global leukocyte DNA methylation in Bipolar disorder. Mol Neuropsychiatry. 2015;1(2):76–81. https://doi.org/10.1159/000430867.

    Article  PubMed  PubMed Central  Google Scholar 

  82. Huzayyin AA, Andreazza AC, Turecki G, Cruceanu C, Rouleau GA, Alda M, et al. Decreased global methylation in patients with bipolar disorder who respond to lithium. Int J Neuropsychopharmacol. 2014;17(4):561–9. https://doi.org/10.1017/S1461145713001569.

    Article  CAS  PubMed  Google Scholar 

  83. Pisanu C, Merkouri Papadima E, Melis C, Congiu D, Loizedda A, Orru N, et al. Whole genome expression analyses of miRNAs and mRNAs suggest the involvement of miR-320a and miR-155-3p and their targeted genes in Lithium response in Bipolar disorder. Int J Mol Sci. 2019;20(23). https://doi.org/10.3390/ijms20236040.

  84. Cui G, Li Z, Li R, Huang J, Wang H, Zhang L, et al. A functional variant in APOA5/A4/C3/A1 gene cluster contributes to elevated triglycerides and severity of CAD by interfering with microRNA 3201 binding efficiency. J Am Coll Cardiol. 2014;64(3):267–77. https://doi.org/10.1016/j.jacc.2014.03.050.

    Article  CAS  PubMed  Google Scholar 

  85. Veit JA, Scheckenbach K, Schuler PJ, Laban S, Wiggenhauser PS, Thierauf J, et al. MicroRNA expression in differentially metastasizing tumors of the head and neck: adenoid cystic versus squamous cell carcinoma. Anticancer Res. 2015;35(3):1271–7.

    CAS  PubMed  Google Scholar 

  86. Ha TY. MicroRNAs in human diseases: from Cancer to cardiovascular disease. Immune Net. 2011;11(3):135–54. https://doi.org/10.4110/in.2011.11.3.135.

    Article  Google Scholar 

  87. Tufekci KU, Meuwissen RL, Genc S. The role of microRNAs in biological processes. Methods Mol Biol. 2014;1107:15–31. https://doi.org/10.1007/978-1-62703-748-8_2.

    Article  CAS  PubMed  Google Scholar 

  88. Jope RS, Nemeroff CB. Lithium to the rescue. Cerebrum. 2016;2016.

  89. Atack JR, Broughton HB, Pollack SJ. Inositol monophosphatase--a putative target for Li+ in the treatment of bipolar disorder. Trends Neurosci. 1995;18(8):343–9. https://doi.org/10.1016/0166-2236(95)93926-o.

    Article  CAS  PubMed  Google Scholar 

  90. Seelan RS, Parthasarathy LK, Parthasarathy RN. Lithium modulation of the human inositol monophosphatase 2 (IMPA2) promoter. Biochem Biophys Res Commun. 2004;324(4):1370–8. https://doi.org/10.1016/j.bbrc.2004.09.199.

    Article  CAS  PubMed  Google Scholar 

  91. Giannini AL, Gao Y, Bijlmakers MJ. T-cell regulator RNF125/TRAC-1 belongs to a novel family of ubiquitin ligases with zinc fingers and a ubiquitin-binding domain. The Biochemical journal. 2008;410(1):101–11. https://doi.org/10.1042/BJ20070995.

    Article  CAS  PubMed  Google Scholar 

  92. Zhao H, Li CC, Pardo J, Chu PC, Liao CX, Huang J, et al. A novel E3 ubiquitin ligase TRAC-1 positively regulates T cell activation. J Immunol. 2005;174(9):5288–97. https://doi.org/10.4049/jimmunol.174.9.5288.

    Article  CAS  PubMed  Google Scholar 

  93. Arimoto K, Takahashi H, Hishiki T, Konishi H, Fujita T, Shimotohno K. Negative regulation of the RIG-I signaling by the ubiquitin ligase RNF125. Proc Natl Acad Sci U S A. 2007;104(18):7500–5. https://doi.org/10.1073/pnas.0611551104.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

We thank the clinicians and nursing staff for the recruitment and clinical characterization of individuals with BD. We thank the Centre de ressources biologiques_site cochin in Paris where LCLs were established.

Funding

N.C. and A.C. received funding from Ricerca Corrente (Italian Ministry of Health, MoH). A.C. has been also supported by PSR (University of Milan).

Author information

Authors and Affiliations

Authors

Contributions

N.C., C.C., E.M. and Cr.M. contributed to the experiments in the in vitro model. B.E. and F.B. contributed to the patients’ recruitment. N.C., E.M. and Cr.M. managed the literature searches. N.C., E.M. and Ca.M. performed the bioinformatics and statistical analyses. N.C. contributed to the first draft of the manuscript. C.C., E.M., Ca.M., Cr.M., B.E., F.B. and C.M-C. contributed to the revision. A.C. revised all the versions of the manuscript and approved the final one. All the authors have contributed to and have approved the final manuscript.

Corresponding author

Correspondence to Annamaria Cattaneo.

Ethics declarations

Ethics approval and consent to participate

This research has been approved by the ethical committee (Comité de Protection des Personnes – La Pitié- Salpétrière Hospital – Paris – France) (reference: P111002-IDRCB2008-AO1465–50). All participants provided written informed consent prior to inclusion. All methods were carried out in accordance with relevant guidelines and regulations.

Consent for publication

Not applicable.

Competing interests

The authors declare no conflict of interest.

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 Supplementary Table 1.

List of 335 differentially expressed mRNAs identified in Li responders (R) versus non-responders (NR) with |FC| ≥ 1.2 and p-value ≤0.05.

Additional file 2 Supplementary Table 2.

List of 57 common transcripts between those differentially expressed in our dataset (BD responders (R) versus non-responders (NR)) and other available transcriptomic studies in Pubmed. *: LCLs from Li-responders versus Li non-responders patients diagnosed with BD. #: Li-treated vs vehicle-treated LCLs from patients diagnosed with BD. Δ: before Li-treatment versus after Li-treatment in healthy controls. +: Li-treated vs vehicle-treated LCLs from Li-responders and Li non-responders patients diagnosed with BD. •: Li-treated LCLs versus vehicle-treated LCLs from Li-responders patients diagnosed with BD. □: Li-treated versus vehicle-treated LCLs from healthy controls.

Additional file 3 Supplementary Table 3.

List of 54 statistically significant pathways modulated by 335 differentially expressed genes identified using IPA. Ratio is calculated as the number of analysis-ready genes in a given pathway, divided by the total number of genes in the reference dataset that makes up that pathway. Z-score predicts the activation state of the canonical pathways, using the gene expression patterns of the genes within the pathway. A negative z value connotes an overall pathway’s inhibition while a positive z value is representative of an overall pathway’s activation. A NaN z-score is assigned to a pathway when gene expression patterns within that pathway are not sufficient to perform the analysis and predict the activation/inhibition state.

Additional file 4 Supplementary Table 4.

List of 21 relevant functional networks identified by performing a network analysis in IPA on 335 differentially expressed mRNAs.

Additional file 5 Supplementary Table 5.

List of 77 differentially expressed mature miRNAs identified in Li responders (R) versus non-responders (NR) with |FC| ≥ 1.2 and p-value ≤0.05.

Additional file 6 Supplementary Table 6.

List of 77 statistically significant pathways targeted by differentially modulated miRNAs identified using Diana Tool software (miRPath v.3). MiRNAs and their targeted mRNAs are also reported.

Additional file 7 Supplementary Table 7.

List of 5450 significant mature miRNA-mRNA interactions with p-value ≤0.05 and Pearson correlation index ρ ≤ − 0.38.

Additional file 8 Supplementary Table 8.

List of 2027 significant mature miRNA-mRNA interactions with p-value ≤0.05 and Pearson correlation index ρ ≤ − 0.5.

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

Cattane, N., Courtin, C., Mombelli, E. et al. Transcriptomics and miRNomics data integration in lymphoblastoid cells highlights the key role of immune-related functions in lithium treatment response in Bipolar disorder. BMC Psychiatry 22, 665 (2022). https://doi.org/10.1186/s12888-022-04286-3

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12888-022-04286-3

Keywords