Population Performance and Individual Agreement of Coronary Artery Disease Polygenic Risk Scores

Study Population

The primary analysis was conducted using data from the All of Us (AOU) Research Program, a National Institutes of Health-funded biobank composed of adult volunteers across the United States who have given written consent for analysis of their deidentified EHR and genetic data.^16,17 At the time of analysis, 237,568 participants had available whole genome sequencing. For each participant, genetic sex was imputed using PLINK’s check_sex function. A binary CAD phenotype was assigned in the presence of at least one of the following codes: 410, 411, 412, 413, 41, V45.81 (ICD-9), I21, I22, I24, Z95.1, Z98.61, I20.0 (ICD-10). Age was calculated using participant birth year and the most recent data release cutoff date. Individuals without EHR data were excluded. All analyses were performed on deidentified data from All of Us Controlled Tier Dataset v7 by authorized researchers. Similar approaches were used for the Penn Medicine Biobank (PMBB) and UCLA ATLAS Precision Health Biobank (UCLA) (eMethods).^18–20

Selection and Creation of Polygenic Risk Scores

Polygenic risk scores for CAD were selected from the Polygenic Score (PGS) Catalog (eMethods).⁹ In addition, two novel CAD PRSs were created by applying Bayesian approaches – PRS-CSx-auto and LDPred2-auto – to summary statistics from a multi-population meta-analysis of genome-wide association studies (GWAS), comprising 380,508 individuals with and 1,836,455 individuals without CAD (eFigure 1); these scores were named “PGS_LDP2Auto” and “PGS_prscsx,” after the methods used to construct them (eMethods).^21–23 None of the selected or created scores utilized AOU, PMBB, or UCLA data for training.

Calculation of Polygenic Risk Scores

Polygenic scores for each individual in the AOU, PMBB, and UCLA cohorts were calculated using pgsc_calc, which combines individual genotypes with the existing PRS weights.^9,24 Scores were calculated for each sample using whole-genome sequence data from AOU and imputed genotypes from PMBB and UCLA, and weights from all PGS Catalog scores meeting described criteria, as well as PGS_LDP2Auto and PGS_prscsx. Match rate, or the percent of variants in a weight file that were successfully identified in the genomic data, was recorded.

Heterogeneity of allele frequencies and linkage disequilibrium patterns across populations can influence raw PRS distributions, limiting interpretation and generalizability.^8,25 To minimize these effects, scores were adjusted using a principal component analysis (PCA)-based method to normalize both mean and variance to the 1000 Genomes + HGDP reference panel (eMethods).²⁶ All downstream analysis utilized PCA-normalized scores. These values were accordingly translated into risk percentiles based on a standard normal distribution.

PRS performance in diverse populations depends on representation of genetically similar individuals in the GWASs used to construct the PRS.^15,27 CAD PRSs vary widely in their use of diverse GWAS data. To characterize the impact of genetic background on PRS performance, population subgroup-stratified sensitivity analyses were conducted. Individuals were assigned to one of six population groups based on their genetic similarity to populations included in the 1000 Genomes + HGDP reference panel: African (AFR), admixed American (AMR), East Asian (EAS), European (EUR), Middle Eastern (MID), and Central/South Asian (CSA). Sensitivity analyses were conducted in the AFR and EUR groups. Population group assignment was implemented using pgsc_calc (eMethods).²⁸

Population-Level Assessment: Identification of equivalently performing polygenic risk scores

The association between each PRS and prevalent CAD was assessed using generalized linear regression models with a logit link, including age and sex as covariates. These covariates were included to minimize ascertainment bias. Genetic principal components were not included in the primary analysis model on the premise that they are already accounted for by PCA-based score normalization. Inclusion of PCs, as well as exclusion of all covariates were both evaluated as sensitivity analyses. The approach to risk model evaluation considered relevant published reporting standards^11,12: risk score distribution was assessed visually, and PRS effect size was quantified with odds ratios; discrimination and calibration were visualized with calibration plots and ROCs, and quantified using Brier score and area under the receiver operator curve (AUROC). Bayesian methods were used to perform between-model comparisons and identify a set of PRSs with statistically and practically equivalent performances relative to the ‘best’ performing score for each metric (eMethods).

Individual-Level Polygenic Risk Score Assessment

United States regulatory guidelines state that new clinical tests should be evaluated by assessing the agreement between the results of new test and that of a reference standard.²⁹ Scores were therefore evaluated through this lens of agreement to determine how consistently scores with comparable population-level performance assigned individual risk.

First, descriptive statistics were calculated to understand the individual-level distributions of genetic risk estimations provided by equivalently performing scores. Each individual’s average risk percentile was calculated, defined as the mean of percentiles assigned by each score meeting the Region of Practical Equivalence (ROPE) of 0.02, as ‘practical equivalence’ criteria. In addition, we calculated the standard deviation and the ratio of the standard deviation to the mean (i.e., the coefficient of variation (CV)) across scores per individual. Distributions of mean risk percentiles and standard deviations across biobank populations were plotted, and we calculated the population-level median of the individual-level average scores, standard deviations, and CV. Bootstrapped 95% confidence intervals were obtained using the “boot” and “simpleboot” R packages.

Agreement across scores was next assessed using inter-rater reliability. If one considers PRSs with equivalent performance as independent, equally qualified testers of the same phenomenon (in this case, genetic liability for CAD), one can apply tests of inter-rater reliability to determine whether scores that are equivalent are also practically interchangeable. The intraclass correlation coefficient (ICC), which is used as a metric of the consistency of quantitative measurements made by different raters measuring the same quantity, was used to compare PRS percentile as a continuous variable.³⁰

Since many tested and proposed clinical applications of PRS use percentile cutoffs to define ‘high-risk’ groups, agreement of binary assignment above/below these thresholds was also assessed by evaluating inter-rater agreement of the categorical classification. This was quantitatively assessed using Light’s Kappa (eMethods).³¹

Pearson’s correlation coefficient (r) of individual risk estimates was calculated for all pairs of scores in the primary cohort.

Statistical Analysis

All analyses were performed using R (version 4.3.0 in PMBB, Version 4.3.1 in AOU and UCLA). Statistically equivalent scores were defined as those with a less than 95% probability of a positive difference of both Brier score and AUROC. Practically equivalent scores were those with a >95% probability of model performance being ‘practically equivalent’ using a ROPE estimate.^32,33 An a priori practical effect size was set at 0.02, and in sensitivity analyses ranges of 0.01 and 0.005 were considered. For individual-level analysis, Light’s Kappa (κ) and ICC were used to compare categorical and continuous agreement of PRS percentiles. Although both ICC and κ are descriptive measures, guidelines exist to aid interpretation of intermediate values. For both, 0 indicates agreement that is no better than random chance, and 1 indicates perfect agreement. A Kappa between 0.41 and 0.60 can be interpreted as fair agreement; an ICC coefficient between 0.50 and 0.75 is generally interpreted as a fair to good relationship. ^34–37 Methods are described in further detail in eMaterials.

Clinical Characteristics

The AOU study population comprised 17,589 (10.3%) participants with and 153,506 (89.7%) participants without CAD (eResults). 97,265 (63.4%) of those without CAD and 7,682 (43.7%) of those with CAD were female. 100,493 (58.7%) of participants were most genetically similar to a European reference population, 35,590 (20.8%) to an African reference population, 29,805 (17.4%) to an Admixed American reference population, and the remaining to Central/South Asian, East Asian, and Middle Eastern reference populations (ST1, eFigure 2). Those with CAD were older on average (69.3 years versus 54.9) (ST1). Both PMBB and UCLA featured a higher rate of prevalent disease. Of 41,193 PMBB participants, 9,215 (22.4%) had CAD (ST2). Of 53,092 UCLA participants, 12,513 (23.5%) had CAD (ST3). These cohorts are representative of national demographic patterns and are well-powered for PRS assessment and replication, allowing for reliable and generalizable population- and individual-level assessment.

Population-Level PRS Performance

We first aimed to identify a set of polygenic scores that predicted prevalent CAD in a statistically equivalent way based on population-based performance prediction metrics (Figure 1). We obtained total of 50 CAD polygenic scores, consisting of 48 from the PGS catalog and two additional scores constructed specifically for this study (PGS_LDP2Auto and PGS_prscsx, Methods, ST4). The majority of genetic variants included in each score were present for AOU, PMBB, and UCLA participants (ST6-8). We observed that two scores had negative associations with prevalent CAD and were subsequently excluded from downstream analysis (ST9, eFigure3, eResults). When added to a model of prevalent disease that included age and sex as covariates in AOU, the remaining 48 scores were significantly associated with CAD (p<0.05) with odds ratios that ranged from 1.10 (PGS000059) to 1.46 (PGS_LDP2Auto) per standard deviation increase in each PRS (eFigure3, ST9), confirming the expected statistical association of these models.

• Download figure
• Open in new tab

Figure 1: Overview of Approach.

50 polygenic risk scores were calculated for all individuals in the All of Us Research Program, Penn Medicine Biobank, and UCLA ATLAS. In each cohort, the scores with superior population-level performance when included in a model of prevalent disease were identified, and performance metrics were compared with all other scores. Using this to define a list of practically equivalent scores, we considered agreement between risk percentiles assigned individuals by different, ‘equivalent’ scores.

We next evaluated population-based performance measures for a range of prediction metrics. We calculated Brier scores and AUROC for all PRS scores (Methods). We found that the model for the PGS_LDP2Auto score had the best calibration measured by Brier score (0.0825, 95% Credible Interval 0.0829-0.0827), whereas the model including PGS003725 had the best discrimination as measured by AUROC (0.777, 95% Credible Interval 0.777, 0.778) (ST10). To determine the set of score that were practically equivalent, we applied a prespecified primary region of practical equivalence effect size of 0.02 (Methods). All 48 scores had practically equivalent population-level performance (Figure 2). Sensitivity analyses that applied more stringent practical effect size margins of 0.01 and 0.005 resulted in 19 and 5 equivalent scores, respectively (ST10). Two scores met our prespecified definition for statistically equivalent performance across both measures of calibration and discrimination. Analyses in PMBB and UCLA revealed similar results, identifying 33 and 48 ROPE 0.02 equivalent scores, respectively (eResults, ST11-14, eFigure 4-7). This result indicates that most CAD PRS are indistinguishable using population-level performance metrics.

• Download figure
• Open in new tab

• Download figure
• Open in new tab

Figure 2: Differences in Polygenic Risk Score AUROCs and Brier Scores in AOU.

The mean of the posterior distribution of the difference between the AUROC/Brier Score of each combination of scores is plotted. The score with the ‘best’ model metric is denoted with bold text. Plus signs indicate scores with practically equivalent metrics (ROPE 0.02). Asterisks indicate scores with statistically equivalent metrics. Scores are ordered by year of publication, with PGS_LDPred2Auto and PGS_prscsx (newest) on the right.

Finally, we performed sensitivity analyses in AOU to assess the robustness of the broad population performances of these measures. When stratifying by genetic similarity to continental reference populations or removing model covariates, we observed broadly similarly hierarchies of score performance, but with unique sets of scores meeting equivalence criteria (eResults, ST15-21, eFigure8-15). This result indicates that the best-performing score at the population-level can differ within as well as between cohorts.

Agreement of Individual PRS Risk Assessments

Having identified a set of CAD PRS that have equivalent population-level metrics, we next evaluated the consistency of individual-level risk estimates. We defined an individual’s average risk percentile as the mean of the PRS percentile risk estimates provided by the 48 scores. Across AOU participants, the median standard deviation of risk percentiles for each individual was 22.94 (22.92, 22.96) and the median coefficient of variation for each individual was 0.504 (0.503, 0.505), providing scale-dependent and -independent evidence of the variability of individual-level estimates (Figure 3, eResults). We observed similar in analyses of PMBB and UCLA ATLAS and across sensitivity analyses (eFigure16-21).

• Download figure
• Open in new tab

Figure 3: Within-Person Score Concordance in AOU.

Concordance of individual score percentiles across all scores meeting practically equivalent ROPE 0.02 criteria in All of Us (48 scores). A) Mean individual risk percentile; median 48.38(48.26, 48.49). B) Standard Deviation of the mean individual risk percentile; median 22.94 (22.92, 22.96). C) Coefficient of Variation; median 0.5039 (0.5028, 0.5052).

Next we examined the distribution of scores in the top / bottom 20^th percentiles of risk across all 48 scores. We observed that 84% of individuals had at least one risk estimate above the top 20^th risk percentile, and one in the bottom 20^th risk percentile (ST22). To help illustrate this individual-level variability, we plotted individual-level risk percentiles for five randomly selected participants across all practically equivalent scores (Figure 4) and designed an interactive web app that allows users to explore CAD PRS percentiles for randomly selected individuals from the 1000 Genomes + HGDP reference population (sabramow.github.io/PRS_Var/). These results empirically highlight how an individual’s estimated risk can differ depending on which score was used.

• Download figure
• Open in new tab

Figure 4: Risk Predictions for Five Randomly Selected Individuals.

48 practically equivalent (defined as ROPE 0.02 in AOU) polygenic risk scores are plotted on the x-axis, in order of year of publication (newest on the right). On the y-axis is polygenic risk score percentile. Color corresponds to a unique randomly selected All of Us participant. Each participant’s risk percentile according to each score is plotted.

We next quantified interrater reliability across scores. Among the 48 practically equivalent scores, ICC was 0.351 (0.349, 0.352) which, according to common interpretation frameworks can be considered poor (Table 1a).^35–37 ICCs for the sets of 5 and 19 scores that met practical equivalence criteria with a ROPE of 0.005 and 0.01 were 0.734 (0.732, 0.736) and 0.555 (0.551, 0.558) respectively. The five ROPE 0.005 equivalent scores included two pairs of scores from the same studies (PGS_prscsx and PGS_LDP2Auto, and PGS003725 and PGS003726). The ICC between the two scores which performed equivalently at the population-level (PGS003725 and PGS_LDP2Auto) of 0.646 (0.643, 0.649) can be considered moderate to good.^35–37 Because thresholding PRSs to assign individuals into “high-risk” categories is commonly applied, Light’s Kappa was then used to assess whether scores agreed in their risk categorization. Across a range of possible thresholds, performance was almost universally poor but tended to improve with more stringent equivalence thresholds (Table 1b, ST22, eResults).

Next, we repeated our approach in our other settings. In each of these analyses, we used the corresponding sets of equivalent and practically equivalent scores determined for each analysis cohort. We found that metrics and patterns of variability and risk percentile congruence were similar we replicated our analysis pipeline in population-stratified AOU analysis, and in PMBB and UCLA (eResults, ST23-24, ST26-27). These results confirmed that low inter-rater reliability was not exclusive to our primary cohort.

Noting that some pairs of scores seemed to produce more highly correlated individual estimates than others, we considered whether improvements in performance associated with improvements in agreement. For our primary cohort, we produced a heatmap of the correlation coefficients between all score pairs, and also plotted this coefficient by the difference between the performance metrics of the two scores. We observed that correlation between risk estimates provided by pairs of scores varied widely. Two main trends emerged: correlation was highest among pairs of scores derived from the same studies, and also tended to be higher among scores that performed similarly overall and better compared to the rest of the rest of the scores (eFigure 22-23). This trend suggests that using a common data set to generate multiple scores can inflate agreement, and that inter-score agreement may be improving with time and score performance.

Limitations

This study should be interpreted in the context of its limitations. First, although we selected scores that were available from the PGS Catalog at the time of analysis, PRSs for CAD are continuously being developed. Second, primary conclusions regarding relative population-level performance of these scores are specific to tested models of prevalent CAD (including population-stratified, and with and without age and sex as covariates). Prevalent disease was chosen as a primary outcome to maximize number of cases and therefore power to detect differences in score performance. Other predictive modeling scenarios have been proposed⁵³ and our results do not preclude the possibility that scores can be more conclusively differentiated from these other models. Third, minor differences in sequencing and imputation methodology (which in turn affects match percentage) and phenotype definitions used in UCLA, PMBB, and AOU may account for some of the variability between score performance metrics obtained from each biobank.⁵⁴ Fourth, tests of inter-rater reliability assume that raters (in our case, PRSs) provide independent estimates of risk, but our tested PRSs cannot be considered truly independent raters because in most cases the scores build on a core set of CAD GWAS data. Although a lack of true independence does limit the formal validity of our inter-rater reliability metrics, the effect should be to bias the test statistic towards a higher degree of agreement/correlation, and thus only strengthens our interpretation of substandard inter-rater reliability.

Conclusions and Relevance

When tested across three diverse biobanks, CAD PRSs that performed similarly at the population level demonstrated highly variable individual-level estimates of risk. Approaches to clinical implementation of CAD PRSs need to consider the potential for incongruent individual risk estimates from otherwise indistinguishable scores.