A nanobody multi-property optimization method based on a protein language model

Abstract

The application relates to the technical field of bioinformatics and computational antibody design, and particularly relates to a nanobody multi-attribute optimization method based on a protein language model. The method acquires a wild-type VHH sequence, determines a complementarity determining region and a framework region, establishes a candidate sequence set and a historical elite sequence set, selects a parent sequence to perform CDR fragment crossing, takes unmasked residues of a recombined sequence as a sequence context, obtains an amino acid conditional probability distribution of a masked position through a protein language model, and forms a candidate VHH sequence set; the evaluation results of the binding free energy change, the melting temperature, the naturalness and the human origin are directionally converted and numerically mapped to form a multi-attribute vector; the candidate set and the historical elite set are combined, the next round of candidate set is updated according to the non-dominated relationship and the hyper-volume contribution value, and the target candidate VHH sequence set is output. The application connects the conditional mutation and the recombined sequence context, and retains the candidate sequence of the multi-attribute combination.

Description

Technical Field

[0001] This invention relates to the fields of bioinformatics and computational antibody design technology, and in particular to a method for optimizing the multiple attributes of nanobodies based on a protein language model. Background Technology

[0002] In the fields of bioinformatics and computational antibody design, computational modification of nanobody sequences typically uses wild-type VHH sequences as a basis. This involves mutating, recombinating, or replacing candidate residues in the complementarity-determining region (CDR), followed by screening candidate VHH sequences using sequence models, structural models, or antibody-antigen binding evaluation models. Some existing methods utilize protein language models to predict the amino acid probability distribution at masked positions, while others employ genetic algorithms to perform crossover and mutation on CDR fragments, selecting subsequent candidates based on evaluation results of binding free energy, stability, human origin, or other attributes.

[0003] Current candidate generation processes typically treat CDR fragment crossover and amino acid mutation as independent operations. After crossover forms a new combination of CDR residues, subsequent mutations may still involve random residue substitutions or selection of candidate amino acids based on probability results at various positions in the original sequence, failing to fully utilize the complete sequence context formed by unmasked residues in the recombinant VHH sequence. When the frame region is preserved and multiple CDR fragments recombine, the residue relationships around the masked positions have changed. If the selection of candidate residues is not sufficiently aligned with this sequence context, it can easily lead to the formation of candidate VHH sequences that do not match the current recombinant sequence.

[0004] Furthermore, considering that free energy change, melting temperature, naturalness, and human origin have different numerical ranges and optimization directions, existing schemes often use attribute thresholds, fixed weights, or comprehensive scores for candidate ranking. This approach compresses multiple attributes into a single evaluation value, making it difficult to retain candidates that differ across different attribute combinations. When selection is based solely on the current round's score, candidates formed in previous iterations are easily lost in subsequent screenings. Therefore, it is necessary to address the insufficient connection between the recombinant sequence context and conditional mutation, and the difficulty in continuously updating candidates across different rounds based on the contributions of multiple attribute sets, focusing on CDR candidate generation and multi-attribute candidate selection within the VHH framework constraint. Summary of the Invention

[0005] To address the aforementioned technical problems, this invention provides a method for optimizing multiple attributes of nanobodies based on a protein language model, comprising:

[0006] S100. Obtain wild-type VHH sequences, determine the complementarity-determining region and the frame region, with the complementarity-determining region position as the variable position and the frame region position as the maintained position; establish a candidate VHH sequence set and a historical elite sequence set based on wild-type VHH sequences.

[0007] S200. Select two parent VHH sequences from the candidate VHH sequence set, cross-pair the CDR fragments, and generate a recombinant VHH sequence. Set a mask in the recombinant VHH sequence, use the unmasked residues as the sequence context input to the protein language model, and generate a candidate VHH sequence set according to the amino acid conditional probability distribution of the mask position.

[0008] S300: Obtain the binding free energy change, melting temperature, naturalness and human origin evaluation results of candidate VHH sequences; perform direction transformation and numerical mapping, and use the mapping results as independent components to generate multi-attribute vectors;

[0009] S400: Merge the candidate VHH sequence set with the historical elite sequence set to remove duplicates, and determine the non-dominated candidate set based on the multi-attribute vector; use the difference in hypervolume between the set before and after removing the non-dominated candidate as the hypervolume contribution value, and update the next round of candidate VHH sequence set and historical elite sequence set; if not terminated, return the next round of candidate VHH sequence set as the candidate VHH sequence set to S200, and output the target candidate VHH sequence set when terminated.

[0010] Furthermore, in S100, the wild-type VHH sequence is numbered using the IMGT numbering system, and FR1, CDR1, FR2, CDR2, FR3, CDR3 and FR4 are identified sequentially; the positions of CDR1, CDR2 and CDR3 are written into the variable position set, and the positions of FR1, FR2, FR3 and FR4 are written into the retention position set; the amino acids at the retention positions remain unchanged during subsequent crossover, masking and residue substitution processes.

[0011] Further, in S200, two parent VHH sequences are selected from the candidate VHH sequence set; a crossover position is determined in the same complementarity determination region of the two parent VHH sequences, and the CDR fragments on one side of the crossover position are exchanged; the exchanged CDR fragments, the unexchanged CDR fragments, and the frame region are combined to generate two recombinant VHH sequences, wherein the frame region and sequence length of the recombinant VHH sequences are consistent with the wild-type VHH sequences.

[0012] Furthermore, masks are set at two variable positions in each recombinant VHH sequence; using unmasked residues in the recombinant VHH sequence as sequence context, the protein language model outputs the conditional probabilities of 20 standard amino acids corresponding to the two mask positions; the top two amino acids with the highest conditional probabilities are extracted at each mask position to generate a candidate amino acid set; the amino acids in the two candidate amino acid sets are combined and the mask is replaced.

[0013] Furthermore, the candidate VHH sequences generated by mask substitution are completely deduplicated; the candidate VHH sequences are compared with the corresponding recombinant VHH sequences to verify the sequence length, standard amino acid characters, and frame region residues that retain their positions; and candidate VHH sequences that are duplicated, have changed sequence length, contain non-standard amino acids, have changed frame region residues, or have not undergone amino acid substitution are deleted.

[0014] Furthermore, in S300, the antibody-antigen complex structures corresponding to the candidate VHH sequence and the wild-type VHH sequence are input into the binding free energy change evaluation model to generate the binding free energy change; the candidate VHH sequence is input into the melting temperature evaluation model to generate the predicted melting temperature; and the candidate VHH sequence is input into the VHH naturalness evaluation model and the human origin evaluation model respectively to generate the naturalness evaluation result and the human origin evaluation result.

[0015] Furthermore, based on the optimization directions and numerical ranges of the combined free energy change, predicted melting temperature, naturalness evaluation results, and human origination evaluation results, the Sigmoid function is used for direction transformation and numerical mapping, mapping the four evaluation results to [0, 1]. In accordance with the order of combined free energy change, predicted melting temperature, naturalness evaluation results, and human origination evaluation results, the four mapping results are used as four independent components of the multi-attribute vector.

[0016] Further, in S400, the dominance state between the merged and deduplicated candidate VHH sequences is determined according to the Pareto dominance relationship, the dominated candidate VHH sequences are deleted, and a non-dominated candidate set is generated; a hypervolume reference point is set where each component is not better than the corresponding attribute component of the non-dominated candidate set; the hypervolume of the set containing a non-dominated candidate and the hypervolume of the set after removing the non-dominated candidate are calculated respectively, and the difference between the two is written into the hypervolume contribution value of the non-dominated candidate.

[0017] Furthermore, the non-dominated candidate set is sorted from largest to smallest according to the supervolume contribution value; non-dominated candidates are filtered according to the candidate set capacity, and the filtered non-dominated candidates are simultaneously written into the next round candidate VHH sequence set and the updated historical elite sequence set; the two parent VHH sequences used in the next round of crossover are both selected from the next round candidate VHH sequence set.

[0018] Further, the difference in supervolume between the non-dominated candidate sets of adjacent iterations is calculated to generate a supervolume increment; the iteration is terminated when the number of iterations reaches a preset number of iterations, or when the supervolume increment of multiple consecutive iterations is not greater than a preset increment threshold; when the candidate VHH sequence set is empty in this round, the historical elite sequence set before the update is returned as the candidate VHH sequence set for execution S200; when the iteration is terminated, the non-dominated candidates in the updated historical elite sequence set are written into the target candidate VHH sequence set.

[0019] The innovative aspects of this invention include:

[0020] (1) The VHH framework region is kept constrained, the CDR fragment crossover is coupled with the protein language model conditional prediction; the recombinant VHH sequence after crossover retains the framework region residues and uses the unmasked residues to form a new sequence context, and the protein language model forms the amino acid conditional probability distribution of the mask position and the candidate VHH sequence set accordingly.

[0021] (2) The binding free energy change, melting temperature, naturalness and human origin of the candidate VHH sequence are transformed and numerically mapped respectively, and retained as independent components in the multi-attribute vector; based on the Pareto dominance relationship between candidates and the difference in set hypervolume before and after removing candidates, the corresponding hypervolume contribution value is formed.

[0022] (3) Merge the current round of candidate VHH sequence set with the historical elite sequence set to remove duplicates, and select the next round of candidate VHH sequence set in a unified multi-attribute space; the next round of candidate VHH sequence set is used again as the source of the parent generation VHH sequence, so that the multi-attribute selection result participates in subsequent CDR fragment crossover and conditional mutation. Attached Figure Description

[0023] Figure 1 This is a flowchart illustrating a multi-attribute optimization method for nanobodies based on a protein language model, provided in an embodiment of this application. Detailed Implementation

[0024] Example 1: Refer to Figure 1 This is a flowchart illustrating a multi-attribute optimization method for nanobodies based on a protein language model, provided in an embodiment of the present invention. The flowchart may include at least steps S100-S400:

[0025] S100. Obtain wild-type VHH sequences, determine the complementarity-determining region and the frame region, with the complementarity-determining region position as the variable position and the frame region position as the maintained position; establish a candidate VHH sequence set and a historical elite sequence set based on wild-type VHH sequences.

[0026] S200. Select two parent VHH sequences from the candidate VHH sequence set, cross-pair the CDR fragments, and generate a recombinant VHH sequence. Set a mask in the recombinant VHH sequence, use the unmasked residues as the sequence context input to the protein language model, and generate a candidate VHH sequence set according to the amino acid conditional probability distribution of the mask position.

[0027] S300: Obtain the binding free energy change, melting temperature, naturalness and human origin evaluation results of candidate VHH sequences; perform direction transformation and numerical mapping, and use the mapping results as independent components to generate multi-attribute vectors;

[0028] S400: Merge the candidate VHH sequence set with the historical elite sequence set to remove duplicates, and determine the non-dominated candidate set based on the multi-attribute vector; use the difference in hypervolume between the set before and after removing the non-dominated candidate as the hypervolume contribution value, and update the next round of candidate VHH sequence set and historical elite sequence set; if not terminated, return the next round of candidate VHH sequence set as the candidate VHH sequence set to S200, and output the target candidate VHH sequence set when terminated.

[0029] S100. Obtain wild-type VHH sequences, determine the complementarity-determining region (MDR) and frame region (FR), using the MMR position as the variable position and the FR region position as the maintained position; establish a candidate VHH sequence set and a historical elite sequence set based on the wild-type VHH sequences.

[0030] This implementation is executed by a computing device equipped with a sequence processing program, a protein language model inference program, a multi-attribute evaluation program, and a multi-target selection program. Upon receiving the nanobody optimization task, the sequence processing program first reads the wild-type VHH sequence. VHH stands for Variable Domain of Heavy Chain of Heavy-Chain Antibody. The wild-type VHH sequence is represented by a single-letter amino acid code, sourced from sequencing result files, antibody sequence database exported files, or manually verified sequence files. The sequence file retains only the VHH body, excluding expression tags, linker peptides, and signal peptides. If spaces, newlines, or position numbers exist in the file, the sequence processing program first removes non-sequence characters and then checks whether the remaining characters belong to the 20 standard amino acid characters.

[0031] Before performing region partitioning, the sequence processing program checks the length and conserved residue positions of the wild-type VHH sequence. If non-standard amino acid characters appear in the sequence, or if consecutive missing residues prevent conserved positions from matching, the complementarity determination region (CDR) identification process is not initiated, and a sequence verification anomaly status is output. For duplicate or invisible characters caused by file formatting, the sequence processing program retains the original file, forming a cleaned sequence to be numbered. For residue positions that cannot be distinguished as sequencing deletions or file errors, the current optimization task is stopped, and S100 is retried using the confirmed sequence. This processing ensures that subsequent comparisons of crossover positions, mask positions, and frame regions are all based on the same residue index sequence.

[0032] After sequence verification, the wild-type VHH sequence was numbered using the International Immunogenetics Information System (IMGT). The numbering process, based on the positional relationships of conserved residues in the VHH, divided the sequence into Framework Region 1, Complementarity Determining Region 1, Framework Region 2, Complementarity Determining Region 2, Framework Region 3, Complementarity Determining Region 3, and Framework Region 4, denoted as FR1, CDR1, FR2, CDR2, FR3, CDR3, and FR4, respectively. Here, CDR stands for Complementarity Determining Region, and FR stands for Framework Region. The region division results simultaneously preserved the start and end positions and corresponding residue sequences of each region; subsequent steps continued to use this positional relationship, and renumbering was not performed during iterations.

[0033] The sequence processing program writes the residue positions covered by CDR1, CDR2, and CDR3 into the variable position set and the residue positions covered by FR1, FR2, FR3, and FR4 into the retention position set. The variable position set defines the range of positions where crossovers, masks, and amino acid substitutions occur; the retention position set serves as a benchmark for candidate sequence verification. In any subsequent round of candidate generation, if the residues at the retention positions of a candidate VHH sequence are inconsistent with the wild-type VHH sequence, that candidate will not proceed to multi-attribute evaluation. Because this implementation only performs amino acid substitutions and not residue insertions or deletions, the variable position set and the retention position set maintain the same sequence index throughout all iterations.

[0034] After the region division is completed, the sequence processing program establishes a candidate VHH sequence set and a historical elite sequence set based on wild-type VHH sequences. The candidate VHH sequence set represents the set of sequences available for selection by the parent generation in the current iteration round; the historical elite sequence set represents the retained set formed by the non-dominant candidates from the current round and previous rounds. Before the first iteration begins, both sets use wild-type VHH sequences as the initial sequences. At this time, the multi-attribute vector of the wild-type VHH sequence has not yet been formed. Its combination of free energy change, predicted melting temperature, naturalness, and anthropogenicity evaluation results is obtained in S300 according to the same evaluation path as the candidate VHH sequences, and then a correspondence is established with the sequence.

[0035] In addition to the complete amino acid sequence, the candidate VHH sequence set is associated with its current iteration round and the corresponding region partitioning result. The region partitioning result is not stored repeatedly with the candidate sequences, but is referenced using wild-type VHH sequence identifiers. The historical elite sequence set retains wild-type VHH sequences before the first round of selection, serving as the benchmark for multi-attribute comparisons of candidates in this round. The candidate VHH sequence set, variable position set, and hold-position set formed in S100 enter S200; among them, the variable position set controls CDR fragment crossover and mask position selection, while the hold-position set participates in recombinant sequence assembly and candidate validity verification.

[0036] S200. Select two parent VHH sequences from the candidate VHH sequence set, cross-pair the CDR fragments, and generate a recombinant VHH sequence; set a mask in the recombinant VHH sequence, use the unmasked residues as sequence context input to the protein language model, and generate a candidate VHH sequence set based on the conditional probability distribution of amino acids at the masked positions:

[0037] S200 receives the candidate VHH sequence set, variable position set, and preserved position set formed by S100, and executes them sequentially by the crossover recombination program and the protein language model inference program. The crossover recombination program first selects two parent VHH sequences from the current candidate VHH sequence set. In one implementation, parent selection uses uniform sampling with replacement; when the first round of candidate VHH sequence set contains only wild-type VHH sequences, these sequences are sampled twice to form two parent inputs. The first round of crossover does not introduce sequence differences between the two identical parents; candidate differences are formed by subsequent mask replacement. In the second and subsequent rounds, the parents are selected from the next round of candidate VHH sequence set formed by the previous round of overvolume contribution screening.

[0038] The two parent VHH sequences were indexed using a unified position index established by S100. Since each candidate VHH sequence involves residue substitutions only at variable positions and has the same sequence length as the wild-type VHH sequence, the CDR boundaries and frame regions of the two parents can be directly aligned. The crossover procedure selects a complementarity-determining region (CDR) from CDR1, CDR2, or CDR3 and determines the crossover position within that CDR. The crossover position does not fall within the frame region and does not alter the length of the CDR. If the selected CDR contains only one residue or cannot form a crossover position, another CDR is selected.

[0039] After determining the crossover position, the crossover recombination procedure swaps the CDR fragments of the two parent VHH sequences located on one side of the crossover position. The first recombinant CDR fragment is formed by combining the CDR fragment before the first parent crossover position and the CDR fragment after the second parent crossover position; the second recombinant CDR fragment is formed using two fragments from opposite sources. The CDR fragments that did not participate in the crossover retain the residue arrangement of their respective parents, while FR1, FR2, FR3, and FR4 retain the frame region residues of the wild-type VHH sequence. The crossover recombination procedure then combines the recombinant CDR fragment, the non-crossover CDR fragment, and the frame region in their original positions to form two recombinant VHH sequences.

[0040] After assembly, the sequence lengths and retention positions of the two recombinant VHH sequences are checked. If the length of the recombinant sequence is inconsistent with that of the wild-type VHH sequence, it indicates a conflict in the boundary index of the cross fragment, and the recombinant sequence is removed from further processing. If the retention position shows that the parent CDR residues are mistakenly written into the frame region, the corresponding frame region residues are restored according to the retention position set of S100, and the cross boundary is rechecked. Only recombinant VHH sequences with the same length, the same retention position, and a clear CDR fragment origin are included in the masking process.

[0041] For each recombinant VHH sequence that passes the assembly check, the masking procedure selects two positions without replacement from the variable position set. The two mask positions can be located in the same complementarity-determining region (CDR) or in two separate CDRs. The selection of mask positions is associated with the corresponding recombinant VHH sequence to avoid mismatches between the model output and the sequence positions. If the randomly selected position shows no residue differences in either parent sequence, the mask position is retained because the protein language model subsequently recalculates the amino acid distribution at that position based on the complete recombinant sequence context, rather than directly using the parent residues.

[0042] The masking process replaces the amino acids at two selected positions with mask characters defined by the protein language model, forming a masked VHH sequence. Except for the two masked positions, the remaining residues in the recombinant VHH sequence remain unchanged; these unmasked residues together constitute the sequence context. After reading the entire masked VHH sequence, the protein language model outputs the conditional probabilities of the 20 standard amino acids corresponding to the two masked positions. These conditional probabilities are not calculated independently based on individual local fragments, but rather within the complete sequence context formed by the crossed CDR fragments, non-crossed CDR fragments, and the frame region. Therefore, when the parent CDR fragment combination changes, the conditional probability distribution at the same masked position also changes with the input sequence.

[0043] The protein language model is loaded into the model inference program before the start of this optimization task, and the model parameters are not changed during iteration. In one implementation, a masked protein language model pre-trained with natural protein sequences is used. The model receives single-letter amino acid sequences and mask characters, and outputs a 20-dimensional probability vector for each mask position. The model inference program removes the output of the start marker, stop marker, and the mask marker itself, retaining only the standard amino acid probabilities; then, the 20 probability values ​​are normalized so that their sum is 1. If there are non-numerical results, missing components, or the normalized sum still does not conform to the model output format, the model inference is re-executed; if an anomaly occurs again, the masked VHH sequence is not considered a candidate result.

[0044] For each mask position, the model inference program arranges 20 standard amino acids in descending order of conditional probability, extracts the top two amino acids, and forms a candidate amino acid set for that position. Each mask position has two candidate amino acids. The combination program performs Cartesian combination on the two candidate amino acid sets, resulting in four residue combinations. Each residue combination replaces the corresponding two mask positions, forming four candidate VHH sequences. When the two recombinant VHH sequences undergo the same processing, the number of candidates formed by a single parent crossover does not exceed eight; the actual number is determined by subsequent deduplication and verification results.

[0045] After candidate sequences are formed, complete deduplication is performed based on the complete amino acid sequence. Complete deduplication is based on position-by-position consistency throughout the entire sequence, and repetitive sequences are not retained due to differences in parental source, crossover positions, or mask positions. For each candidate VHH sequence after deduplication, the sequence processing program compares it with the corresponding recombinant VHH sequence and the wild-type VHH sequence. The comparison with the recombinant VHH sequence verifies whether amino acid substitutions have actually occurred at the two mask positions; the comparison with the wild-type VHH sequence verifies the sequence length, standard amino acid characters, and frame region residues at retention positions.

[0046] A candidate VHH sequence is deleted when its length changes, it contains non-standard amino acids, the residues in the frame region at the retention position change, or the entire sequence is identical to the corresponding recombinant VHH sequence after two mask positions are replaced. If one of the two mask positions retains the original residues while the other is replaced, and the entire candidate sequence still shows actual amino acid changes, the candidate is retained. The verified candidate VHH sequences constitute the current round of candidate VHH sequence sets, retaining their parent sequence identifiers, crossover positions, and mask positions for subsequent verification of sequence origin.

[0047] If the current candidate VHH sequence set is empty, the model inference program does not call the four-attribute evaluation model in S300, but instead passes the empty set state to S400. The iterative control process then uses the historical elite sequence set before the update as the candidate VHH sequence set, and re-executes parent selection, CDR fragment crossover, and mask position selection. In one implementation, the random sampling state of the mask position is changed during re-execution to avoid continuously using the same parent and the same mask position. If the current candidate VHH sequence set is not empty, the complete sequences in the set and their residue differences from wild-type VHH sequences enter S300.

[0048] S300: Obtain the binding free energy change, melting temperature, and naturalness / human origin evaluation results of candidate VHH sequences; perform direction transformation and numerical mapping, and use the mapping results as independent components to generate multi-attribute vectors.

[0049] S300 is executed by a multi-attribute evaluation program. The input includes the set of candidate VHH sequences generated in S200, the wild-type VHH sequence, and the antibody-antigen complex structure corresponding to the wild-type VHH. The antibody-antigen complex structure is represented in a three-dimensional coordinate file, which includes the VHH chain, the antigen chain, and residue numbers. After loading the structure file, the multi-attribute evaluation program first maps the position indices of the wild-type VHH sequence to the VHH residue numbers in the structure file. For missing frame region residues or CDR residues in the structure, if the missing position does not affect the structural mapping of the candidate mutation position, the existing structure is retained for evaluation. If the candidate mutation position itself lacks spatial coordinates, the evaluation status of the candidate's binding free energy change is recorded as incomplete.

[0050] The binding free energy change evaluation uses the wild-type antibody-antigen complex structure as the structural benchmark. A multi-attribute evaluation program reads the mutation positions of the candidate VHH sequence relative to the wild-type VHH sequence and maps the corresponding amino acid substitutions to the VHH chain in the complex structure. A structural mutagenesis program updates the side chain conformation at the mutation position according to the type of the substituted side chain and performs local conformational adjustments on residues adjacent to the mutation position to form the antibody-antigen complex structure corresponding to the candidate VHH sequence. The candidate complex structure and the wild-type complex structure are then input into the binding free energy change evaluation model to obtain the difference in antibody-antigen binding free energy before and after the candidate mutation.

[0051] In this implementation, the change in binding free energy is denoted as ΔΔG, with units of kcal / mol, and is calculated in the order of candidate complex binding free energy minus wild-type complex binding free energy. The smaller the ΔΔG value, the better the binding state of the candidate complex under this evaluation model. If the candidate sequence contains multiple mutation sites, the binding free energy change evaluation model calculates the change in binding free energy as a whole for multiple mutation sites, and does not directly add the evaluation results of multiple unit points. If the structure mapping fails, the local conformation adjustment does not converge, or the model does not return a value, the candidate will not be included in the current round of non-dominated relation calculation; the other three attributes can still be retained in the evaluation state, but will not form a complete multi-attribute vector.

[0052] The predicted melting temperature is calculated by a melting temperature evaluation model. This model receives the complete amino acid sequence of a candidate VHH and outputs the predicted melting temperature in °C. Before inputting the sequence into the model, the sequence processing program checks the length of the candidate sequence against the input range supported by the model during training; candidates exceeding the supported range are not used by the model. Since the candidates in this embodiment only undergo residue substitution, the sequence length is always consistent with the wild-type VHH sequence, and under normal circumstances, there will be no model input conflict due to length changes. The higher the predicted melting temperature value, the better the attribute value of the candidate VHH in this evaluation dimension.

[0053] The VHH naturalness assessment model receives the complete amino acid sequence of a candidate VHH and outputs a naturalness assessment result based on the distribution of natural Camelidae VHH sequences. This result reflects the similarity between the candidate sequence and the distribution of residue combinations in natural VHHs. The humanization assessment model outputs a humanization assessment result based on the distribution of the variable region sequence of human antibodies. The two assessment results are retained separately and not merged into the same homology value. For the same candidate VHH sequence, the input sequences of the two models are the same, but the sequence distributions referenced by the models are different. Therefore, the naturalness assessment result and the humanization assessment result occupy different components in the subsequent multi-attribute vector.

[0054] In the first round of multi-attribute evaluation, wild-type VHH sequences obtained predicted melting temperature, naturalness evaluation results, and human origin evaluation results following the same path as the candidate VHH sequences in this round. The binding free energy change corresponding to wild-type VHH was set to 0, or a value of 0 was obtained by benchmarking the wild-type complex with itself using the binding free energy change evaluation model. A correspondence was established between the four evaluation results and the wild-type VHH sequences, allowing the first round of S400 to use wild-type VHH as historical elite candidates for non-dominance judgment. In subsequent rounds, if a historical elite sequence already has complete attribute evaluation results, the same version of the evaluation model will not be used repeatedly.

[0055] Since the bound free energy change, predicted melting temperature, naturalness evaluation results, and anthropogenicity evaluation results have different numerical ranges and different optimization directions, the multi-attribute evaluation program performs Sigmoid function transformation on the four results respectively. For the bound free energy change x1, which indicates a better attribute with a smaller value, the mapping value y1 is calculated according to y1 = 1 / [1 + exp(a1(x1 - b1))]. For the predicted melting temperature x2, naturalness evaluation result x3, and anthropogenicity evaluation result x4, which indicate a better attribute with a larger value, the mapping value yⱼ is calculated according to yⱼ = 1 / [1 + exp(-aⱼ(xⱼ - bⱼ))], where j takes values ​​of 2, 3, and 4, a1 to a4 are positive scaling parameters, and b1 to b4 are the center parameters of the corresponding attributes.

[0056] The scaling and center parameters are fixed before the start of each optimization task. The center parameter is set based on the calibration sample distribution of the evaluation model or a preset attribute reference value, while the scaling parameter is set based on the effective numerical range of the corresponding attribute. During iteration, the parameters are not refitted based on the current candidate set to avoid the same original attribute value obtaining different mapping results in different rounds. If an original evaluation result exceeds the effective numerical range declared by the evaluation model, the evaluation status of that candidate attribute is recorded as abnormal, and the extrapolated mapping value is not included in the non-dominated judgment.

[0057] After direction transformation and numerical mapping, all four evaluation results fall within the range [0, 1]. The multi-attribute evaluation program, following a fixed order of combining free energy change, predicted melting temperature, naturalness evaluation results, and anthropogenicity evaluation results, treats y1, y2, y3, and y4 as four independent components to form a multi-attribute vector Y = [y1, y2, y3, y4] for candidate VHH sequences. The four components are not summed with fixed weights, nor are they compressed into a single comprehensive score. Each candidate VHH sequence only forms a complete multi-attribute vector after all four evaluations have obtained valid results; candidates with incomplete attribute evaluations are removed from the candidate set merged in this round, but their sequences and abnormal states are retained until the end of this round and do not participate in hypervolume calculations.

[0058] S300 outputs the correspondence between the current round of candidate VHH sequences and multi-attribute vectors, and also passes the existing multi-attribute vectors of wild-type VHH sequences or historical elite sequences to S400. S400 no longer calls the attribute evaluation model, but directly compares the dominance relationships between candidates based on the four independent components. The candidate sequence identifier serves as the correspondence index between the sequence and the multi-attribute vector; parent source, crossover position, and mask position do not participate in the Pareto dominance calculation.

[0059] S400: Merge the candidate VHH sequence set with the historical elite sequence set to remove duplicates, and determine the non-dominated candidate set based on the multi-attribute vector; use the difference in hypervolume between the sets before and after removing non-dominated candidates as the hypervolume contribution value, and update the next round of candidate VHH sequence set and historical elite sequence set; if not terminated, return the next round of candidate VHH sequence set as the candidate VHH sequence set to S200; if terminated, output the target candidate VHH sequence set.

[0060] S400 is executed by a multi-objective selection program and an iterative control program. The multi-objective selection program receives the current round of candidate VHH sequences and their multi-attribute vectors output by S300, and reads candidates that have already undergone multi-attribute evaluation from the historical elite sequence set. The two types of candidates are merged and deduplicated according to the complete amino acid sequence; when the same sequence appears in both the current round of candidates and the historical elite sequence set, only one sequence and its complete multi-attribute vector are retained. If the same sequence produces two different evaluation results under the same model configuration, it will not directly participate in the current round of selection. First, the evaluation model configuration and input structure are checked, and the re-evaluated result is used to replace the conflicting result.

[0061] After merging and deduplication, the multi-objective selection program compares the Pareto dominance relationships between candidate VHH sequences based on the four components of the multi-attribute vector. For candidate A and candidate B, if all four components of candidate A are not less than the corresponding components of candidate B, and one of its components is strictly greater than the corresponding component of candidate B, then candidate A dominates candidate B. If the two candidates dominate on different components, they do not dominate each other. The program iterates through the merged candidate set, deletes sequences dominated by the other candidate, and the remaining sequences form the non-dominated candidate set.

[0062] After the non-dominated candidate set is formed, a hypervolume reference point is set. Since the four mapping results all indicate that the larger the value, the better the attribute, each component of the reference point is set to be no greater than the minimum value of the non-dominated candidate set on the corresponding component. In one implementation, when all four attributes are mapped to [0, 1] and the mapping result of the valid candidate is greater than 0, the hypervolume reference point is set to r = [0, 0, 0, 0]. If the valid mapping result of a certain attribute is equal to 0, a fixed reference value less than 0 is used on that attribute component, and the same reference point is maintained throughout all iterations.

[0063] Using the reference point r and the multi-attribute vectors of the non-dominated candidate p as relative vertices, a four-dimensional axis-aligned region can be formed. The hypervolume of the non-dominated candidate set is the measure of the union of the four-dimensional regions corresponding to each candidate, denoted as HV(P, r), where P represents the non-dominated candidate set. The hypervolume calculation is performed using the same numerical precision and the same reference point, without changing the calculation boundaries between candidates. If a multi-attribute vector contains missing components or non-numerical components, the candidate has been removed from the merged candidate set before the hypervolume calculation.

[0064] For each candidate p in the non-dominated candidate set, the multi-objective selection procedure first calculates the hypervolume HV(P, r) of the complete non-dominated candidate set, then removes candidate p from set P and calculates the hypervolume HV(P - {p}, r) of the remaining set. The difference between the two, C(p) = HV(P, r) - HV(P - {p}, r), is taken as the hypervolume contribution value of candidate p. This value reflects the independent contribution of candidate p to the area covered by the current non-dominated candidate set, rather than the sum of the four attribute components of candidate p. Candidates located within areas already covered by other candidates may have a hypervolume contribution value of 0; candidates located at the boundary of the non-dominated candidate set have their contribution value determined by their newly added coverage area.

[0065] The multi-objective selection program ranks non-dominated candidates from largest to smallest according to their hypervolume contribution value and reads a pre-set candidate set capacity M. When the number of non-dominated candidates is not greater than M, all non-dominated candidates are written into the candidate VHH sequence set for the next round; when the number of non-dominated candidates is greater than M, the candidates ranked in the top M positions are selected. If two candidates have the same hypervolume contribution value, in one implementation, the binding free energy change mapping component, the predicted melting temperature mapping component, the naturalness mapping component, and the anthropogenic mapping component are compared in sequence; if the four components are still the same, the candidate is selected according to a fixed order of candidate sequence identifiers to avoid the formation of different candidate sets when the same input is run repeatedly.

[0066] The non-dominated candidates selected are simultaneously written into both the next round's candidate VHH sequence set and the updated historical elite sequence set. The historical elite sequence set does not retain the dominated candidates from this round, nor does it replace the entire set with the maximum value of a single attribute. Both parent VHH sequences required for the next round of crossover are selected from the next round's candidate VHH sequence set; therefore, the set selection result in S400 directly changes the next round's parent fragment combination. After the parent fragment combination changes, the recombinant VHH sequences in S200 change accordingly, and the unmasked sequence context received by the protein language model also changes.

[0067] After completing the current round of candidate selection, the iterative control program calculates the difference between the current non-dominated candidate set hypervolume and the previous round's non-dominated candidate set hypervolume, forming the hypervolume increment. If there is no hypervolume from the previous round in the first iteration, the first round's hypervolume is used as the baseline value, and the termination comparison for consecutive rounds is not performed. Subsequent rounds save the hypervolume value and hypervolume increment, respectively, and the termination condition is jointly defined by the preset iteration round and the preset increment threshold. When the preset iteration round is reached, or when the hypervolume increment of multiple consecutive iteration rounds is not greater than the preset increment threshold, the iterative control program terminates the candidate generation process.

[0068] If the termination condition has not yet been met, the iterative control program returns the next round of candidate VHH sequence set as the new candidate VHH sequence set to S200. If there is only one sequence in the next round of candidate VHH sequence set, S200 obtains two parent inputs with replacement; if there are multiple sequences in the set, two parent VHH sequences are selected from the set for CDR segment crossover. The historical elite sequence set is retained when returning to S200, and in the next round S400, it is merged with the newly generated candidate VHH sequence set again.

[0069] If S200 forms an empty set of candidate VHH sequences after candidate deduplication and validity verification, the iterative control program does not update the historical elite sequence set, nor does it replace the current candidate state with an empty set. In this case, the historical elite sequence set before the update is returned to S200 as the candidate VHH sequence set, and parent selection, crossover, and masking processes are re-executed. When consecutive empty sets occur, the historical elite sequence set and the previous round's overvolume value are retained; empty sets are not considered normal iteration rounds with a zero overvolume increment.

[0070] Upon termination of the iteration, the non-dominated candidates in the updated historical elite sequence set constitute the target candidate VHH sequence set. Each target candidate VHH sequence maintains a correspondence with the complete amino acid sequence, mutation location, four original evaluation results, multi-attribute vector, and hypervolume contribution value. The target candidate VHH sequence set is output as a sequence file and an evaluation result file for subsequent antibody-antigen complex structure verification, coding nucleic acid design, or experimental validation. Subsequent steps do not alter the non-dominated state and multi-attribute evaluation results already established in this implementation method.

Claims

1. A method for optimizing multiple attributes of nanobodies based on a protein language model, characterized in that, include: S100. Obtain wild-type VHH sequences, determine the complementarity-determining region and the frame region, with the complementarity-determining region position as the variable position and the frame region position as the maintained position; establish a candidate VHH sequence set and a historical elite sequence set based on wild-type VHH sequences. S200. Select two parent VHH sequences from the candidate VHH sequence set, cross-pair the CDR fragments, and generate a recombinant VHH sequence. Set a mask in the recombinant VHH sequence, use the unmasked residues as the sequence context input to the protein language model, and generate a candidate VHH sequence set according to the amino acid conditional probability distribution of the mask position. S300: Obtain the binding free energy change, melting temperature, naturalness and human origin evaluation results of candidate VHH sequences; perform direction transformation and numerical mapping, and use the mapping results as independent components to generate multi-attribute vectors; S400. Merge the candidate VHH sequence set and the historical elite sequence set to remove duplicates, and determine the non-dominated candidate set based on the multi-attribute vector; use the difference in hypervolume between the set before and after removing the non-dominated candidate as the hypervolume contribution value, and update the candidate VHH sequence set and the historical elite sequence set for the next round. If the process is not terminated, the next round of candidate VHH sequence set will be returned as the candidate VHH sequence set to S200. If the process is terminated, the target candidate VHH sequence set will be output.

2. The method according to claim 1, characterized in that, In S100, the wild-type VHH sequence is numbered using the IMGT numbering system, and FR1, CDR1, FR2, CDR2, FR3, CDR3 and FR4 are identified sequentially; the positions of CDR1, CDR2 and CDR3 are written into the variable position set, and the positions of FR1, FR2, FR3 and FR4 are written into the retention position set; the amino acids at the retention positions remain unchanged during subsequent crossover, masking and residue substitution processes.

3. The method according to claim 1, characterized in that, In step S200, two parent VHH sequences are selected from the candidate VHH sequence set; a crossover position is determined in the same complementary decision region of the two parent VHH sequences, and the CDR fragments on one side of the crossover position are swapped; the swapped CDR fragments, the unswapped CDR fragments, and the frame region are combined to generate two recombinant VHH sequences, wherein the frame region and sequence length of the recombinant VHH sequences are consistent with the wild-type VHH sequences.

4. The method according to claim 3, characterized in that, Masks are set at two variable positions in each recombinant VHH sequence; using unmasked residues in the recombinant VHH sequence as sequence context, the protein language model outputs the conditional probabilities of 20 standard amino acids corresponding to the two mask positions; the top two amino acids with the highest conditional probabilities are extracted at each mask position to generate a candidate amino acid set; the amino acids in the two candidate amino acid sets are combined and the mask is replaced.

5. The method according to claim 4, characterized in that, Complete deduplication is performed on the candidate VHH sequences generated by mask substitution; the candidate VHH sequences are compared with the corresponding recombinant VHH sequences to verify the sequence length, standard amino acid characters, and frame region residues that retain their positions; candidate VHH sequences that are duplicated, have changed sequence length, contain non-standard amino acids, have changed frame region residues, or have not undergone amino acid substitution are deleted.

6. The method according to claim 1, characterized in that, In S300, the antibody-antigen complex structures corresponding to the candidate VHH sequence and the wild-type VHH sequence are input into the binding free energy change evaluation model to generate the binding free energy change; the candidate VHH sequence is input into the melting temperature evaluation model to generate the predicted melting temperature; the candidate VHH sequence is input into the VHH naturalness evaluation model and the human origin evaluation model respectively to generate the naturalness evaluation result and the human origin evaluation result.

7. The method according to claim 6, characterized in that, Based on the optimization directions and numerical ranges of the combined free energy change, predicted melting temperature, naturalness evaluation results, and human origin evaluation results, the Sigmoid function is used for direction transformation and numerical mapping, mapping the four evaluation results to [0, 1]. In the order of combined free energy change, predicted melting temperature, naturalness evaluation results, and human origin evaluation results, the four mapping results are used as four independent components of the multi-attribute vector.

8. The method according to claim 7, characterized in that, In S400, the dominance status between candidate VHH sequences after merging and deduplication is determined according to the Pareto dominance relationship, the dominated candidate VHH sequences are deleted, and a non-dominated candidate set is generated. Set a hypervolume reference point where each component is no better than the corresponding attribute component of the non-dominated candidate set; Calculate the set hypervolume when a non-dominated candidate is included and the set hypervolume after removing the non-dominated candidate, and write the difference between the two into the hypervolume contribution value of the non-dominated candidate.

9. The method according to claim 8, characterized in that, The non-dominated candidate set is sorted from largest to smallest according to the supervolume contribution value; non-dominated candidates are filtered according to the capacity of the candidate set, and the filtered non-dominated candidates are simultaneously written into the next round candidate VHH sequence set and the updated historical elite sequence set; the two parent VHH sequences used in the next round of crossover are both selected from the next round candidate VHH sequence set.

10. The method according to claim 9, characterized in that, Calculate the difference in hypervolume between the non-dominated candidate sets in adjacent iterations to generate hypervolume increments; terminate the iteration when the number of iterations reaches a preset number of iterations, or when the hypervolume increments of multiple consecutive iterations are not greater than a preset increment threshold. When the current candidate VHH sequence set is empty, the previous historical elite sequence set is returned as the candidate VHH sequence set for execution S200; when the iteration terminates, the non-dominant candidates in the updated historical elite sequence set are written into the target candidate VHH sequence set.

Images