www.nature.com/npjgenmed ARTICLE OPEN Complement component C4 structural variation and quantitative traits contribute to sex-biased vulnerability in systemic sclerosis Martin Kerick 1,82 ✉, Marialbert Acosta-Herrera 1,2,82 ✉, Carmen Pilar Simeón-Aznar 3, José Luis Callejas4, Shervin Assassi 5, International SSc Group*, Susanna M. Proudman6, Mandana Nikpour7, Australian Scleroderma Interest Group (ASIG)*, PRECISESADS Clinical Consortium*, Nicolas Hunzelmann8, Gianluca Moroncini 9, Jeska K. de Vries-Bouwstra10, Gisela Orozco 11,12, Anne Barton 11,12, Ariane L. Herrick13, Chikashi Terao 14, Yannick Allanore 15, Carmen Fonseca16, Marta Eugenia Alarcón-Riquelme 17, Timothy R. D. J. Radstake18, Lorenzo Beretta19, Christopher P. Denton16, Maureen D. Mayes 5 and Javier Martin 1 ✉ Copy number (CN) polymorphisms of complement C4 play distinct roles in many conditions, including immune-mediated diseases. We investigated the association of C4 CN with systemic sclerosis (SSc) risk. Imputed total C4, C4A, C4B, and HERV-K CN were analyzed in 26,633 individuals and validated in an independent cohort. Our results showed that higher C4 CN confers protection to SSc, and deviations from CN parity of C4A and C4B augmented risk. The protection contributed per copy of C4A and C4B differed by 1234567890():,; sex. Stronger protection was afforded by C4A in men and by C4B in women. C4 CN correlated well with its gene expression and serum protein levels, and less C4 was detected for both in SSc patients. Conditioned analysis suggests that C4 genetics strongly contributes to the SSc association within the major histocompatibility complex locus and highlights classical alleles and amino acid variants of HLA-DRB1 and HLA-DPB1 as C4-independent signals. npj Genomic Medicine (2022)7:57 ; https://doi.org/10.1038/s41525-022-00327-8 INTRODUCTION pathogenesis. In this sense, structural variants not captured by Systemic sclerosis (SSc) is a chronic immune-mediated inflamma- GWAS, such as copy number (CN) polymorphisms, which have tory disease (IMID) more frequently observed in women (female/ been implicated in the etiology of several diseases8,9, could male ratio ~3.8–11.5) that affects the connective tissue and is contribute substantially to the genetic risk of SSc. Several CN associated with considerable morbidity and mortality1–3. The variants in immunological genes have been found to be heterogeneous clinical manifestations of SSc are characterized by associated with autoimmune diseases10–13, although technical functional and structural vasculopathy, fibrosis of the skin and limitations and the complexity of CN polymorphisms have internal organs, in addition to inflammatory and immunological reduced the impact of their analysis in understanding alterations like auto-antibody production1. The individual genetic autoimmunity14,15. background, together with environmental risk factors and The complement system plays an important role in innate epigenetics factors, play an important role in the pathogenesis immunity and forms a bridge to the adaptive immune of the disease4,5. response16–18. Functional abnormalities in the complement A recent genome-wide association study (GWAS) has identified system have been widely described in rheumatic diseases, such new genes and pathways implicated in the development and as rheumatoid arthritis (RA) or systemic lupus erythematosus (SLE), progression of SSc6. Similar to other IMIDs, these genetic and to a lesser extent in SSc19. Furthermore, genetic variability in variations account for a limited portion of estimated heritability7, several complement components may contribute to the develop- making clear that additional genetic variants remain to be found ment of inflammatory and autoimmune diseases20,21. with the potential to bring novel insights into disease etiology and 1 Department of Cell Biology and Immunology, Institute of Parasitology and Biomedicine López-Neyra, CSIC, Granada, Spain. 2Systemic Autoimmune Disease Unit, Hospital Clínico San Cecilio, Instituto de Investigación Biosanitaria Ibs. GRANADA, Granada, Spain. 3Department of Internal Medicine, Valle de Hebrón Hospital, Barcelona, Spain. 4Department of Internal Medicine, Hospital San Cecilio, Granada, Spain. 5Department of Rheumatology, The University of Texas Health Science Center at Houston, Houston, TX, USA. 6 Rheumatology Unit, Royal Adelaide Hospital and University of Adelaide, Adelaide, SA, Australia. 7The University of Melbourne at St. Vincent’s Hospital, Melbourne, VIC, Australia. 8 Department of Dermatology, University of Cologne, Cologne, Germany. 9Department of Clinical and Molecular Science, Università Politecnica delle Marche e Ospedali Riuniti, Ancona, Italy. 10Department of Rheumatology, Leiden University Medical Center, Leiden, The Netherlands. 11Center for Genetics and Genomics Versus Arthritis, Division of Musculoskeletal and Dermatological Sciences, School of Biological Sciences, Faculty of Biology, Medicine and Health, The University of Manchester, Manchester, UK. 12NIHR Manchester Biomedical Research Center, Manchester University NHS Foundation Trust, Manchester, Greater Manchester, UK. 13Division of Musculoskeletal and Dermatological Sciences, The University of Manchester, Northern care Alliance NHS Foundation Trust, Manchester Academic Health Science Centre, Manchester, UK. 14Laboratory for Statistical and Translational Genetics, RIKEN Center for Integrative Medical Sciences, Yokohama, Kanagawa, Japan. 15Department of Rheumatology A, Hospital Cochin, Paris, Île-de-France, France. 16Center for Rheumatology, Royal Free and University College Medical School, London, UK. 17Center for Genomics and Oncological Research (GENYO), Pfizer-University of Granada-Andalusian Regional Government, Granada, Spain. 18Department of Rheumatology and Clinical Immunology, University Medical Center Utrecht, Utrecht, The Netherlands. 19Referral Center for Systemic Autoimmune Diseases, Fondazione IRCCS Ca’ Granda Ospedale Maggiore Policlinico di Milano, Milan, Italy. 82These authors contributed equally: Martin Kerick, Marialbert Acosta-Herrera. *Lists of authors and their affiliations appear at the end of the paper. ✉email: mkerick@ipb.csic.es; m.acostaherrera@ipb.csic.es; javiermartin@ipb.csic.es Published in partnership with CEGMR, King Abdulaziz University M. Kerick et al. 2 Complement component 4 (C4), encoded by two closely linked, The association between C4 copy number variants and SSc is highly polymorphic genes C4A and C4B within the major modified by sex and HERV-K histocompatibility complex (MHC) class III region on chromosome We found a higher C4 CN to be protective for SSc (Fig. 1a and 6, is an important protein in the classical and lectin complement Table 1). Less than four copies of C4 were found in 28.6% of SSc activation pathways, which are major effectors for controlling patients, and 2.8% had less than three copies. Interestingly, less microbial infections and for promoting clearance of apoptotic cells than four copies of C4 were found in 28% of female and 32% of and soluble immune complexes17,19. C4A and C4B encode proteins male patients. In a simple additive model C4A, C4B and HERV-K with distinct affinities for their molecular targets21–23 and present CNs exhibited 1.8-fold variation in their relative risk of SSc (95% variability in genomic copy number (CN)24 and length. The long confidence interval (95% CI), [1.51–2.33]; P = 1.04 × 10−16). form is determined by the presence of a 6.4-kb human Logistic-regression analysis estimated a rather small difference endogenous retrovirus K (HERV-K) element in intron 9 of both between the protection afforded by each copy of C4A (odds ratio genes. C4 CN studies in autoimmunity have been mainly focused (OR) = 0.73; 95% CI = 0.69–0.78) or C4B (OR = 0.82; 95% on SLE providing no definitive results14,22,25 in part due to the CI = 0.77–0.87) (Table 1). We replicated our calculation in the complexity of the genetic variation of C4A and C4B and the high second cohort and performed a meta-analysis, showing consistent linkage disequilibrium (LD) of the MHC locus, which contains the results (Table 1). HLA genes, the strongest genetic associations with most The number of subjects in the first cohort permits the simple autoimmune diseases including SLE, RA, and SSc6,26,27. A recent additive model to be expanded to a more complex one study, however, managed to attribute most of the genetic investigating the predictors that influence each other. We found association of SLE and Sjögren’s syndrome (SjS) of the MHC locus evidence for three two-way interactions: an interaction of C4A to C428. Of note, the same study made the complex genetic with C4B, a sexual dimorphism of C4A and C4B, and an interaction variation of C4A and C4B accessible through the imputation of SNP of HERV-K CN with C4A and C4B (Supplementary Table 3). The full data using a large multi-ancestry panel of 2530 reference model log(risk) ~ b1 Sex + b2 C4A + b3 C4B + b4 HERV-K + b5 haplotypes28,29. Sex:C4A + b6 Sex:C4B + b7 HERV-K:C4A + b8 HERV-K:C4A + b9 CN variations at C4 genes have been implicated as a source of C4A:C4B predicted 7.7-fold variation in the relative risk of SSc. sexual dimorphism in two systemic autoimmune diseases, SLE and The following analyses have been derived from this complex 1234567890():,; SjS28, diseases that have a higher prevalence in women than men. model and its calculated coefficients. In light of the above, we set out to investigate the contribution of The relationship between SSc risk and C4A and C4B gene CNs C4 CN to SSc using data from the largest GWAS cohort of SSc, exhibited consistent, logical patterns across the 18 different CN combined with genetic, RNA-sequencing, and C4 serum data from combinations of C4A and C4B (Fig. 1a, b), which is based on an an independent cohort. interaction between C4A and C4B CNs that was suggested by logistic-regression analysis (bc4a:c4b = −0.14, P = 2.1 × 10−5). While higher total C4 CN provide protection, strong deviations from the RESULTS 1:1 ratio of C4A and C4B are of higher risk like e.g., four copies of Experimental design C4A and zero copies of C4B (Fig. 1b). To investigate the role of C4 genetics in SSc, we analyzed two We found evidence for a sexual dimorphism for C4A and C4B independent cohorts of European descent (see “Methods”, but not for HERV-K CN. Stratified analysis showed C4A to be more Supplementary Fig. 1 and Supplementary Table 1). We used the protective in men (bmale = −0.49, P = 1.7 × 10−9 vs bfemale = −0.27, first cohort (genetic data, N = 26,633, 34% SSc, 10 countries) to P = 6.5 × 10−12) while C4B showed statistical evidence to be determine the association of C4 CN to SSc. The second cohort protective only in women (bfemale = −0.22, P = 1.1 × 10−10 vs (genetic, whole blood expression, and C4 serum data, N = 857, bmale = −0.09, P = 0.23) (Table 1). Indeed, C4B CNs of two or 39% SSc, 9 countries) was analyzed to detect C4 expression higher seem to augment the risk for men (Fig. 1c). Logistic quantitative trait loci (eQTL) and for C4 expression- and C4 protein regression with interaction terms confirmed our results for the level- modeling. Expression-model eQTLs from the second cohort C4A:sex interaction, while the significance of the C4B:sex that explained expression variance were subsequently used in the interaction was only suggestive (bmale:c4a = −0.17, P = 0.0034, first cohort as co-factors to determine residual genetic association bmale:c4b = 0.13, P = 0.065) (Supplementary Table 3). For total C4 to SSc unrelated to C4 genetics in the MHC region on CN, we did not find a significant sex bias (bmale:c4 = −0.1, chromosome 6. Finally, we utilized classical HLA alleles and amino P = 0.074). We calculated the statistical power for C4B association acids (AAs) to model remnant C4-independent association in the in males to be 0.41 and for the interaction of C4B:sex to be 0.33. MHC region. HERV-K copies generally augment the risk for SSc (Fig. 1d and Table 1). The protection afforded by C4A and C4B CNs are affected C4 haplotype diversity and its correlation with classical HLA by HERV-K copies in a slightly different way. Protection associated alleles with each C4B copy is affected more strongly by HERV-K (bHERV- K:C4B = 0.085, P = 0.053) in comparison to C4A (bHERV-K:C4A = 0.042, The complex genetic variation of C4A and C4B, which consists of P = 0.03). A logistic-regression analysis with separate terms for the many haplotypes with different numbers of C4A and C4B genes long and the short forms of C4A and C4B confirmed that the short was recently made accessible for analysis in large cohorts28. We forms confer more protection than the long forms (Table 1). imputed 29 C4 haplotypes independently in both cohorts. Haplotype frequencies were found to be similar in both cohorts and comparable to published results28 but varied substantially C4 copy number affects C4 expression and C4 protein levels in across countries (Supplementary Fig. 2 and Supplementary Table whole blood 2A). We determined nine C4 haplotypes to have a moderate to We investigated the effect of C4 and HERV-K CN on whole blood strong correlation (r2 > 0.4) to at least one classical HLA allele expression levels of C4 (C4 = C4A + C4B, Fig. 2a). Specific analysis (Supplementary Dataset 1). Correlation to classical HLA alleles for C4A and C4B showed significant positive correlation for both associated to SSc on a genome-wide level was found to be small C4A CN (rPearson,C4A = 0.15, P = 1.2 × 10−5) and C4B CN (rPear- −24 (r2 < 0.3) for most haplotypes with the exceptions of HLA-B*08:01 son,C4B = 0.34, P = 1.7 × 10 ) with C4A or C4B expression. HERV- and HLA-C*07:01 which had a strong correlation (r2 = 0.82 and K CN had a weakening effect on C4 expression at all C4 CN levels 0.67, respectively) with the BS haplotype and HLA-C*16:01 with AL- (Fig. 2a). C4A or C4B expression models performed better if they BS (H2111B) of 0.63 (Supplementary Dataset 1). contained separate CNs for long (L) and short (S) versions of C4 npj Genomic Medicine (2022) 57 Published in partnership with CEGMR, King Abdulaziz University M. Kerick et al. 3 a C4 copy number and SSc risk b Interaction of C4A and C4B frequency greater (%) 0.6 C4A copy number 0 1 0 1 5 10 20 4 2 relative SSc risk 0.5 3 4 0.48 0.32 0.15 C4B copy number 3 0.4 SSc risk 2 0.3 1 0.2 0.1 0 2 3 4 5 6 0 1 2 3 4 total C4 copy number C4A copy number c C4B copies and SSc risk in males d HERV−K copies and SSc risk 0.30 C4B copy number 0.45 0 1 0.25 2 SSc risk in male individuals 3 4 0.40 0.20 SSc risk 0.15 0.35 0.10 HERV−K copy number 0.30 0.05 1 2 3 4 0.00 5 2 3 4 5 6 2 3 4 5 6 total C4 copy number total C4 copy number Fig. 1 C4 and HERV-K copy numbers and Systemic Sclerosis risk. a depicts relative systemic sclerosis (SSc) risk vs total C4 copy number stratified by C4A CN. The SSc risk score is calculated per individual as the sum of effect sizes (betas) multiplied with the design matrix. Betas of C4A, C4B and C4A:C4B were taken from the most complex model “d” (see “Methods”). Crosses are calculated as average relative risk per rounded C4 CN + /− 2 standard deviations (y axis). Linear regression lines are colored by C4A CN and drawn to visualize the interaction effect of C4A and C4B. The y axis contains a color code to aid a comparison with (b). b depicts the relative SSc risk of combinations of C4A and C4B CNs. Relative risk is calculated as in (a). Outer circles are drawn according to population frequency ranges of each C4A, C4B combination and highlight more common combinations. Diagonal dotted lines help to identify combinations of equal total C4 CN. c depicts relative SSc risk in male individuals vs total C4 CN stratified by C4B CN. Relative risk is calculated like in (a) using effect sizes of C4A, C4B, C4A:C4B, Sex:C4A, and Sex:C4B. Crosses are calculated as average relative risk per rounded C4 CN + /− 2 standard deviations (y axis). Cubic regression lines are colored by C4B CN and drawn to visualize the interaction effect of C4A and C4B. d depicts relative SSc risk vs total C4 CN stratified by HERV-K CN. Relative risk is calculated like in (a) using effect sizes of C4A, C4B, C4A:C4B, HERV-K:C4A, and HERV-K:C4B. Crosses are calculated as average relative risk per rounded C4 CN + /− 2 standard deviations (y axis). Linear regression lines are colored by HERV-K CN. (AS, AL, BS, BL) instead of C4A plus C4B plus total HERV-K CNs. C4BShort, and C4BLong can be predicted well (r2C4AL = 0.5, r2C4AS = Interestingly, expression models were best if CNs of C4B were 0.58, r2C4BL = 0.54, r2C4BS = 0.67, all P < 2.2 × 10−16) using C4 eQTLs included in the model of C4A and vice versa. About 21% (27%) of forward selected to explain C4A or C4B expression variance. expression variance of C4A (C4B) can be attributed to C4A and C4B Blood serum concentrations of C4 protein were well correlated CNs (Supplementary Table 4B, C). with C4 CNs (rPearson = 0.25, P = 1.3 × 10−12, Fig. 2b). Regression Using C4A and C4B eQTLs alone, we were able to explain up to analysis determined that C4 and HERV-K CNs, sex, age, and SSc 42% (C4B: 38%) of expression variance with 19 (C4B: 15) SNPs explained about 12% of C4 serum concentration variance (Supplementary Table 4B, C). This seems to be the upper bound of (P = 1.4 × 10−20) with HERV-K copies again weakening C4 levels C4 expression variance explained by C4 genetics as C4A and C4B (bHERV-K = −0.03, P = 9 × 10−7) (Fig. 2b and Supplementary Table CN and eQTLs together could not explain more than 40% of 5). CN-corrected serum C4 levels determined men to have more expression variance albeit with only 12 (C4B: 13) additional SNPs C4 protein than women (bmale = 0.04, P = 4.6 × 10−4) independent (Supplementary Tables 4B, C). C4 eQTLs seem to integrate C4 and HERV-K CN information. Indeed, copy numbers of C4AShort, C4ALong, of disease (Fig. 2e and Supplementary Table 5). Individuals with Published in partnership with CEGMR, King Abdulaziz University npj Genomic Medicine (2022) 57 M. Kerick et al. 4 Table 1. Logistic-regression analysis for total C4, C4A, C4B, and HERV-K copy numbers. 1st cohort (N = 26,633) 2nd cohort (N = 857) Meta-analysis Model Model terms Beta s.e. P Beta s.e. P Beta s.e. P a: all total C4 −0.23 0.03 6.3E-17 −0.20 0.13 0.12 −0.23 0.03 1.9E-17 HERV-K 0.12 0.02 2.9E-10 0.16 0.10 0.10 0.12 0.02 7.9E-11 b: all C4A −0.31 0.03 3.7E-19 −0.36 0.18 0.04 −0.31 0.03 4.7E-20 C4B −0.20 0.03 7.0E-11 −0.14 0.14 0.33 −0.19 0.03 4.6E-11 HERV-K 0.16 0.02 2.6E-13 0.23 0.11 0.04 0.16 0.02 3.4E-14 b: female C4A −0.27 0.04 6.5E-12 −0.81 0.33 0.01 −0.28 0.04 2.6E-12 C4B −0.22 0.03 1.1E-10 −0.10 0.27 0.70 −0.22 0.03 2.3E-13 HERV-K 0.15 0.02 1.4E-10 0.41 0.19 0.03 0.15 0.02 1.5E-14 b: male C4A −0.49 0.08 1.7E-09 −1.19 0.48 0.01 −0.51 0.08 1.1E-10 C4B −0.09 0.07 2.3E-01 −0.73 0.43 0.09 −0.11 0.07 1.2E-01 HERV-K 0.18 0.05 2.8E-04 0.66 0.28 0.02 0.19 0.05 7.5E-05 c: all C4AShort −0.29 0.18 1.0E-01 −0.43 0.32 0.18 −0.32 0.16 3.9E-02 C4ALong −0.16 0.02 2.0E-10 −0.11 0.18 0.55 −0.16 0.02 1.1E-15 C4BShort −0.20 0.03 1.4E-10 −0.12 0.28 0.67 −0.20 0.03 2.5E-11 C4BLong −0.04 0.03 1.8E-01 0.09 0.15 0.56 −0.04 0.03 2.3E-01 Depicted are beta values from the logistic-regression analysis of three different models (blocks of rows, see “Methods”). All models contained sex and five genetic principal components as co-variables. Logistic-regression analysis for the first cohort additionally contained cohort as co-variable. Model b was also calculated separately for females and males. Models contain copy numbers as calculated from the imputed C4 alleles per individual as dosages. SSc had less CN-corrected C4 protein than healthy subjects SNPs had a smaller effect with residual association (P < 5 × 10−8) independent of sex (bSSc = −0.02, P = 0.012). to SSc in three regions centered on HLA-DPB1, HLA-DRB1, and HLA- Despite that only about a third of SSc patients have less than four B (Supplementary Fig. 4B, C), which suggests a stronger copies of C4, we found significant downregulation of C4 expression contribution of C4A to SSc. (pfemale = 0.001, pmale = 0.003) and C4 protein levels (P = 0.004) in SSc While the ten independent C4 eQTL SNPs absorb as much SSc patients compared to healthy controls (Fig. 2c, d). association as possible by design (forward selection), they might have been selected due to the extensive genetic linkage in the C4 genetics can explain a part of the SSc association to the MHC region without being implicated in SSc pathogenesis. We MHC region therefore asked if SSc association with the MHC could be We performed conditional association analysis for genetic markers explained by C4 eQTLs selected to explain C4A and C4B expression across the MHC genomic region. Conditioning on C4A, C4B, and variance in the second dataset. Using expression-model eQTL HERV-K CN or on a risk score calculated using the complex C4 CN SNPs, selected in the second cohort to explain expression interaction model derived above, showed an impact on residual variance, as co-factors in the first cohort to determine residual association levels limited to the vicinity of the C4 gene genetic association with SSc again rendered most MHC associa- (Supplementary Fig. 3B, C). In addition, we calculated a C4 risk tion with SSc nonsignificant (P > 5 × 10−8), except for the peaks score recently proposed for SLE and SjS based only on C4A and shown in Fig. 3c centered on HLA-DPB1 and HLA-DRB1. C4B CN (as: risk = (2.3)C4A CN + C4B CN)28. We applied the SLE/SjS score for conditional association analysis in the first dataset and Remaining MHC signal after conditioning on C4 genetics found almost no effect (Supplementary Fig. 3A). highlights HLA-DRB1 and HLA-DPB1 Given the strong association of C4 CN with SSc and its rather Having attributed most SSc association within the MHC to C4 local impact in conditional analysis, we focused on C4 eQTLs as genetics, we investigated which classical HLA alleles and HLA potential modifiers of C4 CN risk as both CN and eQTLs affect C4 amino acids (AA) demonstrated C4-independent association to expression levels. Our analysis above suggests that the predictive SSc. Conditioning on expression derived independent C4 eQTLs power of C4 CN and C4 eQTLs on C4 expression levels is at least results in residual significance (P < 5 × 10−8) for classical alleles partly redundant. In fact, C4 eQTLs alone can explain more C4 and AAs of HLA-DRB1, HLA-DPB1, HLA-DPA1 and HLA-DQA1, HLA- expression variance than C4 CNs and C4 eQTLs are significant DQB1, and HLA-B (Supplementary Dataset 2A, B). Forward predictors for copy numbers of C4. We obtained 10,680 eQTLs of selection to derive independent residual signals marked HLA C4A and C4B from the GTEx v8 database and found that 37% of classical alleles for HLA-DRB1 (*07:01, *11:03, *11:04, *13:01), HLA- these are associated to SSc with pGWAS < 10−5 using the first DPB1 (*13:01, *26:01, *40:01, *06:01), HLA-DQA1 (*04:01), and HLA- dataset (Fig. 3a). This encouraged us to find independent C4 DQB1*05:01. eQTLs using forward selection to explain their contribution to the We tested if classical alleles and AAs of HLA-DRB1 and HLA-DPB1 MHC association to SSc. alone could explain residual C4-independent association. We Conditioned association analysis for genetic markers across the found that 9 HLA classical alleles (HLA-DRB1: *11:04, *08:01, *07:01, MHC genomic region on ten independent C4 eQTL SNPs rendered *13:01, *11:03 and HLA-DBP1: *13:01, *26:01, *40:01, *06:01) most association to SSc nonsignificant (P > 5 × 10−8), except for a together with 16 expression-model-derived independent C4 peak around HLA-DPB1, the initial association of which was eQTLs can explain almost all associations of SSc with the MHC enhanced by conditioning (Fig. 3b). While conditioning on 13 region (Supplementary Fig. 5B). The same is true for AAs of HLA- independent C4A exclusive eQTL SNPs had a similar effect on the DRB1 and HLA-DPB1, which can also explain residual C4- residual association profile, conditioning on 12 C4B exclusive eQTL independent association with eight independent AAs (HLA- npj Genomic Medicine (2022) 57 Published in partnership with CEGMR, King Abdulaziz University M. Kerick et al. 5 a C4 expression in blood b C4 protein in plasma HERV−K HERV−K 3 0.6 0 0 normalized C4 protein level residualized C4 expression 1 1 2 2 2 0.5 3 3 4 4 5 5 1 0.4 0 0.3 −1 0.2 −2 0.1 −3 0.0 2 3 4 5 2 3 4 5 C4 copy number C4 copy number c d e C4 expression C4 protein C4 protein by age 0.6 0.6 normalized C4 protein per C4 copy 4 female female female male male male normalized C4 protein level residualized C4 expression 0.5 0.5 2 0.4 0.4 0.3 0 0.3 0.2 0.2 −2 0.1 0.1 −4 ** 0.0 *** * ** 0.0 Controls SSc Controls SSc Controls SSc 25 30 35 40 45 50 55 60 Age (years) Fig. 2 C4 expression and C4 protein concentrations in whole blood. a depicts residualized total C4 expression levels by total C4 copy number (CN) stratified by HERV-K CN. C4 expression is calculated as the sum of C4A and C4B expression as obtained by RNA-Sequencing. The residualized expression has been calculated by regressing out 20 (18) principal components for controls and cases, respectively. Data has been grouped by rounded C4 and HERV-K CN dosage. b depicts normalized C4 protein levels in plasma by total C4 CN stratified by HERV-K CN. C4 protein levels have been normalized across 10+ laboratory sites. c depicts residualized total C4 expression levels (like in a) for SSc and controls, stratified by sex. Significant comparisons are highlighted by asterisk (*P < 0.05, **P < 0.01, ***P < 0.001). d depicts normalized C4 protein levels (like in b) for SSc and controls, stratified by sex. Significant comparisons are highlighted by asterisk (*P < 0.05, **P < 0.01, ***P< 0.001). e depicts normalized C4 protein levels in blood from 119 adult men (blue) and 447 adult women (red) as a function of age with locally estimated scatterplot smoothing (LOESS). Protein levels are normalized to the number of C4 gene copies in an individual’s genome. All boxplot are drawn with default settings in R 4.0.3: lines are defined as first, second and third quartile (Q1, Q2, Q3), whiskers depict the most extreme data points within Q1–1.5 interquartile range (IQR), and Q3 + 1.5 IQR. Boxplot notches are defined as 95% confidence interval of the median. DRB1: 37SY, 58A, 74RL, 96E, and HLA-DPB1: 8, 76I, 91H, 96x) SSc (Fig. 1b). Our results might differ from recent observations in (Fig. 3d). SLE and SjS where C4A and C4B copies have been described to act To complement our analysis, we repeated our search for in an additive way, but the authors did not describe the residual C4-independent association to SSc, this time conditioning interaction with HERV-K copies in detail28. on the C4 genetic signal which was not derived from the Our results showed a sexual dimorphism with respect to the expression dataset but from the first dataset as ten independent protection afforded by C4A and C4B. While in female individuals, C4 eQTL signals as described. Repeating the forward selection of C4A copies grant slightly more protection than C4B copies, our AAs or classical HLA alleles conditioning on the ten first dataset- data suggest that in male individuals only C4A confers protection derived independent C4 eQTLs signals resulted in four indepen- while we did not observe a strong effect for C4B. In male dent AAs from HLA-DPB1 (9F, 76I, 91H, 96x) or five independent individuals, C4B might therefore function like a null allele with alleles (HLA-DRB1: *11:04 and HLA-DBP1: *13:01, *26:01, *28:01, respect to protection from SSc as higher CNs of C4B are associated *30:01) (Supplementary Fig. 4D, E) supporting the role of HLA- to higher SSc risk (Fig. 1c). However, as the power of our study to DRB1 and HLA-DPB1 as a C4-independent association. detect significant C4B signals in males was limited, and the sex:C4B interaction was only of suggestive significance, we cannot rule out that C4B has a protective effect in males. While C4 alleles DISCUSSION have been described to act more strongly in men, no distinction In this study, we found a strong association of low C4 and high was made between C4A and C4B activity in SLE or SjS in a recent HERV-K CN with SSc in two independent SSc datasets and their study of similar size28. It has been described that activated C4A meta-analysis, supporting the protective role of C4 copies in IMIDs. bonds preferably with protein antigens, such as immune C4A gene copies were slightly more protective than C4B as has complexes, while activated C4B reacts rapidly with carbohydrate been shown in SLE and SjS28 but our data suggest a complex antigens, such as bacterial cell walls30. This could partly explain interaction of C4A and C4B CNs that has to be evaluated in the the greater susceptibility to and severity of infections reported in context of HERV-K copies and sex. We found that in SSc, an equal men and the higher incidence of autoimmune diseases in number of C4A and C4B gene copies grants more protection than women31,32. In addition, low C4B CNs have been associated to (strongly) imbalanced numbers which we found to be a risk for cardiovascular disease33 where the incidence in women is usually Published in partnership with CEGMR, King Abdulaziz University npj Genomic Medicine (2022) 57 M. Kerick et al. 6 Association to SSc a b C4 eQTLs conditioning 60 C4A 60 C4B none both on 10 C4 eQTLs −log10(p[SSc]) −log10(p[SSc]) [GWAS data selected] HLA− C4 DPB1 40 40 20 20 0 0 28 29 30 31 32 33 34 28 29 30 31 32 33 34 Chromosome 6 [mb] Chromosome 6 [mb] c d conditioning conditioning 60 60 none HLA− none on 16 C4 eQTLs DRB1 on 16 C4 eQTLs −log10(p[SSc]) −log10(p[SSc]) [Expr.data selected] HLA− & DRB1 + DPB1 AAs DPB1 40 40 20 20 0 0 28 29 30 31 32 33 34 28 29 30 31 32 33 34 Chromosome 6 [mb] Chromosome 6 [mb] Fig. 3 MHC region conditional association with systemic sclerosis. Association is calculated in the first dataset (N = 26,633) using logistic regression with cohort, genetic background (PC1-5), and sex as covariates and depicted as position (GRCh38) by significance (Manhattan plot) in gray if no additional covariates were used. The dotted line represents the genome-wide significance cutoff P = 5 × 10−8. a Manhattan plot with marked C4 eQTLs obtained from the GTEx v8 database. b Manhattan plot with additional conditioning on ten independent C4 eQTLs, obtained by forward selection in the first dataset, depicted in blue. The arrow marks the position of HLA-DPB1. c Manhattan plot with additional conditioning on 16 independent C4 eQTLs (obtained by forward selection to explain expression variance in the second dataset (N = 857) depicted in blue. The arrows mark the positions of HLA-DPB1 and HLA-DRB1. d Manhattan plot with additional conditioning on 16 independent C4 eQTLs (obtained by forward selection to explain expression variance in the second dataset) and 8 independent amino acids of DRB1 and DPB1 (obtained by forward selection in the first dataset conditioning on 16 independent C4 eQTLs) depicted in blue. lower than in men34. Interestingly, in a recent paper studying the production by bone-marrow-derived monocytes and macro- female-biased expression in human skin, several genes from the phages can restore humoral response in c4 deficient mice42. complement activation pathway were identified as a molecular Interestingly C4A and C4B expression models both profit strongly signature and in genome-wide co-expression networks35. Sexual from the other gene’s CN as a predictor, which supports the dimorphism has been extensively reported in vascular physiology genetic interaction between them suggested in this study. and pathophysiology36, where women more commonly develop Furthermore, the distinction between the long and the short microvascular dysfunction, and in autoimmune-related interstitial forms of C4: AS, AL, BS, and BL as expression predictors instead of lung diseases, where young women are most commonly C4A, C4B, and HERV-K CNs alone, greatly favors the accuracy of the affected37. All of these clinical manifestations, for which the role expression model. This suggests that HERV-K acts specifically on of C4 is yet to be elucidated, are hallmarks of SSc. the gene where it is located, to suppress its expression, which is Our data confirm that C4 and HERV-K CNs are strong predictors supported by studies in brain and serum22,29,43. C4A and C4B CNs of C4 expression levels in blood and other tissues29. While the were able to explain about 20% of C4A and C4B expression major site of C4 expression is the liver38, it has been shown that variance, which is clearly lower than the ability of C4 eQTLs, which whole blood can be used with some caution as a surrogate tissue could explain ~40%. Although we most likely over-fitted the for quantitative trait analysis39–41. In addition, local complement expression data, SNPs seem to be the superior instrument in npj Genomic Medicine (2022) 57 Published in partnership with CEGMR, King Abdulaziz University M. Kerick et al. 7 predicting C4 expression as they can integrate CN as well as Our study has some limitations. First, the number of samples in classical eQTL signals. This finding might help to bring C4 genetics the second dataset is very low in terms of GWAS. While we were to the clinic in the form of simple genetic tests. The interconnect- able to replicate the association of C4A, C4B, and HERV-K CNs, edness of C4 CN and eQTLs is further supported by C4 eQTLs replication of the most complex model was out of reach and being able to predict C4 CNs (AS, AL, BS, BL) with coefficients of needs to be the subject of further study. Second, we did not determination of r2 > 0.5. stratify the patients by clinical or serological subgroups. While our C4 and HERV-K CNs are strong predictors of C4 protein levels in results on HLA-DRB1 AAs, being associated with SSc indepen- blood serum43. Interestingly, it was reported that men had on dently of C4 genetics, point towards anti-topoisomerase auto- average 27% more C4 protein per C4 CN than women and that antibodies and probably diffuse cutaneous SSc, the topic is too this bias is stronger during reproductive years28. While we found vast to explore in this study. Third, unfortunately, we could not the difference between men and women to be smaller; men had distinguish C4A and C4B protein levels in serum, which would on average 17% more C4 protein per C4 CN than women; we think have been very useful, to further investigate the sexual dimorph- our dataset confirms this finding and its timeframe (Fig. 2f) ism described. Fourth, the exact amino acid positions and classical reinforcing the role of C4 in the differential susceptibility between alleles from the models calculated in our study might change in men and women observed in SSc. The deficiency of C4 may the future. New imputation reference panels might provide new trigger an inappropriate clearance of apoptotic debris and associations that could influence the models as forward selection stimulate chronic activation of myeloid cells. Also, it results in a is sensitive to statistical fluctuations. Last, C4 forms a genetic defect to eliminate autoreactive B-cell clones, and a higher module termed RCCX with three genes: serine/threonine nuclear tendency to form self-reactive germinal centers23 and has been protein kinase RP, steroid 21-hydroxylase CYP21, and extracellular previously associated with more severe SLE with earlier disease matrix protein tenascin TNX51. Although we only assessed C4 CNs onset44. We also observed that patients with SSc have lower associated with SSc, we cannot discard the possibility that this C4 serum levels than unaffected individuals even after correcting module plays a role in disease susceptibility. Specifically, TNX is for C4 gene CN, suggesting that hypocomplementemia in SSc is involved in the maintenance of collagen networks and tissue not simply due to C4 genetics but also reflects disease effects on integrity52, as well as in TGFB activation and signaling, typical for background complement levels45. fibrotic conditions such as SSc53,54. C4 expression was clearly down-regulated in SSc patients Many rheumatic diseases, including SSc, could benefit from compared to healthy individuals, as were C4 protein levels, therapies directed toward the complement system. These are although to a lesser extent (Fig. 2c, d). This might be explained by currently under active development and are not only focused on the difficulty in standardizing C4 protein assays across 10+ inhibitory mechanisms, but on activators or downstream activa- laboratory sites but might also point to differential mRNA stability tion fragments23. The inhibition of the complement pathway has adding another layer of complexity yet to be analyzed. Indeed, proven challenging. Eculizumab, a C5 inhibitor, is a complement- while we observed clear disease-independent differences of C4 targeting approved drug for a variety of vascular disorders and has protein levels between men and women, there was no significant recently been approved in kidney diseases55. Moreover, it has differential expression of C4 between healthy men and women, been studied in idiopathic inflammatory myopathies and SSc renal which suggests post-translational effects to play a role. In this line crisis, with promising results56,57. Our data suggest that C4 it has been proposed that IFN-gamma increases the stability of C3 genetics in SSc, by affecting expression and C4 protein levels, and C4 mRNA46 and a recent expression analysis in SSc detected a plays an important role in mediating the genetic association in the strong IFN signature in a subset of patients47. MHC locus and might also be involved in the epidemiological sex More than a third of more than 10,000 C4 eQTLs from the GTEx bias of SSc. This highlights the contribution of the complement to v8 database are associated with SSc (pGWAS < 10−5) and C4 eQTLs the development of SSc and to autoimmune disorders in general, alone can be used to explain most of the association of SSc within which could benefit from therapies directed towards the the MHC region, further supporting their importance in SSc. complement system. Our findings might help to bring C4 genetics Interestingly, C4A-specific eQTLs can explain more SSc association to the clinic in the form of simple genetic tests. than C4B-specific eQTLs (Supplementary Fig. 4B, C), which supports a stronger role for C4A in SSc. While C4 eQTLs could in principle be associated with SSc by the strong linkage structure METHODS present in the MHC, our data suggest that C4 eQTLs, forward Patients selected to explain expression variance in blood, can also explain All patients fulfilled the classification criteria of the 2013 American most genetic association with the MHC. Both analyses raise the College of Rheumatology (ACR) or The European League Against possibility that C4 genetics is indeed the main signal on Rheumatism (EULAR) or the criteria proposed by LeRoy and chromosome 6 for SSc, as has been suggested for SLE and SjS28, Medsger for early SSc58,59. CSIC’s Ethics Committee approved the both rheumatic diseases that can co-occur with SSc. study and written informed consent was obtained in accordance C4-independent genetic association with SSc centers on two with the Declaration of Helsinki. peaks (Fig. 3c), most of which can be explained by four AAs each of HLA-DPB1 and HLA-DRB1 (Fig. 3d). Interestingly, the AA positions for HLA-DPB1 and HLA-DRB1 overlap and all 8 AAs Cohorts and datasets positions can be associated with four of five binding pockets We (re-)analyzed two independent cohorts of European descent: described for class II HLA molecules48 likely interfering with (auto-) First cohort: genome-wide genotyped data from 14 indepen- antigen binding. In addition, three of the HLA-DRB1 AA positions dent epidemiological cohorts comprising a total of 28,179 (37, 58, and 74) are close to sites (30, 60, and 74) which have been unrelated individuals (9846 SSc patients and 18,333 healthy described to play a role in binding the consensus antigenic subjects) from 10 European countries6. To identify ancestry peptide of the topoisomerase I epitope, auto-antibodies to which outliers ~100,000 quality-filtered independent SNPs were selected define the ATA+ subgroup of SSc patients49. Furthermore, we from each case–control GWAS cohort. Principal component (PC) found that the C4-independent genetic association with SSc can analysis was performed using PLINK v1.07. Samples showing be explained by 10 independent classical HLA-Alleles instead of >4 standard deviations from the cluster centroids of each cohort AAs, seven of which overlap with a model of nine independent were considered outliers and removed from further analyses. The HLA-Alleles recently described50, which supports the indepen- presence of relatives and/or duplicates was assessed by comput- dence of C4 and HLA associations with SSc. ing identity-by-descent (IBD) estimation using PLINK v1.07. An Published in partnership with CEGMR, King Abdulaziz University npj Genomic Medicine (2022) 57 M. Kerick et al. 8 individual from each pair of relatives (Pi_Hat > 0.45) or duplicates software imputec429 and https://github.com/freeseek/imputec4 (Pi_Hat > 0.99) was removed. After exclusion of non-European and the reference panel downloaded from the dbGaP study samples, we recalculated genetic PCs using the merge of all accession: phs001992.v1.p128. Weighted imputation accuracy was imputed datasets, selecting ~100,000 independent markers using calculated by multiplying r2Allele by Allele frequency in Supple- PLINK v1.9. Missing data values due to the different platforms mentary Table 2B. used for genotyping were corrected by PLINK v1.9 (parameter –correct_for_missingness). We obtained informative principal C4 copy numbers. Each C4 haplotype carries a specific number of components as the visualization of the first two PCs can be C4 isotypes (C4A, C4B) and HERV-K elements (Supplementary Table interpreted as a “map” of the European continent (Supplementary 2C). We calculated total C4, C4A, C4B, and HERV-K CN dosages by Fig. 1). Second cohort: this cohort included genome-wide multiplying the allele dosages of the structural haplotype by the genotyped data, whole blood expression data and blood serum number of copies of each C4 isotype and HERV-K that the C4 protein concentrations from 333 SSc patients and 524 healthy haplotype contains. For instance, the haplotype AL-BL contains individuals from 9 European countries39. This second cohort is a one C4A gene and one C4B gene and two HERV-K copies. The subset of a larger cohort of seven immune-mediated diseases plus numbers of short and long forms of C4A and C4B (AL, AS, BL, BS) controls described here60. Individuals were excluded on the basis per haplotype are self-evident for 17 of 29 imputed haplotypes. of incorrect sex assignment, high missingness (>10%), non- For the remaining, long and short forms were inferred by the European ancestry (<55% using Frappe61 and high relatedness consensus that ~95% of C4A is present in the long form43,68–70. (PLINK v1.9 Pi_Hat >0.5). In addition, population stratification was The haplotype AL-BS for instance can be coded as 0.95 AL, 0.05 AS, also analyzed by PC analysis selecting ~100,000 independent 0.05 BL, and 0.95 BS. CNs per haplotype can be found in markers using PLINK v1.9. We obtained informative principal Supplementary Table 2C. components as the visualization of the first two PCs can be interpreted as a “map” of the European continent (Supplementary Classical HLA alleles and HLA amino acids (AA). Data for the Fig. 1). classical HLA alleles and AA variants were obtained from the first Basic clinical epidemiological information by cohort can be cohort by imputation using SNP2HLA71 and the reference panel found in Supplementary Table 1. from the Type 1 Diabetes Genetic Consortium72, described in ref. 50. After genotyping QC, all variants were imputed for each Expression data case–control dataset separately in the extended MHC region in chromosome 6. Imputed data were also filtered for 95% success Whole blood expression data was obtained from alpha and beta call rate for alleles and amino acids, deviation from HWE globin depleted (globinCLEAR, Ambion) total RNA. Single end considering a P value of <0.001 for SNPs in controls and 95% 50 bp stranded sequencing was performed on a HiSeq2500 total call rate for individuals50. Illumina within the PRECISESADS consortium60 and processed with bcl2fastq (Illumina), Cutadapt62, STAR v2.5.2 (2-pass default mapping to GRCh1963, and RSEM v1.2.3164 to obtain estimated Pearson correlation of C4 haplotypes and classical HLA Alleles counts per gene. Raw count data were normalized for quantitative Was calculated among the C4 haplotype dosages and the allele trait analyzes39. Briefly: Three genetic principal components (PCs) dosages from the HLA imputation. were regressed out from VSN-normalized65 raw read count data. Potential non-genetic influences were regressed out for SSc and C4 copy number association analysis controls separately by 20 (SSc:18) PCs calculated from inter- Logistic-regression models from simple to complex were calcu- sample expression correlation matrices. lated (using the function glm in R 4.0.3) to assess the association of total C4 dosage and its isotype dosages with the disease. We C4 protein data included cohort, five genome-wide principal components (PCs) Human complement C4 serum data was obtained from the and sex as covariates, assuming their effects were not collinear: PRECISESADS consortium60 from a turbidimetric immunoassay method according to the manufacturer’s recommendations (a) SSc ~ C4 + HERV-K + PC1-5 + cohort + sex (SPAPLUS analyzer)66. A corrective factor was calculated in order (b) SSc ~ C4A + C4B + HERV-K + PC1-5 + cohort + sex to normalize the data between the centers as described60. (c) SSc ~ C4Ashort + C4Along + C4Bshort + C4Blong + PC1-5 + cohort + sex Imputation The number of subjects in our first cohort permits us to expand SNPs. For both cohorts, we imputed SNPs from chromosome 6 the simple additive model to a more complex one investigating using the TOPMed reference panel with default settings at https:// the predictors that influence each other. We included three two- imputation.biodatacatalyst.nhlbi.nih.gov/67. Stringent QC mea- way interaction terms in the logistic-regression model: sures were applied to both cohort’s pre-imputation as follows: SNPs with call rates < 0.98; minor allele frequencies (MAFs) <0.01; (d) SSc ~ C4A + C4B + HERV-K + PC1-5 + cohort + sex + and those that deviated from Hardy–Weinberg equilibrium (HWE; C4A:C4B + C4A: HERV-K + C4B: HERV-K + C4A:sex + C4B:sex P < 0.001 in both case and control subjects) were filtered out from further analysis; samples with call rates <0.95 were removed. Meta-analysis was conducted with Metasoft73 using data from Relatives and/or duplicated samples were removed. Post- model a, b, and c from both datasets. imputation quality control included filtered for imputation quality (r2 > 0.3), MAF > 0.05, and HWE, which 6,39 resulted in 9,068 SSc patients and 17,565 healthy individuals for C4 haplotype imputa- Power calculation tion for the first cohort. Power calculations in CNV studies are problematic because effect sizes and models of the association are based on approximations C4 haplotypes. A set of 7021 SNPs TOPMed imputed SNPs were that may be unrealistic74. selected as they were (a) imputed in all individuals in both cohorts and (b) overlapped the C4 CN reference panel. C4 haplotype C4. Power calculations for C4B in males was carried out using the imputation was carried separately for both cohorts using the GAS Power Calculator [https://csg.sph.umich.edu/abecasis/ npj Genomic Medicine (2022) 57 Published in partnership with CEGMR, King Abdulaziz University M. Kerick et al. 9 gas_power_calculator/]. Here, we calculated the disease allele C4 gene expression analysis in whole blood frequency as (sum(CN C4B < 2)/N) = 0.26 and the genotype Using raw count data, we included disease, blood cell composi- relative risk as e0.09 = 1.09. Using an additive disease model and tion, and effective library size (calculated by EdgeR in R 4.0.3) in 1278 male cases and 6875 male controls, results in a power of the final model. While cell type-specific expression changes 0.406 to detect an association with P < 0.05. between SSc and controls were found significant at a nominal level for most cell types, the direction of expression change C4:sex interaction. We calculated the power to detect C4B:sex coincided for all cell types. We decided to report only whole blood interactions using “powerGWASinteraction”75 in R 4.0.3 which can expression changes controlling for blood cell composition. treat sex as environmental variable. We used: prevalence = 0.00034, pEnv = 0.306, betaC4B = −0.04, beta.sex = −1.35, beta.c4b:sex = 0.13, caseControlRatio = 0.34, ORgeneEnviron- C4 protein blood serum analysis ment = 1.03, alpha = 0.05 and alpha1 = 1. pGene probability was We included disease, sex, age, AS, AL, BS, and BL in the final calculated as (sum(CN C4B > = 2)/N) = 0.73. This results in a power model. The significance for the difference between SSc and to detect a C4B:sex interaction with P values <0.05 of 0.33. controls in men and women was calculated with both the Mann–Whitney test and a t test. Calculation of composite C4 risk score for SSc For each individual (i) a composite C4 risk score can be calculated Residual association of genetic variants across the MHC region as the sum of betas “Sb,i” from the (model-specific) effect sizes to SSc multiplied by the design matrix (of CN dosages, sex, interactions.. We performed conditional association analysis for genetic markers etc.) of the predictors. An individual relative risk score was then across the MHC genomic region. The first dataset was analyzed. In calculated as riski = eSb,i/(1+eSb,i). all models, we included cohort, five genome-wide PCs and sex as To visualize the interaction effect of C4A and C4B CNs on basic covariates. Association analysis of MHC region variants was relative risk, one multiplies the effect sizes (betas) of the most conditioned on the basic covariates plus: complex model “d” of C4A, C4B, and C4A:C4B with the design (1) nothing; matrix and calculates the relative risk score as above. To visualize (2) a risk score: 2.3 × C4A CN + C4B CN as proposed;28 the effect of HERV-K on risk we calculated a composite score with (3) covariates from model “b”: C4A CN + C4B CN + HERV-K CN; model “d” betas and C4A, C4B, HERV-K, C4A:HERV-K, C4B:HERV-K, (4) covariates from the most complex model “d” and C4A:C4B. To visualize the effect of C4B CN in males subjects, described above; we used effect sizes from model “d” for sex, C4A, C4B, sex:C4A, (5) C4 (C4A or C4B or both) eQTLs from GTEx v8 (obtained by sex:C4B, and C4A:C4B. forward selection in the first dataset until no SNP had pSNP < 10−5, see Supplementary Tables 4B, 6 and 4C, 6); Pearson correlation of C4 CN and C4 expression and C4 serum (6) C4A -specific eQTLs from GTEx v8 (obtained by forward levels selection in the first dataset until no SNP had pSNP < 10−5). EQTLs Was calculated with C4 CN dosages, the PC residualized were called C4A-specific if no C4B eQTL was reported in GTEx v8 expression data (see above) and the center corrected C4 protein with P < 0.01 for each SNP; serum concentrations. For visualization, CN dosages were rounded (7) C4B-specific eQTLs from GTEx v8 (obtained by forward to integers. selection in the first dataset until no SNP had pSNP < 10−5); (8) expression-model SNPs (with pGWAS < 10−5) (obtained by forward selection in the second dataset as described above, see C4 expression modeling Supplementary Tables 4B, 7 and 4C, 7). Total C4 expression was calculated as the sum of C4A and C4B expression. We used the linear model function “lm” in R 4.0.3 to calculate the adjusted coefficient of determination (r2) for each Residual, C4-independent, the association of the MHC region with SSc model with C4 CN, C4 CN + C4 eQTLs and C4 eQTLs alone as predictors. Model evolution is noted in Supplementary Table After accounting for the contribution of C4 genetics with models 4A–C. To add C4 eQTLs to the C4 CN model as expression “5” or “8” above, we sought to model residual, C4-independent, modifiers, we used forward selection. In a stepwise manner, we association of MHC SNPs with (a) forward selection of classical HLA selected the SNP to add to the model which had the most alleles; (b) forward selection of classical HLA alleles of HLA-DRB1 significant P value conditioning on all predictors already in the and HLA-DPB1; (c) forward selection of AAs of HLA genes; (d) model until no one more SNP was found with P < 0.01. In the same forward selection of AAs of HLA-DRB1 and HLA-DPB1. Forward way, SNPs were selected for the eQTLs only model until no more selection was carried out until no more HLA allele or AA was found SNPs were found with P < 0.01. To select SNPs in the expression with P < 5 × 10−8. (=second) dataset which were to be used for conditioned analysis of the MHC SNPs in the first dataset, forward selection was applied Reporting summary with SNPs which had pGWAS < 10−5 until no more SNP was found Further information on research design is available in the Nature with P < 0.01. Research Reporting Summary linked to this article. Modeling of C4 copy numbers using eQTLs DATA AVAILABILITY We tested if the eQTLs found to explain C4A or C4B expression Summary statistics of the SSc meta-GWAS are available through the NHGRI-EBI GWAS variance (see C4 expression modeling) can predict copy number Catalog (https://www.ebi.ac.uk/gwas/downloads/summary-statistics): GCST009131. dosages of the long and short forms of C4A and C4B: AS, AL, BS Data from the PRECISESADS consortium are hosted by ELIXIR Luxembourg https:// and BL. We used the linear model function “lm” in R 4.0.3 to elixir-luxembourg.org/ and are available upon request. The access procedure is calculate the adjusted coefficient of determination (r2) for each described on the data landing page (https://doi.org/10.17881/th9v-xt85). All other model with either C4A eQTLs or C4B eQTLs as predictors. data are contained in the article file and its Supplementary Information. Published in partnership with CEGMR, King Abdulaziz University npj Genomic Medicine (2022) 57 M. Kerick et al. 10 CODE AVAILABILITY 28. Kamitaki, N. et al. Complement genes contribute sex-biased vulnerability in All analysis has been performed with either the software described in “Methods” or diverse disorders. Nature 582, 577–581 (2020). within R 4.0.3. For logistic-regression analysis, the glm function of R was used. 29. Sekar, A. et al. Schizophrenia risk from complex variation of complement com- ponent 4. Nature 530, 177–183 (2016). 30. Yu, C. Y. et al. Dancing with complement C4 and the RP-C4-CYP21-TNX (RCCX) Received: 3 March 2022; Accepted: 22 September 2022; modules of the major histocompatibility complex. Prog. Nucleic Acid Res. Mol. Biol. 75, 217–292 (2003). 31. Klein, M. et al. Contribution of CD8+ T cells to inflammatory cytokine production in systemic sclerosis (SSc). Autoimmunity 49, 532–546 (2016). 32. Ingersoll, M. A. Sex differences shape the response to infectious diseases. PLoS REFERENCES Pathog. 13, e1006688 (2017). 1. Denton, C. P. et al. Systemic sclerosis. Lancet 390, 1685–1699 (2017). 33. Blasko, B. et al. Low complement C4B gene copy number predicts short-term 2. Elhai, M. et al. Mapping and predicting mortality from systemic sclerosis. Ann. mortality after acute myocardial infarction. Int. Immunol. 20, 31–37 (2008). Rheum. Dis. 76, 1897–1905 (2017). 34. Walli-Attaei, M. et al. Variations between women and men in risk factors, treat- 3. Bergamasco, A. et al. Epidemiology of systemic sclerosis and systemic sclerosis- ments, cardiovascular disease incidence, and death in 27 high-income, middle- associated interstitial lung disease. Clin. Epidemiol. 11, 257–273 (2019). income, and low-income countries (PURE): a prospective cohort study. Lancet 4. Angiolilli, C. et al. New insights into the genetics and epigenetics of systemic 396, 97–109 (2020). sclerosis. Nat. Rev. Rheumatol. 14, 657–673 (2018). 35. Liang, Y. et al. A gene network regulated by the transcription factor VGLL3 as a 5. Truchetet, M. E. et al. Current concepts on the pathogenesis of systemic sclerosis. promoter of sex-biased autoimmune diseases. Nat. Immunol. 18, 152–160 (2017). Clin. Rev. Allergy Immunol. https://doi.org/10.1007/s12016-021-08889-8 (2021). 36. Boese, A. C. et al. Sex differences in vascular physiology and pathophysiology: 6. Lopez-Isac, E. et al. GWAS for systemic sclerosis identifies multiple risk loci and estrogen and androgen signaling in health and disease. Am. J. Physiol. Heart Circ. highlights fibrotic and vasculopathy pathways. Nat. Commun. 10, 4955 Physiol. 313, H524–H545 (2017). (2019). 37. Han, M. K. et al. Female sex and gender in lung/sleep health and disease. 7. Bossini-Castillo, L. et al. Immunogenetics of systemic sclerosis: defining herit- increased understanding of basic biological, pathophysiological, and behavioral ability, functional variants and shared-autoimmunity pathways. J. Autoimmun. 64, mechanisms leading to better health for female patients with lung disease. Am. J. 53–65 (2015). Respir. Crit. Care Med. 198, 850–858 (2018). 8. Zhang, F. et al. Copy number variation in human health, disease, and evolution. 38. Blanchong, C. A. et al. Genetic, structural and functional diversities of human Annu Rev. Genomics Hum. Genet. 10, 451–481 (2009). complement components C4A and C4B and their mouse homologues, Slp and 9. Henrichsen, C. N. et al. Copy number variants, diseases and gene expression. C4. Int. Immunopharmacol. 1, 365–392 (2001). Hum. Mol. Genet. 18, R1–R8 (2009). 39. Kerick, M. et al. Expression quantitative trait locus analysis in systemic sclerosis 10. Fanciulli, M. et al. FCGR3B copy number variation is associated with susceptibility identifies new candidate genes associated with multiple aspects of disease to systemic, but not organ-specific, autoimmunity. Nat. Genet. 39, 721–723 pathology. Arthritis Rheumatol. 73, 1288–1300 (2021). (2007). 40. Basu, M. et al. Predicting tissue-specific gene expression from whole blood 11. Yang, Y. et al. Gene copy-number variation and associated polymorphisms of transcriptome. Sci. Adv. 7, eabd6991 (2021). complement component C4 in human systemic lupus erythematosus (SLE): low 41. Mu, Z. et al. The impact of cell type and context-dependent regulatory variants copy number is a risk factor for and high copy number is a protective factor on human immune traits. Genome Biol. 22, 122 (2021). against SLE susceptibility in European Americans. Am. J. Hum. Genet. 80, 42. Gadjeva, M. et al. Macrophage-derived complement component C4 can restore 1037–1054 (2007). humoral immunity in C4-deficient mice. J. Immunol. 169, 5489–5495 (2002). 12. de Cid, R. et al. Deletion of the late cornified envelope LCE3B and LCE3C genes as 43. Wouters, D. et al. High-throughput analysis of the C4 polymorphism by a com- a susceptibility factor for psoriasis. Nat. Genet. 41, 211–215 (2009). bination of MLPA and isotype-specific ELISA’s. Mol. Immunol. 46, 592–600 (2009). 13. McKinney, C. et al. Association of variation in Fcgamma receptor 3B gene copy 44. Juptner, M. et al. Low copy numbers of complement C4 and homozygous defi- number with rheumatoid arthritis in Caucasian samples. Ann. Rheum. Dis. 69, ciency of C4A may predispose to severe disease and earlier disease onset in 1711–1716 (2010). patients with systemic lupus erythematosus. Lupus 27, 600–609 (2018). 14. Olsson, L. M. et al. Copy number variation in autoimmunity—importance hidden 45. Esposito, J. et al. The association of low complement with disease activity in in complexity? Eur. J. Immunol. 42, 1969–1976 (2012). systemic sclerosis: a prospective cohort study. Arthritis Res. Ther. 18, 246 (2016). 15. Usher, C. L. et al. Complex and multi-allelic copy number variation in human 46. Mitchell, T. J. et al. IFN-gamma up-regulates expression of the complement disease. Brief. Funct. Genomics 14, 329–338 (2015). components C3 and C4 by stabilization of mRNA. J. Immunol. 156, 4429–4434 16. Carroll, M. C. The role of complement and complement receptors in induction (1996). and regulation of immunity. Annu. Rev. Immunol. 16, 545–568 (1998). 47. Beretta, L. et al. Genome-wide whole blood transcriptome profiling in a large 17. Walport, M. J. Complement. First of two parts. N. Engl. J. Med. 344, 1058–1066 European cohort of systemic sclerosis patients. Ann. Rheum. Dis. 79, 1218–1226 (2001). (2020). 18. West, E. E. et al. Complement and the regulation of T cell responses. Annu. Rev. 48. Agudelo, W. A. et al. Quantum chemical analysis of MHC-peptide interactions for Immunol. 36, 309–338 (2018). vaccine design. Mini Rev. Med. Chem. 10, 746–758 (2010). 19. Chen, M. et al. The complement system in systemic autoimmune disease. J. 49. Kongkaew, S. et al. Interactions of HLA-DR and topoisomerase I epitope modu- Autoimmun. 34, J276–J286 (2010). lated genetic risk for systemic sclerosis. Sci. Rep. 9, 745 (2019). 20. Holers, V. M. Complement and its receptors: new insights into human disease. 50. Acosta-Herrera, M. et al. Comprehensive analysis of the major histocompatibility Annu. Rev. Immunol. 32, 433–459 (2014). complex in systemic sclerosis identifies differential HLA associations by clinical 21. Goicoechea de Jorge, E. et al. Common and rare genetic variants of complement and serological subtypes. Ann. Rheum. Dis. 80, 1040–1047 (2021). components in human disease. Mol. Immunol. 102, 42–57 (2018). 51. Blanchong, C. A. et al. Deficiencies of human complement component C4A and 22. Yang, Y. et al. Diversity in intrinsic strengths of the human complement system: C4B and heterozygosity in length variants of RP-C4-CYP21-TNX (RCCX) modules serum C4 protein concentrations correlate with C4 gene size and polygenic in caucasians. The load of RCCX genetic diversity on major histocompatibility variations, hemolytic activities, and body mass index. J. Immunol. 171, 2734–2745 complex-associated disease. J. Exp. Med. 191, 2183–2196 (2000). (2003). 52. Matsumoto, K. I. et al. The roles of tenascins in cardiovascular, inflammatory, and 23. Wang, H. et al. Complement C4, infections, and autoimmune diseases. Front. heritable connective tissue diseases. Front. Immunol. 11, 609752 (2020). Immunol. 12, 694928 (2021). 53. Kasprzycka, M. et al. Tenascins in fibrotic disorders-from bench to bedside. Cell 24. Banlaki, Z. et al. Fine-tuned characterization of RCCX copy number variants and Adh. Migr. 9, 83–89 (2015). their relationship with extended MHC haplotypes. Genes Immun. 13, 530–535 54. Valcourt, U. et al. Tenascin-X: beyond the architectural function. Cell Adh. Migr. 9, (2012). 154–165 (2015). 25. Wu, Z. et al. Association between complement 4 copy number variation and 55. Ricklin, D. et al. The renaissance of complement therapeutics. Nat. Rev. Nephrol. systemic lupus erythematosus: a meta-analysis. Clin. Exp. Med. 20, 627–634 14, 26–47 (2018). (2020). 56. Gouin, A. et al. Role of C5 inhibition in idiopathic inflammatory myopathies and 26. Okada, Y. et al. Genetics of rheumatoid arthritis contributes to biology and drug scleroderma renal crisis-induced thrombotic microangiopathies. Kidney Int. Rep. discovery. Nature 506, 376–381 (2014). 6, 1015–1021 (2021). 27. Langefeld, C. D. et al. Transancestral mapping and genetic load in systemic lupus 57. Devresse, A. et al. Complement activation and effect of eculizumab in scler- erythematosus. Nat. Commun. 8, 16021 (2017). oderma renal crisis. Medicine 95, e4459 (2016). npj Genomic Medicine (2022) 57 Published in partnership with CEGMR, King Abdulaziz University M. Kerick et al. 11 58. van den Hoogen, F. et al. 2013 classification criteria for systemic sclerosis: an Reumáticas (RIER) from Instituto de Salud Carlos III (RD16/0012/0013). This work has American college of rheumatology/European league against rheumatism colla- received funding from the Innovative Medicines Initiative 1 & 2 Joint Undertaking borative initiative. Ann. Rheum. Dis. 72, 1747–1755 (2013). (JU) under grant agreements No 115565 (PRECISESADS) and No 831434 (3TR). The JU 59. LeRoy, E. C. et al. Criteria for the classification of early systemic sclerosis. J. receives support from the European Union’s FP7 and Horizon 2020 research and Rheumatol. 28, 1573–1576 (2001). innovation programs and EFPIA. MAH was supported by the Juan de la Cierva 60. Barturen, G. et al. Integrative analysis reveals a molecular stratification of systemic Incorporacion program, grant IJC2018-035131-I funded by MCIN/AEI/10.13039/ autoimmune diseases. Arthritis Rheumatol. 73, 1073–1085 (2021). 501100011033. This work is dedicated to the memory of Annette Kerick (1945-2020). 61. Tang, H. et al. Estimation of individual admixture: analytical and study design considerations. Genet Epidemiol. 28, 289–301 (2005). 62. Martin, M. Cutadapt remove adapter sequences from high-throughput sequen- AUTHOR CONTRIBUTIONS cing reads. EMBnet J. 17, 10–12 (2011). M.K., M.A.H., and J.M. contributed to the conception and study design. M.K. and 63. Dobin, A. et al. STAR: ultrafast universal RNA-seq aligner. Bioinformatics 29, 15–21 M.A.H. contributed to data collection, QC, imputation, and data analysis. C.P.S-A., (2013). J.L.C., S.A., S.M.P., M.N., N.H., G.M., J.K.V-B., G.O., A.B., A.H., C.T., Y.A., C.F., M.E.A-R., 64. Li, B. et al. RSEM: accurate transcript quantification from RNA-Seq data with or T.R.D.J.R., L.B., C.P.D., and M.D.M. contributed to GWAS and RNA-Sequencing data without a reference genome. BMC Bioinforma. 12, 323 (2011). collection. All co-authors made substantial contributions to data acquisition, data 65. Huber, W. et al. Variance stabilization applied to microarray data calibration and interpretation, and revised the work critically for important intellectual content. to the quantification of differential expression. Bioinformatics 18, S96–S104 (2002). 66. Capaldo, C. et al. The active immunological profile in patients with primary COMPETING INTERESTS Sjogren’s syndrome is restricted to typically encountered autoantibodies. Clin. Exp. Rheumatol. 34, 722 (2016). The authors declare no competing interests. 67. Taliun, D. et al. Sequencing of 53,831 diverse genomes from the NHLBI TOPMed Program. Nature 590, 290–299 (2021). 68. Mason, M. J. et al. Low HERV-K(C4) copy number is associated with type 1 dia- ADDITIONAL INFORMATION betes. Diabetes 63, 1789–1795 (2014). Supplementary information The online version contains supplementary material 69. Zai, C. C. et al. Association study of the complement component C4 gene in available at https://doi.org/10.1038/s41525-022-00327-8. tardive dyskinesia. Front. Pharm. 10, 1339 (2019). 70. Mariaselvam, C. M. et al. The complement C4 genetic diversity in first episode Correspondence and requests for materials should be addressed to Martin Kerick, psychosis of the OPTiMiSE cohort. Schizophr Bull Open 2, sgab003 (2021). Marialbert Acosta-Herrera or Javier Martin. 71. Jia, X. et al. Imputing amino acid polymorphisms in human leukocyte antigens. PLoS ONE 8, e64683 (2013). Reprints and permission information is available at http://www.nature.com/ 72. Brown, W. M. et al. Overview of the MHC fine mapping data. Diabetes Obes. reprints Metab. 11, 2–7 (2009). 73. Han, B. et al. Random-effects model aimed at discovering associations in meta- Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims analysis of genome-wide association studies. Am. J. Hum. Genet. 88, 586–598 in published maps and institutional affiliations. (2011). 74. Rucker, J. J. et al. Phenotypic association analyses with copy number variation in recurrent depressive disorder. Biol. Psychiatry 79, 329–336 (2016). 75. Dai, J. Y. et al. Two-stage testing procedures with independent filtering for Open Access This article is licensed under a Creative Commons genome-wide gene-environment interaction. Biometrika 99, 929–944 (2012). Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative ACKNOWLEDGEMENTS Commons license, and indicate if changes were made. The images or other third party We would like to thank Guillermo Barturen Briñas and Elena Carnero-Montoro for material in this article are included in the article’s Creative Commons license, unless fruitful discussions and Sofia Vargas and Gema Robledo for their excellent technical indicated otherwise in a credit line to the material. If material is not included in the assistance. We would like to thank Elena López-Isac for organizing all SSc GWAS article’s Creative Commons license and your intended use is not permitted by statutory datasets and all members of the PRECISESADS consortium, especially Ralf Lesche, regulation or exceeds the permitted use, you will need to obtain permission directly Sepideh Babaei, Anne Buttgereit, Suzana Makowska and Martina Runge for preparing from the copyright holder. To view a copy of this license, visit http:// the RNA Seq data and Johan Frostegård and Jacques-Olivier Pers for preparing and creativecommons.org/licenses/by/4.0/. normalizing the serum C4 data. We greatly appreciate the patients and healthy donors for their essential participation in these studies. This work was supported by © The Author(s) 2022 grant RTI2018101332-B-100 funded by MCIN/AEI/10.13039/501100011033 by “ERDF A way of making Europe”, Red de Investigación en Inflamación y Enfermedades INTERNATIONAL SSC GROUP P. Carreira20, I. Castellvi21, José Luis Callejas4, R. Ríos4, R. García Portales22, A. Fernández-Nebro23, F. J. García-Hernández24, M. A. Aguirre25, B. Fernández-Gutiérrez26, L. Rodríguez-Rodríguez26, P. García de la Peña27, E. Vicente28, J. L. Andreu29, M. Fernández de Castro29, F. J. López-Longo30, and Carmen Pilar Simeón-Aznar 3, V. Fonollosa3, A. Guillén3, G. Espinosa31, C. Tolosa32, A. Pros33, E. Beltrán33, M. Rodríguez Carballeira34, F. J. Narváez35, M. Rubio Rivas35, V. Ortiz-Santamaría36, A. B. Madroñero37, M. A. González-Gay38, B. Díaz39, L. Trapiella39, M. V. Egurbide40, P. Fanlo-Mateo41, L. Saez-Comet42, F. Díaz43, J. A. Roman-Ivorra44, J. J. Alegre Sancho45, M. Freire46, F. J. Blanco Garcia47, N. Oreiro47, T. Witte48, A. Kreuter49, G. Riemekasten50, P. Airò51, C. Magro10, A. E. Voskuyl52, M. C. Vonk53, R. Hesselstrand54, A. Nordin55, C. Lunardi56, A. Gabrielli57, A. Hoffmann-Vold58, J. H. W. Distler59, L. Padyukov55 and B. P. C. Koeleman60 20 Department of Rheumatology, 12 de Octubre University Hospital, Madrid, Spain. 21Department of Rheumatology, Santa Creu i Sant Pau University Hospital, Barcelona, Spain. 22 Department of Rheumatology, Virgen de la Victoria Hospital, Málaga, Spain. 23Department of Rheumatology, Carlos Haya Hospital, Málaga, Spain. 24Department of Internal Medicine, Virgen del Rocío Hospital, Sevilla, Spain. 25Department of Rheumatology, Reina Sofía/IMIBIC Hospital, Córdoba, Spain. 26Department of Rheumatology, San Carlos Clinic Hospital, Madrid, Spain. 27Department of Rheumatology, Madrid Norte Sanchinarro Hospital, Madrid, Spain. 28Department of Rheumatology, La Princesa Hospital, Madrid, Spain. 29 Department of Rheumatology, Puerta de Hierro Hospital-Majadahonda, Madrid, Spain. 30Department of Rheumatology, Gregorio Marañón University Hospital, Madrid, Spain. 31 Department of Internal Medicine, Clinic Hospital, Barcelona, Spain. 32Department of Internal Medicine, Parc Tauli Hospital, Sabadell, Spain. 33Department of Rheumatology, Hospital Del Mar, Barcelona, Spain. 34Department of Internal Medicine, Hospital Universitari Mútua Terrasa, Barcelona, Spain. 35Department of Rheumatology, Bellvitge University Published in partnership with CEGMR, King Abdulaziz University npj Genomic Medicine (2022) 57 M. Kerick et al. 12 Hospital, Barcelona, Spain. 36Department of Rheumatology, Granollers Hospital, Granollers, Spain. 37Department of Internal Medicine, Hospital General San Jorge, Huesca, Spain. 38 Epidemiology, Genetics and Atherosclerosis Research Group on Systemic Inflammatory Diseases, DIVAL, University of Cantabria, Santander, Spain. 39Department of Internal Medicine, Hospital Central de Asturias, Oviedo, Spain. 40Department of Internal Medicine, Hospital Universitario Cruces, Barakaldo, Spain. 41Department of Internal Medicine, Hospital Virgen del Camino, Pamplona, Spain. 42Department of Internal Medicine, Hospital Universitario Miguel Servet, Zaragoza, Spain. 43Department of Rheumatology, Hospital Universitario de Canarias, Tenerife, Spain. 44Department of Rheumatology, Hospital Universitari i Politecnic La Fe, Valencia, Spain. 45Department of Rheumatology, Hospital Universitari Doctor Peset, Valencia, Spain. 46Department of Internal Medicine, Thrombosis and Vasculitis Unit, Complexo Hospitalario Universitario de Vigo, Vigo, Spain. 47 Department of Rheumatology, INIBIC-Hospital Universitario A Coruña, La Coruña, Spain. 48Department of Clinical Immunology, Hannover Medical School, Hannover, Germany. 49 Department of Dermatology, Josefs-Hospital, Ruhr University Bochum, Bochum, Germany. 50Clinic of Rheumatology, University of Lübeck, Lübeck, Germany. 51Service of Rheumatology and Clinic Immunology Spedali Civili, Brescia, Italy. 52Department of Rheumatology, VU University Medical Center, Amsterdam, The Netherlands. 53Department of Rheumatology, Radboud University Nijmegen Medical Center, Nijmegen, Netherlands. 54Department of Rheumatology, Lund University, Lund, Sweden. 55Division of Rheumatology, Department of Medicine, Karolinska University Hospital, Karolinska Institute, Stockholm, Sweden. 56Department of Medicine, Università degli Studi di Verona, Verona, Italy. 57Istituto di Clinica Medica Generale, Ematologia ed Immunologia Clinica, Università Politecnica delle Marche, Ancona, Italy. 58Department of Rheumatology, Oslo University Hospital, Oslo, Norway. 59Department of Internal Medicine 3, Institute for Clinical Immunology, University of Erlangen-Nuremberg, Erlangen, Germany. 60 Department of Genetics, University Medical Center Utrecht, Utrecht, The Netherlands AUSTRALIAN SCLERODERMA INTEREST GROUP (ASIG) Susanna M. Proudman6, Mandana Nikpour7, W. Stevens7, J. Zochling61, J. Sahhar62, J. Roddy63, P. Nash64, K. Tymms65, M. Rischmueller66 and S. Lester66 61 Menzies Research Institute Tasmania, University of Tasmania, Hobart, TAS, Australia. 62Department Rheumatology, Monash Medical Centre, Melbourne, VIC, Australia. 63 Rheumatology, Royal Perth Hospital, Perth, WA, Australia. 64Research Unit, Sunshine Coast Rheumatology, Maroochydore, QLD, Australia. 65Canberra Rheumatology, Canberra, ACT, Australia. 66Department Rheumatology, The Queen Elizabeth Hospital, Woodville, SA, Australia PRECISESADS CLINICAL CONSORTIUM Lorenzo Beretta19, Barbara Vigone19, Jacques-Olivier Pers67, Alain Saraux67, Valérie Devauchelle-Pensec67, Divi Cornec67, Sandrine Jousse-Joulin67, Bernard Lauwerys68, Julie Ducreux68, Anne-Lise Maudoux68, Carlos Vasconcelos69, Ana Tavares69, Esmeralda Neves69, Raquel Faria69, Mariana Brandão69, Ana Campar69, António Marinho69, Fátima Farinha69, Isabel Almeida69, Miguel Angel Gonzalez-Gay Mantecón38, Ricardo Blanco Alonso38, Alfonso Corrales Martínez38, Ricard Cervera31, Ignasi Rodríguez-Pintó31, Gerard Espinosa31, Rik Lories70, Ellen De Langhe70, and Nicolas Hunzelmann8, Doreen Belz8, Torsten Witte48, Niklas Baerlecken48, Georg Stummvoll71, Michael Zauner71, Michaela Lehner71, Eduardo Collantes72, Rafaela Ortega-Castro72, Ma Angeles Aguirre-Zamorano72, Alejandro Escudero-Contreras72, Ma Carmen Castro-Villegas73, María Concepción Fernández Roldán73, Norberto Ortego74, Enrique Raya75, Inmaculada Jiménez Moleón75, Enrique de Ramon76, Isabel Díaz Quintero76, Pier Luigi Meroni77, Maria Gerosa77, Tommaso Schioppo77, Carolina Artusi77, Carlo Chizzolini78, Aleksandra Zuber78, Donatienne Wynar78, Laszló Kovács79, Attila Balog79, Magdolna Deák79, Márta Bocskai79, Sonja Dulic79, Gabriella Kádár79, Falk Hiepe80, Velia Gerl80, Silvia Thiel80, Manuel Rodriguez Maresca81, Antonio López-Berrio81, Rocío Aguilar-Quesada81 and Héctor Navarro-Linares81 67 Centre Hospitalier Universitaire de Brest, Hospital de la Cavale Blanche, Brest, France. 68Pôle de pathologies rhumatismales systémiques et inflammatoires, Institut de Recherche Expérimentale et Clinique, Université catholique de Louvain, Brussels, Belgium. 69Centro Hospitalar do Porto, Porto, Portugal. 70Katholieke Universiteit Leuven, Leuven, Belgium. 71 Medical University Vienna, Vienna, Austria. 72Servicio Andaluz de Salud, Hospital Universitario Reina Sofía Córdoba, Córdoba, Spain. 73Servicio Andaluz de Salud, Complejo hospitalario Universitario de Granada (Hospital Universitario San Cecilio), Granada, Spain. 74Department of Medicine, University of Granada, Granada, Spain. 75Servicio Andaluz de Salud, Complejo hospitalario Universitario de Granada (Hospital Virgen de las Nieves), Granada, Spain. 76Servicio Andaluz de Salud, Hospital Regional Universitario de Málaga, Málaga, Spain. 77Università degli studi di Milano, Milan, Italy. 78Hospitaux Universitaires de Genève, Genève, Switzerland. 79University of Szeged, Szeged, Hungary. 80 Charite, Berlin, Germany. 81Andalusian Public Health System Biobank, Granada, Spain. npj Genomic Medicine (2022) 57 Published in partnership with CEGMR, King Abdulaziz University