RESEARCH ARTICLE Protein surface roughness accounts for binding free energy of Plasmepsin II‐ligand complexes Mario E. Valdés‐Tresanco 1 | Mario S. Valdés‐Tresanco 2 | Pedro A. Valiente 1 | Germinal Cocho 3 | Ricardo Mansilla 4 | J.M. Nieto‐Villar 5 1 Computational Biology and Biomolecular Dynamics Laboratory, Center for Proteins Studies, Faculty of Biology, University of Havana, Havana, Cuba 2 Faculty of Biology, University of Havana, Havana, Cuba 3 C3 Complex Systems Institute and UNAM Physics Institute, Mexico 4 Center for Interdisciplinary Investigations of Humanities and Sciences, UNAM, Mexico 5 Department of Chemical‐Physics, Faculty of Chemistry and H. Poincare Group of Complex Systems, Faculty of Physics, University of Havana, Havana, Cuba Correspondence Mario E. Valdés‐Tresanco, Computational Biology and Biomolecular Dynamics Laboratory, Center for Proteins Studies, Faculty of Biology, University of Havana, Havana 10400, Cuba. Email: mevaldes@fbio.uh.cu J. M. Nieto‐Villar, Department of Chemical‐ Physics, Faculty of Chemistry and H. Poincare Group of Complex Systems, Faculty of Physics, University of Havana, Havana 10400, Cuba. Email: nieto@fq.uh.cu Abstract The calculation of absolute binding affinities for protein‐inhibitor complexes remains as one of the main challenges in computational structure‐based ligand design. The present work explored the calculations of surface fractal dimension (as a measure of surface roughness) and the relationship with experimental binding free energies of Plasmepsin II complexes. Plasmepsin II is an attractive target for novel therapeutic compounds to treat malaria. However, the structural flexibility of this enzyme is a drawback when searching for specific inhibitors. Concerning that, we performed separate explicitly solvated molecular dynamics simulations using the available high‐resolution crystal structures of different Plasmepsin II complexes. Molecular dynamics simulations allowed a better approximation to systems dynamics and, therefore, a more reliable estimation of surface roughness. This constitutes a novel approximation in order to obtain more realistic values of fractal dimension, because previous works considered only x‐ray structures. Binding site fractal dimension was calculated considering the ensemble of structures generated at different simulation times. A linear relationship between binding site fractal dimension and experimental binding free energies of the complexes was observed within 20 ns. Previous studies of the subject did not uncover this relationship. Regression model, coined FD model, was built to estimate binding free energies from binding site fractal dimension values. Leave‐ one‐out cross‐validation showed that our model reproduced accurately the absolute binding free energies for our training set (R 2 = 0.76; <|error|> =0.55 kcal/mol; SD error = 0.19 kcal/mol). The fact that such a simple model may be applied raises some questions that are addressed in the article. KEYWORDS drug design, fractal dimension, molecular descriptor, molecular dynamics, plasmepsin II, protein surface roughness 1 | INTRODUCTION The calculation of absolute binding affinities for protein‐ligand com- plexes remains one of the major challenges in computational struc- ture‐based drug design strategies. Several approaches to this problem have been developed. Some of them are based on free ener- gies calculations, such as the rigorous free energy perturbation (FEP) and thermodynamic integration (TI) methods. 1,2 However, FEP and TI approaches are quite time consuming, and they bring along associated sampling and convergence problems. 1,2 To overcome these problems, end point methods such as Linear Interaction Energy (LIE) and Molec- ular Mechanics‐Poisson Boltzmann/Generalized Born‐Surface Area (MM‐PB/GB‐SA) have been used. The main drawback of these methods is the use of several adjustable parameters, thus limiting its application to specific systems. 3 Most rapid methods for estimation of binding free energies, such as, so‐called empirical or knowledge‐based (statistical) scoring functions, have been proposed. These methods are based on very simple energy functions 4-6 or on the frequency of occurrence of different atom‐atom contact pairs in complexes of known structure. 7,8 The simplicity of the energy function along with the lack of conformational sampling and explicit water treatment makes these approaches very fast, but usually at the cost of accuracy. 9 Thus, further development of fast and accurate methods to complement the already existing ones is still needed. Received: 1 February 2017 Revised: 10 August 2017 Accepted: 11 August 2017 DOI: 10.1002/jmr.2661 J Mol Recognit. 2017;e2661. https://doi.org/10.1002/jmr.2661 Copyright © 2017 John Wiley & Sons, Ltd. wileyonlinelibrary.com/journal/jmr 1 of 9