Modeling Blood-Brain Barrier Penetration Using E-sol

Tutorial Created with Software Release: 2026-3
Topics: Small Molecule Drug Discovery
Products Used: Jaguar, MacroModel, Maestro

Tutorial files

59 MB
This tutorial is written for use with a 3-button mouse with a scroll wheel.
Words found in the Glossary of Terms are shown like this: Workspacethe 3D display area in the center of the main window, where molecular structures are displayed

 

Tip: You can hover over a glossary term to display its definition. You can click on an image to expand it in the page.
Abstract:

This tutorial demonstrates how to compute energy of solvation (E-sol) values and analyze the results to assess blood-brain barrier (BBB) permeability for small molecule ligands.

Tutorial Content
  1. Introduction
  2. Creating Projects and Importing Structures
  3. Setting up E-sol Calculations
  4. Analyzing the Results
  5. Conclusion and References
  6. Glossary of Terms

1. Introduction

The hydration free energy (HFE), defined as the energy of transferring a compound from the vacuum gas phase to water, quantifies how favorably a molecule interacts with water. A more negative HFE value indicates a stronger preference for the aqueous phase. There are two approaches available within the Schrödinger suite to calculate HFE:

  • E-sol: A quantum mechanics (QM) approach to calculate solvation free energy to screen candidates for BBB penetration
  • FEP+: A molecular mechanics (MM) approach

In this tutorial, you will use the E-sol method to calculate the HFE. E-sol has been shown to be a good categorical predictor of Kp,uu (unbound brain-to-plasma partition coefficient) and Efflux Ratio (ER), which implies rather than accurately predicting the precise numerical value of these endpoints, it separates the compounds with favorable and unfavorable endpoint values.

Kp,uu and ER are briefly described below:

  • Kp,uu represents brain exposure as a combination of passive transport and active influx/efflux across the BBB. It quantifies how freely a central nervous system drug crosses the BBB and remains unbound in the brain tissue relative to the blood. Higher Kp,uu ensures effective drug concentrations at the brain at minimal doses, and reduced side-effects at other organs.
  • The BBB is enriched with efflux transporters like multidrug resistance protein 1 (P-glycoprotein (P-gp)) and breast cancer resistance protein. ER represents the ratio of drug secretion (how fast a drug is pumped out) to drug absorption (how fast a drug is pumped in). A high ER indicates that the drug is rejected and pumped out into the blood almost as fast as it enters the brain, and hence, for a brain-targeted therapy to be effective, a low ER is aimed.

While E-sol is the recommended first-line approach for computing HFE (as a categorical predictor for Kp,uu and ER), it relies on MacroModel-based conformational sampling and an implicit solvent model, both of which have well-defined limits of applicability. For compounds that are large (>500 MW), possess more than 8-10 rotatable bonds, are macrocycles, or are otherwise peptidic, MacroModel's conformational search may not adequately sample the relevant conformational space, and the resulting E-sol values tend to show weaker trends with Kp,uu or ER endpoints. In these cases, FEP+ is the recommended alternative, since it does not depend on the same implicit-solvent and conformer-sampling assumptions that constrain the domain of applicability of E-sol.

The QM-based E-sol involves the following steps:

  1. Starting from practically any charge state or tautomer of each ligand, the structure and state penalty (SP) of the lowest energy uncharged neutral form are determined using Epik.
  2. The conformers of each ligand are generated using MacroModel. Conformers generated in vacuum, in water or a combination of both may be selected. By default, 12 conformers are generated for each phase. Increasing this number may improve E-sol accuracy while also increasing the overall computational cost of the workflow.
  3. The geometry of each conformer is optimized in the gas phase using Jaguar. A DFT-based approach is recommended (more accurate but slower), but MLFF is an alternative if higher throughput is necessary. Evac is obtained.
  4. The solvation energy is calculated using PBF, a Poisson-Boltzmann implicit water model. Ewat is obtained. Note that if full MLFF preset is chosen, CPCM-X solvation model is used, which is computationally much faster.
  5. A raw E-sol (unmodified by the state penalty) for each compound is calculated using the conformer with lowest energy in vacuum (step 3) and water (step 4) as: E-sol raw = Ewat,min – Evac, min
  6. The final E-sol value is calculated from the hydration energy (E-solraw) and the state penalty (SP) as follows:E-sol = E-solraw – SP

This tutorial uses a small set of compounds from the OpenADMET training set. In this tutorial, you will learn how to compute the energy of solvation and analyze the results to assess the BBB permeability for small molecule ligands. Since these compounds lack Kp,uu data but have known ER, you will perform a retrospective E-sol evaluation using the ER endpoint.

2. Creating Projects and Importing Structures

  1. Open Maestro and create a new project named solvation_energy.prj for this tutorial.
  1. Download the tutorial zip file including input files and reference outputs here: https://www.schrodinger.com/sites/default/files/s3/release/current/Tutorials/zip/solvation_energy.zip
  2. After downloading the zip file, unzip the contents in your Working Directorythe location where files are saved for ease of access throughout the tutorial.

Figure 2-1. Importing the structures.

We have provided a small set of compounds from the OpenADMET training set with known endpoints like Kinetic Solubility (μM), Caco-2 Efflux Ratio, Mouse Plasma Protein Binding (% unbound), and Mouse Brain Protein Binding (% unbound). To proceed with this tutorial, you will import these ligand structures from the tutorial files.

  1. Go to Files > Import Structures.
  2. Navigate to and choose openadme-uncurated.maegz in your Working Directorythe location where files are saved and click Open.
    • The structures are added to the Entriesa simplified view of the Project Table that allows you to perform basic operations such as selection and inclusion and the top entry is included in the Workspacethe 3D display area in the center of the main window, where molecular structures are displayed.

3. Setting up E-sol Calculations

Before submitting the E-sol calculations, we strongly recommend curating the compounds using known experimental data like Kinetic Solubility, Recovery %, and Brain/Plasma Tissue Binding (% unbound) to reduce noise coming from the measurements and get the clearest picture of the correlation of E-sol with Kp,uu and ER endpoints.

Note

To build robust BBB predictive models, we recommend filtering compounds with low solubility (< 1 μM), low recovery (< 50%), and high tissue-binding percentage (unbound < 1%). These compounds likely fail to penetrate the BBB for reasons unrelated to passing the barrier itself.

The following interactive tutorial walks you through curating the compounds against known experimental data (in this case, Solubility (μM) and Unbound Brain/Plasma Tissue %), setting up and running E-sol calculations.

The E-sol Calculation panel offers three presets, each striking a different balance between accuracy and computational cost:

  • MLFF Optimization achieves a categorical accuracy of 77.4% at a moderate speed of 48 CPU-minutes per ligand, making it the recommended choice for production-quality predictions.
  • Full MLFF trades some accuracy (69.4%) for a substantial speed advantage – at only 18 CPU-minutes per ligand, it is well suited for rapid screening of large compound sets where throughput matters more than precision.
  • Full DFT performs geometry optimization and single-point energies entirely at the DFT level without machine-learning force field filtering. It delivers a categorical accuracy of 74.5% but at a much higher cost of roughly 30 CPU-hours per ligand. It is well suited for spot validation of a small number of high-priority compounds (where the highest confidence is needed) rather than routine use.

Please note the reported accuracy metrics are based on a broader validation dataset, not the specific set of compounds used in this tutorial.

We recommend using the MLFF Optimization preset for most E-sol calculations. It uses the MPNICE machine learning force field to pre-filter and geometry optimize the conformational ensemble, then scores only the 6 best conformers with a full DFT single-point energy calculation. (M06-2X/LACVP** with PBF implicit solvation). This two-stage strategy preserves the QM-quality final energies that govern E-sol accuracy while eliminating the most expensive step in the workflow – DFT geometry optimization of every conformer. In practice, MLFF Optimization typically matches the accuracy of the most exhaustive DFT protocol at roughly one-third the computational cost, making it the best balance of speed and reliability.

What happens if you do not have any MLFF licenses?

If you lack the MLFF licenses, the calculation defaults to the xTB level of theory rather than failing. Please note that the runtime and the results obtained using xTB may differ from the recommended MLFF Optimization output.

4. Analyzing the Results

In this section, you will evaluate how well the computed solvation energies discriminate between compounds with favorable and unfavorable ADME endpoints – in this tutorial, Efflux Ratio. 

Instead of performing regression analysis which fits a linear correlation between E-sol and the endpoint (Kp,uu or ER), we recommend performing classification analysis. The recommended classification approach constructs a confusion matrix from categorical cutoffs (e.g., E-sol and endpoint thresholds) and evaluates discriminating power using Receiver Operating Characteristic (ROC) curve statistics, such as the true positive and false positive rates, which better capture E-sol's categorical accuracy.

In this example, the optimal E-sol cutoff was determined to be -21.22 kcal/mol. Compounds with E-sol more negative than this value are predicted to be P-gp substrates with high efflux, while those less negative are predicted to have low efflux. The resulting classification accuracy of 93%, combined with an AUC of 0.96 and a Fisher exact test p-value of 0.005, confirms that E-sol is a statistically significant predictor of efflux liability for this chemical series. Importantly, the specificity of 1.00 means no compounds were falsely predicted to have low efflux.

Two considerations are worth noting:

  • The E-sol cutoff of -21.22 kcal/mol is specific to this dataset and chemical series, but would be useful to prospectively classify new ideas in this project as high-efflux or low-efflux.
  • While R² = 0.42 might appear discouraging when taken at face value, it is entirely consistent with the expected noise in both the E-sol calculation (~ 0.5-0.8 kcal/mol from conformational sampling stochasticity) and the efflux ratio measurement itself. This is precisely why categorical accuracy, not linear correlation, is the appropriate metric for evaluating E-sol performance.

Note

It is important to note that the E-sol cutoff is series-specific: -18 to -16 kcal/mol can be a general starting range for Kp,uu when no project data is available, but as your own experimental data accumulates, the cutoff should be recalibrated. Some chemical series are systematically shifted left or right, and a different cutoff may provide better discrimination for your specific scaffold.

5. Conclusion and References

In this tutorial, you learned how to compute the energy of solvation and analyze the results to assess BBB permeability.

Once you have established and validated an E-sol cutoff on a retrospective dataset from your chemical series, you can apply that same cutoff to newly designed or enumerated compounds to triage them before synthesis (for example, refer to this publication to learn how E-sol predictions guided the design of brain penetrant DLK inhibitors). Compounds with E-sol values more negative than the cutoff are predicted to have unfavorable efflux or poor brain penetration and can be deprioritized or discarded from the synthesis queue – they are predicted to be too hydrophilic to cross the BBB passively and are likely substrates for P-gp efflux. Conversely, compounds with E-sol values less negative than the cutoff are predicted to have favorable permeability and should be prioritized for synthesis and testing.

6. Glossary of Terms

Entries - a simplified view of the Project Table that allows you to perform basic operations such as selection and inclusion

Included - the entry is represented in the Workspace, the circle in the In column is blue

Project Table - displays the contents of a project and is also an interface for performing operations on selected entries, viewing properties, and organizing structures and data

Recent actions - This is a list of your recent actions, which you can use to reopen a panel, displayed below the Browse row. (Right-click to delete.)

Scratch Project - a temporary project in which work is not saved, closing a scratch project removes all current work and begins a new scratch project

Selected - (1) the atoms are chosen in the Workspace. These atoms are referred to as "the selection" or "the atom selection". Workspace operations are performed on the selected atoms. (2) The entry is chosen in the Entries Table (and Project Table) and the row for the entry is highlighted. Project operations are performed on all selected entries

Working Directory - the location where files are saved

Workspace - the 3D display area in the center of the main window, where molecular structures are displayed