BERT2DAb: a pre-trained model for antibody representation based on amino acid sequences and 2D-structure

ABSTRACT

Prior research has generated a vast amount of antibody sequences, which has allowed the pre-training of language models on amino acid sequences to improve the efficiency of antibody screening and optimization. However, compared to those for proteins, there are fewer pre-trained language models available for antibody sequences. Additionally, existing pre-trained models solely rely on embedding representations using amino acids or k-mers, which do not explicitly take into account the role of secondary structure features. Here, we present a new pre-trained model called BERT2DAb. This model incorporates secondary structure information based on self-attention to learn representations of antibody sequences. Our model achieves state-of-the-art performance on three downstream tasks, including two antigen-antibody binding classification tasks (precision: 85.15%/94.86%; recall:87.41%/86.15%) and one antigen-antibody complex mutation binding free energy prediction task (Pearson correlation coefficient: 0.77). Moreover, we propose a novel method to analyze the relationship between attention weights and contact states of pairs of subsequences in tertiary structures. This enhances the interpretability of BERT2DAb. Overall, our model demonstrates strong potential for improving antibody screening and design through downstream applications.

KEYWORDS:

• Affinity prediction
• antibody screening
• antibody sequences
• secondary structure
• pretrained language model

Introduction

The high specificity and strong neutralizing ability of antibodies have made them highly valued in the prevention and treatment of diseases such as tumors and viral infections.^1–3 Compared to standard antibody discovery approaches, computer-based screening methods offer a promising alternative in terms of efficiency and cost-effectiveness.^4–6 Computer-aided screening of neutralizing antibodies typically involves several steps, including structural modeling, affinity prediction, and developability assessment, and this process usually requires multiple iterations.^7,8 Compared to using structures as inputs, the use of sequences can significantly reduce computational requirements and is generally associated with lower data acquisition hurdles for accomplishing the aforementioned tasks. Leveraging an extensive collection of antibody sequences obtained through sequencing efforts performed by other researchers, we can pre-train an antibody language model (ALM) to subsequently enable downstream tasks related to antibody screening based on sequence information.^9–11

Within the realm of protein sequences, several protein language models (PLM) have been developed that leverage general protein datasets (e.g., UniProtKB and BFD) to learn common sequence features and transfer these features to specific tasks such as biophysical property prediction.^12–14 ESM2, ProtTrans, and UniRep are examples of PLM.^13,15,16 Among these models, the first two are based on self-attention and include multiple models with varying parameter sizes, while the last one is a single model based on LSTM. However, the evolution of antibodies relies on gene rearrangement and somatic hypermutation induced by antigens, which is fundamentally distinct from ordinary protein (e.g., enzymes, structural protein, transport protein) evolution.¹⁷ Consequently, antibody sequences may possess unique features, and using PLMs to capture characteristics of antibody sequences may be inappropriate. AbLang, PARA, and AntiBERTy are both ALMs trained on antibody sequence datasets using self-attention.^18–20

While these language models have achieved promising performance in tasks such as structure modeling, there is still room for further optimization.²¹ First, with the exception of AbLang, all other language models use a single model to embed both light and heavy chains. However, significant differences exist between the two chains in terms of sequencing data availability (which affects the amount of available training data), physicochemical properties, and gene expression characteristics (supplementary Fig.S1). As a result, using the same model to learn and represent sequences of both chains may dilute features specific to the light chain. Second, all of the aforementioned language models use amino acids as the fundamental embedding unit for pre-training. This approach may not effectively consider the collective impact of secondary structures on the spatial structure and functionality of antibodies. Additionally, previous studies have revealed the presence of linear motifs of varying lengths within protein sequences, which essentially represent conserved fragments in the sequence.²² By referring to frequency-based vocabulary construction methods used in natural language processing, Asgari and colleagues developed a vocabulary of protein conserved sequence fragments that has been applied to tasks such as motif recognition, toxicity prediction, subcellular localization prediction, and protease prediction, with promising performance outcomes reported.²³

The secondary structure of antibodies is closely related to antibody-specific antigen recognition and stability.^24,25 Firstly, the complementarity determining regions (cCDRs), especially CDR H3, are the most variable regions in antibodies. The secondary structure of these regions is related to the positioning and overall conformation of the antibody binding site. Their variations affect the orientation and exposure of the CDRs, determining whether optimal spatial matching can occur with the antigen, thereby achieving high specificity binding between the antigen and antibody. Secondly, the secondary structure of the framework region is relatively conserved, forming a stable scaffold to ensure the overall stability of the antibody. Past research has shown that mutations can occur in the framework regions of antibodies, resulting in changes in secondary structure and subsequently affecting the stability and functionality of the antibodies.²⁶ In summary, the close association between the secondary structure of antibodies and stability/functionality suggests that we need to consider the impact of the secondary structure as a whole when designing models.

Here, we propose a pair of antibody language models specifically designed for variable (V) regions that utilize self-attention and incorporate secondary structure information. These models undergo separate pre-training on the light and heavy chain sequence data and are collectively referred to as Bidirectional Encoder Representation from Transformers for Antibody Sequences based on Secondary Structure (BERT2DAb), which includes BERT2DAb_H and BERT2DAb_L (Figure 1). We demonstrated the performance of BERT2DAb across several tasks, including classification of binding of mutant trastuzumab to human epidermal growth factor receptor 2 (HER2), classification of binding of multiple antibodies to multiple coronavirus antigens, and prediction of the change in binding free energy (ΔΔG) after antibody mutation (highlighted by the red box in Figure 1e). While deep neural networks can be challenging to interpret due to opaque network parameters, explainability is essential in biomedical research.²⁷ To address this, Vig and colleagues explored the relationship between attention weights and tertiary structure, inspiring our novel approach to explore the interpretability of BERT2DAb.²⁸ We examine the association between the attention weights of pre-trained models and contact states of pairs of subsequences (highlighted by the yellow box in Figure 1e).

Figure 1. Workflow. Steps ‘a’ through ‘e’ correspond to data cleaning (a), sequence splitting (b), subsequence vocabulary construction (c), pre-training (d), and pre-training model evaluation (e), respectively. In step ‘b’, the ‘H’ refers to 310-helix, alpha-helix, and pi-helix, the ‘e’ refers to beta-bridge and beta-strand, and the ‘C’ refers to high curvature loop, beta-turn, and coil. In step ‘C’, the vocabulary training examples are highlighted using a light blue background. The symbol ‘##’ indicates that the subsequence is located within a secondary structure.

Results

BERT2DAb

We trained BERT2DAb_H and BERT2DAb_L separately for 12,1508K and 976.4K steps, respectively. The training process took approximately three weeks. The final losses were 1.09 and 0.69 for BERT2DAb_H and BERT2DAB_L, respectively. Each model consisted of approximately 110 million trainable parameters. The pre-training process was speeded on 8 parallelized A100 GPUs.

We applied t-SNE to perform dimensionality reduction on the antibody sequence embeddings obtained from BERT2DAb and visualized the results. This showcases the model’s ability to extract diverse B-cell features from the sequences (supplementary Fig.S2).²⁹

Results of three downstream tasks

The classifier based on BERT2DAb can distinguish antibodies that bind to HER2 (multiple antibodies and an antigen)

Here, we aimed to evaluate BERT2DAb by using the classifier based on it to screen for trastuzumab mutants that can bind to HER2. Firstly, we trained a classifier using the Tra_dataset and compared it with the CNN classifier used in the original study through testing. We plotted the Receiver Operating Characteristic Curve (ROC) and Precision-Recall Curve (PR Curve) based on the output results (P) of the testing set. As shown in Figure 2, the area under the ROC curve for the classifier based on BERT2DAb was 0.90 and the average precision was 91.0%, while the area under the ROC curve for the classifier used in the original study was 0.91 and the average precision was 83.0%.

Figure 2. ROC and PR curve. The AUC of the ROC and average precision are robust measures of model accuracy and precision. The blue dashed line on the ROC curve represents the ROC for a random classifier.

To compare with other language models, we trained and tested five language models using the balanced_Tra_dataset. As shown in Table 1, the classifier based on BERT2DAb achieved the best performance in terms of F1-score (86.26%), Recall (87.41%), and Accuracy (85.46%), especially Recall, which improved by nearly 1.5% compared to the AbLang model that was also pre-trained using antibody sequences. The classifier based on BERT2DAb achieved the second best performance on precision (85.15%). Moreover, we observed that the results of classifiers based on ALMs (BERT2DAb, AbLang, and AntiBERTy) were superior to those of classifiers based on PLMs (esm2_t30_150M_UR50D and ProtBert).¹⁵

Table 1. Benchmarks results in classifying the binding of trastuzumab mutants to HER2 (%).

	Language Model	F1 score	Precision	Recall	Accuracy
PLM	Esm2_150M	83.86	82.49	85.27	82.87
	ProtBert	84.68	84.09	85.28	83.89
ALM	AbLang	86.03	86.07	85.99	85.42
	AntiBERTy	85.36	84.15	86.61	84.49
	BERT2DAb	86.26	85.15	87.41	85.46

Due to the difficulty in obtaining labeled data, we further evaluated the performance of our model using subsets of balanced_Tra_dataset at 100%, 60%, and 20% (Table 2). The results demonstrated that the classifier based on BERT2DAb achieved high classification accuracy even when trained with limited data, indicating its robustness and efficiency.

Table 2. Performance of BERT2DAb-based classifier with different training data volumes in classifying the binding of trastuzumab mutants to HER2 (%).

Training sequences	F1 score	Precision	Recall	Accuracy
19440	86.33	84.29	88.48	85.37
11664	83.43	81.77	85.16	82.41
3888	81.70	79.01	84.58	80.09

The classifier based on BERT2DAb can screen antibodies that bind to coronaviruses (multiple antibodies and multiple antigens)

To compare with other language models, we also trained and tested five language models using CoV-Ab-Bind (to evaluate the classification performance of the classifier for the binding of multiple antibodies to multiple coronavirus). As shown in Table 3, the classifier based on BERT2DAb achieved the best performance in terms of F1-score (90.30%), Recall (86.15%), and Accuracy (83.54%). Based on F1-score, the performance of the classifier based on ALMs ((88.81%+89.33%+90.30%)/3 = 89.48%) is superior to that of the classifier based on PLMs ((89.78%+87.11%)/2 = 88.45%).

Table 3. Benchmarks results in classifying the binding of antibodies and coronavirus (%).

	Language Model	F1 score	Precision	Recall	Accuracy
PLM	Esm2_150M	89.76	94.11	85.80	82.60
	ProtBert	87.11	96.00	79.73	79.02
ALM	AbLang	88.81	94.00	84.17	81.16
	AntiBERTy	89.33	94.55	84.66	82.04
	BERT2DAb	90.30	94.86	86.15	83.54

The predictor based on BERT2DAb can predict ΔΔG of antigen-antibody complexes

ΔΔG is an important parameter for evaluating the reversibility of chemical reactions, which reflects the energy difference between molecules in different states, such as the free energy change before and after mutation in antigen-antibody complexes.^30,31 As a quantitative indicator, ΔΔG is often used to describe changes in binding free energy in biological macromolecular systems, such as receptor-ligand, enzyme-substrate, and protein complex. Specifically, when ΔΔG value is negative, the reaction is more likely to occur; conversely, when the value is positive, progression of the reaction is difficult and requires external energy input to achieve. In the case of the antigen-antibody complex before and after mutation, if the ΔΔG value is negative, it indicates that the complex is more stable after mutation; conversely, if the ΔΔG value is positive, it indicates that the complex is less stable than the original.

ΔΔG is a more direct indicator to describe the changes in affinity before and after mutation of antigen-antibody complex. Here, we evaluated the performance of the predictor based on BERT2DAb for predicting ΔΔG and compared it with other language models. As shown in Table 4, the ΔΔG prediction based on BERT2DAb demonstrated excellent performance in both datasets (AB-Bind-allMut and AB-Bind-single). For single point mutations, the TopNetTree model,³² which is based on tertiary structure information, achieved an Rp of 0.65 (SOTA), while our model achieved an Rp of 0.77. Moreover, our predictor can handle multiple point mutations (currently few predictors can adapt to this situation) with an Rp of 0.73. Compared with other language-model based predictors that only use sequences as inputs, the BERT2DAb-based predictor has demonstrated significant advantages.

Table 4. Benchmarks results in predicting the ΔΔG.

		LM	RMSE	Rp (p-value)
AB-Bind-SingleMut	PLM	Esm2_150M	4.45	0.55 (1.14× 10^−5)
		ProtBert	2.64	0.52 (4.65× 10^−5)
	ALM	AbLang	4.36	0.61 (7.05× 10^−7)
		AntiBERTy	2.02	0.52 (4.25× 10^−5)
		BERT2DAb	1.14	0.77 (5.51× 10-12)
AB-Bind-allMut	PLM	Esm2_150M	4.44	0.52 (6.42× 10^−6)
		ProtBert	2.15	0.55 (1.27× 10^−6)
	ALM	AbLang	4.48	0.59 (9.99× 10^−8)
		AntiBERTy	6.10	0.56 (5.63× 10^−7)
		BERT2DAb	2.77	0.73 (2.19× 10-12)

The attention weights of BERT2DAb reflect the 3D structure of the antibody

To investigate the relationship between attention weights of token pairs in BERT2DAb and the contact status of antibody sub-sequences, we fitted a binary logistic regression (using linear models to analyze information on attention weights) and further visualized the connection using Contact correlation coefficient.

Firstly, Table 5 shows that the model achieved the highest accuracy of 77.71% ($l_7$) and 73.59% ($l_3$) for heavy chains and light chains, respectively. The performance was better for heavy chains possibly due to the larger pre-training data compared to light chains. Additionally, high hidden layers exhibited higher sensitivity, while low hidden layers demonstrated higher specificity. This could be attributed to the ability of high hidden layers to learn short-range relationships (contacts), whereas the lower hidden layers learned long-range relationships (non-contacts), resulting in different spatial information expression capabilities for each hidden layer. Moreover, we identified the top three attention heads in logistic regression corresponding to each hidden layer, sorted by the size of the regression coefficient (supplementary Table S1).

Table 5. Results of logistic regression(%).

	Layer ($l$)	Accuracy	Sensitivity	Specificity
Heavy chain	1	77.09	51.10	86.81
	2	77.13	49.15	87.59
	3	76.74	46.55	88.03
	4	77.58	46.18	89.32
	5	77.70	47.73	88.90
	6	76.82	47.55	87.76
	7	77.71	50.97	87.71
	8	76.92	50.16	86.92
	9	76.42	49.09	86.64
	10	77.62	51.03	87.57
	11	75.47	54.86	83.17
	12	75.60	52.10	84.39
Light chain	1	63.06	48.86	69.71
	2	73.52	34.25	91.88
	3	73.59	45.21	86.87
	4	73.07	45.38	86.02
	5	72.45	40.25	87.52
	6	73.36	40.61	88.67
	7	71.73	45.45	84.03
	8	71.48	45.77	83.50
	9	71.17	48.95	81.57
	10	68.70	47.15	78.79
	11	62.25	50.73	67.64
	12	56.55	53.12	58.15

Secondly, heatmaps of contact correlation coefficients for different attention heads and hidden layers were generated and displayed in Figure 3. Overall, the contact correlation coefficients of each attention head gradually increased as the threshold increased, with some attention heads exceeding 90% when the threshold was greater than 0.2. This result indicated that in these attention heads, pairs of subsequences were more likely to contact in the tertiary structure if the attention weight of the token pair was large. In order to further demonstrate the relationship between attention weights and subsequence pair contacts, we used Pymol (2.6) to visualize the tertiary structures of three randomly selected antibodies (Figure 4) and annotated the spatial distance and attention weights of pairs of subsequences. The visualization showed that in the same attention head, the attention weight of a contact subsequence pair was greater than that of a non-contact subsequence pair. We believe that the pre-trained model extracts spatial features from sequence pairs of subsequences. In another study, under a threshold of 0.3, the highest “contact correlation coefficient” in TapeBert, ProtAlbert, ProtBert, ProtBert-BFD, and ProtXLNet was 44.7%, 55.7%, 58.5%, 63.2%, and 44.5%, respectively,²⁸ while our study achieved a maximum of 99.99%. We regarded the basic structural unit as a whole, rather than embedded representations of single amino acids, which may better preserve protein topology information or motif features. Consequently, the structural features learned by the pre-training model are more in line with the essence of tertiary structure of antibody molecules.

Figure 3. Contact correlation coefficient. Darker colors correspond to higher proportions in the heatmap, with the horizontal axis representing the attention head and the vertical axis representing the hidden layer. Blank cells indicate no token pair with weight greater than the threshold in the attention head.

Figure 4. Relationship between attention weights and contact states of pairs of subsequences. the subfigures display the structures of heavy chain (upper row) and light chain (lower row) for each antibody. We analyzed BERT2DAb_H through $l_2\mathrm{\_}{h}_8$ and $l_{11}\mathrm{\_}{h}_6$, as well as BERT2DAb_L through $l_2\mathrm{\_}{h}_{12}$ and $l_4\mathrm{\_}{h}_9$. Minimum distances between control, contact, and non-contact subsequences are indicated by red dotted lines, with a distance less than 8 Å considered as contact. In the same attention head, the attention weight (orange box) of a contact subsequence pair is greater than that of a non-contact subsequence pair. The loop highlighted in red corresponds to CDRH3 or CDRL3.

Discussion

We utilized secondary structure information to construct the vocabulary for antibody subsequences and trained the Transformer-based pre-training model BERT2DAb, which comprises the light chain model (BERTA2DAb_L) and heavy chain model (BERTA2DAb_H). The attention weights were statistically connected with the tertiary structure of antibodies, thus demonstrating that our pre-trained model had learned spatial information from sequence data. This demonstrates some interpretability of BERT2DAb, which refers to the information learned by the model from the sequences and preserved in the form of model parameters, capable of being reflected in the structure of antibodies. This indicates that the model indeed learns relevant features, which are similar to the information reflected by structure modeling tools based on self-attention mechanisms, such as AlphaFold2 and IgFold.^21,33 In the future, BERT2DAb could potentially be used for fragment-based structure modeling and fragment-based paratope recognition tasks. Moreover, BERT2DAb accomplished critical tasks in developing neutralizing antibodies, such as antibody-antigen-specific binding classification and ΔΔG prediction.

We found that classifiers or predictors based on ALMs outperformed those based on PLMs in three downstream tasks. We believe that the main reason for this is the significant differences between protein and antibody sequences in terms of variability (antibodies have highly variable CDR sequences for recognizing different antigens, whereas protein sequences are more conserved and stable) and secondary structure composition (antibodies typically use beta sheet and loop structures to form their variable regions for antigen binding, while proteins generally have more regular structures such as alpha helices).³⁴ The incorporation of secondary structure information may enable BERT2DAb to extract additional features (such as the composition of secondary structures) from the sequence, leading to improved performance of the model in downstream tasks.

Moreover, in two binary classification tasks, we found that classifiers based on BERT2DAb achieved higher recall (87.41% and 86.15%) compared to classifiers based on other language models. We attribute this mainly to the incorporation of secondary structure information. The secondary structure of an antibody can affect the conformation of its CDRs and ultimately its binding affinity. Certain types of secondary structures, such as beta sheets and turns, are often present in CDRs and can contribute to antigen recognition and binding. Furthermore, the orientation and flexibility of these secondary structures can affect the geometry of the binding interface and the stability of the antibody-antigen complex.^6,35 By embedding the secondary structure information into the representation, we can highlight its role as a whole.

The superior performance of BERT2DAb in the aforementioned antibody function prediction can be attributed to the incorporation of secondary structure information. To further illustrate this point, we conducted validation experiments on BERT2DAb for Paratope prediction (supplementary Table S2) and secondary structure prediction (supplementary Table S3). For Paratope prediction, BERT2DAb achieved the highest F1-score (69.3%) and, compared to AntiBERTy, Parapred, and ProtBERT, demonstrated a relatively high recall (72.8%) while maintaining a satisfactory precision (66.1%).^15,36,37 The higher recall can be attributed to the fragment embedding representation based on secondary structure, which enables the model to capture spatial features of the fragments and facilitates the identification of segments that can form optimal spatial matches with antigens. These segments serve as the building blocks for the complete antibody binding sites (conformational epitopes). Regarding secondary structure prediction, BERT2DAb outperformed AntiBERTy and AbLang, achieving the highest accuracy (VH: 99.5%, VL: 99.4%).^18,19 These results demonstrate the explicit incorporation of secondary structure information by BERT2DAb, enabling downstream models to effectively capture the features related to variations in secondary structure types. In conclusion, the incorporation of secondary structure information through fragment embedding representations makes BERT2DAb more capable of predicting whether an antibody can tightly bind to an antigen, thereby facilitating antibody function (e.g., binding/non-binding) prediction.

There are some limitations in our study that need to be further investigated in future work. Firstly, although we validated BERT2DAb in various related tasks, it has yet to be applied to tasks such as structure modeling and developability prediction. Therefore, we plan to expand the model’s application scenarios and develop a versatile antibody-specific model. Our model embeds subsequences, making it possible to attempt modular antibody design and three-dimensional structure modeling, or add constraints to complete antibody generation or design. Secondly, we used a simple neural network to accomplish downstream tasks, and further tuning is required to find an output module more suitable for specific tasks. Thirdly, the attention mechanism and relationship that drive productive and functional heavy and light chain pairing are paramount for structural stability, as well as specificity, yet less explored by language models reported so far. This could limit the generalizability of the pre-train models to transfer to diverse and real-world antibody tasks. So, we aim to further pre-training BERT2DAb to elucidate the underlying patterns governing the pairing of light and heavy chains, and learn the inter-chain spatial relationships of light and heavy chains. Fourthly, due to time constraints and limited computational resources, the current model has a parameter count of approximately 110 M. However, in future work, we intend to train even larger models to further enhance model performance.

In conclusion, our pre-training model provides insights into the relationship between amino acid sequences and protein structures, offering a promising approach for antibody optimization and drug development. Our findings suggest that integrating secondary structural information can effectively enhance the performance of pre-training models for antibody-related tasks.

Materials and methods

Data

Antibody sequences for pre-training

We conducted pre-training using human V regions sequences from Observed Antibody Space (OAS), which currently holds the largest collection of antibody sequence data sourced from various studies.³⁸ Following data cleaning, the light (VL, the V region of light chain) and heavy (VH, the V region of heavy chain) chain datasets contained 1.2 billion and 33 million sequences, respectively (Figure 1a).

Datasets for benchmarks

To evaluate the performance of BERT2Db, we fine-tuned the language model on three benchmark datasets related to antibody screening. Firstly, the trastuzumab mutation dataset, which we named Tra_dataset, included 38,839 sequences of heavy chain CDR 3 (CDRH3) mutations of trastuzumab (we supplemented the sequence of heavy chain V region except CDRH3 and the sequence of light chain V region) and their binding outcomes to HER2 in wet lab experiments (‘1’ for binding and ‘0’ for non-binding), with 11,300 (29.09%) of them being binding sequences.⁸ To facilitate baseline comparison with other language models, we randomly sampled a subset of 21,600 sequences (with approximately 50% binding samples) from the Tra_dataset and named it balance_Tra_dataset. Secondly, the Coronavirus antibody dataset was derived from CoV-AbDab and curated as CoV-Ab-Bind, containing 15,912 antibody VH and VL sequences along with their binding outcomes to various strains and subtypes of coronavirus (‘1’ for binding and ‘0’ for non-binding), including SARS-CoV1, SARS-CoV2_Alpha, SARS-CoV2_Omicron-BA1, SARS-CoV2_WT, SARS-CoV2_Omicron-BA2.11, and MERS.³⁹ We represented the antigenic epitopes (S protein RBD or NTD) sequences as a feature vector of length 147 using Pyprotein,⁴⁰ as the antigenic feature representation. Binding samples comprised 88% of CoV-Ab-Bind. Thirdly, the AB-Bind dataset is commonly used as a benchmark for ΔΔG (supplementary Formula 1–2) prediction, consisting of 1,101 data points of measured ΔΔG for antigen-antibody complex mutations.⁴¹We extracted two subsets, namely AB-Bind-allMut (681) and AB-Bind-single (536), where the former contains samples with multi-point mutations while the latter contains only single-point mutations. We extracted the antibody and antigen sequences from the wild-type complex (PDB) and obtained the mutated antibody and antigen sequences based on the mutation sites. Similarly, we represented the antigenic sequences as feature vectors using Pyprotein.⁴⁰ In AB-Bind-allMut and AB-Bind-single, the antibody sequence refers to the VH and VL sequences, while the antigen sequence is the full sequence extracted from PDB (we excluded complexes with multiple peptide chains).

Dataset for explore the relationship between attention weights and subsequence pair contacts

In tertiary structure, the closer the pairs of amino acids are, the more likely they are to make contact.⁴² According to the definition used in the Critical Assessment of protein Structure Prediction (CASP) competition, two amino acids are considered to be in contact if their Cβ atoms are less than 8 Å apart.⁴³ Based on this definition, we assumed that two subsequences were in contact if they contained contacting amino acid pairs (<8 Å).

The Thera-SabDab dataset, which includes 828 antibody sequences, was used to perform an interpretability analysis of BERT2DAb.⁴⁴ Firstly, the light and heavy chains of the antibodies in Thera-SabDab were modeled using IgFold.²¹ The contact relationship between each subsequence and other subsequences was determined based on the predicted atomic coordinates and the previously defined subsequence contact concept. Secondly, we fed the heavy chain and light chain sequences separately into BERT2DAb to obtain attention weights of all token pairs in each attention head of each hidden layer.

Pre-training of BERT2DAb

Introduction of secondary structure information

Secondary structures are represented as continuous segments on the primary sequence, and we split each secondary structure in the antibody sequence by inserting spaces. The following example explains how we introduced the secondary structure information before pre-training (Figure 1b). ‘VSGSPGQSVTISCTGTSSDVGAY’ represents a fragment of an antibody sequence. The proteinUnet can predict the secondary structure for each amino acid in the sequence.⁴⁵ Based on the prediction results, consecutive amino acids belonging to the same secondary structure type were considered as a segment and separated by spaces. The sequence fragment was split into VSGSPGQ, SVTISCT, GTSSDVG, and AY, which were equivalent to ‘words’ in English, so we referred to them as ‘sequence words’. We split the secondary structures for all sequences in the light and heavy chain datasets used for pre-training.

Construction of vocabulary

To make each secondary structure a cohesive token in the embedded representation and also to control the length of the vocabulary, we use the wordpiece algorithm to train separate vocabularies for light and heavy sub-sequences on the segmented secondary structure sequences.^46,47The main steps of wordpiece in this study are described below (Figure 1c).

Initializing the subsequence vocabulary

The initial subsequence vocabulary comprises 20 standard amino acids and identifiers required for pre-training (which are not involved in the iteration process).

Finding the subsequence pair with the highest score

We calculate the score of any two pairs of subsequences in the subsequence vocabulary using formula (1) and add the subsequence pair with the highest score to the vocabulary.

Iterating the second step

We continue iterating the second step until we reach the maximum vocabulary size set in advance. In this study, the subsequence vocabulary size was set to 40,000.

(1) $Score=\frac{freq\mathrm{\_}of\mathrm{\_}pair}{freq\mathrm{\_}of\mathrm{\_}first\mathrm{\_}element\times freq\mathrm{\_}of\mathrm{\_}second\mathrm{\_}element}$(1)

Pre-training

We trained BERT2DAb_H and BERT2DAb_L using Pytorch 3.8 with a masked language model (MLM, as shown in Figure 5) as a self-supervised pre-training task. Both models were based on Bidirectional Encoder Representations from Transformers and had 12 hidden layers, each with 12 attention heads. The hidden layer embedding size was 768, the maximum sequence length was 128 (we only encoded the variable region sequence, which is typically no more than 128 amino acids, and we used subsequences as tokens for embedding. Thus, setting a maximum length of 128 is sufficient). The pre-training hyperparameters are listed in supplementary Table S4.

Figure 5. Masked language model. Each sequence segmented based on the secondary structure was further divided into subsequences (tokens) based on the subsequence vocabulary using the principle of maximum length matching. Additionally, 15% of tokens in each sequence were randomly masked. After passing through the embedding layer, an input vector of shape (128,768) was obtained. Sequence representations (128,768) were generated through 12 transformer layers.

Benchmarks for BERT2DAb

We compared BERT2DAb with four other self-attention-based language models: esm2_t30_150M_UR50D, ProtBert, AbLang, and AntiBERTy, the former two being PLMs and the latter two being ALMs. Their parameters are 150 M, 430 M, 85 M, and 26 M, respectively. The main aim of this study is to explore a new method of incorporating secondary structure information for embedding representations. To control the impact of parameter count on model performance, we selected models for baseline testing that have a parameter count comparable to BERT2DAb. The esm2_t30_150M_UR50D, ProtBert, and AntiBERTy use the same model to represent light chain and heavy chain, so their parameter count will double during fine-tuning. We compared the performance of five language models across three downstream tasks, and in fine-tuning, all language model parameters were set to be trainable.

In the three downstream tasks, the input to the downstream classification or prediction network is as follows: Task one (classification of binding of mutant trastuzumab to HER2): concatenated VH embedding representation vector and VL embedding representation vector. Task 2 (classification of binding of multiple antibodies to multiple coronavirus antigens): concatenated VH embedding representation vector, VL embedding representation vector, and antigen sequence representation vector (147). Task 3 (prediction of ΔΔG after antibody mutation): concatenated wild-type antigen-antibody complex VH embedding representation vector, VL embedding representation vector, and antigen sequence representation vector (147), concatenated with mutated antigen-antibody complex VH embedding representation vector, VL embedding representation vector, and antigen sequence representation vector (147). According to the use of different pre-trained language models, the lengths of the VH and VL representation vectors are different, specifically: 128 (BERT2DAb), 640 (esm2_t30_150M_UR50D), 1024 (ProtBert), 768 (AbLang), and 512 (AntiBERTy). BERT2DAb employs subsequence-based tokenization for embedding, while other language models employ amino acid-based tokenization for embedding. During fine-tuning, all models used the embedding representation of the last hidden layer output.

As our objective was to evaluate the performance of BERT2DAb and compare it with other language models, we used simple neural network layers to process the output of the language models for all three tasks. For task one, we reduced the dimensionality of the antibody sequence embedding representation to 1 using a fully connected layer and used a sigmoid function to output the binding probability (P, where a mutant antibody is considered to bind HER2 when P > .5). For task two, we concatenated the antibody sequence representation with the antigen feature vector and used a fully connected layer to reduce the dimensionality, followed by a Softmax to output the binary classification result. For task three, the embedded representation vector of the wild-type antigen-antibody complex and the embedded representation vector of the mutated antigen-antibody complex were then processed through two fully connected layers to reduce their dimensions, respectively. We computed the prediction of ΔΔG by subtracting the embedding representation of the wild-type complex from that of the mutated complex and summing the resulting value.

In all three downstream tasks, we performed a dataset split (train: 80% train, validation: 10%, test: 10%, supplementary Table S5), and initially conducted hyperparameter search using the train set and validation set. Subsequently, we merged the train set and validation set, trained the models using the best hyperparameter combination obtained, and evaluated their performance on the test set. We used the Adam optimizer to iteratively optimize the parameters of the model. We performed the hyperparameter search with the following hyperparameters being searched: train epoch, batch size, learning rate, Adam weight decay, dropout rate, and warm-up steps ratio. The results of the hyperparameter search are shown in supplementary Table S6. For binary classification tasks (task one and two), we evaluated the models based on F1-score, precision, recall, and accuracy. For regression tasks (task three), we evaluated the models based on Root Mean Square Error (RMSE) and Pearson correlation coefficient (Rp) (supplementary Formula 3–7).

Explore the relationship between attention weights and subsequence pair contacts

We trained BERT2DAb using a multi-head attention mechanism.⁴⁸ For any sequence $S\mathrm{\_}T=\left(t_1,\cdots ,t_L\right)$ with a length of $L$ after tokenization, each attention head h$=\left(h_1,\cdots ,h_{12}\right)$ of each hidden layer $l=\left(l_1,\cdots ,l_{12}\right)$ in the model would generate a set of attention weights for each token, and each $w_{i,j}$（0<$w_{i,j}$ <1，${\sum }_j^L{w}_{ij}=1$） represented the attention weight of $t_i$to $t_j$ ，that was, the degree of correlation between them.

To analyze the relationship between the interaction of pairs of subsequences and the corresponding attention weights of the token pairs, we carried out the steps indicated below.

Binary logistic regression

We conducted a binary logistic regression analysis on the contact relationship of pairs of subsequences (dependent variable) and the attention weights of the corresponding token pairs (independent variable) in 12 attention heads. The purpose was to determine the ability to predict the contact relationship of pairs of subsequences through attention weights, which reflects the model’s ability to learn spatial information of antibodies from sequences during pre-training. We used Accuracy, Sensitivity (predictive ability for contact samples), and Specificity (predictive ability for non-contact samples) to evaluate the performance of logistic regression. The calculation formulas are shown in equations (2)-(4). We obtained the True Positive value (TP), True Negative value (TN), False Positive value (FP), and False Negative value (FN) by comparing the prediction results of the test set with the actual labels.

(2) $Accuracy=\frac{TP+TN}{TP+FN+FP+TN}$(2)

(3) $Sensitivity=\frac{TP}{TP+FN}$(3)

(4) $Specificity=\frac{TN}{TN+FP}$(4)

Contact correlation coefficient

To further investigate the relationship between attention weights and contact states of pairs of subsequences, we defined a ‘contact correlation coefficient’. This coefficient represents the proportion of token pairs with a contact relationship among the token pairs whose attention weight is greater than the threshold (t) in each attention head of each hidden layer. The thresholds used were 0.05, 0.1, 0.2, and 0.3. The contact correlation coefficient was averaged across the sequence level.

Data and code availability

The OAS dataset that support the pretraining of this study is available in https://opig.stats.ox.ac.uk/webapps/oas/.³⁸ The Trastuzumab dataset that support the classification of specific binding of mutant Trastuzumab to HER2 is available from github repository: https://github.com/dahjan/DMS_opt/.⁸ The CoV-AbDab dataset that support the classification of specific binding of multiple antibodies to multiple coronavirus antigens is available in https://opig.stats.ox.ac.uk/webapps/covabdab/.³⁹ The AB-Bind dataset that support the prediction of ΔΔG after antibody mutation is available from github repository: https://github.com/sarahsirin/AB-Bind-Database.⁴¹ The Thera-SabDab dataset that support the analyze the relationship between the attention weights of pre-trained models and the contact states of pairs of subsequences is available in https://opig.stats.ox.ac.uk/webapps/newsabdab/therasabdab/.⁴⁴

The pre-trained model and source data files for downstream task model training and data analyses in this study are provided on https://huggingface.co/w139700701.

The source code and code for analyses in the study are opened on GitHub: https://github.com/Xiaoxiao0606/BERT2DAb.

Author contributions

X.W.L. designed the study, implemented the code, performed the experiments, analyzed the results and wrote the paper. F.T. implemented the code and analyzed the results. W.B.Z. implemented the code and performed the experiments. X.W.Z. implemented the code and analyzed the results. J.Y.L. implemented the code and analyzed the result. J.L. performed the experiments and analyzed the results. D.S.Z. designed and supervised the study, analyzed the results and wrote the paper. All the authors revised the manuscript.

Acknowledgments

Our gratitude goes to the developers of datasets used in this study, including OAS, CoV-AbDab, AB-Bind, Thera-SabDab and the dataset of mutated Trastuzumab. Their excellent work and the public resources enable us to engage in this research.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Supplementary material

Supplemental data for this article can be accessed online at https://doi.org/10.1080/19420862.2023.2285904

Additional information

Funding

The author(s) reported there is no funding associated with the work featured in this article.