Copy Number Variant Risk Scores Associated With Cognition, Psychopathology, and Brain Structure in Youths in the Philadelphia Neurodevelopmental Cohort

Key Points Question How do copy number variants (CNVs) combine with common genetic variants and environmental factors to help explain variability in cognition and psychopathology in a community sample? Findings In this community-based cohort study including 9498 youths in the Philadelphia Neurodevelopmental Cohort, elevated CNV risk scores were associated with lower cognitive ability and more subtly associated with both higher overall psychopathology and higher psychosis-spectrum symptoms. Statistical models of cognitive and psychopathological outcomes were significantly improved when CNV risk scores were combined with polygenic scores and quantitative measures of environmental stress. Meaning It is important to integrate rare genetic, common genetic, and environmental factors in investigations of clinically relevant developmental outcomes.

CNV associations with cognition and psychopathology, including environmental factors in multiancestry sample eFigure 6.Combined models of developmental outcomes and their joint associations with CNV scores, environmental factors, and common variant polygenic scores (PGSs) eFigure 7. Brain imaging (MRI) deviation from normative models and CNV risk scores eFigure 8. Sensitivity analysis excluding participants with known pathogenic CNVs eFigure 9. Sensitivity analysis with random effect for Illumina Infinium Beadchip subversion eFigure 10.Full distributions of individual brain imaging-based centile scores eFigure 11.Combined models of developmental outcomes and their joint associations with CNV scores, environmental factors, and common variant polygenic scores eTable 1. Demographic information for CNV samples after quality control eTable 2. Definitions of terms related to CNV risk scores eTable 3. Logistic regression models of categorical diagnoses eTable 4. Models of cognitive and psychopathological outcomes eTable 5. Exploratory analysis of interaction effects eTable 6. Known pathogenic CNVs identified in the PNC sample eTable 7. Association with CNV risk scores excluding known pathogenic CNVs eTable 8. Associations with CNV risk scores including X chromosome CNVs eTable 9. Associations with CNV risk scores including Affymetrix arrays eTable 10.No association between pHI CNV risk scores and demographic variables Participants were also administered the GOASSESS, a comprehensive computerized tool for structured evaluation of psychopathology domains based on the Kiddie Schedule for Affective Disorders and Schizophrenia (K-SADS) 1,5 .From GOASSESS items, overall psychopathology was calculated using a bifactor model 6 .This single continuous measure, sometimes referred to as the p-factor 7,8 , represents liability to overall psychopathology in youths 9 .Other domains of psychopathology were "Psychosis-Spectrum", "Externalizing", "Fear (phobias)" and "Mood-Anxiety" factors.We also report results from a correlated traits model, which includes the same specific symptom domains without extracting an orthogonal measure of overall psychopathology 6 .Categorical determinations were also made using the GOASSESS to identify individuals meeting criteria for ADHD, depression, anxiety, or psychosis spectrum 1,10 .For consistency with prior work, all psychopathology outcomes were age-normalized, using the residuals of a linear model with age as the predictor variable, and z-scored prior to further analysis.
As previously reported 1,3,6 , cognitive and psychopathological measures correlate non-trivially.Our factor analytic approach was designed to capture cognitive and behavioral domains supported by psychological theory and prior literature, rather than to constrain factors to be orthogonal.We note that it is common to investigate multiple cognitive or psychopathological outcome measures that are not expected to be independent of one another, which reflects the reality of high correlation among domains of cognitive function as well as psychiatric comorbidity 11 .False Discovery Rate (FDR) correction for multiple comparisons across multiple outcome measures that are not independent is expected to be conservative, i.e., this procedure is not expected to inflate the false discovery rate 12 .However, it is important to take these points into consideration when interpreting results; it is not justified to infer that reported associations with different cognitive and psychopathological outcomes are necessarily independent.As discussed above, an exception is the bifactor model of psychopathology domains, which models psychopathology factors orthogonally by accounting for an overall psychopathology factor that "absorbs" the inter-factor correlations.Thus for example, the overall psychopathology factor is orthogonal to the psychosis factor estimated from the bifactor model.This additional step is taken to justify the inference that, for example, an effect on psychosis is not attributable to an eMethods 3. CNV quality control and descriptive statistics CNV processing steps are illustrated in eFigure 1.We applied a series of quality control procedures to CNV output that have been used previously in similar published work 15,16 .Code is available at https://github.com/BGDlab/cnv-studyand https://github.com/MartineauJeanLouis/MIND-GENESPARALLELCNV.We inspected chiplevel attributes and determined if they met stringent quality cutoffs with respect to call rate (≥95%) log r ratio standard deviation (<0.35), B allele standard deviation (<0.08) and absolute value of wave factor (<0.05).All probe coordinates were updated from hg18 to hg19.
CNVs were called using the QuantiSNP 17 and PennCNV 18 algorithms.Both algorithms are freely available and detect CNVs using hidden Markov models (HMMs) of SNP genotype data, and have been validated for Illumina and Affymetrix arrays 19 .Implemented in MATLAB, QuantiSNP uses an objective Bayes approach to detect CNVs based on the LogR ratio and the B allele frequency for each SNP marker 17 .Implemented in PERL and C, PennCNV infers CNVs based on multiple sources of information including the log r ratio and B allele frequency at each SNP marker, the distance between neighboring SNPs, and an adjustment for "genomic waves" using a regression model GC content 20 .To minimize false discoveries, the overlap between CNVs detected from each algorithm was detected using CNVision 21 , and only CNVs with 70% overlap between the algorithms were included in subsequent analysis.
To minimize false discoveries, we also implemented a series of additional series of filtering steps.Following REF 16 , we initially excluded arrays with a suspiciously high number of CNVs detected (defined as ≥50 for low resolution arrays with <1 million probes, and ≥200 for high resolution arrays with ≥1 million probes) (see eFigure 2).CNVs were then filtered to ensure they were greater than 50,000 base pairs in size; they had a confidence score rating ≥ 30 from either QuantiSNP or PennCNV; they did not have >50% overlap with segmental duplications, the major histocompatibility complex (MHC), centromeric regions or telomeric regions.Finally, chips were excluded if, after CNV filtering, they had >30 CNVs or >8Mb total CNV burden.After exclusion criteria based on chip quality control and genetic relatedness within the sample (with a threshold of third degree relatives), information about CNVs was available for N=7,543 participants.
It is important to note that the PNC was genotyped in 15 batches, using ten different types of Affymetrix and Illumina arrays by the CHOP Center for Applied Genomics.Arrays have different probe coverage of given genomic regions, which has the potential to impact CNV detection sensitivity/specificity; this is the principal reason we limited our analysis large CNVs >50,000 bps which are more reliably detected across platforms.We elected to use all available probes in CNV detection to maximize discovery potential.To assess the impact of differences in array technology, main analyses were conducted using subjects genotyped with Illumina arrays (N=7,126, all genotyped with Illumina Infinium Beadchip subversions) with additional subjects included in sensitivity analyses (eMethods 8).We also conducted analyses where Illumina Chip Version was modelled as a random effect in statistical models of the association with outcome variables, as well as selected analyses within Chip Version, which showed convergent associations with outcomes (eMethods 8).
For families with more than one participant in the PNC, determined via kinship analysis in KING, we chose a single random family member to be included in subsequent analyses.For individuals with more than one chip, we selected the highest quality chip by log r standard deviation to be included in subsequent analyses.While systematic re-genotyping to validate CNVs calls was not possible for this study because sample aliquots for the PNC are largely exhausted, a subset of 1,515 participants had been genotyped more than once on different chips by the CHOP Center for Applied Genomics.For these participants, we compared the highest quality chip by log r standard deviation (which was used in all subsequent analyses) with the second highest quality chip by log r standard deviation, in terms of the reliability of the CNV risk scores as indexed by pHI (Pearson's r=0.83, t=58.018,df=1,513, p=4.44e-53).This suggests high reliability for CNV risk scores that are the basis of the results reported in this paper.

eMethods 4. CNV annotation and CNV risk scores
As defined in eTable 2, with respect to terminology we use the phrase "CNV risk score" when referring to quantitative measures of predicted risk based on annotating individual's CNVs.More specifically, we use the term CNV burden when referring to annotations of total size or number of genes encompassed by CNV, the term intolerance when referring to annotations based on intolerance to loss of function (i.e.pLI and LOEUF scores), and the term dosage sensitivity when referring to annotations based on deletion or duplication sensitivity (i.e., pHI and pTS).We use the phrase "known pathogenic" CNVs when referring to CNVs with known associations with cognition or psychopathology from prior literature.
We annotated CNVs using Gencode V19 (hg19) with ENSEMBL (https://grch37.ensembl.org/index.html) to derive CNV risk scores.Specifically, CNVs were scored on multiple measures of CNV burden, including total size (number of base pairs) and the total number of genes (the number of genes fully encompassed by the CNV) (eFigure 3C).Two measures of CNV intolerance were also calculated (eFigure 3D). 1) The probability of loss intolerance (pLI) is the probability that a gene is intolerant to a loss of function mutation, ranging from 0-1 with a relatively bimodal distribution.For a CNV, the cumulative pLI was calculated as the sum of pLI of the genes encompassed by the CNV.2) The loss of function observed/expected upper bound fraction (LOEUF) is a continuous measure ranging from 0 to 2, where values below 0.35 are typically suggestive of intolerance.For a CNV, the cumulative inverse LOEUF (i.e., 1/LOEUF) was calculated as the sum of the inverse LOEUF of the genes encompassed by the CNV 16 .Finally, measures of dosage sensitivity were calculated as the probability of haploinsuffiency (pHI) and probability of triplosensitivity (pTS), which quantify the probability (ranging from 0-1) that a gene is sensitive to copy number loss or copy number gain, respectively (eFigure 3E).The pHI and pTS used in the current analyses were derived from a previously reported statistical model based on 753,994 genomes, which mapped dosage sensitivity for all genes across the genome 22 .In contrast to pHI and LOEUF which were used to annotate both deletions and duplications (although the impact of deletions and duplication was modelled separately), pHI was used to annotate deleted genes while pTS was used to annotate duplicated genes.For a CNV deletion, the cumulative pHI was calculated as the sum of the pHI of genes encompassed by the CNV; for a CNV duplication, the cumulative pTS was calculated as the sum of the pTS of genes encompassed by the CNV.Participant-level CNV risk scores were summed across all of a participant's CNVs, with deletions and duplications considered separately.
As previously reported, pHI and pTS are generally highly correlated per gene (Pearson's r=0.69), though exceptions certainly occur.As would be expected due to the cumulative effect of multiple genes also being related to CNV size, the correlation is even higher when considering cumulative pHI and pTS within large CNVs -e.g.within 118 known pathogenic autosomal CNVs [23][24][25][26][27] (Pearson's r=0.93).Note that simply comparing pHI scores with pTS scores across CNVs in terms of their numeric values does not provide a meaningful indicator of predicted risk.While both pHI and pTS range from 0-1, numeric equivalence does not imply equal predicted risk.For example in REF 22 , different thresholds were used to establish genes that were categorically haploinsufficient (pHI>0.84)and triplosensitive (pTS > 0.993) with comparable effect sizes.This is one reason why, in terms of their effect on outcome measures, deletions and duplications were always modelled with separate terms in our statistical models.As described in the main text, and consistent with the fact that deletions are typically more damaging than duplications, the effect size of deletions was in general markedly higher than that of duplications across outcome measures.
Consistent with the fact that CNVs with large effect sizes on developmental outcomes are rare events, CNV risk scores are positively skewed (eFigures 3-4).It is important to note that modelling assumptions were not violated by skewness in an independent variable 30 , and because of this skewness we also explored categorical and logarithmic transformations of CNV scores.Our main results were convergent when using these nonlinear transformations.Moreover, model fits as assessed by AIC were superior for untransformed CNV scores (e.g., see Table 1), which are also consistent with prior related work 15,16 .We note that regression approaches designed to model rare events, such as Poisson regression and related approaches, would be applicable if modelling CNV scores as the outcome of interest (rather than a predictor variable) 31 .It will be instructive to further explore the impact of model assumptions on CNV associations in future work.

eMethods 5. PGS calculation, empirical determination of ancestry, and kinship
The approach to derive polygenic scores (PGSs) in the PNC has recently been described in detail 32 .Briefly, prior to imputation, PLINK 1.9 33 was used to remove single nucleotide polymorphisms (SNPs) with > 5% missingness, samples with more than 10% missingness, and samples in which genotyped sex was different from reported sex.Genotypes were phased (Eagle v.2.4) and imputed to the 1000 Genomes Project Other/Mixed GRCh37/hg19 reference panel (Phase 3 v.5) using Minimac 4 via the Michigan Imputation Server 45 .Following imputation, only polymorphic sites with imputation quality R2 ≥ 0.7 and MAF ≥ 0.01 were retained.Multidimensional scaling (MDS) of the imputed genotypes was conducted using KING v.2.2.4 34 to identify the top ten ancestry principal components (PCs) for each sample.These PCs were projected onto the 1000 Genomes Project PC space, and genetic ancestry was inferred using the e1071 35 support vector machines package in R. Based on these inferences, a European ancestry cohort was established and other ancestry groups were excluded from PGS analyses based on strong evidence of European bias in existing GWAS and PGS calculations 32 .A second round of unprojected MDS was then performed within the European ancestry group to produce ten PCs that were included as additional covariates in statistical models within the European ancestry cohort.KING was also used to identify all pairwise relationships out to third degree relatives based on estimated kinship coefficients and inferred IBD segments, and a single random family member was selected for subsequent analyses.
PRS-CS (PRS using SNP effect sizes under continuous shrinkage) 36 was used to infer posterior effect sizes of the SNPs based on specific GWAS summary statistics.Like other Bayesian-based techniques such as LDpred 37 and SBayesR 38 , PRS-CS has the substantial advantage of not relying on processes of linkage disequilibrium (LD) clumping and p-value thresholding 36 and these approaches have been shown to yield more predictive PGSs than those produced using older methodologies 39 .To ensure convergence of the underlying Gibbs sampler algorithm, we ran 25,000 Markov chain Monte Carlo (MCMC) iterations and designated the first 10,000 MCMC iterations as burn-in.The PRS-CS global shrinkage parameter was set to 0.01 when the discovery GWAS had a SNP sample size that was less than 200,000; otherwise, it was learned from the data using a fully Bayesian approach.Default settings were used for all other PRS-CS parameters.Raw PGS were produced from the posterior effects using PLINK 1.9, and PGSs were standardized using a z-transformation for use in all subsequent analyses.GWAS summary statistics for specific disorders were based on six GWAS, for attention deficit hyperactivity disorder (ADHD) 40 ; autism spectrum disorder (ASD) 41 ; major depressive disorder (MDD) 42 ; bipolar disorder (BPD) 43 ; schizophrenia (SCZ) PGC wave 3 44,45 ; and intelligence 46 .Sample sizes of the underlying GWAS referenced above were as follows for ADHD (n=20,183 cases, 35,191 controls), ASD (n=18,381 cases and 27,969 controls), MDD (n=246,363 cases and 561,190 controls), BPD (n=20,352 cases and 31,358 controls), SCZ (n= 69,369 cases and 236,642 controls), and intelligence (n = 269,867).The estimated proportion of variance due to common variants based on the above GWAS (h 2 SNP) were as follows for ADHD (h 2 SNP=0.22),ASD (h 2 SNP=0.12),MDD (h 2 SNP=0.089),BPD (h 2 SNP=0.17-0.23 depending on assumed prevalence [0.5-2.0%)]),SCZ (h 2 SNP=0.24), and intelligence (h 2 SNP=0.19).Note that for ADHD, we based our PGS on GWAS excluding participants from CHOP (who may have overlapped with the PNC) 40 .
Given prior evidence of an association between pathogenic CNVs and alterations in brain size 48 , and well-replicated associations between brain size and intelligence 49 and schizophrenia detectable with MRI 50 , we predicted that CNVs would impact volumetric structural brain MRI phenotypes in the PNC.The wide range of neuroanatomical alterations reported by structural neuroimaging studies of specific CNVs, however, prohibited a targeted neuroanatomical hypothesis 48 as well as a unique hypothesis about the directionality of alterations.Indeed, pathogenic CNVs are known to cause increases or decreases (depending on the specific variant) in total brain volume and regional brain structures [49][50][51][52] .We therefore explored a hypothesis based on deviation from normative developmental expectations in brain size, in any direction and across multiple imaging phenotypes.Using a recently described approach 55 , we deployed a normative model of brain development to quantify the imaging phenotypes of individual PNC participants relative to age-and sex-normalized growth charts.Cortical gray matter volume (GMV), subcortical gray matter volume (sGMV) and white matter volume (WMV) were classified as either within the normative range (10%-90%), infranormal (<10%) or supranormal (>90%).Categorically, individuals were classified as 'low deviation' if they were within the normative range for all phenotypes, and 'high deviation' if they were either infranormal or supranormal in any phenotype.The neuroimaging analysis focused on this omnibus test of a deviation from normative trajectories in order to maximize statistical power given the relatively smaller number of participants in the neuroimaging sample.The rationale for this procedure was also based also on its similarity to recent work showing categorical deviations from normative ranges in schizophrenia 56 .
Statistical tests were performed using logistic regression with imaging deviation as the outcome variable and with CNV-status as the predictor variable of interest.Covariates included ancestry, SES, trauma exposure, and a manual rating of image quality 47 .As normative imaging models were stratified for sex and age, these were not included as covariates in baseline models; however, the neuroimaging results reported in the main text and below were unchanged when including sex and age as additional covariates as well as for models without covariates included.
To supplement the results describing neuroimaging deviation associated with CNV risk scores derived from pHI and pTS as presented in the main text (Figure 3), we replicated the analysis using LOEUF-based CNV risk scores.After genetic and imaging quality control, there were 920 multi-ancestry participants with structural imaging.Of this cohort, 38 participants were characterized as having high CNV risk scores, defined as total deletion 1/LOEUF score > 2.94 (corresponding to LOEUF < 0.34, which has been suggested as the threshold for loss-of-function intolerance at the single gene level 16,28 ).Of these participants with high CNV risk scores, 58% were also categorized as high deviation from neuroimaging normative models, compared to only 40% of participants with low CNV risk scores (logistic regression β=0.68, p=0.047) (see eFigure 7A).We performed an additional supplemental analysis where we re-characterized CNV risk scores to incorporate a category of medium CNV risk scores defined as 0 < total 1/LOEUF score ≤ 2.94, which comprised 208 participants.The low CNV risk score group then comprised remaining participants with total 1/LOEUF score = 0.This stratification resulted in an association with neuroimaging-based brain deviation as follows: low CNV risk score group, 38% high neuroimaging deviation; medium CNV risk score group, 46% high neuroimaging deviation; high CNV risk score group, 58% high neuroimaging deviation (high risk score group, β=0.77, p=0.024; medium risk score group, β=0.34, p=0.036; eFigure 7B).Note that these β coefficients represent the effect on the log-odds ratio compared to the low CNV risk score reference group.See eFigure 10 for plots of the full distributions of the neuroimaging centile scores.
In sum, although the overall power to detect neuroimaging alterations are low because of the greatly decreased sample size in the neuroimaging cohort compared the whole PNC, results support the hypothesis of increased brain deviation from normative ranges associated with elevated CNV risk scores.

eMethods 7. Tests of CNV-PGS and CNV-environment interaction effects
It is of considerable theoretical and clinical interest whether CNV risk scores have additive or interactive effects with both common variant PGSs as well as environmental factors.However, compared to modelling additive effects, statistical power to reliably distinguish interaction effects is likely to be substantially lower 57 , a problem which may be compounded when considering relatively rare events such as the presence of CNVs with high risk scores.Without a clear rationale for targeted testing of specific interactions, we performed an exploratory analysis where we systematically tested for interactions between CNV risk scores (indexed by pHI and pTS) and environmental factors (SES and trauma exposures), and between CNV risk scores and the six PGSs queried in main analyses (ASD, ADHD, SCZ, BPD, MDD, and intelligence) for four cognitive outcomes (overall accuracy, executive complex cognition accuracy, social cognition accuracy, memory accuracy) and five psychopathological outcomes (mood, fear, externalizing, psychosis spectrum, and overall psychopathology from bifactor model described above).Not surprisingly given the combinatorial multiple testing burden (144 separate interaction effects), none of these effects was significant after FDR-correction for multiple comparisons.eTable 5 shows exploratory results (uncorrected p<0.05).We conclude that reliably distinguishing interactive effects will require larger sample sizes or samples enriched for genetic mutations -or alternatively, targeted testing of a smaller number of a priori specified hypotheses.This remains an important area for future work.

eMethods 8. Sensitivity analyses
As described in eMethods 4, we identified 130 total participants with 38 known pathogenic CNVs in the sample, and we excluded these participants in sensitivity analyses in order to ascertain whether main results were due entirely to the contribution of known pathogenic CNVs.eTable 6 provides a list of known pathogenic CNVs identified as well as the cumulative dosage sensitivity scores (pHI and pTS) associated with these CNVs.As shown in eTable 7 and eFigure 8, main results -including significant associations between CNV risk scores and cognitive and psychopathological outcomes -were in general robust to this sensitivity analysis.Specifically, the associations between overall accuracy and high-risk-score deletions indexed by pHI (β=-0.049,p=8.853e-06) and high-risk-score duplications indexed by pTS (β=-.047,p=1.26e-05); the association between overall psychopathology and pHI (β=0.031p=4.20e-03); and the association between psychosis spectrum symptoms and pHI (β=0.038,p=1.19e-03) remained statistically significant after FDR correction for multiple comparisons.As also shown in eFigure 8, reported associations between PGS and neurodevelopmental outcomes, and between environmental factors and neurodevelopmental outcomes, were also robust to this sensitivity analysis.We note that the effect size was reduced for the association between CNV risk scores and cognitive outcomes when excluding known pathogenic CNVs from the sample -for deletions but not for duplications.This may be due to diminished statistical power to detect the true effect size when excluding 130 participants who disproportionately have high CNV risk scores -given previous evidence that the effect size of the association between CNV risk scores and general intelligence is similar whether or not known pathogenic CNVs were included in models 15,16 .For the neuroimaging analysis, 17 of the 59 participants in the main analysis with high CNV risk scores (pHI score or pTS score >1) had known pathogenic CNVs.Unsurprisingly, for a sensitivity analyses excluding participants with pathogenic CNVs, the association between CNV risk scores and brain deviation from the normative model dropped to "trend level" in terms of statistical significance (β=0.62,p=0.052) due to the reduced sampled size, though the effect size was undiminished (in fact marginally increased) compared to the observed association when known pathogenic CNVs were included in the model (original β=0.56, original p=0.039), supporting the conclusion that the association between CNV risk scores and brain deviation was similar whether or not CNVs were "pathogenic" in the sense of having known associations with cognition or psychopathology from prior literature.Future work investigating neuroimaging-CNV associations will benefit from larger samples or samples enriched for high CNV risk scores.
As described in the main results, we also performed a series of additional sensitivity analyses to assess the associations between CNV risk scores, cognition and psychopathology.First, an analysis including X chromosomal CNVs yielded similar outcomes compared to our main result (eTable 8).This sensitivity analysis could not include pHI because that score was derived for autosomal genes only, but the results for intolerance scores LOEUF and pLI were convergent with our main analyses.Second, an analysis including Affymetrix chips yielded an additional 442 participants with similar outcomes to the main results (eTable 9).Third, we modelled Illumina Chip Version as a random effect in statistical models of the association with outcome variables, showing highly convergent results with main analyses (e.g.pHI association with Overall Accuracy, β=-0.11,p=5.86e-25; pTS association with Overall Accuracy, β=-0.05,p=1.17-06; pHI association with psychosis spectrum symptoms, β=-0.49,p=3.49e-05) and suggesting that Chip Version did not confound the association between high-risk-score CNVs and developmental outcomes (see eFigure 9).We also performed selected associations with outcome variables within Chip Version to further show that it was unlikely any biases between Chip Versions accounted for the reported results.Results were convergent with main analyses that combined data across Chip Version (e.g.eTable 12 shows results for the 5 most frequent Chip Versions, which together account for >90% of the data, for the association between CNV risk scores and total accuracy).Finally, main results were robust to sensitivity analyses that did not include covariates sex, self-identified race or ancestry PCs in models (eTable 11).All outcomes measures were age-normalized, but all results were also robust to including age as an additional regressor in statistical models.

eMethods 9. Analysis of family history
Mental illness in the parents of the present sample is likely to be related to outcome measures, as well as not only genetic but also environmental factors considered as predictor variables in the present study.For instance, parental history of psychiatric disorder is an indicator of genetic risk but is likely also to be associated with traumatic exposures.While information about parental history of mental illness per se is unfortunately not available in the PNC, a measure of family history of psychosis (in any first degree relative) was collected and described in prior work 58 .As would be predicted by prior work, this measure of family history was associated with traumatic exposures (β=0.11,p=5.96e-10) as well as neighborhood-level SES (β=0.03,p=0.00136).We reran models of environmental and genetic influence on outcomes including this measure of family history, and all results were robust to the inclusion of this measure in statistical models (eFigure 11).Moreover, there was an independent association between family history and overall psychopathology (β=0.06,p=6.39e-05).Speculatively, this supplemental analysis is consistent with the hypothesis that family history is related to but does not fully account for the environmental associations with neurodevelopmental outcomes reported in the present study.Limitations include the fact that this measure was focused on psychosis in particular and also that it was based on self-or caregiver-reported family history rather than independent clinical evaluation of family members.Further work, in studies designed to address these questions, is clearly required to investigate the combined influence of genetic and environmental effects on neurodevelopmental outcomes.

CNV risk scores
Quantitative measures of predicted risk based on annotating CNVs in terms of burden, intolerance or dosage sensitivity.

Burden
Measure of total size or number of genes encompassed by CNV.

Intolerance
Measure based on genes' intolerance to loss of function, for genes encompassed by a CNV (i.e.pLI and LOEUF scores).pLI Probability that a gene is intolerant to a loss of function mutation, where values closer to 1 are suggestive of intolerance 58 .

LOEUF
Loss of function observed/expected upper bound fraction, a measure ranging from 0 to 2, where values below 0.35 are typically suggestive of intolerance 55 .

Dosage sensitivity
Measure based on genes' deletion or duplication sensitivity, for genes encompassed by a CNV.The probability of haploinsufficiency (pHI) is measure of a gene's sensitivity to deletion, and the probability of triplosensitivity (pTS) is a measure of gene's sensitivity to duplication. 22

Pathogenic CNVs
CNVs with specific, known clinical associations with cognition or psychopathology based on prior literature.CNVs with high risk scores are not necessarily pathogenic.

eTable 3. Logistic regression models of categorical diagnoses
Logistic regression models of categorical diagnoses (psychosis spectrum, ADHD, anxiety, and depression) and their association with CNV risk scores.Note that this table reports association between CNV risk scores and categorical diagnoses rather than dimensional psychopathology outcomes, which are the focus of main analyses.Covariates included sex, self-identified race, and 10 ancestry PCs.This table was generated from the multi-ancestry sample of participant genotyped on Illumina arrays that met quality control criteria (n=7101), and p-values were corrected for 8 comparisons using the Benjamini-Hochberg false discovery rate (FDR).

eTable 4. Models of cognitive and psychopathological outcomes
Models of cognitive and psychopathological outcomes and their association with CNV risk scores (cumulative pHI for deletions and pTS for duplications), environmental factors (neighborhood SES and traumatic exposures), and polygenic scores (PGSs for ASD, ADHD, SCZ, BPD, MDD, and intelligence).AIC and adjusted r 2 are shown for models with increasing complexity, from left to right: demographic covariates only (sex, self-identified race and 10 ancestry principal components); CNV risk scores and demographic covariates; environmental factors and demographic covariates; CNV risk scores, environmental factors, and demographic covariates; CNV risk scores, environmental factors, and demographic covariates.This table was generated from the European sample of participant genotyped on Illumina arrays that met quality control criteria (n=4,482).All outcome measures were age-normalized.

eMethods 1 . 2 . 9 . 1 . 2 . 3 .
Cognitive and psychopathological phenotyping and factor analysis eMethods Quantification of environmental factors eMethods 3. CNV quality control and descriptive statistics eMethods 4. CNV annotation and CNV risk scores eMethods 5. PGS calculation, empirical determination of ancestry, and kinship eMethods 6. Neuroimaging eMethods 7. Tests of CNV-PGS and CNV-environment interaction effects eMethods 8. Sensitivity analyses eMethods Analysis of family history eFigure Schematic of CNV data processing pipeline eFigure Exclusion plot of chip frequency by number of CNV and markers eFigure Distribution of CNVs and CNV risk scores in the study sample eFigure 4. CNV frequency generally decreases in association with increases in CNV risk score eFigure 5.

eFigure 2 .
Exclusion plot of chip frequency by number of CNV and markersAs part of chip quality control illustrated in eFigure 1, chips were assessed for their number of CNVs relative to their number markers (probes).Following 16, we excluded arrays with a suspiciously high number of CNVs detected (defined as ≥50 for low resolution arrays with <1 million probes, and ≥200 for high resolution arrays with ≥1 million probes).Red dotted boxes indicate catchment areas for exclusion, meaning that chips inside red boxes were excluded.Distribution of CNVs and CNV risk scores in the study sample for deletions (red) and duplications (blue).All results are after the full quality control procedure described in eMethods 4, presented as histograms where the height of the lines corresponds to the number of CNVs (A-B) or to the number of participants (C-E).A) Total number of CNVs within 5 size bins.B) Total number of CNVs distributed across autosomes, from left to right: all CNVs; genic CNVs (encompassing at least one gene); CNVs with cumulative pHI or pTS score > 0; CNVs with cumulative pHI or pTS score > 1.Note this Panel contains the same information as Main Figure1Apresented in a different fashion.C) Participant burden scores: Left Panel, total size in Mb of all of a participant's CNVs; Right Panel, the number of genes fully encompassed by CNVs.D) Participant intolerance scores: Left Panel, the sum of genes' probability of loss intolerance (pLI) scores; Right Panel, the sum of the inverse of genes' loss of function observed/expected upper bound fraction (LOEUF) scores.E) Participant dosage sensitivity scores: the sum of genes' probability of haploinsuffiency (pHI) and probability of triplosensitivity (pTS) scores for deletions and duplications, respectively.

eFigure 4 .eFigure 5 .
CNV frequency generally decreases in association with increases in CNV risk scoreCNV frequency generally decreases in association with increases in CNV risk score (measured by cumulative probability of haploinsufficiency (pHI) in deletions and cumulative probability of triplosensitivity (pTS) in duplications.After applying quality control criteria outlined in eFigures 1-2 and eMethods 3, CNVs with 50% reciprocal overlap were considered members of the same CNV "cluster."A) CNV deletions shown with cumulative pHI of genes in each CNV cluster.B) CNV duplications shown with cumulative pTS of genes in each CNV cluster.In clusters, pHI or pTS was estimated as the mean of all CNVs in the cluster.Clusters with a mean pHI or pTS greater than 10 are labeled in the hexplot with their chromosome number and the mean start and stop (Mb) of the CNVs within the cluster.The x-axis shows the total number of participants with each CNV.As expected, most CNVs were unique (1-member clusters without 50% reciprocal overlap to another CNV) and had CNV risk scores close to 0. Note that clustering based on reciprocal overlap was performed for visualization purposes for this figure, but in all other reported analyses each participant's CNV risk score was calculated independently.Participantlevel CNV risk scores were summed across all of a participant's deletions and duplications respectively.Note that this Figure shows CNV risk scores including known pathogenic CNVs, which were included in main results but excluded in sensitivity analysis (eMethods 8).CNV associations with cognition and psychopathology, including environmental factors in multiancestry sample Analogous to Main Figure1but including environmental factors in addition to CNV risk scores.CNV risk scores were quantified as the cumulative probability of haploinsufficiency (pHI, a measure of sensitivity to deletion) or probability of triplosensitivity (pTS, a measure of sensitivity to duplication).Dot plots show the effect size (standardized beta coefficient) for the association of pHI and pTS with eight cognitive outcomes (Top Panel) including speed and accuracy scores for specific (e.g."Memory", "Social Cognition") and global ("Overall") cognitive measures."Slow speed" is summarized from items requiring deliberation, while "fast speed" is summarized from items requiring rapid decisions.For psychopathology outcomes, dot plots of the effect size (standardized beta coefficient) for pHI and pTS generated via bifactor models (Middle Panel) and correlated traits factor models (Bottom Panel).All outcome measures were age normalized, and additional covariates included self-identified race, sex and 10 ancestry principal components in all models.This figure was generated based on the multi-ancestry sample of participant genotyped on Illumina arrays that met quality control criteria (n=7,101), and p-values were corrected for 68 comparisons using the Benjamini-Hochberg false discovery rate (FDR).

eFigure 6 .
Combined models of developmental outcomes and their joint associations with CNV scores, environmental factors, and common variant polygenic scores (PGSs) Analogous to Main Figure2but showing additional terms from the model and omitting depiction of 95% confidence intervals for clarity.Each point in the dot plots indicate the value of a given predictor variable's effect size in models of cognition or psychopathology for models of cognition (Top Panel) and models of psychopathology (Bottom Panel).CNV risk scores were quantified by deletion cumulative probability of haploinsufficiency (pHI) and duplication cumulative probability of triplosensitivity (pTS).Environmental factors comprised neighborhood socioeconomic status (SES) and trauma exposure.PGSs were included for ADHD, ASD, BPD, 'g' (general intelligence), MDD and SCZ.All outcome measures were age normalized, and additional covariates included self-identified race, sex and 10 ancestry principal components in all models.This analysis was conducted in the European ancestry sample genotyped with Illumina arrays that met quality control criteria (n=4,482) and p-values were corrected for 90 comparisons using the Benjamini-Hochberg false discovery rate (FDR).

eFigure 7 .eFigure 8 . 9 .eFigure 10 .eFigure 11 .
Brain imaging (MRI) deviation from normative models and CNV risk scoresAnalyses of the association between brain imaging (MRI) deviation from normative models and CNV risk scores using LEOUF-based risk score and incorporating a medium-risk-score category (compare with Main Figure3).A) Proportion of individuals with high-risk-score CNVs (cumulative 1/LOEUF score > 2.85) who are categorized as low brain deviation (2nd to 9th decile in all imaging phenotypes) versus high brain deviation (1st or 10th decile in at least one imaging phenotype).B) Addition of medium-risk-score category (0 < cumulative 1/LOEUF ≤ 2.85).LOEUF, loss-of-function observed/expected upper bound fraction.Sensitivity analysis excluding participants with known pathogenic CNVs Results of sensitivity analysis excluding 130 total participants with 38 known pathogenic CNVs from the sample.The re-analysis of CNV risk score (Top Panel) was conducted in the multiancestry sample (n=7,101), while the analysis of polygenic scores (PGS, middle panel) was conducted in the European Ancestry sample (n=4,482).For consistency with Main Figure2, the analysis of environmental factors (Bottom Panel) was also conducted in the European Ancestry sample (n=4,482).FDR-correction was across 90 total comparisons.All outcome measures were age normalized, and additional covariates included self-identified race, sex and 10 ancestry principal components in all models.eFigure Sensitivity analysis with random effect for Illumina Infinium Beadchip subversion Results of sensitivity analysis incorporating a random effect for Illumina Infinium Beadchip subversion into statistical models (multi-ancestry sample size, n=7,101; FDR correction across 32 associations).Results were highly convergent with main analyses (compare with Main Figure1).All outcome measures were age normalized, and additional covariates included self-identified race, sex and 10 ancestry principal components in all models.Full distributions of individual brain imaging-based centile scores Illustration of full distributions of individual brain imaging based centile scores (age-and sexstratified).Violin plots of individual brain imaging based centile scores: A) cortical gray matter volume (GMV), B) subcortical gray matter volume, C) cerebral white matter volume (WMV).CNV risk scores are determined as for Main Figure3(total pHI > 1 or total pTS > 1).Combined models of developmental outcomes and their joint associations with CNV scores, environmental factors, and common variant polygenic scores Analogous to eFigure 6 but in models that also include a measure of family of history of psychosis.Each point in the dot plots indicate the value of a given predictor variable's effect size in models of cognition or psychopathology for models of cognition (Top Panel) and models of psychopathology (Bottom Panel).CNV risk score is quantified by deletion cumulative probability of haploinsufficiency (pHI) and duplication cumulative probability of triplosensitivity (pTS).Environmental factors comprised neighborhood socioeconomic status (SES) and trauma exposure.PGSs were included for ADHD, ASD, BPD, 'g' (general intelligence), MDD and SCZ.All outcome measures were age normalized, and additional covariates included self-identified race, sex and 10 ancestry principal components in all models.
This analysis was conducted in the European ancestry sample genotyped with Illumina arrays that met quality control criteria (n=4,482) and p-values were corrected for 99 comparisons using the Benjamini-Hochberg false discovery rate (FDR).Demographic information for CNV samples after quality control © 2022 Alexander-Bloch A et al.JAMA Psychiatry.eTable 2. Definitions of terms related to CNV risk scores

.0759,-0.0322) 1.31E-06
Exploratory analysis of interaction effectsInteraction effects were explored between CNV risk scores and environmental factors (neighborhood-level SES and trauma exposure) and between CNV risk scores and common variant risk indexed by polygenic scores.Statistics are shown for a subset of nominal effects (uncorrected p<0.05).No effects remain significant after FDR-correction for 144 multiple comparison.All models also included additive effects for pHI, pTS, PGS and environmental factors, self-reported race, sex and 10 ancestry principal components.This analysis was run on the European ancestry sample genotyped on Illumina chips that met quality control criteria (n=4,482) Known pathogenic CNVs identified in the PNC sample Known pathogenic CNVs defined by systematic literature review identified in the sample.See eMethods 4; these CNVs were excluded in sensitivity analysis.Sorted by frequency (number of participants, n, shown in leftmost column).Cumulative probability of haploinsufficiency (pHI) and triplosensitivity (pTS) are defined as the sum of pHI and pTS of genes encompassed by CNVs for deletions or duplications, respectively.Start and stop are in hg19 coordinates.Association with CNV risk scores excluding known pathogenic CNVs Association between CNV risk scores and overall cognitive accuracy, sorted by lowest AIC, showing the effect on models of excluding pathogenic CNVs.Overall accuracy scores were agenormalized, and additional covariates included sex, self-identified race and 10 ancestry principal components.This table was generated from the multi-ancestry sample of participant genotyped on Illumina arrays that met quality control criteria, and p-values were corrected for 16 comparisons using the Benjamini-Hochberg false discovery rate (FDR).This table is analogous to Main Table1but showing the effect of excluding pathogenic CNVs.To facilitate comparison, the original results including pathogenic CNVs are also appended below the double bar, with rows marked by asterisks (*), while the results for the sensitivity analysis excluding known pathogenic CNVs appear above the double bar.Associations with CNV risk scores including X chromosome CNVs Associations between CNV risk scores and overall cognitive accuracy including X chromosome CNVs.Rows are sorted by lowest AIC.Overall accuracy scores were age-normalized, and additional covariates included sex, self-identified race and 10 ancestry principal components.Analogous to Main Table1but including CNVs on the X chromosome.This table was generated from the multi-ancestry sample of participant genotyped on Illumina arrays that met quality control criteria, and p-values were corrected for 16 comparisons using the Benjamini-Hochberg false discovery rate (FDR).Associations with CNV risk scores including Affymetrix arrays Association between CNV risk scores and overall cognitive accuracy including Affymetrix arrays.With Affymetrix array, total n=7,543.Rows sorted by lowest AIC.Overall accuracy scores were age-normalized, and additional covariates included sex, self-identified race and 10 ancestry principal components.Analogous to Main Table1but including Affymetrix arrays for the multi-ancestry sample of participant that met quality control criteria.p-valueswerecorrectedfor16 comparisons using the Benjamini-Hochberg false discovery rate (FDR).Akaike information criterion; CI, confidence interval; pLI, probability of loss intolerance; 1/LOEUF, inverse of loss-of-function observed/expected upper bound fraction; pHI, probability of haploinsufficiency; pTS, probability of triplosensitivity eTable 11.Associations without demographic covariates included in statistical models Association between CNV risk scores and overall cognitive accuracy without demographic covariates included in statistical models.Overall accuracy scores were age-normalized.Analogous to Main Table1but without demographic covariates included in statistical models.This table was generated from the multi-ancestry sample of participant genotyped on Illumina arrays that met quality control criteria, and p-values were corrected for 16 comparisons using the Benjamini-Hochberg false discovery rate (FDR)., copy number variant; AIC, Akaike information criterion; CI, confidence interval; pLI, probability of loss intolerance; 1/LOEUF, inverse of loss-of-function observed/expected upper bound fraction; pHI, probability of haploinsufficiency; pTS, probability of triplosensitivity eTable 12. Impact of Illumina beadchip subversion Assessing the impact of Illumina Chip Version on the association between CNV risk score and developmental outcomes.Table shows results for overall accuracy within chip for the 5 most frequent chips which together account for >90% of the sample.Overall accuracy scores were age-normalized, and additional covariates included sex, self-identified race and 10 ancestry principal components.Analogous to Main Table1but showing models of data from each chip separately.The row "all chip" shows data from the full sample for comparison.This table was generated from the multi-ancestry sample of participant genotyped on Illumina arrays that met quality control criteria.P-values shown in this table were not corrected for multiple comparisons.