This study is part of the ethical approval 2023/6-1524-984, Highly-immune compatible iPS cells as source for regenerative medicine and cell therapy-oriented applications, from the Vilnius Regional Biomedical Research Ethics Committee (Lithuania) to Vilnius University, and 2023/4-1507-968, Analysis of the distribution of Human Leukocyte Antigen (HLA; Encoding Genes - HLA) alleles and haplotypes in the group of the Lithuanian unrelated bone marrow donor registry, to Vilnius University Hospital Santaros Klinikos. Written consent was obtained from the participants of the study.

For population-based analyses of HLA frequencies, the study included 3,496 individuals from the Lithuanian unrelated bone marrow donor registry, characterized at third-field resolution for HLA-A, HLA-B, and HLA-C. For the isolation of dermal fibroblasts, individuals were healthy adults who provided study-specific informed consent and were selected based on their known HLA class I genotypes. Individuals aged over 55 years, those with known inherited genetic disorders, or those diagnosed with non-environmentally caused diseases were excluded from dermal biopsy collection to ensure that fibroblast samples were free from age-associated mutations or pathogenic genetic variants.

HLA typing for registry donors’ peripheral blood was performed at the EFI-accredited immunogenetics laboratory at Vilnius University Hospital Santaros Klinikos (Vilnius, Lithuania) using sequencing-based typing, and at the ASHI-accredited laboratory HistoGenetics (Ossining, NY, USA) using next-generation sequencing. Exons 2 and 3 for class I HLA were covered.

Skin samples were collected using a 2–3-mm biopsy punch needle and fragmented with a sterile scalpel and needle. Fibroblasts were grown in AmnioPrime Complete Medium (cat. no. APR-B, Capricorn Scientific, Germany), supplemented with amphotericin B (cat. no. AMP-B, Capricorn Scientific, Germany), for 21 to 45 days until fibroblasts migrated from tissue sections and reached 80%–90% confluence. The medium was changed every 3 days to ensure optimal cell growth. Fibroblasts were routinely passaged with 0.25% Trypsin-EDTA at a density of 2 × 10⁵ cells/cm². Genomic DNA from fibroblasts was purified using the DNeasy Blood and Tissue Kit (cat. no. 69504, Qiagen, Germany) and genotyped using the primers HLAA-P1: TCCAGGTGGACAGGTAAGGA, HLAA-P2: GTCACTGCCTGGGGTAGAAC, HLAB-P1: TGCATTCTGGGTTTCTCTACTGG, HLAB-P2: CACGCGAAACATCCCAATCA, HLAC-P1: AGGTAAGGCAAAGGGTGGGA, and HLAC-P2: AGGCCGCCTGTACTTTTCTC. Samples were Sanger sequenced using the primers HLAA-P3: ACCCTCGTCCTGCTACTCTCG, HLAB-P3: ACCCTCCTCCTGCTGCTCTG, and HLAC-P3: CGTTGGGGATTCTCCACTCC at Microsynth, Germany.

Python and R scripts used for data analysis, along with anonymized datasets, are available through the Supplementary Data and or through the open-source GitHub developer platform in the repository https://github.com/Arias-Lab/superdonors.

The total allele count in the dataset was divided by the number of alleles (n = 2) times the number of individuals in this study (n = 3,496), all of whom had at least third-field resolution.

The observed genotypes present in the population were quantified (n = 3,496). The allele frequencies were determined using the sampled genotype count, and the expected genotype frequencies were calculated. The observed and expected genotype counts were compared with a χ² test. The χ² test is reliable for genotypes present more than five times in the population. Genotypes with a count < 5 times were filtered from the Hardy–Weinberg equilibrium (HWE) analyses. The degrees of freedom (df), calculated as (n(n + 1)/2) – n, were estimated based on the number of possible genotypes and the number of alleles identified in the sampled population for each HLA class I gene: 44 for HLA-A, 83 for HLA-B, and 45 for HLA-C.

Allele frequencies were extracted from the publicly available data from European-American (19) and British (20) populations and compared to the allele frequency from our study. Linear regression analyses (y ~ mx + c) were performed using R for pairwise comparison of allele frequencies of HLA-A, HLA-B, and HLA-C. Frequencies are calculated as frequency = allele count(dataset)/n(dataset).

Monte Carlo population haplotypes were simulated based on the published allele frequencies of European-American and British cohort studies. Data were processed with one-hot encoding to convert allele entries per individual into 1 or 0, using the caret library in R (22). Centroids and Euclidean distances were calculated from the principal components. Distances were represented as edges and as heatmaps.

The sequences for all alleles at the protein, transcript, and gene levels were downloaded as FASTA files from the IPD-IMGT-HLA database version 3.58 (23) and analyzed in Python and R. Allele sequences were extracted based on the HLA alleles present in the population. Cas9 binding sites were extracted with Python and analyzed in R using CrisprScore (24). Transmembrane prediction was conducted with DeepTMHMM (25).

The Lithuanian Bone Marrow Donor Registry, located at Vilnius University Hospital Santaros Klinikos, includes 13,884 individuals, with 11,153 characterized at the second field (protein level) for HLA-A, HLA-B, and HLA-C. Of these, 3,496 individuals are characterized in the third field ( Figure 1A ). We found that 858 individuals are at least homozygous for one HLA class I gene. A total of 542 individuals are homozygous for the coding sequence of HLA-A, 233 individuals are homozygous for HLA-B, and 338 individuals are homozygous for HLA-C ( Figure 1B ). The HLA types identified and their prevalence in the population are summarized in Figure 1 and Supplementary Table S1 . The five most frequent HLA-A alleles are A*02:01:01, A*03:01:01, A*24:02:01, A*01:01:01, and A*11:01:01, which together account for 74.6% of the population ( Figure 1C ). Notably, HLA-A*02:01:01 is the most frequent HLA class I allele, representing 31.6% of the population. Similarly, the five most frequent HLA-B alleles are B*07:02:01, B*13:02:01, B*15:01:01, B*44:02:01, and B*40:01:01, which account for 43.2% of the population ( Figure 1D ). HLA-B*07:02:01 alone represents 15.1% of the Lithuanian population. Furthermore, the five most frequent HLA-C alleles are C*07:02:01, C*06:02:01, C*04:01:01, C*02:02:02, and C*07:01:01, with a cumulative frequency of 59.2% in the population ( Figure 1E ). It is important to highlight that the HLA-B gene exhibits the largest diversity of alleles, followed by HLA-A and HLA-C ( Supplementary Table S1 ), as also observed in previous studies (19–21). Of the HLA homozygotes, a total of 153 are double homozygous ( Figure 1B ; Supplementary Table S2 ): 58 for HLA-A and HLA-B, 76 for HLA-A and HLA-C, and 172 for HLA-B and HLA-C. Remarkably, 51 individuals are triple homozygous for HLA-A, HLA-B, and HLA-C ( Figure 1B ; Table 1 ). Haplotype frequencies of the complete dataset (3,496 individuals) are available in the Supplementary Data .

Comparisons of the Lithuanian Class I HLA frequencies with those reported for the European-American and British populations using linear regression models show strong correlations between the three cohorts ( Figures 2A–C ). The linear regression analyses yielded an average slope of 0.914 for HLA-A, 0.827 for HLA-B, and 0.860 for HLA-C. This indicates the populations closely resemble each other in the composition and prevalence of allele variants. Principal component analysis (PCA) was performed on the genotypes of the Lithuanian population and on genotypes reconstructed from published datasets using Monte Carlo analysis based on reported allele frequencies. The results showed that the Lithuanian population clustered in close proximity to the compared populations ( Figure 2D ). The Euclidean distances between the centroids of the populations were quantified and represented in the PCA and as a heatmap ( Figure 2E ). The distance metrics indicate that the centroid of the Lithuanian population is proximal to the European-American and British populations, with distances of 1.00 and 0.79 relative units, respectively. The British and European-American populations closely resemble each other, with a Euclidean distance of 0.27 relative units. Hardy–Weinberg equilibrium analyses show that some genotypes, including the 10 most frequent allele types, occur at higher frequencies than expected ( Supplementary Data ; Supplementary Figure S1 ).

We stochastically arranged the 3,496 donors and interrogated whether the subset of HLA-A, HLA-B, and HLA-C triple homozygous (51 samples) and double homozygous (153 samples) individuals were compatible with the 3,496 patients ( Figure 3 ). We found that our cohort of triple homozygous individuals matches 60.46% of the Lithuanian population ( Figure 3A ). Likewise, the double homozygous cohort matches 33.32% of the Lithuanian population. In comparison, a randomly selected subset of 153 or 51 samples from the dataset could match only 11.84% ( Figure 3B ) and 4.1% ( Figure 3C ) of the Lithuanian population, respectively. We then evaluated the matching provided by our triple homozygous and double homozygous cohorts to the European-American and British populations. We assessed their immune compatibility with Monte Carlo datasets reconstructed from allele frequencies reported for European-American and British individuals. Remarkably, we found that the 51 triple homozygous samples of our cohort match 13.4% of the British population, while the double homozygous cohort matches 5.2% ( Figure 3D ). Additionally, we found that triple homozygous samples match 7.4% of the European-American population, and double homozygous samples match 3.3% ( Figure 3E ).

We extracted the Cas9-binding site sequences from the HLA alleles present in the Lithuanian population. First, we focused on the analysis of target regions encompassing the gene body, from the 5′UTR to the 3′UTR. We found 1,996 unique target sites in HLA-A alleles, 2,342 unique target sites in HLA-B, and 2,300 unique target sites in HLA-C. We calculated the activity prediction score based on the rule set 1 of nuclease catalytic activity (26). We found that, as in non-hyper polymorphic genes, the activity scores of all HLA alleles are centered in the inactive Q4 quadrant. We show this distribution for the five most frequent alleles of HLA-A, HLA-B, and HLA-C ( Figure 4A ). The potential of HLA gene knockout to modulate immune compatibility is well accepted in the literature. Although pairs of guide RNAs can be used in conjunction to create exon-spanning knockouts, we focused on guide RNAs in exon regions. From the guide RNAs present in the gene body, we found 679 unique target sites in the HLA-A exons of Lithuanian alleles, 698 in HLA-B, and 687 in HLA-C ( Figure 4B ). Since HLA-A, HLA-B, and HLA-C are class I single-span transmembrane proteins ( Figure 4C ), only guide RNAs targeting the ectodomain have the capacity to create knockouts that eliminate plasma membrane expression of HLA genes. We predicted the transmembrane spanning region (25) of the allele sequences and focused on guide RNAs directed to the N-terminus, upstream of the predicted transmembrane domain. We found there are 615 unique target sites in Lithuanian alleles on HLA-A ectodomains, 658 on HLA-B, and 613 on HLA-C ( Figure 4B ). Of those useful for ectodomain targeting, a fraction have predicted activity scores greater than 0.5. These include 54 for HLA-A, 75 for HLA-B, and 66 for HLA-C ( Figure 4B ).

Naturally occurring triple and double homozygous samples are particularly useful for gene engineering approaches as they allow bi-allelic targeting with a single programmable nuclease in a one-step intervention. Next, we modelled the impact of HLA-A and HLA-B knockouts on the immune compatibility of the double and triple homozygous samples when matching them to the Lithuanian population and other European datasets ( Figure 5 ). We included all 51 triple homozygous individuals from our cohort ( Figure 5A ). From the 153 double homozygous individuals identified, we focused on those that are HLA-A and HLA-B double homozygous, comprising seven individuals ( Figure 5B ). The 51 triple homozygous samples, when in an HLA-C-retained (HLA-A and HLA-B double knockout) configuration, match a maximum of 0.9799 of the Lithuanian population ( Figure 5A ). These 51 samples achieve a match of 0.9577 in the European-American population ( Figure 5C ) and 0.9556 in the British population ( Figure 5D ).

Since the triple homozygote individuals identified in this study are immune-compatible with a large fraction of the Lithuanian and other European populations, we sampled these volunteers. The collected dermal fibroblast samples were used to establish biobank stocks and cultures. Primary fibroblast cultures were robustly established for 15 triple homozygotes ( Figure 6A ; Supplementary Table S3 ). PCR products of exons 2 and 3 ( Figure 6B ) display single bands, and Sanger sequencing yields clear chromatograms, both of which are characteristic of homozygous samples ( Figures 6C–E ). Sanger sequencing of exons 2 and 3, which code for the ectodomains of HLA-A, HLA-B, and HLA-C, revealed characteristic residues for each allele. Characteristic amino acids p.F33 and p.R121 were confirmed for HLA-A*02:01:01 ( Figure 6C ), p.Y33 and p.W119 for HLA-B*13:02:01 ( Figure 6D ), and p.D33 and p.L119 for HLA-C*06:02:01 ( Figure 6E ). These findings were consistent for both XY (donor SD9) and XX (donor SD6) individuals with the homozygous haplotype HLA-A02:01:01-HLA-B13:02:01-HLA-C06:02:01 ( Figure 6F ).