Overview

Manganese superoxide dismutase (Mn-SOD) is a key antioxidant enzyme but suffers from poor stability, rapid activity loss, and short half-life, limiting its industrial and biomedical applications. Traditional enzyme engineering approaches often fail to simultaneously improve both thermostability and catalytic activity, and experimental screening can be costly and time-consuming. Modeling has played a crucial role in our project. We used transformer-based deep learning (DL) to design Mn-SOD variants with improved thermostability and catalytic performance. DL models predicted beneficial mutations that enhance catalytic function. This strategy enabled efficient, high-confidence identification of promising Mn-SOD variants, providing a critical foundation for determining our experimental planning. Additionally, we also constructed a multi-scale mathematical model from microcatalytic kinetics to macroscopic 0-dimensional homogeneous system and then to 1-dimensional spatial reaction-diffusion, aiming to simulate and optimize the process of Mn-SOD eye drops from molecular mechanism to in vivo spatiotemporal efficacy, thereby significantly improving the efficiency and accuracy of formula screening and drug delivery protocol design.

2 Machine Learning（ML）and Deep Learning (DL)

Machine learning (ML) and deep learning (DL) are transformative approaches in protein engineering, enabling the prediction of mutation effects from sequence, structural, and functional datasets. Unlike traditional methods that rely on manual feature engineering, ML/DL automatically extracts relevant features from raw sequences—avoiding bias and missed information—and models complex, nonlinear, high-dimensional sequence–function relationships. This scalability allows for rapid prediction across entire sequence spaces, far outperforming physics-based simulations in speed and efficiency [1-3].

In our project, we focused on Transformer models, a type of DL architecture that uses self-attention mechanisms to capture long-range dependencies between amino acids. This was critical for identifying interactions between distant residues (e.g., catalytic sites and distal loops) that influence protein function—relationships that traditional sequential models (like RNNs) might miss. By training Transformer models on large-scale Mn-SOD datasets, we could detect patterns associated with thermostability and catalytic activity. Our Transformer model has successfully identified variants, which were promising candidate for further experimental validation. The integration of ML/DL into our pipeline not only accelerated the discovery of beneficial mutations but also enhanced the precision of predicting their impact on SOD performance (Figure1).

Figure1. Transformer model architecture.

3 Mathematical Model

In order to scientifically design and evaluate the efficacy of Mn-SOD eye drops in practical application, we constructed a multi-scale mathematical model to quantitatively simulate the efficacy-time relationship in the eyeball, and further verify the feasibility of Mn-SOD eye drops for early prevention of cataracts. Our mathematical models employ a bottom-up, simplified to complex modeling strategy. Firstly, at the molecular biological mechanism (microscopic) level, the catalytic reaction mechanism of Mn-SOD catalytic scavenging of ROS was characterized, and the enzymatic reaction rate equation was derived based on this. On this basis, at the macroscopic 0D system level, the ocular environment is simplified into a homogeneous reaction chamber, and ordinary differential equations (ODEs) are established to simulate the dynamic process of drug and ROS concentrations with time. Further, at the macroscopic 1D spatial expansion level, reaction-diffusion partial differential equations (PDEs) are introduced to simulate the spatial distribution gradient of drugs and ROS in the eye, so as to achieve higher fidelity prediction of the spatiotemporal evolution of drug efficacy.

3.1 Microscopic Catalytic Kinetic Model

Microscopic Catalytic Kinetic Model is based on the Ping-Pong Mechanism in which manganese superoxide dismutase (Mn-SOD) catalyzes the disproportionation reaction of superoxide anion (O2⁻). This mechanism accurately describes the cyclic catalytic process of manganese ions in the enzemic active center between the oxidation state (Mn³⁺-SOD, Eox) and the reducing state (Mn²⁺-SOD, Ered), with each step and a half of the reaction consuming a superoxide molecule and generating the corresponding enzyme-substrate intermediate complex (C1 and C2). The two-and-a-half-step reaction is shown in Figure 2 below. The model derives the total reaction rate equation from the first principles by applying the Briggs-Haldane steady-state approximation (assuming zero rate change in the concentration of the intermediate complex) combined with the flow equilibrium condition of the two-step reaction. The equation can eventually be organized into a Mies-like form, but its parameters are determined by a combination of microscopic rate constants (k₁, k₂, k₃, k₄), and the effective catalytic constant reflects that the rate of the entire catalytic cycle is limited to slower half-reactions [4].

Figure2. Chemical Reaction Mechanism of Mn-SOD Scavenging ROS.

The main method process of this model mainly includes the following steps. Firstly, the reactions of each element and their rate equations are written according to the ping-pong mechanism. Subsequently, a steady-state approximation (d[C1]/dt ≈ 0, d[C2]/dt ≈ 0 was applied to the intermediate complexes C1 and C2, and the concentrations of oxidation and reducing enzymes were expressed as intermediate complex concentrations using the condition that the two halves of the reaction flow rate in the catalytic cycle were equal (k₂[C1] = k₄[C2]). These expressions are then substituted into the total enzyme conservation equation (Figure3). The expression of [C1] is finally solved by algebraic rearranging, and the explicit equation of total ROS clearance rate is obtained accordingly, which clearly expresses the dependence of clearance rate on the substrate concentration S and the total enzyme concentration Et.

Figure3. The total enzyme conservation equation.

3.2 Macroscopic Zero-Dimensional (0D) Homogeneous System Dynamics Model

The model simplifies the aqueous humor environment in the eyeball into an ideal, well-mixed continuous stirred kettle reactor (CSTR), and its core assumption is spatial uniformity, that is, the concentration of drug and ROS in the entire aqueous humor is uniform at any time [5]. The model describes the changes of drug concentration and ROS concentration over time by establishing a coupled ordinary differential equation (ODE) system (Figure 4). Among them, the change of drug concentration is only determined by the purge process (such as degradation and loss), following the exponential decay law; while the change of ROS concentration is controlled by the dynamic balance of its endogenous generation rate and enzyme-catalytic clearance rate. The key model parameters in the zero-dimensional (0D) homogeneous system dynamics model are shown in Table 1.

Figure4. O**rdinary differential equation (ODE) system in Zero-Dimensional (0D) Homogeneous System Dynamics Model*.*

The model construction contains two main ODEs describing the exponential decay of drug concentrations and the dynamic equilibrium of ROS concentrations, where V(S; Et) is the reaction engine embedded from the microscopic model. After obtaining all the necessary parameters (initial dose E₀, dose, ROS background concentration S₀, clearance constant kclear, generation rate constant k_gen, etc.), the ODE system can be solved by numerical methods (such as Euler's method, Runge-Kutta method) to simulate the complete curve of drug and ROS concentrations after a single dose over time.

Table1*. Key Model Parameters* in Zero-Dimensional (0D) Homogeneous System Dynamics Model*.*

3.3 Macroscopic One-Dimensional (1D) Reaction-Diffusion Model

The target that drugs for preventing cataracts need to act on is the lens. Traditional models cannot answer the core questions of "how much drug will reach the lens" and "how long will it take?" In order to directly and dynamically simulate the penetration process of drugs in the anterior chamber, we designed the Macroscopic One-Dimensional (1D) Reaction-Diffusion Model. This model introduces spatial dimensions and diffusion effects. Its principle is based on Fick's second law, which suggests that changes in the concentration of substances in the eye are not only due to local chemical reactions, but also from their diffusion flux between different spatial locations. The model describes this process through the reaction-diffusion partial differential equation (PDE), in which temporal evolution and spatial distribution are coupled (Figure 5). Tthe spherical coordinate system (Figure 6) can more realistically capture the geometric dilution effect caused by the spherical geometry of the eyeball due to the introduction of the (2/r) (∂C/∂r) term, that is, the natural decrease in concentration due to the increase in the area of the sphere through which the material passes when it diffuses outward [6].

Figure5. The reaction-diffusion partial differential equation.

Figure6. T**he spherical coordinate system equation*.*

The main method process of this model mainly includes the following steps. Firstly, the coordinate system is selected according to the geometric characteristics of the simulation domain and the corresponding reaction-diffusion PDE system is established. For the concentration of drugs and ROS C(x,t) or C(r,t), the control equation is ∂C/∂t = D * ∇²C + R(C), where R(C) represents the chemical reaction term (i.e., microscopic engine V). Subsequently, physiologically reasonable initial conditions (e.g., initial distribution of ROS S₀(x) before dosing, drug distribution at the moment of administration) and boundary conditions (e.g., flux input conditions at the corneal boundary to simulate eye drop administration, and no flux conditions at the lens boundary to simulate its physical barrier) were set. Finally, numerical methods (such as finite difference method or finite element method) were used to discretize and solve the PDE to obtain the complete distribution of drug and ROS concentrations in space and time C(x, t).

4 Analysis and Discussion of Model Results

4.1 Rational design of Mn-SOD based on Transformer deep learning model

The identified Mn/FeSOD protein sequences were first retrieved and obtained from IPR001938 entries in the InterPro database. Considering that Mn/FeSOD functions in the form of dimers under physiological conditions, all extracted sequences are repeated to simulate its natural state of action to improve the accuracy of model training. In the process of sequence adaptive distance modeling, a codec structure neural network based on Transformer encoder is adopted. In the training phase, 40% of amino acid residues are randomly masked to enhance the generalization ability of the model. The model training parameters were set as follows: using the Adam optimizer, the learning rate was 3×10⁻⁴, the training round was 20, and the batch size was 32. After the training was completed, all residue sites were sorted according to the logit, and the top 50 residue sites with the highest average scores were selected as potential key candidates for the adaptive evolution of Mn-SOD (Figure7). In the subsequent wet experimental validation, we selected three predicted mutations, namely D107L, K46S, and G93H, for experiments.

Figure7. Transformer model prediction results*.*

4.2 Microscopic Catalytic Kinetic Model

The core purpose of this microscopic model is to build an accurate "biochemical reaction engine". It quantitatively describes the catalytic efficiency of Mn-SOD from the lowest molecular mechanism, providing fundamental and reliable chemical reaction terms. Without this microscopic engine, any macroscopic simulation would lose its foundation in biochemical authenticity. The steady-state approximation and rate equation derivation process of the model are shown in Figure 8.

Figure8. Derivation process of Microscopic Catalytic Kinetic Model*.*

After successfully deriving the relevant equations of the model, we visually compared the catalytic efficiency of natural Mn-SOD and engineered Mn-SOD in the scavenging of reactive oxygen species (ROS) using the Michaelis-Menten curve. The analysis results showed that the maximum reaction rate (Vmax) of engineered Mn-SOD reached 4000.0 μM/s, which was about 2.67 times higher than that of natural Mn-SOD (1500.0 μM/s). This significant improvement means that our engineered variant can perform the role of scavenger with extreme efficiency when oxidative stress occurs in the eye and ROS concentrations rise sharply, quickly reducing ROS concentrations below safe thresholds, thereby effectively protecting lens proteins from oxidative damage, which will theoretically significantly reduce the time window for the protective effect after eye drop administration. Although its Mie constant (Km) is slightly elevated (12.0 μM → 10.0 μM), indicating a subtle change in substrate affinity, in cataract-related pathological settings (often accompanied by significant increases in ROS concentration), the catalytic capacity gain from extremely high Vmax far outweighs the small effects that Km change may have. In summary, the excellent catalytic efficiency of engineered Mn-SOD has laid a solid theoretical foundation for its potential to achieve low-dose and high-efficacy in the development of subsequent ophthalmic drop formulations, and greatly enhanced the confidence of this project from molecular design to practical application (Figure 9).

Figure9. Plot of the comparative Michaelis-Menten kinetics for Natural vs. Engineered *Mn-*SOD.

4.3 Macroscopic Zero-Dimensional (0D) Homogeneous System Dynamics Model

Macroscopic Zero-Dimensional (0D) Homogeneous System Dynamics Model is designed to quickly assess the overall pharmacodynamic profile of dosing regimens and serves as a powerful tool for conducting efficient parameter scanning and preliminary screening. It visualizes the intensity (how low ROS can be suppressed) and duration (how long ROS remains low) for comparing the strengths and weaknesses of different formulations (corresponding to different kclear values), analyzing dose-response relationships, and initially evaluating the dosage adjustment strategies required for different levels of oxidative stress (simulated by adjustment kgen). Its simple and fast calculations make it ideal for high-throughput virtual screening in the early stages of a project.

We systematically evaluated the efficacy boundaries of engineered Mn-SOD eye drops under different pathological conditions. The heat map results visually showed that the final ROS residual concentration (color from purple to yellow indicates high to low) formed a clear efficacy demarcation between different initial Mn-SOD doses (horizontal axis, 0–1 μM) and different oxidative stress levels (vertical axis, 10–70 μM/h). The results showed that ROS could be inhibited at a safe level (yellow region) with very low doses (about 0.2–0.3 μM) of engineered Mn-SOD under low oxidative stress conditions (e.g., physiological aging). In environments with high oxidative stress (e.g., pathological states caused by metabolic diseases), the dose needs to be significantly increased to more than 0.6 μM to achieve effective clearance (breakthrough to the yellow-green region). This simulation result not only quantitatively characterizes the efficacy advantages of engineered Mn-SOD, but more importantly, provides clear dosage window guidance for clinical precision drug use, that is, the need to adopt differentiated drug delivery strategies for different risk groups, thus providing a key theoretical basis for the personalized prevention of cataracts (Figure 10).

Figure10. Heatmap showing efficacy (ROS remaining) for different doses and stress levels.

4.4 Macroscopic One-Dimensional (1D)Reaction-Diffusion Model

The core purpose of this model is to make a leap from "temporal dynamic" prediction to "spatiotemporal dynamic" prediction to obtain simulation results with higher biological realism and prediction accuracy. It can accurately simulate the spatiotemporal process of how drugs penetrate aqueous humor and finally reach the target (lens) after dripping from the cornea, predicting its actual exposure concentration and history on the target surface, so as to provide the final and most reliable decision-making basis for designing precise dosing regimens (dosage and frequency) that can truly ensure the effectiveness of the target. At the same time, its visualizations, such as heat maps, provide an excellent representation of the spatial dynamics of drug penetration and ROS clearance.

We simulated the spatiotemporal dynamics of engineered Mn-SOD eye drops after dropping into the eye and its scavenging effect on reactive oxygen species (ROS). As shown in Figure 11, the model simulation results reveal changes in SOD and ROS concentrations over time (hours) within the eye (from the cornea to the lens, expressed as normalized distance): the distribution of SOD (Figure 11a) shows that its concentration is relatively uniform spatially, with only slight attenuation over time, indicating that the drug can effectively diffuse and maintain therapeutic concentrations; The dynamics of ROS (Figure 11b) showed a clear spatiotemporal gradient, with a higher initial concentration in the near lens region (bottom), but gradually decreased with the continuous effect of SOD, especially in the later stage of administration. This result confirms that engineered Mn-SOD can efficiently penetrate into the target area and achieve spatiotemporal specific clearance of ROS, which provides a key basis for the optimization of the frequency and dosage of eye drops, and verifies its feasibility and reliability as a cataract prevention strategy.

Figure11. Generates a conceptual plot of the spatiotemporal dynamics from a 1D-PDE model.

5 References

​ [1] Alley E C, et al. Unified rational protein engineering with sequence-based deep representation learning. Nat Methods, 2019, 16(12): 1315-1322.

​ [2] 蒋迎迎, 曲戈, 孙周通. 机器学习助力酶定向进化. 生物学杂志, 2020, 37(4): 1-11.

​ [3] Cui Y, Chen Y, Sun J, et al. Computational redesign of a hydrolase for nearly complete PET depolymerization at industrially relevant high-solids loading. Nat Commun, 2024, 15(1): 1417.

​ [4] Truscott RJ. Age-related nuclear cataract: a lens transport problem. Ophthalmic Res. 2000 Sep-Oct;32(5):185-94.

​ [5] Davies NM. Biopharmaceutical considerations in topical ocular drug delivery. Clin Exp Pharmacol Physiol. 2000 Jul;27(7):558-62.

​ [6] Guo X, Wen H, Hao H, Zhao Y, Meng Y, Liu J, Zheng Y, Chen W, Zhao Y. Randomness-Restricted Diffusion Model for Ocular Surface Structure Segmentation. IEEE Trans Med Imaging. 2025 Mar;44(3):1359-1372.