Charged ion ionization free energy acquisition method, system, equipment and medium

Abstract

The invention provides a charged ion ionization free energy obtaining method, system and device and a medium, and belongs to the field of molecular dynamics simulation computation.The method comprises the following steps that physical parameters of charged ions are obtained, and a charged ion physical model is constructed according to the physical parameters; constructing a hydration environment of charged ions according to the charged ion physical model; based on the hydration environment of the charged ions, gradually increasing the layering number of the hydration process, performing free energy disturbance molecular dynamics simulation on the hydration process of the charged ions every time the layering number is increased, and performing bidirectional sampling simulation by adopting different sampling step numbers in the free energy disturbance molecular dynamics simulation process, analyzing a bidirectional sampling simulation result by using a BAR method, and determining an optimal sampling step number and an optimal layering number; and performing free energy disturbance molecular dynamics simulation on the hydration process of the charged ions by using the optimal sampling step number and the optimal layering number, and calculating the hydration free energy of the charged ions.

Description

Technical Field

[0001] The present invention belongs to the field of molecular dynamics simulation calculations, and particularly relates to a method, system, device and medium for obtaining the hydration free energy of charged ions. Background Art

[0002] Molecular Dynamics (MD) simulation is a computational method widely used in the fields of biophysics and computational chemistry to study the dynamic behavior of molecular structure systems. When calculating the hydration free energy of charged ions, MD simulation calculates the interaction energy between ions and water molecules by simulating the movement of ions in aqueous solution. To obtain reliable statistical results, traditional MD simulation methods require long-time simulation, which leads to high consumption of computing resources and high computing costs. Especially when dealing with large-scale and complex biomolecular systems, the computing time may be as long as several weeks or even months. In addition, in complex systems, especially for weak-affinity binders, traditional MD simulation may not be able to fully sample all conformational spaces, thus affecting the result accuracy. For example, during the protein-ligand binding process, some low-energy conformations may be difficult to be sampled due to energy barriers. At the same time, the MD simulation results highly depend on the selected force field parameters, and different force fields may produce significantly different results, which poses new challenges to the computing reliability and repeatability.

[0003] The Monte Carlo (MC) method is another commonly used method for calculating free energy. The MC method explores the conformational space through random sampling and calculates the free energy in different states. However, although the MC method can effectively explore the conformational space, its sampling efficiency is still limited when dealing with complex systems, especially in low-probability regions. For example, during the protein folding process, some intermediate states are difficult to be sampled due to energy barriers.

[0004] In summary, the existing methods for calculating the hydration free energy of charged ions have problems such as large computational amount, low accuracy of calculation results and low computational efficiency. Summary of the Invention

[0005] In order to overcome the deficiencies of low accuracy of calculation results and low computational efficiency existing in the above-mentioned prior art, the present invention provides a method for obtaining the hydration free energy of charged ions, including the following steps:

[0006] Obtain the physical parameters of the charged ions, and construct a physical model of the charged ions according to the physical parameters; construct the hydration environment of the charged ions according to the physical model of the charged ions;

[0007] Based on the hydration environment of the charged ions, the number of layers in the hydration process is successively increased. Each time the number of layers increases, a free energy perturbation molecular dynamics simulation is performed on the hydration process of the charged ions. During the free energy perturbation molecular dynamics simulation, two-way sampling simulations are carried out with different numbers of sampling steps, and then the BAR method is used to analyze the results of the two-way sampling simulations to determine the optimal number of sampling steps and the optimal number of layers;

[0008] Perform a free energy perturbation molecular dynamics simulation on the hydration process of the charged ions using the optimal number of sampling steps and the optimal number of layers, and calculate the hydration free energy of the charged ions.

[0009] Preferably, after calculating the hydration free energy of the charged ions, it further includes correcting the obtained hydration free energy of the charged ions. Specifically: using a three-dimensional periodic simulation box to correct the free energy results of the finite-size molecular dynamics simulation box to obtain the corrected hydration free energy of the charged ions.

[0010] Preferably, the free energy perturbation molecular dynamics simulation further includes setting the periodic boundary conditions, equilibrium time, sampling step size, cutoff radius, and initial temperature of the free energy perturbation molecular dynamics simulation. The periodic boundary is set in the X, Y, and Z directions.

[0011] Preferably, the physical parameters of the charged ions include: the mass, charge, van der Waals force, and electrostatic force parameters of the charged ions.

[0012] Preferably, the construction of the physical model of the charged ions according to the physical parameters includes the following steps:

[0013] Construct a PDB file according to the physical parameters of the charged ions;

[0014] Import the PDB file into the Molecule File Browser module in VMD software, and use the AutoPSF function of the Molecule File Browser module to generate a PSF file;

[0015] Import the PDB file and the PSF file into VMD software simultaneously to obtain a separate physical model of the charged ions.

[0016] Preferably, the construction of the hydration environment of the charged ions is specifically: in VMD software, select the TK Console menu item in Extensions, enter package require solvate in the VMD TK Console window, and wait for VMD software to load the solvate module; after loading is completed, enter the command solvate Li.psfLi.pdb -t 50 -o Li-solvate to obtain A cubic water box.

[0017] The present invention also provides a charged ion hydration free energy acquisition system, including:

[0018] A model construction model for obtaining physical parameters of charged ions, constructing a physical model of charged ions according to the physical parameters; constructing a hydration environment of charged ions according to the physical model of charged ions;

[0019] An optimal parameter acquisition module for successively increasing the number of layers in the hydration process based on the hydration environment of the charged ions. Each time the number of layers increases, a free energy perturbation molecular dynamics simulation is performed on the hydration process of the charged ions. During the free energy perturbation molecular dynamics simulation, bidirectional sampling simulations are performed with different numbers of sampling steps, and then the BAR method is used to analyze the results of the bidirectional sampling simulations to determine the optimal number of sampling steps and the optimal number of layers;

[0020] A hydration free energy calculation module for performing a free energy perturbation molecular dynamics simulation on the hydration process of charged ions using the optimal number of sampling steps and the optimal number of layers, and calculating the hydration free energy of the charged ions.

[0021] The present invention also provides a computer device, including a memory and a processor; the memory stores a computer program, and the processor is used to run the computer program in the memory to execute the method for obtaining the hydration free energy of charged ions.

[0022] The present invention also provides a computer-readable storage medium, which stores a computer program, and the computer program is suitable for being loaded by a processor to execute the method for obtaining the hydration free energy of charged ions.

[0023] The method for obtaining the hydration free energy of charged ions provided by the present invention has the following beneficial effects:

[0024] In the process of free energy perturbation molecular dynamics simulation of the present invention, by successively increasing the number of layers in the hydration process, and performing a free energy perturbation molecular dynamics simulation on the hydration process of charged ions each time the number of layers increases, the number of layers can be optimized, so as to obtain a more optimal number of layers. Performing a free energy perturbation molecular dynamics simulation at this number of layers and calculating the hydration free energy of charged ions can minimize the consumed computing resources while ensuring the calculation accuracy, and can significantly reduce the calculation time cost; by optimizing the number of sampling steps and using bidirectional sampling simulation, it can ensure that all possible states are fully sampled, significantly reduce the calculation error, and improve the calculation accuracy. Description of the Drawings

[0025] To more clearly illustrate the embodiments of the present invention and their design concepts, the accompanying drawings required for this embodiment will be briefly introduced below. The drawings in the following description are only partial embodiments of the present invention. For those of ordinary skill in the art, without creative efforts, other drawings can be obtained based on these drawings.

[0026] Figure 1 It is a schematic flow chart of the present invention;

[0027] Figure 2 It is a structural diagram of the initial solution constructed by the present invention;

[0028] Figure 3 It is the influence of different numbers of layers on the free energy results;

[0029] Figure 4 It is the influence of different sampling steps on the free energy change;

[0030] Figure 5 It is the variation trend of the free energy of lithium ion hydration with λ; Specific Embodiments

[0031] To enable those skilled in the art to better understand the technical solutions of the present invention and implement them, the present invention will be described in detail below with reference to the accompanying drawings and specific embodiments. The following embodiments are only used to more clearly illustrate the technical solutions of the present invention and cannot be used to limit the protection scope of the present invention.

[0032] In the description of the present invention, it should be understood that the terms "center", "longitudinal", "transverse", "length", "width", "thickness", "upper", "lower", "front", "rear", "left", "right", "vertical", "horizontal", "top", "bottom", "inner", "outer", "axial", "radial", "circumferential", etc. indicate the orientation or positional relationship based on the orientation or positional relationship shown in the drawings. It is only for the convenience of describing the technical solutions of the present invention and simplifying the description, rather than indicating or implying that the device or element referred to must have a specific orientation, be constructed and operated in a specific orientation, and therefore cannot be understood as a limitation of the present invention.

[0033] In addition, terms such as "first" and "second" are for descriptive purposes only and should not be construed as indicating or implying relative importance. In the description of the present invention, it should be noted that unless otherwise clearly specified or defined, the terms "connected" and "connected to" should be understood in a broad sense. For example, it can be a fixed connection, a detachable connection, or an integral connection; it can be a mechanical connection or an electrical connection; it can be directly connected or indirectly connected through an intermediate medium. For those of ordinary skill in the art, the specific meanings of the above terms in the present invention can be understood according to specific circumstances. In the description of the present invention, unless otherwise stated, the meaning of "a plurality of" is two or more, which will not be elaborated here.

[0034] Embodiment

[0035] The present invention provides a method for obtaining the hydration free energy of charged ions, specifically as Figure 1 shown, including the following steps:

[0036] Step 1: Obtain the physical parameters of the charged ions, construct a physical model of the charged ions according to the physical parameters; construct a hydration environment of the charged ions according to the physical model of the charged ions.

[0037] In this embodiment, the selected charged ion is a lithium ion, and the physical parameters of the lithium ion include the mass charge, van der Waals force, and electrostatic force parameters of the lithium ion.

[0038] Construct a lithium ion structure file and a coordinate file according to the physical parameters, so as to construct a physical model of the charged ions, specifically including the following steps:

[0039] (1) Download the CHARMM36 force field file, find the lithium ion-related parameters in the toppar_ions_won.str file under toppar / stream / misc, and construct a PDB (Protein Data Bank) file.

[0040] (2) Import the constructed PDB file in the Molecule File Browser module of the VMD (Visual Molecular Dynamics) software, and use the AutoPSF function to generate a PSF (Protein Structure File) file.

[0041] (3) Import the PDB file and the PSF file into the VMD software together, and a separate lithium ion model can be seen.

[0042] The visualization of the lithium ion physical model can be realized by importing the PDB file and the PSF file of the charged ions into the visualization software VMD software.

[0043] To simulate a more realistic environment for the presence of lithium ions and determine the size of the water box to be generated, the present invention constructs in the Solvate module of the VMD software a cubic water box, that is, constructs a hydration environment. The specific process is as follows: Select the TK Console menu item in Extensions, and enter in the VMD TK Console window: package require solvate, and wait for VMD to load the solvate module; after loading is completed, enter the command: solvate Li.psf Li.pdb -t 50 -o Li-solvate, and then obtain Li-solvate.psf and Li-solvate.pdb, generating the solvated PDB file and PSF file. The hydration model is as Figure 2 shown. The lithium ions are in the water box. The entire water box system has 96598 atoms, 64401 chemical bonds, and 32199 water molecules.

[0044] Step 2: Based on the hydration environment of the charged ions, gradually increase the number of layers in the hydration process. Each time the number of layers increases, perform a free energy perturbation molecular dynamics simulation on the hydration process of the charged ions. And during the free energy perturbation molecular dynamics simulation, perform two-way sampling simulations with different sampling steps, and then use the BAR (Bennett Acceptance Ratio) method to analyze the results of the two-way sampling simulations to determine the optimal sampling steps and the optimal number of layers.

[0045] (1) Prepare the required force field parameter files, including periodic boundary conditions, equilibration time, simulation sampling step size, and hydration environment temperature, etc. During the simulation, the lithium ions are not fixed and will move out of the boundary under the action of force, resulting in an out-of-bounds phenomenon. Therefore, the present invention sets periodic boundaries in all three directions of X, Y, and Z; the equilibration time is 1 ns, the sampling step size is 1 fs, the cutoff radius is The initial temperature is 298.15 K.

[0046] (2) Run the NAMD (Nanoscale Molecular Dynamics) software for energy minimization, relaxation, and forward and backward free energy molecular dynamics simulations. Free energy perturbation molecular dynamics simulation FEP (Free Energy Perturbation) is a method for calculating the free energy difference between different states. In FEP, the free energy difference between two states is obtained by calculating the free energy between the reference state λ = 0 and the target state λ = 1. According to perturbation theory, the change in free energy is expressed as:

[0047]

[0048] where kB where \(k\) is the Boltzmann constant, \(T\) is the temperature, and the Hamiltonian \(H\) can be expressed as: \(H(\lambda)=H_0 + V\), where \(H_0\) is the Hamiltonian of the reference state system; \(V\) is the perturbation term; \(\lambda\) is a perturbation parameter between 0 and 1 used to control the strength of the perturbation term. When \(\lambda = 0\), the system is in the reference state Hamiltonian \(H_0\); when \(\lambda = 1\), the system is in the target state Hamiltonian \(H_1\).

[0049] The present invention simulates and statistically calculates the two-way (forward and reverse) free energy change process, and uses the BAR (Bennett Acceptance Ratio) method to estimate the free energy of the data obtained from the two-way sampling simulation respectively, which is expressed by the formula as follows:

[0050] Forward: \(\exp(-\beta\Delta A)=\langle\exp(-\beta\Delta U)\rangle_0\);

[0051] Reverse: \(\exp(+\beta\Delta A)=\langle\exp(+\beta\Delta U)\rangle_1\);

[0052] where \(\beta = -1 / (k\) B T), \(k\) B is the Boltzmann constant, \(T\) is the temperature, and \(\Delta U\) is the potential energy difference of the system in the layering process.

[0053] The free energy difference is determined by comparing the energy difference between the reference state and the target state during the perturbation process. In addition, the results of the forward and backward simulation changes are combined, and the solution to the two-way sampling problem can be found by considering the intermediate state between the reference state and the target state, which is expressed by the formula as follows:

[0054]

[0055] where \(Q_1\) is the energy value of the target state; \(Q_0\) is the energy value of the reference state; \(Q\) i is the energy value of the system when the perturbation parameter is \(i\); \(U\) i is the potential energy of the system when the perturbation parameter is \(i\); \(U_0\) is the potential energy of the system in the reference state; \(U_1\) is the potential energy of the system in the target state.

[0056] By changing the perturbation process to simulate and calculate the accuracy of the analysis results of the number of layers, gradually increase the number of layers, such as 4, 8, 16, 32 layers, etc. After each increase, perform a complete FEP simulation and record the free energy change. Evaluate the optimal solution of the selected number of layers through the calculation results of the BAR method. Increasing the number of layers will cause the final free energy results to tend to be stable, and the error gradually decreases, improving the accuracy of the calculation results. In addition, as the number of layers increases, the simulation calculation time will increase significantly. Therefore, it is necessary to find a balance point in the layering strategy to reduce the calculation time cost while ensuring the calculation accuracy. After the free energy perturbation molecular dynamics simulation of the present invention is completed, the free energy change output during the entire simulation time can be obtained. Through data processing with Origin software, qualitative analysis of the free energy obtained after sampling is performed, such as Figure 3 shown. Combining the free energy results processed by the BAR estimator, it is found that when a 16-layer layering strategy is adopted, not only is the system the most stable, but also the computational resources consumed are the least, which can significantly reduce the calculation time cost. Therefore, it is most appropriate to select a 16-layer layering strategy.

[0057] Reasonably selecting the simulation sampling steps is crucial for obtaining a reliable free energy result estimate. Sufficient sampling steps can ensure that all possible states within the system are fully sampled, significantly reducing statistical errors. By studying the free energy change trend under different sampling steps (as Figure 4 shown), it is found that when the number of steps reaches approximately 25,000 steps, the free energy value has tended to be stable and has a very small difference from the experimental value. Therefore, the present invention finally selects 25,000 steps as the standard simulation sampling steps, and the obtained lithium ion hydration free energy is -110.99 kcal / mol, with an error of 0.17 kcal / mol.

[0058] The most important method for applying the free energy perturbation molecular dynamics simulation to the molecular dynamics simulation is the layering strategy. The convergence degree of the finally obtained free energy data is related to the number of layers. During the perturbation process, the reference state and the target state are stratified in a series of n equally wide layers, and independent trajectories generated in both directions are used to evaluate the lithium ion hydration free energy. When the layering strategy is insufficient, it will lead to poor convergence of the calculation results. In order to obtain higher convergence, the present invention performs simulation calculations on the lithium ion hydration free energy simulation using 8, 16, and 32 stratifications to observe the optimal number of layers. The obtained results are as Figure 3 shown. It can be seen from the figure that when n = 16, the coincidence rate of the hydration free energy is the highest, the result is the most stable, and the computational resources consumed at this time are the least. The free energy for forward and backward conversion and the predicted values of the BAR estimator are represented by black, red, and blue solid lines respectively.

[0059] Step 3: Use the optimal sampling steps and the optimal number of layers to perform free energy perturbation molecular dynamics simulation on the hydration process of charged ions, calculate the hydration free energy of charged ions, and perform size correction.

[0060] At the optimal number of sampling steps and the optimal number of layers, the hydration free energy of lithium ions is -111.07 kcal / mol. The variation trend of the lithium ion hydration free energy with λ is as Figure 5 shown. In addition, the above results ignore the influence of the size of the system (cubic water box) on the results, that is, ΔG will vary with the size of the simulation region. Therefore, the results need to be corrected. In the case of Ewald type simulation, a three-dimensional periodic simulation box is used, and the size-related correction formula is +1 / 2ξ Ewald (q1 2 -q0 2 ), where q1 and q0 are the reference state and the target state, that is, +1 and 0, and ξ Ewald is the size correction parameter, and ξ Ewald = -2.837 / L, where L is the cell length. After calculating the size correction, the hydration free energy of lithium ions is -115.86 kcal / mol, and this result is very close to the experimental measurement value.

[0061] The present invention also provides a system for obtaining the hydration free energy of charged ions, including:

[0062] A model construction model for obtaining the physical parameters of charged ions, constructing a physical model of charged ions according to the physical parameters; and constructing a hydration environment of charged ions according to the physical model of charged ions;

[0063] An optimal parameter acquisition module for successively increasing the number of layers in the hydration process based on the hydration environment of charged ions. Each time the number of layers increases, a free energy perturbation molecular dynamics simulation is performed on the hydration process of charged ions. During the free energy perturbation molecular dynamics simulation, bidirectional sampling simulations are performed using different sampling steps, and then the BAR method is used to analyze the results of the bidirectional sampling simulations to determine the optimal sampling steps and the optimal number of layers;

[0064] A hydration free energy calculation module for performing a free energy perturbation molecular dynamics simulation on the hydration process of charged ions using the optimal sampling steps and the optimal number of layers, and calculating the hydration free energy of charged ions.

[0065] The present invention also provides a computer device, including a memory and a processor; the memory stores a computer program, and the processor is used to run the computer program in the memory to execute the method for obtaining the hydration free energy of charged ions.

[0066] The present invention also provides a computer-readable storage medium storing a computer program, and the computer program is suitable for being loaded by a processor to execute the method for obtaining the hydration free energy of charged ions.

[0067] The above-described embodiments are only preferred specific embodiments of the present invention, and the protection scope of the present invention is not limited thereto. Any simple variations or equivalent replacements of technical solutions that can be obviously obtained by those skilled in the art within the technical scope disclosed by the present invention shall fall within the protection scope of the present invention.

Claims

1. A method for obtaining hydration free energy of charged ions, characterized in that: The steps include: Acquire the physical parameters of the charged ions, and construct a physical model of the charged ions according to the physical parameters; and construct a hydration environment of the charged ions according to the physical model of the charged ions; Based on the hydration environment of the charged ions, the number of layers in the hydration process is gradually increased, and each time the number of layers is increased, a free energy perturbation molecular dynamics simulation is performed on the hydration process of the charged ions, and in the free energy perturbation molecular dynamics simulation process, different sampling steps are used to perform bidirectional sampling simulation, and then the BAR method is used to analyze the bidirectional sampling simulation results to determine the optimal sampling steps and the optimal number of layers; The free energy perturbation molecular dynamics simulation of the hydration process of charged ions was carried out using the optimal sampling steps and the optimal number of layers, and the hydration free energy of the charged ions was calculated.

2. The method for obtaining hydration free energy of charged ions according to claim 1, characterized in that: After calculating the hydration free energy of the charged ions, the method further includes correcting the obtained hydration free energy of the charged ions, specifically, using a three-dimensional periodic simulation box to correct the free energy result of the finite-size molecular dynamics simulation box to obtain the corrected hydration free energy of the charged ions.

3. The method for obtaining hydration free energy of charged ions according to claim 1, characterized in that: The free energy perturbation molecular dynamics simulation also includes setting periodic boundary conditions, equilibrium time, sampling step, cutoff radius and initial temperature of the free energy perturbation molecular dynamics simulation, and the periodic boundary is to set boundaries in three directions of X, Y and Z.

4. The method for obtaining hydration free energy of charged ions according to claim 1, characterized in that: The physical parameters of the charged ions include: mass charge, van der Waals force and electrostatic force parameters of the charged ions.

5. The method for obtaining hydration free energy of charged ions according to claim 1, characterized in that: The method of constructing a charged ion physical model according to physical parameters comprises the following steps: Construct PDB files based on the physical parameters of charged ions; Import the PDB file into the Molecule File Browser module in VMD software, and use the AutoPSF function of the Molecule File Browser module to generate a PSF file; Import the PDB file and PSF file into VMD software at the same time to obtain a separate charged ion physical model.

6. The method for obtaining hydration free energy of charged ions according to claim 1, characterized in that: The construction of the hydration environment of charged ions is specifically as follows: in the VMD software, select the TK Console menu item in Extensions, enter package require solvate in the VMD TKConsole window, and wait for the VMD software to load the solvate module; after loading, enter the command solvate Li.psfLi.pdb -t 50 -oLi-solvate to obtain Cube water box.

7. A charged ion hydration free energy acquisition system, characterized in that: include: A model building model is used to obtain the physical parameters of the charged ions and build a physical model of the charged ions according to the physical parameters; Construct the hydration environment of charged ions based on the charged ion physical model; An optimal parameter acquisition module is used to gradually increase the number of layers in the hydration process based on the hydration environment of the charged ions. Each time the number of layers increases, a free energy perturbation molecular dynamics simulation is performed on the hydration process of the charged ions. In the free energy perturbation molecular dynamics simulation process, different sampling steps are used to perform a bidirectional sampling simulation, and then the BAR method is used to analyze the bidirectional sampling simulation results to determine the optimal sampling steps and the optimal number of layers. The hydration free energy calculation module is used to perform free energy perturbation molecular dynamics simulation on the hydration process of charged ions using the optimal sampling step number and the optimal number of layers, and calculate the hydration free energy of charged ions.

8. A computer device, characterized in that: It comprises a memory and a processor; the memory stores a computer program, and the processor is used to run the computer program in the memory to execute the method for obtaining the hydration free energy of charged ions as described in any one of claims 1-6.

9. A computer-readable storage medium, characterized in that: The computer-readable storage medium stores a computer program, and the computer program is suitable for being loaded by a processor to execute the method for obtaining hydration free energy of charged ions according to any one of claims 1 to 6.

Images