Ensemble molecular dynamics for predicting binding energies of DNA intercalators

Ensemble molecular dynamics (MD) simulations have been widely used to improve binding energy predictions in protein-ligand systems, but their application to DNA-ligand systems is limited. Following previous studies on DNA-doxorubicin and DNA-proflavine systems, the present study extends the ensemble MD simulations for a total of 16 intercalators spanning different molecular sizes, intercalating scaffolds, flexibilities, and charges. For each DNA-intercalator complex, 25 independent 10 ns MD replicas were simulated using randomized initial velocities. The uncorrected Molecular Mechanics/Poisson-Boltzmann Surface Area (MM/PBSA) energies significantly overestimated the binding energies and had a poor correlation with the corresponding experimental values (with a correlation coefficient R ² of 0.45). Including entropy and deformation energy correction terms to the MM/PBSA energies reduced the discrepancy between experimental and predicted binding energies of the cationic intercalators to within ca. 2 kcal mol^-1 (with an R ² of 0.70). Two additional monocationic intercalators were included for validation. Their predicted binding energies were within 1.2 kcal mol^-1 of the corresponding experimental values. Bootstrap analyses showed that around 10 independent replicas are generally sufficient to reduce the uncertainty in binding energy predictions for cationic intercalators to below 1 kcal mol^-1, although up to 20 replicas could be needed for more flexible ligands that have relatively high degrees of freedom. Overall, MM/PBSA energies estimated from an ensemble of short MD simulations, and corrected for entropy and deformation, are reliable and computationally efficient for ranking and predicting the binding affinity of the cationic intercalators in this study. This method provides a valuable alternative when experimental binding measurements are impractical or costly.