Full viral proteome sequences for SARS-CoV-2 were downloaded from UniProt (14). FASTA-formatted protein sequence data for each protein and protein class were examined individually and in combination. We produced potential peptides of 15 amino acids length (15mers), using sliding windows over the entire proteome. Proteins of fewer than 15 amino acids in length were not examined. Analyses were performed considering proteins individually, some protein domains individually, and in groupings of proteins (both the whole proteome, and all structural, non-structural, and accessory proteins).

HLA population frequencies were obtained from US unrelated stem cell donor registry National Marrow Donor Program (NMDP)/Be The Match. High resolution HLA Class II haplotype frequencies (DRB1, DRB3/4/5, DQA1, DQB1, DPA1 and DPB1 loci) were estimated using an expectation-maximization algorithm [as described by Gragert et al., 2013 (15)] utilising a cohort of 8.9 million US volunteer donors (NMDP/Be The Match registry snapshot 29/05/2020) HLA typed by molecular methods ( Table S1 showing number of individuals genotyped). US population categories were developed based on a race/ethnicity questionnaire included on the donor consent form. There were four broad population categories (European American, African American, Asian or Pacific Islander and Hispanic) and 21 detailed populations subcategories (African American, African Black, South Asian Indian, American Indian - South or Central American, Alaska native of Aleut, North American Indian, Caribbean Black, Caribbean Hispanic, Caribbean Indian, European Caucasian, Filipino, Hawaiian or other Pacific Islander, Japanese, Korean, Middle Eastern or North Coast of Africa, Mexican or Chicano, Chinese, Hispanic - South or Central American, Black - South or Central American, Southeast Asian, Vietnamese). Some individuals were assigned only membership of a broad population category.

Following Hardy-Weinberg equilibrium proportions, multi-locus HLA Class II genotypes were generated by randomly sampling two haplotypes from the same population HLA haplotype frequency distribution. Simulated genotypes were generated for each of the four broad and 21 detailed population groups, with ten replicates at each population of size 1,000, 5,000 and 10,000.

For population-based analyses at the haplotype level, we analysed haplotypes up to 99% cumulative coverage within each population. HLA alleles, when examined individually, included all HLA alleles present in any of these selected haplotypes for the four broad and 21 detailed population groups.

Peptide binding affinity was assessed for all Class II HLA that featured in haplotypes in any of the ethnic populations studied using NetMHCIIpan v4.0 (16). Peptides were examined both for their predicted binding affinity (nM) and their percentage rank (compared to a pool of representative peptides for the corresponding HLA). Peptides were defined as binders if their binding affinity was equal to or less than 500nM and their percentage rank was equal to or greater than 2% (default parameter for strong HLA class II peptide binders). Within the 9,590 15-mer peptides in the entire SARS-CoV-2 proteome, 4,289 peptides were predicted to bind strongly to one or more HLA Class II molecule. Analyses were also performed using alternative binding threshold cut-offs of ≤50nM peptide binding affinity threshold and ≤0.5% percentage rank in combination with a ≤500nM peptide binding affinity threshold. The NetMHCIIpan-4.0 programme also identifies the predicted 9-mer binding core. Of all HLA observed, on average 0.52% of haplotypes contained an HLA that was not in the list of 5,620 HLA available to run on NetMHCIIpan-4.0. Total predicted peptide counts for an individual HLA per protein or per whole proteome were calculated by counting peptides only if they were both unique to one another (i.e., a unique 15-mer), but also if the predicted binding core (9-mer) was also unique. This prevented the count from appearing falsely elevated due to sequential overlapping peptides which were presenting the same core, by preventing them from being counted more than once. Total peptide counts (in the context of an HLA haplotype or genotype) were also enumerated from peptides with both unique 15-mers and 9-mer cores, from binding peptides presented in any Class II HLA within the haplotype or genotype.

Peptide-HLA scoring models were assessed using scikit-learn (17) in Python, specifically using sklearn.metrics.roc_auc_score (AUROC), sklearn.metrics.average_precision_score (Average Precision), sklearn.metrics.accuracy_score (Accuracy), and sklearn.metrics.classification_report (Sensitivity and Specificity) functions (see Table S3 , with n=93 for peptide-HLA combinations). In order to assess AUROC for two predictor variables, a composite predictor was calculated consisting of the normalised sum of both. Population comparisons of peptide scores were performed by calculating the mean and standard deviation using NumPy in Python (see Table S4 , N=10,000 individuals in simulated population groups). Population simulations were assessed for normality using statistics.shaipro (Shapiro-Wilk test) and following demonstration of non-normality for all distributions examined, with statistics.kruskal (Kruskal–Wallis one-way analysis of variance) to assess differences between replicates of population simulations (all N=10,000), both implemented in Python.

The principal objective of our study was to examine whether genetic variation at the HLA loci may influence immune responses, and therefore COVID-19 clinical outcomes or response to vaccination, among patients from different ethnic populations. HLA class II molecules present peptides from exogenous antigens for CD4⁺ T-cell recognition and are therefore critical components for an effective adaptive immune response that incorporates humoral (B-cell) and cytotoxic (CD8⁺ T-cell) arms.

We analysed peptide binding for all classical HLA class II loci (HLA-DRB1, -DRB3/4/5, -DQA1, -DQB1, and -DPA1, -DPB1). The main analysis focused on four broad population categories, African Americans (AFA), European Caucasians (CAU), Hispanics (HIS) and Asian/Pacific Islanders (API) and sub-analysis included 21 detailed population subgroups of these broad categories. To account for linkage disequilibrium between HLA class II loci, we analysed HLA haplotypes rather than assuming allele frequencies were independent at each locus. HLA haplotype frequencies were estimated utilising HLA genotype data obtained from a cohort of approximately 8.9 million donors (NMDP/Be The Match registry; Table S1 ). We analysed a total of 25,128 unique haplotypes to enable 99% coverage of each population which in turn required the analysis of 803 HLA class II molecules (see Supplemental Information ). Figure S1 depicts the number of distinct HLA-DRB1, -DRB345 alleles and HLA-DQA1/DQB1, DPA1/DPB1 heterodimers examined and associated haplotype coverage in each population. To confirm the robustness of our observations at the individual level, we analysed multiple replicate genotype samples and demonstrate with replicate sets of simulated genotypes of 10,000 individuals that HLA diversity, as measured by the alpha parameter fits to the power law distribution, was stable, and the number of unique haplotypes needed to reach 95% cumulative frequency was also stable across replicates ( Table S6 ).

Immunogenic T-cell epitopes were identified using the full-length reference SARS-CoV-2 sequence (18) with 9,744 amino acids, which was subdivided into the four structural proteins (E, M, N, and S) and 7 additional open reading frames encoding non-structural proteins (NSP1-16) and accessory proteins (proteins 3a, 6, 7a, 7b, 8 and 10). A sliding window of 15 amino acids length was used and, to minimise redundancy, peptides were only counted towards totals if the HLA class II binding core was unique. To evaluate peptide binding and stable display by HLA class II molecules we employed NetMHCIIpan-4.0 which outputs predicted peptide-HLA binding affinity (IC50) in nanomolar units and percentile rank of binding affinity compared to a set of 100,000 random natural peptides. The percentile rank enables incorporation of information from the biological antigen presentation pathway in addition to the peptide-MHC binding event (16, 19). We aimed to detect strongly binding peptides and, therefore, elected to use a threshold of ≤2% percentile rank (default parameter for HLA class II strong peptide binding) combined with an affinity threshold ≤500nM (16). To further increase the stringency of criteria for prediction of peptide binding, such that we maximise precision and reduce the rate of false positive peptides (accepting we may not recall all possible peptides that could be displayed), we also explored a model involving a threshold of ≤0.5% percentage rank and ≤500nM binding affinity. Finally, we examined a binding affinity threshold of ≤50nM, as previously suggested (13). To validate the computational models we analysed publicly available datasets of experimentally determined, immunogenic SARS-CoV-2 peptides. We focused on the two largest datasets recently described by Snyder et al. (20, 21) and Nolan et al. (22) (250 HLA class II peptides) and by Mateus et al. (23) (135 HLA class II peptides) that contain peptides from the entire SARS-CoV-2 proteome. We also examined two relatively small datasets encompassing 9 nucleocapsid and 25 structural protein-derived (Spike, nucleocapsid or membrane) peptides (24, 25). In the aforementioned datasets, the HLA restriction was not known and validation focused on positive identification of experimental peptides by the computational models ( Supplemental Information ). We also used an HLA-DRB1*04:01 restricted peptide dataset determined using an in vitro peptide-HLA stability assay (26). Our analyses showed that scoring peptide binding based on a combination of ≤2% percentile rank and ≤500nM binding affinity achieved the best true positive rate (sensitivity) for predicting experimentally derived SARS-CoV-2 peptides ( Table S2 ). Similarly, in the HLA restricted dataset by Prachar et al, the combined ≤2% percentile rank and ≤500nM binding affinity threshold had an AUROC of 0.85 and provided the best combination of precision and specificity in classifying stable peptide binders compared to alternative scoring methods ( Table S3A ). Finally, we validated our approach using a recently published dataset of experimentally determined CD4⁺ T-cell epitopes (27) and demonstrated that our approach does not introduce systematic bias across the studied populations and does not over-predict the number of immunogenic epitopes ( Table S3B ). It should be noted that using a 50nM threshold, without accounting for binding characteristics of natural peptides, resulted in 46% of HLA alleles examined lacking presentation of any SARS-CoV-2 peptides and a bias towards HLA-DR as the major SARS-CoV-2 peptide presenting molecules ( Figure S2 ). Nevertheless, we have confirmed key aspects of our analyses using both a ≤0.5% percentage rank in combination with a ≤500nM peptide binding affinity threshold, and using a ≤50nM peptide binding affinity threshold alone. Overall, a total of 9,590 15-mer peptides, derived from all 11 SARS-CoV-2 genes, were examined and 4,289 peptides were predicted to bind strongly according to our defined criteria to at least one HLA class II molecule.

We first examined presentation of viral epitopes by all HLA class II molecules contained in haplotypes representing 99% of the four major broad ethnic populations. Assessment of presentation capacity at the entire viral proteome level ( Figure 1A ) showed that the majority of HLA class II molecules are capable of presenting SARS-CoV-2 peptides, albeit with significant variability. HLA-DR alleles have the highest viral peptide presentation capacity followed by HLA-DP and, to a significantly lower extent, HLA-DQ molecules. Notably, certain common individual HLA molecules were predicted to have very limited ability to present viral peptides, including DQA1*03:01~DQB1*02:01 (no peptide presentation from any protein within the SARS-CoV-2 proteome; with frequency of 3.2% within AFA haplotypes, <0.01% within API haplotypes, 0.1% within CAU haplotypes and 0.5% within HIS haplotypes), DRB1*03:02 (two peptides presented in total from the entire proteome; with frequency of 6.3% within AFA haplotypes, 0.01% within API haplotypes, 0.04% within CAU haplotypes and 1.0% within HIS haplotypes). Consistent with its known high immunogenicity (10, 27–31), Spike protein derived peptides showed strong binding for the majority of HLA class II molecules although presentation capacity again varied and was lowest (no peptides presented) for relatively common alleles such as DQA1*01:01~DQB1*05:03, with a frequency of 1.9-5.4% in each of the four broad population groups ( Figure 1B ). This observation was more prominent for Nucleocapsid derived peptides where strong binding was predicted to be absent in 117 out of 306 HLA molecules found in 99% of haplotypes in the four broad ethnic populations (mostly reflecting HLA-DP molecules and to some extent -DQ molecules; Figure 1C ). Similar findings were noted for the relatively small Membrane and Envelope structural proteins, with the latter predicted to be non-immunogenic for the majority of common HLA class II molecules ( Figures 1D, E ). For non-structural proteins, HLA class II presentation capacity was variable and dependent upon protein amino acid sequence length (data not shown).

The capacity of individuals to present viral peptides for recognition by CD4⁺ T-cells depends on the composition of their HLA class II alleles in their inherited haplotypes. Given variation in population-specific HLA haplotype frequencies, we hypothesised that potential differences in SARS-CoV-2 proteome immunogenicity (as reflected by HLA class II peptide presentation) at the population level may reflect disparities in capacity for effective anti-SARS-CoV-2 immunity which could in turn influence response to vaccination or underpin observed variability in COVID-19 clinical outcomes among different ethnic populations. To investigate this hypothesis, we examined viral proteome presentation by HLA class II, accounting for the distribution of HLA haplotypes in a population. This analysis showed that the overall capacity for HLA peptide presentation, at the whole SARS-CoV-2 proteome level, among the four broad ethnic populations examined was remarkably similar ( Figure 2A ). This was also the case considering T-cell epitopes from each SARS-CoV-2 protein ( Table S4 ), suggesting that polymorphism at the HLA genomic region is unlikely to underpin potential differences in immune responses and, thus, in clinical outcomes among ethnic groups. It was notable, however, that within ethnic populations the capacity of individual HLA haplotypes to present viral peptides varied widely (data on HLA haplotypes examined and SARS-CoV-2 peptide presentation for the four broad ethnic populations is provided in Supplemental Information ), with the top 5% of haplotypes predicted to present between 497 and 591 peptides from the entire viral proteome (such as DRB1*04:01~DRB4*01:01~DQA1*03:01~DQB1*03:01~DPA1* 01:03~DPB1*04:01 which represents 1% of AFA haplotypes, 0.14% of API haplotypes, 3.9% of CAU haplotypes and 0.9% of HIS haplotypes, predicted to present 506 SARS-CoV-2 peptides) as opposed to 5% of haplotypes at the opposite end of the spectrum predicted to present approximately 316 viral peptides or fewer (e.g. DRB1*03:01~DRB3*01:01~DQA1*05:01~DQB1* 02:01~DPA1*02:01~DPB1*01:01 which represents 0.25% of AFA haplotypes, 0.02% of API haplotypes, 3.5% of CAU haplotypes and 1.4% of HIS haplotypes, predicted to present 258 peptides). This observation suggests that individual capacity to mount CD4⁺ T-cell immune responses against SARS-CoV-2 is not uniform and is likely dependent on HLA phenotype. Similar inter-individual variability was noted for peptide presentation derived from structural and from non-structural proteins as well as for distinct SARS-CoV-2 proteins examined ( Figures 2B–E and S3 ). Recent experimental studies suggested that up to 70% of CD4⁺ T-cell responses against SARS-CoV-2 target the Spike, Membrane and Nucleocapsid antigens (10); our analysis showed the predicted CD4⁺ T-cell response to these structural proteins is highly variable at the haplotype level ( Figures 2B–D ). In agreement with the observed immunogenicity of Spike protein in experimental studies (28–31), relatively high numbers of Spike-derived peptides were predicted to be recognised both within and across ethnic populations ( Figure 2B ). In contrast, our analysis suggests that on average 13.9% of HLA haplotypes within each population have low capacity (≤2 peptides) to present Nucleocapsid-derived peptides (e.g. DRB1* 03:02~DRB3*01:01~DQA1*04:01~DQB1*04:02 ~DPA1*02:02~DPB1*01:01 which presents one peptide, and accounts for 4.98% of haplotypes found in AFA populations and 0.44% in HIS populations, although it is rare in API and CAU populations; Figure 2C ). This inter-individual variability may, in part, account for the heterogeneity in the presence and magnitude of CD4⁺ T-cell and antibody responses against the Nucleocapsid protein noted in recent COVID-19 studies (10, 32, 33). Among non-structural viral proteins, our analysis suggested NSP3, NSP4 and NSP12 as the most immunogenic, in part reflecting their size, in each population ( Figure S3 ). The above noted similarity in overall capacity for SARS-CoV-2 peptide presentation among different populations was also observed at different (more stringent) thresholds for HLA-peptide binding, albeit with even higher individual variability within ethnic populations (500nM binding affinity and ≤0.5% percentage rank or ≤50nM binding affinity; Figure S4 ).

To further explore the consequences of the above observations at the individual level, and given that every individual expresses two HLA haplotypes, we generated 10 replicates each of random genotype datasets encompassing 1,000, 5,000 and 10,000 simulated individuals for each population, as described in the methods. Populations encompassing 10,000 individuals achieved >95% cumulative HLA haplotype coverage in every population. As shown in Figure 3 , this analysis confirmed equivalent HLA class II presentation of SARS-CoV-2 peptides across all four broad ethnic populations, both at the entire viral proteome level and for individual proteins ( Table S4 ). As noted for the HLA haplotype analysis, there was wide inter-individual variability in predicted potential for CD4⁺ T-cell immune responses, according to HLA genotype. We noted significant, but variable, capacity for T-cell reactivity against the entire Spike glycoprotein across individuals, whereas reactivity against the Nucleocapsid protein was predicted to be weaker for 10% of individuals in each population, on average (range 7.8-13.3%, as defined by HLA presentation of less than 5 nucleocapsid peptides). These observations were confirmed (but, again, inter-individual variability was higher) using more stringent thresholds for defining HLA class II peptide presentation (data not shown). In further analysis, we considered immune reactivity against the Receptor Binding Domain (RBD) of Spike glycoprotein, as it represents a proposed target of coronavirus subunit vaccines currently in clinical trials (34–36). Although, overall, there was no significant difference in predicted CD4⁺ T-cell reactivity at the population level, there were notable differences at the individual level (both based on HLA haplotype and on genotype analyses) with wide variation in predicted RBD specific peptide presentation ( Figure 3F ). The above analyses were consistent irrespective of the size of the population sampled and among the 10 replicates at each population size (Kruskal-Wallis test p-value >0.05 for peptide comparisons at the population level).

We next examined HLA class II presentation of SARS-CoV-2 peptides for a further 21 detailed population subgroups, as described above ( Figure S5 and Table S4 ). Overall, analysis of the entire viral proteome identified similar capacity for HLA class II viral peptide presentation across population subgroups. This was also the case for structural proteins, including Spike and RBD, mirroring the findings above for the broad population groups. Although we did not identify SARS-CoV-2 vulnerability of particular populations at the HLA level, again we observed inter-individual variation in predicted cellular immunity within ethnic groups. This was reflected in the range of predicted viral peptide presentation within simulated populations of 10,000 individuals including 24-112 for Spike, 4-25 for RBD, 1-24 for Nucleocapsid and 208-854 for the entire proteome ( Table S4 ).

The effectiveness of peptide and subunit vaccine formulations against SARS-CoV-2 depends on robust presentation by individual HLA class II molecules and, therefore, investigated vaccines should account for linkage disequilibrium and HLA haplotype frequencies in different ethnic populations to achieve universal coverage. Capitalising on extensive HLA haplotype frequency data from the NMDP/Be The Match registry, we generated maps of SARS-CoV-2 immunogenicity for the entire viral proteome. This analysis showed that each viral protein contains immunogenic peptide segments with variable degree of population coverage which is, on the whole, similar across different populations for a given protein region, as well as peptide segments of variable length that are non-immunogenic in any population ( Figure 4 ). Similar observations were made for coronavirus subunit components that are being investigated as potential vaccines with several immunogenic peptides predicted to achieve universal coverage across population groups ( Figure 4 ). This was particularly evident for the Spike protein underlying its inherent immunogenicity and its potential as vaccination target. Table S5 , depicts SARS-CoV-2 immunogenic peptide segments predicted to cover over 90% of HLA genotypic variation in every broad ethnic group examined.