Subscriber access provided by University of Sunderland
B: Biomaterials and Membranes
Membrane Permeability of Terpenoids Explored with Molecular Simulation Josh V. Vermaas, Gayle J Bentley, Gregg T Beckham, and Michael F. Crowley J. Phys. Chem. B, Just Accepted Manuscript • DOI: 10.1021/acs.jpcb.8b08688 • Publication Date (Web): 23 Oct 2018 Downloaded from http://pubs.acs.org on October 27, 2018
Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.
is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.
Page 1 of 45 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
Membrane Permeability of Terpenoids Explored with Molecular Simulation Josh V. Vermaas,† Gayle J. Bentley,‡ Gregg T. Beckham,∗,‡ and Michael F. Crowley∗,† †Biosciences Center, National Renewable Energy Laboratory, Golden, CO 80401 ‡National Bioenergy Center, National Renewable Energy Laboratory, Golden, CO 80401 E-mail:
[email protected];
[email protected] 1
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Abstract Terpenoids constitute a class of compounds with remarkable potential for pharmaceutical, fragrance, specialty chemical, and biofuel applications. However, their industrial production is limited by their rarity within their native plant hosts, creating considerable interest in microbial hosts capable of manufacturing terpenoids. To reduce production costs, non-destructive product recovery from these microbial hosts is preferred, and is achievable using a hydrophobic organic overlay. Our prior research has indicated that oxidized fatty acyl products may permeate faster through host membranes, increasing overall biorefinery productivity. To test this hypothesis for terpenoids, we computed membrane permeabilities of conventional terpenoid target products (e.g. limonene, bisabolene, farnesene), and related oxidized compounds through molecular dynamics simulations. These simulations indicate that terpenoid product permeabilities from cytosol to overlay are oxidation independent, as increases in membrane extraction efficiency due to product oxidation are proportionally offset by decreases in the membrane crossing rate if the membrane and organic phase are in close contact. However, if aqueous extraction is required, oxidation will accelerate the slow product extraction from the membrane. Experimental toxicity assays performed indicated that most terpenoids tested were tolerated by microbial hosts, although exposure to oxidized terpenes often retarded microbial growth compared with conventional terpenes. Thus, terpenoid oxidation is not expected to significantly increase or decrease the extraction productivity in an industrial setting where cells are in close contact, unlike the previously studied fatty acyl products.
2
ACS Paragon Plus Environment
Page 2 of 45
Page 3 of 45 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
Introduction Terpenes and terpenoids, naturally-occurring compounds primarily extracted from plants, display diverse functionality which allow for their use in a wide variety of chemical applications. 1,2 Still, it remains difficult to rely on plant-based extraction for the large-scale production of terpenes naturally produced in low quantities from plants or animals that are difficult to cultivate, or are derived from endangered species. Compounded with these limitations is crop yield variability, further affecting the yield and types of terpenes extracted by plant-based methods. 3 In response to these limitations, many research groups have successfully engineered microbes for the production of an increasingly wide range of terpenes. 4–10 Additionally, engineering microorganisms for terpene production increases the space of terpene chemistries that are accessible for overproduction due to the abundance of terpene synthases that enable unique chemistries. 11 From a process perspective, toxicity and ease of extraction of the terpene product are critical factors that affect the feasibility of scalable and cost effective production. Terpenes can function as antimicrobial agents to protect their natural hosts, with antibacterial activity occurring via disruption of the lipid membrane. 12–14 The use of a hydrophobic organic overlay such as dodecane is often used at bench scale during microbial production of terpenes to eliminate loss of a volatile product and to reduce the immediate toxicity of the terpene product. 9,10,15–17 Here, we examine the effect of terpene oxidation on cellular toxicity and extraction efficiency of the terpene or terpenoid into an organic overlay. Since terpenoid extraction faces the same challenges as extraction of other biological products, it may be informative to apply findings from other systems to terpenoids. Using molecular simulation, our recent results have indicated that hydroxyl groups maximize fatty acid extraction into a hydrophobic organic phase. 18 Fatty alcohol and fatty aldehyde products were shown to efficiently cross lipid bilayers in addition to reducing the barrier to extraction of these typically lipophilic molecules. Terpenoid target molecules are generally devoid of oxygenation and are highly lipophilic, although recent work has also shown pathways towards 3
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 4 of 45
creating oxidized terpenoids, 19–24 and so there is a significant potential for improvement if the trend observed for fatty acyl products 18 is similar for terpenoids.
OH
OH
Pinene
O
Verbenol Verbenone Farnesene HO
Bisabolene
Bisabolol
Farnesol OH
Limonene
Terpineol
OH
Perillyl alcohol
Figure 1: Terpenoid compounds and their oxidized equivalents considered in this study. The colored boundaries group together equivalent chemical functionalities, which are used consistently in other visualizations throughout the text. Red highlights bisabolene and its predominant oxidized form, bisabolol, blue highlights pinene and its derivatives verbenol and verbenone, green highlights limonene with derivatives terpineol and perillyl alcohol, and violet surrounds farnesene and its oxidation product farnesol. To test whether hydroxylation of terpenoids improves the extraction rate from cell membranes to an organic overlay, in this study we have conducted molecular dynamics (MD) simulations to determine the membrane permeabilities and expected membrane extraction rates of both oxidized and non-oxidized terpenoid compounds (Fig. 1). The oxidized terpenoids studied here have biological synthetic routes, with examples of metabolic engineering to produce farnesol, 19,20 bisabolol, 21 verbenol, 22 verbenone, 23 terpineol, 24 and perillyl alcohol. 25 The MD findings reflect the previous trend observed for fatty acyl products, in that adding a hydroxyl group improves membrane extraction rates while retarding membrane 4
ACS Paragon Plus Environment
Page 5 of 45 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
crossing. Since the determined crossing and extraction rates are so similar, we predict that terpenoid oxidation does not significantly improve the theoretical extraction rate of organic phase extraction if the cell and organic phase are in close contact. However, if cells avoid organic phases through an active mechanism, as has been suggested previously based on microscopy measurements, 26 our results demonstrate that oxidized terpenoids may improve extraction rates if the cell’s immediate environment is predominantly aqueous. In addition, we find through microbial growth experiments that most terpenoid products are tolerated at current experimental titers of up to 3%, although some oxidized products exhibited greater toxicity, resulting in decreased relative to their lipophilic equivalents. Thus, targeting oxidized terpenoids may be a viable path forward in optimizing microbial production of these species, depending on the environmental conditions within industrial settings.
Methods Both computational and experimental assays were conducted. The molecular simulations were used to determine the permeability of these compounds. Experimental assays measure the impact of alternate terpenoids on growth. The computational methods used to assess the bilayer crossing and organic phase partitioning of terpenoids are largely derived from prior studies on fatty acyl products 18 and the wider membrane permeability literature, 27 which details the theory behind computing permeability coefficients separately for extraction and crossing steps. In brief, equilibrium MD simulations were performed to monitor spontaneous membrane crossing of the terpenoid compound present. These equilibrated systems formed the basis for Replica Exchange Umbrella Sampling (REUS) 28 calculations used to determine local diffusivity and free energy profile for membrane crossing and extraction processes for all of the terpenoids indicated in Fig. 1. The diffusivity and free energy profiles are used to
5
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 6 of 45
compute membrane permeability (P) using the inhomogeneous solubility diffusion model: 29,30
P
−1
Z
ξu
= ξl
exp (∆G (ξ) β) dξ D (ξ)
(1)
In this formalism, the free energy change (∆G) and diffusivity (D) are known along the reaction coordinate (RC) ξ, which contributes to the permeability along the given RC. The details of the simulation system preparation are described in the subsequent subsections. To understand what the practical effect these compounds might have on their toxicity and growth of the host organism, which is difficult to assess computationally, growth assays were also carried out, with their methods detailed at the end of this section.
Simulation System Setup Two types of simulation systems are required to accurately describe the terpenoid extraction process, one with the product in aqueous solution (Fig. 2A) representing the terpenoid adsorption to the membrane, and a second (Fig. 2B) representing terpenoid extraction from the membrane into a proximal organic phase, in this case dodecane. We consider this proximal organic phase to be a better descriptor of how dodecane interacts with cells at experimental dodecane concentrations, as scenarios where the dodecane either is sandwiched between the membrane leaflets or extracts significant lipids from the membrane would likely lead to cell death by disrupting the membrane. Both the aqueous and organic systems share a common membrane composition with the comparison study of fatty acyl product extraction to aid comparison, 18 with each leaflet having a 30:14:12:9:4:1 molecular composition of phosphatidyl choline (PC), phosphatidyl ethanolamine (PE), phosphatidyl inositol (PI), ergosterol, phosphatidyl serine (PS), and phosphatidic acid (PA) headgroups and sterols and a 2:1:1 mixture of palmityl oleyl, stearyl-linoleyl and dioleyl lipids. This composition was initially chosen as a stand-in for an oleaginous yeast membrane based on the available experimental evidence for both headgroups 31–33 as well as tail compositions. 34,35 The membrane
6
ACS Paragon Plus Environment
Page 7 of 45 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
B
A
Figure 2: Example states near the start of simulation, detailing the relative positioning and abundance of the modeled dodecane overlay (tan) compared with the membrane (atomicallycolored spheres), ions (yellow for Na+ , cyan for Cl – ) and the single terpenoid product per simulation system (green). Hydrogen atoms are omitted for clarity, and the surrounding water is represented as a transparent blue surface where solution is present. (A) The overlayfree system is used as a control. (B) A high-dodecane overlay, representative of a lab-scale overlay in close proximity with the membrane. In addition to the main findings regarding membrane permeability, three of these high-dodecane overlays were considered to predict the membrane structural changes depending on the degree of bilayer hydration. The three initial spacings between the membrane and the organic phase were 5, 10, and 15 ˚ A, with water to lipid ratios of approximately 12, 27, and 38 as the result. Here, we only show the 5˚ A case.
7
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
was assembled using CHARMM-GUI. 36–38 For each terpenoid in Fig. 1, both aqueous and organic systems were initially constructed with a single terpenoid molecule added outside of the bilayer, followed by solvation and ionization using the Solvate and Autoionize plugins of VMD 39 to add water and dodecane as appropriate for the system (Fig. 2). Due to the observed ordering of the lipids under conventional constant pressure simulation that are discussed later, aqueous solvation layers between the membrane and organic phase (Fig. 2B) of different thicknesses were tested, resulting in water to lipid ratios of 12, 27, and 38, depending on the degree of hydration. These three systems cover the range of lipid ratios used in previous simulation studies of lipid structural changes in response to hydration. 40–42 In all instances, the total aqueous salt concentration is 150 mM NaCl, which includes the 30 ion sodium bias that results from the net charge of the bilayer (e.g. there are 34 Na+ and 4 Cl – ions in the organic system with the least water, shown in Fig. 2B). All systems have an initial total size of approximately 65 ˚ A×65 ˚ A×100 ˚ A.
Simulation Details Equilibration for each system and biased simulations used to compute membrane permeability share specific simulation parameters, which are described here. Simulations were performed with NAMD 2.12, 43 using the CHARMM36 force field for lipids 42 and sterols. 44 The dodecane overlay and the terpenoids used parameters from the CHARMM general force field 45 (CGenFF) as determined by the ParamChem webserver. 46 Consistent with the standard parameterization methodology for CHARMM force fields, 45 simulations were performed using a 12 ˚ A nonbonded cutoff and the TIP3P water model. 47 Long-range electrostatic interactions were treated using particle mesh Ewald 48,49 with a 1.2 ˚ A grid. A Langevin thermostat using γ=1 ps−1 maintained the system temperature 50 at either 300 K (equilibration) or 310 K (biased simulation). Anisotropic pressure coupling was maintained via the Langevin piston method 51,52 to a pressure of 1 atm, with periodic cell growth along the membrane normal 8
ACS Paragon Plus Environment
Page 8 of 45
Page 9 of 45 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
axis decoupled from growth along the membrane parallel axes. In most cases, the membrane parallel axes were allowed to move freely; however, to combat observed membrane ordering, one set of REUS calculations was carried out in a fixed area ensemble, where the only dimension controlled by the barostat is along the membrane normal. Each simulation timestep was 2 fs, enabled by using the SETTLE algorithm 53 to fix bond lengths to hydrogen. Equilibrium Simulation Prior to biased simulation, five independent 100 ns equilibrium simulations were performed starting from the structures shown in Fig. 2 to increase sampling of the single terpenoid such that the terpenoids within the subsequent steered MD (SMD) trajectories would not experience the same membrane environment. The 100 ns timescale permits the membrane to adopt its preferred geometry given the imposed simulation conditions and for the terpenoid positions to be thoroughly decorrelated, minimizing the bias from the initial geometry when seeding the initial positions for the REUS calculations. However, it must be emphasized that this timescale is significantly shorter than what is required for significant mixing, clustering, or leaflet exchange of lipids, which have been previously demonstrated to occur on up to the multi-microsecond timescale. 54–56 Further exploration of lipid ordering within our simulated membranes was conducted through additional sets of five independent simulations of each terpenoid for 20 ns each where additional hydration is present, as described in Fig. 2B. Bilayer Crossing and Extraction Preparation For many of the terpenoid compounds shown in Fig. 1, there are comparatively few rotatable bonds, as the ring structures present constrain the conformational space accessible to each compound. Thus, for pinene, bisabolene, limonene, and their derivatives, the center of mass of the molecule accurately describes the progress along the extraction and crossing RC. For these molecules, we prepared a complete trajectory where the compound center of mass ranged from the membrane midplane to the plane 45 ˚ Aaway from the bilayer center. This
9
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 10 of 45
was accomplished by taking the five independent positions resulting from the final snapshot of the equilibrium trajectory, and using SMD to move the molecules to the end points of the RC in two separate simulations (“forward” along the RC and “backward”). The force constant used for this step was 7 kcal mol−1 ˚ A−2 . The total time for the “forward” and “backward” simulations was 10 ns, with the length of each simulation dependent on the distance between the initial state and the final target for the RC. For farnesene and farnesol, however, there is considerable conformational flexibility through rotation of the single bonds, which makes the center of mass an unreliable estimator of extraction progress. 57,58 Similar to the approach taken when studying fatty acyl compounds, 18 we split the RC into two steps, creating independent RC for bilayer adsorption/desorption and bilayer crossing. For bilayer crossing, the assumed RC, measuring the position the end of the molecule where the alcohol would be on farnesol, ranged from -25 ˚ A to 25 ˚ A along the membrane normal, with the zero-point located at the membrane midplane. This RC intentionally overestimates the 40 ˚ A thickness of the bilayer 59 so that we can be sure the free energy minimum is captured by this RC. As was done for the other compounds, “forward” and “backward” simulations were performed to take the end state from the equilibrium trajectories to the boundaries of the RC. These pulling trajectories were performed over 10 ns ˚−2 . in total using a force constant of 2 kcal mol−1 A For extraction of farnesene and farnesol, we compose a RC based on the number of contacts between the single terpenoid compound and the membrane and, if present, the dodecane overlay. This RC is implemented as a coordination number collective variable from the colvars 60 module of NAMD, 43 using the same cutoffs and parameters as were used for fatty acyl products: 18 4 XX 1 − |xi − xj | /8˚ A C (g1 , g2 ) = ˚ 10 i∈g1 j∈g2 1 − |xi − xj | /8A
(2)
The coordination number is defined over the farnesene heavy atoms (g1 ), with g2 defined
10
ACS Paragon Plus Environment
Page 11 of 45 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
as either the even-numbered acyl tail carbons of the lipids, creating a lipid-farnesyl contact number (Clipid ), or as a subset of the dodecane heavy atoms, creating a dodecane-farnesyl contact number (Cdodecane ). For extraction into aqueous solution, we drive Clipid to 0. For extraction into an organic phase, we drive the sum of Clipid − Cdodecane to be as negative as possible. Due to greater direct access to the heavy atoms offered by the unsaturations present in farnesene and farnesol, the coordination numbers during equilibrium simulation tend to be somewhat larger than they were in the fatty acyl products (Fig. S1). Thus, the extraction RC from membrane to aqueous solution ranges from [0,650], and from [-850,650] for extraction to organic solution. To construct continuous pulling trajectories to seed the REUS calculations, each of the 5 end states from the equilibrium trajectory were pulled with “forward” and “backward” simulations with total duration of 10 ns for each replicate. The extraction pulls used a 0.01 kcal mol−1 contact−2 force constant. The five independent contiguous trajectories created by the “forward” and “backward” simulations for each compound were used to seed the initial coordinates for the REUS calculations, with each trajectory contributing an equal number of initial positions randomly distributed across each RC considered. In this way, we also sample from different places within the heterogeneous membrane environment created by the mixed bilayer, rather than simply one potential path of bilayer translocation, increasing the accuracy of the result. Replica Exchange Umbrella Sampling Parameters For each of the three distinct RC sampled (center of mass, farnesene translocation, and farnesene extraction), the REUS simulations themselves are all 20 ns long. Exchanges between alternating replica neighbors were attempted every picosecond, yielding a 2 ps overall exchange rate, consistent with reported best practices suggesting frequent exchange attempts. 61 Other parameters changed depending on the RC under study. For center of mass translocations, the [0,45] ˚ A RC was divided into 60 equal umbrellas with a 5 kcal mol−1 ˚ A−2 force constant ensuring consistent sampling, resulting in a 12% ex-
11
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
change acceptance probability. For farnesene and farnesol bilayer crossing simulations, which were carried out in systems without dodecane solvent present, the [-25,25] ˚ A range of the RC ˚−2 force constant, yielding was covered by 70 equally spaced umbrellas with a 5 kcal mol−1 A a 15% exchange acceptance probability. The spacing and force constants are consistent with other umbrella sampling studies of bilayer translocation, 18,27,62–64 although they do result in somewhat lower exchange acceptance probabilities than the 20 % optimum determined from temperature replica exchange studies. 65 Extracting farnesene or farnesol into either aqueous solution or a dodecane overlay proceeds in much the same manner, with 55 equally spaced replicas along the [0,650] RC into aqueous solution, and 130 equally spaced replicas along the [-850,650] RC describing the extraction into a dodecane overlay. For both extraction RCs, the applied umbrella used a 0.01 kcal mol−1 contact−2 force constant. For farnesene and farnesol, to convert the eventual permeability coefficient into units that are relevant to experiment, we determine the conversion rate between contact number and distance by a functional fit of the relationship between the biased contact RC and the observed displacement, as was done for the fatty acyl products. 18 Since the membrane was observed to become more ordered due to compression of the periodic cell when the organic dodecane phase was present, the center of mass translocations with organic solvent were repeated a second time using a fixed area barostat. During this second set of REUS simulations, the cell dimensions of the system in the membrane plane were fixed to be 60.2 ˚ A, consistent with the periodic cell dimensions seen in purely aqueous membrane solvation environments (Fig. S2). Analysis Analysis of the simulation trajectories was primarily carried out using a combination of builtin and purpose built VMD 39 scripts and features, leveraging the numpy and scipy libraries 66 to handle and manipulate trajectory data. Individual dataplots were created with the matplotlib library, 67 and include particle tracking during equilibration, lipid order parameters,
12
ACS Paragon Plus Environment
Page 12 of 45
Page 13 of 45 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
and naturally also the free energies and diffusivities used to compute permeabilities. Estimation of the free energy profile required for Eq. 1 was performed using a modified version of BayesWHAM, 68 which has been accelerated by using Gibbs sampling of the known Dirichlet prior 69 rather than the slower Metropolis-Hastings originally implemented. 68 The local diffusivity is estimated from the variance and time autocorrelation along the RC, as first formulated by Hummer 70 , and used for fatty acyl products previously. 18 To account for the exchanges introduced by REUS sampling, the time autocorrelation is determined through an exponential fit to the autocorrelation function at timescales less than the exchange frequency. This approach was also used in determining diffusivity for the fatty acyl products, 18 and yields diffusivities in broad agreement with experimentally determined small-molecule diffusivities of 10−6 cm2 s−1 for low molecular weight hydrocarbons in an aqueous environment. 71,72 The estimates for free energy and diffusivity are then numerically integrated to determine the permeability using the formalism from Eq. 1, with the zero level on the free energy surface taken to be the free energy minima nearest to the level of lipid carbonyls. This means that the permeabilities computed are relative to an already inserted compound in one of the two membrane leaflets.
Strains, Media and Chemicals Bacterial Growth Bacterial strains Pseudomonas putida KT2440 and Escherichia coli MG1655 were used for toxicity growth assays. P. putida and E. coli were cultivated in Luria-Bertani (LB) medium (Lennox) containing 10 gL−1 tryptone, 5 gL−1 yeast extract, and 5 gL−1 NaCl for precultures. For growth assays, P. putida and E. coli were cultivated in M9 minimal medium containing 6.78 gL−1 Na2 HPO4 , 3.00 gL−1 K2 HPO4 , 0.50 gL−1 NaCl, 1.66 gL−1 NH4 Cl, 0.24 gL−1 MgSO4 , 0.01 gL−1 CaCl2 , and 0.002 gL−1 FeSO4 , supplemented with 3.60 gL−1 glucose and/or different concentrations of terpenes. Stock solutions of each terpenoid were prepared as 20% v/v in H20 . All the chemicals used for the study were obtained from 13
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Sigma-Aldrich (St. Louis, MO, USA) or Alfa Aesar (Haverhill, MA, USA- bisabolene only). Yeast Growth The yeast strain Saccharomyces cerevisiae Y294 (ATCC 201160) was used for toxicity growth assays. Pre-cultures of S. cerevisiae were grown in yeast-peptone-dextrose (YPD) liquid medium containing 10 gL−1 yeast extract powder, 20 gL−1 peptone powder, 20 gL−1 glucose. For growth assays, S. cerevisiae was cultivated in yeast nitrogen base (YNB) medium containing 2% (v/v) glucose at pH 4.5 (Sigma Aldrich, St. Louis. MO, USA) containing proper auxotrophic supplements (leucine, uracil, histidine, and tryptophan) at 1 gL−1 and/or different concentrations of terpenes.
Growth Assay Toxicity Assays Toxicity of the terpene compounds was evaluated in microplate growth assays performed in a Bioscreen C MBR microplate reader (Growth Curves US, Piscataway, NJ). Pre-cultures were inoculated by scraping strains preserved in 25% v/v glycerol into 5 mL LB (bacteria) or YPD (yeast) medium in 15 mL disposable culture tubes. Pre-culture tubes were incubated shaking at 225 rpm, 30 ◦ C. When strains reached OD600 1–2, cells were harvested by centrifugation at 4,500 rpm, and the cell pellets were washed and resuspended in M9 salts without a carbon source. OD600 were measured using a Beckman DU640 spectrophotometer (Beckman Coulter, Brea CA). These resuspended cells were used to inoculate microplate wells containing 150 µL of M9 (bacteria) medium or YNB (yeast) to an OD600 of 0.1. Medium volume of each well was adjusted based on the volume of terpenoid compound required to reach the tested concentration. After strain inoculation, each of the three microbes was grown three times under each of the different conditions tested. These conditions include a control without terpenoid, as well as 0.5%, 1%, and 3% concentrations of bisabolene, bisabolol, farnesene, farnesol, 14
ACS Paragon Plus Environment
Page 14 of 45
Page 15 of 45 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
limonene, terpineol, and perillyl alcohol. Pinene and its derivatives were not tested due to their comparatively high price. To ensure homogeneity in the inoculum, each terpenoid was vortexed vigorously immediately prior to addition to the microplate well. Microplates were then incubated at 30 ◦ C with maximum shaking and growth was measured by reading the absorbance (OD420−580 ) every 15 min. Growth on Perillyl Alcohol and Terpineol as a Sole Carbon Source We examined the ability for P. putida KT2440 to grow on perillyl alcohol and terpineol as sole carbon sources. Pre-cultures were prepared as above and harvested, washed cells were used to inoculate microplate wells containing 150 µL of M9 medium prepared without glucose to an OD600 of 0.1. Concentrated terpenoids were added to each microplate well as described above. Microplates were then incubated at 30 ◦ C with maximum shaking and growth was measured by reading the absorbance (OD420−580 ) every 15 min.
Results System Evolution at Equilibrium In the prior study of fatty acyl products, individual compounds were initially placed within the membrane. 18 However, as indicated by Fig. 2, here the compounds were initially placed in solution or at the aqueous/organic interface. The behavior of the individual compounds depends on their chemistry and environment. When initially placed in aqueous solution, as in Fig. 2A, the compounds tend to insert into the bilayer (Fig. 3). The insertion process is not uniformly fast. In some equilibrium trajectories, insertion is not realized prior to the end of 100 ns of equilibration (Fig. S3), whereas others insert rapidly, within the 50 ns timescale previously observed for lipid insertion into a bilayer. 57 The observed trend is that the species with the greatest lipophillicity, such as farnesene, limonene, and pinene, are the fastest to insert, with all five 100 ns trajectories leading to a 15
ACS Paragon Plus Environment
The Journal of Physical Chemistry
Bisabolene Bisabolol Limonene Terpineol Perillyl alcohol Pinene Verbenol Verbenone Farnesene Farnesol
Probability
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
0
0 10 20 30 40 50 |Compound-Membrane Center-Center Distance| (˚ A)
Figure 3: Distribution of the position the center of mass of each compound relative to the membrane center over all equilibration trajectories. For each compound, there are two distributions, a solid line for the aqueous equilibration trajectories (Fig. 2A), and a dashed line for the organic equilibration trajectories (Fig. 2B). Similar compounds are grouped by color as in Fig. 1, with more hydrophilic compounds given lighter colors. membrane-embedded compound. The hydroxylated counterparts insert more slowly, with always at least one trajectory failing to insert within 100 ns and remaining in aqueous solution. The compounds also differ in their preferred insertion depth. Hydrogen bonding with lipid carbonyl groups causes the hydroxylated compounds to have a probability peak in the range of 10-15˚ A, with the precise value depending on the size of the molecule, consistent with prior simulation at both atomic 18,73 and coarse-grained resolution. 74 The unmodified terpenoids are all found to be enriched at the membrane center after insertion from aqueous solution (Fig. S3), similar to previous findings for benzene 75 and alkanes. 18 During equilibration trajectories where a dodecane phase was present, similar to Fig. 2B, spontaneous insertion into the bilayer was rarely observed (Fig. S3). Where insertion into the bilayer from the organic phase does occur, as in bisabolene, limonene, and farnesene, insertion appears to be reversible (Fig. S3). The reversibility of the insertion and the overall population balance between organic and membrane phases suggests that the organic phase may be the preferred partition for these terpenoid compounds. Within the organic phase, compounds are observed to partition based on hydrophilicity. Compounds bearing a hydroxyl
16
ACS Paragon Plus Environment
Page 16 of 45
Page 17 of 45 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
group prefer to remain on the interface between water and the organic phase, and compounds without potential hydrogen bond interaction sites are nearly uniformly distributed within the dodecane phase (Fig. 3). From these aggregate results, we can already begin to make predictions regarding the overall shape of membrane crossing and extraction free energy surfaces. Oxidized compounds will cross the bilayer more slowly than do the reduced terpenoids, based on the hydrogen bonds created with the phospholipid head groups. Similarly, the oxidized compounds will extract only as far as the dodecane-water interface, and further extraction would bear an energetic penalty. Based on the free diffusion of reduced terpenoids within the organic phase, no additional free energy cost to bury these compounds entirely within the dodecane phase would be expected.
Free Energy and Local Diffusivity To test these expectations, as well as to quantify the local free energy and diffusivity such that a permeability can be computed using Eq. 1, we analyze the REUS calculations described previously. However, since farnesene and farnesol use a different RC than the other compounds tested due to the additional rotateable bonds present in their molecular structure (Fig. 1), we must present these results separately. Ring Terpenoids We begin with the majority of the compounds, where the center of mass of the compound relative to the membrane center was used as the RC due to their rings restricting the available degrees of freedom (Fig. 4). The composition of Fig. 4 depends on both aqueous and organic extraction REUS calculations. The free energy and local diffusivity profiles are fused together at the midpoint of the membrane under the assumption that the environment experienced by these small molecules is unaffected by the solvent environment at the membrane periphery. This assumption is not violated so long as the organic extraction is performed in a fixed-area 17
ACS Paragon Plus Environment
The Journal of Physical Chemistry
2 0
∆G (kcal/mol)
−2 −4 −6 −8 −10 −12
Bisabolene Bisabolol Limonene Terpineol Perillyl alcohol Pinene Verbenol Verbenone
45 40 35 30 D (˚ A2 /ns)
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 18 of 45
25 20 15 10 5 0
−40
−30 −20 −10 0 10 20 30 Position Relative to Membrane Center (˚ A)
40
Figure 4: Free energy profile and local diffusivity for compact terpenoid compounds. For the free energy profile, the reference point was chosen to be the compound in aqueous solution (left), which transits the bilayer (center), and is eventually extracted into the organic phase (right). Thus, the RC shown here is a composite of the two REUS calculations, joined at the membrane center (z=0), specifically of the extraction to aqueous solution (left), and the extraction to the organic phase in a fixed area ensemble (right). As a visual aide, a simplified atomistic representation of the respective components is used as a background, with water oxygens shown in blue, dodecane heavy atoms in green, and lipid heavy atoms following a standard color scheme (grey for carbons, red for oxygens, blue for nitrogen, and brown for phosphorus).
18
ACS Paragon Plus Environment
Page 19 of 45 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
ensemble, as it was for the results presented in Fig. 4. The most straightforward way to test if the presence of the organic solvent interferes with bilayer structure in the membrane interior is to check the profiles for symmetry around the membrane center. When compared with Fig. S4, where the free energy profile is determined from REUS calculations carried out under a semi-isotropic ensemble typically used with membrane simulations, Fig. 4 has significantly more symmetry in the region near the membrane center. The asymmetry with the semi-isotropic barostat comes from lipid ordering induced by the presence of the organic solvent, the origins of which are detailed in a subsequent section.
A
B
Figure 5: Simulation snapshots capturing interfacial behavior of small compounds (green) with both the membrane, water, and dodecane (brown). (A) Bisabolene is relative large, and part of the molecule can interact directly with dodecane while the center of mass is around 20 ˚ A. (B) Pinene is relatively small, and cannot directly interact with dodecane. Instead, in this case the dodecane acts via perturbing water and membrane structure to lower the desorption barrier. The observed barriers for compound extraction from the membrane into either aqueous solution (left) or dodecane (right) are unequal due to differential interaction between the compounds and either water or dodecane. Since the terpenoid compounds are not point
19
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
objects, but rather are extended molecules, larger compounds such as bisabolene may interact directly with dodecane even as their center of mass is within the membrane or near the interface (Fig. 5A). The favorable interactions formed between the dodecane and the terpenoid would naturally yield a lower energetic barrier to desorption than the equivalent extraction into water, which lacks these hydrophobic interactions. However, these smaller barriers also persist in smaller molecules that are too small to directly interact with dodecane while within the membrane headgroup region, such as pinene (Fig. 5B). The reduced organic extraction barrier in small molecules can be ascribed to dodecane ordering the water structure. Aqueous solution is largely bulk, disordered water, where the individual molecules are free to rotate and diffuse, and can find an interaction partner in every direction in space. Near the interface, water molecules preferentially orient to satisfy their hydrogen bonding requirements with nearby lipids, and the impact of this can still be easily observed to 30 ˚ A from the bilayer center. 76 From our own equilibrium trajectories, we determine the geometric orientation of water molecules as a function of distance from the bilayer center (Fig. S5). Unlike the aqueous simulations, where the ordering breaks down further from the membrane as in previous studies, 76 simulations with an organic phase present demonstrate ordering near the water-organic phase interface in addition to the ordering near the membrane. The ordering of the waters near the organic phase reduces the number of easily accessible states, and thereby their entropy with respect to bulk water. The reduced entropy directly lowers the energetic barriers to moving a hydrophobic molecules into solution, since the typical entropy loss associated with creating an ordered water “cage” around the hydrophobic molecule has already happened. The other notable trend from Fig. 4 is the relative extraction barriers between hydroxylated and non-hydroxylated forms of the terpenoids. Hydroxylation uniformly lowers the membrane extraction barrier, with the typical barrier to membrane desorption being reduced by approximately 1 kcal mol−1 relative to the non-hydroxylated product. This would reflect favorable hydrogen bonding between the added hydroxyl and the surrounding wa-
20
ACS Paragon Plus Environment
Page 20 of 45
Page 21 of 45 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
ter molecules as the terpenoids enter solution. The trade-off resulting from these hydrogen bonds is that the hydroxylated terpenoids face a 1-4kcal mol−1 larger barrier to crossing the membrane. One final point on the free energies themselves is that the transfer of terpenoid product from the lipid bilayer to the organic phase is always favorable (Fig. 4). This is dissimilar from the fatty acyl results, where transfer from the bilayer to the dodecane phase was typically an energetically unfavorable process, albeit barely. 18 The origin of the discrepancy lies in the orientation imposed by the lipid tails. Fatty acyl products are highly similar to lipid acyl tails in their conformational flexibility, and therefore preferentially segregate to the membrane environment rather than the disordered organic phase. By contrast, the rigid terpenoids are rotationally confined within the membrane, but are free to explore these rotational degrees of freedom in the dodecane phase, lowering the free energy there. Similar mechanisms have been previously proposed when aromatic amino acid analogs are transferred into organic solvents within a bilayer. 63 The local diffusivities reported in Fig. 4 are all similar in magnitude and shape, reflecting the chemical similarity of the individual compounds. Based on the local variation in the diffusivity (Fig. 4), we estimate the uncertainty of individual values to be below 10 %, which will have only a small effect on the computed permeability. The shape of this shared profile reflects the local structure of the medium being traversed. In regions of disorder, such as in solution or the organic phase, the local diffusion is highest, with aqueous solution demonstrating higher diffusion than the organic phase due to the faster diffusion of water relative to dodecane. Similarly, near the membrane center where the tails are the least ordered, there is another increase in diffusivity (Fig. 4) similar to that seen in some small molecules, 77–79 as well as the comparison fatty acid simulations. 18 The mechanism for these regions of high diffusion may be voids at the interface between leaflets formed by tail movement. The reduced pressure environment of the microscopic voids would draw in terpenoids, raising the diffusion relative to the remainder of the membrane, where molecular collisions with the
21
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
dense lipid tails retards motion. 79 Farnesene and Farnesol
Figure 6: Free energy profile and local diffusivity for farnesene and farnesol membrane crossing, as well as their nearest equivalent fatty acyl products previously studied 18 for comparison.
Farnesene and farnesol were treated separately from the other terpenoids considered due to the different RCs required to describe membrane crossing and extraction. For membrane crossing, only the position of the C1 carbon was used to drive the transition from one membrane leaflet to the other (Fig. 6), analogous to the RC used for prior fatty acid studies, 18 making these previous results the natural point of comparison. The free energy profiles and diffusivities for farnesene and farnesol show the same chemistry-dependent trends as the
22
ACS Paragon Plus Environment
Page 22 of 45
Page 23 of 45
equivalent fatty acyl products. Farnesene prefers to be at the membrane center, similar to the other hydrophobic terpenoids, whereas farnesol exhibits a small crossing barrier similar to the other terpene alcohols (Fig. 4). The measured diffusivity for farnesene and farnesol crossings are also in line with the other terpenoids, with a diffusion rate near the bilayer center of around 20 ˚ A2 ns−1 , which is less than the previously reported diffusivities for comparable fatty acid products (Fig. 6). Mechanistically, the extra unsaturation sites in the farnesene rigidify the molecule, preventing farnesene from easily filling in voids within the membrane that do not have the right shape, slowing diffusion relative to the comparatively flexible saturated fatty acids. A
B Pentadecane Cetyl alcohol Farnesene Farnesol
∆G (kcal/mol)
15 10 5 0 −5 D (Contact Number2/fs)
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
Contact Number
Contact Number
8 6 4 2 0 0
100
200
300 400 Contact Number
500
600
−800
−600
−400
−200 0 200 Contact Number
400
600
Figure 7: Free energy profile and local diffusivity for farnesene and farnesol membrane extraction into aqueous solution (A) and into an organic phase (B). For comparison, the nearest equivalent fatty acyl products previously studied, 18 cetyl alcohol and pentadecane, are also included. The definition for the contact number is provided in Eq. 2, keeping in mind that extraction from the membrane proceeds from right to left as the contact number is decreased. Extraction of farnesene and farnesol into both aqueous (Fig. 7A) and organic (Fig. 7B) environments also exhibit existing trends observed for the other terpenoid compounds tested. The barrier to extraction into both environments is lower for farnesol, just as it is for fatty 23
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
alcohols or hydroxylated terpenoids more generally. Following the same trend as the other terpenoids, farnesene and farnesol have a lower free energy in dodecane where rotational degrees of freedom for the compound as a whole are not constrained as they would be within the bilayer. The diffusivities are broadly in line with what was observed for the fatty acyl products, although no direct comparison can be made due to the different units that the diffusion was measured in relative to the crossing step or that of the other terpenoids more broadly. To handle the different units for diffusivity in Eq. 1, we determine the relationship between contact number and displacement in Fig. S6.
Induced Lipid Ordering by Water/Organic Interfaces Before determining the permeabilities, we revisit the impact of the water-dodecane interface on membrane dynamics, which necessitated moving to a constant area simulation barostat to mitigate asymmetries in the internal bilayer free energy profile (Fig. S4). In the original design of the simulations, the significant surface tension generated at the water-dodecane interface 80 and its coupling to membrane structure was not appreciated. Typical surface tensions for membrane components and water are between 1-5 mN m−1 , 81 while those for the dodecane water interface are around 50 mN m−1 . 80 This means that the preferred contact area between water and dodecane will be smaller than the contact area between water and membrane. In a periodic simulation system, the shrinking water-dodecane interface will therefore compress the membrane-water interface. The coupling of the two interfaces in a conventional barostat will therefore shrink both the membrane and doecane-water interfacial surface areas (Fig. S2), resulting in a thicker membrane. We highlight here that the membrane thickening is unnatural. In an experimental setting, the two interfacial tensions would not be coupled as they are in a periodic membrane simulation, and so there is no reason to suspect membranes thicken in the presence of an organic solvent. The thickening of the bilayer changed the behavior of the individual membrane lipids as well. Specifically, the thicker bilayer packed together individual lipid tails more closely, 24
ACS Paragon Plus Environment
Page 24 of 45
Page 25 of 45 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
Figure 8: Deuterium order parameters for each of the prominent lipid tails compared between aqueous (black) and organic (gray) systems during unbiased equilibrium simulation in a typical constant ratio barostat.
25
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 26 of 45
resulting in a more ordered bilayer (Fig. 8). Small molecules have been previously shown to partition between coexisting ordered and disordered membrane phases depending on their specific chemistry. 82–84 Thus, the free energy profile changes between Figs. 4 and S4 are purely an artifact of lipid ordering induced by the dodecane-water interfacial surface tension if the barostat is allowed to control all dimensions of the periodic cell. Since the permeabilities are exponentially dependent on the free energies (Eq. 1), these free energy differences will significantly impact the result. In the previous fatty acyl product study, 18 as well as the aqueous crossings performed in this work, no additional surface tension was present that would thicken the bilayer, since these calculations were performed in the absence of an organic layer. To mimic this, we use a fixed area barostat for the organic simulations involved in generating Fig. 4 to fix the membrane width, and thereby preserve membrane internal structure.
Permeability Having determined all the requirements to compute a permeability based on Eq. 1, we can now tabulate the permeabilities of each step (Table 1), which in turn are proportional to the anticipated flux for a given concentration differential. 18 The overall permeability for the total crossing and extraction process (Pt ) is limited by the slower of the crossing (Pc ) and extraction (Pe ) permeabilities due to the additive nature of Eq. 1:
Pt =
Pc
Pc−1
Pc Pe 1 = ' Pe −1 + Pe Pc + Pe 1P 2
if Pc