Genome-wide analysis suggests the importance of vascular processes and neuroinflammation in late-life antidepressant response

Antidepressant outcomes in older adults with depression is poor, possibly because of comorbidities such as cerebrovascular disease. Therefore, we leveraged multiple genome-wide approaches to understand the genetic architecture of antidepressant response. Our sample included 307 older adults (≥60 years) with current major depression, treated with venlafaxine extended-release for 12 weeks. A standard genome-wide association study (GWAS) was conducted for post-treatment remission status, followed by in silico biological characterization of associated genes, as well as polygenic risk scoring for depression, neurodegenerative and cerebrovascular disease. The top-associated variants for remission status and percentage symptom improvement were PIEZO1 rs12597726 (OR = 0.33 [0.21, 0.51], p = 1.42 × 10−6) and intergenic rs6916777 (Beta = 14.03 [8.47, 19.59], p = 1.25 × 10−6), respectively. Pathway analysis revealed significant contributions from genes involved in the ubiquitin-proteasome system, which regulates intracellular protein degradation with has implications for inflammation, as well as atherosclerotic cardiovascular disease (n = 25 of 190 genes, p = 8.03 × 10−6, FDR-corrected p = 0.01). Given the polygenicity of complex outcomes such as antidepressant response, we also explored 11 polygenic risk scores associated with risk for Alzheimer’s disease and stroke. Of the 11 scores, risk for cardioembolic stroke was the second-best predictor of non-remission, after being male (Accuracy = 0.70 [0.59, 0.79], Sensitivity = 0.72, Specificity = 0.67; p = 2.45 × 10−4). Although our findings did not reach genome-wide significance, they point to previously-implicated mechanisms and provide support for the roles of vascular and inflammatory pathways in LLD. Overall, significant enrichment of genes involved in protein degradation pathways that may be impaired, as well as the predictive capacity of risk for cardioembolic stroke, support a link between late-life depression remission and risk for vascular dysfunction.


Supplementary Methods
Genome-wide analysis suggests the importance of vascular processes and neuroinflammation in late-life antidepressant response 1. Quality Control

Genotyping and Quality Control
Endorsed by the Psychiatric Genomics Consortium, the PsychArray contains 265,000 tag SNPs from the Illumina HumanCore BeadChip, 245,000 markers from the HumanExome BeadChip and 50,000 SNPs which have previously been associated with common psychiatric disorders.
From our initial sample of 453 patients, we excluded 91 individuals who were withdrawn for several reasons (see Supplementary Figure 1).Thirteen individuals were removed due to excessive heterozygosity based on an inbreeding coefficient greater than two standard deviations from the sample mean.We checked individuals for discordance between self-reported and genetic sex, during which one individual was excluded due to Ychromosome abnormalities.Lastly, one individual showed excessive relatedness ( ̂ > 0.185, i.e., second cousins), and two individuals had excessive missing genotypes (more than 10%).Overall, we excluded 17 individuals who failed genetic quality control based on one or more criteria, resulting in a final sample of 345 individuals who entered imputation.
Genetic ancestry was assessed using multidimensionality scaling in PLINK v.1.9. 1 First, we pruned SNPs based on linkage disequilibrium (LD) using a 50 SNP window and shifting by five SNPs with an r 2 threshold of 0.2.We also removed predefined high-LD regions. 2 Outliers were defined as individuals' principal components 1 and 2 loadings beyond the six standard deviations from the centre of the ancestral cluster. 3Divergent individuals (i.e., discrepant between self-reported and genetic ancestry) that were visibly clustering well with other continental populations were reclassified appropriately (e.g., African-ancestry vs. European-ancestry, N=2); however, those with ambiguous ancestries were reclassified as "admixed" (N=5, see Supplementary Figures).

Imputation
Per marker, we exclude variants based on violations of Hardy-Weinberg equilibrium at p < 10 -7 , low genotyping call rate < 95%, and low minor allele frequency < 1%.Whole-genome imputation was conducted using the genipe pipeline, 4 which uses IMPUTE2 v2.2 5 in 5-Mb segments per chromosome after pre-phasing with SHAPEIT2 6 and the 1000 Genomes reference panel (Phase 3). 7We filtered for biallelic SNPs and retained those with an imputation score 0.7, completion rate > 90%, and a minor allele frequency ≥ 5%.Hard genotype calls were made using a probability threshold of 90%.We included 4,471,676 bi-allelic variants with a genotyping rate of 99.1% in 329 individuals (307 European-ancestry, 22 African-ancestry).

Secondary Analyses
Linear mixed-effects models For top-associated variants, we constructed linear mixed-effects models across the eight time-points of treatment.Our outcome of interest was the MADRS score at the end of treatment (week 12).We included fixed effects from age, sex, current depressive episode duration, baseline MADRS score, and additive SNP genotype (i.e., 0, 1, 2), as well as random effects from individual ID and site of recruitment (denoted below by "|").

Time-to remission analysis
First, we obtained the Kaplan-Meier survival curve for the entire population, as well as stratified by the SNP.We used time-to remission as our 'failure event' and assumed that censoring time is independent of failure time.The estimated survival curves across SNP genotypes or alleles were compared curves using a Mantel-Haenszel,  2 test with one degree of freedom at an α=0.05.Subsequently, to assess the effects of baseline covariates on time-to remission, we fit a Cox proportional hazards regression.For this, we fit two models: Model 1:   () =  0 ()  1 + +  2 ()+  3 1+  4 (2)+ 5  +  6  + 7 1+ 8 2 Model 2:   () =  0 ()  1 + +  2 ()+  3 1+  4 (2)+ 5  +  6  + 7 1+ 8 2+ 9  These two models were then compared using a likelihood ratio test (LRT, α=0.05) to assess whether the inclusion additive SNP genotypes in estimating the survival curves improves model fit.For the final Cox model (i.e., Model 2), the overall model significance was assessed using an LRT compared to a model only including the intercept.Also, ANOVA was used to assess the significance of the individual coefficients.
Given that the Cox regression assumes proportional hazards over time for validity, this assumption was tested for all predictors in the model by assessing predictor interaction with time.We obtained the Pearson productmoment correlation (⍴) between the scaled Schoenfeld residuals and log(time) for all variables, including additive SNP genotype.In addition, we also considered the global test for all interactions.These tests were considered at α=0.05, which would indicate violators of the proportionality assumption.In addition, we plotted the scaled Schoenfeld residuals against transformed time with a smoothed line, including ± two standard deviations to inspect any non-proportional effects in the model across the exposure and covariates.Any observed systematic deviations from a horizontal line were interpreted as an indication of non-proportional hazards.Lastly, we explored model outlier using deviance residuals (±3 standard deviations) given their assumed normality.