Molecular Simulation and Statistical Learning Methods toward Predicting Drug–Polymer Amorphous Solid Dispersion Miscibility, Stability, and Formulation Design
1.28.2021
Originally Published in Molecules
Amorphous solid dispersions (ASDs) have emerged as widespread formulations for drug delivery of poorly soluble active pharmaceutical ingredients (APIs). Predicting the API solubility with various carriers in the API–carrier mixture and the principal API–carrier non-bonding interactions are critical factors for rational drug development and formulation decisions. Experimental determination of these interactions, solubility, and dissolution mechanisms is time-consuming, costly, and reliant on trial and error. To that end, molecular modeling has been applied to simulate ASD properties and mechanisms. Quantum mechanical methods elucidate the strength of API–carrier non-bonding interactions, while molecular dynamics simulations model and predict ASD physical stability, solubility, and dissolution mechanisms. Statistical learning models have been recently applied to the prediction of a variety of drug formulation properties and show immense potential for continued application in the understanding and prediction of ASD solubility. Continued theoretical progress and computational applications will accelerate lead compound development before clinical trials. This article reviews in silico research for the rational formulation design of low-solubility drugs. Pertinent theoretical groundwork is presented, modeling applications and limitations are discussed, and the prospective clinical benefits of accelerated ASD formulation are envisioned.
Keywords: solubility; bioavailability; drug development; amorphous solid dispersions; molecular modeling; molecular dynamics; machine learning
Drug solubility plays a vital role in the drug discovery and development pipeline, with around 40 to 70% of drugs displaying poor aqueous solubility due to highly stable solid lattice arrangements and/or pronounced hydrophobicity [1,2,3]. Decreases in bioavailability due to poor solubility are prevalent for oral dosages in which therapeutic compounds must be in solution to reach systemic circulation and achieve clinically relevant concentrations at the site of action. Considerable attention during discovery and lead optimization is therefore given to compound structural characteristics and the resulting physicochemical properties that affect solubility [4,5], with notable parameters, among many others, being a drug’s molecular weight, octanol:water partition coefficient (logP) [6], and possible acid-base ionization equilibria (pKa) [7]. There is complex interplay and tradeoff between the optimal chemical scaffold and functional groups for physiological response at the target and the ensuing effect on solubility and oral bioavailability [8]. Overcoming low drug solubility while maintaining critical molecular components for therapeutic activity often necessitates precise formulation strategies during the development stage [9,10,11].
A number of approaches exist for addressing the poor solubility of BCS (Biopharmaceutics Classification System) class II (low solubility, high permeability) and IV (low solubility, low permeability) active pharmaceutical ingredients (APIs) [1], including particle size reduction (micronization) [12,13], salt formation and co-crystallization [14,15], lipid-based delivery [16,17], and cyclodextrin complexation [18]. Among the variety of strategies, amorphous solid dispersions (ASDs) have seen consistent use as formulation and drug delivery systems, especially for BCS class II APIs [19,20,21,22]. In this system, lipophilic API is dispersed with a carrier compound, often a water-soluble amphiphilic polymer excipient, forming a multi-component matrix (Figure 1A). Compared to the more thermodynamically stable API crystalline form, once dispersed within the carrier matrix, the API occupies an elevated energetic state, and solubility and dissolution rate enhancement of the mixed ASD phase relative to crystalline API results from increased solvent-exposed surface area [23]. While this instability is beneficial due to the increased dissolution rate, relative ASD instability incurs a practical challenge of crystallization out of the amorphous state during shelf life [24]. An ideal ASD must strike a delicate balance between sufficient API saturation within the carrier and recrystallization [25]. Predicting the stability and miscibility of a given API within a carrier and the resulting solubility and dissolution rate of the mixture (Figure 1B) is of paramount importance to guide pharmaceutical formulation and production decisions.
Several experimental methods have been developed and applied towards understanding and characterizing ASD solid-state nanostructure and properties [26,27]. Key solid-state endeavors include probing for regions of API crystallinity and the presence of impurities, both of which can serve as further crystal nucleation sites [28]. Measurement of the glass transition temperature (Tg) routinely guides API–carrier miscibility, homogeneity, and physical stability evaluations [29]. Experimental deduction of intermolecular interactions occurring between API and carrier yields fundamental indicators of miscibility and overall ASD stability [30]. While experimental methods have yielded considerable understanding of the nanostructure and solid-state properties of ASD formulations, the metastable ASD system can also provide a significant challenge for experimental characterization. Appropriately matching a given excipient carrier to an API remains a chance approach in terms of concurrently optimizing API–carrier miscibility, ASD rate of drug release, and long-term ASD storage stability. The difficulty of rational excipient screening and lengthened development cycles due to ASD characterization experiments both contribute to the relatively small number of drugs prepared as solid dispersions ultimately reaching the market, as shown in Table 1 [31,32].
Computational approaches allow further understanding of ASD phenomena, especially when working in tandem to complement and inform experiments, and to assist in the development of predictive models for rational formulation design (Figure 1C). Molecular modeling and simulation techniques assist in elucidating important stabilizing intermolecular interactions between API and carrier, predict solubility parameters, simulate ASD formation and dissolution mechanisms, and generate descriptors for quantitative structure–property relationships (QSPR) [33,34]. Furthermore, machine learning methods are only just recently showing prominent applications to solid dispersion modeling and prediction. In this review, recent examples of molecular and statistical modeling and prediction as applied to the study of ASD intermolecular interactions, miscibility, formation, and stability are highlighted.
Two widely applied theoretical approaches for predicting the solubility of solid dispersions are the Hildebrand and Hansen solubility parameters and the Flory–Huggins (FH) interaction parameter [35,36,37,38,39]. Both models were originally developed to probe thermodynamic solute–solvent miscibility. In the context of API–polymer carrier miscibility, the theories are extended for polymer carrier systems by assuming the polymer as the solvent and the API as solute in a mixed system. The principle ultimately governing favorable mixing is the Gibbs free energy change upon mixing, given by Equation (1):
where T is the absolute temperature, and ∆Gmix, ∆Hmix, and ∆Smix represent the change in free energy, enthalpy, and entropy between the unmixed and mixed components at constant pressure.
The Hildebrand solubility parameter, δ, is given as the square root of the cohesive energy density (CED) [38,40], defined as the per volume energy required to break intermolecular bonds and vaporize one mole of molecules from the liquid phase (Equation (2)):
where Evap is the energy of vaporization, ∆Hvap is the change in the enthalpy of vaporization associated with the phase change, V is the molar volume of the liquid at the temperature of vaporization, R is the ideal gas constant, and T is the absolute temperature. Since the cohesive energy density is proportional to the enthalpy of vaporization, the magnitude of δ gives an indication of the strength of the intermolecular forces in the condensed state. The change in enthalpy upon mixing can be defined in terms of solubility parameters, as shown in Equation (3):
where Vt is the total volume of the mixture, δ and Φ are the solubility parameters and volume fractions of the individual mixture components, respectively. The entropic contribution to the free energy of mixing, ∆Smix, is positive as more ordered states disperse during mixing. Hence, favorable mixing occurs, as indicated by a negative ∆Gmix, when the individual solubility parameters of the mixture components (δ1 and δ2) are close in value and ∆Hmix approaches zero. As applied to API–polymer carrier systems, a heuristic for estimating miscibility based on differences between component solubility parameters indicates ∆δAPI-polymer < 7 MPa1/2 as miscible and ∆δAPI-polymer > 10 MPa1/2 as immiscible, and values in between the thresholds as borderline miscible [37].
Hansen extended the Hildebrand solubility parameter to a three-dimensional space based on intermolecular attractive forces in order to more accurately model solutions containing polar compounds, as Hildebrand’s initial theory was developed based on hydrocarbon solvents [39]. In Hansen space, the solubility parameter is decomposed into a sum of contributions from three non-bonding interactions, as shown in Equation (4):
where δh, δp, and δd are the hydrogen-bonding, polar, and dispersion contributions, respectively.
Group contribution methods calculate solubility parameters from molecular structure and can provide a rapid method for miscibility estimation if the structural fragment constants are available for all components of the compound of interest [41]. Two commonly applied group contribution methods include that of Fedors [42] and the Hoftyzer-Van Krevelen method [43], the latter shown below (Equation (5)) for estimating the Hansen solubility parameter:
where Ehi, Fpi, and Fdi are the molar group attraction constants due to hydrogen-bonding, polar, and dispersion components, respectively, and V is the molar volume.
Another common theory for estimating miscibility between API and polymers is Flory–Huggins theory. One disadvantage of solubility parameters addressed by FH theory is the entropic contribution to ∆Gmix, which potentially becomes important to incorporate due to the structural complexity and flexibility of polymeric carriers as the mixture deviates from an ideal solution [36]. In FH theory, the API–polymer system is modeled as a mixture where the average energy of interaction between each component depends on the FH interaction parameter, χ, and is related to the free energy change upon mixing via Equation (6) [11]:
where r1 and r2 represent the number of monomer units in each polymer, and χ is the unitless FH interaction parameter that relates the strength of the API–polymer interaction. Equation (6) describes the mixing of two polymers, and in the case of API–polymer mixing, components 1 and 2 are drug and carrier, respectively, and r1 = 1. The first two terms of Equation (6) (right-hand side) relate to the combinatorial entropy of mixing, while the last term corresponds to the enthalpic contribution. Since mixing generally increases entropy, it is assumed that the dominance of the interaction parameter χ within the enthalpic term decides miscibility. The key threshold is χ = 0.5, with values of χ > 0.5 indicating likely phase separation between the two components, while χ < 0.5 predicts API–polymer miscibility. Negative values of χ further indicate that drug–polymer forces are greater than self-association forces (API–API, polymer–polymer).
Upon dispersion of the API within the polymer carrier, hydrogen-bonding, acid–base/ionic, dipole–dipole, and van der Waals non-bonding interactions discourage self-association, increase API–carrier interaction energy, and promote ASD physical stability over time [24,30,44]. Understanding miscibility between a given API and carrier requires the determination of the extent and types of interactions taking place between APIs and carriers in the ASD matrix. While this is a challenging task from a modeling perspective, as the carriers are large, highly flexible, multi-functional polymers, there are several molecular modeling techniques and software available for studying ASDs in detail (Table 2). The following sections describe quantum mechanics (QM) and molecular mechanics (MM)/molecular dynamics (MD) approaches for investigating ASD intermolecular interactions, solubility and FH miscibility parameters, formation, and stability.
Table 2. Molecular modeling software applied in this review to studying ASD systems.
QM methods such as density functional theory (DFT) provide insight into the non-bonding interactions taking place between API and polymer carriers through direct construction of the hypothesized molecular complexes and subsequent calculation of the associated interaction energy. Recent examples of ASD non-bonding interaction modeling highlight DFT as a relatively rapid (as compared to molecular dynamics; see Section 3.2.1) yet powerful computational methodology for predicting the interaction energy of model complexes and the identification of ionic and strong hydrogen-bonding pairs. As a computational modeling tool for studying solid-state pharmaceutical compounds, periodic DFT computations have seen increased applications in the literature [82]. Furthermore, combined DFT/MD methods may also be applied (e.g., Car–Parrinello MD) in which standard empirical force fields (see Section 3.2.1) are replaced with an ab initio description of the electronic structure [83,84,85,86]. Pure QM approaches commonly seen in ASD modeling produce a static snapshot and therefore require adequate, often manual, sampling of the conformational space as a prerequisite. Additionally, in order to facilitate DFT calculations of the large chemical systems, polymer molecules can be truncated to their monomeric or dimeric forms to model systems, which drastically reduces the time of the calculation. This reduction of calculation runtime as a result of the truncated model system allows for more systematic and detailed structural investigations but ultimately presents a tradeoff between a chemically realistic system and computational cost. Simulating a larger and more structurally realistic system does not necessarily result in an increase in accuracy but remains a desirable objective in ASD molecular modeling due to the ability to predict a higher number of properties and more diverse chemical behavior.
Miscibility and stability both depend on strong interactions between API and carrier, and understanding the critical interactions controlling ASD stability and solubility is a crucial step toward rational excipient screening and ASD formulation [87]. DFT applications in this context probe hydrogen-bonding and charge interactions. Maniruzzaman et al. investigated the interactions between cationic propranolol hydrochloride and diphenhydramine hydrochloride APIs with anionic Eudragit L100 and Eudragit L100-55 polymer carriers [46]. The polymer compounds were represented by their monomer forms, with multiple arrangements between monomer and polymer constructed, considering both the absence and presence of an explicit chloride anion within the complex. DFT computations found that the hydroxyl and amine groups of the API formed favorable hydrogen bonds with the ester and hydroxyl groups of the monomer units. Inclusion of the chloride ion in the complex was found to disrupt hydroxyl group interactions, but not amine group hydrogen-bonding, suggesting that the optimal interaction is between the amine groups of propranolol and diphenhydramine with the carbonyl groups of the L100 and L100-55 polymers. The proposed hydrogen-bonding motif was confirmed using nuclear magnetic resonance (NMR) and X-ray photoelectron spectroscopy experiments. Maniruzzaman et al. built upon the monomeric findings with a systematic study of dimeric API–carrier hydrogen-bonding interactions [47]. DFT computations showed 5 to 8 kcal/mol increase on average in binding energy between API and polymer for each additional hydrogen-bonding interaction, with the most stabilizing interaction being tertiary amine and carboxyl hydrogen-bonding (20 to 28 kcal/mol). In similar work, Nie et al. probed the molecular interactions between the antileprosy therapeutic clofazimine with hypromellose phthalate as carriers within an ASD system [48]. Partial charge analysis identified potential strong donors and acceptors, and a truncated model system with acetic acid as a structural proxy for the hypromellose phthalate carboxylic acid was used for DFT geometry and energy calculations. Based on these calculations and spectroscopic studies, formation of an ion pair complex was proposed as critical for API–carrier miscibility.
In an example of DFT calculations providing insight into ASD storage stability, Wang et al. studied the ability of polyvinylpyrrolidone (PVP) to inhibit crystallization of PVP–resveratrol and PVP–griseofulvin ASDs [88]. DFT calculations and Fourier-transform infrared (FTIR) spectroscopy investigated relative API–API and API–polymer interaction strengths, revealing greater stabilizing hydrogen-bonding interactions within the resveratrol complex than the griseofulvin–PVP system. The stronger interaction energy translated to higher stability, measured as a continued high dissolution rate after 90 days of storage of the PVP–resveratrol system compared to the lower dissolution rate of the griseofulvin–PVP ASD.
With the advent of computing power over the last decade, MD simulations of amorphous API–API [89,90], API–micelle [91], polymer–polymer [92,93,94,95], polymer–membrane [96], and polymer–plasticizer [97] are prominent throughout the literature. The scope of the following sections focuses specifically on MD simulations of ASD systems composed of APIs interacting with polymer carriers. As applied to ASDs, MM-based docking and MD simulations predict a range of phenomena and properties. Recent literature examples are detailed pertaining to docking simulations probing API–carrier interactions, MD simulations to generate solubility and FH interaction parameters as predictors of API–carrier miscibility, and MD simulations probing ASD formation, polymer interactions with API crystals, and ASD hydration and dissolution mechanisms.
Studying systems as large and complex as ASDs often necessitates the use of MD in order to approach a more physically realistic view of the phenomena under study. Application of MD requires that all components of the system are parametrized by a MM–based force field dictating the total energy as a sum of empirical potential energy functions (potentials) [98]. Common force fields applied to API–excipient ASD modeling (Table 3) are Condensed-phase Optimized Molecular Potentials for Atomistic Simulation Studies (COMPASS) [99,100], Polymer Consistent Force Field (PCFF) [101], CHARMM [102,103], Amber [104,105], and Optimized Potential for Liquid Simulations (OPLS) [106]. The force field describes energy as dependent on atomic position, and its subcomponent potentials can be divided into two types: (1) those for bonded atoms, describing bond lengths, angles, torsions, and (2) those for non-bonded atoms, accounting for Coulombic and van der Waals interactions. MD adds a temporal component to a molecular system by probing chemical phenomena over time. Chemical phenomena critical to predicting API–carrier miscibility and molecular mobility [107,108], such as hydrogen bond formation or motion of flexible chains, can emerge and recede over the course of the simulation and can be taken as average values over the simulation length (typically nanoseconds to microseconds) and even as averages from multiple runs originating from different initial states. Classical MD follows Newtonian laws of motion, tracking atomic position, velocity, acceleration, force, and energy across all atoms for each time step of the simulation. The total energy at a given time is then a function of atomic positions, where the positions are updated for each time step of the simulation in an iterative fashion based on the forces (and therefore acceleration) encountered by the atoms.
Table 3. Recent studies of ASD API–polymer miscibility from MD-derived solubility and interaction parameters.
Docking provides an algorithmic method for rapidly sampling and scoring complexes between API and polymer carriers and, within this context, is often used to generate favorable preliminary binding conformations, which are then subjected to more in-depth analysis and energy refinement through MD simulations. Barmpalexis et al. investigated ternary ASD systems of ibuprofen and carbamazepine APIs, Soluplus (SOL) polymer, and polyethylene glycol (PEG) plasticizer [50]. As with QM approaches, monomeric or dimeric model systems of the polymer are often employed to probe specific interactions and reduce computational cost. Fule et al. undertook docking simulations to investigate intermolecular interactions between the antimalarial lumefantrin with monomers of SOL, Kollidon VA-64, and vinylpyrrolidone-vinyl acetate polymers [78]. Strong hydrogen-bonding was identified between SOL and vinylpyrrolidone-vinyl acetate monomers and lumefantrin hydroxyl, chlorine, and amine groups. Docking can be used not only to identify potential non-bonding interactions but also to rapidly screen the stability of potential API tautomeric forms prior to MD simulations. Gangurde et al. applied docking to assess the binding energy of diketo and keto-enol forms of curcumin with Eudragit EPO monomers and dimers [79], finding that the diketo tautomer was overall favored, with binding controlled primarily by van der Waals interactions. These initial docking studies were followed by MD simulations to refine the docking energies and provide more detailed structural insight into curcumin–Eudragit EPO interactions. In a comprehensive docking study, Macháčková et al. detail systematic docking simulations with full-size polymers [60]. In this study, the anchoring ability of cyclosporin A was investigated with short (~7 nm), medium (13–14 nm), and long (~20 nm) chain polymers of l/d–polylactide, chitosan, polyglycolic acid, PEG, and cellulose. For each type of polymer modeled at a given length, one million complexes were generated between stationary cyclosporin A and flexible polymer, and the resulting interaction energies predicted chitosan and cellulose as most miscible with cyclosporin A.
In the context of MD atomistic simulations, cohesive energy density and solubility parameters are derived by comparing the total energy of the system in a vacuum (infinitely separated) to the interacting, modeled amorphous state (bulk). This difference in energy relates to the energy of vaporization, and the MD-derived Hildebrand solubility parameter is taken as Equation (7) [109]:
where Evacuum is the potential energy of a single system component (API or polymer) in vacuum, Ebulk is the potential energy of the same component in the simulated system, and V is the volume of the periodic cell. Since MD force field potentials explicitly account for Coulombic and van der Waals forces, the predicted solubility parameter can be further decomposed into its non-bonding electrostatic and dispersion terms, yielding the Hansen solubility parameter. From MD-derived electrostatic (accounting for both polar and hydrogen-bonding interactions) and dispersion contributions to the total potential energy, the Hansen solubility parameter is given by Equation (8):
where Ee and Ed are the electrostatic and dispersion contributions to the potential energy, respectively, between infinitely separated molecular components and the simulated interacting molecular assembly for all constituents of the system (∆E = Evacuum − Ebulk) [68,74]. MD simulations can also be used to predict FH interaction parameters. The computed cohesive energy densities (CED = ∆E/V) are related to the energy of mixing, ∆Emix, and finally to the FH interaction parameter, χ, via Equation (9):
and Equation (10):
where the Φ terms represent the volume fractions of the individual components within the mixed system. It should be noted that the Flory–Huggins model does not account for entropic effects related to the structural properties of the polymer. Entropic energetic penalties such as restricted polymer movement within the ASD matrix from API–carrier intermolecular hydrogen-bonding are not explicitly accounted for. Accurate calculation of the entropic contribution is considered a significant cause of the observed discrepancies between predicted and experimental miscibility evaluations (see end of Section 3.2.4) [110].
MD methodology provides an attractive and versatile approach for calculating solubility and FH interaction parameters (Table 3). Gupta et al. applied MD simulations to compute solubility parameters and compared the MD-derived values with solubility parameters calculated using the Hoftyzer–Van Krevelen and Hoy group contribution methods [57]. Calculated solubility parameters predicted that indomethacin API was miscible with polyethylene oxide (PEO), borderline miscible with sucrose, and immiscible with glucose. These predictions agreed with differential scanning calorimetry (DSC) experiments, which showed indomethacin miscibility with PEO and general immiscibility with sucrose and glucose carriers. Yani et al. investigated ibuprofen, fenofibrate, and alprazolam API miscibility with poly(vinylpyrrolidone-co-vinyl acetate) 64 (PVP–VA 64), hydroxypropyl methylcellulose (HPMC), and Eudragit EPO polymer carriers using Hildebrand and Hansen solubility parameters [63]. MD simulations were used to derive Hansen solubility parameters and to calculate the average hydrogen bond lifetime and total hydrogen bond interaction energy per monomer between API and polymer functional groups. Based on the predicted solubility, hot melt extrusion was used to produce ASDs for ibuprofen/PVP–VA 64, ibuprofen/Eudragit EPO, and fenofibrate/PVP-VA 64. Weaker predicted relative intermolecular interactions in the fenofibrate/PVP–VA 64 ASD coincided with its observed recrystallization. Barmpalexis et al. took the most stable poses from docking simulations between ibuprofen and carbamazepine binding with SOL/PEG polymer/plasticizer blends (Section 3.2.2) and evaluated miscibility using Hoftyzer–Van Krevelen and Hildebrand solubility parameters [50]. Both ternary systems were shown as miscible from calculated solubility parameters and DSC experiments. FTIR spectroscopy was applied to probe the intermolecular interaction within the ASD and confirmed strong hydrogen-bonding observed from docking and MD between the carbamazepine primary amine group and the SOL carbonyl and amide groups. MD studies by Iesavand et al. calculated the enthalpy of mixing (∆Hmix) and Hildebrand solubility parameters of 6-Mercaptopurine API with polylactic acid (PLA) and PEG-modified PLA [58].
The FH interaction parameter has also seen recent application for probing API–carrier miscibility in combination with MD simulation. Macháčková et al. used insights from docking (Section 3.2.2) and MD–derived FH interaction parameters to identify chitosan and cellulose as the most miscible polymers with cyclosporin A [60]. Additionally, MD probed the polymer flexibility as measured by radius of gyration and elucidated all hydrogen-bonding interactions within the API–polymer complexes, revealing that miscibility primarily correlated with polymer chain length, cyclosporin A–polymer interaction energy, and the number of hydrogen bonds between cyclosporin A and polymer. Xiang et al. computed both solubility and FH interaction parameters between amorphous systems of indomethacin and PVP, varying API, polymer, and water composition in each MD simulation [68]. While differences in the solubility parameters of indomethacin and PVP predicted close to borderline miscibility (δIMC-PVP = 6.5 MPa1/2), FH interaction parameters predicted complete miscibility (χIMC-PVP = − 0.61). In later work, Xiang et al. similarly calculated FH parameters for ASDs of felodipine API and HPMC polymer for several compositions of API, polymer, and water, with both predicted solubility and FH parameters showing complete miscibility across all modeled compositions [69]. Erlebach et al. investigated mixtures of indomethacin with PEG and PLA polymers, studying ASD formation from simulated annealing, API–polymer miscibility from MD–predicted FH interaction parameters, and encapsulation efficiency of the polymer carrier [56]. Eslami et al. calculated both solubility and FH interaction parameters for tacrine interacting with chitosan and polybutylcyanoacrylate (PBCA) polymers [74]. Tacrine showed higher miscibility with PBCA with both approaches, and the MD simulations further revealed that increases in polymer chain length result in higher MD–derived interaction energies. In a recent study, Kapourani et al. disclosed a computational and experimental evaluation of ASDs formed between simvastatin API and PVP [52]. MD-computed Hansen solubility and FH interaction parameters predicted miscibility, which was verified using DSC experiments.
While MD–derived miscibility assessments provide utility and often complement experimental results, the accuracy of these methods remains largely untested as direct comparisons between MD–predicted solubility and FH interaction parameters and experimental values are uncommon. Rather, comparisons to experimental results typically rely on more qualitative observations. Examples include FTIR confirmation of strong intermolecular interactions and single observed Tg values for ASDs predicted as miscible, or the onset of rapid ASD recrystallization for a system predicted as immiscible (or even complete lack of ASD yield upon synthesis). Within a small scale of studied APIs and carriers, predicted and observed trends tend to agree (Table 3). With respect to quantitative experimental comparisons, Turpin et al. disclosed a systematic screen of API–carrier ASDs, allowing direct evaluation of predicted and experimental miscibility parameters [110]. High-throughput assays measured the miscibility limits for combinations of nine APIs and six polymer carriers (poly (glycerol–adipate) (PGA), three PGA derivatives, PVP, and PVP–VA 64). Both MD-predicted solubility and FH interaction parameters showed poor correlation with experimental miscibility. For FH interaction parameters predicted with PGA carrier, six of the nine API/PGA mixtures were predicted as having complete miscibility (χ < 0), in contrast to experimental results.
Knowledge pertaining to the mechanisms and rates of API dissolution [111], ASD formation, and ultimately ASD dissolution and drug release [112] from the carrier matrix is of great value to guide and streamline candidate excipients for a given poorly water-soluble therapeutic, resulting in lower development times and faster progression to clinical trials. Dissolution and drug release are highly mechanistically complex processes, with three underpinning hypotheses [19,20]: (1) rapid dissolution forming a supersaturated solution that promotes the generation of amorphous or crystalline API nanoclusters within the polymer carrier, (2) gradual release of amorphous API–polymer complexes from the greater bulk ASD complex, and (3) gradual release of API–polymer, but plasticization by water causes crystalline API–polymer clusters to form at the ASD–solvent interface. Which cases apply to a given ASD depends on the properties of the API, carrier, and the ASD formulation. In tandem with more time-consuming and resource-intense experimental mechanistic studies, recent literature examples show that computational methods can provide detailed insight into mechanistic aspects of ASD formation and dissolution [27]. Chan et al. applied MD simulations to the investigation of formation and dissolution of ASDs composed of ibuprofen and PEG, PVP, and poloxamer 188 (P188) [66]. ASD formation was modeled using the simulated annealing method followed by water immersion to investigate dissolution behavior, with small solid dispersion particles displaying rapid ASD dissolution. For larger particles, the dissolution behavior was dependent on the polymer. Higher API mobility within PEG and PVP polymers tended towards API aggregation, a precursor for API recrystallization (Figure 2), while P188 slowed API aggregation. A combined experimental and molecular modeling approach by Han et al. investigated a ternary ASD system containing glipizide API, PEG carrier, and sodium dodecyl sulfate/Tween 80 (Tween 80 is derived from polyethoxylated sorbitan and oleic acid) surfactants [67]. MD simulations studied ASD formation (simulated annealing), molecular mobility, and ASD dissolution.
Maintaining the supersaturated high energy state of ASDs increases solubility for poorly water-soluble APIs, but ASDs are not an equilibrium system, and there is a constant pull towards the thermodynamic stability of the API crystalline phase (Figure 2). The kinetic interplay between API dissolution and crystallization requires careful tuning of the properties of polymer carriers that promote extremely rapid dissolution, as fast API dissolution can lead to crystallization due to rapid supersaturation [113,114].
ASD formulation design decisions seek API–excipient mixtures that will inhibit or prolong crystallization, providing optimal long-term storage and physical stability, but modeling and predicting crystallization is highly challenging. The methodologies outlined in Section 3.1 and Section 3.2, related to quantifying the relative strength of API–polymer interaction and miscibility prediction of solubility and interaction parameters, are thermodynamic evaluations, giving no kinetic information on crystallization [115]. Furthermore, MD simulations that truly capture phase separation events when modeling API–excipient mixtures would need to run for exceedingly long simulation time scales in order to increase the likelihood of observing such phenomena. A recent example by Brunsteiner et al. discloses MD simulations of eight different API–polymer combinations of phenacetin and flufenamic acid APIs with Eudragit, polystyrene sulfonic acid, poly acrylic acid, and PVP polymers [72]. Thermodynamic mixing energy, kinetic diffusivity, and mobility were both investigated in detail, and this work highlights the importance of obtaining thermodynamic and kinetic information for confidence in ASD stability over time. As the 20 microsecond simulations of ~20,000 atoms took months to complete on standard computing machines, adapting and scaling such an approach for industrial and pharmaceutical applications such as API–excipient screening and long-term physical stability assessments currently remains limited by computational power. API crystalline lattices may also be constructed and modeled directly as an approach to furnish insight into ASD crystallization. Gao et al. studied interactions of diblock copolymer poly(ethylene glycol)-block-poly(lactic acid) (PEG-b-PLA) at the interface of tolazamide crystal slabs in aqueous conditions [76]. MD simulations found that PEG-b-PLA interacted most rapidly and strongly with the (001) face, followed by the (010) face, with minimal interactions at the interface of the (100) crystal face.
Of significant importance to ASD development and formulation are the effects of water on ASD stability and dissolution rate, with direct implications for both understanding moisture uptake during storage and in vivo solubility and absorption. Jha et al. used in-depth MD simulations to probe the intermolecular interactions and dissolution behavior between poorly soluble phenytoin API and cellulose-type polymer excipients [59]. ASD dissolution behavior was investigated by varying water weight percentages and tracking the rate of phenytoin diffusion, where diffusion constants were determined to have an exponential relationship with increases in water weight. Razmimanesh et al. investigated the effects of hydration on API loading efficiency from MD simulations of gemcitabine API and chitosan polymer carrier with varying concentrations of API [61]. Both the maximum drug loading efficiency and the strongest intermolecular interactions between gemcitabine and chitosan were predicted for the intermediate API concentration (40%) in simulated aqueous conditions. In addition to FH interaction parameters (Section 3.2.4), MD simulations by Xiang et al. studied the effects of hydration on hydroxypropyl methylcellulose polymer physical properties, finding that increasing water fractions lead to decreased Tg and increased structural relaxation and molecular mobility [69]. Recent work has also pointed to the role of solvation free energy in inhibiting API crystallization. Mosquera-Giraldo et al. investigated ASDs formed from telaprevir API and a range of synthesized cellulose polymer derivatives [54]. QM free energies of hydration were computed for cellulose monomers from application of implicit solvation using the polarized continuum model for water. Monomers with more favorable free energies of hydration corresponded to cellulose polymers yielding more stable ASDs with telaprevir, highlighting that polymer ionization is important not only for intermolecular API–polymer interactions but also for solubilizing the carrier itself.
Many of the complexities seen experimentally and via molecular modeling when delineating ASD formation, stability, and dissolution mechanisms can potentially be overcome using statistical learning models. Machine learning (ML), a subset of artificial intelligence (AI), has been at the forefront of major speed and efficiency improvements among in silico drug development approaches over the last decade. By enabling the extraction of complex and often non-linear relationships between the input features and the target feature, ML has been increasingly popular in diverse areas of healthcare and pharmaceutical [116,117,118] and chemical research applications [119,120]. With increasing amounts of data available experimentally, data-driven supervised ML and AI algorithms have shown progress in developing effective drug formulation models [121,122,123]. Even with sparse data, there has been a rise in ML algorithms such as transfer learning, one-shot, zero-shot learning, and Bayesian-based optimization approaches to improve the model performance given small amounts of data [124,125,126].
A subfield of ML is deep learning (DL) [127], exemplified by artificial neural networks (ANNs). ANNs mimic the neural connectivity of the brain, with all nodes (neurons) connected to every other node in the network either directly or indirectly through the layers. ANNs take in information from the input layer, which is then processed through one or more connected hidden layers and finally sent as responses to the output layer. ANNs are especially advantageous for deciphering non-linear and complex unknown relationships between input and output variables. Amongst the many available ML methodologies, ANNs have shown increasing and diverse applications to drug development and process optimization [128,129] and appear in the majority of recent disclosed ML models for predicting and optimizing ASD composition, stability, and dissolution rate (Table 4).
Table 4. Recent ML prediction models for API–carrier ASD properties, stability, and formulation.
Often, many ML algorithms are applied to the same problem, allowing prediction accuracy of the same target trained on identical data to be compared across multiple algorithms, as seen by Han et al. in Table 4 [130]. In addition to ANNs, other ML methods surveyed in this review include genetic algorithms (GA) [131], multiple linear regression (MLR), logistic regression (LR) [132], decision tree (DT) [133], random forest (RF) [134], k-nearest neighbors (kNN) [135], Naïve Bayes (NB) [136], and light gradient boosting machine (LGBM) [137].
There are relatively few examples of statistical learning models predicting properties and phenomena of ASDs formed between small-molecule APIs and polymer carriers, likely due to limited experimental data for model training. ML models tend to focus on predicting target features of amorphous small-molecule compounds. Nevertheless, ML models developed for predicting amorphous API properties provide important understanding and analysis towards rational API–polymer ASD development and formulation [142]. As an example, Nurzyńska et al. developed MLR models for predicting the long-term physical stability of the amorphous forms of 25 poorly water-soluble compounds (no carrier, only the amorphous form of the drug compound) using physicochemical properties directly from two-dimensional structure as input features in addition to measured thermodynamic and kinetic solid properties such as melting point, glass transition temperature, enthalpy of fusion, relaxation time, and configurational free energy [143]. Features that correlate with predicting the amorphous behavior of pure APIs may also show importance in prediction models of API–carrier ASD systems. In another model solely focused on small molecule solubility properties, Przybyłek et al. constructed a model to predict Hansen solubility parameters from a dataset of 130 compounds for which measured solvent solubility parameter data were available [144]. A large collection of connectivity features, indices, and physicochemical properties were generated directly from SMILES (simplified molecular-input line-entry system) data and used as input features to train multivariate adaptive regression splines for solvent solubility parameter prediction. The model was further extended to solvent–polymer miscibility predictions via a binary classification model and predicting solvent-dependent drug-like solid dissolution.
In addition to thermodynamics and molecular modeling approaches, ML has also been applied to evaluating the likelihood of ASD formation between a given API and polymer carrier. Moore et al. developed a LR model to predict the potential of ASD formation using data from twelve small-molecule therapeutics co-solidified with PVP/vinyl acetate copolymer as a stabilizing agent [139]. The target feature, ASD formation potential, was determined by experimental characterization of miscibility and stability, with mixtures that exhibited rapid crystallization or phase separation deemed immiscible. Topological and molecular indices were generated from atom connectivity and three-dimensional structures of the small-molecule compounds. The most significant descriptors identified for predicting dispersion potential within the dataset were the atomic mass-weighted third-order R autocorrelation index and the topological distance (i.e., number of path-connected bonds) between oxygen and chlorine atoms.
The optimal ratio for all components within the ASD is critical for rational formulation design and relates to ASD physical properties, stability, and drug release behavior. Furthermore, ASD forms and subsequent properties are sensitive to small changes in starting material ratios during synthesis and formulation. Predicting such properties and behavior is challenging for binary mixtures and becomes increasingly challenging for systems with even more components. Barmpalexis et al. approached this optimization problem using the statistical learning methodologies of ANNs and GAs. ASDs of nimodipine API with PEG polymer carrier were characterized experimentally, and the ANN and GA models trained to optimize physical formulation and drug release properties of effervescent controlled release floating tablets [138]. The ANN input layer consisted of PEG, PVP, HPMC, effervescent, and nimodipine proportions within the formulation, and an output layer with responses for percent drug released at 60 min, time to reach 90% dissolved, float strength, and float duration. Inclusion of trigonometric and exponential/logarithmic functions when deriving the GA equation accounted for complex non-linearity. Input feature analysis indicated the proportion of effervescent agents as most important, followed by percent HMPC, while PEG and PVP proportions as features showed minor importance in the model. Applying a similar approach, Medarević et al. constructed an ANN model to predict the dissolution rate of carbamazepine API within ternary ASDs of SOL and P188 carriers [141]. Carbamazepine release rate was found to correspond most strongly with ASD compositions lowest in carbamazepine and highest in P188. ASD dissolution rate was also modeled by Barmpalexis et al. by constructing an ANN to predict the rate of dissolution (measured as percent API dissolved after thirty minutes) of ASDs formed from mixtures of tibolone and different molecular weights of PEG polymers [140]. Dissolution rate was determined to primarily correlate with percent tibolone in the mixture, molecular weight of the PEG polymer, and the mixing temperature.
The long-term physical stability of ASD formulations remains a significant challenge, and one that is often not sufficiently addressed by solubility and interaction parameters alone [110]. Given specific data related to ASD stability over time, machine learning models can be trained to directly predict ASD stability and dissolution. In a recent study, Han et al. combined molecular modeling techniques and eight different machine learning algorithms towards the development of models for ASD physical stability prediction [130]. The collected ASD formulation data encompassed formulation, process parameters, stability measurement conditions, and experimental stability measurements (“stable” vs. “unstable” labels after both three and six months). The random forest (RF) algorithm gave the highest accuracy models, and RF feature importance analysis determined drug loading ratio, relative humidity, temperature (storage and preparation), and molecular weight (polymer carrier and API) as key features for the prediction of three- and six-month physical ASD stability.
Formulating poorly water-soluble drugs as ASDs with polymeric carriers is an effective method for increasing solubility yet is often hindered by complex characterization experiments and lengthy development cycles. It is likely that the representation of approved ASDs coming to market will remain relatively low until significant advances are achieved in predicting compatible carriers and excipients for novel therapeutic compounds with poor aqueous solubility. Much attention has been given to the use of solubility and FH interaction parameters for estimating API–polymer thermodynamic miscibility. Molecular modeling via MM, QM, and MD approaches for studying API–polymer miscibility and intermolecular non-bonding interactions, in addition to investigating mechanisms of ASD formation, physical stability, dissolution, and drug release, have increased alongside advances in computational power. These studies have provided a breadth of information but are often limited to a small number of APIs and carriers. More large-scale studies with quantitative comparisons between predicted and experimental results will assist in theoretical progress, in which a remaining major challenge is the accurate calculation of the entropic contribution to miscibility. Furthermore, combined DFT–MD methods providing an electronic description of molecular systems remain highly underexplored in studies of ASD API–polymer interactions and may also yield improved model predictions. Finally, even today’s computational power remains a limiting factor in the simulation of important ASD kinetic behavior such as molecular mobility and mechanisms of phase separation and crystallization over desired storage and physiologically relevant timescales, all of which dictate long-term ASD physical stability. These limitations must be overcome prior to ASD computational modeling tools arriving as more mainstay applications in process chemistry and pharmaceutical settings.
ML methods have the potential to provide the next transformative leap forward toward rapid polymer screening and formulation design, with recent examples already showcasing the ability of ML to predict API–polymer miscibility, optimize ASD formulation for desired physical properties, and predict dissolution rates. Of special note for future applications and expanded research is close integration of molecular modeling and machine learning pipelines [145]. Prediction accuracy and model robustness will only improve over time as more and more data become available, where data encompass both experimental targets and more representative chemical and physical input descriptors of ASD molecular systems. Many experimental and molecular modeling studies test less than ten polymers for a given API, with a fraction of the combinations yielding promising results, yet an advantage of ML is that negative outcomes remain valuable as the training set should avoid skew towards only positive experimental results. In this sense, almost all characterization data and excipient screening results are pertinent to ML model development, and the limiting factors become data organization, labeling, and public accessibility.
The overarching driving force of ASD modeling is ultimately to bring more poorly soluble therapeutics to market in order to rapidly combat much needed disease areas. Conventional approaches to excipient screening and ASD development have relied on trial and error experiments, both to match a given API with a compatible polymer and to further optimize the formulation properties once a suitable match is located. There remains only around thirty commercial ASD formulation products available on the market (Table 1). Experimentally, it is intractable to survey even half of all available polymers and excipients, with potential synthesized novel carriers expanding the chemical space further. This bottleneck necessitates that computational approaches continue to play a vital role. Advances in predictive ASD models will facilitate more widespread rational excipient screening, promote de novo excipient construction, and accelerate ASD development, thus expediting the progression from early-stage assessment and characterization to stable formulations with optimal physical properties primed for preclinical and clinical testing.
Writing—original draft preparation, D.M.W. and Y.B.; writing—review and editing, A.J.; writing—review and editing, project administration, V.A. and K.C.; writing—review and editing, project administration, funding acquisition, J.V. All authors have read and agreed to the published version of the manuscript.
This research was funded by VeriSIM Life: Inc.
The authors declare no conflict of interest.
Learn more about VeriSIM Life’s BIOiSIM platform and unique Translational Index™️ technology.
AtlasGEN Novel Drug Designer efficiently identifies new molecules based on favorable chemical properties and simultaneously evaluates them for clinical safety and effectiveness based on biology using VeriSIM Life’s groundbreaking Translational Index™ technology.
Iterative ranking of target engaging hits by their Translational Index values reduces the time and cost associated with other computational hit discovery approaches, while compressing the hit-to-lead refinement and optimization process by pre-validating hits for translatability. This unique integration dramatically reduces the number of compounds to evaluate in subsequent experimental research.
The ability to generate target engaging compound hits depends largely on the size of the molecular search library. The AtlasGEN search space is supported by more than 1 trillion compounds. Generative AI and deep learning techniques deliver unprecedented discovery capability.
AtlasGEN accurately predicts protein-lygand binding affinity by combining geometric conformer analysis with machine learning based on experimental data 72% more efficiently than molecular docking-based approaches.
.png)
