Elucidating mechanisms of genetic cross-disease associations at the PROCR vascular disease locus

Many individual genetic risk loci have been associated with multiple common human diseases. However, the molecular basis of this pleiotropy often remains unclear. We present an integrative approach to reveal the molecular mechanism underlying the PROCR locus, associated with lower coronary artery disease (CAD) risk but higher venous thromboembolism (VTE) risk. We identify PROCR-p.Ser219Gly as the likely causal variant at the locus and protein C as a causal factor. Using genetic analyses, human recall-by-genotype and in vitro experimentation, we demonstrate that PROCR-219Gly increases plasma levels of (activated) protein C through endothelial protein C receptor (EPCR) ectodomain shedding in endothelial cells, attenuating leukocyte–endothelial cell adhesion and vascular inflammation. We also associate PROCR-219Gly with an increased pro-thrombotic state via coagulation factor VII, a ligand of EPCR. Our study, which links PROCR-219Gly to CAD through anti-inflammatory mechanisms and to VTE through pro-thrombotic mechanisms, provides a framework to reveal the mechanisms underlying similar cross-phenotype associations. Many individual genetic risk loci associate with multiple diseases, but the molecular basis of these loci often remains unclear. Here, the authors provide a framework to reveal the genetic cross-disease associations at the PROCR vascular disease locus.

Many individual genetic risk loci have been associated with multiple common human diseases. However, the molecular basis of this pleiotropy often remains unclear. We present an integrative approach to reveal the molecular mechanism underlying the PROCR locus, associated with lower coronary artery disease (CAD) risk but higher venous thromboembolism (VTE) risk. We identify PROCR-p.Ser219Gly as the likely causal variant at the locus and protein C as a causal factor. Using genetic analyses, human recall-by-genotype and in vitro experimentation, we demonstrate that PROCR-219Gly increases plasma levels of (activated) protein C through endothelial protein C receptor (EPCR) ectodomain shedding in endothelial cells, attenuating leukocyte-endothelial cell adhesion and vascular inflammation. We also associate PROCR-219Gly with an increased pro-thrombotic state via coagulation factor VII, a ligand of EPCR. Our study, which links PROCR-219Gly to CAD through anti-inflammatory mechanisms and to VTE through pro-thrombotic mechanisms, provides a framework to reveal the mechanisms underlying similar cross-phenotype associations. https://doi.org/10.1038/s41467-022-28729-3

OPEN
A full list of author affiliations appears at the end of the paper. G enome-wide association studies (GWAS) have revealed widespread pleiotropy of disease-associated genetic variants. A recent study of cross-phenotype genetic association data in the UK Biobank has shown that 96% of traitassociated variants (minor allele frequency (MAF) ≥ 1%) are associated with more than one ICD-10 code, with some showing associations with more than 50 codes 1 . The vast majority of these pleiotropic variants were found to impact the risk of multiple diseases in a directionally consistent manner, but 1.9% of loci (excluding the major histocompatibility complex) showed evidence of both higher and lower risk effects attributable to the same allele 1 . One such example is rs9349379 A > G, a wellcharacterized regulatory variant at the PHACTR1-EDN1 locus, which is associated with a higher risk of coronary artery disease but a lower risk of four other vascular diseases including migraine headache and hypertension 2 .
Another example of a pleiotropic variant is p.Ser219Gly (rs867186 A > G) in the PROCR gene, which encodes the endothelial protein C receptor (EPCR), a key regulator of the protein C (PC) pathway. The minor G allele of this variant has been shown to correlate with a lower risk of CAD 3,4 and myocardial infarction 5 , but a higher risk of venous thromboembolism (VTE) [6][7][8] . This pattern of opposing associations seems paradoxical because several conventional cardiovascular risk factors (e.g., measures of adiposity) show directionally concordant associations for CAD and VTE 9 . Further, GWAS of cardiovascular intermediate traits have reported associations between rs867186-G and components of the coagulation cascade, including higher plasma levels of PC 10 and coagulation factor VII 11,12 . However, the causal relevance of these intermediate traits to cardiovascular diseases remains uncertain.
The thrombomodulin-protein C pathway serves as a key mediator of the cross-talk between coagulation and inflammatory processes. It comprises molecular components that can respond to a range of pathophysiological environments in different vascular beds [13][14][15] . At the vascular endothelium, thrombomodulin binds to thrombin, directly inhibiting its clotting and cell activation potential and converting PC to activated PC (APC) (reviewed in 15,16 ). The activation of PC by the thrombin-thrombomodulin (TM) complex is markedly enhanced when PC is presented by EPCR 17 , a type I transmembrane protein that is mainly expressed on the endothelium of large blood vessels. 18,19 Once APC dissociates from EPCR, it binds to protein S to inactivate the coagulation factors Va and VIIIa, thereby inhibiting further thrombin generation. In addition, APC promotes fibrinolysis by decreasing the levels of plasminogen activator inhibitor type 1 (PAI-1), and reduces inflammation by inhibiting the production of tumor necrosis factor (TNF)-α and interleukin(IL)-1β (reviewed in 15,16 ).
A soluble form of EPCR (sEPCR) is present in plasma, which is generated by ectodomain shedding of EPCR from the endothelium. Plasma sEPCR levels in healthy individuals display a bimodal distribution, with higher levels being associated with one of the four frequent haplotypes at the PROCR locus 8,[20][21][22][23] . This haplotype (denoted A3 or H3) is tagged by the minor allele of the p.S219G variant. Functional studies showed that the variant results in increased shedding of EPCR from the endothelial surface by rendering the receptor more sensitive to cleavage by metalloprotease 21 and by forming an alternatively spliced, truncated transcript 24 . The shedding is effectively regulated by TNF-α and IL-1β 25 . sEPCR retains its ability to bind both PC and APC but does not enhance PC activation 26,27 . However, the precise molecular mechanism underlying the PROCR-p.S219G functional variant and its influence on the cardiovascular intermediate phenotypes that may mediate the risk of CAD and VTE is incompletely understood.
In this study, we aim (1) to systematically assess the association of the PROCR-p.S219G variant with a range of cardiometabolic outcomes and relevant risk factors; (2) to evaluate causality of individual components of the protein C pathway on cardiovascular diseases; and (3) to help uncover the molecular and cellular chain-of-events that connect the PROCR-219Gly allele to a lower risk of CAD but a higher risk of VTE. The results of our integrative epidemiological and functional analyses ( Fig. 1) reveal new insights underlying the PROCR association locus for arterial and venous diseases and have potential implications for the development of therapeutic strategies targeting components of the protein C pathway.

Results
Association of PROCR-p.S219G with cardiovascular diseases and risk factors. To search for associations of PROCR-p.S219G with a broad range of human diseases, we conducted a phenomewide association analysis across 1402 electronic health recordderived ICD-codes from the UK Biobank. The association of PROCR-p.S219G with each of these codes was tested using SAIGE 28 , a generalized mixed model association test that accounts for case-control imbalance and sample relatedness, as implemented in PheWeb (Methods). The data implicated diseases of the circulatory system, e.g., phlebitis/thrombophlebitis (Phe-WAS code 451; P = 4.2 × 10 −8 ) and coronary atherosclerosis (PheWAS code 411.4; P = 2.9 × 10 −5 ) (Fig. 2a).
We  binding activity of the PC and APC SOMAmer reagents to a range of relevant proteins, specifically, PC, APC, sEPCR, thrombin, FV, FVIIa, protein S and thrombomodulin (Methods). We confirmed that the APC SOMAmers bind the proteins in a specific manner. However, we found that the PC SOMAmer binds to both the zymogenic and activated form of protein C (Supplementary Table 1), which may contribute to the observed difference in the magnitude of effect sizes observed for the immuno-and SomaScan assays (Fig. 2c). Additionally, we confirmed that the presence of relevant binding partners of PC and APC do not interfere with SOMAmer binding (Supplementary Table 1).
Identification of shared genetic etiology at the PROCR locus.
Despite data showing associations of the rs867186 variant at the PROCR locus with CAD and VTE, it has been uncertain whether they reflect a shared causal variant and mechanism. To address this, we performed statistical colocalization analyses. We applied a Bayesian algorithm, Hypothesis Prioritization in multi-trait Colocalization (HyPrColoc) 35 , which allows for the assessment of colocalization across multiple complex traits simultaneously (Methods). We found colocalization of the genetic association data of CAD and DVT as well as factor VII, PC and APC levels at the PROCR locus, with a posterior probability of colocalization of 99.37% (Fig. 3a). The variant rs867186 was found to be the likely causal variant at the locus explaining 99.31% of the posterior probability (Fig. 3b). Thus, these data provide support for a common genetic mechanism underlying the PROCR locus.
Causal evaluation of protein C in arterial and venous diseases.
The association data suggest that genetic variants at the PROCR locus influence PC and APC abundance, FVII activity and susceptibility to CAD and DVT. However, these data do not necessarily imply that these molecular traits have a causal relationship with the disease phenotypes. To help define this relationship, we conducted Mendelian randomization (MR) analyses, using genetic variants as instrumental variables to avoid confounding and reverse causation 36 . We constructed a multi-allelic genetic score to estimate the causal associations between the putative risk factors and cardiovascular outcomes (Methods). The score comprised of approximately independent (r 2 < 0.1) SNPs at the PROCR region with P value ≤5 × 10 −8 (Methods; Supplementary Table 2). Our data showed that every geneticallypredicted increment (per 1 SD) in PC levels is associated with a lower risk of CAD  Supplementary Fig. 1). We also performed these analyses with APC, resulting in similar effect sizes and association P values (Supplementary Table 3). Findings were robust to the use of a range of different MR approaches, i.e., inverse-variance weighting (IVW) method, median-based methods (simple and weighted) and MR-Egger regression (Methods). We conducted further sensitivity analyses, confirming the validity of our results to potential violations of the MR assumptions (Methods; Supplementary Fig. 2). We applied reverse MR to evaluate evidence for causal effects in the reverse direction by modeling disease phenotypes as the exposure and PC or APC level as the outcome using genome-wide significant predictors of disease (Methods). These analyses revealed no reverse causality of CAD or DVT/VTE on the levels of PC or APC (Table 1; Supplementary Table 3). Taken together, these analyses provide evidence of causal relationships between the levels of zymogenic and activated protein C and CAD and VTE outcomes, in opposite directions.
Validation of 'focal' phenotype associations in the protein C pathway. To determine the molecular and cellular effects of rs867186, the causal variant at the PROCR locus, we performed a recall-by-genotype study. Such recall-studies allow for the strict control of experimental conditions (e.g., identical processing of blood samples), statistical efficiency (i.e., balanced recruitment based on genotype independent of MAF) and deep-phenotypic characterization of the collected samples (e.g., in vitro challenge experiments) (reviewed in 37 ). From a genotyped panel of healthy volunteers, we selected 52 individuals stratified by rs867186 genotype and matched for sex and age (Methods). In these individuals, we measured four biomarkers in plasma representing focal phenotypes that describe the functional state of the protein C pathway, i.e., levels of protein C (inferred from a chromogenic assay measuring PC activation in response to an exogenous stimulus), APC, sEPCR and thrombin-antithrombin (TAT) complex (Methods). We found that the minor (G) allele of rs867186 associated with higher plasma levels of sEPCR (β = 1.10, P = 3.29 × 10 −22 ) (Fig. 4a). This finding is consistent with previous reports 8,23,[38][39][40][41][42] . We also found that the G allele associated with elevated PC activity, a marker for PC levels (β = 0.59, P = 2.05 × 10 −4 ) (Fig. 4a). These data are concordant with and complementary to the data that we report from the epidemiological studies above (Fig. 2c), in that the chromogenic assay used here is not affected by potential binding-affinity effects of proteinaltering variants often detected in protein-binding assays. We did not observe genotypic effects on plasma levels of either APC (β = −0.14, P = 0.42) or TAT complex (β = −0.13, P = 0.46) (Fig. 4a). Together, these data provide a direct comparison of the genotypic effect of the PROCR causal variant on the functional PC pathway, and independent confirmation that both sEPCR and PC are higher in carriers of PROCR-rs867186-G.
Quantification of EPCR expression and shedding in endothelial cells. We next aimed to identify the direct upstream molecular determinants of elevated sEPCR levels due to the rs867186-G genotype. Using transcriptomic data across 27 mature hematopoietic cell types from the BLUEPRINT Blood Atlas 43 , we found that PROCR is highly expressed in human umbilical vein endothelial cells (HUVECs) and modestly expressed in macrophages, but not expressed in any other cell type analyzed ( Supplementary Fig. 3). Consistent with these data, in flow cytometric analyses, we determined high expression of membrane-bound EPCR in HUVECs (Supplementary Figs. 4,5). We detected 1.9-fold lower levels of EPCR in untreated HUVECs obtained from homozygotes of the rs867186-G-allele compared to homozygotes of the A-allele (P = 0.0051) (Fig. 4b). In HUVECs treated with phorbol 12myristate 13-acetate (PMA), a potent agent to enhance ectodomain shedding, we found lower levels of EPCR compared to HUVECs treated with vehicle control in both homozygote groups, i.e., 4.1-fold (P = 0.0124) and 5.3-fold (P = 0.0046) for rs867186-G-allele and -A-allele homozygotes, respectively (Fig. 4c). Taken together, these findings are consistent with increased EPCR shedding from endothelial cells in carriers of PROCR-rs867186-G. We also performed flow cytometric analyses in a monocytic cell line (U937 cells), for which we observed modest levels of EPCR expression (Supplementary Fig. 6). We then differentiated these cells into macrophagelike cells using PMA and showed a~30% reduction in EPCR expression relative to the undifferentiated cells ( Supplementary  Fig. 6). Finally, we sought to determine whether rs867186-G also affects EPCR shedding on primary neutrophils and monocytes purified from the individuals from our recall-study. However, we did not detect the presence of EPCR on the surface of either of these cell types by flow cytometry (Supplementary Figs. 7, 8).
Effect of sEPCR on leukocyte-endothelial cell adhesion.
Leukocyte-endothelial cell adhesion is a critical step in atherosclerosis that triggers vascular infiltration of monocytes and subsequently leads to microvascular inflammation 44 . Previous in vitro studies have highlighted EPCR as a potential modulator of the leukocyte-endothelial cell interaction. Specifically, sEPCR is a binding partner for the integrin macrophage-1 antigen (Mac-1) 45 , which is expressed on the surface of activated leukocytes and is a key mediator of adhesion to the endothelium. Consequently, we investigated the effects of increasing concentrations of recombinant human sEPCR on leukocyte-endothelial cell adhesion using an in vitro static adhesion model. In brief, U937 cells were differentiated into macrophage-like cells using PMA and then dispensed onto a monolayer of TNF-α-activated HUVECs (Methods). Cell adhesion events were quantified following incubation with increasing concentrations of anti-Mac-1 antibody (positive control) and recombinant sEPCR (Methods). We found that increasing concentrations of anti-Mac-1 antibody (compared to an IgG control; P = 0.029), but not sEPCR (P > 0.05) (Supplementary Fig. 9) led to a reduction of adhesion events.
Effect of APC on cell adhesion molecule expression and leukocyte-endothelial cell adhesion. Inflammatory cytokines such as TNF-α activate endothelial cells by increasing the expression of cellular adhesion molecules. We investigated whether human plasma-derived APC mitigates the TNF-α-associated increase in gene expression of cellular adhesion molecules, such as intercellular adhesion molecule 1 (ICAM-1) and vascular cell adhesion molecule 1 (VCAM-1). Using reverse transcription quantitative PCR (RT-qPCR), we showed that increasing concentrations of APC attenuate the TNF-α-induced increase in ICAM1 mRNA levels in both HUVECs (P = 0.0003) and human coronary artery endothelial cells (HCAECs) (P = 0.0009) but not b a VCAM1 and CCL2 mRNA levels (P > 0.05) (Fig. 5a). Notably, APC exposure also reduced PROCR gene expression in HUVECs (P = 0.0066) (Fig. 5a). Finally, in static leukocyte-endothelial cell adhesion assays, we showed that APC treatment leads to a reduction of leukocyte-endothelial cell adhesion events in HUVECs (P = 0.0011) and HCAECs (P = 0.0246) (Fig. 5b). Together, these data suggest that in carriers of the PROCR-219Gly genotype, who exhibit elevated APC levels (as measured on the SomaScan platform), the lower genetic susceptibility to arterial disease may be due to a reduced number of leukocyte-endothelial cell adhesion events at sites of vascular inflammation.

Discussion
Elucidation of the molecular basis of cross-disease associations affords a major opportunity to advance understanding of disease etiology. Leveraging recent advances in population biobanks, statistical genomics and translational epidemiology, we illustrate an integrative, multi-modal approach to address this challenge. We applied this approach to two vascular diseases oppositely associated with the missense variant p.S219G (rs867186) in PROCR. We showed that PROCR-219Gly protects against CAD but increases susceptibility to VTE through distinct chains of molecular events, summarized in Fig. 6.
The data from our study show that PROCR-219Gly leads to a perturbed PC pathway, which acts focally to modulate the circulating levels of APC and has downstream effects on the biological mechanisms of associations with VTE and CAD.
We found that PROCR-219Gly is associated with higher circulating plasma sEPCR and lower EPCR levels on endothelial cells (Fig. 4), which is consistent with an increase in membrane shedding of EPCR and confirms findings from previous studies 8,21,38,42 . As only the membrane-bound form of EPCR is capable of activating PC 26 , we anticipated that this reduction in EPCR would result in increased PC but reduced APC levels. Accordingly, we (Fig. 4a) and others 46,47 have observed higher plasma PC levels in PROCR-219Gly carriers. However, in our phenome-scan, we observed an unexpected increase in APC levels as measured on the SomaScan platform (Fig. 2c). We performed extensive testing to confirm the specificity of the APC SOMAmer (Supplementary Table 1), indicating that this finding is not due to cross-reactivity with PC or other coagulation factors.
Since the primary driver of PC activation in vivo is the TM complex 48 , not EPCR, the higher levels of APC observed in PROCR-219Gly carriers may be due to an upregulation of TM activity in these individuals. In this scenario, an increase in APC would represent a homeostatic mechanism attempting to compensate for the increased thrombotic potential in PROCR-219Gly carriers and may be indicative of an acquired APC resistance. Indeed, APC resistance in the absence of Factor V Leiden is estimated to be prevalent in the general population (10-15%) 49 . Alternatively, given both PC and APC bind to sEPCR with the same affinity as EPCR 26 , it is conceivable that the higher sEPCR levels observed in PROCR-219Gly carriers may serve to stabilize and prolong the presence of PC/APC in the circulation. This would be particularly salient for APC given its short half-life of 15 minutes 50 . Furthermore, when bound to sEPCR, APC is unable to inactivate FV or FVIII 26,51 . Therefore, the sequestering of APC by sEPCR in PROCR-219Gly carriers may inhibit the anticoagulant activity of APC, resulting in APC resistance and increased thrombotic potential in these individuals. However, we did not observe statistically significant associations of PROCR-219Gly with FV or FVIII levels in our phenome scan (Fig. 2c).
In addition to its well-known role as an anticoagulant, APC has also been shown to function as a cytoprotective and antiinflammatory agent via the protease-activated receptor 1 (PAR-1) 52 . Indeed, TNF-treated endothelial cells exposed to APC have reduced mRNA and surface protein levels of key intercellular adhesion molecules, such as intercellular adhesion molecule-1 (ICAM-1) (Fig. 5a), which regulate the adhesion of leukocytes to the endothelium [53][54][55][56][57] . By performing static adhesion assays, we provided evidence that increasing concentrations of APC reduce the adhesion of activated monocytes to endothelial cells (Fig. 5b). Based on these data, we propose that the APC/PAR-1 signaling pathway may be critical in protecting against CAD by reducing leukocyte-endothelium adhesion and vascular inflammation in the coronary arteries of PROCR-219Gly carriers (Fig. 6).
Our findings have implications for therapeutic strategies targeting the PC pathway for vascular diseases. Despite early positive clinical data that proposed the use of Drotrecogin alfa (Xigris ® ; a recombinant form of APC) as a therapeutic intervention for sepsis and septic shock, the medicine was withdrawn due to the lack of replication in subsequent trials and its associated risk of bleeding 58 . However, APC has since emerged as a potential candidate for the treatment of stroke. Clinical trials are currently ongoing to test in patients with acute ischemic stroke the efficacy of 3K3A-APC, a recombinant form of APC that lacks its anticoagulant activity but retains its PAR-1 cell-signaling activities 59 . Preliminary results showed that patients receiving 3K3A-APC had reduced hemorrhage volume and hemorrhage incidence on day 30 following initial drug infusion, relative to a placebo group 60 . The findings from these clinical studies are consistent with the results of our wide-angled genetic association scan (Fig. 2), and provide a rationale to define and catalogue the disease relationships of pleiotropic variants on a genome-wide level to inform the development of new medicines. The presented phenome-scan also showed a significant association of PROCR-219Gly with higher plasma levels of FVII (Fig. 2c). Recently, FVII has been identified as a ligand for EPCR and shown to bind EPCR with the same affinity as PC 61 . Although EPCR does not affect the activation of FVII, the interaction of EPCR with FVII leads to the clearance of FVII/ FVIIa from the circulation through endocytosis 61 . Our data are consistent with this observation, as PROCR-219Gly is not only associated with higher plasma levels of FVII but also reduced levels of EPCR (Fig. 4b, c). Thus, the reduced availability of EPCR could directly contribute to the reduced internalization of FVII/ FVIIa and increased accumulation in the circulation, which in turn may increase thrombotic potential. Further research is necessary to confirm that PROCR-219Gly is indeed associated with reduced FVII/FVIIa internalization, for example, through performing endocytosis assays in genotype-specific or CRISPR/ Cas9-edited endothelial cell lines. Nevertheless, this proposed mechanism is consistent with the suggestive genetic association signals observed in PROCR-219Gly carriers for higher levels of D-dimer 62 , a marker of blood clot degradation, and shorter prothrombin time 63 (P = 3.70 × 10 −6 and P = 9.98 × 10 −8 , respectively). The association with shorter prothrombin time was replicated at genome-wide significance in the Japanese population (P = 5.64 × 10 −24 ) 64 .
We acknowledge that our study has limitations. First, many hemostatic factors have short half-lives 50,65 , which presents a technical challenge for studies seeking to quantify accurately these markers. Second, contrary to the statistically significant association between rs867186 and APC levels as measured using the SomaScan assay, in our recall-study, we found no evidence of an association. This is likely due to the difference in statistical power between the two experiments, with sample sizes of 3,301 and 52 individuals, respectively ( Supplementary Fig. 10). Indeed, previous studies that aimed to ascertain an association between rs867186 and APC have been hampered by relatively small sample sizes, yielding mixed findings 8,39,40 . Replication in independent large cohorts is needed. Third, EPCR is expressed on the  surface of platelets 66 . As the blood processing in our study likely resulted in platelet-poor as opposed to platelet-free plasma, it is possible that some of the sEPCR signal observed in our recall experiment (Fig. 4a) may have originated from platelet-associated EPCR. However, PROCR mRNA levels are very low in human platelets ( Supplementary Fig. 3), suggesting that any plateletassociated signal is likely to be negligible. Fourth, further studies are required to elucidate the complex interactions between EPCR and its ligands PC (APC) and FVII (FVIIa), as well as the downstream consequences of these interactions on hemostasis. Several aspects of our approach are generalizable to the study of other cross-disease associations (Fig. 1). First, the availability of large, disease-agnostic population biobanks with linked genomic, molecular phenotype and health record data, such as UK Biobank, provides an opportunity to systematically characterize the molecular underpinnings of health outcomes. Second, publicly available bioinformatics tools, including SAIGE 28 and PhenoScanner 67 , allow for the mining of these data and the generation of specific hypotheses about the underlying biological mechanisms at individual genetic association loci. Third, the use of freely available software for statistical colocalization 35 and Mendelian randomization 68 analyses enables evaluation of the extent to which associated phenotypes share the same causal variant and the causal relationship of molecular biomarkers with a disease outcome. Fourth, the increasing availability of volunteers in bioresources (e.g., UK National Institute for Health Research BioResource) who have agreed to participate in biomedical studies on the basis of their genetic and/or phenotypic characteristics enables targeted mechanistic studies tailored to specific hypotheses. This includes recall-by-genotype studies, which afford an efficient approach to detailed phenotyping that can be applied to different study designs, biological samples and experimental techniques 37 .
Taken together, our study provides new insights into the role of the PC pathway in arterial and venous diseases. We demonstrate that the combination of population biobank data and advanced statistical methods can help identify causal biomarkers and pathways, and that recall-by-genotype is a powerful experimental approach that can yield informative mechanistic insights. Overall, our study provides a framework for mapping molecular mechanisms that underlie cross-phenotype associations.

Methods
PROCR-rs867186 phenome-scan. The phenome-wide association scan of PROCR-p.S219G (rs867186) across electronic health record-derived ICD-codes from the UK Biobank was conducted using PheWeb v1.1.17 (http:// pheweb.sph.umich.edu/SAIGE-UKB/variant/20:33764554-A-G). To assess the effects of PROCR-rs867186 genotype on cardiovascular intermediate traits and outcomes, we collated data from the latest available genome-wide association studies using PhenoScanner v2, a database of human genotype-phenotype associations 67 . To allow comparative analyses, we considered data from individuals of European ancestry where possible. We focused our analyses on cardiometabolic traits and outcomes; thus, not all genome-wide significant associations are reported. The following association statistics were retrieved: stroke outcomes from the MEGASTROKE consortium 69 ; venous thromboembolism outcomes from the INVENT consortium 7 or UK Biobank; hypertension and aortic aneurysm from UK Biobank; coronary artery disease from van der Harst et al. 4 ; type 2 diabetes from Mahajan et al. 70 ; blood lipids from the Global Lipids Genetics consortium 71 ; hematological traits from Astle et al. 72 ; and plasma proteins of the coagulation cascade and protein C pathway from the ARIC study 10 , CHARGE consortium 34,73,74 , Sun et al. 30 and Suhre et al. 75 . To enable a comparison of the magnitude of the effect sizes, we conducted analyses with standardized units of measurement for each quantitative trait. Supplementary Data 1 provides further details of all data used in our analyses. The data availability section provides further information on the results of the complete data query. Supplementary Data 2 shows an overview of the associations of PROCR-rs867186 with protein levels measured by the SomaScan platform.
Determination of equilibrium binding constants. Equilibrium binding constants (K d values) of modified aptamers were determined by filter binding assay. K d values of modified aptamers were measured in SB18T buffer (40 mM Hepes pH 7.5, 102 mM NaCl, 5 mM KCl, 5 mM MgCl 2 , 0.01% Tween-20). Modified aptamers were 5ʹ end-labeled using T4 polynucleotide kinase (New England Biolabs) and γ-[ 32 P]ATP (Perkin-Elmer). Commercially available proteins to be used in the filter binding assay (protein C, APC, sEPCR, thrombin, factor V, factor VIIa, protein S and thrombomodulin) were biotinylated by covalent coupling of EZ-Link NHS-PEG4 -Biotin (Thermo Scientific) following the manufacturer's protocol. Briefly, proteins were combined with a 10-fold molar excess of EZ-Link NHS-PEG4 -Biotin in SB18T buffer and incubated at room temperature for 30 min. Free biotin was removed via YM-3 filtration (Millipore). Following biotinylation, protein concentrations were determined using a Micro BCA Protein Assay kit (Thermo Fisher). Radiolabeled aptamers (~20,000 CPM, 0.03 nM) were mixed with biotinylated proteins at concentrations ranging from 10 −7 to 10 −12 M and incubated at 37°C for 40 min. Bound complexes were partitioned on MyOne streptavidin beads (Invitrogen) and captured on Durapore filter plates (EMD Millipore). The fraction of bound aptamer was quantified with a phosphorimager (Typhoon FLA 9500, GE) and data were analyzed in ImageQuant TL (GE). To determine binding affinity, data were fit using the equation: y = (max − min)(Protein)/(K d + Protein) + min.
Competition binding assays. Competition binding assays were performed to test whether sEPCR and thrombomodulin interfere with the SOMAmer reagent 2961-1_2 binding to protein C or whether sEPCR, protein S and factor VIIa interfere with the binding of SOMAmer reagents 2961-1_2, 3758-63_3 and 3758-68_3 to APC. These experiments were performed by pre-incubating equal volumes of biotinylated protein C (80 nM) or biotinylated APC (48 nM or 80 nM) with competitor protein concentrations ranging from 10 -5 to 10 -10 M at 37°C for 30 min in SB18T buffer in the presence of 2 µM polyanionic competitor Z-block (a 30-mer modified DNA sequence, [AC(BndU) 2 ] 7 AC) 23 to allow protein complexes to form. Following the 30-min incubation, the reaction was diluted in half with radiolabeled SOMAmer reagent (20,000-60,000 CPM, 0.03 nM) and returned to 37°C for an additional 30 min. Bound complexes were partitioned on MyOne streptavidin beads and captured on Durapore filter plates. The amount of bound aptamer was quantified with a phosphorimager and data were analyzed in ImageQuant. The fraction of SOMAmer bound at each competitor concentration was normalized to the signal in the no competitor control well.
Multi-trait colocalization. We performed colocalization analysis at the PROCR gene locus (chr20: 31,916,110-35,505,723 bp; hg19), as defined based on recombination rates 3 . Details about the GWAS summary statistics used for this analysis are provided in Supplementary Data 1. Variants with both imputation (INFO)score <0.7 and MAF < 0.01, or variants with INFO-score <0.3 and MAF > 0.01 were removed. The remaining 4,264 SNPs shared across each of the datasets were aligned to the DNA plus-strand (hg19) prior to colocalization analyses. We used a Bayesian algorithm, implemented in the Hypothesis Prioritization in multi-trait Colocalization (HyPrColoc) v1.0 method 35 , to perform colocalization across all traits simultaneously. HyPrColoc extends the established colocalization methodology 76 by approximating the true posterior probability of colocalization with the posterior probability of colocalization at a single causal variant and a small number of related hypotheses 35 . If all traits do not share a causal variant, HyPr-Coloc employs a novel branch-and-bound selection algorithm to identify subsets of traits that colocalize at distinct causal variants at the locus. We used uniform priors as primary analysis and set strong bounds for the regional and alignment probabilities as default, i.e., the P R * (regional probability threshold) = P A * (alignment probability threshold) = 0.8, so that the algorithm identified a cluster of traits only if the posterior probability of full colocalization (PPFC) = P R P A > 0.64. We also performed sensitivity analysis with non-uniform priors to assess the choice of priors, which used a conservative variant-level prior structure with P = 1 × 10 −4 (prior probability of a SNP being associated with one trait) and P c = 1 − γ = 0.02 (P c is the conditional colocalization prior that a SNP is causal for an additional trait given that it is causal for one trait), i.e., 1 in 500,000 variants is expected to be causal for two traits.
Selection of instrumental variables for MR analysis. We obtained regional association statistics at the PROCR region for plasma PC levels from the ARIC study and plasma APC levels from the INTERVAL study to assess the causal effects of PC (APC) on cardiovascular outcomes. Details about the GWAS data on cardiovascular outcomes are provided in Supplementary Data 1. To select genetic variants as instrumental variables for PC levels, we first removed SNPs with MAF < 0.01 and INFO-score <0.8. Next, we performed LD clumping to obtain approximately independent SNPs. In brief, the algorithm groups SNPs in LD (r 2 ≥ 0.1 in 4,994 participants from the INTERVAL study 77 ) within ±1 MB of an index SNP (i.e., SNPs with association P value ≤5 × 10 −8 ). The algorithm tests all index SNPs, beginning with the smallest P value and only allowing each SNP to appear in one clump. Thus, the final output contains the most significant proteinassociated SNPs for each LD-based clump across the genomic region. An overview of the instrumental variables is provided in Supplementary Table 2. This analysis was performed using PLINK v1.90 78 .
NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-022-28729-3 ARTICLE Mendelian Randomization analyses. We used two-sample Mendelian randomization (MR) 68,79 to estimate the causal associations between PC and cardiovascular outcomes. The MR approach was based on the following assumptions: (i) the genetic variants used as instrumental variables are associated with PC levels; (ii) the genetic variants are not associated with any confounders of the exposure-outcome relationship; and (iii) the genetic variants are associated with the outcome only through changes in PC levels, i.e., a lack of horizontal pleiotropy. We applied the inverse-variance weighting (IVW) method in a multiplicative random-effect metaanalysis framework 79 , MR median-based method (simple and weighted) 80 , MR-Egger regression 81 and MR-PRESSO 82 to estimate the causal effects. We also performed several sensitivity analyses to assess the robustness of our results to potential violations of the MR assumptions, given these analyses have different assumptions for validity: (i) heterogeneity was estimated using the MR-IVWQstatistic; (ii) horizontal pleiotropy was estimated using MR-Egger's intercept; (iii) the median-based methods have greater robustness to individual genetic variants with strongly outlying causal estimates compared with the inverse-variance weighted and MR-Egger methods; and (iv) influential outlier instrumental variables due to pleiotropy were identified using MR-PRESSO and (v) MR-Steiger filtering 83 was used to eliminate spurious results due to reverse causation. We also applied reverse MR 84 to evaluate evidence for causal effects in the reverse direction by modeling disease phenotypes as the exposure and PC or APC levels as the outcome. Instrumental variants for phenotypes of interest (i.e., CAD, DVT/VTE) were selected from their original GWAS data (Supplementary Data 1). The effects of these GWAS SNPs on PC levels were derived from Sun et al. 30 . The power and strength of the instrumental variables was assessed using the variance explained (R 2 ) and F-statistics (F = β 2 /se 2 ) 85  Recall-by-genotype study. The study was approved by the Leicester Central Research Ethics Committee and Health Research Authority (Reference: 17/EM/ 0028). Healthy volunteers were recruited from the NIHR Cambridge BioResource with informed consent. Participants who were older than 18 years of age and of European ancestry were selected based on PROCR-rs867186 genotype and homozygosity of the major allele for both F5-R506Q (rs6025; Factor V Leiden) and F2-G20210A (rs1799963; Factor II). Participants across the three rs867186 genotype groups were matched at the end of the study with respect to sex and age (within 10 years). Study participants were excluded that had a diagnosis of (i) a chronic disease; (ii) hypertension (or history of consistently high blood pressure readings, i.e., >140/90 mmHg); and/or (iii) hypercholesterolemia (or history of consistently high cholesterol levels, i.e., >6 mmol/l). Participants agreed to fast and abstain from caffeinated drinks for at least four hours prior to the study visit and to not receive any vasoactive medication for up to seven days prior to procedures. The study design was informed by a power calculation (Supplementary Fig. 10).
Assessment of baseline characteristics of recall-by-genotype study participants. Participants reported past medical conditions, demographic factors (e.g., ethnicity) and lifestyle factors (e.g., smoking and alcohol consumption). Height and weight/body fat were measured using a stadiometer and bioelectrical impedance (i.e., Tanita scale), respectively. Blood pressure and heart rate were assessed in onemin intervals using a validated, automated device while seated and again after 3-5 min standing. All measurements were done in triplicate using the same arm. An overview of the characteristics of the study participants is provided in Supplementary Table 4. These characteristics are presented as mean and standard deviation or percentage. Continuous and categorical variables between homozygous groups were compared using the two-sample t-test and chi-square test, respectively.
Blood sample collection and processing. A total of 46 ml of peripheral blood was collected from each donor in our recall-study using a 21 gauge needle unless clinically contraindicated. We collected blood in two S-Monovette 7.5-ml K3 EDTA tubes and two S-Monovette 10-ml sodium citrate 3.2% (1:10) 9NC tubes (Sarstedt). Samples were immediately centrifuged at 4°C and 1000 × g for 15 min. Multiple aliquots of the top phases were stored at −80°C within 30 min of blood draw. A full blood count for all donors was obtained from blood collected in a S-Monovette 1.2-ml K3 EDTA tube using a Sysmex Hematological analyzer.
Quantification of plasma biomarkers. Samples were thawed at 37°C for 15 min, mixed and then centrifuged at room temperature and 3000 × g for 10 min immediately prior to assay. Soluble EPCR levels were determined using an Asserachrom sEPCR kit (00264; Diagnostica Stago); thrombin/antithrombin III complex levels using an Enzygnost TAT micro immunoassay (OWMG15; Siemens Healthcare Diagnostics Limited); APC levels using an Activated Protein C assay kit (CSB-E09909H; Cusabio Biotech); and PC levels using a HemosIL Protein C chromogenic assay (0020300500; Instrumentation Laboratory). All assays were performed according to the manufacturer's instructions. Samples were analyzed in random order and laboratory staff were blinded to genotype status. Participants with biomarker levels (or activity levels) 3 standard deviations above or below the population mean were removed.
Quantification of EPCR levels on genotyped HUVECs. Cultured HUVECs at baseline or treated with Phorbol myristate acetate (PMA, Sigma) or control DMSO (Sigma) were collected by trypsinization and then re-suspended in 1% BSA/PBS to a final concentration of 1 × 10 5 cells/500 μl. Rat Anti-Human EPCR monoclonal antibodies and isotype controls (BD Biosciences) were added as appropriate at a final concentration of 0.125 μg/500 μl and incubated at room temperature for 20 min in the dark. Cells were washed once with cold 1% BSA/PBS and resuspended in 1 ml ice-cold PBS prior to flow cytometric analysis using Gallios Flow Cytometer (Beckman Coulter) with Cytomics CXP software v2.2. Results were recorded as median fluorescence intensity and raw data were analyzed using Kaluza Analysis v1.3 (Beckman Coulter). We used one-tailedt-tests to test for differences in mean fluorescence intensities between the specific groups. We applied paired tests when comparing PMA vs vehicle control and unpaired when testing between genotypes.
Quantification of EPCR levels on human monocytes and neutrophils. We lysed 100-μl citrated whole blood samples at room temperature for 10 min using Lysing Solution 10X Concentrate (349202; BD Biosciences). Lysed blood was then centrifuged at 4°C and 600 × g for 6 min, and the pellet re-suspended in HEPES buffered saline (Sigma-Aldrich). Cultured cells were also re-suspended in HEPES buffered saline, to a final concentration of 10 5 cells/100 μl. Rat Anti-Human EPCR monoclonal antibodies and isotype controls were added as appropriate at a final concentration of 0.125 μg/100 μl and incubated at room temperature for 20 min in the dark. Samples were diluted in 0.5 ml ice-cold HEPES buffered saline prior to flow cytometric analysis using either a Cytomics FC500 with Cytomics CXP software v2.2 or a CytoFLEX S Flow Cytometer with CytExpert Acquisition and Analysis software v2.3 (Beckman Coulter). CD14 + Monocytes and CD16 + neutrophils from blood lysates were gated using forward and side light scatter, enabling discrimination by cell size and granularity, respectively. The gating strategy was validated using mouse anti-human CD14 and CD16 antibodies ( Supplementary  Fig. 7) added to blood lysates at a final concentration of 0.125 μg/100 μl. Results were recorded as median fluorescence intensity.
Reverse transcription quantitative PCR (RT-qPCR). HUVECs were seeded at a density of 31,250 cells/cm 2 (3 × 10 5 cells per well of a 6-well plate) in 2 ml medium and left to attach overnight. Cells were then co-incubated with 1 ng/ml TNF-α and varying concentrations of APC (0, 0.1, 1, 10, 100 nM) for a further 24 hr prior to cell lysis and RNA extraction using a Quick-RNA Microprep kit (Zymo Research). RNA was quantified using a NanoDrop Lite Spectrophotometer (Thermo Scientific). 1 μg of RNA was reverse transcribed using a Maxima H Minus First-Strand cDNA Synthesis kit with dsDNase (Thermo Scientific), and cDNA was diluted 1:20 in ddH 2 O. Quantitative PCR (qPCR) reactions were performed in solution containing 10 μl SYBR Green PCR Master Mix (Thermo Scientific), 70 nM of each forward and reverse primer, 4 μl cDNA and ddH 2 O to a total volume of 20 μl. The sequences of all primers used in this study are as follows: ACTB: forward 5ʹ-CCC TGG AGA AGA GCT ACG AG-3ʹ, reverse 5ʹ-GGA TGC CAC AGG ACT CCA T-3ʹ; GAPDH: forward 5ʹ-CCC ACT CCT CCA CCT TTG AC-3ʹ, reverse 5ʹ-CCA CCA CCC TGT TGC TGT A-3ʹ; RPLP0: forward 5ʹ-GCA TCT ACA ACC CTG  AAG TGC-3ʹ, reverse 5ʹ-TTG GGT AGC CAA TCT GCA GA-3ʹ; GUSB: forward  5ʹ-ACG TGG TTG GAG AGC TCA TT-3ʹ, reverse 5ʹ-TCT GCC GAG TGA AGA  TCC C-3ʹ; ICAM1: forward 5ʹ-TGA TGG GCA GTC AAC AGC TA-3ʹ, reverse 5ʹ-GCG TAG GGT AAG GTT CTT GC-3ʹ; VCAM1: forward 5ʹ-TGT GAA GGA  ATT AAC CAG GCT G-3ʹ, reverse 5ʹ-TGA CAC TCT CAG AAG GAA AAG C-3ʹ;  CCL2: forward 5ʹ-CAT GAA AGT CTC TGC CGC C-3ʹ, reverse 5ʹ-GGT GAC  TGG GGC ATT GAT TG-3ʹ; PROCR: forward 5ʹ-CGG TAT GAA CTG CGG GAA  TT-3ʹ, reverse 5ʹ-GTG TAG GAG CGG CTT GTT TG-3ʹ. qPCR reactions were run using a QuantStudio 6 Flex Real-Time PCR instrument with QuantStudio software v1.3 (Thermo Scientific). After an initial step of 15 min at 95°C, samples were subjected to 40 cycles of 30 sec at 95°C and 30 sec at 59°C, followed by dissociation curve analysis. Target C t -values were normalized using the arithmetic mean of four endogenous control genes (ACTB, GAPDH, RPLP0, GUSB) and results were analyzed using the Delta-Delta C t method. We applied a linear regression model for gene expression level~log(APC concentration). To test for significance of the observed APC effects, we used the F-test of the linear regression model. We tested the residuals for normality using the Shapiro-Wilk test and for equal variance using the Bartlett test.
In vitro static adhesion assay to assess the effects of recombinant sEPCR. To quantify U937-HUVEC interactions, we used an in vitro static adhesion assay, as previously described 87 . U937 cells were seeded at a density of 1 × 10 5 cell/ml in T25 flasks and differentiated into macrophages in the presence of 100 ng/mL phorbol 12-myristate 13-acetate (PMA) for 48 h. HUVECs were seeded at a density of 27,174 cells/cm 2 (i.e., 1 × 10 5 cells per well of a 12-well plate) in 1 ml medium and left to attach overnight. HUVECs were then treated with 10 ng/ml TNF-α or vehicle control for 4 h. U937 cells were collected and re-suspended in fresh medium at a concentration of 1 × 10 5 cells cell/ml, and then incubated with various concentrations of anti-Mac-1 or recombinant human sEPCR (i.e., 0, 3, 6, 12 ng/ml). HUVEC monolayers (at ≥90% confluence) were rinsed in Phosphate NaCl (PBSA) buffer and incubated with 1 ml U937 cell suspension comprising 1 × 10 5 cells (±anti-Mac-1/sEPCR) at 37°C for 5 min. After aspirating the U937 suspension, the HUVECs and any adherent U937 cells were gently rinsed four times in PBSA, and then a further 2 ml PBSA was added to the well. Using a phase-contrast videomicroscope (Leica Microsystems, DMI3000B), pictures at 10-fold magnification were taken, choosing four different fields at random. Quantification of cell adhesion events was performed using the ImagePro v6.3 software. We applied a linear regression model for adhesion events~treatment concentration. To test for significance of the observed sEPCR and anti-Mac-1 effects, we used the F-test of the linear regression model. To test for significance of the IgG effect, we used a nonparametric linear model. We assessed the residuals for normality using the Shapiro-Wilk test and for equal variance using the Bartlett test. We also tested for a difference of slope coefficients between the IgG and anti-Mac-1 conditions by fitting the regression model with an interaction term (adhesion events~treatment concentration*condition). The P value was calculated using an ANOVA of the linear regression model.
In vitro static adhesion assay to assess the effects of APC. U937 cells were seeded and treated with PMA as described above. HUVECs and HCAECs were seeded at a density of 27,174 cell/cm 2 (i.e., 1 × 10 5 cells per well of a 12-well plate) in 1 ml medium and left to attach overnight. Endothelial cells were treated for 24 h with either: (i) 1 ng/ml TNF alone, (ii) 1 ng/ml TNF and 100 nM APC, or (iii) vehicle. U937 cells were then collected and re-suspended as above, but with no further treatments. Endothelial cell monolayers (at ≥90% confluence) were rinsed in PBSA buffer and incubated with 1 ml U937 cell suspension comprising 1 × 10 5 cells at 37°C for 5 min. The monolayers were then rinsed, and adhesion events recorded and quantified as outlined above. We used paired one-tailedt-tests to test for differences in adhesion events between the TNF and TNF + APC conditions.
Reporting summary. Further information on research design is available in the Nature Research Reporting Summary linked to this article.