Task control method and system for calculating FEP parameters in calculation of predicted binding free energy

Abstract

The application relates to the technical field of drug design, and discloses a task control method and system for calculating FEP parameters in calculation of predicted binding free energy calculation, which comprises the following specific steps: according to a topology and a coordinate file of an input system, selecting a molecular dynamics engine, and constructing initial FEP parameters; performing molecular dynamics simulation, calculating potential energy difference between a current window and a next window; calculating actual dissipation work of the potential energy difference and variance of free energy of the window; comparing the actual dissipation work and the variance with a maximum theoretical dissipation work and an error threshold value derived from an error estimation formula, adjusting a current window scheme or predicting a step length and time of next simulation according to a comparison result, and performing next window simulation; and running an analysis program to calculate binding free energy and using a convergence judgment tool to perform statistical analysis. The application solves the problem that the existing technology needs a long time to calculate binding free energy, and has the characteristics that a suitable simulation scheme can be formulated according to the disturbance size of different systems.

Description

Technical Field

[0001] This invention relates to the field of drug design technology, and more specifically, to a task control method and system for combining computational prediction with FEP parameter calculation in free energy calculation. Background Technology

[0002] The key to a drug's efficacy lies in its binding affinity to its target. A precise and efficient method for predicting drug-target affinity can greatly accelerate the drug development process, thereby significantly reducing the cycle and cost of new drug development.

[0003] To date, numerous techniques have been developed for calculating drug-target binding free energy, such as scoring functions built into docking procedures, and methods based on implicit solvent models like MM / PBSA (Molecular Mechanics / Poisson-Boltzmann Surface) and MM / GBSA (Molecular Mechanics / Generalized Berne Surface). However, because these methods typically rely on experimental data and require the introduction of empirical constants, the calculated results often differ significantly from actual experimental data. Methods based on free energy perturbations for calculating binding free energy have a rigorous theoretical foundation, do not require empirical parameters, and can achieve accurate predictions of drug-target affinity. However, this method lacks an optimal simulation scheme. Using a fixed simulation window and simulation duration may be insufficient to achieve convergence for systems with large perturbations, while for systems with small perturbations, some computational effort may be wasted.

[0004] The error in free energy calculation can be used as a criterion for judging whether the free energy calculation has converged. The error in free energy calculation is divided into two parts: systematic error (bias) and statistical error (variance). When the sampling size is large, the error mainly comes from the statistical error. Different perturbation systems have different perturbation magnitudes. By controlling the statistical error and developing a task control method that combines computational prediction with the FEP parameter in free energy calculation, reasonable molecular dynamics simulation schemes can be formulated for different perturbation systems, solving the technical shortcomings of long calculation time and low computational efficiency.

[0005] An existing technology provides a perturbation method for absolute free energy to predict drug-target binding strength, comprising: constructing a Lig system and a Rec-Lig system; performing molecular dynamics simulations on the parameters of each state in the Lig and Rec-Lig systems to obtain the molecular dynamic trajectory of each state; calculating the potential energy difference between each adjacent state by performing single-point energy calculations using force field files of adjacent states; statistically analyzing the probability distribution of the potential energy difference between each adjacent state and fitting it using a linear combination of multiple Gaussian functions; regenerating the potential energy difference of each state based on the fitted probability distribution, and calculating the binding free energy between the drug and the target.

[0006] However, existing technologies have the problem of long calculation time required to calculate the combined free energy. Therefore, how to invent a time-saving task control method and system for calculating and predicting FEP parameters in the combined free energy calculation is a technical problem that urgently needs to be solved in this field. Summary of the Invention

[0007] To address the problem of long calculation times required for the combined free energy in existing technologies, this invention provides a task control method and system for calculating and predicting FEP parameters in the combined free energy calculation. It has the characteristic of being able to formulate appropriate simulation schemes according to the magnitude of disturbances in different systems.

[0008] To achieve the above-mentioned objectives of this invention, the technical solution adopted is as follows:

[0009] The task control method that combines FEP parameters in free energy calculation with prediction includes the following specific steps:

[0010] Step S1: Based on the input system topology and coordinate file, select the molecular dynamics engine and set the initial file and initial FEP parameters required for molecular dynamics simulation;

[0011] Step S2: Perform molecular dynamics simulation based on the initial FEP parameters; sample the current window using molecular dynamics simulation, and calculate the potential energy difference ΔU between the current window and the next window at an interval of Δλ using the force field files of adjacent windows; calculate the actual dissipation work W of the potential energy difference ΔU. ΔU The variance of the free energy of that window

[0012] Step S3: Derive the theoretical maximum value of dissipated work W using the error estimation formula. Δλ ; Calculate the error threshold The maximum theoretical dissipation work W Δλ and error threshold With actual dissipated work W ΔU and variance The comparison is performed, and the prediction combined with the FEP parameter in the free energy calculation is adjusted according to the comparison results. The task control method adjusts the current window scheme or predicts the step size and time of the next simulation to perform the next window simulation, until all window simulations are completed.

[0013] Step S4: Run the analysis program to calculate the binding free energy and perform statistical analysis using the convergence judgment tool.

[0014] Step S1 specifically includes the following steps:

[0015] Based on the input system's topology and coordinate files, select the molecular dynamics engine and set the initial files and initial FEP parameters required for the molecular dynamics simulation. Specifically, this includes setting the initial window λ, initial step size Δλ, the next adjacent window λ′ = λ + Δλ, the initial simulation time, and the overall error threshold Var for a thermodynamic process.

[0016] Furthermore, step S2 specifically includes the following steps:

[0017] Step S2.1: Perform molecular dynamics simulation on the current window λ and sample the molecular dynamics trajectory of the current window;

[0018] Step S2.2: Based on the molecular dynamics trajectory of the current window, calculate the potential energy difference ΔU between adjacent windows with a spacing of Δλ using the force field files of adjacent windows, and calculate the actual dissipation work W of the potential energy difference ΔU. ΔU Specifically:

[0019]

[0020] Step S2.3: Calculate the variance of the free energy during the window. Specifically: Set the forward step size and the length of the sliding window, move the window from the first frame to the last frame, and use the FEP formula to calculate ΔU in each sliding window segment, and calculate ΔG accordingly.

[0021] ΔG=-k B Tln <e -β(ΔU) >

[0022] Where, k B Let be the Boltzmann constant, and T be the ambient temperature; calculate the variance of a series of ΔG as...

[0023]

[0024] Furthermore, step S3 specifically includes the following steps:

[0025] Step S3.1: Derive the theoretical maximum dissipation work W using the error estimation formula. Δλ The specific steps are as follows:

[0026] Introducing the error estimation formula of Jarzyski's equation:

[0027]

[0028] Where N is the number of sampled trajectories, W dis Dissipated work is the work done during state transitions.

[0029] By introducing parameter 'a' to correct N in the error estimation formula of Jarzynski's equation, we obtain the error formula for FEP calculation:

[0030]

[0031] Among them, e a The physical meaning of a is the interval between independent samples. a is independent of N but related to λ. The parameter a is determined by fitting the actual error obtained from the simulation.

[0032] Let the total error of a thermodynamic process calculation be Var. Using the error formula for FEP calculation, we can calculate the work dissipated by energy sampled in the current window with a step size of Δλ under a specified total error Var. Δλ Maximum value:

[0033]

[0034] Step S3.2: Calculate the error threshold of ΔG calculated in the current window with a step size of Δλ under the specified total error Var. Specifically:

[0035]

[0036] Step S3.3: Calculate the theoretical maximum dissipation work W. Δλ and error threshold With actual dissipated work W ΔU and variance The comparison is performed, and the current window scheme is adjusted or the step size and time of the next simulation are predicted based on the comparison results, and then the next window simulation is performed.

[0037] Specifically, step S3.3 includes the following steps:

[0038] like And W ΔU ≤W Δλ If the current window energy with a step size of Δλ has converged, then the Δλ′ of the next simulation can be appropriately increased. Specifically, the value of the next step size Δλ′ can be predicted using a formula.

[0039] like or W ΔU >W Δλ If Δλ and Δλ′ are reduced, the potential energy difference ΔU and the variance of ΔU are recalculated from the saved potential energy file of adjacent windows using the reduced Δλ′. With dissipation work W ΔU ; Re-apply the theoretical maximum value of dissipated work W Δλ Combined with the set error threshold With actual dissipated work W ΔU and variance Compare;

[0040] If the current Δλ is reduced to 0.005, the problem still exists. or W ΔU >W Δλ Then, based on the existing trajectory, the sampling time is extended, and the theoretical maximum dissipation work W is recalculated. Δλ Combined with the set error threshold With actual dissipated work W ΔU and variance Compare them.

[0041] Furthermore, the specific steps for predicting the value of the next step size Δλ′ using the formula are as follows:

[0042] Derivation of the actual dissipation work W ΔU Relationship with λ:

[0043]

[0044] The relationship between the potential energy difference and the total potential energy difference is obtained:

[0045]

[0046] The relationship between the potential energy difference and Δλ between two specific states λ1 and λ2 is obtained:

[0047]

[0048] Simultaneously consider:

[0049]

[0050] Then ΔU and Δλ, With Δλ 2 Proportional; the dissipation work W under the next window λ′ with step size Δλ′ is obtained. dis 'for:

[0051]

[0052] To obtain the value of Δλ′ in the next step:

[0053]

[0054] Furthermore, step S4 specifically includes the following steps: running the analysis program to calculate the binding free energy, and performing statistical analysis using a convergence judgment tool, specifically:

[0055] Run the analysis program to calculate the binding free energy using the Bennet Acceptance Ratio algorithm; use convergence judgment tools to perform analysis, construct time series convergence plots, ΔG / Δλ plots, and use the reweighting algorithm to generate a heatmap of the difference in ΔG between adjacent windows.

[0056] The task control system, which combines computational prediction with FEP parameter calculation in free energy calculation, includes a cascaded task control module, dynamic simulation module, FEP parameter prediction module, and statistical analysis module.

[0057] The task control module is used to control all steps of the task control method for combining the computational prediction with the FEP parameters in the free energy calculation; including constructing the initial files and initial FEP parameters required for molecular dynamics simulation; and controlling the execution of other modules as a workflow.

[0058] The aforementioned dynamics simulation module is used to perform molecular dynamics simulations based on FEP parameters; it samples the current window for molecular dynamics simulation and calculates the potential energy difference ΔU between the current window and the next window at an interval of Δλ using the force field files of adjacent windows; it then calculates the actual dissipation work W of the potential energy difference ΔU. ΔU The variance of the free energy of that window

[0059] The FEP parameter prediction module is used to derive the theoretical maximum dissipation work W using the error estimation formula. Δλ ; Calculate the error threshold The maximum theoretical dissipation work W Δλ and error threshold With actual dissipated work W ΔU and variance The comparison is performed, and the current window scheme is adjusted or the step size and time of the next simulation are predicted based on the comparison results, and then the next window simulation is performed.

[0060] The statistical analysis module is used to run analysis programs to calculate the binding free energy and to perform statistical analysis using convergence judgment tools.

[0061] The beneficial effects of this invention are as follows:

[0062] This invention discloses a task control method for calculating and predicting FEP parameters in free energy calculation. It calculates the potential energy difference ΔU between the current window and the next window at an interval of Δλ from the force field files of adjacent windows by sampling the current window using molecular dynamics simulation; and calculates the actual dissipation work W of ΔU. ΔU The variance of the free energy of that window The maximum theoretical dissipation work W Δλ and error threshold With actual dissipated work WΔU and variance The simulation scheme is compared, the current window scheme is adjusted according to the comparison results, and the next window simulation is performed until all window simulations are completed. Statistical analysis is then performed. Thus, this invention formulates appropriate simulation schemes based on the magnitude of perturbations in different systems, enabling the system energy to converge in the shortest simulation time. This significantly reduces the time required for calculating the binding free energy, providing guidance for the optimization and modification of small molecule drug structures. Attached Figure Description

[0063] Figure 1 This is a flowchart illustrating the task control method for calculating and predicting FEP parameters in the free energy calculation of this invention.

[0064] Figure 2 The embodiment is a flowchart of the task control method for calculating and predicting FEP parameters in the free energy calculation of the present invention.

[0065] Figure 3 This is a convergence plot of ΔΔG over time for the perturbation 23-46 of the tyk2 system in the task control method of the present invention, which combines the calculation and prediction of FEP parameters in the calculation of free energy.

[0066] Figure 4 This is the ΔΔG / Δλ diagram of the perturbation 23-46 of the tyk2 system in the task control method of the present invention, which combines the calculation and prediction of FEP parameters in the calculation of free energy.

[0067] Figure 5 This is a heatmap of the difference between ΔΔG between adjacent windows generated by the reweighting algorithm in the task control method of the present invention, which calculates and predicts FEP parameters in the free energy calculation.

[0068] Figure 6 This is a simulation time distribution diagram of eight protein systems in the task control method of the present invention, which combines the calculation and prediction of FEP parameters in the calculation of free energy. Detailed Implementation

[0069] The present invention will now be described in detail with reference to the accompanying drawings and specific embodiments.

[0070] Example 1

[0071] like Figure 1 As shown, the task control method that combines FEP parameter calculation with free energy calculation includes the following specific steps:

[0072] Step S1: Based on the input system topology and coordinate file, select the molecular dynamics engine and set the initial file and initial FEP parameters required for molecular dynamics simulation;

[0073] Step S2: Perform molecular dynamics simulation based on the initial FEP parameters; sample the current window using molecular dynamics simulation, and calculate the potential energy difference ΔU between the current window and the next window at an interval of Δλ using the force field files of adjacent windows; calculate the actual dissipation work W of the potential energy difference ΔU. ΔU The variance of the free energy of that window

[0074] Step S3: Derive the theoretical maximum value of dissipated work W using the error estimation formula. Δλ ; Calculate the error threshold The maximum theoretical dissipation work W Δλ and error threshold With actual dissipated work W ΔU and variance The comparison is performed, and the prediction combined with the FEP parameter in the free energy calculation is adjusted according to the comparison results. The task control method adjusts the current window scheme or predicts the step size and time of the next simulation to perform the next window simulation, until all window simulations are completed.

[0075] Step S4: Run the analysis program to calculate the binding free energy and perform statistical analysis using the convergence judgment tool.

[0076] In this embodiment, as Figure 2 As shown, selectable molecular dynamics engines include, but are not limited to, Amber, Gromacs, and Open-MM.

[0077] Example 2

[0078] In this embodiment, the relative binding free energy was calculated using a test set of 8 targets and 190 ligands;

[0079] In one specific embodiment, the initial files and initial FEP parameters required for molecular dynamics simulation include: setting the initial window λ, the initial step size Δλ, the next adjacent window λ′ = λ + Δλ, the initial simulation time, and the overall error threshold Var for a thermodynamic process.

[0080] In this embodiment, the input files are the topology and coordinate files of the receptor protein-ligand system and the individual small molecule ligand system. The molecular dynamics simulation engine selected is AMBER, the force field is AMBER FF14SB / GAFF1.8, and a thermodynamic cycle for calculating the relative binding free energy is constructed. The initial molecular dynamics simulation parameters include an initial window λ = 0.0, an initial step size Δλ = 0.01, the next adjacent window λ′ = 0.01, an initial simulation time t = 10 ps, ​​and a global standard deviation threshold Var = 0.25 kcal / mol for a thermodynamic process.

[0081] In one specific embodiment, a molecular dynamics simulation is performed on the current window λ to obtain the molecular dynamic trajectory of the current window. Based on the molecular dynamic trajectory of the current window, the potential energy difference ΔU between adjacent windows with a spacing of Δλ is calculated using the force field files of adjacent windows. Simultaneously, the potential energy of each window with a spacing of Δλ = 0.005 between two adjacent windows is also saved for subsequent adjustments to the current simulation window, avoiding resimulation. The actual dissipation work W of the potential energy difference ΔU is calculated. ΔU Specifically:

[0082]

[0083] In one specific embodiment, the variance of the free energy of the window is calculated. Specifically: Set the forward step size and the length of the sliding window, move the window from the first frame to the last frame, and use the FEP formula to calculate ΔU in each sliding window segment, and calculate ΔG accordingly.

[0084] ΔG=-k B Tln <e -β(ΔU) >

[0085] Where, k B Let be the Boltzmann constant, and T be the ambient temperature; calculate the variance of a series of ΔG as...

[0086]

[0087] In one specific embodiment, the theoretical maximum dissipation work W is derived using an error estimation formula. Δλ The specific steps are as follows:

[0088] Introducing the error estimation formula of Jarzyski's equation:

[0089]

[0090] Where N is the number of sampled trajectories, W dis Dissipated work is the work done during state transitions.

[0091] By introducing parameter 'a' to correct N in the error estimation formula of Jarzynski's equation, we obtain the error formula for FEP calculation:

[0092]

[0093] Among them, e a The physical meaning of a is the interval between independent samples. a is independent of N but related to λ. The parameter a is determined by fitting the actual error obtained from the simulation.

[0094] Let the total error of a thermodynamic process calculation be Var. Using the error formula for FEP calculation, we can calculate the work dissipated by energy sampled in the current window with a step size of Δλ under a specified total error Var. Δλ Maximum value:

[0095]

[0096] In one specific embodiment, the error threshold is calculated. Specifically:

[0097] Calculate the error threshold of ΔG calculated in the current window with a step size of Δλ under a specified total error Var.

[0098]

[0099] In one specific embodiment, the theoretical maximum dissipation work W Δλ and error threshold With actual dissipated work W ΔU and variance The comparison is performed, and the current window scheme is adjusted based on the comparison results. The specific steps are as follows:

[0100] like And W ΔU ≤W Δλ If the current window energy with a step size of Δλ has converged, then the Δλ′ of the next simulation can be appropriately increased. Specifically, the value of the next step size Δλ′ can be predicted using a formula.

[0101] like or W ΔU >W Δλ If Δλ and Δλ′ are reduced, the potential energy difference ΔU and the variance of ΔU are recalculated from the saved potential energy file of adjacent windows using the reduced Δλ′. With dissipation work W ΔU ; Re-apply the theoretical maximum value of dissipated work W Δλ Combined with the set error threshold With actual dissipated work W ΔU and variance Compare;

[0102] If the current Δλ is reduced to 0.005, the problem still exists. or W ΔU >W Δλ Then, based on the existing trajectory, the sampling time is extended by 50 ps each time, and the theoretical maximum dissipation work W is recalculated. Δλ Combined with the set error threshold With actual dissipated work W ΔU and variance Compare them.

[0103] In one specific embodiment, the value of the next step size Δλ′ is predicted using a formula, and the specific steps are as follows:

[0104] Derivation of the actual dissipation work W ΔU Relationship with λ:

[0105]

[0106] The relationship between the potential energy difference and the total potential energy difference is obtained:

[0107]

[0108] The relationship between the potential energy difference and Δλ between two specific states λ1 and λ2 is obtained:

[0109]

[0110] Simultaneously consider:

[0111]

[0112] Then ΔU and Δλ, With Δλ 2 Proportional; the dissipation work W under the next window λ′ with step size Δλ′ is obtained. dis 'for:

[0113]

[0114] To find the value of Δλ′ in the next step:

[0115]

[0116] In one specific embodiment, if Δλ′ is greater than 0.05, it is set to 0.05 to prevent the calculation from being difficult to converge due to an excessively large step size.

[0117] In one specific embodiment, the analysis program is run to calculate the binding free energy, and statistical analysis is performed using a convergence assessment tool, specifically:

[0118] Once the λ=1.0 window is completed, the molecular dynamics simulation program is stopped, and the binding free energy is calculated using the Benneet Acceptance Ratio algorithm. Convergence assessment tools are used for analysis, constructing a time-series convergence plot, a ΔG / Δλ plot, and a heatmap of the difference in ΔG between adjacent windows using a reweighting algorithm. The time-series convergence plot is used to determine the energy convergence of each simulation window; the ΔG / Δλ plot is used to observe the repeatability of each calculation under repeated computations; and the heatmap of the difference in ΔG between adjacent windows is used to determine if any simulation anomalies occur in a particular window.

[0119] In this embodiment, RBFE calculations were performed five times in parallel for each pair of ligands. Outliers were detected and removed using the K-Means clustering algorithm. The average of the removed data was used to obtain the ΔG value, and a cyclic closure convergence tool was employed to obtain the final ΔG value. The comparison of the obtained test results with FEP+ and AMBER-TI is shown in Table 1.

[0120]

[0121] Table 1 shows the mean absolute error (MAE), root mean square error (RMSE), and correlation coefficient R for predicting ΔG of eight protein systems using FEP+, AMBER-TI, and the method described in this patent (autoFEP). 2 In summary with Kendall, Table 1 shows that, compared to FEP+, although enhanced sampling was not used and the simulation time was very short, the results are close to those of FEP+. Compared to AmberTI's results, all indicators are better. This demonstrates that the method of this invention has high accuracy.

[0122] A perturbation (23-46) in the Tyk2 system is selected, and statistical analysis is performed using convergence testing tools to illustrate its effectiveness. The convergence plot of ΔG over time series is shown below. Figure 3 As shown, the energy curves of the red forward time series estimation and the blue reverse time series estimation converge rapidly to a certain value from both ends, while the green ΔG fluctuation curve remains stable over time.

[0123] ΔΔG / Δλ diagram as shown Figure 4 As shown, the trends of ΔΔG / Δλ in the five parallel calculations are consistent and the fluctuations are small.

[0124] The heatmap of the difference ΔG between adjacent windows generated using the reweighting algorithm is shown below. Figure 5 As shown, each small square is a lighter color, indicating that the energy between adjacent windows converges and is self-consistent.

[0125] In this embodiment, the simulation duration of five parallel tests for each perturbation is directly averaged, and the distribution of the simulation duration of all perturbations is as follows: Figure 6 As shown, 95% of the perturbations are completed within 20 ns, and 74% of the perturbations are completed within 10 ns. autoFEP can rationally allocate simulation time according to the magnitude of perturbations in different systems. Compared with FEP+ and AmberTI, which have simulation times exceeding 70 ns per system, autoFEP reduces the simulation time by several to tens of times, significantly reducing the required computational resources.

[0126] Therefore, using a test set of 8 targets and 190 ligands as examples, the proposed method was used to simulate the system. Under conditions of no enhanced sampling and very short simulation time, the mean absolute error (MAE), root mean square error (RMSE), correlation coefficient R², and Kendall's τ of the eight protein systems obtained were close to the FEP+ method and better than the AmberTI method, while the simulation time was reduced by several to tens of times compared to FEP+ and AmberTI. This indicates that the task control method for predicting binding free energy calculations using FEP parameters described in this invention can accurately and quickly calculate binding free energy based on different system perturbation magnitudes, greatly reducing the computational load. It holds promise for application in drug design to guide the design and modification of lead compounds.

[0127] Example 3

[0128] The task control system, which combines computational prediction with FEP parameter calculation in free energy calculation, includes a cascaded task control module, dynamic simulation module, FEP parameter prediction module, and statistical analysis module.

[0129] The task control module is used to control all steps of the task control method for combining the computational prediction with the FEP parameters in the free energy calculation; including constructing the initial files and initial FEP parameters required for molecular dynamics simulation; and controlling the execution of other modules as a workflow.

[0130] The aforementioned dynamics simulation module is used to perform molecular dynamics simulations based on FEP parameters; it samples the current window for molecular dynamics simulation and calculates the potential energy difference ΔU between the current window and the next window at an interval of Δλ using the force field files of adjacent windows; it then calculates the actual dissipation work W of the potential energy difference ΔU. ΔU The variance of the free energy of that window

[0131] The FEP parameter prediction module is used to derive the theoretical maximum dissipation work W using the error estimation formula. Δλ ; Calculate the error threshold The maximum theoretical dissipation work W Δλ and error threshold With actual dissipation work W ΔU and variance The comparison is performed, and the current window scheme is adjusted or the step size and time of the next simulation are predicted based on the comparison results, and then the next window simulation is performed.

[0132] The statistical analysis module is used to run analysis programs to calculate the binding free energy and to perform statistical analysis using convergence judgment tools.

[0133] Obviously, the above embodiments of the present invention are merely examples for clearly illustrating the present invention, and are not intended to limit the implementation of the present invention. Any modifications, equivalent substitutions, and improvements made within the spirit and principles of the present invention should be included within the protection scope of the claims of the present invention.

Claims

1. A task control method that combines calculation and prediction with FEP parameter calculation in free energy calculation, characterized by: The specific steps include the following: Step S1: Based on the input system topology and coordinate file, select the molecular dynamics engine and set the initial file and initial FEP parameters required for molecular dynamics simulation; Step S2: Perform molecular dynamics simulation based on the initial FEP parameters; sample the current window for molecular dynamics simulation, and calculate the force field relationship between the current window and the interval Δ using the force field files of adjacent windows. λ The potential energy difference ΔU between the next windows; calculate the actual dissipation work of the potential energy difference ΔU. The variance of the free energy of that window ; Step S3: Derive the theoretical maximum value of dissipated work using the error estimation formula. ; Calculate the error threshold The maximum theoretical dissipation work and error threshold With actual dissipation work and variance The comparison is performed, and the prediction combined with the FEP parameters in the free energy calculation is adjusted according to the comparison results to adjust the current window scheme or the step size and time of the next simulation, and the next window simulation is performed until all window simulations are completed. Step S3 specifically includes the following steps: Step S3.1: Derive the theoretical maximum value of dissipated work using the error estimation formula. The specific steps are as follows: Introducing the error estimation formula of Jarzyski's equation: Where N is the number of sampled trajectories, Dissipated work is the work done during state transitions. By introducing parameter 'a' to correct N in the error estimation formula of Jarzynski's equation, we obtain the error formula for FEP calculation: in, The physical meaning is the interval between independent samples, a is independent of N, and is related to The parameter 'a' is determined by fitting the actual error obtained from the simulation. Let the total error of a thermodynamic process calculation be Var. Using the error formula from FEP calculations, we can calculate the step size under a specified total error Var. The dissipation work of energy obtained from the current window sampling Maximum value: Step S3.2: Calculate the next step size of the specified total error Var. The current window is calculated Error threshold Specifically: ； Step S3.3: Calculate the theoretical maximum dissipation work. and error threshold With actual dissipation work and variance The comparison is performed, and the current window scheme is adjusted or the step size and time of the next simulation are predicted based on the comparison results, and then the next window simulation is performed. Specifically, step S3.3 includes the following steps: like and This indicates that... To ensure the energy convergence of the current window size, appropriately increase the step size for the next simulation. Specifically, it involves using a formula to determine the next step size. Values ​​are predicted; like or Then decrease and From the saved adjacent window potential energy file, the reduced Recalculate the potential energy difference ΔU and the variance of ΔU. With dissipation power ; Re-apply the maximum theoretical dissipation work Combined with the set error threshold With actual dissipation work and variance Compare; If the current It still exists even after being reduced to 0.

005. or Then, based on the existing trajectory, the sampling time is extended, and the theoretical maximum dissipation work is recalculated. Combined with the set error threshold With actual dissipation work and variance Compare; Step S4: Run the analysis program to calculate the binding free energy and perform statistical analysis using the convergence judgment tool.

2. The task control method for predicting and combining FEP parameters in free energy calculation according to claim 1, characterized in that: Step S1 specifically includes the following steps: Based on the input system's topology and coordinate files, select the molecular dynamics engine and set the initial files and initial FEP parameters required for the molecular dynamics simulation. Specifically, this includes setting the initial window. λ Initial step size Δ λ The next adjacent window = λ + Δ λ Initial simulation time, and the overall error threshold Var for a thermodynamic process.

3. The task control method for predicting and combining FEP parameters in free energy calculation according to claim 2, characterized in that: Step S2 specifically includes the following steps: Step S2.1: Perform molecular dynamics simulation on the current window λ and sample the molecular dynamics trajectory of the current window; Step S2.2: Based on the molecular dynamics trajectory of the current window, calculate the potential energy difference ΔU between adjacent windows with a spacing of Δλ using the force field files of adjacent windows, and calculate the actual dissipation work of the potential energy difference ΔU. Specifically: Step S2.3: Calculate the variance of the free energy during the window. Specifically, this involves setting the forward step size and the length of the sliding window, moving the window from the first frame to the last frame, and using the FEP formula to calculate ΔU for each segment of the sliding window. : in, Boltzmann's constant, For ambient temperature; calculate a series of variance as : 。 4. The task control method for predicting and combining FEP parameters in free energy calculation according to claim 1, characterized in that: The formula is used to determine the next step size. The specific steps for predicting the value are as follows: Derivation of actual dissipation work Relationship with λ: The relationship between the potential energy difference and the total potential energy difference is obtained: Two specific states are obtained and The potential energy difference between them and Relationship: Simultaneously consider: but and , and Proportional; obtained at step size The next window The dissipation of the work for: Seeking the next step Value: 。 5. The task control method for predicting and combining FEP parameters in free energy calculation according to claim 1, characterized in that: Step S4 specifically includes the following steps: running the analysis program to calculate the binding free energy, and performing statistical analysis using a convergence judgment tool, specifically: Run the analysis program to calculate the binding free energy using the Bennet Acceptance Ratio algorithm; perform analysis using convergence assessment tools, including constructing a time series convergence plot. Figure 1. Using the Reweighting algorithm to generate the spacing between adjacent windows. The difference in heatmaps.

6. A task control system that combines computational prediction with FEP parameter calculation in free energy calculation, characterized in that: It includes a cascaded task control module, a dynamics simulation module, an FEP parameter prediction module, and a statistical analysis module; The task control module is used to control all steps of the task control method for combining the computational prediction with the FEP parameters in the free energy calculation; including constructing the initial files and initial FEP parameters required for molecular dynamics simulation; and controlling the execution of other modules as a workflow. The aforementioned dynamics simulation module is used to perform molecular dynamics simulations based on FEP parameters; it samples the current window for molecular dynamics simulations and calculates the relationship between the current window and the interval Δ using the force field files of adjacent windows. λ The potential energy difference ΔU between the next windows; calculate the actual dissipation work of the potential energy difference ΔU. The variance of the free energy of that window ; The FEP parameter prediction module is used to derive the theoretical maximum dissipation work value through the error estimation formula. ; Calculate the error threshold ; The maximum theoretical dissipation work and error threshold With actual dissipation work and variance The comparison is performed, and the current window scheme is adjusted or the step size and time of the next simulation are predicted based on the comparison results, and then the next window simulation is performed. The FEP parameter prediction module includes the following steps: The maximum theoretical dissipation work is derived using the error estimation formula. The specific steps are as follows: Introducing the error estimation formula of Jarzyski's equation: Where N is the number of sampled trajectories, Dissipated work is the work done during state transitions. By introducing parameter 'a' to correct N in the error estimation formula of Jarzynski's equation, we obtain the error formula for FEP calculation: in, The physical meaning is the interval between independent samples, a is independent of N, and is related to The parameter 'a' is determined by fitting the actual error obtained from the simulation. Let the total error of a thermodynamic process calculation be Var. Using the error formula from FEP calculations, we can calculate the step size under a specified total error Var. The dissipation work of energy obtained from the current window sampling Maximum value: Calculate the next step size of the specified total error Var. The current window is calculated Error threshold Specifically: ； The maximum theoretical dissipation work and error threshold With actual dissipation work and variance The comparison is performed, and the current window scheme is adjusted or the step size and time of the next simulation are predicted based on the comparison results before proceeding to the next window simulation. Specifically, this includes the following steps: like and This indicates that... To ensure the energy convergence of the current window size, appropriately increase the step size for the next simulation. Specifically, it involves using a formula to determine the next step size. Values ​​are predicted; like or Then decrease and From the saved adjacent window potential energy file, the reduced Recalculate the potential energy difference ΔU and the variance of ΔU. With dissipation power ; Re-apply the maximum theoretical dissipation work Combined with the set error threshold With actual dissipation work and variance Compare; If the current It still exists even after being reduced to 0.

005. or Then, based on the existing trajectory, the sampling time is extended, and the theoretical maximum dissipation work is recalculated. Combined with the set error threshold With actual dissipation work and variance Compare; Step S4: Run the analysis program to calculate the binding free energy and perform statistical analysis using the convergence judgment tool; The statistical analysis module is used to run analysis programs to calculate the binding free energy and to perform statistical analysis using convergence judgment tools.

Images