The B cell reservoir sequencing data were downloaded from the Observed Antibody Spaces (OAS) as of Oct 05 2022 (23). We filtered the OAS data to include datasets that are from human samples and characterize both the light chain and heavy chain of the B cell receptors. To profile the bulk B cell reservoir, we excluded datasets that only profile certain isotypes (e.g. IgG or IgM). We also excluded samples with low numbers of reads (a total number of reads less than 2000). In total, the filtered dataset contains 3300237 reads (2091286 light chain reads and 1208951 heavy chain reads) from 199 subjects.

The amino acid level sequences of the 93 therapeutic antibodies were obtained from the IMGT monoclonal antibody database (mAb-DB) (24). The clinical ADA rate information was obtained from the immunogenicity section of FDA labels. In cases where multiple ADA rates are available from different clinical trials, the maximum rates reported by monotherapy trials were used. The ADA rate of other antibodies (antibodies withdrawn from the market or investigational antibody drugs) were obtained from published literature, previous FDA labels or the Genentech internal database. Other information about the therapeutic antibodies, including the approval status, the year of the first US approval, the first approved indication and the origin (human or humanized) was obtained from the Antibody Society website.

For the OAS dataset, the V and J alleles were identified by the OAS team based on previously published workflow and were included as part of the OAS dataset (23). Briefly, The nucleotide sequences of the B cell repertoire sequencing data were aligned to germline and translated using IgBLASTn 1.17.1 (25). The VJ germline databases for humans were created using ImMunoGeneTics (IMGT) germline sequences derived from IMGT (24). For the therapeutic antibodies, the V and J alleles were identified by aligning the amino acid sequences of therapeutic antibodies to the amino acid-level sequences of the germline V or J alleles. The amino acid-level sequences of the germline V or J alleles were obtained from the gene database (GENE-DB) of IMGT (24). The germline sequence with the smallest number of mismatches with the antibody sequence was selected. In some cases, multiple alleles could be selected due to a tie in the number of mismatches. These alleles were all considered as the matching alleles in the downstream analysis.

The allele usage of V and J segments were measured by dividing the number of sequences that use the V or J allele by the total number of sequences in the B cell repertoire sequencing dataset. The allele usage was first measured in the B cell repertoire of each individual. The median percentage of the V and J alleles was used to represent the allele usage in the OAS dataset. For therapeutic antibodies whose sequence only has one closest allele, the usage of the allele was used to correlate with the immunogenicity. For therapeutic antibodies whose sequence has more than one closest allele with the same number of mismatches, the percentages of the alleles were summed up.

The gene usage of V and J segments were measured by dividing the number of sequences that use the V or J gene by the total number of sequences in the B cell repertoire sequencing dataset. The gene usage was first measured in the B cell repertoire of each individual. The median percentage of the V and J gene was used to represent the allele usage in the OAS dataset. For therapeutic antibodies whose sequence only has one closest gene, the usage of the gene was used to correlate with the immunogenicity. For therapeutic antibodies whose sequence has more than one closest gene with the same number of mismatches, the percentages of the genes were summed up.

The deletion of V genes are detected using a previously published criteria (26). Briefly, V genes with percentage below 0.1% were set as candidates for deletion in the sample. A binomial test was performed to test if the percentage of the candidate gene was significantly lower than the percentage of the gene in non-deletion-candidate samples. For a given gene, candidate samples with a FDR below 0.01 were marked as deleted. The frequency of deletion was defined as the number of individuals in which a V gene was deleted divided by the total number of individuals in the OAS dataset.

For each non-deleted-gene, an individual can either have one (homozygous) or two (heterozygous) alleles of the gene. However, because of the introduced random mutations, a sequence in B cell repertoire has a low but non-zero probability to be mapped to an allele even if the allele is not present in the genome. We reason that, despite the existence of mutations, the majority of BCR sequences will be matched to its originating alleles. For a given gene, we first calculate the abundance of the alleles by dividing the number of sequences with a certain allele by the number of sequences with the gene that the allele belongs to. The allele is set as the candidate allele if its abundance is greater than 20%. For each individual, If there is only one candidate allele, the person is annotated as homozygous. If two or more alleles are set as candidates, the individual is annotated as heterozygous of the top two alleles. The frequency of the allele is measured as the number of individuals with a certain allele divided by the total number of individuals who carry the gene in the OAS dataset.

Spearman correlation was used to measure the association between clinical ADA rate and variables of interest, including the number of mismatches, the allele usage, gene usage, gene deletion frequencies, and allele frequencies. Spearman correlation was used to measure the correlation of allele usage between different disease cohorts. T-test was used to compare the differences of mismatch numbers between human and humanized antibodies. As a rank-based correlation metrics, spearman correlation does not require data transformation to reduce the effect of outliers. However, we used log transformation to variables that are not normally distributed for data visualization purposes. For data that contain zeros, pseudo counts are added to allow log transformation.

Multivariable linear regression was used to test the association between ADA rate and V gene usage among human-origin antibodies while controlling for other factors that could impact the immunogenicity, including the class of drug target, the disease category (based on drug’s first approved indication) and the year of drug’s first approval. The V gene usage was measured as the log2(percent of sequences with light chain V * percent of sequences with heavy chain V). The purpose of the log transformation was to reduce the effect of outliers and make the data normally distributed.

R (version R-4.1.1) was used for all statistical analysis (27).

It has been long known that the usages of different V and J genes in the B cell receptor (BCR) repertoire are not equal (28, 29). Some V and J genes are more frequently used, both in terms of the number of sequences and in terms of the number of unique clones, than other V and J genes. In addition, V and J genes are polymorphic in the human population (30, 31), which could increase the chance of immunogenicity even if the antibody is humanized or has human-origin. For example, there are 14 amino acid differences between the IGHV4-34*10 and IGHV4-34*11 alleles, a 14.4% sequence difference. An antibody that uses the IGHV4-34*10 allele could potentially be seen as a foreign antigen in patients that only carry the IGHV4-34*11 allele.

We leveraged the BCR repertoire sequencing data from the observed antibody space (OAS) database to quantify the usage of different V and J alleles in the BCR repertoires of 199 individuals (23). We identified the closest somatic V and J alleles of each sequence in the BCR sequencing data and quantified the usage of the alleles as the percent of total sequences. V and J alleles from both heavy and light chains show imbalanced usage, with some alleles used in > 10% sequences and other alleles rarely being detected ( Figures 1A, B , see the accompanying dataset (32)).

The OAS dataset includes not only samples from healthy donors but also samples from patients with different diseases, including COVID-19, POEMS syndrome, celiac disease, HIV infection, Light Chain Amyloidosis and Myasthenia Gravis. Our analysis shows that the usage of V and J genes are largely conserved across different disease cohorts ( Figure 1C ).

We then examined the V and J alleles that are used in the therapeutic antibodies. We analyzed the sequences of 93 human or humanized antibody drugs. Among them, 50 antibodies are humanized and 43 antibodies have human-origin. 71 antibodies are commercially available and 22 antibodies are either investigational or withdrawn from the market. Bispecific antibodies, antibody-drug conjugates and chimeric antibodies are excluded from the current analysis.

We performed sequence alignment between antibody therapeutic sequence and somatic V/J alleles. We identified the closest somatic V and J alleles to each antibody therapeutic sequence. We found that the therapeutic antibodies use a diverse set of V and J alleles (see the accompanying dataset (32)). The sequence alignment also allows us to quantify the number of mismatches between the antibody therapeutic sequence and somatic V/J alleles ( Figure 2A ). Interestingly, we see a clear bi-modal distribution of the total number of mismatches, which is the sum of mismatches in heavy chain J gene, heavy chain V gene, light chain J gene and light chain V gene. Further investigation reveals that the bimodal distribution represents the difference between human-origin antibodies and the humanized antibodies ( Figure 2B ). The total number of mismatches are significantly correlated with the clinical ADA rate, suggesting that the introduced mutations increase the risk of immunogenicity of the therapeutic antibodies ( Figure 2C ).

We examined the relationship between the usage of V/J alleles and the immunogenicity of the therapeutic antibodies. We found that the clinical ADA rate is positively correlated with the usage (measured as the percent of total sequences in OAS dataset) of V alleles in both light chain and heavy chains of the therapeutic antibody ( Figure 3A ). We did not see significant correlation between the ADA rate and J gene usage ( Figure 3A ).

We used the product of V allele usages in light chain and heavy chain to quantify the overall usage of V alleles of a therapeutic antibody. We again see a significant correlation between the V gene usage and the clinical ADA ( Figure 3B ). Most antibodies that have high clinical ADA rates use rare V alleles (alleles with low percent of total sequences in OAS dataset), such as donanemab and alemtuzumab ( Figure 3C ) The antibodies that used popular V alleles (alleles with high percent of total sequences in OAS dataset) tend to have low clinical ADA rate, such as cemiplimab, denosumab and daratumumab ( Figure 3C ).

The humanized antibodies have a higher number of mutations than the human-origin antibodies ( Figure 2B ). The higher number of mutations could affect our analysis in two ways. First, the number of mutations is directly correlated with ADA ( Figure 2C ). Second, the higher number of mutations makes it more difficult to identify the originating somatic V alleles used by the therapeutic antibody. Therefore, we next examined the relationship between ADA rate and V gene usage separately in humanized and human-origin antibodies.

Among human-origin antibodies, the usages of V alleles in light chain and heavy chain, as well as their product, are significantly correlated with clinical ADA ( Figures 4A, C ). We did not see significant correlation between the ADA rate and J gene usage ( Figure 4A ). Among humanized antibodies, we did not see significant correlation between V/J usage and the clinical ADA ( Figures 4B, D ).

We further performed multivariable regression analysis to test the association between ADA rate and V gene usage among human-origin antibodies while controlling for other factors that could impact the immunogenicity, including the class of drug target, the disease category (based on drug’s first approved indication) and the year of drug’s first approval. The result of the multivariable regression again shows a significant association between V gene usage and clinical ADA rate (p value = 0.0038).

To summarize, our analysis showed that V alleles are associated with ADA in human-origin antibodies, but not in humanized antibodies. J alleles are not associated with ADA of human-origin and humanized antibodies.

The usage of V alleles in a patient population is determined by both its availability in the genomes and the gene-level usage in the BCR repertoires. Specifically, the usage of V alleles can be decomposed into three factors: allele variation (30), gene deletion (26) and gene usage (28, 29) ( Figure 5A ). The deletion of V genes refers to the phenomenon that an entire V gene is missing in the genome of some individuals. For example, the IGHV2-70D gene is missing in 57% of individuals and IGHV4-30-2 gene is missing in 17% of individuals (26).

We used the OAS dataset to estimate the usage of the different V genes in BCR repertoires, the frequency of specific alleles of each V gene, and the frequency of deletions of the V genes ( Figures 5B–D , see Methods for details). We next explored how the three factors contribute to the association between ADA rate and V allele usage. We found that the allele frequency and the gene deletion are significantly associated with the ADA rate of therapeutic antibodies ( Figures 5E, F ). In contrast, ADA is not significantly associated with the usage of V genes ( Figure 5G ). The results suggest that the genomic variations (allelic variation and gene deletion) in the population contribute to the immunogenicity of the therapeutic antibodies.