Replicable and Coupled Changes in Innate and Adaptive Immune Gene Expression in Two Case-Control Studies of Blood Microarrays in Major Depressive Disorder

Background Peripheral inflammation is often associated with major depressive disorder (MDD), and immunological biomarkers of depression remain a focus of investigation. Methods We used microarray data on whole blood from two independent case-control studies of MDD: the GlaxoSmithKline–High-Throughput Disease-specific target Identification Program [GSK-HiTDiP] study (113 patients and 57 healthy control subjects) and the Janssen–Brain Resource Company study (94 patients and 100 control subjects). Genome-wide differential gene expression analysis (18,863 probes) resulted in a p value for each gene in each study. A Bayesian method identified the largest p-value threshold (q = .025) associated with twice the number of genes differentially expressed in both studies compared with the number of coincidental case-control differences expected by chance. Results A total of 165 genes were differentially expressed in both studies with concordant direction of fold change. The 90 genes overexpressed (or UP genes) in MDD were significantly enriched for immune response to infection, were concentrated in a module of the gene coexpression network associated with innate immunity, and included clusters of genes with correlated expression in monocytes, monocyte-derived dendritic cells, and neutrophils. In contrast, the 75 genes underexpressed (or DOWN genes) in MDD were associated with the adaptive immune response and included clusters of genes with correlated expression in T cells, natural killer cells, and erythroblasts. Consistently, the MDD patients with overexpression of UP genes also had underexpression of DOWN genes (correlation > .70 in both studies). Conclusions MDD was replicably associated with proinflammatory activation of the peripheral innate immune system, coupled with relative inactivation of the adaptive immune system, indicating the potential of transcriptional biomarkers for immunological stratification of patients with depression.

Depression and inflammation are often associated with one another. Depressive symptoms in a large population sample were significantly related to blood concentrations of C-reactive protein (CRP; odds ratio w1.8 for depressive symptoms in people with CRP . 3 mg/L vs. CRP , 1 mg/L) (1). Multiple case-control studies of major depressive disorder (MDD) have reported increased peripheral blood concentrations of CRP (Cohen's d w0. 50) and proinflammatory cytokines such as interleukin 6 (d w0. 50) and tumor necrosis factor (d w0. 40) in MDD (2,3). The prevalence of comorbid depression is increased in many nonpsychiatric inflammatory disorders (4).
There is growing evidence for a causal effect of inflammation on depression. Peripheral inflammation precedes the emergence of depressive symptoms in longitudinal epidemiological studies (5) and in about 30% of patients receiving proinflammatory interferon-a treatment for hepatitis C (6,7). Experimental challenge with peripheral proinflammatory stimuli in animals robustly induces a syndrome of illness behavior and anhedonia that approximates depressive symptoms (8). Peripheral immune cells and cytokines are known to mediate signals across the blood-brain barrier by several mechanisms (9). Activation of microglia can locally amplify the effects of even a weak peripheral proinflammatory signal on neuronal function and behavior (10).
These observations suggest that pharmacological disruption of peripheral proinflammatory signals could be therapeutically effective, at least for a subgroup of patients with depression. However, it is most unlikely that any single anti-inflammatory drug will prove to be superior to all existing treatments for all patients (11). Only about a third of patients with depression have biological evidence for peripheral inflammation [e.g., CRP . 3 mg/L (12)], and an antiinflammatory drug seems likely a priori to be most effective for an inflamed subgroup of patients with depression. There are many markers of the peripheral immune system that can be conveniently measured in venous blood samples from patients with depression, including cytokines, CRP, and other proteins; cell counts from flow cytometry; and gene transcription (13).
Transcriptional (messenger RNA [mRNA]) biomarkers have the potential advantages of assay stability (compared with cytokines) and target specificity (compared with cell counts). However, previous case-control studies of peripheral blood gene expression and MDD have been inconsistent (14)(15)(16)(17)(18) (Supplemental Table S1). To date, the Netherlands Study of Depression and Anxiety (NESDA) (14) is the largest single case-control study of the peripheral blood transcriptome in MDD (882 current patients with MDD, 635 remitted patients with MDD, and 331 control subjects). The NESDA reported statistically significant (false discovery rate [FDR] = 10%) differential expression of 129 genes enriched for interleukin 6 and natural killer (NK) cell signaling pathways (19).
We were primarily motivated by the hypothesis that MDD is associated with peripheral blood transcriptional markers of innate immune system activation (8,20). We were also concerned to focus on results that were more likely to replicate across case-control studies of gene expression in MDD. We report Affymetrix microarray data on 18,863 probes from two independently designed and conducted case-control studies of MDD. We used a Bayesian method to identify genes that were differentially expressed in both studies. Focusing on a consensus set of 165 genes, we investigated the functional significance of the genes that were differentially (over-or under-) expressed in cases compared with controls. We also explored the secondary hypothesis that innate immune system activation is coupled to relative inactivation of the adaptive immune system in patients with MDD (21)(22)(23).

Samples
We analyzed data from two case-control studies of depression: the GlaxoSmithKline-High-Throughput Disease-specific target Identification Program (GSK-HiTDiP) study and the Janssen-Brain Resource Company (Janssen-BRC) study. Other aspects of these studies have been previously reported (24)(25)(26)(27)(28); demographic and clinical details on the samples are provided in Supplemental Table S2.
GSK-HiTDiP. This study was designed primarily to test for an association between genetic (DNA) variation and diagnosis of depression. Minimal sociodemographic and clinical data were collected, but microarray data were available for analysis (after quality control) from whole-blood samples stored for less than 6 years on a sample comprising 113 patients with MDD prospectively balanced for comorbid anxiety disorder [57 with generalized anxiety disorder diagnosed by the Mini-International Neuropsychiatric Interview (29) and 56 without anxiety disorder] and 57 healthy control subjects. All participants provided informed consent in writing. The study was approved by an independent ethics review board.
Janssen-BRC. This study was designed primarily for biomarker discovery. Microarray data were available for analysis (after quality control) from whole-blood samples stored for less than 1 year on a sample comprising 94 patients with MDD (40 with generalized anxiety disorder diagnosed post hoc by the Mini-International Neuropsychiatric Interview and 54 without anxiety disorder) and 100 healthy control subjects. Additional data on melancholic symptom severity, anxiety, substance use, and body mass index (BMI) were available for patients with MDD. All participants provided informed consent in writing. The study was approved by an independent ethics review board.
Whole-blood samples from both studies were analyzed using the Affymetrix Human Genome U133 plus 2.0 array. We applied identical quality control, normalizing, and annotation algorithms to the two datasets, resulting in the estimation of mRNA expression at each of 18,863 unique probes for each participant.

Differential Expression Analysis
To determine differential gene expression between cases and controls, we adopted the same strategy for both studies; see Supplemental Figure S1 for a schematic overview of the data analysis strategy. For each gene i = 1, ., 18,863, we fitted a linear regression model that included group (coded Gr; twolevel factor, case/control), batch (B; two-level factor, 1/2), gender (Ge; two-level factor, male/female), age (Ag; continuous), and presence/absence of anxiety (An; two-level factor, 0/1) as covariates. Denoting samples by j = 1, ., n and the design matrix by X, the model is where ε i w N(0, s i 2 ). For identifiability, we imposed the con- The model was fitted using the R package limma (30). Subsequently, we tested the null hypothesis H i0 : b i2 Gr = 0; that is, there is no difference in expression of the ith gene between the two groups using the moderated t statistic (31). For each of the 18,863 probes in each study, the p value was generated using the asymptotic approximate distribution (and also nonparametrically by a permutation algorithm; see Supplemental  Tables S3 and S4).

Combining p Values for Differential Gene Expression From Two Studies
To identify MDD-related genes that were replicated in both the GSK-HiTDiP and Janssen-BRC datasets, we set the p-value threshold for significance of differential expression of each gene in each study to optimize in some sense the number of Peripheral Immune Transcripts and Depression Biological Psychiatry January 1, 2018; 83:70-80 www.sobp.org/journal genes that were differentially expressed in both studies compared with the number of coincidental differences expected by chance. To do this, we used the method of Blangiardo and Richardson (32), implemented in the R package sdef (33), and specified the p-value threshold as q 2 , which represents the largest (most lenient) p value , .05 for which there are at least twice as many significant case-control differences in common between the two studies as expected by chance.

Gene Coexpression Network Analysis
We used weighted gene correlation network analysis (34) to construct a normative coexpression network in which nodes represent genes and weighted edges represent correlations between the expression of pairs of genes in healthy control subjects (34). To maximize the amount of data available for this estimation, and to ensure that the resulting network was representative of both studies, we included the healthy control data from both the GSK-HiTDiP and Janssen-BRC studies (N = 157 in total). We used the consensus weighted gene correlation network analysis method to construct a weighted undirected graph that could be decomposed into a set of modules of strongly coexpressed genes (34). Eigengenes were used to test for case-control differences in the gene expression profile on average within each module of the normative weighted gene correlation network analysis transcriptome (Supplemental Table S11).

Gene Ontology Enrichment Analysis
We performed enrichment analysis for various gene lists such as the list of genes differentially expressed in both studies and the list of genes comprising each module of the coexpression network. We used the R package topGO (35) to determine whether these lists were significantly enriched for specific Gene Ontology (GO) terms (36) using a stringent Bonferroni correction (q , .05) for multiple comparisons across all 10,124 terms, resulting in a p-value threshold of 4.94 3 10 26 .

Protein-Protein Interaction Mapping
We used the Search Tool for the Retrieval of Interacting Genes/Proteins [http://string-db.org (37)] to determine the network of known protein-protein interactions among the genes differentially expressed in both studies (38).

Cell Type-Specific Gene Expression Analysis
We investigated the genes differentially expressed in both studies in relation to a third independent microarray (Affymetrix) dataset designed to examine cell type-specific expression patterns (39). In particular, we focused on the expression of MDD-related genes across 37 sorted cell samples of the following immune cell types: erythroblasts (8 samples), monocytes (6 samples), monocyte-derived dendritic cells (5 samples), neutrophils (3 samples), B cells (4 samples), CD41 T cells (5 samples), CD81 T cells (5 samples), and NK cells (6 samples). We used the Bio-Layout Express 3D software (40,41) to visualize clusters or subgroups of genes that shared similar patterns of expression across different cell types.

RESULTS
Combining Differential Expression Across Two Case-Control Studies of MDD Separate differential expression analyses of the GSK-HiTDiP and Janssen-BRC datasets yielded two lists of p values for the same set of 18,863 gene expression probes (Supplemental Tables S3 and S4). The Bayesian method identified q 2 = .0246 as the largest p value threshold associated with twice the number of genes differentially expressed in both studies compared with the number of coincidental case-control differences expected by chance. At this threshold, 173 genes were differentially expressed in both studies. We further refined this gene list by restricting attention to the 165 consensus genes (95%) that showed the same direction of overexpression (UP) or underexpression (DOWN) in the two casecontrol studies (Supplemental Table S5).
These findings were corroborated by Fisher's chi-square test for combining p values, which identified 393 genes that were differentially expressed in both studies with FDR = 10%, of which 146 were included in the list of replicably, concordantly, differentially expressed genes defined by Bayesian analysis (Supplemental Table S6).
Furthermore, the set of 165 replicably concordant genes (henceforth the MDD-165 consensus set) was partially corroborated by the prior results of the NESDA case-control study (14). Considering the 15,830 genes that were measured by the microarrays used in all three studies (NESDA, Janssen-BRC, and GSK-HiTDiP), 150 of the MDD-165 consensus set were included, of which 7 genes (CD247, PRKCH, PGLYRP1, NFATC2, ST6GAL1, MAPK14, and MTSS1) were also differentially expressed, with the same sign of fold change, in the NESDA study at FDR = 5%. When the FDR threshold for the NESDA study was relaxed to 10% and 20%, 18 and 45 of the 150 genes differentially expressed in both the GSK-HiTDiP and Janssen-BRC studies, respectively, were also differentially expressed in the NESDA study (Supplemental Table S5). By comparison, we expect by chance to find 1.95, 8.30, and 30.40 genes in common with the NESDA study at FDR thresholds of 5%, 10%, and 20%, respectively. In contrast, there was less overlap between the MDD-165 consensus set and the top 29 most significantly differentially expressed genes reported in another large casecontrol study of MDD (16); only SRSF5 was underexpressed in both of these lists.
To assess the robustness of our results to key modeling assumptions, we conducted two sensitivity analyses. First, we repeated the analysis using a model for differential expression that did not include comorbid anxiety as an explanatory factor. Only 24 genes were identified by Bayesian analysis as significantly coexpressed; of these 24 genes, 21 were concordant for sign of fold change, and 15 of these 21 genes were also included in the MDD-165 consensus set (Supplement). The reduced sensitivity when not adjusting for anxiety likely reflects the importance of refining or homogenizing the clinical phenotype in the search for biomarkers.

Peripheral Immune Transcripts and Depression
Second, we repeated the analysis using a model for differential expression that included BMI as a covariate. BMI was not measured in the GSK-HiTDiP study, so this analysis was restricted to 74 patients and 80 control subjects from the Janssen-BRC dataset. Inclusion of BMI in the model tended to increase the p values for MDD-related differences, and only 73 genes in the MDD-165 consensus set were differentially expressed at the probability threshold q 2 = .0246 The top 20 GO terms for biological processes are ranked according to their p values by Fisher's exact test for significant enrichment. Bonferroni correction specifies a p-value threshold of p = 4.94 3 10 26 to achieve significance at q , 0.05. Such a correction is well known to be too stringent given that many of the GO terms are correlated; however, the top 17 terms were significantly enriched in the overexpressed UP gene set after (Supplemental Figure S3). The reduced sensitivity of the analysis including BMI as a covariate likely reflects the known proinflammatory effect of obesity and the reduced sample size available.

Gene Set Enrichment Analyses Across GO Terms
We conducted a series of enrichment analyses to functionally characterize the MDD-165 consensus set. The top 10 most significantly enriched GO terms all were related to the immune response to infection (Bonferroni q , .05 for the top 9 terms) (Supplemental Table S7). Among the MDD-165 consensus set, 90 genes were overexpressed in cases compared with controls. The UP gene list was significantly enriched for GO terms related to the innate immune system and response to infection (Bonferroni q , .05 for the top 17 terms) ( Table 1). In addition, 75 genes were significantly underexpressed in cases compared with controls. As we hypothesized based on prior literature (19), this gene list was enriched for GO terms related to the adaptive immune system (p # 10 23 for top 20 GO terms; Table 1), but there was no significant enrichment of the DOWN genes at the more stringent Bonferroni threshold.

Protein Interaction and Whole-Genome Transcriptional Network Analyses
We used Search Tool for the Retrieval of Interacting Genes/ Proteins analysis (37) to visualize the network of proteinprotein interactions between the proteins coded by the MDD-165 consensus set. These genes were significantly enriched for protein-protein interactions (permutation test, p = 7 3 10 24 ) that were concentrated around mitogen-activated protein kinase 14 (MAPK14) and matrix metalloproteinase 9 (MMP9), which thus can be regarded as the most highly interactive hubs of this immune signaling network (Figure 1). Many of the proteins in this network are coded by genes that were independently identified by the NESDA study as differentially expressed in MDD; see Figure 1.
We constructed a graph representing significant coexpression of a pair of genes (nodes) as an edge drawn between them. As previously reported (42), this whole-genome transcriptional network or transcriptome had complex topological properties, including a community structure comprising modules of coexpressed genes enriched for specific GO terms (Figure 2). The UP genes were concentrated in the normative module (red) significantly enriched for innate immune response GO terms (e.g., myeloid cell activation involved in immune response, q = 5 3 10 25 , Bonferroni corrected), whereas the DOWN genes were concentrated in a normative module (pink) specialized for translation-based terms ( Figure 2). We confirmed in a hypothesis-driven analysis that the module associated with DOWN gene expression was also significantly enriched for adaptive immune response terms (e.g., T-cell differentiation, p = .0002). There were significant case-control differences in eigengene scores summarizing expression of all genes within modules enriched for innate immune response (red, p = .031; yellow, p = .015), translation (pink, p = .002), and one additional module with no significant enrichment terms (tan, p = .025); see Supplemental  Table S11.

Cell Type-Specific Expression Patterns for the MDD-165 Consensus Set
To assess the cellular specificity of the MDD-165 consensus set, we used microarray data from an independent study of cell type-specific gene expression across eight major classes of immune cells in healthy control subjects (39). Remarkably, 89 of the 90 UP genes formed four clusters of gene coexpression, three of which represented strong overexpression of a subset of UP genes in one or two myeloid cell classes ( Figure 3). Likewise, 71 of the 75 DOWN genes formed two clusters of gene coexpression, each representing cell-specific overexpression of a subset of DOWN genes in one or two lymphoid or erythroblast classes ( Figure 3).  (38). Only high-confidence (.0.7) links are retained, and disconnected genes are not shown. Red (blue) nodes correspond to genes over-(under-) expressed in major depressive disorder in both the GlaxoSmithKline-High-Throughput Disease-specific target Identification Program and Janssen-Brain Resource Company datasets. Smaller inner circles highlight proteins that are coded by genes also differentially expressed and with the same sign of fold change in a third large independent case-control study of major depressive disorder (Netherlands Study of Depression and Anxiety [NESDA]) (14) thresholded to control false discovery rate (FDR) at 20% (white circles), 10% (gray circles), and 5% (black circles).

Coupled UP and DOWN Gene Expression in Patients With MDD
To explore the extent to which UP gene overexpression was related to DOWN gene underexpression, we estimated the mean UP gene expression (over all 90 differentially overexpressed genes) and the mean DOWN gene expression (over all 75 differentially underexpressed genes) for each of 113 patients with MDD in the GSK-HiTDiP study and each of 94 patients with MDD in the Janssen-BRC study. In both studies, UP and DOWN expression scores were strongly negatively correlated (Figure 4), indicating that these case-control differences were coupled at the level of individual patients (Janssen-BRC: r = 2.82, p , 10 24 ; GSK-HiTDiP: r = 2.74, p , 10 24 ). Given the strong linear relationship between UP and DOWN gene expression in patients with MDD, we mapped each patient's data to a point on the fitted regression line, thereby characterizing each patient by a bioscalar that locates it on a single axis or dimension of coupled UP/DOWN gene expression; see Figure 4A (inset). As expected, there were significant case-control differences in the MDD-165 bioscalar in both studies (Cohen's d = 0.68 for GSK-HiTDiP and 0.58 for Janssen-BRC, corresponding to a medium effect size; Figure 4C, D), with patients on average having more positive values, indicating a greater shift in the direction of coupled innate activation and adaptive inactivation. Receiver operating characteristic analysis indicated that cases and controls were classified with an area under the curve of 0.71 for the GSK-HiTDiP study and an area under the curve of 0.67 for the Janssen-BRC study ( Figure 4C, D). This moderately accurate classification performance reflects the fact that the bioscalar distributions of cases and controls are overlapping. We also assessed the performance of the bioscalar not to discriminate cases from controls (which is ultimately a clinical diagnostic decision) but rather to identify the top third "most inflamed" patients with MDD. To do this, we defined a cutoff value for the bioscalar corresponding to the top tertile critical value of its distribution in each study (0.57). Participants with MDD-165 . 0.57 were classified correctly as belonging to such an inflamed MDD subgroup (sensitivity = 100% by definition) with a high specificity (93% and 88% in GSK-HiTDiP and Janssen-BRC studies, respectively).
We provisionally explored correlations between the MDD-165 bioscalar and sociodemographic or clinical differences among the patients with MDD in the Janssen-BRC study. BMI (r 72 = .20, p = .099), the CORE total score for melancholic symptom severity (r 88 = .18, p = .097), and the presence or absence of comorbid substance abuse disorder (two-tailed t 63 = 1.9, p = .07) were not significantly associated with bioscalar scores. There was a significant difference in bioscalar scores between subgroups of patients with MDD with or without comorbid anxiety disorder (t 83 = 22.4, p = .02).

DISCUSSION
We have reported the differential expression analysis of wholegenome microarray data from two independent case-control studies of patients with MDD compared with healthy control subjects. Using Bayesian methods, we identified a set of 165

Peripheral Immune Transcripts and Depression
Biological Psychiatry January 1, 2018; 83:70-80 www.sobp.org/journal consensus genes that were replicably associated with MDD, sharing the same direction of fold change in the cases compared with controls in both studies. The robustness of these results was further supported by comparison with the results of the largest prior case-control microarray study of depression, the NESDA study (14); several genes differentially expressed in the NESDA study were likewise differentially expressed in these data, for example, MAPK14 and MMP9. Both of these genes were consistently overexpressed in patients with MDD and both code proteins (MAPK14 and matrix metallopeptidase 9) that were hubs of a network of interactions between immune signaling proteins coded by many of the other differentially expressed genes. A MAPK14 inhibitor has been tried for treatment of major depression but did not demonstrate consistently significant effects on symptom rating scales compared with placebo at the single dose tested (43).
The consensus set of 165 genes associated with MDD was divided into two approximately equal-sized subsets; here, 90 so-called UP genes were overexpressed in patients and 75 socalled DOWN genes were underexpressed in patients. The overexpressed UP genes were significantly enriched for GO terms related to the response to infection and the innate immune system. These gene transcripts were not functionally unrelated to or independent of each other. Most UP genes were affiliated with the module of the normative gene transcriptional network specialized for innate immune response. Smaller clusters of highly correlated UP genes were typically enriched in neutrophils, monocytes, and monocyte-derived dendritic cells. These results are compatible with prior data implicating activation of the innate immune system and increased proinflammatory cytokine signaling in (some people with) depression (8,20).
In contrast, we found that the underexpressed DOWN genes were enriched for GO terms related to T-cell function and adaptive immunity. The DOWN genes were affiliated with modules of the transcriptome partially specialized for adaptive immune function, and smaller clusters of strongly coexpressed DOWN genes were enriched in T cells, B cells, and NK cells. These results are compatible with prior data implicating relative suppression of the cellular immune system in (some people with) depression, including evidence for decreased NK cell cytotoxicity and decreased proliferation of lymphocytes challenged with mitogens in vitro (44), and with recent data suggesting that deficiencies of NK and T cells co-occur with inflammatory monocyte activation as related phenomena in the same patients with MDD (21,22,45,46).
In both studies, we also observed a strong statistical association between mean UP gene overexpression and mean DOWN gene underexpression in each patient. This provided robust evidence that these complementary immunophenotypes are indeed coupled at the level of an individual patient. It also motivated the concept that each individual patient may be located by a single number (bioscalar) on a spectrum of coupled change in innate and adaptive immune system function. Patients with MDD had significantly higher scores on this bioscalar, indicating a mean shift to relatively increased UP gene expression and decreased DOWN gene expression. The subgroup of most-inflamed patients with MDD, defined as the top tertile of the MDD-165 bioscalar distribution, were identified with high sensitivity and specificity, suggesting that this may in the future prove to be a useful biomarker for defining an abnormally inflamed subgroup of patients with MDD.
So far, we have used the terms overexpression and underexpression simply to describe higher and lower levels of measured mRNA in cases compared with controls. At least three (mutually nonexclusive) explanatory factors are plausible: cellular, genetic, and environmental. First, case-control differences in expression could reflect differences in cell counts. For example, increased monocyte counts have been reported in MDD (47) and could cause apparent overexpression of innate immune system genes measured in whole-blood samples (and relative underexpression of adaptive immune genes). It will be important in the future to combine cytometry and transcriptional measurements in the same patients and to measure case-control expression differences in sorted cell samples. Second, mRNA changes could be quantitative traits determined by DNA variation at expression quantitative trait loci. This genetic explanation would require that (at least some) patients with MDD had a consistent profile of allelic variation in inflammation-related genes (48), but large genome-wide association studies have so far failed to identify genetic variants robustly linked to risk for major depression (49). The absence of significant genome-wide association study findings in MDD could be regarded as problematic for an expression quantitative trait locus interpretation of expression changes, or it could be discounted on the grounds that classically designed association studies were underpowered to detect DNA variations occurring in only a subgroup of immunologically dysfunctional patients. Third, gene expression changes could have been induced in patients by exposure to some shared environmental stimulus. The risk for depression is associated with adverse events in the biological and social environments such as infection, childhood abuse, and bereavement (50). In addition, psychosocial stress has been linked to peripheral immune state changes, such as increased proinflammatory cytokines and monocyte activation (51), that are compatible with overexpression of innate immune system genes.
It is a methodological limitation that the GSK-HiTDiP study, designed primarily to identify risk genes for MDD, did not provide data on severity of depressive symptoms or BMI. The case-control comparisons are not controlled for potentially confounding effects of cigarette smoking, race, comorbid medical disease, or socioeconomic status on peripheral immune status. Future studies of transcriptional biomarkers in MDD, with more detailed clinical and immunological phenotyping and more complete control of potential confounds, will be required to evaluate the generalizability of these results.
In short, we have reported replicable new evidence in support of peripheral immune gene expression markers for MDD. It may be fruitful in the future to investigate coupled overexpression of innate immune genes and underexpression of adaptive immune genes as a predictor of antidepressant response to novel immunotherapeutics. We thank David Hume for comments on the manuscript; Fiona Kelly and Stewart Bates for supporting access to the GSK-HiTDiP dataset; Pieter Peeters, Wilhelmus Drinkenberg, and William Talloan for supporting access to the Janssen-BRC and HiTDiP datasets; and Evian Gordon from the BRC.

ACKNOWLEDGMENTS AND DISCLOSURES
The microarray and clinical data from the Janssen-BRC study are accessible to scientists through the independent YODA Project process (http:// yoda.yale.edu). The GSK-HiTDiP data are accessible on GEO (https://www. ncbi.nlm.nih.gov/geo) as GSE98793. For more information on combining the datasets please contact the YODA Project (yodap@yale.edu).
GMW and WCD are employees and stockholders of Janssen Research & Development of Johnson & Johnson. SK and RH are employees and stockholders of GlaxoSmithKline (GSK). ETB is employed half-time by GSK and half-time by the University of Cambridge; he holds stock in GSK. ETB is also a shareholder in the BRC. CMP has received research funding from Johnson & Johnson and speaker fees from Lundbeck. The remaining authors all report no biomedical financial interests or potential conflicts of interest.