Method for simulating movement patterns of a protein molecular machine
By combining molecular dynamics simulations and deep learning with VR technology, the binding of the S protein in the RBD region of the Omicron variant SARS-CoV-2 virus to ACE2 was accurately simulated and predicted. This solved the problem of insufficient simulation in existing technologies, revealed the mechanism by which mutations lead to enhanced infectivity, and supported medical diagnosis and treatment.
Patent Information
- Application Number
- CN202211404471.X
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2022-11-10
- Publication Date
- 2026-05-15
- Estimated Expiration
- 2042-11-10
AI Technical Summary
Existing technologies are insufficient to accurately simulate and predict the binding of the S protein in the RBD region of the Omicron variant SARS-CoV-2 virus to the human receptor ACE2, resulting in inadequate analysis of its infectivity and antibody resistance.
By combining molecular dynamics simulation with deep neural networks, structural data of protein molecular machines are acquired, molecular dynamics simulation and binding free energy calculation are performed, NRI models are built to predict motion trajectories, and VR devices are used to observe the three-dimensional structure and motion trajectory.
This study achieved accurate simulation and prediction of the binding ability of the Omicron variant S protein to ACE2, revealing the mechanism by which mutations in the RBD region lead to enhanced infectivity, and providing a basis for medical diagnosis and treatment.
Smart Images
Figure CN115662496B_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of vehicle control technology, and in particular to a method and apparatus for online estimation of tire-road adhesion coefficient. Background Technology
[0002] Protein molecular machines are machines composed of molecular-scale substances capable of performing certain processing functions. Their main components are biomolecules such as proteins, and they participate in numerous important life processes. Research on the structure, function, and movement patterns of protein machines is fundamental to understanding these high-dimensional, complex molecular machines and is of great significance for revealing the essence of life phenomena. The latest SARS-CoV-2 variant, Omicron, has caused global panic due to its infectivity and vaccine escape mutations. The infectivity and antibody resistance of SARS-CoV-2 variants depend on mutations in the receptor-binding domain (RBD) of the spike protein. Omicron exhibits numerous mutations in the RBD of the spike protein, which has attracted significant attention from the scientific community and the public. Exploring the binding of the Omicron variant's S protein to the human receptor ACE2 has thus become particularly important. Summary of the Invention
[0003] To address the aforementioned problems, this invention provides a method for simulating the motion patterns of protein molecular machines, which combines molecular dynamics simulation and deep neural networks to achieve accurate simulation and prediction of the motion patterns of protein molecular machines.
[0004] The technical solution provided by this invention is as follows:
[0005] A method for simulating the motion patterns of a protein molecular machine, comprising:
[0006] S10 acquires a protein molecular machine structure data file, wherein the protein molecular machine includes a complex of wild-type ligand and the RBD region S protein of Omicron Ba.2 with receptor ACE2.
[0007] S20 simulates the molecular dynamics of the protein molecular machine to obtain the topology and coordinate files of the complex and its receptor and ligand motion, and calculates its binding free energy.
[0008] S30 predicts the motion trajectory of a protein molecular machine based on the constructed NRI model and coordinate file. The NRI model includes an encoder for predicting the trajectory interaction of a given dynamic system and a decoder for predicting the trajectory of a given dynamic system in the interaction graph. In the NRI model, the input is the motion coordinates of the complex and the output is the predicted motion coordinates of the complex.
[0009] S40 uses VR equipment to observe and record the three-dimensional structure and motion trajectory of the protein molecular machine.
[0010] The present invention provides a method for simulating the motion patterns of protein molecular machines. This method employs molecular dynamics simulations to explore the dynamic interaction between the S protein in the RBD region of the ligand Omicron Ba.2 and the receptor ACE2, comparing it with wild-type and ACE2 systems. Simulations revealed that Omicron exhibits a stronger binding affinity to human cells, and many important mutations occur in the RBD region in contact with ACE2, indicating that mutations in the RBD may lead to enhanced infectivity, providing some basis for subsequent medical diagnosis and treatment. Attached Figure Description
[0011] The preferred embodiments will now be described in a clear and easy-to-understand manner, with reference to the accompanying drawings, to further explain the above-mentioned characteristics, technical features, advantages, and implementation methods.
[0012] Figure 1 This is a schematic diagram of the motion mode simulation method for protein molecular machines in this invention. Detailed Implementation
[0013] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the specific implementation methods of the present invention will be described below with reference to the accompanying drawings. Obviously, the drawings described below are merely some embodiments of the present invention. For those skilled in the art, other drawings and other implementation methods can be obtained based on these drawings without any creative effort.
[0014] Molecular dynamics is a computer simulation method based on molecular mechanics, primarily used to simulate and observe the motion and changes of molecules under different environments. Molecular dynamics simulation is the closest simulation method to experimental conditions among molecular simulations, capable of providing the microscopic evolution of a system at the atomic level, intuitively demonstrating the mechanisms and laws governing experimental phenomena. This drives research towards greater efficiency, economy, and predictability; therefore, molecular dynamics simulation plays an increasingly important role in research in biology, pharmacy, chemistry, and materials science.
[0015] Deep learning has made remarkable progress in recent years, achieving impressive results in models, algorithms, and large-scale applications. The emergence of deep learning represents a significant revolution in machine learning and a powerful driving force for the development of artificial intelligence. Traditional neural networks represent a shallow form of machine learning, while deep learning is a new generation of neural networks developed from traditional ones. Deep learning achieves feature extraction from external input data, from low-level to high-level, by building and simulating the information processing neural structures of the human brain, thereby enabling machines to understand and learn from the data and acquire information.
[0016] VR, relying on computer platforms and advanced electronic technology, creates a highly realistic simulated environment. By applying various sensing devices, it further enhances the user experience, allowing users to interact with objects in the virtual space in a more intuitive and tangible way, thus optimizing the simulation experience. This experiment utilizes VR equipment to more intuitively view the three-dimensional structure of the protein machine and better understand its movement patterns.
[0017] Based on this, the present invention provides a method for simulating the motion patterns of protein molecular machines, such as... Figure 1 As shown, it includes:
[0018] S10 acquires the structural data file of the protein molecular machine, which includes the complex of the wild-type ligand and the RBD region of the OmicronBa.2 S protein with the receptor ACE2.
[0019] S20 simulates the molecular dynamics of protein molecular machines, obtaining topological and coordinate files of the complex and its receptor and ligand motion, and calculates their binding free energy.
[0020] S30 predicts the motion trajectory of protein molecular machines based on the constructed NRI model and coordinate files. The NRI model includes an encoder for predicting the trajectory interaction of a given dynamic system and a decoder for predicting the trajectory of a given dynamic system in the interaction graph. In the NRI model, the input is the motion coordinates of the complex and the output is the predicted motion coordinates of the complex.
[0021] The S40 uses VR equipment to observe and record the three-dimensional structure and motion trajectory of protein molecular machines.
[0022] Specifically, in step S10, complexes of the wild-type (pre-mutation) and Omicron Ba.2 (Omicron variant) RBD region S proteins with ACE2 need to be prepared separately. The structural files of these complexes are all from the Protein Data Bank (PDB) website. The complex formed with Omicron Ba.2 is from pdbID 7ZF7, while the complex with wild-type is from pdbID 6LZG. Furthermore, Pymol software is needed to split the PDB files of the complexes into receptor and ligand components to facilitate subsequent calculations of binding free energy (in the two complexes mentioned above, the wild-type and Omicron Ba.2 RBD region S proteins are ligands, and ACE2 is the receptor). Before performing molecular dynamics simulations, the complexes, receptors, and ligands need to be preprocessed using UltraEdit, retaining only the Atom field information and processing their disulfide bonds.
[0023] In step S20, molecular dynamics simulations of the protein molecular machine are performed using Amber to obtain the trajectory file of the complex. Specifically, the software Amber16 and the force field of leapc.ff03ua are used to simulate the wild-type and Omicron Ba.2 / ACE2 complex systems, respectively: First, molecular construction is performed under the Linux system to obtain the topology and coordinate files of the complex, receptor, and ligand (as a reference for subsequent trajectory prediction). During this process, the complex is centered in a cubic box, and appropriate ions are added to the system to electrostatically neutralize the molecular system. The complex is then dissolved using the tip3pbox water model to obtain the topology and coordinate files of the solvent-based complex. Afterward, energy minimization is performed on the solvent-based complex, and after a 50 ps temperature rise, the complex is subjected to 50 ps density equilibrium under weak confinement conditions, followed by 500 ps isobaric equilibrium at 300 K. Throughout the simulation, SHAKE constraints are applied to hydrogen atoms, and a 2 fs time step and Langevin kinetics are used to control the temperature. The root mean square deviation (RMSD) of the protein backbone relative to the energy-minimized structure is calculated to determine whether the conformation has stabilized during equilibrium. Here, whether the conformation is stable during equilibrium is determined by whether the calculated RMSD fluctuates within a certain range. For example, if the value starts from 0 and then fluctuates within a certain range after a period of time... arrive Fluctuations within a certain range indicate that the conformation is stable during equilibrium. Sudden increases or decreases indicate that the conformation is not yet stable.
[0024] After assessing whether the molecular dynamics simulation brought the solvent complex to equilibrium based on RMSD properties, the energy term of the S protein-ACE2 interaction was further calculated using a script. This method is based on molecular mechanics with Poisson-Boltzmann surface area (MM-PBSA).
[0025] In ideal real-world scenarios, the binding free energy can be directly calculated using the binding of acceptors and ligands in a solvent. However, in simulations of these solvated states, most of the energy contribution comes from solvent-solvent interactions, and the total energy fluctuation will be an order of magnitude larger than the binding energy. Therefore, the calculation of the binding energy requires a very long time to converge. A more efficient method is to calculate the various parts of the thermodynamic cycle separately. A simplified process involves first binding the acceptor and ligand in the gas phase, and then placing the complex in a solvent environment. Since the initial and final states of both paths are the same, from an energy conservation perspective, the binding free energy ΔG is much larger. Solv,bind For equation (1):
[0026] ΔG Solv,bind =ΔG Gas,bind +ΔG Solv,comp -(ΔG Solv,Lig +ΔG Solv,Rec(1)
[0027] Wherein, ΔG Gas,bind ΔG represents the energy of ligand-receptor binding in a solvent. Solv,comp ΔG represents the energy required for the conversion of a complex between its gaseous and solvent states. Solv,Lig ΔG represents the energy of ligand conversion between gaseous and solvent states. Solv,Rec This represents the energy required for the acceptor to switch between gaseous and solvent states.
[0028] From the perspective of combined energy, combined with free energy ΔG bind This can be further expressed as equation (2):
[0029] ΔG bind =ΔH-TΔS(2)
[0030] Here, ΔH represents the enthalpy change, which can be decomposed into gas phase energy ΔE. MM and solvation energy ΔG sol Ignoring the effect of entropy (-TΔS) on energy ΔG bind The contribution of enthalpy change ΔH is expressed by equations (3) to (5):
[0031] ΔH=ΔE MM +ΔG sol (3)
[0032] ΔE MM =ΔE ele +ΔE vdW +ΔE int (4)
[0033] ΔG sol =ΔG pb +ΔG np (5)
[0034] Where, ΔE ele Represents the electrostatic term, ΔE vdW Denotes the van der Waals term, ΔE int Represents the internal energy term (ΔE) int Internal energy includes bond energy, bond angle energy, and dihedral angle energy, ΔG pb Represents the polar solvation energy, ΔG np This represents the nonpolar solvation energy.
[0035] In addition to calculating the binding free energy, after obtaining the trajectory file based on the above kinetic simulation process, the root mean square fluctuation (RMSF), solvent reachable surface area (SASA), and minimum residue distance (Mindist) are further calculated and analyzed using Amber software. RMSF represents the distance a residue has moved from its original position; a smaller movement indicates better binding, a stronger binding free energy, and a higher solvent reachable surface area (SASA). Similarly, a smaller minimum residue distance (Mindist) indicates a stronger binding free energy. The residues are considered contact residues. Furthermore, hydrogen bonding was analyzed using VMD tools to determine the hydrogen bonding situation, where the distance between the donor and acceptor atoms was less than [a certain value]. Furthermore, the standard is that the hydrogen atom formation angle of the donor-acceptor-donor connection is less than 45 degrees. The ligplot software is used to analyze the protein-ligand 2D structure to help explain the conclusion that the binding free energy is stronger after mutation.
[0036] After calculating the binding free energy, in step S30, the motion trajectory of the protein molecular machine is predicted based on the constructed NRI model and coordinate file. The NRI model includes an encoder for predicting the trajectory interaction of a given dynamic system and a decoder for predicting the trajectory of a given dynamic system in the interaction graph. In the NRI model, the input is the motion coordinates of the complex, and the output is the predicted motion coordinates of the complex. Specifically, the input of the NRI model consists of N nodes, and the feature vector (position and velocity in the x, y, and z dimensions) of node i is represented as X at time t. i t The feature set of all N nodes is represented as X. t ={X1 t ,...,X t N The trajectory of node i is represented as X. i ={X i 1 ,...,X i t}, where t represents the time step. Finally, all trajectory data are recorded as X = {X} 1 ,...,X t The NRI model, based on an unknown graph z, learns edge values and reconstructs the future trajectory of a dynamic system in an unsupervised manner. The interaction between nodes i and j is represented by the latent variable Zj. i The interaction types are modeled in the form {1,...,K}, where K is the number of interaction types. These interaction types do not have any predefined meanings, but the model learns to assign a meaning to each type.
[0037] After obtaining the predicted trajectory file in step S30, place it in a subfolder of the Assets folder, and then categorize and place the resources into these folders. Here, the FBX format model file is used, which allows for the simultaneous export and use of animation, materials, and cell attributes along with the model. The 3D scene and character models are optimized to achieve real-time rendering of 3D meshes in the interactive system. Specifically, the Xfrog Organic modeling tool is used to create the entire process of wild-type and Omicron S proteins binding to human ACE2, improving the efficiency of 3D animation production.
[0038] The VR device includes: a physical space positioning module, used to locate the user's relative position in real space through a head-mounted sensor to obtain the user's spatial position data; a VR environment interaction module, used to generate interaction information based on the spatial position obtained by the physical space positioning module and the interaction actions captured by the handheld sensor; a message processing module, used to process the interaction information generated by the VR environment interaction module; and a VR virtual environment simulation module, used to update the parameters of scene objects according to the interaction information generated by the message processing module, complete the rendering output of a new frame, and feed back the virtual environment interaction information to the user's sensors and the physical space positioning module.
[0039] During operation, the physical space positioning module uses a head-mounted sensor to locate the user's relative position in real space, acquiring the user's spatial location data. After acquiring spatial location and interaction behavior data, this data is passed to the interaction system for interpretation. The location information in real space, after coordinate transformation and conversion, is mapped to different location nodes in the corresponding virtual environment. Interaction action data captured by the handheld sensor is uniformly converted into different types of message structures and sent to various message processing modules. Finally, the interaction behavior feedback module transmits the corresponding interaction response data and processing events to the VR virtual environment simulation system to update the parameters of scene objects, complete the rendering output of a new frame, and feed back the virtual environment interaction information to the user's sensors and spatial positioning device, thereby allowing for a clearer observation and perception of the protein motion trajectory after the completion of molecular dynamics simulation.
[0040] In this invention, molecular dynamics simulations were used to explore the kinetic interaction between ACE2 and the RBD, and compared with wild-type and ACE2 systems. The results showed that Omicron exhibited a stronger binding ability to human cells, and many important mutations occurred in the RBD in contact with ACE2, indicating that mutations in the RBD may lead to enhanced infectivity. Simultaneously, the binding properties of the receptor ACE2 with the wild-type and Omicron Ba.2 RBD region S proteins were studied using molecular dynamics methods. Analysis was performed on hydrogen bonds, the number of contact residues, and the solvent-accessible surface area, and the binding free energy of their RBD region S proteins to ACE2 was calculated.
[0041] In predicting motion trajectories, an NRI model was built, and the optimal prediction model was selected based on the accuracy of the training and test sets. Specifically, Anaconda software was used, and a deep learning system was built using the Python language. Python is an object-oriented interpreted computer programming language with cross-platform compatibility. The most significant characteristic of a BP neural network is that the signal propagates forward, while the error propagates backward. This can be understood as continuously adding input data as a forward propagation signal, and observing and feeding back the results to better optimize the model; this is understood as the backward propagation of the error. Its learning rule uses the steepest descent method, continuously adjusting the network's weights and thresholds through backpropagation to minimize the sum of squared errors.
[0042] Protein allostericity is a biological process facilitated by long-range intraprotein communication in space, where ligand binding or changes in distal amino acids remotely affect active sites. Molecular dynamics (MD) simulations provide a powerful computational approach for exploring allosteric effects. However, current MD simulations cannot reach the timescale of the entire allosteric process. The advent of deep learning has made it possible to assess both short-range and long-range communication in space to understand allostericity. To this end, a neural relational inference model based on graph neural networks was applied. This model employs an encoder-decoder architecture to infer potential interactions, probing the protein allosteric process as a dynamic network of interacting residues. From MD trajectories, the model successfully learned long-range interactions and pathways. Furthermore, the model can detect allosteric-related interactions earlier in MD simulation trajectories and predict relative free energy changes at mutations more accurately than other methods.
[0043] Molecular dynamics simulations can directly probe the motion of biomolecules, but due to the finite timescale of the simulations and the high dimensionality and complexity of 3D trajectory data, they may fail to capture meaningful functional information. Furthermore, many challenging molecular dynamics (MD) analysis problems lack suitable methods for probing long-range communication. Computational techniques for simulating protein allosteric communication rely on graph-theoretic metrics to identify long-range couplings between two distal active sites. Typically, proteins can be mapped to a graph where each node represents a residue and each weighted edge represents an interaction between two nodes. The shortest path between allosteric sites and active sites in a protein can be important for signal propagation in allosteric communication. Early graphical models used static crystal structures to calculate the shortest path between one residue and other residues, which may not account for all potential contacts and associated allosteric behavior in dynamic proteins. A later, well-known allosteric approach, Perturbation Response Scan (PRS), uses a Hessian-based Elastic Network Model (ENM) to obtain the relevant dynamics of positions. This model studies how a perturbation on a single residue triggers a cascade of perturbations (signals) to other nodes in the elastic network, thus enabling allosteric communication. To more accurately simulate the response to ligand binding or mutation, the reciprocal of the Hessian model was replaced with a covariance matrix that incorporates the dynamic characteristics of the system.
[0044] The NRI model is suitable for learning simulated trajectories of biomolecules in molecular dynamics (MD) simulations, where biomolecules are formed by chemically bonded atoms whose motion is described by Newtonian mechanics. The model uses a generative neural network (GNN) to learn the network's dynamic embeddings by minimizing the reconstruction error between the reconstructed and simulated trajectories; the NRI model then infers the edges between the residuals, represented by latent variables. The learned embeddings essentially abstract the fundamental roles of key residues in conformational transitions, which helps decipher the mechanisms of protein allosteric changes.
[0045] Unity3D is a professional cross-platform game development and virtual reality engine developed by Unity Technologies, supporting multiple scripting languages such as C# and JavaScript. This invention uses the widely acclaimed HTC Vive VR device, combined with the Unity engine as the development platform, to develop a VR simulation platform for biological control. Ultimately, VR is used to observe cell structure and its movement trajectory.
[0046] It should be noted that the above embodiments can be freely combined as needed. The above are merely preferred embodiments of the present invention. It should be pointed out that for those skilled in the art, several improvements and modifications can be made without departing from the principle of the present invention, and these improvements and modifications should also be considered within the scope of protection of the present invention.
Claims
1. A method for simulating the motion patterns of a protein molecular machine, characterized in that, include: S10 Obtain the protein molecular machine structure data file, wherein the protein molecular machine includes a complex of wild-type ligand and the RBD region S protein of OmicronBa.2 with receptor ACE2. S20 simulates the molecular dynamics of the protein molecular machine to obtain the topology and coordinate files of the complex and its receptor and ligand motion, and calculates its binding free energy; S30 predicts the motion trajectory of a protein molecular machine based on the constructed NRI model and coordinate file. The NRI model includes an encoder for predicting the trajectory interaction of a given dynamic system and a decoder for predicting the trajectory of a given dynamic system in the interaction graph. In the NRI model, the input is the motion coordinates of the complex and the output is the predicted motion coordinates of the complex. S40 Observe and record the three-dimensional structure and motion trajectory of the protein molecular machine based on VR equipment; In step S20, Amber is used to simulate the molecular dynamics of the protein molecular machine to obtain the trajectory file of the complex; The software Amber16 and the force field of leaprc.ff03ua were used to simulate the wild-type and Omicron Ba.2 and ACE2 complex systems, respectively; In step S40, the VR device generates the corresponding fbx file based on the pdb trajectory file generated by the NRI model, and displays the entire process of the binding of the wild-type and Omicron Ba.2 RBD region S protein to the receptor ACE2.
2. The method for simulating the motion patterns of protein molecular machines as described in claim 1, characterized in that, In step S20, the free energy ΔG is combined with Solv,bind It can be represented as: ΔG Solv,bind =ΔG Gas,bind +ΔG Solv,comp -(ΔG Solv,Lig +ΔG Solv,Rec ) Wherein, ΔG Gas,bind ΔG represents the energy of ligand-receptor binding in a solvent. Solv,comp ΔG represents the energy required for the conversion of a complex between its gaseous and solvent states. Solv,Lig ΔG represents the energy of ligand conversion between gaseous and solvent states. Solv,Rec This represents the energy required for the acceptor to switch between gaseous and solvent states.
3. The method for simulating the motion patterns of protein molecular machines as described in claim 1, characterized in that, In step S20, the binding free energy ΔG bind It can also be expressed as: ΔG bind = ΔH-TΔS Here, ΔH represents the enthalpy change, which can be decomposed into gas phase energy ΔE. MM and solvation energy ΔG sol Ignoring the effect of entropy TΔS on energy ΔG bind The contribution of enthalpy change ΔH is expressed by the formula: ΔH = ΔE MM + ΔG sol ΔE MM = ΔE ele + ΔE vdW + ΔE int ΔG sol = ΔG pb + ΔG np Where, ΔE ele Represents the electrostatic term, ΔE vdW Denotes the van der Waals term, ΔE int Let ΔG represent the internal energy term. pb Represents the polar solvation energy, ΔG np This represents the nonpolar solvation energy.
4. The method for simulating the motion patterns of protein molecular machines as described in claim 1, characterized in that, The VR device includes: The physical space positioning module is used to locate the user's relative position in real space through a head-mounted sensing device and obtain the user's spatial position data. The VR environment interaction module is used to generate interactive information based on the spatial position obtained by the physical space positioning module and the interactive actions captured by the handheld sensor. The message processing module is used to process the interactive information generated by the VR environment interaction module; The VR virtual environment simulation module is used to update the parameters of scene objects based on the interaction information generated by the message processing module, complete the rendering output of a new frame, and feed back the virtual environment interaction information to the user's sensors and physical space positioning module.