Introduction

Interleukin-10 (IL-10) is a critical anti-inflammatory cytokine for treating inflammatory bowel disease (IBD). However, in the inflamed gut, extreme pH and higher local temperatures reduce IL-10’s activity. The structure of IL-10 has been highly studied and many studies report that the biologically active form of IL-10 can be thermally unstable in its homodimer form, as it has a short half-life and can be easily degraded in vivo.

To address this, our Dry Lab aimed to engineer IL-10 variants with improved thermostability. We designed two optimized IL-10 variants (M1 and M2) using a computational pipeline that integrates ThermoMPNN, PyMOL, AlphaFold3, and ESM modeling. This approach allowed us to target specific point mutations predicted to stabilize IL-10’s fold without disrupting its dimeric structure or receptor binding.

ThermoMPNN

ThermoMPNN is a deep neural network trained on vast datasets of experimental stability measurements, allowing it to learn the complex relationships between sequence, structure, and stability. Based on the well documented ProteinMPNN framework, ThermoMPNN analyzes a protein's 3D structure to predict how changes in its amino acid sequence (point mutations) will impact a protein’s entire thermal stability as denoted through the point mutation’s resulting change in -ΔΔG.

We began with the high-resolution crystal structure of IL-10 (PDB 1ILK) and ran ThermoMPNN via Google Colab.

Figure 1. Workflow of the ThermoMPNN model used to identify stabilizing point mutations. The pipeline combines structural features of IL-10 (e.g., backbone distances) with sequence embedding processed by an encoder–decoder network. The output is passed through a stability-prediction module (LA and MLP blocks) to calculate the predicted change in folding free energy (ΔΔG) for each possible amino-acid substitution. Negative ΔΔG values (blue on the heat-map) indicate mutations predicted to increase protein stability, whereas positive ΔΔG values (red) suggest destabilizing effects. This method builds on the design principles described in Mega-scale experimental analysis of protein folding stability in biology and design, https://www.nature.com/articles/s41586-023-06328-6

Further detail is provided below:

Starting with a Research Collaboratory for Structural Bioinformatics protein data bank (RCSB PDB) file of our desired protein IL-10. We chose the crystal structure PDB 1ilk-IL10, of IL-10, derived using X-Ray diffraction with a high PDB validation report. To process the ThermoMPNN software, we used Google Collab’s implementation.

The app is organized into easy to follow steps used for running the software. After confirming each step, the user hits the play key on the left side of the platform.

Following, the user uploads a desired PDB file and chooses which protein chains for the software to focus on.

Standard ThermoMPNN mode was used for our dry lab. Importantly, running the “Include Cysteines” option determines whether the model can introduce cysteine residues into designed sequences. Allowing cysteines takes into account potential disulfide bonds, which can greatly improve protein thermostability but also risk misfolding or unwanted aggregation strain on the structure. For IL-10 engineering, our group consulted with Maggie Horst postdoc and chose to toggle this setting off as a cysteine-free sequence is safer for expression and folding.

Next, the multichain interference was experimented with in our dry lab outputs. The multichain interference allows the model to account for the interface between the two chains. For IL-10 this prevents designs that stabilize one monomer but destabilize the dimer as a whole. Our results when toggling this option on versus off did not differ much or introduce new mutations that were significant in improving thermostability. This is logical as the IL-10 chains mirror one another.

Lastly, post-processing cleans up the raw sequence output by removing duplicates and low-confidence mutations. It also filters against unwanted motifs, like buried polar residues or problem cysteines, that could cause instability or misfolding. For example sample output could list “Identified the following disulfide engaged residues: [].” Finally, it ranks the remaining candidates by predicted thermostability so only the most promising variants move forward.

Raw Data

ThermoMPNN produced a heat-map of ~7,000 candidate substitutions, showing stabilizing (blue) and destabilizing (red) mutations.

Figure 2. Thermostability heat map of IL-10 single-site mutations generated using ThermoMPNN.

The heat map illustrates predicted changes in folding free energy (ΔΔG, kcal/mol) across all amino-acid substitutions and positions of IL-10. Negative ΔΔG values (blue) correspond to mutations predicted to enhance thermostability, while positive ΔΔG values (red) suggest destabilizing effects. Although most substitutions showed little improvement, we observed that mutations clustered in positions 32–38 and 133–139 displayed a distinct blue gradient.

Although these effects are subtle, we believe they may contribute to incremental increases in protein stability. Notably, because IL-10 is a homodimer, these regions exhibit symmetry in their interactions, enabling us to identify pairs of mutations across the dimer interface that could act synergistically to enhance thermostability.

Exported Data (to Google Sheets)

We exported the raw data to Google Sheets for further curation. Only ~1,018 mutations had negative ΔΔG values. Among the viable hits, R97V and R97I ranked highest for stabilizing effect and became focus points for deeper analysis.

Figure 3. ThermoMPNN’s data results were organized into Google Sheets. The ddG (kcal/mol) value represents the predicted change in protein stability when a mutation is introduced compared to the wild type. A negative ddG suggests the mutation is stabilizing the protein, while a positive value indicates it may destabilize the structure.

In total, close to 7 thousand point mutations were made, however, only the first 1018 mutations had a stabilizing ddG. These 1018 became the focus of our selection. Here, many mutations that could yield high thermostability results were ruled out by marking “X.” The “X” denotes that ThermoMPNN could not map that mutation to a resolved residue position in the structure used. This usually happens when the residue exists in the input sequence but is missing from the PDB file, or when the alignment between sequence and structure does not line up cleanly due to gaps or numbering mismatches. Another possibility is that the residue was filtered out during preprocessing because no backbone coordinates were available to calculate stability. After establishing this, we focused on the available mutations that we would be able to clearly see the PDB and thermoMPNN mutation position it corresponded to. Mutations at the top of the list such as R97V and R97I immediately caught our attention. These had the most significant stabilizing impact on IL-10’s thermostability and needed more confirmations to decide on its viability.

Our dry lab addresses the problem of wild-type IL-10 being naturally unstable with a ΔG of positive folding value of ~4.36 kcal/mol, predicted by ThermoMPNN. This explains its limited half-life and limited therapeutic potential. Screening of nearly 7,000 mutations, our group utilizedThermoMPNN to predict over 1,000 candidates that have thermally stabilizing influences on the structure of IL-10, with clusters between residues 32–38 and 133–139 showing uniform improvements.

Adding to further confirmations, we modeled all suggested mutations with high thermostability using PyMOL for visualization and biochemical interaction analysis. This allowed us to rule out mutations that imposed steric strain, destroyed salt bridge interactions, or decreased dimer geometry. Moreover, this analysis revealed variants that combined positive rotamers with strong ΔΔG scores. Various combinations of mutations with these affirmed confluences were tested for their ability to function as whole through ESM-generative modeling. These results were successful in showing lower ΔG scores from the original sequence. AlphaFold3 modeling further showed that high-priority mutations preserved IL-10's key receptor interactions with IL-10RA and IL-10RB, evidenced by stable ipTM and pTM values across mutant complexes. Together, these findings suggest that targeted point mutations can greatly enhance IL-10 thermostability without compromising structural integrity or receptor affinity, allowing the feasibility of more stable therapeutic variants.

ESM Model

Baseline Stability of IL-10

When we analyzed the wild-type (unmutated) IL-10 monomer, the ESM predictor estimated a folding free energy of approximately ΔG ≈ +4.36 kcal/mol. A positive ΔG implies that the isolated monomer, on its own, would not be strongly stable in the folded state. However, IL-10 is naturally a homodimer, meaning two identical subunits pair together. This dimerization greatly improves the effective stability of the protein in vivo, so IL-10 as a whole is functional despite its intrinsically less stable monomer.

Figure 4. Baseline ESM-predicted folding free energy (ΔG) for wild-type IL-10. The ESM model estimated the folding free energy of unmutated IL-10 to be ΔG ≈ +4.36 kcal/mol, indicating that the isolated monomer is relatively unstable and tends to unfold at equilibrium. This baseline value served as the reference for evaluating whether our engineered mutations could lower ΔG and thus improve the protein’s predicted thermostability.

Application of ESM

Using an ESM-based predictor, unmutated IL-10 shows ΔG fold ​ ≈ +4.36 kcal/mol, which indicates that the isolated monomer would be predominantly unfolded at equilibrium, signifying weakness in thermostability. IL-10 functions as a homodimer in vivo, and the dimerization free energy contributes substantially to stabilization such that the dimer can be thermodynamically favored despite the weak intrinsic monomer stability. ESM estimates are derived from statistical patterns in large protein sequence datasets rather than physics-based simulations, so we report values together with the sign convention and interpret them as overall thermostability predictions.

To validate whether our collection of point mutations addresses these problems in functional stability, we ran our point mutations through the software to find the combinations of them that would decrease the predicted ΔG the most.

These single-residue changes may push this ΔG downward, making the protein more likely to hold its folded shape under heat or stress. In practice, decreasing the ΔG means the protein should resist unfolding better, which is the goal when trying to engineer higher thermostability.

Using the ThermoMPNN predicted results, we used visualization softwares to filter for practical mutations that would still be impactful for thermostability improvements (see protein folding). Then, running through the ESM-predictive model, our best combination mutations 1 and mutation 2 of IL-10 yielded results of 4.2645 and 4.25949 ΔG, a significant decrease and improvement from the original IL-10 value of 4.357622. Because ESM stability scores are statistical proxies, not calorimetric measurements it is difficult to estimate the improvement in degrees of these predictions. Hence, wet lab confirmation of this order-of-magnitude prediction is required.

Figure 5. ESM-predicted folding free energy (ΔG) for selected IL-10 mutation set. The mutation combination N148M, R32K, E45M, N126M, T26Q, and Q33D was evaluated using the ESM predictor. Chain A: ΔG ≈ 4.309 kcal/mol, showing a notable reduction compared to the wild-type IL-10 baseline of 4.36 kcal/mol; Chain B: ΔG ≈ 4.378 kcal/mol, which is comparable to the wild type.

These results suggest that the engineered substitutions can incrementally improve the stability of Chain A, which may be beneficial in the context of IL-10’s dimeric structure, although experimental validation is required to confirm these computational predictions.

Figure 6. ESM-predicted folding free energy (ΔG) for alternative IL-10 mutation set. The mutation combination N36C, E45M, S2M, N139M, Q33D, and D19E was evaluated using the ESM predictor. The predicted folding free energy was ΔG ≈ 4.2595 kcal/mol, which represents a greater reduction compared to the wild-type IL-10 baseline of 4.36 kcal/mol than the first mutation set achieved. This result suggests that this alternative combination may offer a more pronounced predicted gain in thermostability, reinforcing its potential as a candidate for experimental validation in improving IL-10’s stability.

Summary

Given above are the mutation lists and their respective ESM-generative model results. These mutations are combinations of those that passed through pyMOL visualization and biochemical analysis of R-group interactions. As shown above, both sets of mutations show a reduction in the deltaG predicted value from 4.36 in the original amino acid sequence down to 4.30895 and 4.2594 on chain A respectively from both groups. These reductions, although numerically small, are meaningful because even fractional shifts in ΔG can translate to noticeable changes in folding stability under physiological conditions. In practice, such improvements often correlate with higher melting temperatures (Tm), allowing the protein to remain functional and resist unfolding at elevated temperatures. This kind of stability is especially valuable for therapeutic proteins like IL-10. While experimental validation is necessary, the computational trend points toward a practical gain in thermostability that could enhance the reliability of IL-10 in clinical applications.