In silico prediction of post-translational modifications in therapeutic antibodies

ABSTRACT

Monoclonal antibodies are susceptible to chemical and enzymatic modifications during manufacturing, storage, and shipping. Deamidation, isomerization, and oxidation can compromise the potency, efficacy, and safety of therapeutic antibodies. Recently, in silico tools have been used to identify liable residues and engineer antibodies with better chemical stability. Computational approaches for predicting deamidation, isomerization, oxidation, glycation, carbonylation, sulfation, and hydroxylation are reviewed here. Although liable motifs have been used to improve the chemical stability of antibodies, the accuracy of in silico predictions can be improved using machine learning and molecular dynamic simulations. In addition, there are opportunities to improve predictions for specific stress conditions, develop in silico prediction of novel modifications in antibodies, and predict the impact of modifications on physical stability and antigen-binding.

KEYWORDS:

• In silico prediction
• post-translational modifications
• chemical stability
• developability
• therapeutic antibody development

Introduction

Monoclonal antibodies (mAbs) have become one of the predominant classes of therapeutic proteins and they are used to treat various diseases, such as breast cancer, multiple sclerosis, and asthma.¹ Therapeutic antibodies can be engineered to have high specificity and affinity for their targets.² However, physical and chemical instabilities can have a negative impact on the manufacturability, safety, and efficacy of therapeutic antibodies.³ Like all proteins, mAbs are susceptible to chemical degradation (e.g., oxidation)⁴ and enzymatic modifications (e.g., sulfation) during cell culture.⁵ Chemical and enzymatic modifications contribute to heterogeneity. For example, asparagine (Asn) deamidation can generate charge variants; tryptophan (Trp) oxidation can generate hydrophilic or hydrophobic variants.⁶ In addition, chemical modifications can affect the physical stability and biological activity of antibodies. For example, isomerization within the Fab region can reduce conformational stability,⁷ whereas deamidation within complementary-determining region (CDR) loops can reduce binding affinity.⁸ Therefore, identifying liable sites for post-translational modifications (PTMs) has become a critical step in assessing the developability of therapeutic candidates.

Developability assessments aim to identify candidates with long-term stability, manufacturability, and low heterogeneity.⁹ Forced degradation with thermal, pH, and light stress has been used to accelerate chemical degradation and identify liable residues.¹⁰ Peptide mapping can identify the specific sites for chemical modifications after forced degradation. In contrast, chromatographic techniques such as cation exchange chromatography and hydrophobic interaction chromatography can monitor the overall change in charge and hydrophilic variants, respectively.¹¹ However, experimental approaches for identifying PTM liabilities are time-consuming and require high quantities of the purified protein.¹² The sample preparation and data analysis for peptide mapping is incredibly labor-intensive.¹³ At earlier stages of drug development, the number of forced degradation conditions is limited by the low availability of the purified protein.¹⁰ Computational tools are becoming more common during developability assessments due to the low cost, lack of sample consumption, and high speed. In the past decade, computational tools have been used to predict PTM liable sites and engineer antibodies with better chemical stability.¹⁴

Computational approaches for predicting PTMs

Computational approaches for predicting PTMs can be divided into three categories: sequence-based, structure-based, and physics-based. Sequence-based approaches either flag individual residues prone to chemical degradation (e.g., methionine oxidation) or liable motifs (e.g., NG, NS, and NT for deamidation).¹⁵ Liable motifs for deamidation and isomerization were identified by investigating the effects of protein sequence on deamidation¹⁶ and isomerization¹⁷ rates for model peptides. Model peptides are suitable for assessing the effects of protein sequence on chemical degradation due to the generation of substantial chemical degradation in a short time.¹⁷

Structure-based approaches predict PTM liabilities by using structural features correlated with enzymatic and chemical modifications.¹⁸ Common structural features include, but are not limited to, secondary structure, water coordination number (WCN), and solvent-accessible surface area (SASA).⁹ WCN represents the average number of water molecules within the radius of an atom;¹⁹ SASA represents the surface area of the protein that interacts with the solvent.²⁰ Structural features such as SASA are typically extracted from crystal²¹ or predicted structures for antibodies.²² Alternatively, machine learning algorithms (e.g., NetSurfP) can also be used to predict structural features such as secondary structure, solvent exposure, and structural disorder for protein sequences.²³

Physics-based approaches are based on physical principles: for example, molecular dynamic (MD) simulations use Newtonian physics to predict the spatial position of atoms over time.²⁴ Physics-based approaches can predict the free energy barriers for chemical modifications^25,26 and probe protein dynamics.²⁷ MD simulations can be used to estimate averaged SASA, which captures changes in solvent exposure due to conformational changes. In addition, MD simulation can provide the root-mean-square fluctuations (RMSF) for Cα atoms, which captures structural flexibility.²⁸

Comparison of computational approaches

Sequence-based approaches are simple and easy to implement. Once the protein sequence is available, candidates with poor chemical stability can be eliminated by checking for liable motifs.¹⁵ However, using liable motifs alone can overestimate the number of liable residues or miss potential degradation hotspots.¹² The rate of chemical degradation for model peptides does not represent chemical degradation in native proteins, where the tertiary structure plays a role in chemical and enzymatic modifications.¹⁶ For example, high solvent exposure can increase the risk of chemical modifications due to the increased exposure to degrading agents such as water, peroxides, and light.²⁷

Structure-based approaches have better accuracy than sequence-based approaches because the in silico prediction includes the impact of secondary and tertiary structure on PTMs.²⁹ The accuracy of structure-based approaches is dependent on the quality of crystal or predicted structures and the selection of structural features. Knowledge of the chemical mechanisms for modifications (e.g., deamidation) can improve the accuracy of structure- and physics-based approaches. For example, some studies use the nucleophilic attack (Cγ – N_n+1) distance to predict the risk of deamidation due to the nucleophilic attack on the side chain by the backbone nitrogen during succinimide formation.^30,31

Although physics-based approaches provide valuable insight into the dynamics and energetics of chemical degradation, physics-based approaches are more computationally demanding compared to sequence- and structure-based approaches. It is challenging and time-consuming to predict free energy barriers and conduct MD simulations for larger proteins like antibodies.³² Improvements in the speed and accuracy of MD simulations and quantum chemical calculations are required for the broader adoption of physics-based approaches. Moreover, coarse-grained models can potentially improve the efficiency of MD simulations for larger proteins such as antibodies.³³

More recent studies combine these three computational approaches for in silico prediction of PTMs. Protein sequence and structural features have been incorporated into machine learning algorithms (e.g., random forest) to classify cold and hot spots for chemical modifications.³⁴ This review discusses current computational approaches for predicting isomerization, deamidation, oxidation, and lysine glycation in mAbs. In silico prediction of novel PTMs in mAbs, such as carbonylation, sulfation, and hydroxylation, is also considered.

Aspartic acid isomerization

Aspartic acid (Asp) residues in mAbs can undergo isomerization to form IsoAsp. During isomerization, there is a nucleophilic attack on the carbonyl group of Asp by the ionized amine group at the n + 1 residue. The Asp side chain requires a proton (H⁺) to form the succinimide intermediate and leaving water molecule, and thus isomerization is more favorable at low pH.³⁵ Finally, the succinimide intermediate is hydrolyzed at the n + 1 peptide bond to form isoAsp.¹⁷ At pH 4–7, Asp isomerization forms succinimide, which can cause a change in net charge. In contrast, IsoAsp formation can cause conformational changes because of the methyl group added to the peptide backbone. These conformational changes can cause changes in surface charge distribution or surface hydrophobicity.⁵ Depending on the location of IsoAsp formation, isomerization can generate basic,⁸ acidic,³⁶ hydrophilic, or hydrophobic variants.³⁷ Aspartate isomerization at the CDR loops can decrease the antigen-binding affinity.³⁸

Asp isomerization is sensitive to temperature and the dielectric constant of the solvent. At neutral pH, a high dielectric constant for solvents increases the pKa of Asp, which increases the concentration of the carboxylic acid (COOH) form of the Asp side chain. The COOH form is more reactive and prone to isomerization than the carboxylate form (COO⁻).³⁹ High temperatures can also accelerate the rate of isomerization reactions. Moreover, flanking residues, ionization state, and higher-order structure also influences isomerization.⁴⁰ For example, the risk of isomerization is highest for Asp residues within random coils due to structural flexibility and higher solvent exposure.⁴¹ The risk of isomerization also differs for different CDR loops. Asp isomerization at liable motifs is more likely with the CDR H3, H2, and L1 loops.⁴²

Sequence-based approaches

Sequence-based approaches for predicting Asp isomerization flag liable motifs for Asp isomerization such as DG, DS, DD, DT, and DH.⁴² Flanking glycine (Gly) residues make succinimide formation more favorable due to the lower steric hindrance. Serine and threonine residues can act as proton donors during isomerization.⁴⁰ Positively charged residues at n + 1 or n-1 position can accelerate Asp isomerization due to electrostatic effects.³⁵ Moreover, flanking residues that act as proton donors and acceptors can make Asp residues more reactive and prone to isomerization. For example, histidine residues can act as proton donors (H⁺) to OH⁻ ions at low pH during succinimide formation.⁴³

During the developability stage, protein engineering can be used to remove PTM liabilities and improve chemical stability.¹⁴ The risk of isomerization can be reduced by substituting Asp with glutamic acid, which is less prone to isomerization. Patel et al. identified a liable DS motif within the CDR loop that was prone to isomerization: substituting Asp with glutamic acid led to a loss in bioactivity. However, engineering the n + 1 residue (i.e., serine) retained the bioactivity while reducing the risk of Asp isomerization.⁴⁴

Structure-based approaches

The higher-order structure of antibodies also plays a crucial role in the risk for Asp isomerization. Dihedral angles can influence the reactivity of Asp residues and affect the risk of isomerization.²⁶ Diepold and colleagues observed isomerization for liable motifs within the solvent-exposed and flexible CDR loops. However, isomerization did not occur at the DG motif in the conserved region of the antibody.⁴⁵ Asp isomerization is the least favorable within beta-sheets due to structural rigidity. Sreedhara and colleagues observed succinimide formation for Asp74, which was solvent-exposed and part of a loop, whereas the Asp73 within beta-helix was less prone to isomerization.⁴⁰ In addition, hydrogen bonding for peptide backbone and side chains can make succinimide formation less favorable.⁴⁶

Structure-based approaches for predicting Asp isomerization incorporate solvent exposure, secondary structure, nucleophilic attack (Cγ – N_n+1) distance, and hydrogen bonding.⁴⁶ Sydow and colleagues used the root mean square deviation of Asp C^α atoms, size of the n + 1 residue, and secondary structure to develop a decision tree for classifying liable Asp residues.⁴⁶ They reported that Asp residues with high conformational flexibility and small size for n + 1 residue were more prone to isomerization. Although, a smaller side chain leads to lower steric hindrance for forming succinimide;²⁶ using the size of n + 1 residue for predicting isomerization neglects potential electrostatic effects from flanking residues on isomerization.³⁵

Physics-based approaches

MD simulations can supplement structure-based approaches for predicting Asp isomerization. Sharma et al. developed a logistic regression model to predict Asp isomerization at liable motifs (i.e., DG, DS, DT, DD, DH). ²⁹ The averaged SASA for Asp residues, RMSF for Cα atoms, and the averaged SASA for the n + 1 residue correlated well with Asp isomerization liability. Sharma and colleagues’ regression model successfully predicted five out of six liable Asp residues and all non-liable sites. Most of the parameters used in the study were derived from MD simulations. Thus, MD simulations can provide insight into how conformational flexibility plays a role in chemical degradation.

Recent studies have used MD simulations to predict chemical modifications; however, few studies attempt to predict the free energy barriers for chemical modifications. Plotnikov et al. used quantum mechanical/molecular mechanics (QM/MM) MD simulations to identify isomerization and deamidation hotspots. They estimated the free-energy barriers at different steps for isomerization and deamidation. The steric hindrance for succinimide formation is captured by the free energy barriers for adopting more reactive conformations to form succinimide.²⁶ In addition, Plotnikov et al. included long- and short-range electrostatic interactions in their QM/MM model. Compared to structure-based approaches, free energy predictions can incorporate electrostatic effects on Asp isomerization.²⁶

Therapeutic candidates with high and low risk for isomerization can be classified based on liable motifs and structural predictors such as SASA. In addition, knowledge of the antigen-binding sites, either from crystal structures or in silico predictions, can guide protein engineering to improve chemical stability. For example, if a liable Asp residue is involved in antigen binding, then the n + 1 residue can be substituted to improve chemical stability without affecting the bioactivity of the antibody.⁴⁴

Asparagine deamidation

Asn residues can undergo deamidation to form Asp or IsoAsp. There are three possible pathways for Asn deamidation. A nucleophilic attack on the backbone nitrogen at the n + 1 residue forms succinimide, which is hydrolyzed to form Asp. ⁴⁶ Alternatively, the nucleophilic attack by the backbone carbonyl group can form isoimide, which is hydrolyzed to IsoAsp and Asp (Figure 1). At acidic pH (pH ≤ 4), direct hydrolysis of Asn can also form Asp (Figure 1).⁴⁷ The direct hydrolysis pathway competes with succinimide-mediated deamidation.⁴⁸ Overall, deamidation is more favorable at pH ≥ 6. Forced degradation by thermal and high pH stress are used to confirm liable Asn residues.¹⁰ Flanking residues, secondary and tertiary structure, solvent exposure, and structural flexibility can affect deamidation.⁴⁹

Figure 1. Asparagine deamidation and aspartic acid isomerisation pathway. Reproduced from Sydow JF et al. (2014).⁴⁶ (Creative Commons Attribution License).

Sequence-based approaches

Liable motifs for deamidation include NG, NS, NN, NT, and NH:⁴² the NG motif has the highest deamidation rate, followed by NS, NT, and NH.¹⁶ The high deamidation rate for the NG motif can be attributed to lower steric hindrance for forming succinimide. Alternatively, Gly residues at the n + 1 position can provide more conformational flexibility, allowing Asn residues to adopt more reactive conformations to form succinimide.⁴¹ The n + 1 flanking residues can also play a role in deamidation. Chelius et al. reported the highest level of deamidation at SNG, ENN, LNG, and LNN.⁵⁰ Lu et al. investigated the deamidation for 131 clinical-stage antibodies and reported deamidation at liable motifs was more common at CDR-H2 and CDR-L1 loops.⁴² Deamidation at NT, NF, and NY only occurred within the CDR-H3 loop,¹⁴ which could be due to the length and flexibility of the CDR-H3 loop, allowing it to adopt more conformations.⁵¹ Excluding the CDR-H3 loop, the remaining CDR loops are clustered into canonical structures, where a few critical amino acids at specific positions can determine the CDR conformation.⁵² Non-H3 CDR loops are clustered based on structural similarity: dihedral angles and root mean square deviation have been used to compare and cluster structures for the same type of CDR loops (e.g., CDR-L1).⁵³ For example, the L1-11-A cluster refers to CDR-L1 loops with a length of 11 residues, which fall into the first cluster of CDR-L1 structures. ⁵³ Future studies could investigate how differences in conformations for different CDR-L1 and H2 clusters affect the risk of deamidation.

Structure-based approaches

The higher-order structure also influences the rate of Asn deamidation.⁴⁵ The risk of deamidation is low for Asn within beta-sheets due to structural rigidity and an extensive network of hydrogen bonds. Structure-based approaches for predicting deamidation incorporate secondary structure, hydrogen bonding, SASA, dihedral angles, and nucleophilic attack (Cγ – N) distance.³⁰ Random forest and decision trees are the most common machine learning algorithms for classifying hot and cold spots (Table 1). Delmar et al. developed a random forest model to identify liable Asn residue and predict the deamidation rate.³⁰ The backbone and side-chain dihedral angles are used to assess if the Asn side chain and backbone alignment are favorable for forming succinimide. However, the backbone dihedral angles are more critical than side-chain dihedral for predicting deamidation hot spots.^30,54

Table 1. Summary of parameters used for machine learning methods for predicting asparagine deamidation in mAbs. Abbreviations: RMSD, root mean square deviation; SASA, solvent-accessible surface area

Authors	Machine learning algorithm	Parameters	Data Set	Accuracy and precision^1	Reference
Sydow, J. F. et al. (2014).	Decision tree	RMSD of asparagine Cα-atoms, C-terminal amino acid size (>102.7 Å), SASA (>89.4 Å^2), and phi angles (>75.2º)	Data set consisted of 55 deamidation hotspots and 940 non-hotspots from forced degradation data for 37 mAbs.	Accuracy: 95.7 % Precision: 100.0 %	^Citation46
Yan, Q., Huang, M., Lewis, M. J. & Hu, P. (2018)^2	Decision tree	Secondary structure, SASA (> 10 Å^2), Cγ-Nn+1distance (> 3.4 Å), motifs (NS, NT, NG).	Data set consisted of deamidation sites for five IgG1 and five IgG4 antibodies.	Accuracy: 83.8% Precision: 41.7%	^Citation54
Delmar, J. A., Wang, J., Choi, S. W., Martins, J. A. & Mikhail, J. P. (2019).^2	Random forest	Residue at n + 1 position, half-life of pentapeptide, backbone dihedral angles, Cγ-Nn+1 distance, side-chain dihedral angles, and hydrogen bonding within side chains and secondary structures.	Combination of in-house data and deamidation data from Lu et al.^Citation39	Accuracy: 95.6 % Precision: 100.0 %	^Citation30

¹Accuracy is determined by the following formula: $\frac{\mathit{TP}+\mathit{TN}}{\mathit{TP}+\mathit{TN}+\mathit{FN}+\mathit{FP}}$ and precision was determined by the following formula: $\frac{\mathit{TP}}{\mathit{TP}+\mathit{FP}}$. TP = true positive, TN = true negative, FN = false negative, and FP = false positive.

²The accuracy and precision is reported from Delmar et al. (2019)³⁰

There are contradicting results for the impact of Cγ – N_n+1 distance on Asn deamidation. Yan et al. reported that a Cγ – N_n+1 distance > 3.4 Å was unfavorable for Asn deamidation, ⁵⁴ whereas Sydow et al. found no correlation between deamidation risk and Cγ – N_n+1 distance.⁴⁶ The discrepancy in results could be due to different stress conditions used to identify liable Asn residues for both studies. Yan et al. used high pH stress to identify liable residues, while Sydow et al. thermally stressed their samples at the formulation pH of 6.0. At basic pH, deamidation by forming succinimide is more favorable, and the Cγ – N_n+1 distance is a useful predictor for succinimide formation.⁴⁸ Therefore, the impact of structural features such as Cγ – N_n+1 distance on deamidation can differ between thermal stress at acidic and basic pH. Future studies could select different structural features to develop machine learning models to predict the risk of deamidation for different types of stress conditions.

Antibodies can also change conformations at different pH, which can cause a change in the solvent exposure for Asn residues. At high pH, heavy chain (HC)-Asn325 at the FG loop of the C_H2 domain is buried between the FG and BC loop (Figure 2). There is a conformational change in the BC loop at low pH, which causes HC-Asn325 to become exposed and susceptible to deamidation.⁵⁴ The CDR-H3 loop can adopt multiple conformations due to its length and flexibility. Lan et al. observed different conformations for the CDR-H3 loop at different pH. At neutral pH, the CDR-H3 loop is extended and solvent-exposed. In contrast, at low pH, the CDR-H3 loop is bent. ⁵⁵ Future studies could use MD simulations to identify more reactive conformations for deamidation at different pH.

Figure 2. Ribbon diagram representation for HC-Asn325 between BC and FG loop (PDB: 4BYH).

Oxidation

Oxidation of antibodies can occur due to peroxides within formulations, exposure to trace metals during manufacturing, and light exposure. Aromatic amino acids and the sulfur groups of methionine and cysteine residues are the most susceptible to oxidation.⁵⁶ Methionine (Met) oxidation generates sulfoxides and sulfone. Histidine (His) oxidation generates oxo-His. For tryptophan (Trp) oxidation, singlet oxygen reacts with the indole group to form kynurenine and N-formylkynurenine.^4,49 Met is the most prone to oxidation due to the reactive sulfur atom.⁵⁷

Sequence-based approaches flag all Met and Trp residues in the variable domains as liable residues.⁵⁸ However, sequence-based approaches can overestimate the risk of oxidation due to the impact of higher-order structure on oxidation. Exposed Met, Trp, and His are more prone to oxidation. SASA is a common predictor for identifying liable Met and Trp residues.⁹ More recent studies have successfully incorporated MD simulations and machine learning to predict Met and Trp oxidation.^29,59However, specific factors that affect oxidation risk can vary depending on the type of liable residue and the type of oxidative stress.⁶⁰ The subsequent sections outline in silico prediction for Met, Trp, and His oxidation in more depth.

Methionine oxidation

Met is oxidized by peroxides, light stress, and trace metal ions. Oxidation by peroxides involves transferring oxygen from the peroxide to Met through nucleophilic substitution.⁶¹ Water molecules stabilize the transition state near the sulfur atom.⁶² Alternatively, Met can also react with singlet oxygen to form Met sulfoxide. Met oxidation in mAbs can reduce conformational stability,⁶³ generate hydrophilic variants,⁵ cause structural changes,⁶⁴ and affect antigen binding.⁶⁵ Met residues within CDR loops are more prone to oxidation than Met residues within framework regions. For example, Yan et al. reported oxidation of Met-56 in the solvent-exposed CDR-H2 loops.⁷

SASA and WCN are useful predictors when water molecules play a role in stabilizing the transition state during oxidation by H₂O₂.⁶² WCN represents the average number of water molecules within the radius of an atom.¹⁹ The accuracy of WCN as a predictor depends on the selected radius. If the radius is too small, the risk of Met oxidation is underestimated; if the radius is too big, the risk of Met oxidation is overestimated. A radius of 5 Å- 6 Å for WCN is sufficient for predicting Met oxidation. ^19,66 Recent studies use the averaged SASA values of Met residues from MD simulations. MD simulations capture changes in solvent exposure for Met residues due to the flexibility of flanking residues and conformational changes.²³ In addition, the SASA of the sulfur atom can be a better predictor than overall SASA for Met residues.¹⁹

Machine learning models have been used to identify liable Met residues.^22,59 Yang and colleagues used the SASA of Met side chains derived from a random forest model to identify liable Met residues for 121 antibodies.⁶⁷ They observed a strong correlation between side-chain solvent exposure and experimental Met oxidation, but there were some false negatives. The study used SASA values from static crystal structures, which fails to capture changes in solvent exposure of Met residues due to conformational changes.⁶⁷ Sankar et al. developed a random forest model incorporating SASA, WCN, and distance between sulfur and aromatic residues, to predict Met oxidation by 2,2ʹ-azobis(2-amidinopropane) dihydrochloride (AAPH). Most of the incorrect predictions for Met oxidation were for the CDR-H3 loop. ²² Predicting the structure of CDR-H3 is challenging due to the longer length and high structural variation of the CDR-H3 loops.⁶⁸ Thus, the accuracy of SASA values is limited by the quality of the crystal or predicted structures.

Although MD simulations have provided more accurate SASA values, most studies conduct MD simulations with water. Excipients and buffers impact the risk of Met oxidation and the hydrodynamic radius of antibodies.⁶⁹ Changes in hydrodynamic radius of IgG can lead to changes in solvent exposure for Met residues. Conducting MD simulations with the formulation buffer can further bridge the gap between in silico predictions and experimental data. In addition, in silico tools can be used to identify formulations with a lower risk of Met oxidation.⁷⁰

In silico prediction of Met oxidation by peroxides is very well established. More recent studies have explored in silico prediction for Met oxidation under light stress and forced oxidation by AAPH. ^22,59 The mechanism for Met oxidation by peroxides and photolysis is different. Hydrogen peroxide oxidizes Met residues through nucleophilic substitution, whereas light exposure oxidizes Met by generating free radicals.⁴ In addition to using solvent exposure, Delmar and colleagues also used the closest atomic distance between phenylalanine and Met residues as a predictor for Met oxidation under light stress. Neighboring aromatic residues can prevent Met oxidation by acting as scavengers for free radicals.⁵⁹ Forced oxidation by H₂O₂ mimics exposure to trace peroxides from polysorbates in formulations, whereas photolysis is meant to mimic degradation due to light exposure during shipping and storage.⁷¹ Therefore, different predictors for Met oxidation (e.g., WCN) should be selected for different types of forced oxidation (e.g., peroxides and photolysis).¹⁹

Tryptophan oxidation

Trp oxidation can cause color changes due to kynurenine formation,⁷² reduce physical stability,⁷³ and oxidation of Trp residues in CDR loops can reduce binding affinity.⁶⁵ Structural flexibility, solvent exposure, secondary structure, and side-chain and backbone conformations affect Trp oxidation.⁵⁹ Most aromatic residues get buried within the hydrophobic core during protein folding;⁷⁴ however, the CDR loops are enriched with tyrosine and Trp residues, which play a role in antigen binding.⁷⁵

Trp oxidation is more common in the solvent-exposed and flexible CDR loops:⁹ Trp oxidation in CDR-H3,⁷⁶ CDR-L1,⁷⁷ CDR-L3,^78,79 CDR-H2, and CDR-L2 loops⁶⁰ has been reported. Neighboring aromatic residues and disulfide bonds can prevent photooxidation of Trp residues by acting as a light sink.⁸⁰ IgG1 antibodies have a conserved Cys-Cys-Trp triad, where the disulfide bond protects the neighboring Trp residue from photooxidation. Electron transfer from the Trp residue to the disulfide bond occurs as the disulfide bond dissipates the UV energy absorbed by Trp residues. ⁸⁰ Solvent exposure is the predominant predictor for Trp oxidation. SASA from static structures does not capture fluctuations in the exposure of Trp residues due to conformational changes. An averaged SASA sampled from MD simulations can provide a more accurate prediction for Trp oxidation. Solvent exposure for Trp side chains is a good predictor because the oxidation of the indole ring is a critical step for forming kynurenine.²⁷ In addition, SASA for Trp side chains can be an accurate predictor for oxidation in cases where the Trp backbone is buried, but the side chain is exposed. Sharma et al. observed that a Trp side chain SASA of > 80 Ų correlated well with Trp oxidation after incubation with AAPH.²⁹

Other structural factors also influence the oxidation rate and liability for different Trp residues (e.g., side chain and backbone dihedral angles) and the electrostatic environment.⁶⁰ Recently, Delmar and colleagues developed a random forest model to predict Trp oxidation under light stress, which incorporated secondary structure and psi and phi angles.⁵⁹ Future studies could investigate the role of other structural factors on Trp oxidation for a broader panel of antibodies.

The risk of Trp oxidation can also differ under long-term thermal stress compared to forced oxidation by AAPH and photolysis. For example, Jacobitz et al. reported that solvent exposure of Trp residues correlated better with Trp oxidation under long-term thermal stress than oxidation by AAPH and light stress.⁶⁰ Future studies can continue to tailor in silico prediction of Trp oxidation to specific stress conditions, such as thermal stress, photolysis, and forced oxidation by AAPH and peroxides combined with trace metals.

Histidine oxidation

Photo-oxidation and metal-catalyzed oxidation (MCO) are the primary mechanisms for His oxidation. During light stress or MCO, a singlet oxygen reacts with the imidazole ring of His to form peroxide intermediates that then yield oxo-His or cross-linked products.⁴ Recently, His oxidation has been reported for mAbs: Luo et al. observed oxidation at His-304 and His-428 for IgG2 antibodies after incubation with Cu²⁺/ascorbate,⁸¹ and Amano et al. reported His oxidation in the C_H2 domain of an IgG1 antibody under light stress.⁸² Studies have also reported cross-linking of mAbs caused by His oxidation, with cross-linking being more favorable at higher pH.⁸³ For IgG1 antibodies, oxidized His residues in the Fc and hinge region can cross-link with His, cysteine, and lysine residues after light stress.⁸⁴ Interestingly, His-His cross-linking has also occurred between the His side chain of an IgG4 antibody and the free His in the formulation.⁸⁵ Therefore, in silico tools can be used to select formulations for antibodies with a higher risk of His oxidation.

Solvent exposure,⁸⁵ conformational flexibility,⁸⁴ surrounding residues, and pKa influence His oxidation. Not all solvent-exposed His residues will get oxidized,⁸² but solvent exposure is still a prerequisite for His oxidation. His residues within highly flexible and solvent-exposed regions are more prone to oxidation. For example, the oxidized HC-His440 and HC-His292 were located at the C_H2-C_H3 interface within loops,⁸⁴ whereas the oxidized His231⁸⁴ and His220⁸³ were within the hinge region. Surrounding residues can also enhance or inhibit His oxidation. Neighboring Trp residues can act as photosensitizers and increase the risk of His oxidation. Under light stress, Trp residues can generate singlet oxygen or superoxide anions, which can oxidize neighboring His residues.^86,87

Studies on the impact of pKa on the risk of His oxidation have yielded inconsistent results. For example, Amano et al. reported no correlation between pKa of His residues and His oxidation,⁸² whereas Miyahara et al. reported His residues with lower pKa values were more prone to photooxidation.⁸⁸ The discrepancy in the reported role of pKa on His oxidation could be due to the different stress conditions (i.e., MCO versus UVC light exposure) used in these studies and the different by-products formed after the stress condition (i.e., oxo-His versus Asn and Asp residues). Compared to Met and Trp oxidation, in silico prediction of His oxidation is understudied.⁴⁹ Future studies could investigate the impact of pKa on His oxidation for different types of oxidative stress (e.g., photooxidation versus MCO).

Lysine glycation

Glycation is a non-enzymatic modification, where amino groups at lysine (Lys) and arginine residues or the N-terminal are glycated by reducing sugars such as glucose.⁸⁹ The Schiff base is formed after a condensation reaction between the aldehydes of the reducing sugars and the amine groups of Lys residues (Figure 3(a)). Schiff base formation is reversible; however, the multistep Amadori rearrangement can generate more stable ketoamines (Figure 3(b)).⁹⁰ Glycated Lys residues can further degrade to form advanced glycation end products,⁹¹ which can cause an immunogenic response. Lys glycation typically occurs during long-term storage in formulations, cell culture, and in vivo.⁹² Disaccharides in formulations can degrade to form reducing sugars, leading to glycation during storage.⁹² During cell culture, a high glucose concentration can accelerate Lys glycation and reduce protein yield.⁸⁹ The effects of Lys glycation on antigen binding varies from decreased binding affinity⁹³ to minimal change in binding affinity.⁹⁴ Lys glycation can increase protein aggregation by affecting the net surface charge, reducing the electrostatic repulsion between antibodies.⁹⁵ Forced glycation is used to identify liable Lys residues by incubating the IgG with high concentrations of reducing sugars (e.g., glucose) at high temperatures.⁹² Forced glycation in citrate buffers correlates well with glycation during cell culture.⁹⁶

Figure 3. (a) Mechanism for Schiff base formation. (b) Mechanism for Amadori rearrangement on glycated lysine.

The use of SASA to identify liable residues is a popular approach.²⁷ High solvent exposure is a prerequisite for Lys glycation, but not all exposed Lys residues are glycated.⁹⁶ Lys glycation is also influenced by flanking residues (i.e., polar, ⁹⁷ acidic, and basic residues⁹²), three-dimensional structure, and pKa of Lys residues.⁹⁰ Neighboring Asp and His residues can slow Schiff base formation; however, flanking His and Asp residues can promote Amadori rearrangement by acting as proton donors and acceptors (Figure 3(b)).^90,97 Liable motifs for glycation include KD, KXD, KXK, and KXE.⁹⁶ Structure-based approaches can identify acidic residues that are not part of motifs, but are close enough to catalyze glycation.⁹⁸ For example, Zhang et al. reported glycation at light chain (LC)-Lys49, which was 11 Å apart from LC-Asp31.⁹⁷ They suggested that the flexibility of Lys and Asp side-chains promoted Amadori rearrangement by reducing the distance between LC-Lys49 and LC-Asp31.⁹⁷

Basic residues such as arginine can promote Schiff base formation by lowering the pKa of the Lys side chain.^90,95 Lys side chains with lower pKa values are more likely to adopt the deprotonated form (NH₂), which is more favorable for forming the Schiff base. The deprotonated form of the Lys side chain can perform a nucleophilic attack on the carbonyl carbon of the sugar to form the Schiff base (Figure 3(a)). Miller et al. suggested that neighboring arginine and His residues increased reactivity of HC-Lys98 by influencing the pKa of HC-Lys98.⁹⁵ However, Zhang et al. found no correlation between Lys side chain pKa and glycation risk.⁹⁷ They inferred that Amadori rearrangement (Figure 3(b)) was the rate-limiting step for the glycation of a recombinant humanized mAb. Thus, the impact of the Lys side chain pKa on glycation may depend on the rate-limiting step for Lys glycation. QM/MM simulations of free energy calculations could be used to help identify the rate-limiting step for the glycation of different antibodies. Moreover, future studies could investigate if different backbone and side-chain conformations are more favorable for glycation.

Carbonylation

Metal ions can catalyze oxidative carbonylation of arginine, Lys, proline (Pro), and threonine residues. Transition metals such as iron and copper can convert oxygen (O₂) to superoxide radical anions (O₂^–·).⁹⁹ During carbonylation, free radicals attack the side chain and add an amine or ketone group. Arginine and Pro are converted to glutamic semialdehyde; Lys is converted to aminoadipic semialdehyde, and threonine is converted to 2-amino-3-ketobutyric acid (Figure 4).¹⁰⁰ Exposure to trace metals from stainless steel surfaces and glass vials can cause carbonylation of mAbs during manufacturing and storage.¹⁰¹ Carbonylation of arginine and Lys residues leads to a loss in positive charge, which generates acidic variants. For example, Yang et al. reported increased acidic variants after forced oxidation with ferrous sulfate and hydrogen peroxide.¹⁰² Oxidative carbonylation of mAbs can also increase protein aggregation.¹⁰³

Figure 4. Carbonylation pathway for (a) Lysine (b) Proline (c) Arginine and (d) Threonine.

Lys is most susceptible to carbonylation, followed by Pro as the next most susceptible residue, and arginine and threonine are the least susceptible residues. Pro residues tend to form kinks that make Pro residues more exposed and increases the risk of carbonylation.¹⁰⁴ Carbonylation sites are more likely to occur in regions enriched with arginine, Lys, Pro, and threonine residues. Regions at the N- and C-terminal of the carbonylation sites tend to be enriched with positively charged residues.^105,106 Solvent exposure and structural flexibility also increase the risk of carbonylation.¹⁰⁵ Currently, there is no clear association between carbonylation sites and different types of secondary structures.^104,107 Multiple carbonylation predictors such as CarSPred,¹⁰⁶ iCar-PseCp,¹⁰⁸ and iCarPS,¹⁰⁹ have been developed using random forest and support vector machine (SVM) algorithms (Table 2). The most common features of machine algorithms for predicting carbonylation include the position-specific propensity of amino acid, hydrophilicity, SASA, side-chain interaction parameter, and positive charge (Table 2). ^106,108,109

Table 2. Comparison of online carbonylation predictors used for machine learning methods. Abbreviations: KNN, k-nearest neighbor; ROC, Receiver Operating Characteristic

Carbonylation predictor	Machine learning algorithm	Parameters	Data Set	Area under the ROC^1	Reference
CarSPred	Weighted support vector machine	Position-specific propensity of amino acid, k-spaced amino acid pair, KNN scores, physicochemical properties (electric properties, hydrophobicity, alpha and turn propensities, etc.)	331 lysine, 131 arginine, 128 threonine, and 129 proline carbonylation sites were extracted from 230 carbonylated human proteins. In addition, 22 lysine, 3 arginine, 6 threonine, and 15 proline carbonylation sites were extracted from carbonylated mouse, rabbit and bovine proteins.	Lysine: 0.6704 Arginine: 0.5345 Threonine: 0.6800 Proline: 0.7873	^Citation106
iCar-PseCp	Random forest	pseudo amino acid composition	Data was derived from 230 human carbonylated protein sequences and 20 carbonylated proteins from Photobacterium and Escherichia coli.	Lysine: 0.8728 Arginine: 0.8668 Threonine: 0.8603 Proline: 0.8484	^Citation108
iCarPS	Random forest	3-D conical coordinates and physicochemical properties (hydrophobicity, hydrophilicity, mass, pK1, pK2, pI, rigidity, flexibility, and irreplaceability)	Same benchmark dataset as Lv, et al. (2014).^Citation109	Lysine: 0.789 Arginine: 0.726 Threonine: 0.790 Proline: 0.814	^Citation109

¹Area under the curve was derived for the ROC. The ROC plots the sensitivity (i.e., true positive rate) versus selectivity.

More precise, sensitive, and site-specific analytical methods for detecting carbonylation in antibodies have been developed in the past decade.^5,99,100,110 Joshi et al. identified 27 carbonylation sites within trastuzumab. Most carbonylation sites were within the framework region; they reported no carbonylation sites within the C_H3 domain and CDR loops.¹⁰¹ The abundance of carbonylation sites within the structurally rigid framework regions^100,101 is surprising because previous studies have suggested that carbonylation is more common in structurally flexible or disordered regions. ^104,107 However, residues within the framework region can have a high solvent exposure and undergo carbonylation. ¹⁰⁰ Future studies could investigate the impact of the structure of antibodies on carbonylation. In addition, carbonylation sites from a wide panel of antibodies could be used to develop machine learning models for specifically predicting carbonylation sites for therapeutic antibodies.

Tyrosine sulfation

Tyrosine sulfation is an enzymatic modification that is catalyzed by tyrosylprotein sulfotransferases (TPSTs). During sulfation, a sulfate group is attached to the hydroxyl group of tyrosine residues.¹¹¹ Tyrosine sulfation has been reported for a few monoclonal and bispecific antibodies.^111–114 Sulfation of mAbs can occur in Chinese hamster ovary (CHO) cells during the cell culture. The degree of tyrosine sulfation in CHO cells varies due to differential expression of phosphoadenosine-5ʹ-phosphosulfate (PAPS) synthetase and TPST. PAPS synthetase converts ATP to PAPS; TPST transfers the sulfo group from PAPs to tyrosine residues.¹¹⁵ Sulfation in mAbs generates acidic variants,¹¹¹ but the impact of tyrosine sulfation on the safety and efficacy of therapeutic antibodies has not yet been established.⁵

Common predictors of tyrosine sulfation include secondary structure, flanking residues, structural flexibility,^101, and SASA.¹¹⁶ Sulfation is more favorable for structurally flexible or disordered regions that can fit into the cleft of TPSTs. Tyrosine residues flanked by acidic residues are also more likely to get sulfated (Figure 5(a)).¹¹¹ During sulfation, the flanking acidic residues within the IgG form electrotactic interactions with arginine and Lys residues of TPST.¹¹⁷ The role of structural flexibility in sulfation suggests that sulfation is more likely to occur in the solvent-exposed and flexible CDR loops than in the framework regions.⁵¹ Zhao et al. reported tyrosine sulfation in the CDR-L1 loop,¹¹¹ and Tyshchuk et al. reported sulfation in the CDR-L2 loop.¹¹²

Figure 5. Sequence logos showing the most conserved amino acids around (a) sulfotyrosine (sTyr). (b) Sequence logos showing the most conserved amino acids around 4-hydroxyproline (4Hyp) residues. Reproduced from Tyshchuk, O. et al.¹¹² (Creative Commons Attribution License).

There are multiple online tools for predicting sulfation, such as Sulfinator,¹¹⁸ PredSulSite,¹¹⁹ and Sulfotyrosine.¹²⁰ However, the accuracy of online tools is limited by the training set used to develop machine learning algorithms to predict PTMs.³⁴ The training sets for most online predictors include experimentally identified tyrosine sulfation sites from databases like Swiss Prot. However, mAbs are underrepresented, and an overtrained machine learning algorithm may fail to identify sulfation sites in mAbs. Hence, sulfation predictors have to be used with caution for assessing the risk of tyrosine sulfation in mAbs. Tyrosine sulfation needs to be identified in a broader panel of antibodies to improve tyrosine sulfation prediction for mAbs. Future studies could also explore the impact of structural features such as dihedral angles, pKa of surrounding residues, and position within IgG (i.e., CDR or framework region) on tyrosine sulfation.

Hydroxylation

Hydroxylation is an enzymatic modification that is catalyzed by hydroxylases. Hydroxylation can occur at arginine, tyrosine, Trp, and phenylalanine, but it is more common for Pro and Lys.¹²¹ During hydroxylation, a hydroxyl (OH) group is added to Pro or Lys residues. As a result, Pro is converted to 3-hydroxyPro or 4-hydroxyPro, whereas Lys is converted to 5-hydroxyLys.¹²² Hydroxylation of Lys and Pro residues at consensus motifs (Xaa-Lys-Gly or Xaa-Pro-Gly) is common in collagen, where hydroxylation helps stabilize the collagen triple-helix.¹²³ Common parameters for designing hydroxylation predictors include solvent exposure, intrinsic disorder, hydrophilicity, and sequence.^122,123 Hydroxylation sites tend to be disordered, exposed, and enriched with Pro and Gly residues.¹²⁴

There are multiple online predictors for hydroxylation sites such as PredHydroxy, HydPred, and ModPred.¹²³ Online hydroxylation predictors are suitable for collagen. However, most hydroxylation predictors have limited training sets, which leads to a high number of false positives and false negatives for non-collagen proteins and uncommon hydroxylation motifs. To date, two studies have reported hydroxylation in antibodies. Xie et al. reported Lys hydroxylation for an IgG1 antibody,¹²⁵ and Tyshchu et al. observed Pro hydroxylation for a bispecific antibody.¹¹² Using consensus motifs can lead to the overestimation of hydroxylation sites. For example, Xie et al. did not report Lys hydroxylation at all of the consensus motifs.¹²⁵ Most hydroxylation predictors are not specific to mammalian and human proteins, limiting their utility for therapeutic antibodies.¹²³ Tyshchuk and colleagues used hydroxylation sites for mammalian proteins to generate sequence logos (Figure 5(b)) and develop a kNN model:¹¹² the in silico prediction correlated well with the experimental data for the bispecific antibody. Huang et al. used a similar approach to create HydLoc, which had a training set composed of human proteins and had better accuracy for predicting hydroxylation of human proteins.¹²⁴ Machine learning methods for predicting hydroxylation can be improved with more experimental data for uncommon hydroxylation sites, as well as hydroxylation sites for mammalian proteins, human proteins, and antibodies.¹²³

Perspectives and future directions

In the past decade, in silico prediction of PTMs has evolved from simple sequence-based approaches to more sophisticated approaches incorporating antibody modeling, MD, and machine learning.¹² Although the structure- and physics-based approaches have better accuracy than sequence-based approaches, there is still a gap between the in silico prediction and experimental data from forced degradation studies. External factors such as formulation and cell culture conditions can contribute to the gap between in silico predictions and forced degradation results.³ The quality of predicted structures is a limiting factor for structure- and physics-based approaches. Homology modeling is suitable for modeling the conserved framework regions. Recent advances in loop modeling have provided more accurate predictions for the L1, L2, L3, H1, and H2 loops.¹²⁶ However, modeling the CDR-H3 loop is still challenging due to the length and structural variation.¹²⁷ Future advancements in antibody modeling will improve the accuracy in silico prediction of PTMs. In addition, MD simulations can provide more representative SASA values for predicted and crystal structures and refine the quality of predicted structures.

Machine learning algorithms such as random forest and decision trees have been used to predict the risk of deamidation, isomerization, and oxidation. However, machine learning models are limited by the training set and selection of parameters. Building a training set can require a lot of in-house experimental data. Experimentally confirmed liable sites for publicly available antibody sequences are essential for developing better machine learning models in the future. Selecting the best parameters is complicated by contradictory findings for the impact of structural features (e.g., pKa values) on PTMs. Studies investigating the effect of different structural features on PTMs for a broad panel of antibodies are crucial for improving in silico prediction of PTMs. ^42,67 Moreover, the selection parameters for machine learning models should account for different mechanisms for chemical modifications (e.g., oxidation) under different types of stress conditions, such as forced oxidation by peroxides⁶⁶ and photolysis.

Therapeutic antibodies can also undergo in vivo chemical modifications (e.g., deamidation and oxidation) due to rapid temperature and pH changes after drug administration.¹²⁸ In vivo chemical modifications can impact the safety and efficacy of therapeutic antibodies. For example, in vivo deamidation can lead to a loss in activity and cause immunogenicity. ⁵⁷ Few studies have investigated in vivo modifications by monitoring PTMs for IgG in serum or by using in vitro conditions to mimic in vivo degradation.⁵⁷ Monitoring in vivo modifications is difficult because the human serum is a mixture of endogenous proteins and the administered therapeutic antibodies. It is challenging to separate and distinguish therapeutic antibodies from endogenous proteins. Fluorescence labeling techniques or affinity purification with anti-Fc antibodies are required for mass spectrometry (MS) analysis of in vivo samples. The harsh purification conditions (e.g., low pH elution) can introduce artifacts in liquid chromatography-MS analysis.¹²⁹ Computational approaches for improving the stability of therapeutic antibodies focus on preventing chemical degradation under the types of stresses (e.g., thermal stress) encountered during manufacturing and storage.²⁷ In silico prediction of in vivo modifications remains relatively unexplored. Developing computational approaches to predict in vivo chemical modifications provides opportunities for improving in vivo stability of therapeutic candidates.

The impact of PTMs on heterogeneity, safety, and efficacy can vary depending on the location of the liable residue. A more holistic approach to engineering mAbs could identify PTM hotspots and predict the impact of PTMs on binding affinity and physical stability. Removing liable sites is not always feasible, especially if the amino acid residue plays a role in antigen binding. Paratope prediction and molecular docking of IgG-antigen complex can guide antibody engineering for developability.¹³⁰ The surface charge distribution and surface hydrophobicity¹³¹ could be used to predict the impact of modifications of liable residues on colloidal stability.

Substantial strides have been made in silico prediction of deamidation, isomerization, and Met oxidation.^22,30,46 However, the in silico prediction of glycation, sulfation, hydroxylation, and carbonylation mAbs is underdeveloped compared to other PTMs. Another understudied area is in silico prediction of in vivo chemical modifications.⁵⁷ Advances in analytical characterization of novel modifications (e.g., sulfation) and in vivo modifications are crucial for developing computational approaches for predicting PTMs.

After two decades of establishing a solid foundation for in silico prediction of chemical modifications, future studies can focus on tailoring in silico predictions to more specific stress conditions and formulations.

Abbreviations

Table

AAPH	2,2ʹ-Azobis(2-amidinopropane) dihydrochloride
Asn	asparagine
Asp	aspartic acid
CDR	complementary-determining region
CHO	Chinese hamster ovary
Gly	glycine
HC	heavy chain
His	histidine
kNN	K-nearest neighbor
LC	light chain
Lys	Lysine
mAbs	monoclonal antibodies
MCO	metal-catalyzed oxidation
Met	methionine
MD	molecular dynamics
PAPS	phosphoadenosine-5ʹ-phosphosulfate
Pro	proline
PTM	post-translational modification
QM/MM	quantum mechanics/molecular mechanics
RMSF	root-mean-square fluctuations
SASA	solvent-accessible surface area
SVM	support vector machine
TPST	tyrosylprotein sulfotransferases
Trp	tryptophan
WCN	water coordination number.

Disclosure statement

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

Additional information

Funding

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