Testing High Concentrations of Membrane Active ... - ACS Publications

Apr 20, 2017 - Department of Physics, Brock University, St. Catharines, Ontario L2S 3A1, Canada. ‡. Institute of Molecular Biosciences, University o...
0 downloads 0 Views 5MB Size
Article pubs.acs.org/JPCB

Testing High Concentrations of Membrane Active Antibiotic Chlorhexidine Via Computational Titration and Calorimetry Brad Van Oosten,† Drew Marquardt,‡ and Thad A. Harroun*,† †

Department of Physics, Brock University, St. Catharines, Ontario L2S 3A1, Canada Institute of Molecular Biosciences, University of Graz, 8010 Graz, Austria



S Supporting Information *

ABSTRACT: Coarse-grained strategies for membrane simulations are designed to increase efficiency for larger and more complex molecular dynamics simulations. For membrane active antibiotics, the concentration dependence of their action presents a tremendous challenge in simulation scale. In this study, we examine the effects of concentration for the popular membrane active antibacterial drug chlorhexidine. It presents an interesting biophysical modeling test, where from experimental experience we know that model membranes of 1,2dimyristoyl-sn-glycero-3-phosphocholine (DMPC) can absorb very high quantities of the drug without disruption. We construct a coarse-grained model of chlorhexidine in three different charged states and compare to previous all-atom simulations and new experiments. Using large, long-time, and unbiased simulations of chlorhexidine inserting into the lipid bilayer, we find little changes to the functional structure of a DMPC membrane up through concentrations of 15:100 drug:lipid, where the slowing rate of continued insertion tests the capabilities of even this coarse-grained approach. We validate our simulations with computational calorimetry measurements, and show that they agree with new experimental data from differential scanning calorimetry.



INTRODUCTION

significant way. This seemingly contradicts previous observations of CHX’s lytic properties. Our previous simulation studies were fully atomistic and therefore limited in size to a small patch of membrane, as well as a small number of CHX molecules. Unsatisfied with this type of molecular dynamics, we identified several serious drawbacks. First, with 128 lipids in an 8 × 8× 2 formation, the system was much too small to allow for any significant membrane undulations or defects to be observed with the addition of CHX. Second, we decided on an unbiased simulation, whereby the CHX molecules were introduced into the aqueous phase and allowed to make their own room among the lipids in their own time. Removing lipids prior to simulation, in order to accommodate such a small molecule, would be questionable and defeat the purpose of the simulation. Furthermore, adding a significant number of CHX resulted in a huge expansion of the computational time, as many CHX molecules spent a great deal of time interacting with each other in the aqueous environment before interacting with the lipids. CHX inserting into the membrane was shown to lower the system energy, however the energy difference was not large, and thus the system was slowly driven to equilibrium. We concluded that titrating small numbers of CHX molecules into the system should reduce the time of insertion and equilibrium. This would be true at least until the

Chlorhexidine (CHX) is an effective antibacterial drug that attacks the membrane of many bacteria, yeasts, and certain viruses, making it a leading component in many mouth rinses, industrial cleansers, and medical cleaning solutions.1−4 At low concentrations of 0.5 mg/L, it is able to stop bacteria from reproducing, limiting the spread of contamination, while a concentration of 10 mg/L is enough to lyse the cell, eliminating most of the bacterial culture.5 Although the practical uses of CHX have been widely explored, the molecular mechanism of CHX uptake by cells is still under investigation. CHX has been shown to cause damage to the outer layers of bacteria of yeast cells;6 however, this was not shown to be the major cause of cell death. It was found that CHX partitions into the plasma membrane and cytoplasm of cells, where cell leakage of intracellular components followed, leading to cell lysis.7 At high concentrations of CHX, the intracellular components of the bacteria are seen to coagulate, and the amount of leakage is reduced.8 In our earlier studies, we tested CHX in a model 1,2dimyristoyl-sn-glycero-3-phosphocholine (14:0−14:0-PC, DMPC) membrane9,10 in order to investigate the interaction of CHX with model lipid bilayers. We showed through experimental neutron scattering and Langmuir monolayers9 as well as molecular dynamics simulations10 that CHX has a high affinity for DMPC. However, we also noted that the membrane can absorb a remarkable amount of CHX molecules, as high as 30:100 CHX:DMPC, without losing its bilayer structure in any © 2017 American Chemical Society

Received: December 12, 2016 Revised: April 19, 2017 Published: April 20, 2017 4657

DOI: 10.1021/acs.jpcb.6b12510 J. Phys. Chem. B 2017, 121, 4657−4668

The Journal of Physical Chemistry B



Article

METHODS CHX Parametrization. The coarse-graining (CG) method utilized to create the CHX models was based on our previous study of lipopolysaccharides which showed success utilizing a bottom-up approach, building a coarse-grained model based on all-atom (AA) simulations.10,13 An AA coordinate file of CHX in a water box was converted to CG by mapping real atoms to virtual particles at the location of the center of mass of the atom collection they represent, as shown in the schematic of Figure 1. The mapping was chosen to be similar to the virtual particle

concentration of CHX in the bilayer was sufficiently high that the free energy gained by further insertion is reduced, slowing the time of insertion once again. We also reasoned that a much larger system should be used so that any disruption or instability in the membrane may be observed. Both of these ideas are realistically out of reach for an atomistic simulation in a reasonable computational time frame. The purpose of the MARTINI force field has been to overcome such problems by reducing the computational complexity of the system using a mapping of approximately 4 heavy atoms to 1 virtual particle.11,12 This reduces the total amount of all atoms by an approximate factor of 10 and also allows for a time step 10 times larger, corresponding to an approximate 100× efficiency gain. The reduced degrees of freedom simplify the energy landscape, with the advantage of gaining physical scale with the same computational resources, but at the cost of near total loss of certainty in kinetic and dynamic events. As a result, a large, unbiased simulation of concentrations approaching the expected minimum inhibitory concentration (20:100) can be performed. Beyond this point, the time for additional CHX to enter the bilayer grows, preventing us from reporting extreme concentrations at this time. The MARTINI force field contains a large library of lipid types as well as amino acids and an ever growing array of small molecules, with simulated properties aiming to reproduce many key large-scale properties of all-atom simulations as well as experiments. However, since most small molecules are unique and can not be created out of building blocks such as amino acids, there exists very few available parameters for these. In this study, we introduce a MARTINI-compatible force field for three CHX models (CHX0+, CHX1+, and CHX2+) in order to expand our study on the mechanism of CHX. The MARTINI atom types were chosen to reproduce the known partitioning of CHX into DMPC molecules from neutron diffraction experiments. How does CHX reach high concentrations at the location of the membrane we know it to reside? We perform a unique titration-like experiment where we introduce small amounts of the CHX molecule into a model DMPC membrane at set intervals. From this we are able to investigate the concentration effect of CHX on a DMPC membrane, as well as create equilibrated systems of varying CHX concentration. Finally, we validate our model by comparing the concentration effects with experimental results of (1) previously published neutron diffraction data that established the penetration depth of CHX into DMPC bilayers,9 and (2) new differential scanning calorimetry (DSC) data of the effect of concentration on the main chain melting transition temperature. The overarching goal of this study is to test the ability of molecular dynamics to correctly handle large numbers of CHX molecules in a model DMPC membrane as verified by experimentally observed effects of concentration on the lipid main transition. We show that the amount of energy gained by adding CHX into the membrane decreases as we saturate the membrane with CHX, and look specifically for membrane conformational changes and CHX-localized aggregation, which could lead to membrane instability. We also make observations on the lowered affinity between the molecule and the membrane as well as increased simulation time necessary to reach an equilibrium state, with increased saturation, and thus the importance of using a coarse grained model and the MARTINI force field.

Figure 1. CG mapping and naming schematic of CHX sections and virtual particles (top). Simulation snapshot of AA CHX molecule with virtual particles based on the center of mass (red dot) of the atoms they represent (bottom).

choice of the amino acids Tyr, Arg, and Lys used in the MARTINI force field version 2.2P, as these amino acids guided our original choices for the chlorphenol (CPL), buiguanide (BGU), and hexane (HEX) residues of our AA CHX model.12,14 A 100 ns AA simulation of CHX in a water-filled box using our previously published topology was used to construct histograms of the probability distributions of bond lengths, angles, and dihedrals of the CG-mapped centers of mass. CG simulations of CHX in water were then performed, and the resulting CG histograms of bonds, angles, and dihedrals were directly compared to the respective AA histograms. The spring strength and harmonic equilibrium parameters were iteratively adjusted using 100 ns CG simulations until all distributions for bonds, angles, and dihedrals matched those from the AA simulations. The model was then placed in a DMPC membrane where the spring strength and harmonic equilibrium parameters were iteratively adjusted again to match those of an AA simulation in DMPC. At this point, the atom types were also adjusted to bring the CHX sections into the correct depth in the membrane as compared to the AA simulations and experiment (see Figure 2). Greater detail on this method and results can be found in the Supporting Information. Often CG parametrization is based on matching the measured free energy of partitioning in bulk phases, such as the octanol/water system as a simplistic model of the membrane’s interfaces. However, the partitioning of solutes into bilayer membranes is of a fundamentally different nature. For example, hydrophilicity, due to acceptance of a hydrogen bond, is not reflected by partitioning in octanol.15 Basic amines show higher partition into model membranes than predicted by bulk octanol partitioning, as lipid order plays more of a role. An ordered phospholipid bilayer can support the partitioning of 4658

DOI: 10.1021/acs.jpcb.6b12510 J. Phys. Chem. B 2017, 121, 4657−4668

Article

The Journal of Physical Chemistry B

depth in the bilayer should determine the location of the minimum in the free energy profile. The result is a test to see whether CHX, so situated, can reach the high concentrations established for it without membrane disruption. Membrane Creation and Simulated Titration. A pure DMPC bilayer consisting of 900 CG lipids was constructed using the INSANE membrane building tool and equilibrated for 250 ns.19 Titration of CHX into this model membrane was then performed. A single titration step included adding 9 CHX1+ molecules (1% concentration) to the aqueous solution in a 3 × 3 configuration, arranged as far apart as geometrically possible to minimize their interactions with each other. The CHX were also equally distant from the bilayer and its periodic boundary image, meaning each CHX could freely enter the bilayer from either side. The system was then simulated in 500 ns intervals until all 9 CHX1+ molecules were fully inserted into the membrane. This titration step was repeated 15 times so that it created 16 model membranes with equilibrated concentrations of CHX1+ increasing to 15:100 CHX1+:DMPC. A second membrane consisting of 800 DMPC lipids and 80 cholesterol molecules was also created using the INSANE membrane building tool, and equilibrated for 250 ns. For this membrane, eight CHX molecules were added for each titration step. This process was repeated four more times to bring the membrane to a 5:100 CHX1+:DMPC concentration. For subsequent lower cholesterol concentrations, an equal amount of cholesterol molecules were removed at random from each of the various equilibrated X:1:10 CHX:CHOL:DMPC membranes, with equal amounts removed from both membrane leaflets. The membranes were then equilibrated for 250 ns to allow the membrane to fill in the gaps of the missing cholesterol molecules. Simulation Details. All simulations were performed using GROMACS 5.0.4.20 The coarse-grained simulations were performed using a time step of 10 fs in a periodic xyz orthorhombic box. The temperature was maintained for the lipids, cholesterol, CHX, and the solvent independently with a coupling constant of 10 ps using a Nose Hoover thermostat.21,22 Semi-isotropic coupling was used to maintain a pressure of 1.0 bar using the Parrinello−Rahman barostat23 with a 1.0 ps coupling constant and compressibility of 4.5 × 10−5 bar−1. The Verlet group scheme was used for pair lists which were updated every 20 steps. Long-ranged electrostatics were computed via the particle mesh Ewald scheme, while both short ranged Coulomb and LJ potentials are shifted to zero at an exact cutoff value of 1.2. The polarized water of MARTINI 2.2P24 was used with a screening dielectric of ϵr = 2.5. All atomistic simulations utilized the Slipids force field,25−27 including a pure DMPC membrane consisting of 200 lipids and 6000 TIP3P water molecules built by CHARMM-GUI.28 A time step of 2 fs was used with the system in a periodic xyz orthorhombic box. The temperature was maintained for each component (DMPC, CHX, CHOL, and water) independently with a coupling constant of 1.0 ps using a Nose−Hoover thermostat.21,22 Semi-isotropic coupling was used to maintain a pressure of 1.013 bar using the Parrinello−Rahman barostat23 with a 10 ps coupling constant and compressibility of 4.5 × 10−5 bar−1. The Verlet group scheme was used for pair lists, which were updated every 20 steps. Long-ranged electrostatics were computed via the particle mesh Ewald scheme, while both short-ranged Coulomb and LJ potentials are shifted to zero at an exact cutoff value of 1.2

Figure 2. Top: Simulation snapshots of AA and CG simulations with CHX in pink, lipid tails in green, lipid headgroups in orange, and water in blue. Center: Volume probability of AA DMPC system with DB (width of membrane with less than 50% volume probability of water) and DHH (headgroup to headgroup width) thickness values marked. Bottom: Volume probability of each of the 3 CHX sections: CPL (teal), BGU (red), and HEX (black) for the CHX0+, CHX1+, and CHX2+ AA (solid lines) and CG (dashed lines) models. Bold, black arrows above the plot show the experimentally determined location of the HEX.9

both neutral and positively charged forms of small molecule amines, whereas octanol best supports neutral forms.15,16 This seems to extend to small cyclic hydrocarbons such as benzene as well.17 Even when octanol reproduces the free energy of adsorption of amino acids into the membrane surface, it does not capture the free energy for partitioning into the center of the bilayer at all.18 Thus, we lean more heavily toward what we know of CHX partitioning into DMPC bilayer from experiment, where its 4659

DOI: 10.1021/acs.jpcb.6b12510 J. Phys. Chem. B 2017, 121, 4657−4668

Article

The Journal of Physical Chemistry B Computational Calorimetry. Each membrane was equilibrated at 36 °C for 250 ns. The temperature was then reduced by steps of 2 °C and simulated for 25 ns at each temperature, until the membrane exhibited signs of a transition from liquid to gel state. This point was defined as a sudden change in the area per lipid, in which it no longer followed the linear decrease. At this point, the coldest temperature without any signs of a transition was simulated for 100 ns. If the membrane still went through the transition at this temperature, the next warmest simulation was repeated for 100 ns until no transition occurred. This was set as the new starting point for increased sampling and cooling with a smaller temperature step. The temperature step was reduced to cooling by steps of 0.5 °C and each temperature was simulated for 100 ns. Several other temperature steps and simulation times were tested as well (see the Results section for details). The cooling was continued until the system reached 0 °C or until the system was completely in the gel phase, whichever occurred latest. An approximate measure of specific heat was calculated at each temperature point to simulate a differential scanning calorimetry experiment. The specific heat at constant pressure (Cp) calculated is defined in terms of the fluctuations in energy (E) as

Cp =

molecule in place with a harmonic restraining force on the center of mass equal to 3000 kJ mol−1 nm−1 in the direction normal to the bilayer. PMFs were then calculated using the weighted histogram analysis method29,30 utilizing the gromacs g_wham tool.31 The first 25 ns of each trajectory were not included in the PMF calculation but rather used as further equilibration time for each window. Therefore, the PMF was calculated using a total of 3000 ns of simulation time (75 ns ×40 windows). Experimental Calorimetry. Materials. 1,2-dimyristoyl-snglycero-3-phosphocholine (14:0−14:0-PC, DMPC) was purchased from Avanti Polar Lipids (Alabaster, AL), cholesterol and chlorhexidine dihydrochloride (in powder form) were purchased from Sigma-Aldrich (St. Louis, MO), and all materials were used as received. Lipid stock solutions were prepared by dissolving dry lipid powder in HPLC-grade chloroform. Ultrapure H2O was obtained from a High-Q purification system (Wilmette, IL). Sample Preparation. Lipid:cholesterol:chlorhexidine films were prepared by mixing the appropriate volumes of lipid/ chloroform solution with CHX/chloroform solution to yield the desired ratios. The organic solvent was removed with an N2 stream and gentle heating, followed by drying in vacuo (∼2 h) to ensure that all residual chloroform was removed. Multilamellar vesicles were prepared by hydrating a preheated (30 °C) film with ultrapure H2O to a concentration of 30 mg/mL, followed by vigorous vortexing to disperse the lipid. The resulting MLV suspension was subjected to five freeze/thaw cycles between −78 and 40 °C. For each sample, 30 μL of the MLV suspension was loaded into a Shimadzu aluminum hermetic pan (Mandel Scientific, Guelph, ON) and crimped sealed. DSC Measurements. Differential scanning calorimetry (DSC) thermograms were collected using a Shimadzu DSC60 and a TA-60WS Thermal Analyzer from Mandel Scientific (Guelph, ON). Data was collected as a series of heating and cooling scans at a rate of 2 °C/min using a constant supply of liquid nitrogen as the coolant. The raw DSC data were visualized and analyzed using TA-60 version 2.11.

⟨E2⟩ − ⟨E⟩2 (kBT )2

Where kB is Boltzmann’s constant and T is the average temperature of the simulation. As shown below, where a liquid to gel transition occurs during a simulation run, the decrease in average system energy is included as part of its variance in the calculation. Weighted Histogram Analysis Method (WHAM). Twelve sets of umbrella sampling simulations were conducted, six AA and six CG, to determine the potential of mean force (PMF) and free energy experienced by CHX or solvent as it is pulled through the lipid bilayer.29 The six sets include pulling each of the 3 CHX charge models into a pure DMPC membrane, pulling one CHX1+ into a 10:100 CHX1+:DMPC membrane, pulling water into a pure DMPC membrane, and finally, pulling water into a 10:100 CHX1+:DMPC membrane. For the simulations with a CHX molecule as the pull group, the system was built with the membrane at the center of the box, and a single CHX molecule in the water, equidistant to either side of the membrane, taking periodic boundary conditions into account. The CHX molecule was held in place with a harmonic constraint for a 1 ns equilibration period for the water and CHX to settle. The set of simulations with a water molecule as the pull group used a pre-existing equilibrated membrane, and a water molecule equidistant from either side of the membrane was chosen at random as the pull group. After equilibration, the pulled molecular group (either water or CHX) was brought into the membrane from its starting position in the water to the center of the membrane, with a harmonic force constant of 1000 kJ mol−1 nm−1 and at a rate of 0.001 nm ps−1. From this pulling simulation, 40 configurations were isolated, each with the z distance separating the center of mass of the pulled molecule with the center of the membrane, spaced 0.1 nm apart, ranging from 0 to 4.0 nm. Each of these configurations formed the starting point for a 100 ns simulation at 309 K, with a time step of 10 and 2 fs for CG and AA simulations, respectively. Each simulation held the test



RESULTS The CG CHX Model. Model Validation. Initial validation of our CG model rests on reproducing the structural results our AA model, as well as correct partitioning into a DMPC membrane as compared to our previous neutron diffraction data.9,10 Aside from the overall correct structure of the CG bilayer shown in the top and middle parts of Figure 2, we also looked at the penetration depth of each of the three sections of our CHX models in DMPC. The bottom panel of Figure 2 shows the distribution depth of the CHX sections, and the bold arrows above them indicate the experimentally determined location of the HEX.9 Not only do the AA (solid lines) and CG (dashed lines) overlap quite well, but the charged models agree quite well with experiment. As before, we find that increased charge of CHX brings it closer to the aqueous interface, but it firmly remains inside the bilayer. With correct positioning of each component, we gained assurance in the particle type choices. Further confidence in the particle choice and nonbonded energetics is further explored in the following Energetics section. Bonded interactions as well as the conformation of CHX as either in an extended or collapsed state are compared and discussed in the supplemental results. 4660

DOI: 10.1021/acs.jpcb.6b12510 J. Phys. Chem. B 2017, 121, 4657−4668

Article

The Journal of Physical Chemistry B Energetics. To compare the free energy of CHX insertion into the bilayer between the CG and AA simulations, we used the WHAM method to calculate the PMF along the reaction coordinate of distance from the bilayer center. These PMF plots were used to observe changes in the energy of insertion (ΔGi), the energy of passing through the center of the membrane, or flip-flop, (ΔGff), and the equilibrium distance from the bilayer center corresponding to the location of lowest energy (Z0) . PMF plots generated using different force fields will generally result in different values for ΔGi and ΔGff which can substantially disagree with each other.12,14,32 This is an inherent problem arising from the different nonbonded parameters used by different force fields. Therefore, the reproduction of the absolute value of these energy values was not held at as high of an importance as reproducing qualitative energy changes caused by altering the system. Furthermore, to our knowledge, there is currently no known experimental value for either the ΔGi or ΔGff of CHX in DMPC. We first compare the three charge models of the CHX molecule in a pure DMPC membrane. Figure 3 shows the

Table 1. PMF Calculated Values for the Energy of Insertion (ΔGi), Energy Barrier to Cross the Membrane Center (ΔGff), and the Equilibrium Distance from the Membrane Center (Z0) pulled molecule AA CHX AA CHX1+ AA CHX2+ CG CHX CG CHX1+ CG CHX2+ AA CHX1+ CG CHX1+ AA Water AA Water CG Water CG Water

membrane DMPC DMPC DMPC DMPC DMPC DMPC 10:100 CHX1+:DMPC 10:100 CHX1+:DMPC DMPC 10:100 CHX1+:DMPC DMPC 10:100 CHX1+:DMPC

ΔGi (kJ/mol)

ΔGff (kJ/mol)

Z0 (nm)

43.34 32.94 17.06 52.96 58.24 45.25 23.32

10.43 36.84 48.28 11.54 37.42 55.72 31.69

0.94 1.14 1.15 0.95 1.14 1.19 1.07

47.05

14.53

1.15

-

25.72 29.48

∞ ∞

-

28.6 24.61

∞ ∞

minimum location is still located inside the bilayer for the CHX2+ model, the driving force from the aqueous phase into the bilayer is greatly reduced. Quantitatively, the CG simulations overestimate ΔG i compared to AA simulations. The value of ΔGi decreases with each increase in charge in the AA simulations with values of 43.34, 32.94, and 17.06 kJ/mol for the CHX, CHX1+, and CHX2+ models, respectively. The CG simulations, however, did not have a significant change, having values of 52.96, 58.24, and 45.25 kJ/mol for the CHX, CHX1+, and CHX2+ models, respectively. The higher ΔGi values of the CG models could consequently make the affinity of DMPC for CHX molecules artificially high. PMF plots were also calculated to observe any difference in the energy of inserting a single CHX1+ molecule into a pure DMPC membrane as compared to one already containing a large amount of CHX, at 10:100 CHX1+:DMPC. The central frame of Figure 3 shows that in both the AA and CG simulations, the PMF of the high-concentration CHX bilayer (in red) is broader and more shallow, showing a decrease in both ΔGi and ΔGff. The 10:100 CHX1+:DMPC membrane therefore has a reduced affinity to absorb more CHX molecules. The addition of CHX1+ molecules to the membrane also has the effect of lowering the barrier of a CHX1+ molecule crossing through the membrane. ΔGi decreases by 9.62 kJ/mol, a 29% decrease in energy in the AA simulation, and by 11.19 kJ/mol, a 19% decrease in the CG simulation. The barrier ΔGff decreases by 5.15 kJ/mol, a 14% decrease in energy in the AA simulation, and by 22.89 kJ/mol, a 61% decrease in the CG simulation. Although there is a clear effect on the energy of insertion and barrier at the center of the membrane, the equilibrium position of CHX remains unchanged at 1.15 nm in the CG simulation and slightly closer to the membrane center at 1.07 nm in the AA. Since CHX is putatively a membrane-disrupting molecule, the permittivity of water was also investigated in order to look for evidence of an onset of pore formations. A single water molecule was pulled through a pure DMPC membrane as well as the 10:100 CHX1+:DMPC membrane. Figure 4 demonstrates very little alteration in the PMF plots for the insertion of

Figure 3. PMF of the three CHX charge models pulled through AA (solid lines) and CG (dashed lines) bilayers: neutral CHX0+ (top), CHX1+ (middle), CHX2+ (bottom). The black curves use a pure DMPC membrane. The red curve in the central frame uses a membrane preloaded with CHX1+ to a concentration of 1:10 CHX1+:DMPC.

resulting PMF plots from both AA and CG models, and Table 1 lists the extracted numerical values. The location of lowest energy, Z0, of the CHX molecule moves further from the membrane center with increasing charge. For the AA simulations, the lowest energy position of CHX increased from 0.94 nm, to 1.14 nm, and 1.15 nm for the CHX, CHX1+, and CHX2+ models, respectively. Similarly, the CG models also increase their equilibrium position with increased charge with values of 0.95, 1.14, and 1.19 nm. Both the CG and AA simulations qualitatively agree that increasing the charge also increases the energy barrier for CHX to flip-flop from one leaflet to the other. The neutral CHX0+ has the lowest ΔGff of the three models, being much smaller than ΔGi, by a factor of 4. The neutral CHX0+ molecule is therefore much more likely to flip-flop from one side of the bilayer to the other rather than to leave the membrane entirely. Adding a + 1 charge to the CHX raises ΔGff to be comparable, but slightly smaller in size to ΔGi. The CHX2+ model has the highest ΔGff, notably larger than ΔGi. Although the energy 4661

DOI: 10.1021/acs.jpcb.6b12510 J. Phys. Chem. B 2017, 121, 4657−4668

Article

The Journal of Physical Chemistry B

thickness: the overall bilayer thickness, DB, also known as Luzzati thickness in neutron scattering,33 and distance between lipid headgroups, DHH, commonly used in X-ray scattering (see Figure 2).34,35 The CG membrane has a constant DB thickness value of approximately 2.95 nm until 6:100 CHX:DMPC, where there is a sudden drop in the thickness to 2.82 nm. The DHH thickness of the membrane is less steady than the DB thickness, with values of approximately 1.13 nm with no upward or downward trend visible. The smaller thickness values of the CG simulation than that of the AA simulation arise from the particular atom-groupings of the DMPC model of the MARTINI force field, which groups four acyl-chain carbons into each virtual particle. DMPC having 14 carbons can not evenly be divided by 4, and is therefore modeled using the tails of three virtual particles, representing the same thickness as a lipid with length of 12 carbons (DLPC). The alternate choice to represent DMPC with four particles, which represents a 16 carbon length (DPPC), is no more accurate. Since the MARTINI DMPC lipids consist of lipid tails represented by three particles, the order parameter along the acyl chain is ill defined by usual means. However, a second rank order parameter (Sn) can be calculated using

Figure 4. PMF of one water molecule (solid line AA and dashed line CG) passing through a pure DMPC membrane (top frame) and 10:100 CHX1+:DMPC membrane (bottom frame).

water due to an increase in concentration of CHX for both of our AA and CG models. Therefore, there is no substantial evidence from our simulations that a CHX:DMPC ratio of 10:100 or lower will increase the permittivity of water and thus, no increase in the chances of pore formation in a model DMPC bilayer. Structural Effects. Figure 5 shows the change in key structural features of the DMPC membrane as a function of the

Sn =

3 1 cos2(θ ) − 2 2

where θ is the angle made between the bond between two particles and the bilayer normal. The CG order parameter values are overall higher due to redefining the torsion angles to be between groups of atoms and not between carbon atoms themselves, however, the similar characteristic of higher ordering at the headgroups and lower ordering in the tails, is still present. The upper panel in Figure 6 shows the CG order parameter alongside the order parameter from AA simulations. With or without CHX1+, both the CG and AA simulations show a similar ordering in the headgroups, and a decrease in order in the tails. Although the change in Sn due to increased CHX

Figure 5. Structural dimensions for AA (solid lines) and CG (dashed lines) as a function of CHX concentration. Area per lipid (top), total membrane thickness DB (middle) and membrane backbone to backbone thickness DHH (bottom).

CHX1+ concentration. The data is compared to AA simulations done in a previous study of CHX in DMPC at ratios of 0:100, 1:100, 3:100 and 10:100 CHX:DMPC (see Van Oosten et al. for simulation details10). A steady linear increase in the area per lipid is seen as the CHX1+ concentration is increased. Conversely, there is little change in the lipid thickness, showing no obvious trend. This suggests that the lipids are moving aside, and the membrane expands as the CHX1+ molecule wedges itself in. The AA simulations used as a comparison also show a linear increase in the area per lipid as well as a constant thickness. Two different types of the membrane thickness values were calculated to investigate the changes to the membrane

Figure 6. Top: Deuterium order parameter averaged over both lipid tails for the AA simulation of the pure DMPC membrane (solid black) and the 10:100 CHX1+:DMPC membrane(solid red), as well as the second rank order parameter calculated from the angle between the bond connecting particles in the CG simulations of the pure DMPC membrane (dashed black) and the 10:100 CHX1+:DMPC membrane (dashed red). Bottom: Boxplot representing the maximum, minimum, and average order parameter value at each of the CHX concentrations for the AA simulations (solid color) and CG simulations (hatched color). 4662

DOI: 10.1021/acs.jpcb.6b12510 J. Phys. Chem. B 2017, 121, 4657−4668

Article

The Journal of Physical Chemistry B concentration is small, the lower panel of Figure 6 shows the maximum, minimum, and average values of Sn as a function of CHX concentration, where a small and steady decrease in order with increasing CHX1+ concentration can be observed. Mechanism of Insertion. CHX, with any charge, is mainly hydrophobic and as a result are seen to cluster in aqueous solution and enter the membrane together. After entering the bilayer, they dissociate almost immediately. Figure 7A shows

Figure 8. Top: Two-dimensional, in-plane, CHX1+−CHX1+ radial distribution function for the 10:100 CHX1+:DMPC simulation in the top leaflet (black) and bottom leaflet (red). Bottom: Two-dimensional, in-plane, DMPC phosphate−DMPC phosphate radial distribution of pure DMPC membrane (black) and 10:100 CHX1+:DMPC membrane (green).

membrane, at the depth noted before. Figure 8B shows the two-dimensional radial distribution function for the DMPC phosphate groups. There are nearest-neighbor and nextnearest-neighbor correlation distances of ∼0.9 nm and ∼1.45 nm, respectively, that remain unchanged with the addition of all the CHX1+. This suggests that there is no differential lateral ordering of the lipids in the membrane due to the addition CHX. Since the area per lipid has increased, the extra area is being taken up by CHX and the DMPC acyl tails, while the CHX prevents a thinning of the membrane. The amount of time taken for each of the CHX molecules to insert itself into the membrane from the beginning of each titration step is provided in Table 2. The CHX molecules are numbered by the order in which they insert into the membrane. There is general lengthening in the amount of time taken for each molecule, or group of molecules, to insert for each subsequent titration step. This was expected, as the ΔGi was seen to decrease for CHX in a 10:100 CHX:DMPC membrane as compared to the pure DMPC membrane. Directly analyzing the insertion times of each molecule is problematic due to the aggregation of the molecules in the water. Table 2 shows that, during each of the titration steps, there are groups of CHX molecules that aggregate in the water before inserting into the membrane at the same time. It is also noticeable that the larger the group of CHX, the longer it took for that group to be inserted. Therefore, it is possible that the insertion time is a function of both the CHX concentration in the membrane and the cluster size of CHX trying to insert into the membrane. To deconvolve these two variables, the aggregation would have to be controlled to observe the increase in insertion time for an aggregate of a set size. Our experiment does not have the proper statistical confidence due to the simulations having aggregates of differing size from 1 to 7 CHX molecules, which therefore does not give enough evidence for any one aggregate size. Computational Calorimetry. The specific heat at constant pressure (Cp) was calculated theoretically through “computational calorimetry”, where, by progressively reducing the temperature, we observed anomalously large values for the calculation of Cp at the phase transition. The calculated Cp for a

Figure 7. Time evolution of the position of five CHX1+ molecules from the initial titration step. Top: The distance in the bilayer normal direction from the CHX1+ center of mass to the bilayer center, showing different insertion times of CHX aggregates. Bottom: Intermolecular CHX−CHX distance of CHX1+ molecules showing aggregate breakup after membrane insertion.

the distance in the bilayer-normal direction between the center of the membrane and two clustered groups of CHX1+ molecules numbered 5 (aqua), 6 (blue), 7 (purple), 8 (red), and 9 (orange) from the first titration step. CHX1+ 5, 6, and 7 are seen to form a cluster of three, and CHX1+ 8 and 9 form a cluster of two, moving in tandem. The first cluster forms during the first 20 ns, and around 75 ns they all insert simultaneously, while the second cluster is in the water from 20 to 340 ns before inserting. Figure 7B follows the minimum inter-CHX1+ distance to see how tightly they are clustered. For the first group, they are less than 0.5 nm of one another while in the water, but soon after they enter the membrane at 75 ns, they separate from each other. The same effect is seen with the second cluster, which separate immediately after insertion at 340 ns. The absolute time frame of insertion should not be relied upon, due to how CG effectively samples the energy landscape more frequently during the simulated time period. This speedup is not easily knowable and is molecule and water-model dependent. Regardless of the true effective time scale, the real reduced insertion rate (lengthening insertion time) is simply a reflection of changes to the insertion potential of mean force as the concentration of inserted CHX increases. After the insertion of all 90 CHX molecules, there was an asymmetric distribution of 41 CHX in the lower leaflet and 49 CHX molecules in the upper leaflet. Following an equilibrating period, the two-dimensional CHX1+−CHX1+ radial distribution function in the plane of the bilayer was calculated and is shown in Figure 8A. The radial distribution function rises gradually to 1 and does not show any correlation in the nearest neighbor distance. The CHX molecules are spread evenly through the 4663

DOI: 10.1021/acs.jpcb.6b12510 J. Phys. Chem. B 2017, 121, 4657−4668

Article

The Journal of Physical Chemistry B

Table 2. Time (in ns) for Each CHX Molecule to Insert into the Membrane, Listed by Titration Number (T-1 to T-15)a CHX CHX CHX CHX CHX CHX CHX CHX CHX a

1 2 3 4 5 6 7 8 9

T-1

T-2

T-3

T-4

T-5

T-6

T-7

T-8

T-9

T-10

T-11

T-12

T-13

T-14

T-15

2 26 26 67 75 75 75 342 342

6 13 588 588 588 588 588 588 588

5 35 35 35 35 394 394 394 394

5 10 190 387 387 387 387 387 387

11 41 41 41 41 43 46 364 364

3 9 18 30 50 840 840 840 840

9 13 16 22 597 597 597 597 597

86 86 86 3100 3110 3120 3120 3120 3120

10 997 1010 1010 1010 1010 1010 1010 1010

12 16 16 42 640 1020 1020 1020 1020

1320 1325 1325 1330 1330 1330 1330 1330 1330

25 30 340 340 2320 2320 2320 2320 2320

40 70 150 360 1540 1540 1570 1570 1570

4 20 462 462 462 1132 1132 1132 1132

12 12 2240 2240 2240 2240 2240 2240 2240

CHX molecules are numbered according to the order in which they inserted.

given simulation at any temperature generated a consistent Cp value. However, we observed that the simulation that contained the main lipid transformation contained a large decrease in energy which creates an anomalously high Cp value. Figure 9

Figure 9. Top: System energy as a function of time for the pure DMPC membrane at temperatures of 10.5 °C (black), 10 °C (red), and 9.5 °C (green), as it is cooled through the main transition from liquid to gel. Middle: Specific heat versus simulated temperature. Arrows show the three points of the graph that are calculated from the three energy sections shown above.

Figure 10. Top: Energy versus time for different simulation time lengths at a constant temperature window of 25, 50, and 100 ns (shifted vertically for clarity). Each temperature range begins at 10.5 °C and is reduced by steps of 0.5 °C. Each temperature step is colored differently beginning at 10.5 °C in black at t = 0. Bottom: The area per lipid calculated from the 25 (black), 50 (red), and 100 ns (green) length windows.

shows the resultant Cp versus temperature calculated for the pure DMPC membrane at three different temperatures. The potential energy fluctuations within each temperature window give rise to the calculated Cp values around the peak at the phase transition. The calculated Tm for the pure CG DMPC membrane is ∼10 °C, well below the experimental value, and a well-known drawback of CG models. Nevertheless, the shift of Tm as a function of concentration of CHX is exactly correct with experiments (see below). We investigated the effects of simulation length at each temperature on Cp. Figure 10A shows the energy versus time profile of the cooling using simulation lengths of 25, 50, and 100 ns per simulation temperature window. First, we note that the transition does not always happen at the same temperature, suggesting a degree of error in the temperature of the transition of ±1 temperature window (±0.5 °C). This is likely due to the statistical nature of generating a nucleation site, with the probability increasing with decreased temperature. Second we note that the drop in system energy takes approximately 30 ns. The 25 ns temperature windows are not long enough to capture the transition, and thus the transition was observed to be spread over a large range of windows as the energy decreases as a decaying exponential. The 50 and 100 ns temperature windows are long enough to contain nearly the whole transition

and create a sigmoidal decrease in the system energy, contained in a single simulation window. Figure 10B shows the area per lipid as a function of temperature for each of these three time lengths. For the 25 ns windows, we see the transition from high to low area per lipid takes approximately 5 simulation windows before leveling off at 48 nm2. Similarly we see it takes 2 windows at 50 ns windows, and one point is slightly short enough for 100 ns. This tells us that, although the main transition takes 30 ns, all three time steps take 100−125 ns of equilibration for the system to settle into a steady area per lipid value corresponding to the gel phase. To verify these results experimentally, we performed exothermic DSC scans on CHX:DMPC concentrations of 0:100 up to 50:100. The main transition Tm upon cooling was observed at 20 °C for the pure DMPC sample. This was then chosen as the reference temperature point to compare the rest of the data in which the main transition is measured in terms of T − Tm. A lowering of the main transition temperature as well as a decreasing peak height is observed with increased CHX concentration (see Figure 11A). While a steady lowering of the transition temperature is observed, a jump of approximately 5 °C was seen to occur between the concentrations of 5:100 and 10:100. 4664

DOI: 10.1021/acs.jpcb.6b12510 J. Phys. Chem. B 2017, 121, 4657−4668

Article

The Journal of Physical Chemistry B

MARTINI force field used in this study could not be expected to exactly reproduce exact temperatures.11,12 In Figure 11B, the transition temperature shifts to lower values with increasing CHX concentrations. The trend at which the transition temperature lowers seems to have two distinct regions. Between CHX concentrations of 0:100 to 5:100, we see a linear decrease in Tm of approximately 1 °C per step increase in CHX concentration. This trend is amplified at concentrations of 10:100 and higher, where we see a lowering of the transition temperature of approximately 4 °C per increase of CHX. The height of each peak also decreases with increasing CHX concentrations from 0:100 to 5:100. At concentrations of 10:100 and higher, the peak begins to span more than one temperature simulation and shrinks in size, comparable to background noise. Curiously enough, these two regions seem to be divided around the same CHX concentrations in both simulations and experimental data, and correspond to the previously noted change in area per lipid (cf. Figure 5). The mixture of cholesterol and CHX in the membranes has resulted in nontrivial results. The transition temperature and peak height are plotted in Figure 12 with the experimental data shown with open symbols with connecting lines, and calculated data with filled symbols. The addition of cholesterol does not change the downward trend of the transition temperature and peak height. We do, however, note instances where cholesterol can either reduce the effects of CHX or enhance them. The effects, however, are not straightforward. At 2% cholesterol, for example, we see the effects of CHX hindered at low CHX concentrations and yet aided at high concentrations. Regardless, the relative shifts in transition temperature and change in enthalpy are very accurately modeled by the CG simulations, out to 10:100 CHX:DMPC. Beyond this point, two issues come into play. First, the transition of liquid to gel phase has nearly been erased by the CHX in the CG simulation. Its not clear whether there would be any gel phase with additional CHX. This may be a moot point, as, second, the calculation time required to add another titration step of CHX becomes very prohibitive.

Figure 11. Exothermic plot for differential scanning calorimetry of DMPC:CHX systems with 0% Chol (A,B) and 2% Chol (C,D). A and C were obtained experimentally for DMPC:CHX concentrations of 0:100, 0.3:100, 0.5:100, 1:100, 2:100, 3:100, 5:100, 10:100, 15:100, 20:100, and 50:100 (0% Chol exclusively). B and D were calculated through simulation with 7:100 and 10:100 DMPC:CHX concentration multiplied by a factor of 5.

The concentration scans of CHX:DMPC were repeated in the presence of 1%, 2%, and 5% cholesterol. The effects of cholesterol were minor in that it did not change the overall trend of CHX lowering the transition temperature. Figure 11C shows the sample set for 2% cholesterol in which we note less shift in transition temperature from the 0% cholesterol counterpart before the larger jump between of 5:100 and 10:100 CHX:DMPC. The recorded temperature and height of each peak are plotted in Figure 12 (solid symbols and lines).



DISCUSSION CHX is an effective antibacterial agent, shown to attack, disrupt, and cause cell lysis.5,36−38 We therefore decided to observe CHX in a model DMPC membrane in which we have prior experimental data,9 where the presence of smectic-multilayers at concentrations as high as 30:100 CHX:DMPC indicate the structural integrity of the bilayer remains intact. Thus, our hypothesis was such that we would not necessarily observe signs of membrane instability at concentrations approaching that of the minimum inhibitory concentration (MIC) but perhaps observe a correlation to the intermolecular positioning or interactions between the CHX molecules. Although we have not seen any major disruption in bilayer integrity at these high concentrations of CHX, the CG and AA models agree with experiments in the remarkable accommodation of CHX. This study demonstrates the difficulty in simulating realistic conditions of membrane-active antibiotics. For example, antimicrobial peptides (AMP) exhibit a concentration-dependent effect for their activity that is prohibitive to simulate atomically or as coarse grained, except for a few peptide at a time.39 Nevertheless, the end state of AMP action (typically trans-membrane helix) is sufficiently distinct from its membrane interaction at low concentrations (typically absorbed parallel to the membrane) that a model of

Figure 12. Top: The main phase transition temperature of all CHX:CHOL:DMPC systems in reference to pure DMPC (Tm). Bottom: Relative peak height from the exothermic DSC data with respect to the height of pure DMPC. Experimental values for 0, 1, 2, and 5% cholesterol are presented as solid lines with open symbols, while simulated data appear as solid symbols.

We have reproduced the DSC experiments through simulation as described in the Methods section. The results obtained follow a strikingly similar pattern to the experimental results. The main transition temperature of the DMPC reference simulation occurs at 10 °C. This absolute temperature is of little relevance to us since it has been noted that the 4665

DOI: 10.1021/acs.jpcb.6b12510 J. Phys. Chem. B 2017, 121, 4657−4668

Article

The Journal of Physical Chemistry B action can be determined.40−42 For CHX and related small molecules whose inhibitory effects are also at high concentrations, but molecular interactions are not known, it is insufficient to study bilayer−drug interactions one molecule at a time. However, as we have seen, reaching concentrations in simulation approaching the MIC is a difficult proposition. Assuming a typical rod-shaped Pseudomonas aeruginosa cell to be a cylinder of length 2.5 μs and radius 0.5 μs, with its surface completely composed of lipopolysaccharide (LPS) molecules of area per lipid 1.37 μs2,43 there are roughly 1.4 × 106 LPS molecules covering the surface. This becomes 3 × 106 if such bacteria were covered in DMPC at 0.6 nm2 area per lipid. The actual number lies between theses values since the outer membrane is only composed of approximately 75% LPS molecules with additional, smaller phospholipids.44 Studies to establish MIC values start with a biofilm of as much as 1012 colony forming units (cfu)/L, and if we assume each cfu is established by one cell, then LPS has a concentration in the biofilm of ∼2.3 × 10−6 mol/L (∼5 × 10−6 mol/L for DMPC only).45 This will form a lower limit of lipid concentration, as colonies are more often established by groups of cells. We also take the MIC of of CHX to be 0.5 mg/L, which of course is species dependent, and the molecular mass to be 505 g/mol.5,36,38 Thus, the concentration of CHX at MIC is 1 × 10−6 mol/L, and thus the CHX:lipid ratio at the MIC is roughly in the range of 1:2−8, or as we report in this paper, 12−50:100, assuming all the CHX is taken up by the biofilm. Our results are on the low end of this range. Even at concentrations approaching the MIC, we see no noticeable disruption of the DMPC model membrane, in agreement with our earlier experimental results. The small changes in membrane thickness, and little change in membrane fluidity and order, all point to a cell membrane that remains functionally intact. Unless there is some unidentified target of CHX that is cell membrane bound, there is no sign that even charged CHX would pass through the cell membrane to an interior compartment. We find it unlikely the interaction of CHX with LPS would change these results, and should CHX pass through the outer membrane of Gram-negative bacteria such as P. aeruginosa, it would certainly be easily incorporated in the inner plasma membrane as shown here. One further possibility is that the observed lytic properties of CHX simply cannot be reproduced easily with a simplistic model membrane such as DMPC. In our study we see an increase in the amount of time needed for the insertion of large aggregates of CHX molecules. For a single CHX molecule, we see little to no energy barrier for the insertion from water into the membrane (Figure 3). This suggest that perhaps there is an energy barrier for aggregate sizes greater than one molecule that could act as a function of aggregate size. This energy barrier would also become harder to overcome with the lowering affinity as we saturated our membrane. Unfortunately, we lack the statistical accuracy to make any valuable conclusions on the matter from the study in this work. Lin et al.46 have studied the fusion of micelles of differing sizes with the umbrella sampling technique in order to better study the fusion of micelles of different size. Since we also note a large effect of aggregation of CHX molecules, a combination of these two techniques would be a useful extension to this study. This would include finding the most energetically favorable aggregate size and inserting different aggregate sizes into a model membrane allowing for the observation of the differences in the energy barrier for different

sizes of aggregate inserting into a membrane. A further extension for this study would be to test the energy barrier of these aggregates in membranes of different lipid compositions, which better mimic both a mammalian and bacterial cell membrane, which may lead to insight into why CHX is effective at killing bacteria, but is not as effective in killing mammalian cells.



CONCLUSION In this paper, we presented an extension to the MARTINI force field for the small molecule chlorhexidine (CHX). This allows for an approximate 2 orders of magnitude computational savings when looking for large scale equilibrium structural differences due to high concentrations. This was vital to this study, where we needed to speed-up the ever increasing insertion time, and which contains a total 47.5 μs of production simulation time (36 μs of umbrella sampling and 11.5 μs of simulated titration), on a system of 900 lipids and up to 135 additional CHX molecules. The CHX model was built from our previous all atom (AA) CHX model. The MARTINI particle types were chosen so that the depth in the membrane of each section of the CG CHX molecule matched that of the AA. Bonded parameters determined systematically by adjusting the parameters until the distribution of bonds, angles, and dihedrals of the CG model matched those of the AA. Three different charged CHX models were made: CHX0+, CHX1+, and CHX2+. We then used the CHX1+ model to simulate a titration adding 1% concentration of CHX into a DMPC membrane 15 times so that we had 16 equilibrated systems with CHX concentrations ranging from 0 to 15:100. With these systems, we were able to observe the concentration effects of CHX on the DMPC membrane. The DMPC membrane was easily able to absorb all the CHX molecules, linearly increasing the area per lipid while the thickness remained constant, with the exception of a sudden drop at 6% concentration. Although 15:100 concentration of CHX was not enough to saturate the membrane, a decrease in the affinity between the drug and lipid was observed, as well as a lowered energy barrier for additional molecules to penetrate through the membrane. No internal ordering was found for either the CHX or DMPC molecules. The CHX molecules would aggregate in the aqueous solution, and separate from each other once inserted into the membrane.



ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcb.6b12510. Further parametrization details and model validation are provided, including probability distributions of bonded parameters as well as conformal information on chlorhexidine in both aqueous and membrane environments. (PDF)



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Brad Van Oosten: 0000-0002-2165-2572 Drew Marquardt: 0000-0001-6848-2497 Thad A. Harroun: 0000-0001-9816-2590 4666

DOI: 10.1021/acs.jpcb.6b12510 J. Phys. Chem. B 2017, 121, 4657−4668

Article

The Journal of Physical Chemistry B Notes

from Computer Simulations and Comparison to Experiment. J. Gen. Physiol. 2007, 129, 371−377. (19) Wassenaar, T.; Ingolfsson, H.; Bockmann, R.; Tieleman, P.; Marrink, S. J. Computational Lipidomics With Insane: A Versatile Tool for Generating Custom Membranes for Molecular Simulations. J. Chem. Theory Comput. 2015, 11, 2144−2155. (20) Bekker, H.; Berendsen, H.; Dijkstra, E.; Achterop, S.; Vondrumen, R.; Vanderspoel, D.; Sijbers, A.; Keegstra, H.; Renardus, M. In Physics Computing ’92; DeGroot, R., Nadrchal, J., Eds.; World Scientific Publishing: Singapore, 1993; pp 52−256. (21) Nose, S. A Unified Formulation of the Constant Temperature Molecular Dynamics Methods. J. Chem. Phys. 1984, 81, 511−519. (22) Hoover, W. Canonical Dynamics: Equilibrium Phase-Space Distributions. Phys. Rev. A: At., Mol., Opt. Phys. 1985, 31, 1695−1697. (23) Parrinello, M.; Rahman, A. Polymorphic Transitions in Single Crystals: A New Molecular Dynamics Method. J. Appl. Phys. 1981, 52, 7182−7190. (24) Yesylevskyy, S.; Schafer, L.; Sengupta, D.; Marrink, S. Polarizable Water Model for the Coarse-Grained MARTINI Force Field. PLoS Comput. Biol. 2010, 6, e1000810. (25) Jambeck, J.; Lyubartsev, A. Derivation and Systematic Validation of a Refined All-Atom Force Field for Phosphatidylcholine Lipids. J. Phys. Chem. B 2012, 116, 3164−3179. (26) Jambeck, J.; Lyubartsev, A. An Extension and Further Validation of an All-Atomistic Force Field for Biological Membranes. J. Chem. Theory Comput. 2012, 8, 2938−2948. (27) Jambeck, J.; Lyubartsev, A. Another Piece of the Membrane Puzzle: Extending Slipids Further. J. Chem. Theory Comput. 2013, 9, 774−784. (28) Jo, S.; Kim, T.; Iyer, V.; Im, W. CHARMM-GUI: A Web-Based Graphical User Interface for CHARMM. J. Comput. Chem. 2008, 29, 1859−1865. (29) Kumar, S.; Rosenberg, J. M.; Bouzida, D.; Swendsen, R. H.; Kollman, P. A. The Weighted Histogram Analysis Method for FreeEnergy Calculations on Biomolecules. I. The Method. J. Comput. Chem. 1992, 13, 1011−1021. (30) Roux, B. The Calculation of the Potential of Mean Force Using Computer Simulations. Comput. Phys. Commun. 1995, 91, 275−282. (31) Hub, J. S.; De Groot, B. L.; Van Der Spoel, D. g_wham - A Free Weighted Histogram Analysis Implementation Including Robust Error and Autocorrelation Estimates. J. Chem. Theory Comput. 2010, 6, 3713−3720. (32) Paloncyova, M.; Fabre, G.; DeVane, R. H.; Trouillas, P.; Berka, K.; Otyepka, M. Benchmarking of Force Fields for MoleculeMembrane Interactions. J. Chem. Theory Comput. 2014, 10, 4143− 4151. (33) Luzzati, V.; Husson, F. The Structure of the Liquid-Crystalline Phases of Lipid-Water Systems. J. Cell Biol. 1962, 12, 207−219. (34) Kučerka, N.; Nieh, M.-p.; Katsaras, J. Fluid Phase Lipid Areas and Bilayer Thicknesses of Commonly Used Phosphatidylcholines as a Function of Temperature. Biochim. Biophys. Acta, Biomembr. 2011, 1808, 2761−2771. (35) Kučerka, N.; Van Oosten, B.; Pan, J.; Heberle, F.; Harroun, T.; Katsaras, J. Molecular Structures of Fluid Phosphatidylethanolamine Bilayers Obtained from Simulation-to-Experiment Comparisons and Experimental Scattering Density Profiles. J. Phys. Chem. B 2015, 119, 1947−1956. (36) Thomas, L.; Russell, A. D.; Maillard, J. Antimicrobial Activity of Chlorhexidine Diacetate and Benzalkonium Chloride Against Pseudomonas Aeruginosa and its Response to Biocide Residues. J. Appl. Microbiol. 2005, 98, 533−543. (37) Athanassiadis, B.; Abbott, P.; Walsh, L. The Use of Calcium Hydroxide, Antibiotics and Biocides as Antimicrobial Medicaments in Endodontics. Aust. Dent. J. 2007, 52, S64−S82. (38) Solmaz, G.; Korachi, M. Inhibition and Disruption Properties of Chlorhexidine Gluconate on Single and Multispecies Oral Biofilms. Jundishapur J. Microbiol. 2012, 6, 61−66.

The authors declare no competing financial interest.



ACKNOWLEDGMENTS



REFERENCES

This work was supported by the Engineering Research Council of Hierarchical Academic Research (SHARCNET: www.sharcnet.ca), and da.

Natural Sciences and Canada, the Shared Computing Network Compute/Calcul Cana-

(1) Gilbert, P.; Moore, L. E. Cationic Antiseptics: Diversity of Action Under a Common Epithet. J. Appl. Microbiol. 2005, 99, 703−715. (2) Sreenivasan, P.; Gaffar, A. Antiplaque Biocides and Bacterial Resistance: A Review. J. Clin. Periodontol. 2002, 29, 965−974. (3) Chen, L.; Bromberg, L.; Hatton, T. A.; Rutledge, G. C. Electrospun Cellulose Acetate Fibers Containing Chlorhexidine as a Bactericide. Polymer 2008, 49, 1266−1275. (4) Sinisterra, R. D.; Faris, C. B.; Yue, I. C.; Poff, J.; Cortes, E.; Hildgen, P.; Langer, R.; Shastri, V. P. A Novel Polymeric Chlorhexidine Delivery Device for the Treatment of Periodontal Disease. Biomaterials 2004, 25, 3743−3750. (5) Fraise, A. P., Jean-Yves Maillard, S. S. Principles and Practices of Disinfection Preservation & Sterilization, 4th ed.; Blackwell Publishing Ltd: Malden, MA, 2004. (6) Moug, T. E.; Rogers, D. T.; Furr, J.; El-Falaha, B.; Russell, A. Antiseptic-Induced Changes in the Cell Surface of a ChlorhexidineSensitive and a Chlorhexidine-Resistant Strain of Providencia Stuartii. J. Antimicrob. Chemother. 1986, 16, 685−689. (7) Hiom, S.; Furr, J.; Russell, A.; Dickenson, J. Effects of Chlorhexidine Diacetate and Cetylpyridinium Chloride on Whole Cells and Protoplasts of Saccharomyces Cerevisiae. Microbios 1993, 111−20. (8) Barrett-Bee, K.; Newboult, L.; Edwards, S. The Membrane Destabilising Action of the Antibacterial Agent Chlorhexidine. FEMS Microbiol. Lett. 1994, 119, 249−253. (9) Komljenović, I.; Marquardt, D.; Harroun, T. a.; Sternin, E. Location of Chlorhexidine in DMPC Model Membranes: A Neutron Diffraction Study. Chem. Phys. Lipids 2010, 163, 480−487. (10) Van Oosten, B.; Marquardt, D.; Komljenović, I.; Bradshaw, J. P.; Sternin, E.; Harroun, T. A. Small Molecule Interaction With Lipid Bilayers: A Molecular Dynamics Study of Chlorhexidine. J. Mol. Graphics Modell. 2014, 48, 96−104. (11) Marrink, S.; de Vries, A.; Mark, A. Coarse Grained Model for Semiquantitative Lipid Simulations. J. Phys. Chem. B 2004, 108, 750− 760. (12) Marrink, S. J.; Risselada, H. J.; Yefimov, S.; Tieleman, D. P.; de Vries, A. H. D. The MARTINI Force Field: Coarse Grained Model for Biomolecular Simulations. J. Phys. Chem. B 2007, 111, 7812−7824. (13) Van Oosten, B.; Harroun, T. A MARTINI Extension for Pseudomonas Aeruginosa PAO1 Lipopolysaccharide. J. Mol. Graphics Modell. 2016, 63, 125−133. (14) Monticelli, L.; Kandasamy, S. K.; Periole, X.; Larson, R. G.; Tieleman, D. P.; Marrink, S. J. The MARTINI Coarse-Grained Force Field: Extension to Proteins. J. Chem. Theory Comput. 2008, 4, 819− 834. (15) Smith, D. A.; Jones, B. C.; Walker, D. K. Design of Drugs Involving the Concepts and Theories of Drug Metabolism and Pharmacokinetics. Med. Res. Rev. 1996, 16, 243−266. (16) Barton, P.; Davis, A. M.; McCarthy, D. J.; Webborn, P. J. H. Drug-Phospholipid Interactions. 2. Predicting the Sites of Drug Distribution Using n-Octanol/Water and Membrane/Water Distribution Coefficients. J. Pharm. Sci. 1997, 86, 1034−1039. (17) De Young, L. R.; Dill, K. A. Solute Partitioning into Lipid Bilayer Membranes. Biochemistry 1988, 27, 5281−5289. (18) MacCallum, J. L.; Bennett, W. F. D.; Tieleman, D. P. Partitioning of Amino Acid Side Chains into Lipid Bilayers: Results 4667

DOI: 10.1021/acs.jpcb.6b12510 J. Phys. Chem. B 2017, 121, 4657−4668

Article

The Journal of Physical Chemistry B (39) Perrin, B. S., Jr.; Pastor, R. W. Simulations of MembraneDisrupting Peptides I: Alamethicin Pore Stability and Spontaneous Insertion. Biophys. J. 2016, 111, 1248−1257. (40) Ludtke, S. J.; He, K.; Heller, W. T.; Harroun, T. A.; Yang, L.; Huang, H. W. Membrane Pores Induced by Magainin. Biochemistry 1996, 35, 13723−13728. (41) Heller, W. T.; He, K.; Ludtke, S. J.; Harroun, T. A.; Huang, H. W. Effect of Changing the Size of Lipid Headgroup on Peptide Insertion Into Membranes. Biophys. J. 1997, 73, 239−244. (42) Yang, L.; Harroun, T. A.; Weiss, T. M.; Ding, L.; Huang, H. W. Barrel-Stave Model or Toroidal Model? A Case Study on Melittin Pores. Biophys. J. 2001, 81, 1475−1485. (43) Kirschner, K. N.; Lins, R. D.; Maass, A.; Soares, T. A. A GlycamBased Force Field for Simulations of Lipopolysaccharide Membranes: Parametrization and Validation. J. Chem. Theory Comput. 2012, 8, 4719−4731. (44) Caroff, M.; Karibian, D. Structure of Bacterial Lipopolysaccharides. Carbohydr. Res. 2003, 338, 2431−2447. (45) Fraise, A.; Maillard, J.; Sattar, S. Russell, Hugo and Ayliffeas Principles and Practice of Disinfection, Preservation and Sterilization [electronic resource]; Wiley, Hoboken, NJ, 2012. (46) Lin, D.; Grossfield, A. Thermodynamics of Micelle Formation and Membrane Fusion Modulate Antimicrobial Lipopeptide Activity. Biophys. J. 2015, 109, 750−759.

4668

DOI: 10.1021/acs.jpcb.6b12510 J. Phys. Chem. B 2017, 121, 4657−4668