TSPAN5, ERICH3 and selective serotonin reuptake inhibitors in major depressive disorder: pharmacometabolomics-informed pharmacogenomics.

Millions of patients suffer from major depressive disorder (MDD), but many do not respond to selective serotonin reuptake inhibitor (SSRI) therapy. We used a pharmacometabolomics-informed pharmacogenomics research strategy to identify genes associated with metabolites that were related to SSRI response. Specifically, 306 MDD patients were treated with citalopram or escitalopram and blood was drawn at baseline, 4 and 8 weeks for blood drug levels, genome-wide single nucleotide polymorphism (SNP) genotyping and metabolomic analyses. SSRI treatment decreased plasma serotonin concentrations (P<0.0001). Baseline and plasma serotonin concentration changes were associated with clinical outcomes (P<0.05). Therefore, baseline and serotonin concentration changes were used as phenotypes for genome-wide association studies (GWAS). GWAS for baseline plasma serotonin concentrations revealed a genome-wide significant (P=7.84E-09) SNP cluster on chromosome four 5' of TSPAN5 and a cluster across ERICH3 on chromosome one (P=9.28E-08) that were also observed during GWAS for change in serotonin at 4 (P=5.6E-08 and P=7.54E-07, respectively) and 8 weeks (P=1.25E-06 and P=3.99E-07, respectively). The SNPs on chromosome four were expression quantitative trait loci for TSPAN5. Knockdown (KD) and overexpression (OE) of TSPAN5 in a neuroblastoma cell line significantly altered the expression of serotonin pathway genes (TPH1, TPH2, DDC and MAOA). Chromosome one SNPs included two ERICH3 nonsynonymous SNPs that resulted in accelerated proteasome-mediated degradation. In addition, ERICH3 and TSPAN5 KD and OE altered media serotonin concentrations. Application of a pharmacometabolomics-informed pharmacogenomic research strategy, followed by functional validation, indicated that TSPAN5 and ERICH3 are associated with plasma serotonin concentrations and may have a role in SSRI treatment outcomes.

Introduction

Major depressive disorder (MDD) is the most common psychiatric disorder worldwide, with a lifetime prevalence of approximately 13%.[1, 2] MDD is associated with marked morbidity and premature mortality.[3] Although the causes of MDD are not fully understood, relative deficiency of the neurotransmitter serotonin appears to have a role in the pathophysiology of MDD and, as a result, drugs that enhance serotonergic neurotransmission are used to treat MDD.[4, 5] Selective serotonin reuptake inhibitors (SSRIs), drugs that increase serotonin signaling in the central nervous system by blocking its presynaptic reuptake, are first-line pharmacologic therapy for MDD.[6, 7, 8] However, response to SSRIs is highly variable, with less than half of MDD patients achieving remission during therapy with these drugs.[9, 10, 11]

Twin and other genetic studies suggest that inheritance contributes both to MDD risk[12, 13, 14, 15] and to variation in SSRI response.[16, 17, 18] We[10, 11] and others[19, 20, 21, 22, 23, 24] have performed antidepressant response candidate gene and genome-wide association studies (GWAS), but with only limited success and with few replicated findings.[17, 25, 26, 27] Relative lack of power, variation in study design and phenotypic heterogeneity may all contribute to this state of affairs. The addition of other ‘omics' to genomics might make it possible to achieve enhanced patient subclassification, thus making it possible to identify novel genetic factors that contribute to variation in SSRI response.

We have previously used pharmacometabolomics to help guide and inform genomic studies of SSRI clinical response.[28, 29] Metabolomics is being used increasingly to identify ‘biosignatures' for disease subclassification and/or drug response phenotype(s).[30, 31, 32] Pharmacometabolomics is an emerging field that uses ‘metabolic profiles' to characterize biological response to drug treatment.[28, 29, 33, 34, 35] In the present study, 306 MDD patients were randomly selected from the Mayo Clinic Pharmacogenomics Research Network Antidepressant Medication Pharmacogenomics Study (PGRN-AMPS) SSRI trial who were included in our ‘Clinical SSRI Response' and ‘Citalopram and Escitalopram Metabolism' GWA studies.[11, 36, 37] Plasma samples from those patients were used to perform metabolomic studies through the Pharmacometabolomics Research Network at baseline and after 4 and 8 weeks of SSRI therapy, for a total of 918 samples assayed. Among the metabolites analyzed, plasma serotonin concentrations and changes in plasma serotonin concentrations were associated with the largest number of SSRI treatment outcome measures. Specifically, patients with higher baseline plasma serotonin concentrations and/or greater decreases in plasma serotonin concentrations responded better to SSRI therapy. We then moved from metabolomics to genomics by performing GWAS to identify genes associated with variation in plasma serotonin concentrations or changes in serotonin concentrations during SSRI therapy, followed by the functional pursuit of those genes in neuronal cell models.

Specifically, when GWAS was performed with baseline plasma serotonin concentrations as the phenotype, a genome-wide significant (P=7.84E-09) single nucleotide polymorphism (SNP) signal that was 5' of the Tetraspanin 5 (TSPAN5) gene on chromosome four and a cluster of SNPs across the Glutamate-rich 3 (ERICH3) gene on chromosome one (P=9.28E-08) were identified. Those same SNP signals were identified during GWAS for change in plasma serotonin concentrations after 4 and 8 weeks of SSRI therapy. In addition, the Genome Tissue Expression (GTEx) Database[38] showed that both of those genes were highly expressed in the brain. The SNPs 5' of TSPAN5 were cis-expression quantitative trait loci (eQTLs) for that gene. Follow-up functional genomic experiments performed by knocking down or overexpressing TSPAN5 in a neuroblastoma cell line resulted in significant alterations in the expression of genes encoding serotonin pathway enzymes as well as changes in the concentration of serotonin in the cell culture media. Two of the SNPs in the ERICH3 SNP cluster encoded nonsynonymous variants (ns) that were associated with accelerated proteasome-mediated degradation of ERICH3. In addition, changes in ERICH3 expression significantly altered media serotonin concentrations but did not influence serotonin pathway gene expression. Finally, one of the ERICH3 nsSNPs (rs11580409, P=1.12E-07) was associated with clinical SSRI response in the International SSRI Pharmacogenomics Consortium (ISPC), an observation that was replicated in the Sequenced Treatment Alternatives to Relieve Depression (STAR*D) study. In summary, the application of a ‘pharmacometabolomics-informed pharmacogenomic' research strategy made it possible to identify two novel genes related to plasma serotonin concentration—a phenotype that was associated with SSRI clinical response.

Materials and methods

Trial design, samples and metabolomic assays

Patient selection, treatment outcomes and blood sample collection for the Pharmacogenomics Research Network Antidepressant Medication Pharmacogenomics Study (PGRN-AMPS) SSRI trial have been described in detail elsewhere.[11, 36, 37] Plasma metabolite concentrations were assayed using samples from 306 randomly selected MDD patients at baseline and after 4 and 8 weeks of SSRI therapy using a high-performance liquid chromatography electrochemical coulometric array metabolomics platform.[31, 39] See Supplementary Text for details.

Genotyping and statistical analyses

DNA from PGRN-AMPS SSRI trial patients was genotyped at the RIKEN Center for Genomic Medicine (Yokohama, Japan) using Illumina human 610-Quad BeadChips (Illumina, San Diego, CA, USA), as described previously.[11, 37] GWAS were performed using approximately 7.5 million SNPs. Patients were removed from the analysis for non-compliance or non-Caucasian heritage. Baseline analyses were adjusted for age and sex. Metabolite concentrations and changes in metabolite concentrations after SSRI treatment were tested for association with QIDS-C16 percent change, response and remission. See Supplementary Text for details.

TSPAN5, ERICH3 and SNP function

Lymphoblastoid cell lines (LCLs) were selected from the 'Human Variation Panel' based on TSPAN5 or ERICH3 SNP genotypes to determine whether the SNPs were eQTLs for those genes. The 300 LCLs (100 European-American, 100 African-American and 100 Han Chinese-American subjects) in the ‘Human Variation Panel' that had been SNP genotyped previously have been utilized repeatedly to generate and test pharmacogenomic hypotheses.[40, 41, 42, 43, 44] TSPAN5 SNP function was assessed using electrophoretic mobility shift assays and dual luciferase reporter gene assays. Expression constructs for ERICH3 that encoded wild type (WT) as well as one or both nsSNPs (rs11580409 or rs11210490) were expressed with or without the proteasome inhibitor MG132 or the autophagy inhibitor 3-methyladenosine. See Supplementary Text for details.

TSPAN5 and ERICH3 expression and the serotonin pathway

After TSPAN5 or ERICH3 knockdown (KD) or overexpression (OE) in neurally derived cell lines, serotonin pathway enzyme expression was assessed by quantitative real-time polymerase chain reaction (qRT-PCR) and quantitative western blot. Cell culture media serotonin concentrations were measured by Bioanalytical Systems (BASi, West Lafayette, IN, USA). See Supplementary Text for details.

Results

Plasma metabolite concentrations and their association with clinical outcomes

We set out to use plasma metabolomic profiles of MDD patients being treated with SSRIs to identify metabolites that were correlated with SSRI clinical outcomes and, subsequently, SNPs/genes associated with those metabolite concentrations for functional study in neuronal cell lines. This approach made it possible to move from peripheral plasma metabolomics to genomics and then to test genomic candidates in neural cells—addressing concerns with regard to the relevance of peripheral biomarkers for neuronal function.

Specifically, a liquid chromatography electrochemical coulometric array metabolomics platform was used to quantify 31 known plasma metabolites (Supplementary Table 1), primarily metabolites in the tryptophan, tyrosine, purine and tocopherol pathways, at three time points—baseline and after 4 and 8 weeks of SSRI therapy. We then determined the association of those metabolites with measures of clinical response (remission, response and percent change in QIDS-C16) after 4 and 8 weeks of SSRI therapy. Plasma serotonin concentrations at baseline as well as their change after 4 and 8 weeks of SSRI treatment were more highly associated with SSRI response phenotypes than those for any other metabolite (Table 1 and Supplementary Table 2). The associations listed in Table 1 are ‘nominal' and have not been corrected for multiple comparisons because the purpose was to identify metabolites to use for GWAS. Plasma serotonin concentrations decreased significantly after SSRI treatment at both 4 (P<0.0001) and 8 weeks (P<0.0001) (Figure 1). The odds ratios (OR) and correlation coefficients (r) listed in Table 1 indicated that higher baseline plasma serotonin concentrations as well as larger decreases in plasma serotonin concentrations between baseline and 4 or 8 weeks of therapy were both associated with better clinical outcomes. We then performed GWAS using baseline plasma serotonin concentrations and change in plasma serotonin concentrations at 4 and 8 weeks of SSRI therapy as phenotypes.

Table 1

Association of plasma serotonin concentrations with clinical outcomes

Clinical outcomes	Remission at 4 weeks	Remission at 8 weeks	Response at 4 weeks	Response at 8 weeks	% change at 4 weeks	% change at 8 weeks
Baseline	P=0.012	P=0.028	P=0.007	P=0.047	P=0.015	P=0.019
	OR=1.41	OR=1.31	OR=1.40	OR=1.30	r=−0.14	r=−0.14
Change after 4 weeks	P=0.011	P=0.041	P=0.026	P=0.060	P=0.021	P=0.024
	OR=1.40	OR=1.27	OR=1.31	OR=1.27	r=−0.13	r=−0.13
Change after 8 weeks	P=0.069	P=0.147	P=0.037	P=0.130	P=0.041	P=0.06
	OR=1.27	OR=1.19	OR=1.29	OR=1.21	r=−0.12	r=−0.11

Abbreviation: OR, odds ratio; SSRI, selective serotonin reuptake inhibitor. Plasma serotonin concentrations at baseline and decreases in plasma serotonin concentrations after 4 weeks of SSRI treatment were nominally associated with remission, response and percent change in QIDS-C16 score. The decrease in plasma serotonin between baseline and 8 weeks of SSRI treatment was associated only with the response at 4 weeks and percent change in QIDS-C16 at 4 weeks. OR>1 indicates improvement (associated with higher baseline and larger changes in plasma serotonin concentrations) and the negative r values indicate a decrease in QIDS-C16 scores, i.e., improvement. P-values<0.05 have been bolded.

Figure 1

Patient plasma serotonin concentrations. Relative plasma serotonin concentrations (expressed as a ratio of the standard) in major depressive disorder patient samples were decreased significantly after 4 and 8 weeks of selective serotonin reuptake inhibitor treatment when compared with baseline. ***P<0.0001.

GWAS for plasma serotonin and change in serotonin concentrations

The Manhattan plot of the GWAS for baseline plasma serotonin concentrations showed a genome-wide significant SNP cluster on chromosome four that consisted of 15 SNPs in tight linkage disequilibrium that mapped 15–25 kilobases (kb) 5' of the Tetraspanin 5 (TSPAN5) gene, with the lowest P-value (7.84E-09) for the rs11947402 SNP (Figures 2a and c and Supplementary Table 3). In addition, there was a SNP cluster on chromosome one across the Glutamate-rich 3 (ERICH3) gene that included two nsSNPs (rs11580409 and rs11210490), with the lowest P-value (9.28E-08) for the rs696692 SNP (Figures 2a and b). These same two signals were observed in Manhattan plots of GWAS data for change in plasma serotonin concentrations after SSRI therapy for 4 (TSPAN5: rs11947402, P=5.6E-08; ERICH3: rs696692, P=7.54E-07) and 8 weeks (TSPAN5: rs11947402, P=1.25E-06; ERICH3: rs699848, P=3.99E-07) (Supplementary Figure 1). QQ plots for these GWAS are shown in Supplementary Figure 2.

Figure 2

Baseline serotonin concentration GWAS. (a) GWAS for baseline plasma serotonin concentrations revealed a genome-wide significant signal on chromosome 4 as well as a suggestive SNP cluster on chromosome 1. (b) The locus zoom shows that the SNPs on chromosome 1 are across ERICH3. The SNP most highly associated with baseline plasma serotonin concentration from this cluster was rs696692 (P=9.28E-08). (c) The locus zoom for the genome-wide significant SNP cluster on chromosome 4 shows that the SNPs are approximately 15–25 kb 5' of TSPAN5, with rs11947402 as the most highly associated with baseline plasma serotonin concentration (P=7.84E-09). GWAS, genome-wide association study; SNP, single nucleotide polymorphism.

Variant TSPAN5 SNP genotypes were associated with higher baseline plasma serotonin concentrations and greater decreases in plasma serotonin concentration during SSRI therapy (Supplementary Figure 3a-c). Conversely, ERICH3 variant allele genotypes were associated with lower baseline plasma serotonin concentrations and smaller decreases in plasma serotonin concentrations during SSRI therapy (Supplementary Figure 3d-f).

The minor allele frequency for the SNPs 5' of TSPAN5 was 7% in our European-American MDD patients, consistent with the 6.7% value reported for European populations by the 1000 Genomes Project.[45] The relatively low minor allele frequency value complicated efforts to examine the effect of homozygosity for the variant allele in the functional genomic studies described subsequently. The SNPs across the ERICH3 gene had a minor allele frequency of 35% in the PGRN-AMPS patients, similar to the 34.1% figure reported for European populations by the 1000 Genomes Project.[45] As the SNPs near TSPAN5 were genome-wide significant and the SNPs across ERICH3 were highly suggestive, we pursued the possible functional implications of both signals.

TSPAN5 and ERICH3 SNPs as eQTLs

The initial question that we asked with regard to the SNPs 5' of TSPAN5 and across ERICH3 was whether they might be cis-eQTLs for those genes. Specifically, we selected LCLs from the Human Variation Panel that were either homozygous WT or variant (V) for the SNPs 5' of TSPAN5 or those across ERICH3 and performed qRT-PCR analysis. TSPAN5 mRNA was decreased in LCLs homozygous for the variant genotype as compared with LCLs homozygous for the WT genotype (P<0.05; Figure 3a). However, there was no difference in ERICH3 expression between LCLs homozygous WT and homozygous variant for the ERICH3 SNPs (data not shown).

Figure 3

TSPAN5 and ERICH3 SNP function. TSPAN5 expression is decreased for cells or tissues homozygous for variant (V/V) or with heterozygous (WT/V) SNP genotypes as compared with homozygous wild type (WT/WT) in (a) LCLs, (b) cerebral cortex and (c) frontal cortex. (d) Luciferase assay results comparing WT and variant SNP genotypes (rs1918743, rs59961429 and rs56095565) effects on transcriptional activities indicate decreased transcription for the variant TSPAN5 SNPs in SK-N-BE(2) neuroblastoma cells; (e) ERICH3 plasmids that were WT or contained one or both of the nonsynonymous SNPs (rs11580409 and rs11210490) were expressed in HEK-293T/17 cells. Both P264A (rs11210490) and L1056V (rs11580409) were associated with decreased protein levels as compared with WT, but L1056V was associated with a much greater decrease in protein level. (f) Quantification of ERICH3 protein relative to the GAPDH control for the ERICH3 western blots shown in (e). (g) Plasmids encoding ERICH3 allozymes that were WT or contained one or both of the amino acid substitutions (P264A and L1056V) were expressed in HEK-293T/17 cells with and without a protease inhibitor (MG132) or an autophagy inhibitor (3MA). MG132 prevented ERICH3 SNP-dependent protein degradation but 3MA did not. (h) Quantification of proteasome and autophagy inhibition of the ERICH3 allozyme degradation studies shown in (g). EV, empty vector; NS, non-significant; *P<0.05; **P<0.01; ***P<0.0001. LCL, lymphoblastoid cell lines; SNP, single nucleotide polymorphism.

We next consulted eQTL databases to determine whether the SNPs 5' of TSPAN5 or across ERICH3 might be cis-eQTLs. The GTEx Database[38] showed that brain displayed high expression for both TSPAN5 and ERICH3 (Supplementary Figures 4 and 5). Previous studies in mice had reported that the tissue with the highest TSPAN5 protein level was brain.[46, 47, 48, 49] Unfortunately, GTEx had too few samples to provide reliable eQTL data for TSPAN5. However, the Brain eQTL Almanac (BRAINEAC) database included data for 134 human brain samples.[50] We found that, in the brain areas with the highest TSPAN5 expression (cerebral and frontal cortex), the SNPs were once again cis-eQTLs for TSPAN5 (P=6.9E-05 and 0.027, respectively), with lower expression for homozygous variant or heterozygous genotypes—just as we found for LCLs (Figures 3b and c). Finally, we used the Blood eQTL browser that contains data for 5311 individual blood samples.[51] Once again, the SNPs were identified as eQTLs for TSPAN5 (P=1.36E-14) but, in this case, the Z-score of 7.8 indicated that the variant allele was associated with higher TSPAN5 expression—opposite to what we found using LCLs or BRAINEAC. Therefore, these SNPs are cis-eQTLs for TSPAN5, but there appears to be tissue-specific directionality of associations with the SNPs—an important factor for the interpretation of subsequent functional studies. The SNPs across ERICH3 were not cis-eQTLs in any of these databases.

We next attempted to determine which of the SNPs 5' of TSPAN5 might influence expression. Specifically, we used the TRANSFAC 6.0 database to identify transcription factors that might bind to DNA sequences that contained the SNPs. Eight of the SNPs were predicted to potentially disrupt or create transcription factor-binding sites. Electrophoretic mobility shift assays performed using nuclear protein extracts from neuroblastoma SK-N-BE(2) cells revealed differences in nuclear protein-binding patterns between WT and variant SNP sequences for three of the eight SNPs (rs1918743, rs59961429 and rs56095565) (Supplementary Figure 6). In an attempt to more directly determine the possible role of these three SNPs in transcription, we performed luciferase reporter assays by transfecting luciferase reporter gene constructs containing each of the SNPs into SK-N-BE(2) cells. Each of the variant SNP genotypes significantly decreased luciferase activity when compared with the WT genotype, indicating decreased transcriptional activity in SK-N-BE(2) neuroblastoma cells (Figure 3d)—compatible with the results of our eQTL analyses for brain and LCLs.

ERICH3 SNPs and proteasome-mediated degradation

As the SNPs across ERICH3 were not cis-eQTLs for that gene, we tested the possibility that the two nsSNPs might affect ERICH3 protein concentrations. Proteasome-mediated degradation is a common functional mechanism for the effect of nsSNPs.[52, 53, 54, 55] ERICH3 cDNA constructs that were WT or contained one or both of the nsSNPs (rs11580409 and rs11210490) were transfected into HEK-293T/17 cells. The rs11210490 SNP (Pro264Ala) was associated with a small (28%), but significant (P<0.05) reduction in ERICH3 protein, while the rs11580409 SNP (Leu1056Val) was associated with an 80% decrease of ERICH protein (P<0.001). Constructs with both nsSNPs were associated with a 93% reduction in ERICH3 protein (P<0.001) (Figures 3e and f). Furthermore, the proteasome inhibitor (MG132) increased ERICH3 variant allozyme concentrations, but the autophagy inhibitor (3-methyladenosine) did not (Figures 3g and h), indicating that degradation of the variant ERICH3 allozymes was proteasome-mediated. We next attempted to identify the potential functional relationship of TSPAN5 and ERICH3 with baseline and change in plasma serotonin concentrations.

TSPAN5 and ERICH3 expression and serotonin pathway enzyme gene expression

SK-N-BE(2) neuroblastoma cells were used to perform TSPAN5 functional genomic studies because they are derived from neural cells and express TSPAN5 and serotonin pathway enzymes (Figure 4a). When TSPAN5 was knocked down more than 70% in SK-N-BE(2) cells, there was a significant decrease of mRNA and protein levels for the serotonin pathway enzymes TPH1, TPH2, DDC, MAOA, (Figure 4a) as well as the serotonin transporter SLC6A4 (Figures 4b and d). Furthermore, OE of TSPAN5 was associated with increased expression of TPH1, TPH2, DDC and MAOA (Figures 4b and d). Protein levels were not induced to the same extent as mRNA after TSPAN5 OE, perhaps because TSPAN5 is highly expressed in those cells.

Figure 4

TSPAN5 and ERICH3 association with serotonin biosynthesis and metabolism. (a) Serotonin biosynthesis and metabolism pathway. (b) mRNA expression of genes encoding serotonin pathway enzymes as measured by qRT-PCR were decreased after TSPAN5 KD (black) and increased after TSPAN5 OE (cross-hatched) in SK-N-BE(2) neuroblastoma cells. (c, d) Western blot analysis indicated decreased serotonin enzyme protein levels after TSPAN5 KD but no significant change after TSPAN5 OE, as quantified in (d). (e, f) Culture media serotonin concentrations after TSPAN5 (e) KD and (f) OE. (g, h) Cell culture media serotonin concentrations after ERICH3 (g) KD and (h) OE. *P<0.05; **P<0.01; TPH1/2, tryptophan hydroxylase 1/2; DDC, dopa decarboxylase; KD, knockdown; MAOA/B, monoamine oxidase A/B; OE, overexpression; SLC6A4, serotonin transporter; EV, empty vector.

Similar ERICH3 KD and OE experiments were performed using neurally derived cells including neuroblastoma cells (SK-N-SH and SK-N-BE(2)), human neural progenitor-derived neurons and glioblastoma cells (U251). However, ERICH3 KD and OE did not alter the expression of serotonin pathway enzymes (data not shown).

RBPJ-κ expression and serotonin pathway enzyme gene expression

TSPAN5 has been reported to be involved in Notch signaling[56] through ADAM10 recruitment[57, 58, 59, 60] and the Notch-Recombination Signal Binding Protein for Immunoglobulin Kappa J Region (RBPJ-κ) has been implicated in regulation of the expression of serotonin pathway genes.[61] Therefore, we knocked down RBPJ-κ in SK-N-BE(2) cells, and observed increased expression of TPH1, TPH2, DDC and SLC6A4 (Supplementary Figure 7). This may be one mechanism by which TSPAN5, an integral membrane protein, may influence serotonin biosynthesis, as described in more detail in the Discussion.

TSPAN5 and ERICH3 and serotonin concentrations in cell culture media

Serotonin concentrations in cell culture media decreased significantly after TSPAN5 KD, but TSPAN5 OE did not result in significant changes in cell culture media serotonin concentrations (Figures 4e and f)—consistent with the changes observed in TPH1, TPH2, DDC, MAOA and SLC6A4 protein levels observed after TSPAN5 KD and OE (Figures 4c and d). ERICH3 KD in human neural progenitor cell-derived neurons and OE in SK-N-BE(2) neuroblastoma cells were associated with significantly altered cell culture media serotonin concentrations (Figures 4g and h), perhaps indicating that ERICH3 influences plasma serotonin concentrations through a mechanism that does not directly involve the expression of serotonin biosynthesis and metabolism enzyme genes. These two cell lines were chosen for study because both express the serotonin biosynthesis and metabolism enzymes and because ERICH3 is highly expressed in human neural progenitor cell-derived neurons, but is not expressed in SK-N-BE(2) neuroblastoma cells.

SNP combinations and plasma serotonin concentrations

We next analyzed the possible association of baseline and change in plasma serotonin concentrations in our patients with combinations of genotypes for the top TSPAN5 SNP (rs11947402) and the ERICH3 rs11580409 nsSNP. As anticipated, patients who were homozygous WT for the TSPAN5 SNP (rs11947402) and who were homozygous variant for the ERICH3 SNP (rs11580409) had lower average baseline (P=1.76E-12) and smaller average decreases in serotonin concentrations after 4 (P=6.09E-11) and 8 weeks (P=1.84E-09) of SSRI treatment as compared with patients carrying the TSPAN5 variant SNP allele and/or the ERICH3 nsSNP WT allele (rs11580409). The R-squared values indicate that the TSPAN5 SNP (rs11947402) and ERICH3 nsSNP (rs11580409) account for 18.8% of the baseline variation of serotonin concentrations and 15.4% and 13% of the variation in the change in serotonin concentrations after 4 and 8 weeks of SSRI treatment in this MDD population, respectively (Supplementary Figure 8). The numbers of patients who had each SNP genotype combination are shown in Supplementary Table 4.

ERICH3 and TSPAN5 SNPs and clinical phenotypes

Finally, we attempted to determine whether these SNPs might be associated with SSRI clinical response in the PGRN-AMPS GWAS[11] and/or in two independent SSRI response GWAS (STAR*D[24] and ISPC[10]). The nsSNP (rs11580409) in ERICH3 that displayed striking proteasome-mediated degradation was associated with response at 4 weeks in the ISPC population (P=0.022, OR=1.25) and response at 6 weeks in the STAR*D population (P=0.041, OR=1.17). The SNPs in the cluster 5' of TSPAN5 were not significantly associated with clinical response in any of these SSRI GWAS (Supplementary Table 5).

Discussion

A major goal of molecular psychiatry is to develop a molecular subclassification of psychiatric disease. In theory, that might allow a rational selection of optimal therapy for each patient—that is, ‘Precision' pharmacotherapy. In the present study, we used the most commonly prescribed antidepressant medication—SSRIs—as probes for molecular mechanisms associated with drug response. Specifically, we applied a pharmacometabolomics-informed pharmacogenomic research strategy during which we utilized plasma samples from 306 MDD patients enrolled in the Mayo PGRN-AMPS SSRI trial[11, 37] to perform metabolomic assays for 31 metabolites, primarily metabolites from pathways related to monoamine neurotransmitters. The goal was to associate individual variation in these plasma metabolites with SSRI treatment outcomes—with the understanding that this represented only one step toward determining whether molecular mechanisms identified in the periphery might also have a role in neurotransmitter function in the brain. For that reason, we used GWAS to identify novel genes that might influence concentrations of the metabolite(s) identified in the periphery and then determined whether those same genes might also influence neuronal cell phenotypes.

We found that plasma serotonin concentrations in MDD patients decreased dramatically after SSRI treatment (Figure 1). The effect of SSRIs on plasma serotonin concentrations is not well understood, and published results are contradictory—with some studies reporting decreased plasma serotonin after SSRI treatment,[32, 62, 63, 64, 65] while others report increased concentrations,[66, 67, 68, 69] which may be due, in part, to the plasma collection method or the platform used to assay serotonin. However, we observed a clear decrease in those concentrations in PGRN-AMPS MDD patients as measured with a liquid chromatography electrochemical coulometric array platform—a highly sensitive, quantitative method.[31, 39] We also observed that both higher baseline plasma serotonin concentrations and greater decreases in plasma serotonin concentrations after 4 and 8 weeks of SSRI therapy were associated with better clinical outcomes (response, remission and QIDS-C16 percent change) than were observed for any of the other metabolites assayed (Table 1 and Supplementary Table 1). In an attempt to identify genes associated with individual variation in plasma serotonin concentrations and/or changes in those concentrations during SSRI therapy, we performed a GWAS for both phenotypes.

Unlike previous SSRI response GWAS in which no genome-wide significant SNPs were identified, our GWAS for baseline plasma serotonin concentrations included a genome-wide significant (7.84E-09) SNP cluster 5' of the TSPAN5 gene on chromosome four (Figures 2a and c). We also observed a SNP cluster across the ERICH3 gene on chromosome one (Figures 2a and b). The same SNP signals were observed when GWA studies were performed for change in plasma serotonin concentrations after 4 or 8 weeks of SSRI therapy (Supplementary Figure 1). These observations raised the possibility that TSPAN5 and/or ERICH3 might be involved in the regulation of genes encoding enzymes in the serotonin metabolic pathway—a hypothesis that we tested in neurally derived cell lines.

In an attempt to understand the possible role of the SNPs 5' of TSPAN5 in the regulation of serotonin biosynthesis, metabolism or transport, we first determined that those SNPs were cis-eQTLs for TSPAN5 in LCLs, brain tissue and blood samples (Figures 3d and g and Supplementary Figure 6). However, the variant allele was associated with lower TSPAN5 expression in LCLs and brain tissue but higher expression in blood samples. We then showed that TSPAN5 KD and OE in neuroblastoma cells were associated with changes in the expression of serotonin pathway genes. In addition, TSPAN5 KD in neuroblastoma cells was associated with a significant decrease in serotonin concentration in the cell culture media (Figures 4a and f). These results indicated that the SNPs 5' of TSPAN5 could influence its expression, which, in turn appeared to have a role in the regulation of serotonin-related pathways.

TSPAN5 is a member of the tetraspanin superfamily, a family of proteins characterized by four hydrophobic transmembrane domains.[49, 70] Tetraspanins form molecular complexes within the plasma membrane that can modulate cellular signaling.[71, 72] TSPAN5 function has not been investigated extensively and has not previously been implicated in the regulation of serotonin or variation in SSRI response. However, several recent studies reported that TSPAN5 may promote Notch signaling[56] by facilitating the transport of ADAM10, an α-secretase involved in cleaving the Notch receptor, to the cell membrane.[57, 58, 59, 60] The Notch intracellular domain is then transported to the nucleus where it binds to transcription factors on gene promoters, inducing changes in gene expression.[73, 74] A recent study reported that RBPJ-κ may have a role in the expression of DDC and MAOA,[61] and we showed that RBPJ-κ KD resulted in increased expression of TPH1/2, DDC and SLC6A4 in SK-N-BE(2) neuroblastoma cells (Supplementary Figure 7). Obviously, future studies will be required to clarify the possible functional relationships among ADAM10, TSPAN5, Notch and RBPJ-κ.

We also pursued the function of the chromosome one SNP signal across ERICH3 that included two nsSNPs. These same SNPs were also associated with response in both the ISPC (P=0.022, OR=1.3) and STAR*D (P=0.041, OR=1.2) studies. The SNPs across ERICH3 were not cis-eQTLs, but rather, the nsSNPs were associated with the proteasome-mediated degradation of ERICH3 protein. KD and OE of ERICH3 did not alter the expression of genes encoding serotonin pathway enzymes but were associated with significant changes of serotonin concentrations in the cell culture media. The functional mechanism by which ERICH3 influences serotonin concentrations is unclear.

The present study has shown that baseline plasma serotonin concentrations and decreases in plasma serotonin concentrations after 4 and 8 weeks of SSRI therapy were associated with clinical outcomes in our MDD patients (Table 1). Furthermore, TSPAN5 variant and WT ERICH3 SNP genotypes (Figure 2 and Supplementary Figure 1) were associated with higher baseline plasma serotonin concentrations and larger decreases in plasma serotonin concentrations after SSRI therapy (Supplementary Figure 3). The SNPs across ERICH3 included a nsSNP that resulted in ERICH3 proteasome-mediated degradation. The SNPs 5' of TSPAN5 were eQTLs for that gene in LCLs (Figure 3a), brain (Figures 3b and c) and blood—although with tissue-specific differences in directionality. Higher TSPAN5 expression was associated with the variant SNP genotype in blood, which, based on our functional studies, would suggest higher expression of serotonin pathway genes and elevated serotonin synthesis (Figure 4)—consistent with the higher baseline plasma serotonin concentrations observed in our clinical data—although that hypothesis will require future validation. Finally, the ERICH3 nsSNP was associated with clinical response in two independent SSRI studies.

In summary, the present series of experiments have demonstrated that metabolomics can be a useful tool to help identify novel biology—especially when it is used to guide and inform subsequent genomic studies. By integrating pharmacometabolomic and pharmacogenomic data related to SSRI treatment response, we identified SNPs that are cis-eQTLs for TSPAN5—a gene not previously known to be involved in either SSRI response or the regulation of serotonin-related pathways and nsSNPs in ERICH3 that altered the quantity of ERICH3 protein. Virtually nothing was known about ERICH3 prior to the studies reported here. Future experiments will be required to pursue these observations in depressed and non-depressed populations as will additional pharmacometabolomics-informed pharmacogenomic studies to help us move toward the goal of enhanced molecular subclassification of psychiatric disease and its response to drug therapy.