Multi-Phase Solvation Model for Biological Membranes: Molecular

May 29, 2017 - Video: Taking a Closer Look at Crystal Engineering of Self-Assembled Porous Protein Materials in Living Cells. Crystalline porous mater...
0 downloads 8 Views 2MB Size
Article pubs.acs.org/JCTC

Multi-Phase Solvation Model for Biological Membranes: Molecular Action Mechanism of Amphotericin B J. M. Falcón-González,†,‡ G. Jiménez-Domínguez,† I. Ortega-Blake,∥ and M. Carrillo-Tripp*,†,§ †

Laboratorio de la Diversidad Biomolecular, Centro de Investigación y de Estudios Avanzados Unidad Monterrey, Vía del Conocimiento 201, Parque PIIT, C.P. 66600, Apodaca, Nuevo León, México ‡ Unidad Profesional Interdisciplinaria de Ingeniería Campus Guanajuato, Instituto Politécnico Nacional, Av. Mineral de Valenciana No. 200, Col. Fraccionamiento Industrial Puerto Interior, C.P. 36275, Silao de la Victoria, Guanajuato, México ∥ Instituto de Ciencias Físicas, Universidad Nacional Autónoma de México, Apartado Postal 48-3, C.P. 62251, Cuernavaca, Morelos, México ABSTRACT: Amphotericin B (AmB) is still the most effective drug for the treatment of systemic fungal infections in humans. Despite significant theoretical and experimental efforts trying to understand its molecular mechanism of action, the answer has remained elusive. In this work, we present a computational methodology to test the current membrane related hypotheses, namely, transmembrane ion channel, adsorption, and sterol sponge. We use a thermodynamic approach in which we represent the membrane by a multiphase solvation model with atomic detail (MMPSM) and calculate the free energy of transferring the drug between phases with different dielectric properties. Furthermore, we compare AmB to a chemical analogue with increased safety, an L-histidine methyl ester of AmB. Our findings reveal that both drugs dimerize in all solvents studied here. Also, it is energetically unfavorable for the drugs to penetrate into the hydrophobic core of the membrane, unless their concentration is high. Finally, it is thermodynamically possible that the sterols migrate from the membrane into a drug droplet adsorbed at the surface of the bilayer. In light of our results, several effects could take place in the complex antibiotic process. We suggest a molecular mechanism that connects all three hypotheses through a drug concentration dependence and propose that the drug promotes the formation of membrane toroidal pores. Because MMPSM is of general interest, we made it available at http://tripplab.com/tools/mmpsm. causing an abnormal increase in permeability, leaking of K+ ions and other small organic compounds vital for cell function, leading to eventual cell death.10,11,13,19,27−29 It has also been shown that AmB presents a higher affinity for ergosterol than cholesterol, the main sterols in fungal or mammalian cells, respectively.18,30 It is believed that the AmB preference for ergosterol-containing membranes is related to the stability of the putative channel31 or the self-association of the drug.14 Electrical conductivity is detected in ergosterol-containing membranes at low to high drug concentrations, indicating that AmB has an effect on the bilayer, either in a monomer or aggregated state.14,18 On the other hand, ion permeability in cholesterol-containing membranes is only seen at high AmB concentrations.7,12,14,15 Following such results, it is believed that the low solubility of AmB could be related to its toxicity and side effects on mammals.17 Different approaches have been proposed to reduce the collateral toxicity of AmB. For example, the drug has been incorporated into liposomes, chemically conjugated to

1. INTRODUCTION Amphotericin B (AmB) is the most effective drug for the treatment of systemic fungal infections in humans.1−6 The AmB molecule is amphoteric, owing to the presence of a carboxyl and an amino group, and amphiphilic, due to both nonpolar conjugated heptaene on one side of the macrolactone ring and a polar hydroxyl chain on the other side. Hence, AmB is poorly soluble in water at neutral pH, and even less in salt solutions, making its administration difficult.7 However, AmB has high effectiveness, broad antifungal spectrum, and low occurrence of resistant clinical strains.8 Unfortunately, AmB also acts against mammalian cells, being one of the most toxic drugs used in clinics.9 In the past decades, theoretical and experimental research groups have focused their efforts on understanding the mechanisms of action of AmB.10−26 The design of a new drug with the same efficiency as AmB but higher solubility and lower toxicity than AmB has been one of the main goals for many years.16 Several ideas on how the drug works have been proposed, but despite ample contradictory evidence the “standard” hypothesis has prevailed. This hypothesis states that AmB forms transmembrane barrel-stave ion channels, © 2017 American Chemical Society

Received: March 29, 2017 Published: May 29, 2017 3388

DOI: 10.1021/acs.jctc.7b00337 J. Chem. Theory Comput. 2017, 13, 3388−3397

Article

Journal of Chemical Theory and Computation polymers, or mixed with lipids (see Barrat et al.32 and references therein). More recently, a new scheme involving nanosuspensions of AmB and nanoparticles has been developed.33 Nowadays, the solubilization method of choice consists of the introduction of sodium deoxycholate detergent into the suspension of AmB molecules to form a colloidal system.34 Hence, slow perfusion is the conventional route for AmB application.3 None of these methods have been able to remove the AmB toxicity and maintain the effectiveness of the drug.35 As an alternative approach, chemical substitutions have been performed producing an improved selectivity.14,31,36,37 Although the physicochemical properties of AmB are known and its structure and behavior have been explored by chemical and biophysical studies,11,13,29 the mechanism of action of the drug has not been fully elucidated. In fact, two ideas have recently emerged, in remarkable contrast to the ion channel “standard” hypothesis (H1). The first alternative suggests that AmB kills fungal cells by binding ergosterol on the surface of the membrane and that channel formation is a subsequent process.24 Isothermal titration calorimetry, potassium efflux, and minimum inhibitory concentration experiments support the observation of such ergosterol-binding, membrane permeabilization, and antifungal activities. In this scheme, known as the surface adsorption hypothesis (H2), AmB stays in the polar headgroup region of the membrane, followed by sequestration of ergosterol to the bilayer surface. The second alternative, known as the sterol sponge hypothesis (H3), proposes that AmB is present as large aggregates outside the membrane. There is evidence that such droplets can extract the sterol from inside the bilayer.38 Here, AmB does not penetrate into the membrane. The sterol molecules move to the membrane surface attracted by the drug and subsequently removed from the membrane through a strong interaction. The sponge hypothesis is supported by solid-state NMR, transmission electron microscopy, and cell-based experiments. Previously, we showed that the stabilization and distribution of lateral stress within the hydrophobic core of the membrane is sensitively dependent on the presence of cholesterol.39 Removing cholesterol from the membrane shifts the repulsive lateral chain pressure away from the bilayer interior toward the lipid−water interface. Such observation should remain true for other types of sterols. This may affect conformational equilibrium, and proper function, of integral membrane proteins. The mycosamine moiety in the AmB molecule40 plays a predominant role in the binding of sterols in both scenarios, which would be strictly required for channel formation and antifungal activities.23,24,38 Recently, we have proposed an AmB chemical derivative with increased safety (lower colateral toxicity in therapeutical use).26 Based on the idea that selectivity is also related to the membrane structure,41 an amide substitution from a hydroxil group of AmB was performed in an attempt to force the mycosamine ring to go into the membrane. Specifically, an (L)-histidine methyl ester substitute leads to the AmB analogue, known as A21 (Figure 1). The chemical derivative is consistent with the notion that selectivity of AmB to sterols could be modulated by the modification of this chemical group,31,36 with the mycosamine ring unchanged. We analyzed the AmB and its analogue A21 through molecular thermodynamics, electrophysiological, pharmacological, and preclinical studies. The results show that A21 has a reduced activity against cholesterol-containing membranes, hence, a significant toxicity decrease toward mammals,

Figure 1. Schematic representation of molecules studied in this work: (A) Amphotericin B, (B) AmB chemical analogue A21, (C) DMSO, (D) octanol, (E) AMP, (F) hexadecane, (G) cholesterol, (H) ergosterol, and (I) the two membrane models, M1 and M2, using the solvation approximation MMPSM.

but maintains its effectiveness against fungal cells. In fact, the A21 derivative increases its selectivity toward ergosterol containing membranes vs cholesterol containing ones by 4.33-fold and toward micotic cells vs mammal cells by 7-fold. An additional advantage of A21 over AmB is its increased solubility in water. A21 aggregation starts to occur at a higher concentration than AmB in aqueous solution. The most favorable dimerization free energy corresponds to configurations with an antiparallel electric dipole for both drugs. Analysis of the molecular structure suggests that A21 might form weaker interactions because the polyene rings are perturbed by the presence of the L-histidine that replaces the AmB hydroxyl group. In this work, we test the currently most accepted hypotheses of the mechanism of action of the AmB antimycotic, previously reviewed.25 We take into account the A21 analogue to do a comparative analysis. We use a thermodynamic approach in which we represent the membrane using a two-phase or a three-phase model with atomic detail. Under such theoretical framework, we estimate the solvation and transfer free energy of various compounds to different solvents, mimicking the pass of the drugs through the membrane. 3389

DOI: 10.1021/acs.jctc.7b00337 J. Chem. Theory Comput. 2017, 13, 3388−3397

Article

Journal of Chemical Theory and Computation

2. METHODS 2.1. Molecular Systems. Octanol has been widely used as a rough approximation to model the lipophilic membrane and is still the most common solvent used to evaluate a drug’s partition coefficient in relation to the aqueous solution.42,43 Despite this, octanol cannot appropriately model the hydrophilic and hydrophobic regions of a phospholipid bilayer as two independent phases.42,44 Hexadecane has been used as an alternative molecule to model the membrane interior adequately, given its dynamic, electrostatic, and viscosity properties.45,46 However, neither approximation takes into account the existence of a well-defined polar zone present at the membrane’s surface that directly interacts with the aqueous environment. Hence, it is necessary to identify a solvent that mimics the polar nature of the membrane. Here, we propose the use of acetyl-methyl-phosphate (AMP) since its chemical structure highly resembles that of the phospholipid headgroups (Figure 1). In this work, we use two membrane models to estimate the effect of different approximations on the thermodynamic analysis (Figure 1). The first model (M1) consists only of octanol. The second model (M2) explicitly differentiates the membrane interior, represented by hexadecane, from the membrane surface, represented by AMP. As a whole, we call this the MMPSM approach. It is well known that AmB solubility in water is low (leading to difficulties for its oral administration), moderately soluble in alcohols, and soluble in organic solvents such as methanol, ethanol, and dimethyl sulfoxide (DMSO). Absorption and circular dichroism results reveal that AmB remains in a monomeric form when dissolved in DMSO.7,12,15 We included DMSO in this study because it is the solvent of choice for doing experimental studies,47 in particular those that try to understand the AmB aggregation process. The molecular drug aggregation in solution is a critical step in the mechanism of action of the AmB and the A21. Here, we analyzed the effect of dimerization on the solvation free energy. We can consider four dimer molecular configurations, namely, head-to-head (HH) or head-to-tail (HT), according to the relative arrangement of the polar head and the tail of the molecules, having parallel (P) or antiparallel (A) orientations of the electric dipole (Figure 2). Even though all configurations can occur naturally (drug concentration dependent process), we previously showed that HTA is the most likely configuration found in water for both AmB and A21, with the lowest energy,26 in accordance to previous quantum chemical calculations.48 In fact, the HTA configuration is maintained in other solvents as well, as we show in the following sections. Therefore, we only used the HTA configuration for the solvation free energy calculations of AmB and A21 dimers. Thus, we studied the solvation thermodynamics of six compounds, namely, the AmB in monomer or HTA dimer state, the A21 in monomer or HTA dimer state, the cholesterol in monomer state, and the ergosterol in monomer state. The solvents we used were water, octanol, hexadecane, and AMP. To assess the validity of our theoretical approach, we included DMSO to solvate AmB and A21. Figure 1 shows the chemical structure of all these molecules. Furthermore, we also considered the AmB and A21 to be sterol solvents, following the idea proposed by the sponge hypothesis (H3).38 2.2. Molecular Dynamics. We built the initial configuration of the molecular systems by placing each compound at the center of an ∼8 × 8 × 8 nm3 cubic box filled up with a

Figure 2. Four possible dimer configurations of AmB, namely, headto-head (HH) or head-to-tail (HT), according to the relative arrangement of the polar head and the tail of the molecules, with parallel (P) or antiparallel (A) orientations of the electric dipole (depicted with an arrow). The long axis of the molecule is shown as a dotted line. θ is the angle between the long axis of each molecule in the dimer. The value of cos(θ) equals 1 for the HH configuration and equals −1 for the HT configuration. On average, the linear dimensions of a dimer are 1 × 1 × 3 nm.

solvent. We use the Packmol package49 to obtain the desired density in each case. In this way, the numbers of molecules used to solvate each drug compound were 16 900 (water), 4200 (DMSO), 1970 (octanol), 1190 (hexadecane), and 1150 (AMP). The numbers of molecules used for solvating the sterol monomers were 24 144 (water), 1995 (octanol), 1198 (hexadecane), 1145 (AMP), 504 (AmB), and 347 (A21). All the force fields we used are a united-atom approximation, fully consistent with the GROMOS53A6 set of parameters.50 In the case of the AmB and A21 molecules, a detailed description is found in ref 26. For ergosterol, cholesterol, hexadecane, and AMP, the force fields were taken “as is” from the Automated force field Topology Builder web accessible server (ATB).51 The ATB uses a multistep process in which results from quantum mechanical calculations (B3LYP/6-31G**, in this case) are combined with a knowledge-based approach to ensure compatibility with the parameter set of the specified GROMOS force field (53A6, in this case). DMSO and octanol were represented by the set of parameters already included in the GROMOS53A6 force field. The simple point charge (SPC) model was used to describe the water molecules.52 We performed a series of molecular dynamics (MD) simulations using GROMACS 5.0.53 Energy minimization was accomplished applying the low memory Broyden−Fletcher− Goldfarb−Shanno approximation,54 followed by a steep algorithm53 to achieve greater convergence toward the energy minimum. System equilibration was performed in the canonical ensemble (NVT) for 0.5 ns, followed by the isothermal− isobaric ensemble (NPT) for another 0.5 ns. The Berendsen barostat55 was used to keep the pressure at 1 bar and compressibility at 5 × 10−5 bar−1, with a coupling time of 0.5 ps. The temperature was kept constant at 298 K using the Berendsen thermostat,56 with a coupling time of 1 ps. Longrange electrostatic interactions were handled with the particlemesh Ewald method57 with a cutoff of 1 nm. van der Waals 3390

DOI: 10.1021/acs.jctc.7b00337 J. Chem. Theory Comput. 2017, 13, 3388−3397

Article

Journal of Chemical Theory and Computation

them. ΔG solv is therefore calculated according to the thermodynamic cycle shown in Figure 3. Here, we follow the

interactions were considered within a spherical cutoff radius of 1.4 nm. Covalent bonds involving hydrogen atoms were constrained using the LINCS algorithm.58 A time step of Δt = 0.002 ps was set to integrate the equations of motion. Such protocol ensures a rigorous study of biological systems in thermodynamic equilibrium. At the end, all the calculations amounted to a total of 4400 ns of MD simulations. On average, 10 ns/day were produced at one 24-core node on a high performance computing cluster. This is equivalent to more than 250 000 CPU hours. The equilibrated systems of pure solvents needed in the MMPSM framework can be downloaded at http://tripplab.com/tools/ mmpsm. 2.3. Free Energy. Free energy calculations of a complex molecular system are practically impossible to achieve due to the large number of degrees of freedom present therein. However, it is possible to calculate the difference in free energy between two related states of a system, A and B, using a coupling parameter, λ, in the Hamiltonian of each system,59 represented by H(r, p; λ), where r and p are the respectively coordinates and momenta of the N particles of the system. In this way, the partition function Z can be writen as Z(λ) = (h3N N! )

∫ ∫ e−H(r,p;λ)/k T dr dp

Figure 3. Schematic representation of the thermodynamic cycle used on the MMPSM approximation.

methodology described in ref 63. In short, ΔG1 is the work required to remove all the intramolecular nonbonded interactions of the solute in vacuo; i.e., the Coulomb and Lennard-Jones parameters (in two steps, in that order) are gradually set to zero (state B, a “ghost molecule composed of dummy atoms’”). The intramolecular bonded parameters and masses remain unchanged. ΔG2 is the work required to transfer the ghost molecule from vacuo to the solvent. Because the ghost does not interact with the system and the available volume is the same, ΔG2 is zero. ΔG3 is the work required to remove the intermolecular interactions between the molecule and the solvent, and the molecule’s intramolecular nonbonded interactions. Therefore,

B

(1)

and the free energy G as G(λ) = −kBT ln Z(λ)

(2)

ΔGsolv = ΔG1 − ΔG3

Both are also λ-dependent functions, where h is the Planck constant, kB is the Boltzmann constant, and T is the absolute temperature of the system. The difference in free energy between state A and state B can be expressed as ΔGBA ≡ G(λB) − G(λA) =

∫λ

λB

G′(λ) dλ

A

The logic followed to estimate the free energy cost to transfer a solute from one solvent to another is similar. However, one of the solvents is taken instead of vacuo. Now, ΔG1 is the work required to create a ghost molecule in the first solvent (gradually removing all intra- and intermolecular nonbonded interactions), ΔG2 is the work required to transfer the ghost molecule from the first to the second solvent (zero), and ΔG3 is the work required to create the ghost molecule in the second solvent. Hence, ΔGtransf = ΔG1 − ΔG3. As mentioned earlier, the removal of the nonbonded interactions in the thermodynamic cycle is performed gradually, passing the molecule from its initial state (A) to a second state in which the molecule has totally disappeared (ghost or state B). Here, we first switched off the Coulomb interactions and then the van der Waals interactions, using a soft-core potential to avoid singularities and large fluctuations in energy, as previously suggested.61 Technical details on how to implement this in GROMACS can be found in ref 53. Twenty-five λ-values were sampled between 0 (state A) and 1 (state B). The free energy change associated with each λpoint was calculated using the TI procedure. At each λ-value the systems were minimized and equilibrated, and the average of the observables were taken after 5 ns of simulation. In this stage, we used the V-rescale thermostat62 to keep the temperature at 298 K in all the cases, except those in vacuo for which the Nosé−Hoover thermostat was used.63 All bonds were constrained using the LINCS algorithm.58 The graph

(3)

where G′(λ) ≡ dG/dλ. Differentiating eq 2 with respect to λ, we obtain ∂H(r, p; λ) −H(r, p; λ)/ kBT e ∂(λ)

∫∫ G′(λ) =

∫∫ e

−H(r, p; λ)/ kBT

dr dp

dr dp

which gives the ensemble average of

∂H(r, p; λ) . ∂(λ)

(4)

Therefore, the

free energy given by eq 3 can be rewritten as ΔGBA =

∫λ

λB

A

∂H(λ) ∂λ

dλ λ

(6)

(5)

an expression known as the thermodynamic integration (TI) formula.60 The λ-dependence of the Hamiltonian defines a pathway that connects two states of the system denoted by A and B. One approach to solve eq 5 is to evaluate the ensemble average n times, where n is the number of chosen λ’s, by performing separate simulations for each λ value. Then, the integral can be numerically calculated. 2.3.1. Solvation and Transfer Free Energy. The solvation free energy, ΔGsolv, is the work required to transfer a molecule (solute) from vacuo into a solution (solvent). However, calculations of such physical property are not practical to realize, so one has to deal with thermodynamic cycles used for computational efficiency reasons.59 Since the free energy is a state function, the difference in free energy between two states of a system is independent of the path used to go between

connecting the

∂H(λ) ∂λ λ

values was integrated using the

trapezoidal rule to obtain ΔGBA given by Equation 5. The error in the average was estimated by the block averaging procedure.64 The individual errors were also integrated using the trapezoidal method from λ = 0 to λ = 1 to obtain an estimated error for ΔGBA. 3391

DOI: 10.1021/acs.jctc.7b00337 J. Chem. Theory Comput. 2017, 13, 3388−3397

Article

Journal of Chemical Theory and Computation 2.3.2. Partition Coefficient. In general, the partition coefficient, or partition constant, is a parameter that indicates the hydro- or lipophilicity of a molecule, i.e., it provides information on the preference of the solute for either a polar or organic solvent. The partition coefficient has a great significance in the pharmaceutical industry because it gives a measure of the amount of drug that enters, leaves, or moves through a lipid membrane, and it is directly correlated with the drug effectiveness or toxicity. Such a parameter has been traditionally measured at the equilibrium of a two-phase system such as an octanol−water mixture. A solute with a positive octanol/water partition coefficient is mainly distributed to hydrophobic phases, such as lipid bilayers of cells. Conversely, solutes with negative octanol/water partition coefficients are found primarily in the aqueous phase. Knowing the transfer free energy values, the partition coefficient, log P, can be computed according to the relation, ΔGtransf = −2.303RT log P

(7)

where R is the gas constant and T the absolute temperature.

3. RESULTS AmB aggregation in solution has been the limiting step for its pharmacological use. Here, we made a thermodynamic study of AmB and its chemical derivative, A21, on different environments. We estimated the solvation free energy, ΔGsolv, of the drugs in monomer or dimer configurations immersed in different solutions. As we previously showed, the HTA is the most stable dimer configuration found in water for both drugs.26 However, it is not necessarily true that the HTA should be the most stable dimer configuration found in a solvent other than water, in particular, those we use here to represent the membrane. We performed a molecular dynamics analysis of the AmB and A21 dimer in DMSO, octanol (M1), AMP, and hexadecane (M2). Each simulation started with a dimer in the HTA configuration, running for 50 ns, i.e., long enough to appreciate any relative rearrangement that could occur in the initial configuration. We analyzed the structure of the molecular arrangement of the dimer over time by monitoring the value of cos(θ), where θ is the angle formed between the long axis of each molecule in the dimer (Figure 4). We found that the HTA dimer configuration for AmB and A21 is stable in all solvents but DMSO. This result is concomitant with the known fact that AmB favors the monomeric configuration when solubilized in DMSO.36 In the case of the AMP solvent, there is a structural equilibration within the first 15 ns for both drugs, where a rapid variation of the value of cos(θ) is observed. This is due to a transient dissolution of the dimer; the initial configuration was not the optimal, and the distance between the center of mass of the monomers increases enough for the molecules to rotate. This behavior is more pronounced for A21 because of the larger (L)-histidine methyl ester substitute. After this, the monomers come closer and the dimer stabilizes in the head−tail configuration. As expected, both AmB and A21 have lower ΔGsolv values in DMSO than water (Table 1). Hence, both drugs would migrate from water to DMSO, as shown by the favorable transfer free energy value, ΔGtransf (Water → DMSO). Noteworthy, the drug’s dimerization makes the solvation energy of AmB and A21 significantly more negative in all solvents (Table 2). The partition coefficient values also reflect a greater solvation of the drugs in DMSO than water and its increase by dimerization.

Figure 4. Structural analysis of the dimer configuration for (A) AmB and (B) A21, in different solvents: water (black), DMSO (red), octanol (green), hexadecane (blue), and AMP (brown). The value of cos(θ) is −1 in the HTA configuration.

Table 1. Solvation (ΔGsolv) and Transfer (ΔGtransf) Free Energy Calculations for Monomers AmB and A21, Expressed in kcal/mola monomers

AMP DMSO Octanol Water Hexadecane

AmB

A21

ΔGsolv

ΔGsolv

−115.2 −88.4 −83.7 −75.6 −42.0

± ± ± ± ±

−112.6 −76.8 −65.0 −53.4 −40.6

6.1 5.9 8.1 5.2 5.9

± ± ± ± ±

4.9 4.5 6.9 4.3 4.6

Monomers AmB ΔGtransf Water → DMSO Water → Octanol Water → AMP AMP → Hexadecane Water → Hexadecane

−12.8 −8.1 −39.6 73.2 33.6

± ± ± ± ±

4.1 6.3 4.3 4.9 4.0

A21 log P 9.4 5.9 29.1 −53.8 −24.7

ΔGtransf −23.4 −11.6 −59.2 71.9 12.7

± ± ± ± ±

4.0 6.3 4.3 4.6 4.0

log P 17.2 8.5 43.5 −52.8 −9.4

a

Partition coefficients, log P, are also shown for all the systems (dimensionless), calculated directly from the transfer free energies as indicated in eq 7.

Next, we followed the idea of testing the proposed hypotheses of the mechanism of action of the drugs (H1−3) from a thermodynamic standpoint using two membrane models. We estimated the solvation free energy of AmB and A21, in a monomer or dimer configuration, in octanol (model M1), or AMP and hexadecane (model M2). The transfer free energy represents the probability that a solute spontaneously 3392

DOI: 10.1021/acs.jctc.7b00337 J. Chem. Theory Comput. 2017, 13, 3388−3397

Article

Journal of Chemical Theory and Computation

water is a poor solvent for both sterols (Table 3). In fact, models M1 and M2 correctly predict that the sterols will stay

Table 2. Solvation (ΔGsolv) and Transfer (ΔGtransf) Free Energy Calculations for Dimers AmB and A21, Expressed in kcal/mola

Table 3. Solvation (ΔGsolv) and Transfer (ΔGtransf) Free Energy Calculations for Cholesterol and Ergosterol Expressed in kcal/mola

dimers AmB

A21

ΔGsolv

ΔGsolv

−192.3 −142.9 −137.4 −111.1 −76.0

AMP DMSO Octanol Water Hexadecane

± ± ± ± ±

−212.6 −131.6 −121.0 −83.6 −83.4

13.4 15.4 23.4 11.9 12.9

± ± ± ± ±

cholesterol

10.6 9.4 16.6 9.8 12.0

AmB A21 AMP Octanol Hexadecane Water

dimers AmB ΔGtransf Water → DMSO Water → Octanol Water → AMP AMP → Hexadecane Water → Hexadecane

−31.8 −26.3 −81.2 116.3

± ± ± ±

12.1 20.1 10.2 11.2

35.1 ± 9.6

A21 log P 23.4 19.3 59.6 −85.4 −25.8

ΔGtransf −48.0 −37.5 −129.0 129.2

± ± ± ±

log P

8.1 15.2 9.3 11.4

35.3 27.5 94.7 −94.9

0.2 ± 10.7

−0.1

ΔGsolv

−36.8 ± 2.1 −37.6 ± 3.3 −23.6 ± 1.1 −19.9 ± 2.5 −17.7 ± 0.6 −2.5 ± 0.7 cholesterol

−39.1 ± 2.3 −33.8 ± 3.6 −24.1 ± 0.9 −20.1 ± 1.8 −17.9 ± 1.2 −2.9 ± 0.6 ergosterol

ΔGtransf Octanol → Water Octanol → AmB Octanol → A21 Hexadecane → Water Hexadecane → AMP Hexadecane → AmB Hexadecane → A21 AMP → Water AMP → AmB AMP → A21 AmB → Water A21 → Water

a Partition coefficients, log P, are also shown for all the systems (dimensionless), calculated directly from the transfer free energies as indicated in eq 7.

migrates from one solvent to another. In the approximation of model M1, our results indicate that both drugs would go from the aqueous solution into the membrane, i.e., water → octanol (Table 1). The partition coefficient values also reflect the lipophilic nature of the drugs’ molecular structure, increasing when they dimerize (Table 2). In other words, the drugs will aggregate, and their concentration will be greater in the nonaqueous phase. As mentioned earlier, model M1 is a rough approximation to represent the membrane. Here, we propose a better membrane representation, M2, in which we use AMP to mimic the polar surface region formed by the phospholipids’ headgroups and hexadecane as the low dielectric medium found in the interior of the membrane. M2 predicts that the drugs would migrate from the aqueous solution into the polar phase (Water → AMP) but not go into the membrane’s interior (AMP Hexadecane). Again, we see that the dimerization of the drugs significantly increases the favorable change in the transfer energy. In this scheme, the drug will aggregate and adsorb to the membrane’s surface. Interestingly, we found that the ΔGtransf and the log P values are practically zero in the case of the A21 dimer when going from the aqueous solution to the nonpolar membrane’s core (Water → Hexadecane). This result indicates that this particular molecule’s configuration would be found in the same proportion in a biphase system composed of such solvents. However, we have to keep in mind that in a real system, i.e., the membrane, the drug would have to go through a polar phase, in our case, the AMP. According to the sponge hypothesis (H3), once the drug is adsorbed on the surface of the membrane, the interaction with the sterols would be strong enough to drive the extraction process. We investigate the thermodynamic viability for this to occur by estimating the ΔGsolv values of cholesterol and ergosterol in water, octanol (M1), AMP, and hexadecane (M2) and AmB or A21 (H3). As expected, our results show that

ergosterol

ΔGsolv

17.3 −17.0 −17.7 15.2 −5.9 −19.1 −19.9 21.1 −12.6 −13.4 34.3 35.1

± ± ± ± ± ± ± ± ± ± ± ±

2.5 3.9 5.1 0.6 1.0 2.0 3.2 1.1 2.6 3.8 2.1 3.3

log P −12.7 12.5 13.0 −11.2 4.3 14.1 14.6 −15.5 9.3 9.8 −25.2 −25.8

ΔGtransf 20.1 −19.0 −13.7 15.0 −6.2 −21.2 −15.9 21.2 −14.6 −9.3 36.2 30.9

± ± ± ± ± ± ± ± ± ± ± ±

1.8 3.6 4.9 1.4 1.6 3.0 4.3 1.1 2.7 4.1 2.5 3.8

log P −14.8 14.0 10.1 −11.0 4.5 15.5 11.6 −15.6 10.7 6.8 −26.6 −22.7

a Partition coefficients, log P, are also shown for both sterols (dimensionless), calculated directly from the transfer free energies as indicated in eq 7.

inside the membrane and will not migrate to the aqueous solution (Octanol Water -M1-, Hexadecane Water and AMP Water -M2-). Remarkably, both drugs AmB and A21 are much better solvents for cholesterol and ergosterol than any of the other organic solvents studied here. In this scheme, when a droplet of a drug is adsorbed at the membrane’s surface, the extraction ΔGtransf is favorable (Octanol → Drug -M1-, Hexadecane → Drug and AMP → Drug -M2-), and the sterols should migrate into it. Once absorbed by the droplet, the sterols would remain inside it, as predicted by the ΔGtransf unfavorable values of going from the drug’s environment to the aqueous solution (Drug Water).

4. DISCUSSION It is well accepted that the antibiotic molecular mechanism of action of AmB involves the interaction of the drug with the membrane. Previous theoretical studies make the assumption that AmB can readily penetrate the membrane. Such studies have performed structural and energy analyses with the drug molecules already inside the hydrophobic core of a lipid bilayer (e.g., ref 19). To the best of our knowledge, there is no clear evidence on the drug’s permeation molecular mechanism. Computational methods have been increasingly used to estimate thermodynamic properties of biological systems. In general, they are based on either explicit or implicit solvents. Even though the latter approach is appealing, we might be far from having a robust description of implicit solvent models for the calculation of the solvation free energy of biomolecules in 3393

DOI: 10.1021/acs.jctc.7b00337 J. Chem. Theory Comput. 2017, 13, 3388−3397

Article

Journal of Chemical Theory and Computation

mechanism has been identified in other biological systems. For example, in a similar way, a two-state model has been proposed to explain the antimicrobial action of some peptides (AMPs), where the susceptibility of a cell depends on a critical peptide-to-lipid ratio.72 Experimental evidence shows that AMPs adsorb to the surface of the membrane and disrupt the cholesterol distribution through strong interactions.73 In addition, theoretical studies show that AMPs can also disrupt the structure of the bilayer, even in low concentrations, and form toroidal pores, where peptide aggregation is a prerequisite to pore formation.74 Hence, we speculate that there is a high similarity in the molecular mechanism of action between the antimicotic AmB and antimicrobial peptides. We recognize MMPSM limitations to come from two main sources. Given the minimalistic abstraction used to represent the membrane, the model cannot reproduce any mechanical properties intrinsic to a lipid bilayer. Also, because of the approximations of the current force fields, this approach is restricted to give a qualitative rather than quantitative view in most cases. However, the energy estimations will become more accurate and the error will decrease, as further refinements are made to the molecular force fields. Within these limitations, MMPSM proved to be consistent with previous observations. For instance, we have shown that an increase in drug solubility (reduced dimerization) could lead to an increase in selectivity.26 However, solubility is a kinetic process, distinct from solvation. Hence, the solvation model presented here cannot explain the observed increase in selectivity of A21 over AmB. Nonetheless, MMPSM correctly predicts the ergosterol over cholesterol preference of AmB (solvation free energy shown in Table 3) previously seen experimentally23 and theoretically.22 Furthermore, likely shortcomings of our method should be taken into account regarding the thickness and heterogeneous composition of the different phases in a membrane. Density profiles show that the phospholipid heads cover around 1.5 nm in depth, with some amount of water present. The hydrophobic core of the membrane is around 2 nm deep, with no water present.39,75 MMPSM should be a good approximation when dealing with molecules (solutes) with a size in this range. In particular, the drugs studied in this work are within such limits (Figure 2), especially when there is evidence that the AmB lies flat on the membrane surface at low concentrations.71 In this view, the thickness of the membrane’s head groups and the molecule’s dimensions could be a natural structural restriction for the latter to maintain an orientation such that would cover the walls of a toroidal pore. Regarding the composition of the hydrophilic phase, MMPSM could be improved by adding water to the AMP solvent, using this mixture to calculate the solvation free energy. We plan to do this in a future work. However, we expect that the results presented here would not change significantly. The sponge hypothesis work38 does not mention the size of the drug droplet formed at the surface of the membrane, but we assume it is big enough to consider that the sterols are completely immersed, such that our model is a good approximation to estimate the solvation free energy of a sterol−drug system.

organic solvents. It has been argued that optimization of such methods is still needed to be able to treat solvents other than water.65 On the other hand, computational methods using an explicit system to study the process of membrane permeation are expensive and very sensitive to initial conditions, prone to numerical errors, and difficult to converge.66 A third approach was proposed for the main purpose of studying peripheral and trans-membrane protein systems. In the highly mobile membrane mimetic model, the hydrophobic core of the membrane in an explicit system was replaced by a fluid organic solvent, where the lipid−peptide interactions are described by an all-atom calibrated force field.67 Here, we present a method that focuses on the study of the interaction of small organic compounds with a membrane in a simplified way from a robust thermodynamic standpoint. The membrane multiphase solvation model (MMPSM) we propose takes into account the different changes in dielectric properties when going from an aqueous environment (ϵ = 80), passing through the polar head groups of the phospholipids (ϵ = 8), into the hydrophobic interior of the membrane (ϵ = 2), without sacrificing atomic detail (as opposed to an implicit membrane model68). We compare a two-phase model (M1) to a threephase model (M2). Even though M2 is more expensive computationally, it has the clear advantage over M1 to better represent a biological membrane, due to the explicit polar vs nonpolar phase separation. M1 mixes the two phases into the same solvent. With the use of M2, we could identify the location of the drugs when interacting with the membrane (compare data of Water → Octanol -M1- with data of Water → AMP and AMP → Hexadecane -M2- in Table 2). When applied to the AmB−membrane system, MMPSM predicts that the drugs will stay at the polar phase of the membrane (AMP in our case), in accordance with hypothesis H2. An infrared spectroscopic study showed that the AmB molecules bound to the membrane were located in the polar head groups.69 This was also seen in an explicit membrane molecular dynamics study when a very low AmB concentration was used.70 Here, MMPSM also shows that, at higher concentrations, it is possible that the aggregation and adsorption of the drug at the membrane’s polar phase could act as a tunnel, connecting the aqueous phase with the nonpolar membrane’s interior (see the A21 dimer Water → Hexadecane data in Table 2). In such a case, the drug could spontaneously permeate from the solution into the hydrophobic membrane’s core, in line with the “standard” hypothesis H1. This concentration dependent mechanism was observed before by electron spin resonance spectroscopy.71 At low concentrations, AmB molecules lie flat (membrane adsorption), whereas at higher concentrations, reorientation to the vertical position occurs (membrane penetration). Our results also seem to agree with hypothesis H3, which states that the AmB will extract sterols present in the membrane, with a preference for ergosterol over cholesterol. Such results imply that AmB is a better solvent for the sterol molecules than the hydrophobic interior of the membrane, in accordance with the MMPSM prediction. MMPSM offers a molecular mechanism that connects all three hypotheses through a drug concentration dependence, where several effects could be taking place in the complex antibiotic process. In light of our results, we also propose that the drug promotes the formation of membrane toroidal pores, instead of the putative transmembrane channel-like barrel-stave structures suggested by the “standard” hypothesis. This

5. CONCLUSIONS In this work, we present a computational methodology of general interest to study the thermodynamic interaction of biological membranes with small molecules. When applied to the analysis of the molecular action mechanism of AmB and a 3394

DOI: 10.1021/acs.jctc.7b00337 J. Chem. Theory Comput. 2017, 13, 3388−3397

Article

Journal of Chemical Theory and Computation

III. Molecular structure of polyene antibiotic−cholesterol complexes. Biochim. Biophys. Acta, Biomembr. 1974, 339, 57−70. (12) Mazerski, J.; Bolard, J.; Borowski, E. Self-association of some polyene macrolide antibiotics in aqueous media. Biochim. Biophys. Acta, Gen. Subj. 1982, 719, 11−17. (13) Bolard, J. How do the polyene macrolide antibiotics affect the cellular membrane-properties. Biochim. Biophys. Acta, Rev. Biomembr. 1986, 864, 257−304. (14) Bolard, J.; Legrand, P.; Heitz, F.; Cybulska, B. One-sided Action of Amphotericin B on Cholesterol-Containing Membranes Is Determined by Its Self-Association in the Medium. Biochemistry 1991, 30, 5707−5715. (15) Legrand, P.; Romero, E. A.; Cohen, E.; Bolard, J. Effects of Aggregation and Solvent on the Toxicity of Amphotericin B to Human Erythrocytes. Antimicrob. Agents Chemother. 1992, 36, 2518−2522. (16) Borowski, E. Novel approaches in the rational design of antifungal agents of low toxicity. Farmaco 2000, 55, 206−208. (17) Hartsel, S. C.; Baas, B.; Bauer, E.; Foree, L. T. J.; Kindt, K.; Preis, H.; Scott, A.; Kwong, E. H.; Ramaswamy, M.; Wasan, M. Heatinduced superaggregation of amphotericin B modifies its interaction with serum proteins and lipoproteins and stimulation of TNF-alpha. J. Pharm. Sci. 2001, 90, 124−133. (18) Huang, W.; Zhang, Z.; Han, X.; Tang, J.; Wang, J.; Dong, S.; Wang, E. Ion Channel Behavior of Amphotericin B in Sterol-Free and Cholesterol- or Ergosterol-Containing Supported Phosphatidylcholine Bilayer Model Membranes Investigated by Electrochemistry and Spectroscopy. Biophys. J. 2002, 83, 3245−3255. (19) Baginski, M.; Resat, H.; Borowski, E. Comparative molecular dynamics simulations of amphotericin B cholesterol/ergosterol membrane channels. Biochim. Biophys. Acta, Biomembr. 2002, 1567, 63−78. (20) Récamier, K. S.; Hernández-Gómez, A.; González-Damián, J.; Ortega-Blake, I. Effect of Membrane Structure on the Action of Polyenes: I. Nystatin Action in Cholesterol- and ErgosterolContaining Membranes. J. Membr. Biol. 2010, 237, 31−40. (21) González-Damián, J.; Ortega-Blake, I. Effect of Membrane Structure on the Action of Polyenes II: Nystatin Activity along the Phase Diagram of Cholesterol- and Ergosterol-Containing Membranes. J. Membr. Biol. 2010, 237, 41−49. (22) Neumann, A.; Baginski, M.; Czub, J. How Do Sterols Determine the Antifungal Activity of Amphotericin B? Free Energy of Binding between the Drug and Its Membrane Targets. J. Am. Chem. Soc. 2010, 132, 18266−18272. (23) Palacios, D. S.; Dailey, I.; Siebert, D. M.; Wilcock, B. C.; Burke, M. D. Synthesis-enabled functional group deletions reveal key underpinnings of amphotericin B ion channel and antifungal activities. Proc. Natl. Acad. Sci. U. S. A. 2011, 108, 6733−6738. (24) Gray, K. C.; Palacios, D. S.; Dailey, I.; Endo, M. M.; Uno, B. E.; Wilcock, B. C.; Burke, M. D. Amphotericin primarily kills yeast by simply binding ergosterol. Proc. Natl. Acad. Sci. U. S. A. 2012, 109, 2234−2239. (25) Kamiński, D. M. Recent progress in the study of the interactions of amphotericin B with cholesterol and ergosterol in lipid environments. Eur. Biophys. J. 2014, 43, 453−467. (26) Antillón, A.; de Vries, A. H.; Espinosa-Caballero, M.; FalcónGonzález, J. M.; Flores-Romero, D.; González-Damián, J.; JiménezMontejo, F. E.; León-Buitimea, A.; López-Ortiz, M.; Magana, R.; Marrink, S. J.; Morales-Nava, R.; Periole, X.; Reyes-Esparza, J.; Rodríguez-Lozada, J.; Santiago-Angelino, T. M.; Vargas-González, M. C.; Regla, I.; Carrillo-Tripp, M.; Fernández-Zertuche, M.; RodríguezFragoso, L.; Ortega-Blake, I. An Amphotericin B derivative equally potent to Amphotericin B and with increased safety. PLoS One 2016, 11, e0162171. (27) Ermishkin, L. N.; Kasumov, K. M.; Potzeluyev, V. Single ionic channels induced in lipid bilayers by polyene antibiotics amphotericinB and nystatine. Nature 1976, 262, 698−699. (28) Milhaud, J.; Ponsinet, V.; Takashi, M.; Michels, B. Interactions of the drug amphotericin B with phospholipid membranes containing

chemical analogue with increased safety, our results provide significant physical insight into such a relevant process. We found that (i) both drugs dimerize in all solvents studied here, an early stage in the mechanism of action, (ii) it is energetically unfavorable for the drug to penetrate the hydrophobic core of the membrane at very low concentration, (iii) when the drug’s concentration is increased, it could surpass the polar phase energy barrier and permeate into the membrane’s hydrophobic core, and (iv) it is thermodynamically possible that the sterols migrate from the membrane into a drug droplet adsorbed at the surface of the bilayer. This thermodynamic study integrates all current views regarding the complex molecular action mechanism of AmB.



AUTHOR INFORMATION

Corresponding Author

*(M.C.-T.) E-mail: [email protected]. ORCID

M. Carrillo-Tripp: 0000-0003-0060-6486 Present Address §

(M.C.-T.) Ciencias de la Computación, Centro de Investigación en Matemáticas, A.C., Jalisco S/N, Col. Valenciana, C.P. 36023, Guanajuato, Guanajuato, México.

Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This study was supported by the Consejo Nacional de Ciencia y Tecnologiá México (Grant Number 128575 to I.O.-B, Grant Number 132376 to M.C.-T., and postdoctoral scholarship to J.M.F.-G.). The calculations were performed in the high performance computing center Insurgente at the Centro de Investigación en Matemáticas and Abacus at CINVESTAV, México.



REFERENCES

(1) Rust, D. M.; Jameson, G. The novel lipid delivery system of amphotericin B: drug profile and relevance to clinical practice. Oncol. Nurs. Forum 1998, 25, 35−48. (2) Hann, I. M.; Prentice, H. G. Lipid-based amphotericin B: a review of the last 10 years of use. Int. J. Antimicrob. Agents 2001, 17, 161−169. (3) Torrado, J. J.; Espada, R.; Ballesteros, M. P.; Torrado-Santiago, S. Amphotericin B formulations and drug targeting. J. Pharm. Sci. 2008, 97, 2405−2425. (4) Thornton, S. J.; Wasan, K. M. The reformulation of amphotericin B for oral administration to treat systemic fungal infections and visceral leishmaniasis. Expert Opin. Drug Delivery 2009, 6, 271−284. (5) Hamill, R. J. Amphotericin B formulations: a comparative review of efficacy and toxicity. Drugs 2013, 73, 919−934. (6) Liu, M.; Chen, M.; Yang, Z. Design of amphotericin B oral formulation for antifungal therapy. Drug Delivery 2017, 24, 1−9. (7) Mazerski, J.; Grzybowska, J.; Borowski, E. Influence of net charge on the aggregation and solubility behaviour of amphotericin B and its derivatives in aqueous media. Eur. Biophys. J. 1990, 18, 159−164. (8) Cannon, R. D.; Lamping, E.; Holmes, A. R.; Niimi, K.; Tanabe, K.; Niimi, M.; Monk, B. C. Candida albicans drug resistance-another way to cope with stress. Microbiology 2007, 153, 3211−3217. (9) Sabra, R.; Branch, R. Amphotericin B nephrotoxicity. Drug Saf. 1990, 5, 94−108. (10) Finkelstein, A.; Holz, R. Aqueous pores created in thin lipid membranes by the polyene antibiotics nystatin and amphotericin B. Membranes 1973, 2, 377−408. (11) de Kruijff, B.; Demel, R. A. Polyene antibiotic-sterol interactions in membranes of Acholeplasma laidlawii cells and lecithin liposomes. 3395

DOI: 10.1021/acs.jctc.7b00337 J. Chem. Theory Comput. 2017, 13, 3388−3397

Article

Journal of Chemical Theory and Computation or not ergosterol: new insight, into them role of ergosterol. Biochim. Biophys. Acta, Biomembr. 2002, 1558, 95−108. (29) Volmer, A. A.; Szpilman, A. M.; Carreira, E. M. Synthesis and biological evaluation of amphotericin B derivatives. Nat. Prod. Rep. 2010, 27, 1329−1349. (30) Readio, J. D.; Bittman, R. Equilibrium binding of amphotericin B and its methyl ester and borate complex to sterols. Biochim. Biophys. Acta, Biomembr. 1982, 685, 219−224. (31) Szlinder-Richert, J.; Mazerski, J.; Cybulska, B.; Grzybowska, J.; Borowski, E. MFAME, N-methyl-N-D-fructosyl amphotericin B methyl ester, a new amphotericin B derivative of low toxicity: relationship between self-association and effects on red blood cells. Biochim. Biophys. Acta, Gen. Subj. 2001, 1528, 15−24. (32) Barratt, G.; Bretagne, S. Optimizing efficacy of Amphotericin B through nanomodification. Int. J. Nanomed. 2007, 2, 301−313. (33) Van de Ven, H. V.; Paulussen, C.; Feijens, P. B.; Matheeussen, A.; Rombaut, P.; Kayaert, P.; Van den Mooter, G.; Weyenberg, W.; Cos, P.; Maes, L.; et al. PLGA nanoparticles and nanosuspensions with amphotericin B: Potent in vitro and in vivo alternatives to Fungizone and AmBisome. J. Controlled Release 2012, 161, 795−803. (34) Sokolov, L. B.; Kholodova, G. V.; Naumchik, G. N.; et al., Water-Soluble Derivatives of Polyenic Antibiotics. Advances in the Field of Study and Production of Antibiotics No. 4; Moscow, 1978; pp 4−7. (35) Andes, D.; Safdar, N.; Marchillo, K.; Conklin, R. Pharmacokinetic-Pharmacodynamic Comparison of Amphotericin B (AMB) and Two Lipid-Associated AMB Preparations, Liposomal AMB and AMB Lipid Complex, in Murine Candidiasis Models. Antimicrob. Agents Chemother. 2006, 50, 674−684. (36) Szlinder-Richert, J.; Cybulska, B.; Grzybowska, J.; Bolard, J.; Borowski, E. Interaction of amphotericin B and its low toxic derivative, N-methyl-N-D-fructosyl amphotericin B methyl ester, with fungal, mammalian and bacterial cells measured by the energy transfer method. Farmaco 2004, 59, 289−296. (37) Slisz, M.; Cybulska, B.; Mazerski, J.; Grzybowska, J.; Borowski, E. Studies of the effects of antifungal cationic derivatives of amphotericin B on human erythrocytes. J. Antibiot. 2004, 57, 669−678. (38) Anderson, T. M.; Clay, M. C.; Cioffi, A. G.; Diaz, K. A.; Hisao, G. S.; Tuttle, M. D.; Nieuwkoop, A. J.; Comellas, G.; Maryum, N.; Wang, S. N.; et al. Amphotericin forms an extramembranous and fungicidal sterol sponge. Nat. Chem. Biol. 2014, 10, 400−406. (39) Carrillo-Tripp, M.; Feller, S. E. Evidence for a mechanism by which w-3 polyunsaturated lipids may affect membrane protein function. Biochemistry 2005, 44, 10164−10169. (40) Walters, D. R.; Dutcher, J. D.; Wintersteiner, O. The structure of mycosamine. J. Am. Chem. Soc. 1957, 79, 5076−5077. (41) Ruckwardt, T.; Scott, A.; Scott, J.; Mikulecky, P.; Hartsel, S. C. Lipid and stress dependence of amphotericin B ion selective channels in sterol-free membranes. Biochim. Biophys. Acta, Biomembr. 1998, 1372, 283−288. (42) Feller, S. E. Computational Modeling of Membrane Bilayers; Academic Press: 2008. (43) Ekins, S. Computational Toxicology. Risk Assessement for Pharmaceutical and Environmental Chemicals; John Wiley and Sons Inc.: 2007. (44) Peetla, C.; Stine, A.; Labhasetwar, V. Biophysical interactions with model lipid membranes: applications in drug discovery and drug delivery. Mol. Pharmaceutics 2009, 6, 1264−1276. (45) Venable, R. M.; Zhang, Y.; Hardy, J. H.; Pastor, R. W. Molecular Dynamics Simulations of a Lipid Bilayer Membrane Fluidity. Science 1993, 262, 223−226. (46) Natesan, S.; Lukacova, V.; Peng, M.; Subramaniam, R.; Lynch, S.; Wang, Z.; Tandlich, R.; Balaz, S. Structure-Based Prediction of Drug Distribution Across the Headgroup and Core Strata of a Phospholipid Bilayer Using Surrogate Phases. Mol. Pharmaceutics 2014, 11, 3577−3595. (47) Starzyk, J.; Gruszecki, M.; Tutaj, K.; Luchowski, R.; Szlazak, R.; Wasko, P.; Grudzinski, W.; Czub, J.; Gruszecki, W. I. Self-Association of Amphotericin B: Spontaneous Formation of Molecular Structures

Responsible for the Toxic Side Effects of the Antibiotic. J. Phys. Chem. B 2014, 118, 13821−13832. (48) Jarzembska, K. N.; Kamiński, D.; Hoser, A. A.; Malinska, M.; Senczyna, B.; Wozniak, K.; Gagoś, M. Controlled crystallization, structure, and molecular properties of iodoacetylamphotericin B. Cryst. Growth Des. 2012, 12, 2336−2345. (49) Martínez, L.; Andrade, R.; Birgin, E. G.; Martínez, J. M. Packmol: A package for building initial configurations for molecular dynamics simulations. J. Comput. Chem. 2009, 30, 2157−2164. (50) Oostenbrink, C.; Villa, A.; Mark, A. E.; van Gunsteren, W. F. A. Biomolecular Force Field Based on the Free Enthalpy of Hydration and Solvation: The GROMOS Force-Field Parameter Sets 53A5 and 53A6. J. Comput. Chem. 2004, 25, 1656−1676. (51) Malde, A. K.; Zuo, L.; Breeze, M.; Stroet, M.; Poger, D.; Nair, P. C.; Oostenbrink, C.; Mark, A. E. An Automated Force Field Topology Builder (ATB) and repository: version 1.0. J. Chem. Theory Comput. 2011, 7, 4026−4037. (http://compbio.biosci.uq.edu.au/atb). (52) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; Hermans, J. Interaction Models for Water in Relation to Protein Hydration. In Intermolecular Forces; Pullman, B., Ed.; D. Reidel Publishing Company, 1981; p 331. (53) Abraham, M. J.; van der Spoel, D.; Lindahl, E.; Hess, B.; the GROMACS development team. GROMACS User Manual, version 5.1.2; 2016. www.gromacs.org. (54) (a) Byrd, R. H.; Lu, P.; Nocedal, J. A.; Zhu, C. A limited memory algorithm for bound constrained optimization. SIAM J. Scientif. Statistic. Comput. 1995, 16, 1190−1208. (b) Zhu, C.; Byrd, R. H.; Nocedal, J.; Lu, P. L-BFGS-B: Algorithm 778: L-BFGS-B, FORTRAN routines for large scale bound constrained optimization. ACM Trans. Math. Softw. 1997, 23, 550−560. (55) Berendsen, H. J. C.; Postma, J. P. M.; DiNola, A.; Haak, J. R.; van Gunsteren, W. F. Molecular dynamics with coupling to an external bath. J. Chem. Phys. 1984, 81, 3684−3690. (56) Berendsen, H. J. C. Transport properties computed by linear response through weak coupling to a bath. In Computer Simulations in Material Science; Meyer, M., Pontikis, V., Eds.; Kluwer: 1991; pp 139− 155. (57) Darden, T.; York, D.; Pedersen, L. Particle mesh Ewald: An Nlog(N) method for Ewald sums in large systems. J. Chem. Phys. 1993, 98, 10089−10092. (58) Hess, B.; Bekker, H.; Berendsen, H. J. C.; Fraaije, J. G. E. M. LINCS: A linear constraint solver for molecular simulations. J. Comput. Chem. 1997, 18, 1463−1472. (59) Mezei, M.; Beveridge, D. L. Free energy simulations. In Computer Simulations and Biomolecular Systems; Beveridge, D. L., Jorgenson, W. L., Eds.; Annals of the New York Academy of Science; John Wiley and Sons: 1986; pp 1−23. (60) Mark, A. E. In Encyclopedia of Computational Chemistry; van Rague Schleyer, P., Ed.; John Wiley & Sons: 1998; Vol. 2, p 1211. (61) Villa, A.; Mark, A. E. Calculation of the Free Energy of Solvation for Neutral Analogs of Amino Acid Side Chains. J. Comput. Chem. 2002, 23, 548−553. (62) Bussi, G.; Donadio, D.; Parrinello, M. Canonical sampling through velocity rescaling. J. Chem. Phys. 2007, 126, 014101. (63) (a) Nosé, S. A molecular dynamics method for simulations in the canonical ensemble. Mol. Phys. 1984, 52, 255−268. (b) Hoover, W. G. Canonical dynamics: equilibrium phase-space distributions. Phys. Rev. A: At., Mol., Opt. Phys. 1985, 31, 1695−1697. (64) Hess, B. Determining the shear viscosity of model liquids from molecular simulation. J. Chem. Phys. 2002, 116, 209−217. (65) Zhang, J.; Zhang, H.; Wu, T.; Wang, Q.; van der Spoel, D. Comparison of Implicit and Explicit Solvent Models for the Calculation of Solvation Free Energy in Organic Solvents. J. Chem. Theory Comput. 2017, 13, 1034−1043. (66) Zocher, F.; van der Spoel, D.; Pohl, P.; Hub, J. S. Local Partition Coefficients Govern Solute Permeability of Cholesterol- Containing Membranes. Biophys. J. 2013, 105, 2760−2770. (67) Ohkubo, Y. Z.; Pogorelov, T. V.; Arcario, M. J.; Christensen, G. A.; Tajkhorshid, E. Accelerating Membrane Insertion of Peripheral 3396

DOI: 10.1021/acs.jctc.7b00337 J. Chem. Theory Comput. 2017, 13, 3388−3397

Article

Journal of Chemical Theory and Computation Proteins with a Novel Membrane Mimetic Model. Biophys. J. 2012, 102, 2130−2139. (68) Panahi, A.; Feig, M. Dynamic Heterogeneous Dielectric Generalized Born (DHDGB): An implicit membrane model with a dynamically varying bilayer thickness. J. Chem. Theory Comput. 2013, 9, 1709−1719. (69) Gabrielska, J.; Gagoś, M.; Gubernator, J.; Gruszecki, W. I. Binding of antibiotic amphotericin B to lipid membranes: a 1H NMR study. FEBS Lett. 2006, 580, 2677−2685. (70) Sternal, K.; Czub, J.; Baginski, M. Molecular aspects of the interaction between amphotericin B and a phospholipid bilayer: molecular dynamics studies. J. Mol. Model. 2004, 10, 223−232. (71) Man, D.; Olchawa, R. Two-step impact of Amphotericin B (AmB) on lipid membranes: eSR experiment and computer simulations. J. Liposome Res. 2013, 23, 327−335. (72) Huang, H. W. Action of Antimicrobial Peptides: Two-State Model. Biochemistry 2000, 39, 8347−8352. (73) Qian, S.; Rai, D.; Heller, W. T. Alamethicin Disrupts the Cholesterol Distribution in Dimyristoyl Phosphatidylcholine−Cholesterol Lipid Bilayers. J. Phys. Chem. B 2014, 118, 11200−11208. (74) Sengupta, D.; Leontiadou, H.; Mark, A. E.; Marrink, S.-J. Toroidal Pores Formed by Antimicrobial Peptides Show Significant Disorder. Biochim. Biophys. Acta, Biomembr. 2008, 1778, 2308−2317. (75) Tristram-Nagle, S.; Nagle, J. F. Lipid bilayers: thermodynamics, structure, fluctuations, and interactions. Chem. Phys. Lipids 2004, 127, 3−14.

3397

DOI: 10.1021/acs.jctc.7b00337 J. Chem. Theory Comput. 2017, 13, 3388−3397