Association of Novel ALX4 Gene Polymorphisms with Antidepressant Treatment Response: Findings from the CO-MED Trial

Major depressive disorder (MDD) is a serious and highly prevalent psychiatric illness that requires early detection, diagnosis, treatment, and comprehensive understanding of treatment outcome [1]. Although pharmacological treatments demonstrate efficacy, only about 50% of individuals respond to their first-administered antidepressant medication, and even fewer experience complete remission (i.e., minimal or no symptoms) [2, 3]. The landmark STAR*D (Sequential Treatment Alternatives to Relieve Depression) study [4] examined the effectiveness of different pharmacological treatments in adults with MDD who did not respond to initial antidepressant treatment. Following up to 14 weeks of acute treatment with citalopram in 4,041 subjects, about half (47%) achieved response (i.e., 50% improvement), although only about one-third achieved complete remission [5]. Even after receiving second-line treatment strategies, including both switch and augmentation strategies, only an additional 25–35% achieved remission of symptoms [6, 7]. Rush et al.[8] subsequently reported on the predictors of remission in the STAR*D sample, noting that remission was more likely among individuals who were Caucasian, employed, married or cohabitating, privately insured, had no prior suicide attempts, and had prior intolerance or at least some response to citalopram during the first phase of the study. Many other studies have also reported on predictors of response to antidepressant treatment [9, 10], and while some nonspecific predictors have been identified, there remain few validated clinical or biological markers which may guide clinicians toward personalized, effective first-line treatment [11]. Thus, a trial and error method of antidepressant selection remains necessary.

With the advancement of genetic and molecular technologies, precision medicine for depression is becoming a distinct possibility, however [12]. Identifying common pharmacogenetic variants that correspond with treatment outcome would enable personalized psychiatric treatment and significantly improve treatment efficacy. Currently, there is not a robust genomic biomarker capable of translation into clinical practice, although the search for novel biomarkers in large cohort studies is ongoing. To address this issue, the Combining Medications to Enhance Depression Outcomes (CO-MED) trial [13] was carried out to test different treatment methods (escitalopram [ESCIT] monotherapy or one of two combination therapies: ESCIT/bupropion [BUP] or venlafaxine [VEN]/mirtaz-apine [MIRT]) and assess subsequent treatment outcomes in over 600 participants. While neither medication combination outperformed monotherapy [13], genetic factors plausibly contribute to the differential drug response. Thus, identification of individual genetic profiles might provide guidance for effective medication selection.

Pharmacogenetic studies focus on candidate genes implicated in the mechanisms of action or the pharmacokinetics of antidepressant drugs [14, 15]. Although there are several million single nucleotide polymorphisms (SNPs) in the human genome, identifying a key genetic variation is difficult. The GWAS (genome-wide association studies) approach provides researchers with the potential to identify key genetic variations that may ultimately improve our ability to predict antidepressant treatment response [16]. To date, the most promising genetic variants identified in the literature are insertions/deletions in the promoter region of the serotonin transporter (SLC6A4) gene and an alteration in the ABCB1 gene coding for a P-glycoprotein [17, 18]. Both genes seem to predict treatment response to selective serotonin reuptake inhibitors (SSRIs). In addition, reports suggest that variants of (1) a gene coding for FKBP5 [19-21], a co-chaperone involved in stress hormone signaling, (2) cytochrome P450 (CYP) enzymes (CYP2D6, CYP2C19, CYP2C9, CYP2B6) [22], and (3) 5HT2A [23] and 5-HTTLPR [24] are important, the latter two of which are currently the main candidate genetic variants for prediction of treatment response. However, all of the genetic variants known so far do not effectively guide treatment selection. Further associations based on antidepressant treatment response have been reported for monoamine metabolic enzymes, monoamine oxidase, catechol-O-methyl transferase, tryptophan hydroxylase (MAO-A, COMT, TPH), the glutamatergic receptor gene GRIK4 [25], the enzymatic gene PDER11A [26], inflammation-related genes (CD3E, PRKCH, PSMD9, and STAT3) [27], neurotrophic factor BDNF with the polymorphism of 66 Val/Met variant and for UCN318 expressing a ligand of the CRF2 receptor [27].

The mechanisms by which antidepressants alleviate symptoms of depression are yet not fully understood, and GWAS provide an agnostic approach for identifying genetic variants that may explain the differential responsiveness among patients. GWAS on STAR*D (n = 1,953) [25], MARS (n = 339) [28], GENDEP (n = 706) [9], PGRN-AMPS (n = 529) [29], and ISPC (n = 865) [30] have shown moderate effects in identifying significant SNPs in association with antidepressant response [31]. In the present study, we investigated potential genetic polymorphisms or genetic variants which may predict or moderate treatment outcomes for antidepressant medications. We examined data from the CO-MED study, which compared differential response and remission rates following treatment with a single antidepressant (ESCIT+placebo [PBO]) or one of two combination medication options (BUP-SR plus ESCIT, BUP+ESCIT, or VEN+MIRT) [13]. In the CO-MED study, at week 12, the response and remission rates with most secondary outcomes were not different among treatment groups. The remission rates were 38.8% for ESCIT+PBO, 38.9% for BUP+ESCIT, and 37.7% for VEN+MIRT, and the response rates were even more similar, ranging from 51.6 to 52.4%.

In the present study, we utilized a robust and high density SNP genotyping assay to screen 2.5 million common genetic variants in the human genome to identify genetic variations associated with antidepressant response. Following methodology from STAR*D and other genetic studies, we investigated CO-MED samples to identify pharmacogenetic associations at two time points: an early response at week 6 (QIDS-SR ≤5) and remission at week 12 (QIDS-SR ≤6 and ≤8 at the last two study visits) comparing monotherapy versus combined treatment [13].

CO-MED was a single-blind, randomized, PBO-controlled trial for first-step MDD treatment and includes an acute (12 weeks) and long-term continuation (7 months) component. CO-MED recruited 665 participants who were stratified by clinical sites and assigned to one of the three treatment groups in a 1: 1: 1 ratio. Among them, 459 provided biospecimens at baseline for DNA extraction and hence constitute the analytic sample of this report. The treatment groups in this study were: (1) ESCIT+PBO, (2) sustained-release BUP-SR plus ESCIT, and (3) extended-release VEN plus MIRT. Study visits were conducted at baseline and weeks 1, 2, 4, 6, 8, 10, and 12 for acute phase and weeks 16, 20, 24, and 28 for continuation phase. Participation in continuation-phase treatment (after week 12) was restricted to participants who had either received an acceptable benefit (defined as score of 9 or less on the 16-item Quick Inventory of Depressive Symptomatology – Clinician rated [QIDS-C] by week 12) or had reached a score of 10–13 on QIDS-C in addition to both the physician and participant deciding to continue treatment because of substantial benefit. The methods of CO-MED have been described in greater detail by Rush et al. [13].

For the CO-MED study, broad inclusion and minimal exclusion criteria were used to ensure a representative participant group. Patients seeking treatment at participating clinical sites and planning to continue living in the area of that clinical site for the duration of the study were eligible to enroll. These clinical sites included both primary and psychiatric care clinics and were selected to ensure adequate minority representation. To be included in the study, participants had to meet clinical criteria for nonpsychotic MDD, recurrent (greater than 1 previous episode) or chronic (current episode greater than 2 years) as defined by a clinical interview and confirmed by the MINI International Neuropsychiatric Interview (MINI), had at least a 2-month duration of the current episode, and scored ≥16 on the 17-item Hamilton Rating Scale. Only participants who provided samples for DNA investigation (n = 459) were analyzed. Exclusion criteria for this study are fully listed on the clinicaltrials.gov website [32].

The institutional review boards at UT Southwestern Medical Center at Dallas, the University of Pittsburgh Data Coordinating Center, each participating regional center, and all relevant clinical sites reviewed and approved the study protocol and all consent documents and study procedures. An independent data safety and monitoring board also monitored the study. This trial is registered with Clinicaltrials.gov (Clinicaltrials.gov identifier NCT00590863). The details on antidepressant medications, clinical assessments and research outcomes are discussed in the supplementary methods section.

Study samples were genotyped with Illumina Quad, Human Omni 2.5 bead chip, which contains ∼2.5 million markers. Genotyping was performed using the Illumina protocol (Human Omni 2.5 bead chip assay). About 250 ng of nanodrop-measured genomic DNA was used for genotyping. The markers on the chip were mapped using NCBI build 36 (hg18) (Illumina manifest file Immuno_BeadChip_11419691_B.bpm). Raw image files were imported, and SNPs were identified using a genotyping module of Illumina Genome Studio (GS version 1.9.4). Then, raw genotype outputs from Genome Studio were imported into SNP and Variation suite (SVS version 7.6.8 win64) of Golden Helix for quality control (QC) and downstream analysis. QC was conducted to exclude SNPs and samples with call rate less than 90%, admixture, cryptic relatedness, abnormal heterozygosity, as well as degraded or otherwise unsuccessful genetic samples. Similarly, monomorphic markers and duplicate samples were removed. QC assessments included overall genotype concordance rates based on duplicate sample genotyping and Mendelian inheritance checks based on genotyping of a CEPH trio of two parents and their child. For each SNP, the minor allele frequency (MAF), call rate, and departure from Hardy-Weinberg equilibrium were evaluated.

Principal component analysis (PCA) was done using a total of 571,707 SNPs with MAF ≥0.05 in study population. Genotype-phenotype association analyses were restricted to White Non-Hispanic subjects because Caucasian populations are well characterized and have lot of annotated data resources available in the literature as compared to other ethnic groups. The potential ethnic outliers identified by PCA were analyzed separately. Samples that clustered with CEU HapMap population were considered for downstream analysis. QC filters were applied to both samples and markers. All of the duplicates and those with call rates < 85% were excluded from the analysis. Principal component and quality control analyses were done on 256,890 genotyped SNPs. The genotypes had high concordance rates in duplicate samples (> 99.99%) and low rates of Mendelian inheritance errors (< 0.01%). 30,700 SNPs failed genotyping. Another 28,610 SNPs had allele frequencies below the preset threshold of 0.05 and, consequently, were excluded from analysis due to low power. An additional 49 SNPs were excluded due to significant departures from Hardy-Weinberg Equilibrium (p < 10−⁶). The remaining SNPs (n = 232,877) had over 90% call rates. Genome-wide heterozygosity was within the expected range (online suppl. Fig. 1A; for all online suppl. material, see www.karger.com/doi/10.1159/000487321).

Quality pass samples and markers in each group were used for genetic association tests. A basic allelic association test was performed on common frequency markers (MAF ≥0.05). The association test was controlled for genomic inflation. χ²p values, allele frequencies, and odds ratio were calculated for minor and major alleles using SVS, Golden Helix software. The three CO-MED treatment groups were combined, and genetic association tests were performed on the genotyped data. Outcomes included both response and remission at weeks 6 and 12.

The linkage disequilibrium (LD) structure analysis was done in Haploview v4.2 software. Haplotype analysis was done in SVS, Golden Helix. Regulatory haplotypes were defined based on the following criteria: a coding variation, an UTR or splice variant, an ENCODE’s histone mark, a transcription factor binding site, DNase I hypersensitivity clusters, or an expressed quantitative trait loci (eQTL) in strong LD with peak SNP.

Several public databases and the published literature were used to analyze our study data, including the Encyclopedia of DNA elements (ENCODE), literature on eQTL, recent phase 3 release from 1000 Genome Project, and the PolyPhen/SIFT database [33-39]. The identification and potentially functional impact of an allele on gene transcription was derived from these sources.

We used SNP and Variation suite (SVS) of Golden Helix (version 7.6.8 win64, Golden Helix Inc., Bozeman, MT, USA, www.goldenhelix.com) for genetic analysis and Haploview v4.2 software for visualization of LD plots and haplotype analysis [40]. GraphPad Prism 6.0 software was used for statistical analyses and graphics. Means and standard deviations are presented for continuous outcomes and number and percent for binary outcomes (0, 1). Correlations between continuous variables were determined using Pearson’s r in GraphPad Prism 6.0. Discontinuous variables were compared by Fisher’s exact test. p values < 0.05 were considered significant. With power analysis, our goal was to determine whether any common genetic variant predicts a clinically significant outcome difference among treatment groups in CO-MED (either all groups combined or each group individually). A clinically meaningful difference (defined as a spread of at least three percentage points) [13] was detected on four outcomes: week 6 response, week 12 response, week 6 remission, and week 12 remission. Power calculation was conducted using the CaTS power calculator for GWAS. Applying a two-stage design with genome-wide scan, we calculated a power of 80% to detect genetic effects (allelic model) with a relative risk of 1.60 for SNPs with a MAF of at least 15% and under the assumption of 33% favorable treatment outcome.

Response at week 6 or week 12 was defined as a 50% or greater decrease in QIDS-SR total score from baseline to week 6 or week 12, respectively. Remission at week 6 was defined by two subsequent visits in which the patient had a QIDS-SR score of < 6 and/or < 8. Remission at week 12 was defined in a similar manner. Means and standard deviations are presented for continuous outcomes, and number and percentages for binary outcomes. The associations between treatment group, clinical outcomes, subject characteristics, and haplotype status were tested by χ² tests.

In the CO-MED Study, 459 subjects of the 665 enrolled provided DNA samples and were used in these studies. The age of participants was 43.37 ± 12.5 years. At baseline, approximately one-third of the sample received each treatment combination, as expected due to randomization. By week 6, 416 subjects remained in the study, again with approximately one-third in each group. By week 12, the number of participants dropped to 371 and elicited a slightly unbalanced distribution across groups: 34.50% in the ESCIT+PBO group, 31.27% in the BUP+ESCIT group, and 34.23% in the VEN+MIRT group. Demographic information, baseline clinical characteristics, and treatment outcomes at the 6- and 12-week time points are summarized in Table 1a. The primary “six-week response” GWAS analyses were conducted in the Caucasian group (due to the higher availability of genomic database information) using meta-analyses of responders and nonresponders (n = 210) who were adherent to the protocol. The secondary “twelve-week response” GWAS analyses were conducted using a subgroup of 119 subjects who were adherent to the protocol with remitters versus nonremitters. For the CO-MED GWAS study design, see Figure 1. A list of samples with responders versus nonresponders in Caucasian ethnicity groups are depicted in Table 1b.

Combined GWAS analyses were performed on samples corresponding to response at week 6 and response/remission at week 12, independent of treatment assignment. We observed moderate genetic association in the combined analysis of all three treatment groups for response at week 6 (online suppl. Figure 2A, B) and remission at week 12 (online suppl. Fig. 3A–C). The SNPs identified trended toward statistical significance (p value log^–5 and log^–6) with 6- and 12-week outcome status, respectively. Seventeen SNPs were associated with BUP+ESCIT, and 10 SNPs were associated with VEN+MIRT treatment groups. However, no SNP achieved a genome-wide level of significance (p value ≤10^–8) in these groups.

When separated by treatment group, however, strong and moderate genetic associations were seen with monotherapy and combination therapy, respectively. For example, ESCIT+PBO monotherapy revealed a potentially regulatory novel SNP, rs10769025, in the ALX4 gene on chromosome 11 that was strongly associated (p = 9.85925E^–08) with week 6 treatment response (Fig. 2a, b). The GWAS Manhattan plots evidence that ALX4 has the strongest association in the ESCIT+PBO group as compared to the combinations of medications. Other genes that reached statistically significant (p value of log^–5) level of genetic association with antidepressant treatment response in the three groups are summarized in online supplementary Table 1 but not further described in this report. We did observe a similar trend of ALX4 variants associating with monotherapy outcome at week 12 as well, although the strength of the association was weaker, plausibly due to attrition in responder subjects (online suppl. Fig. 1D). Results from genetic analysis with remission at week 12 across the three treatment groups are provided in online supplementary Table 2.

As described, the ALX4 locus is strongly associated with ESCIT+PBO response in Caucasians at weeks 6 and 12. However, since multiple ALX4 SNPs demonstrated association, we carried out haplotype analysis on potentially regulatory variants in strong LD with peak SNP. The inherent hypothesis for this analysis is that in a region of high linkage, several variants could show association with trait due to their LD with actual causal mutation. So, haplotypes can include all the mutations that may be causal or direct cause of association signal.

The ALX4 gene encodes a paired-like homeodomain transcription factor expressed in the mesenchyme of developing bones, limbs, hair, teeth, and mammary tissue [41, 42]. Mutations in this gene are known to play a role in craniofacial development, mesenchymal-epithelial communication, and hair follicle development. Deletion of a segment of chromosome 11 containing this gene, del (11) (p11p12), causes Potocki-Shaffer syndrome [42], characterized by craniofacial anomalies and intellectual disability [41]. Though there are limited reports regarding ALX4 gene involvement in MDD and antidepressant treatment response, the Psychiatric Genomics Consortium (PGC) data show that there is an ALX4 peak SNP involved in MDD versus controls in Caucasian ancestry (n = 9,240 MDD and n = 9,519 controls). Interestingly, in this data set a baseline association of ALX4 variants with depression disorder has also been observed (http://psychiatry.som.jhmi.edu/metamoodics/gb/mdgwas.php?gene = ALX4). They have shown that ALX4 SNP rs897005 show significant association in MDD cases. We also confirmed association of rs897005 SNP in our data. Furthermore, it is to be noted that index SNP rs897005 show good LD with other variants in the ALX4 gene, suggesting that actual causal variant could be forming a haplotype with the index and other potentially functional variants associated with trait. Our study peak SNP rs10769025 is also listed in their dataset. Although, it does not show the significant association in their analysis, but that could be due to many reasons such as different study design (inclusion, exclusion criteria), treatment heterogeneity among cases, and different statistical correction methods between two studies. Overall, some evidence of association of ALX4 variants with major depression disorder in Psychiatric GWAS consortium does support our observations in the present study [43].

We identified 15 SNPs in the ALX4 gene that showed strong and consistent associations with ESCIT monotherapy treatment outcome. To understand the potential functional impact of these variants on gene expression of ALX4 or other local genes, we annotated these variants with functional data from UCSC genome browser and ENCODE project. We found that SNP rs10769025 and rs7949067 have a strong regulatory score and are thus suggested to impact interactions of transcriptional complex molecules, i.e. GATA2, c-FOS, POLR2A, c-JUN, c-MYC, SMARCA4, and GATA1 with local DNA (Fig. 3a). These two ALX4 variants show strong eQTL effects for neighboring genes ACCS and EXT2 which are downstream of ALX4 (Fig. 3b). This means that responders with the C allele of rs10769025 have reduced expression of ACCS as compared to the nonresponders who carry the G allele. These eQTL effects were observed by Westra et al.[35], in PBMCs, and Fairfax et al.[34] in monocytes. Data from 1000 Genome Project also showed consistent eQTL results with the SNP rs10769025. Overall, these data suggest that higher ACCS expression in ALX4 is associated with lower response, and lower ACCS expression is associated with higher response in the ESCIT monotherapy group as compared to the other two combination groups (Fig. 3c–f). However, the mechanism of action needs to be further studied and validated.

Instead of an individual SNP, we generated haplotypes based on 7 potentially functional SNPs that showed strong association with response and were in strong linkage disequilibrium with each other. Haplotype association analysis showed that the CAAACTG haplotype is significantly associated (odds ratio = 3.4, p = 2.00E–04) with monotherapy responsiveness. The haplotype analysis is shown in Figure 2c, d.

Ingenuity pathway analysis (IPA) may provide insight regarding which genes are upregulated or downregulated in different signaling cascades that may be contributors in MDD. IPA was carried out to better understand the role of the ALX4 network in treatment response. IPA identified a functional pathway of ALX4 from 20 out of 24 associated genes (Fig. 4) that are already known to be significantly associated with MDD antidepressant treatment. ALX4 is known to indirectly affect CYPs, COMT, SLC6A4, PLG, TPH1, and HTR2A, which have been previously implicated in antidepressant responsiveness in MDD studies [44, 45]. In addition, ALX4 is also functionally related to CTNNB1, Lef1, Myc, PI3k, NFκB, and IL-1 signaling cascades. Acting though LEF-1, the signaling cascade that is centered on ALX4 involves key stress-signaling molecules, such as Myc, PI3K, NFκB, and IL-1.

Lastly, subjects were divided into three haplotype groups: responder, nonresponder, and heterozygous (responder haplotype refers to HAP1 sequence; nonresponder haplotype, HAP2 sequence; heterozygous haplotype are from either HAP3, 4, 5, 6, or 7). The number of subjects-by-haplotype status is shown for various subgroups in online supplementary Table 5A. Very few subjects were of the nonresponder haplotype, and they were not included in the clinical characteristics. Overall, week 6 response did not vary significantly by treatment group or haplotype status (p = 0.094) (online suppl. Table 5B). However, a higher response rate at week 6 was seen in the ESCIT+PBO group for those with the responder and heterozygous haplotype. Table 2 shows the distribution by haplotype status for several clinical outcomes and demographic and clinical characteristics for subjects with European ancestry in the ESCIT+PBO group. Week 6 response varied significantly by haplotype status (82.35 in those with responder haplotype versus 34.00% in those with heterozygous haplotype, p = 0.001) and were not attributable to differences in any subject characteristics.

GWAS conducted in STAR*D (Sequenced Treatment Alternatives to Relieve Depression), GENDEP (Genome-Based Therapeutic Drugs for Depression), MARS (Munich Antidepressant Response Signature), PGRN-AMPS (Mayo Clinic Pharmacogenomics Research Network-Antidepressant Medication Pharmacogenomics Study), and IPSC (International SSRI Pharmacogenomics Consortium) depression studies have demonstrated moderate effects of genetic variants on antidepressant treatment outcomes [9, 25, 28-30]. The reports support the prospect that examination of potential SNPs which impact genetic expression and/or protein function may help predict personalized treatment outcome. Here, we demonstrated that GWAS analysis using a subset of CO-MED participants identified a novel ALX4 SNP that was tightly associated with enhanced treatment outcome by ESCIT monotherapy in Caucasian patients.

While ALX4 is a novel GWAS hit in depression studies, there is clinical evidence for its role in moderating treatment outcome in MDD patients. ALX4 is most strongly known for its association with bone development, and MDD patients frequently have significant decreases in bone density [46]. Recently, Kadriu et al. [46] showed that acute ketamine administration corrects abnormal inflammatory bone markers in MDD patients. Peripheral bone inflammatory markers were studied in 28 treatment-resistant depression patients with ketamine infusion and 16 healthy controls. The markers were measured at four time points: baseline, postinfusion, day 1, and day 3. Ketamine significantly increased the plasma osteoprotegrin levels compared in depressed patients relative to healthy controls, suggesting a role of ALX4 in MDD. While drug administration differed, Kadriu’s findings demonstrate the potential for antidepressants affecting ALX4.

In addition to ALX4, the response-associated genetic variants of ALX4 showed cis-expression quantitative stress loci (cis-eQTLs) for 2 additional genes: ACCS (1-Aminocyclopropane-1-Carboxylate Synthase) and EXT2 (Exostosin Glycosyltransferase 2), which are downstream of ALX4. Both demonstrate an indirect association with MDD, as ACCS acts like PLP (pyridoxal phosphatase), an essential enzyme which serves as a cofactor in tyrosine and vitamin B₆ metabolism, the latter of which is associated with the tryptophan-serotonin pathway [47] and depression or depressive symptomology. Based on the published eQTL data on ACCS, higher expression of the enzyme is associated with the GG allele (nonresponder) and lower expression with CC allele (responder). EXT2, another gene regulated by ALX4 variants, has also previously been implicated in major depression through genome-wide linkage analysis of 656 families and 1,494 affected relative pairs genotyped at 418 microsatellite markers [48]. ALX4 is also known to be associated with several regulatory transcription factors that are highly relevant in MDD like GATA1, GATA2, c-FOS, c-JUN, and c-MYC [49-51]. Thus, this report in combination with previous findings lends support to exploration of genetic polymorphisms in ALX4 and associated genes for potentially regulatory effects on ESCIT responsiveness. We further investigated to check if this ALX4 is expressed in the brain using human brain transcriptome, datasets available on the public domain on the Web (http://hbatlas.org). Interestingly, we found that ALX4 is expressed in different regions of the brain such as cerebellar cortex, mediodorsal nucleus of the thalamus, striatum, amygdala, hippocampus, and 11 areas of neocortex (courtesy, http://hbatlas.org/hbtd/images/wholeBrain/ALX4.pdf). The graph shows the expression pattern of ALX4 (online suppl. Fig. 4), in the brain from the time of birth until 82 years. However, at birth, the expression is slightly higher in neocortex compared to other brain regions, but at later stages in lifetime, it almost remains the same. Hence, antidepressant treatment response-associated regulatory variants might be impacting the currently unknown function of ALX4 in the brain. Future functional studies will characterize it hopefully. So, we do not think that there is a lack of biological plausibility for ALX4.

ALX4 hits have not been identified previously in depression studies, although this discrepancy may simply be due to variation among study design or GWAS chips used. For example, in level 1 of STAR*D, citalopram (as opposed to its s-enantiomer used in CO-MED, ESCIT) was used, and the different drugs may account for varying alterations in gene expression levels. Further, the genetic analyses were run on different SNP chips with potentially variable sensitivity. Lastly, our current report was conducted in a defined Caucasian population, due to the relatively expansive genome-wide database available for this ethnicity.

As indicated by IPA, ALX4 is indirectly associated with the activity of several inflammatory and stress factors. Importantly, ALX4 acts through LEF-1 transcription binding, which itself is associated with proteins implicated in mechanisms of depression. For example, overexpression of LEF1 increases the promoter activity of the Indoleamine 2,3-dioxygenase (IDO1) gene. IDO is the catalyst for the rate-limiting step of tryptophan breakdown into kynurenine, a largely regarded detrimental metabolite involved with immunomodulatory processes, neuromodulatory signaling, and serotonin production [52]. Further, the IDO1 gene regulates IFN-γ and TNF-α through its interaction with cis-acting regulatory elements of the promoter region [52]. Thus, while ALX4 is most commonly associated with bone development, IPA suggests that altered ALX4 may act indirectly through different signaling cascades to decrease serotonin availability and increase neuroinflammation and aberrant neuronal signaling, all of which are implicated in the pathogenesis of MDD and responsiveness to treatment.

In the ESCIT monotherapy group, ALX4 haplotypes showed a strong association with week 6 responders (based on the QIDS-SR; p = 0.001). However, only a moderate effect was observed at the week 12 time point (70.59% responder haplotype, 55.77% heterozygous haplotype, p = 0.032). Thus, while ALX4 may predict early responsiveness to ESCIT monotherapy, it may not be ideal for identifying patients likely to reach complete, long-term remission. Further, in the BUP+ESCIT group, ALX4 association had only a mild effect and no strong genome-wide association, so the benefit of ESCIT augmentation strategies still need to be elucidated. Lastly, like most other GWAS reports in MDD, data were analyzed post hoc, thereby necessitating validation of these results in a novel population, which may be difficult provided the heterogeneity among studies and reagents, as described previously in the discussion.

As the search for moderators and predictors of treatment response to antidepressants progresses, elucidation of the contributions of noncoding regulatory genetic variants in the human genome is increasingly important. The present report used genetic analyses among CO-MED participants and identified a regulatory haplotype in the ALX4 gene which has cis-eQTL effects and is strongly associated with treatment response to ESCIT monotherapy in Caucasians. Future clinical trials, such as EMBARC (Establishing Moderators and Biosignatures of Antidepressant Response for Clinical Care) [53], will further investigate ALX4 as a moderator of treatment outcome, with the benefit of a larger sample size and functional validation through mechanistic studies and targeted deep sequencing. In addition, functional studies based on SSRI treatment versus PBO in Caucasian ethnicity with larger sample size will give us a wider swath for understanding the functional role of ALX4, and transcriptome analyses using postmortem tissue brains samples followed up with medications might help us to characterize the functional role of ALX4 in the brain based on treatment response.

This research was funded, in part, by NIMH (N01 MH-90003) and The Hersh Foundation. We wish to thank Dr. Li’s group at the Microarray and Genomics core for performing the SNP genotyping analyses. We thank Dr. Carol Tamminga, Chairman, and Department of Psychiatry for the support and encouragement. The authors are grateful to Jennifer Furman, PhD, and Augustus John Rush, MD, for their substantial help in editing the manuscript and providing valuable feedback, and Jeremy A. Kee, MA, for his editorial and administrative assistance.

The authors have no ethical conflicts to disclose.

Dr. Trivedi has in the past 24 months consulted or served on the advisory board of Alkermes Inc., Akili Interactive Inc., Allergan, Arcadia Pharmaceuticals Inc., AstraZeneca, Brintellix, BMS, Cerecor, Global Medical Education Inc., Health Research Associates, Lundbeck, Medscape, MSI Methylation Sciences Inc., Merck, Naurex Inc., Navitor, Nestle Health Science – Pamlab Inc., One Carbon Therapeutics, Otsuka America Pharmaceuticals Inc., PamLab, Pfizer Inc., Roche, SHIRE Development, Takeda Pharmaceuticals Inc.; has conducted research activities with NIMH, NIDA, J and J, Janssen Research and Development LLC; and has claimed royalties from Janssen Research and Development LLC. Neither Ms. Mayes, nor Drs. Gadad, Raj, Jha, Carmody, Dozmorov, and Wakeland have any conflicts to report.

Drs. Trivedi, Gadad, Raj, Carmody, Jha, Wakeland had full access to all of the data in the study and take responsibility for the integrity of the data and the accuracy of the data analysis.

Study concept and design: Trivedi; acquisition of data: Trivedi, Gadad, Raj, Jha, Carmody, Wakeland; analysis and interpretation of data: Trivedi, Gadad, Raj, Carmody, Jha, Wakeland; drafting of the manuscript: Trivedi, Gadad, Raj; critical revision of the manuscript for important intellectual content: Drs. Trivedi, Gadad, Raj, Jha, Carmody, Wakeland; statistical analysis: Carmody, Raj, Gadad, Trivedi.