A key consideration in the use of binding free energy calculations in real world (industrial or clinical) settings is the reproducibility of the results. Other considerations include computational cost and calculation stability. ESMACS protocols offer advantages in both these regards as they make use of relatively simple and fast classical MD simulations compared to many parallel simulations of intermediate states as required in alchemical calculations of absolute binding free energies . We have shown that the results obtained in this study are robust to changing the ligand charge generation protocol (to use AM1-BCC instead of Gaussian/RESP) and the forcefield used to parameterize the protein (from Amber ff14SB to ff99SBildn). The use of ensemble simulations is the key to obtaining this reproducibility as individual replicas in ensembles varied by as much as 15 kcal mol -1 (which is in line with our own and other groups previous results ). Despite this, performance differences were found for all initial protocols when simulations were initiated from different crystal structures. This observation, along with the fact that some ligands which have very similar experimental binding energies were widely separated even using protocols which accounted for receptor flexibility (1traj-ar and 2traj-ar), prompted us to investigate potential enhancements of the pure MMPBSA protocol. Specifically, we looked at the inclusion of an explicit ligand hydration shell in the receptor and variational entropy which had previously been investigated for single replica simulations by Aldeghi et al. (though they only investigated what we would term "1traj" calculations). These additional components capture chemical and physical features of the system neglected by MMPBSA but at minimal computational cost, a key consideration for practical binding affinity calculation applications. The entropy term reduced extreme outliers but at the expense of decreased overall ranking performance. This observation replicates that obtained by Aldeghi et al. for the DIV dataset bound to BRD4, although they found the term improved results for sensitivity based datasets including multiple proteins. When less than five water molecules were incorporated into the receptor our rankings were improved with the best ranking across the full dataset obtained using this in combination with the 2traj-ar protocol. The most important observation of our work, however, is that the inclusion of these bound water molecules considerably reduced the performance difference between simulations initiated from models based on different crystal structures. A criticism of continuum based methods is that they will be incapable of capturing the effect of crucial water molecules, possible activity cliffs, etc, that are now a well understood feature of structure-activity relationship (SAR) landscapes and medicinal chemistry lead optimization. Here it is shown again how this challenge can be met, with the simple inclusion of explicit water molecules. Future work should address how to consider this in prospective application scenarios and in a wider range of protein targets.