Article Outline

Article Information

PubMed

Related Articles

  • …more

Copyright © 2008 The Biophysical Society. All rights reserved.
Biophysical Journal, Volume 94, Issue 8, 2938-2954, 15 April 2008

doi:10.1529/biophysj.107.118380

Biophysical Theory and Modeling

Biophysical Regulation of Lipid Biosynthesis in the Plasma Membrane

Stephen H. Alley*Oscar CesRichard H. Templer and Mauricio Barahona*Go To Corresponding Author 

* Department of Bioengineering, Imperial College London, South Kensington Campus, London, United Kingdom
Department of Chemistry, Imperial College London, South Kensington Campus, London, United Kingdom
Institute for Mathematical Sciences, Imperial College London, South Kensington Campus, London, United Kingdom

Address reprint requests to M. Barahona.

Abstract

We present a cellular model of lipid biosynthesis in the plasma membrane that couples biochemical and biophysical features of the enzymatic network of the cell-wall-less Mycoplasma Acholeplasma laidlawii. In particular, we formulate how the stored elastic energy of the lipid bilayer can modify the activity of curvature-sensitive enzymes through the binding of amphipathic α-helices. As the binding depends on lipid composition, this results in a biophysical feedback mechanism for the regulation of the stored elastic energy. The model shows that the presence of feedback increases the robustness of the steady state of the system, in the sense that biologically inviable nonbilayer states are less likely. We also show that the biophysical and biochemical features of the network have implications as to which enzymes are most efficient at implementing the regulation. The network imposes restrictions on the steady-state balance between bilayer and nonbilayer lipids and on the concentrations of particular lipids. Finally, we consider the influence of the length of the amphipathic α-helix on the efficacy of the feedback and propose experimental measurements and extensions of the modeling framework.

Introduction

The primary function of the lipids in the plasma membrane is to form a bilayer that provides a permeability barrier between the cytoplasm and the environment. However, whereas lipids were once considered purely passive components, it is now clear that lipids play an active role in a variety of dynamic processes involving the membranes that compartmentalize the cell 1. To achieve this dual role of the membrane as a dynamic boundary and a continuous barrier, the cell must regulate the mechanical properties of the membrane and does so partly by controlling its lipid composition.

Membrane lipids are chemically diverse 2 but they can be classified into the broad categories of bilayer and nonbilayer lipids, depending on their (in)ability to self-assemble into bilayers. Bilayer formation is the result of a thermodynamic equilibrium in which the physicochemical properties of the lipids, such as the chemical structure of the headgroup and hydrocarbon chains, play a crucial role. The cell can therefore regulate the mechanical properties of the bilayer by modifying its lipid composition through lipid biosynthesis. The balance between bilayer and nonbilayer lipids in the plasma membrane has been the subject of many reviews 2,3,4. Experiments have shown that organisms change the lipid composition of their membranes in response to external variations in diet, pressure, and temperature (5–7; see also 18). Moreover, many of the lipids found in biological membranes do not form bilayers under physiological conditions. Subsequent studies have confirmed that most organisms contain significant amounts of at least one nonbilayer lipid 8,9.

The underlying biophysical question is the relationship between the chemical diversity and variability of membrane lipid composition, the mechanical properties of the membrane, and the associated protein functions 10,11. A large experimental effort has been devoted to mapping lipid biosynthetic pathways by characterizing and mutating particular enzymes. There is also an increasing body of experiments that measure the relationship between the biophysical properties of lipids and enzyme activity 10,12. However, there have been few attempts 13 to consider theoretically the interdependence of these two phenomena by modeling the lipid biosynthetic network as an integrated system in which the biochemical and the biophysical descriptions of the metabolic network are fundamentally linked. The focus of the model presented here is to provide a set of tools to understand the interplay between the enzymes and lipids involved in lipid metabolism in relation to the biophysical properties of the bilayer.

Fig. 1 depicts a simplified representation of the connection between the chemical structure of lipids and the mechanical properties of a lipid monolayer. A lipid monolayer consists of conformationally flexible lipids, whose amphiphilic nature leads to a nonuniform pressure distribution across the monolayer. The lateral pressure profile π(z) determines the average “molecular shape” that a lipid adopts and, more importantly, Js, the monolayer spontaneous curvature. Js is an intrinsic property of a lipid species that corresponds to the monolayer curvature in which a lipid in the monolayer is at the conformation with minimum free energy 14. The spontaneous curvature reflects the desire of a lipid monolayer to either curve away from, or curve toward, the membrane-water interface and whether a lipid is a bilayer or nonbilayer lipid.

Display large version of this figure
Display high quality version of this figure
Figure 1
Forces that act between lipids at different depths include electrostatic and hydrogen bond interactions at the headgroup, interfacial tension at the hydrophilic-hydrophobic interface, and the packing of the hydrocarbon chains. The lateral pressure profile π(z) depends crucially on the chemical nature of the lipid head group and the length and saturation of the lipid hydrocarbon chains. The lateral pressure profile determines the spontaneous curvature Js of a lipid monolayer 14.

A lipid bilayer is formed by two monolayers back-to-back. This arrangement means that the monolayers may not be able to adopt their preferred curvature, Js, since the monolayers in the bilayer are held together by the hydrophobic effect. This leads to a difference between the actual curvature of a monolayer, as given by the principal curvatures c1 and c2, and its spontaneous curvature, Js15. Based on this physical picture, Helfrich 16 formulated the stored elastic energy per unit area, g, of a lipid monolayer that is constrained to have principal curvatures c1 and c2:

where F is the Helmholtz free energy, A is the area, and κM is the bending rigidity of the monolayer. The lipids are at the free energy minimum, when the total curvature, c1+c2, is at the value of the spontaneous curvature Js. Because at equilibrium g is minimized, this means it is more difficult for lipids with large spontaneous curvatures to form a bilayer, which is a flat conformation with small c1 and c2. Indeed, it has been suggested that lipids with Js<−1/6nm−1 (the negative sign is a convention to denote that the monolayer curves toward water in an interface) do not form bilayers. Instead, they form curved mesophases, such as the inverse hexagonal phase, which are porous 17.

From a biological perspective, porous mesophases would have severe consequences for cellular function and survival. Gruner 3 hypothesized that the average spontaneous curvature of the lipids in the plasma membrane must be tightly regulated to ensure that the membrane lipids form a (nonporous) bilayer and that the cell is able to control 15 by modifying its lipid composition through the biochemical networks of lipid metabolism.

This insight has been confirmed experimentally. Lipid extracts from the cell-wall-less Mycoplasma Acholeplasma laidlawii grown under different conditions have average spontaneous curvatures in the small range between −1/6.6nm−1 and −1/8.1nm−1 even though the membrane contains lipids with Js outside of this range 18. To achieve this robust regulation, A. laidlawii alters the ratio of its two main glucolipids in response to the length and saturation of exogenously fed fatty acids 5, thus maintaining in a biologically viable “bilayer range” that ensures membrane integrity yet with enough stored elastic energy to allow for its dynamical behavior. Remarkably, although the average spontaneous curvature is controlled, the lipid concentrations exhibit wide variations. This suggests that the control of is not achieved by targeting specific lipid compositions. These observations also apply to Escherichia coli lipid extracts, which begin to form nonbilayer structures close to physiological conditions 7.

Biophysical control mechanisms integrated into lipid biosynthetic networks have been the subject of intense experimental study. An example is given by cytidine triphosphate/phosphocholine cytidyltransferase (CCT), an enzyme involved in the biosynthesis of the ubiquitous lipid phosphatidylcholine. CCT is inactive in the cytoplasm, but becomes active when membrane-bound. It has been shown that its activity is affected by the stored elastic energy in the membrane 12. The biophysical control mechanism arises from the presence of an amphipathic α-helix that affects enzyme activity by regulating the binding of CCT to lipid bilayers. In a broad sense, the amphipathic α-helix can be viewed as a “sensor” of the spontaneous curvature since it binds preferentially to lipid bilayers with large negative thus modulating the activity of the lipid biosynthetic enzyme. This biophysical control mechanism is chemically nonspecific, as it is based on a biophysical interaction between the enzyme and the membrane, and appears to be generic to a number of enzymes present in lipid biosynthetic pathways 12,19, including those present in A. laidlawii, which is the focus of this study.

We have developed a modeling framework for the lipid biosynthetic pathways in A. laidlawii. Building upon the A. laidlawii biochemical network studied in detail by the groups of Lindblom, Rilfors, and Wieslander 5,6,19,20, we formulate a biophysical mechanism, based upon some of the conceptual foundations established in CCT 12, that couples the activity of lipid biosynthetic enzymes to the membrane composition. Our results show that the presence of feedback increases the robustness of the steady state of the system to parameter variations, in the sense that it decreases the probability of inviable values of that would lead to porous phases. From a sensitivity analysis, we identify the enzymes that are most efficient in implementing the control of the network. We also study the restrictions that the network imposes on the steady-state concentrations of particular lipids and show that the system keeps a balance between bilayer and nonbilayer lipids. Finally, we consider the influence of the length of the amphipathic α-helix on the efficacy of the feedback.


The lipid biosynthetic network of a. laidlawii

We take the cell-wall-less Mycoplasma A. laidlawii as our system for the study of cellular models of lipid biosynthesis. This simple organism, which has been studied in great detail 5,6, has two features that make it ideal to showcase our modeling framework. First, virtually all the lipids in A. laidlawii are in the plasma membrane 21. This simplifies the model to a single lipid bilayer, avoiding the complexity of cell walls and intracellular compartments. Second, A. laidlawii cannot synthesize unsaturated fatty acids and is very limited in its synthesis of saturated fatty acids 5. Therefore, A. laidlawii exhibits a significantly reduced number of chemical species in the plasma membrane as it relies on exogenously fed fatty acids for lipid biosynthesis.

The limited fatty acid synthesis implies that the only response of A. laidlawii to variations in its fatty acid diet is to alter the composition of the headgroups of the lipids in the membrane through the network of enzymatic reactions represented in Fig. 2. Indeed, experiments show that the membrane lipid composition of A. laidlawii depends strongly on the length and saturation of exogenously fed fatty acids 5. When A. laidlawii is fed palmitic acid (a short, saturated fatty acid), monoglucosyldiacylglycerol (MGlcDAG) is the most abundant lipid; whereas when A. laidlawii is fed oleic acid (a long, unsaturated fatty acid), diglucosyldiacylglycerol (DGlcDAG) dominates. Central to our study is the observation that although the variation in the lipid composition can be large, the cell maintains the average monolayer spontaneous curvature of the plasma membrane within a “window” in which the bilayer phase is stable 7 (Table 1). The lipid biosynthetic network is able to adjust the lipid composition to achieve a >−1/6nm−1, thus maintaining a dynamic, yet impermeable plasma membrane.

Display large version of this figure
Display high quality version of this figure
Figure 2
A. laidlawii lipid biosynthetic network. (A) The biochemical network. The main lipids in the plasma membrane of A. laidlawii A-EF22 are 6: phosphatidylglycerol (PG), diacylglycerol (DAG), monoglucosyl-DAG (MGlcDAG), diglucosyl-DAG (DGlcDAG), and glycerophosphoryl-DGlcDAG (GPDGlcDAG). Phosphatidic acid (PA), the liponuleotide CDP-DAG, and PG-phosphate (PGP) are lipid intermediates. The top branch is the PG pathway and the bottom branch is the glucolipid pathway. The abbreviated soluble reactants are glucose (Glc) and UDP-Glc, glycerol-3-phosphate (G3P), the inorganic phosphate ions Pi and PPi, and the nucleotide CTP. R indicates an acyl chain. Six of the seven enzymes have irreversible rate equations. A. laidlawii also synthesizes three monoacyl derivatives of the glucolipids: monoacyl-MGlcDAG (MAMGlcDAG), monoacyl-DGlcDAG (MADGlcDAG), and monoacyl-bisglycerophosphoryl-DGlcDAG (MABGPDGlcDAG) 6. However, these lipids have been excluded from the model as they are not always synthesized 5 and their biosynthetic pathways have been postulated, but are not known 59. (B) A biophysical picture of the network. Lipids are color coded according to their Js, which is linked to their molecular shape as shown. The same color code is used to show that the activity of MGS and DGS increases when the plasma membrane has a large negative It can be seen, for instance, that in the case of the lower pathway the effect of MGS and DGS is to increase the effective size of the headgroup of the lipid upon which they are acting and therefore systematically increase the value of Js among DAG, MGlcDAG, and DGlcDAG. By controlling the rate of the steps between DAG/MGlcDAG and MGlcDAG/DGlcDAG, the system is capable of regulating A. laidlawii also synthesizes three monoacyl-derivatives of the glucolipids: monoacyl-MGlcDAG (MAMGlcDAG), monoacyl-DGlcDAG (MADGlcDAG), and monoacyl-bisglycerophosphoryl-DGlcDAG (MABGPDGlcDAG) 6. However, these lipids have been excluded from the model as they are not always synthesized 5 and their biosynthetic pathways have been postulated, but are not known 59.

The lipid biosynthetic network: biochemical and biophysical descriptions

The biochemical description of the lipid biosynthetic network of A. laidlawii is presented in Figure 2A. The first step in the metabolic network is, as in other organisms, the acylation of soluble glycerol-3-phosphate (G3P) to form phosphatidic acid (PA) 22. The network then branches out into two pathways.

The upper branch is the phosphatidylglycerol (PG) pathway, well-studied in bacteria, in which PA is converted into PG through the intermediates cytidine diphosphate diacylglycerol (CDP-DAG) and phosphatidylglycerolphosphate (PGP). The corresponding enzymes CDP-DAG synthase (CDS), PGP synthase (PGPS), and PGP phosphatase (PGPP) have been characterized in E. coli23,24 and in Clostridium perfringens25,26.

The lower branch is a specific pathway in A. laidlawii, deduced from the discovery and purification of the PA phosphatase (PAP) 27 and the two consecutive glucosyltransferases, MGlcDAG synthase (MGS) 28 and DGlcDAG synthase (DGS) 29. The final enzymatic reaction is yet to be characterized since the glycerophosphoryl-DGlcDAG (GPDGlcDAG) synthase (GPDGS) that catalyzes the production of GPDGlcDAG has not been purified yet. However, the genetic similarity of MGS and DGS to the enzymes of Gram-positive bacteria 30 suggests that GPDGlcDAG could be synthesized by the transfer of G3P from PG to DGlcDAG, a reaction that occurs in the synthesis of lipoteichoic acids in the cell walls of Gram-positive bacteria 31.

Figure 2B presents a biophysical interpretation of the network, showing how the molecular shape of each lipid is reflected in its monolayer spontaneous curvature. This physical picture shows that the position of nonbilayer lipids (Js<−1/6nm−1) and bilayer lipids (Js>−1/6nm−1) within the network has an effect on which enzymes can exercise effective control of the of the plasma membrane. By inspection, MGS and DGS are good candidates for the control of since MGS and DGS catalyze the reactions that lead from the lipid with the most negative Js (DAG) to the lipid with the least negative Js (DGlcDAG) (Figure 2B). This intuition is reinforced by a structural feature of these enzymes. MGS and DGS are both peripheral membrane proteins that translocate between the cytoplasm and the membrane. It is postulated that they are only active when they are inserted into the membrane, as suggested by the increased activity of both MGS 29 and DGS 32 in the presence of lipids with large negative Js. This picture leads to a biophysical, intrinsic mechanism for MGS and DGS to control lipid biosynthesis as a function of the average monolayer spontaneous curvature of the membrane. Our model is a mathematical formulation of these ideas.


Cellular model of lipid biosynthesis

The biochemical constituents of our cellular model of lipid biosynthesis are the membrane lipids, the lipid biosynthetic enzymes, and the soluble cytoplasmic reactants. The membrane lipids are assumed to be homogenously distributed over both monolayers of the plasma membrane. Although labeling studies in E. coli show that lipid biosynthesis occurs mainly at the inner leaflet of the plasma membrane 33, we will assume that lipid transport from the inner to the outer leaflet of the membrane maintains a symmetric bilayer. Our assumption of spatial homogeneity for the lipids is based on the fast lateral diffusion of lipids in bilayers 34,35 and leads to a description in terms of ordinary differential equations. The soluble reactants (such as the nucleotide CTP or the inorganic phosphate ions Pi and PPi) are assumed to have constant, regulated cytoplasmic concentrations, due to their involvement in general cellular processes. Therefore, they are only parameters (not variables) of the model.

The A. laidlawii lipid biosynthetic network is modeled as a system of nonlinear differential equations for eight lipids with seven enzymatic reactions. The variables of the model are compiled into the vector of lipid surface concentrations expressed in molar fraction: LT=[{PA} {CDP-PAG} {PGP} {PG} {DAG} {MGlcDAG} {DGlcDAG} {GPDGlcDAG}]. The sum of the lipid molar fractions is 1 at all times: 1TL=1.

Each enzymatic reaction has a nonlinear rate equation of the Michaelis-Menten type, modified using surface dilution kinetics, as explained below, to account for the fact that the reactions take place on the membrane. The enzyme rate equations are compiled into a vector vT=[vCDSvPGPSvPGPPvPAPvMGSvDGSvGPDGS]. The modulation of the enzyme activity due to the biophysical interaction with the mechanical properties of the membrane is introduced through a diagonal matrix Ka=diag([Ka,CDSKa,PGPSKa,PGPPKa,PAPKa,MGSKa,DGSKa,GPDGS]), which incorporates the possibility that some of the enzymatic rates, specifically those of MGS and DGS, could depend on If the enzyme is curvature sensitive, its association constant Ka,Enzyme will depend on Otherwise, the corresponding Ka,Enzyme=1. This is the basis of the biophysical feedback mechanism, which will be introduced in the following section.

The topology of the reaction network is encoded in a stoichiometric matrix N, where Nij is the number of lipid species i consumed (negative) or produced (positive) in reaction j:

The matrix N accounts for the enzymatic reactions and ensures mass conservation. However, our cellular model must also include both the lipid degradation into soluble products and the lipid insertion that enables a growing cell to double the number of lipids before cell division. These processes are incorporated through the transport vector t and the normalization vector n. The transport vector t encapsulates the balance of lipids inserted and extracted. In our model, only PA is inserted at a constant cellular rate V+,PA and lipid degradation is assumed not to play a significant role in A. laidlawii lipid metabolism 36. Therefore, tT=[ V+,PA 0 0 0 0 0 0 0].

The normalization vector n, given by

reduces the surface concentration of each lipid in proportion to its molar fraction while at the same time maintaining the sum of the molar fractions equal to 1.

Combining all the terms, the model can be written compactly as

(1)

This system has stationary points L*.

Finally, to close the system, we need to relate to the lipid concentration. Our underlying, linear assumption is that is well approximated by the weighted average of the spontaneous curvatures of the individual lipids:

(2)

This linear assumption has been shown experimentally to lead to an accurate approximation of the phase behavior of lipid mixtures 37. This linear assumption is also used in many of the experiments that measure the Js of neutral and anionic lipids 38,39,40,41. Using Eq. (2) with the Js values and experimental lipid composition Lexp in Table 1 leads to a calculated of −1/7.9nm−1, which lies within the measured range from −1/6.6nm−1 to −1/8.1nm−1 of A. laidlawii lipid extracts.


Functional form of the lipid biosynthetic enzyme kinetic rates, v(L)

Before considering the biophysical mechanism that couples the biochemical reactions to the mechanical properties of the membrane, we state first some specific features of the enzyme kinetic rate equations of the membrane lipid network. The functional form of the rate equations v(L) in the model differs from standard enzyme kinetics 42 in two respects. First, our cellular model must take into account the number of copies of the enzyme in the cell. Second, we must account for the fact that lipid biosynthetic enzymes have soluble, cytoplasmic reactants that diffuse in three dimensions, whereas their lipid reactants diffuse within the two-dimensional membrane.

Kinetic studies 27,28,29 have fitted the rates of A. laidlawii lipid biosynthetic enzymes to surface-dilution kinetics, in which soluble reactants have a bulk concentration in units of molarity and membrane reactants have a surface concentration in (dimensionless) molar fraction 43. All of the enzymatic reactions in the network, except for the reaction catalyzed by CDS, can be assumed to be irreversible. There is experimental evidence that supports this assumption, e.g., the hydrolyses of the phosphoanhydride bonds in PA and PGP are irreversible 44. Therefore, the rate equations for vPGPS, vPGPP, vPAP, vMGS, and vDGS are of the form

(3)

where {L}i is the surface concentration of the lipid substrate (in molar fraction) and [S] is the bulk concentration of the soluble substrate (in units of molarity). Similarly, KmL is the Michaelis constant of {L}i (in molar fraction) and KmS is the Michaelis constant of [S] (in units of molarity). Experimental values of the enzyme kinetic constants are listed in Appendix B .

Note that Vcell is the rate for all copies of the enzyme in the cell (in units of molar fraction/min):

(4)

where MEnzyme is the total mass of each enzyme in the cell and NLipid is the number of moles of lipid in the cell. Vmax is the standard Michaelis-Menten limiting rate, which typically has units of moles of product synthesized per milligrams of enzyme per minute 42. It is assumed that the ratio MEnzyme/NLipid is kept constant in a growing cell over the cell cycle. In Appendix B we show how we have estimated these parameters.

Two of the enzymatic reactions have slightly different functional forms. The final reaction of the lower path, catalyzed by GPDGS, although irreversible, involves two lipid substrates. As mentioned above, the CDS reaction is modeled reversibly since the equilibrium constant is much less than 1 23. The rate equations of these reactions are listed in Appendix B .


Spontaneous-curvature-sensitive enzymes

We now introduce the terms in the model that describe how the activity of an enzyme with an amphipathic α-helix is modulated as a function of spontaneous curvature, which is in turn a function of the lipid composition. As mentioned above, there is extensive evidence that supports the theory that the activity and function of many proteins, both integral and peripheral, are regulated by the biophysical properties of biological membranes 10,35. Such phenomena differ markedly from specific protein-lipid interactions. Although our model deals with the binding of an amphipathic α-helix to the membrane, the mechanism could be extended to describe the binding of other amphipathic motifs to the membrane.

Enzyme kinetic studies have shown that lipids with large negative Js increase the activity of both MGS 29 and DGS 32. MGS has an amphipathic α-helix between residues 67 and 85 30,45, that shares 5 of its first 8 residues with an α-helix of the E. coli division-site-selection protein MinD 46 that targets heterologous proteins to the membrane 47. Since it has been shown that MGS 19,29, the MGS amphipathic α-helix 20, and the MinD amphipathic α-helix 48 all preferentially bind to membranes with large negative we hypothesize that the curvature-sensitive activity of MGS is a result of the membrane binding of this α-helix.

Through surface plasmon resonance (SPR) experiments, it has been concluded that liposomes bind to MGS through a two-step process 19. The first binding step is independent of lipid composition and has a dissociation constant of ∼10 nM. The second binding step has a large dependence on lipid composition, as its dissociation constant decreases from 10mM to 100nM when the liposomes have large negative 19. Since liposomes with large increase both the activity of MGS and the strength of the second binding step, it follows that MGS is only active after the second binding step.

Amphipathic peptides form random coils in solution. The first binding step corresponds to surface adhesion induced by the electrostatic attraction of exposed basic residues to acidic membrane lipids. The second, subsequent step is the insertion of the hydrophobic residues into the membrane coupled with the emergence of the α-helix, which is entropically favored by the hydrophobic membrane environment. This two-step membrane binding 49 can be summarized through a simple kinetic mechanism,

where Kd1 and Kd2 are the dissociation constants of the binding steps. At steady state, the fraction of membrane-inserted amphipathic α-helices that result in active enzymes is given by the association constant

(5)

We now derive expressions for Kd1 and Kd2 from biophysical considerations.

First binding step, Kd1

From SPR studies, Kd1 is measured to be ∼10 nM 19. It is proposed that this first (irreversible) binding is a result of electrostatic attraction. Structurally, the presence of eight positively charged residues on the 19-residue amphipathic α-helix (see Figure 4B) will produce a strong electrostatic attraction. Indeed, there is ample evidence that negatively charged anionic lipids are essential for the binding and activity of MGS. For instance, it is known that shielding the anionic lipids with 0.75M NaCl prevents the binding of MGS 19.

From simple electrostatic considerations, Kd1 is given by the Boltzmann relation,

(6)

where zpe= +8e is the net charge of the amphipathic α-helix and ψ0 is the membrane surface potential. A dissociation constant of 10 nM would imply ψ0≅−60mV at 40°C, which is comparable to the measured membrane surface potentials of bacterial lipid bilayers 50. This simple estimate reinforces the plausibility of the interpretation of the first binding step in terms of electrostatic interactions. Clearly, the biophysical picture will be complex, including the shielding of charges on the peptide to give an effective valence 51 and the likely involvement of other positively charged enzyme domains.


Second binding step, Kd2

The second binding step involves at least three energetic processes: membrane insertion of the hydrophobic residues; peptide folding to form the α-helix; and lipids bending to accommodate the inserted α-helix. It has been observed that Kd2 decreases dramatically along the lipid sequence dioleoylphosphatidylglycerol (DOPG)>cardiolipin (CL)>dioleoylphosphatidylethanolamine (DOPE)>dioleoylglycerol (DOG) 19, i.e., as Js becomes more negative 39,52. This is the basis for our assumption that the second binding step is dominated by the energy of lipids bending to accommodate the helix.

We can understand this process through the following simplified biophysical picture. Consider a locally flat bilayer with average monolayer spontaneous curvature The diameter of the α-helix is comparable to that of a lipid. Consequently, the insertion of an amphipathic α-helix into a flat membrane does not result in a change of the monolayer curvature, yet it leads to a change in the molecular shape of the lipids alongside the α-helix 12 (Figure 3A). This would translate into a monolayer curvature, cbound≠0, for a monolayer formed entirely by lipids like those surrounding the amphipathic α-helix.

Display large version of this figure
Display high quality version of this figure
Figure 3
(A) Geometric argument used to calculate cbound, the curvature of a lipid monolayer consisting entirely of lipids that lie alongside an amphipathic α-helix. (B) Association constant Ka as a function of for a 19-residue α-helix (dark solid line) and, for a 58-residue α-helix, such as that of CCT, plotted for comparison (dark dashed line). The dashed vertical line is Js=−1/6nm−1 and the solid vertical lines give the measured range of of lipid extracts 5. (Inset) Helical-wheel projection of residues 67–85 of MGS. The bar gives the Eisenberg consensus normalized hydrophobicity scale 76.

Figure 3A sketches a very simple geometrical argument to obtain a first-order estimate of cbound:

(7)

where r=0.45nm is the radius of the α-helix 53; t is the monolayer thickness (the distance between the middle of the α-helix and the bilayer midpoint), which is measured to be 1.71nm 53; and A is the average interfacial surface area of the A. laidlawii lipids, which is measured to be 0.65nm25. A is assumed to be square and the pivotal plane is assumed to coincide with the middle of the α-helix. Equation (7) gives an estimated cbound ≈ −1/3.2nm−1, which is significantly nonflat.

The cylindrical deformation of the α-helix ensures that one of the principal curvatures is zero, c2=0. Therefore, the change in the stored elastic energy in Eq. (1) due to the bending of the lipids alongside the amphipathic α-helix is

(8)

where κM is the bending rigidity of lipids, which we take to be 10 kBT38, and Nh=7.1 is the number of lipids that adopt curvature cbound along both sides of the amphipathic α-helix. Nh is calculated for a 19-residue α-helix of length 2.85nm with 3.5 lipids of length 0.651/2 nm along each side of its long axis. For the range of Js values in Table 1, the free energy is between −1.7kBT/lipid and 0.2kBT/lipid. These energies are not large enough to cause the lateral sequestration of lipids around the α-helix 54, thus justifying the use of

Equation (9) provides an estimate for the binding energy if we assume that the main energetic contribution to this process comes from lipid bending. The dissociation constant of the second binding step would then be given by the Boltzmann relation:

(9)

Note that for the range of Js values in Table 1, the modeled KD2 ranges between 4M and 5 μM. The difference with the experimental values of KD2 may be explained by the enhanced electrostatic attraction due to the absence of divalent cations and the use of zwitterionic lipids in the SPR experiments 19,20. Note that when the binding energy is zero and MGS is equally likely to be bound or unbound. Reassuringly, this bound-to-unbound transition is centered at a value of that lies between the formation of nonbilayer structures and the lower bound of the experimental curvature of lipid extracts in A. laidlawii: −1/6<cbound/2<−1/6.6 (Table 1).

Equations (6) and (9) provide the biophysical feedback for the system in Eq. (1), as the association constant Ka,MGS multiplies the rate vMGS. Given the individual lipid spontaneous curvatures Js in the system (Table 1), and assuming a constant Kd1=10nM, the association constant is constrained to be in the range 1≥Ka,MGS≥0.23 (Figure 3B), and the 19-residue α-helix provides a fourfold regulation of MGS activity. At large negative almost all MGS is active and the synthesis of MGlcDAG increases The opposite effect is produced when is less negative. Clearly, a longer amphipathic α-helix would produce significantly stronger regulation of activity (Figure 3B).



Parameter estimation for the model

The time evolution of the system in Eq. (1) and its corresponding stationary point depend on the model parameters. Most of these parameters have been collected from an extensive survey of the literature, or have been estimated or measured directly. As is usual in the literature, some of the parameters carry substantial uncertainty. In addition, there is an absence of kinetic parameters for some of the enzymes of the A. laidlawii lipid biosynthetic network.

To complete our parametric description, we carry out a constrained nonlinear parameter estimation in which we search for the positive parameter set p that reproduces the experimentally observed lipid concentrations Lexp as close as possible while minimizing the distance to the reliable literature values. Our method of choice to solve this constrained optimization is the Stochastic Ranking Evolutionary Strategy (SRES) 55, an evolutionary strategy with stochastic ranking, which has been shown in a recent survey 56 to be successful in finding feasible parameters in nonlinear biochemical pathways.

For such an underdetermined system, a multiobjective optimization is pursued. The primary objective is to minimize the difference between L*(p), the stationary point of the model in Eq. (1) for the parameter set p, and the experimentally observed lipid composition Lexp given in Table 1,

The secondary objective is to minimize the difference between the estimated parameter set p and the literature parameters plit. Instead of minimizing the Euclidean norm ||pplit|| in our case we minimize a more appropriate measure of the relative distance between the parameter sets, previously introduced to quantify the robustness of dynamical systems 57,

i.e., the sum of the absolute logarithmic errors between p and plit weighted by cj, the confidence in the jth literature parameter. In particular, parameters obtained from A. laidlawii experiments have been assigned cj=1, whereas cj=1/2 for parameters taken from experiments on other organisms, such as E. coli. As a check that the combination of our model and this multiobjective estimation procedure produces plausible parameter sets, we have verified that we can obtain parameters for which the model reproduces all the different lipid compositions that have been observed experimentally 5.

There are a total of 25 kinetic parameters in the model, of which only 18 have literature values. We run our multiobjective optimization by varying 15 parameters (the 7 unknown and 8 parameters with uncertain literature values) and obtain the estimated parameter set p0 presented in Table 2. The distance between each estimated parameter and its literature value has been constrained to be at most two orders of magnitude. The estimated parameter set reflects the sensitivity of the steady state to particular parameters, specifically those appearing in the numerator of the rate equations: Vcell,Enzyme and the soluble substrate concentrations. This emphasizes the importance of measuring intracellular metabolite concentrations.

Our model in Eq. (1) with the estimated parameter set p0 in Table 2 will be our reference system henceforth. This system has a stationary point at the experimental lipid composition shown in Table 1, with . In the next section, we explore the effect of the biophysical regulation mechanism introduced above.



Results

We now investigate the behavior of the cellular model of lipid biosynthesis through the numerical integration of Eq. (1) under a variety of conditions. The first robust feature of the model is that it evolves to a steady state that is independent of the initial condition. Although we have not proved global stability explicitly, numerical integrations from more than 105 randomly generated initial conditions all converge to the same stationary lipid composition. This is strong evidence that the steady state is globally attracting. Our numerical investigation of the model therefore translates into an evaluation of how the fixed point L*(p) changes in response to variations in the parameters or in the presence of feedback.

Which enzyme rates most affect

Experiments indicating that the activities of both MGS and DGS are curvature-sensitive 28,29 have led to the hypothesis that these two enzymes are responsible for the control of If this is true, the rates of MGS and DGS must have a large effect on the steady-state In this section, we use our model to determine which enzymes of the A. laidlawii lipid biosynthetic network have the largest effect on in the absence of feedback. This is an initial step before we introduce the biophysical feedback explicitly in the next section, i.e., in this section the matrix Ka does not depend on the curvature

This question can be posed in the well-known framework of sensitivity analysis, which has been used to characterize, e.g., the robustness of bacterial chemotaxis networks 57. For the reference parameters p0 (Table 2), the model evolves to the experimental lipid composition Lexp (Table 1), i.e., L*(p0)=Lexp. However, the parameters p are inherently noisy. Specifically, the seven Vcell,Enzyme have the most influence on the steady state while at the same time having large variability. Therefore, our sensitivity analysis investigates how changes in the different Vcell,Enzyme translate into variations of the steady-state

The sensitivity analysis is performed through Monte Carlo sampling, one enzyme at a time. The Vcell,Enzyme of the enzyme under study is fixed at the reference value in Table 2. We then produce 106 parameter sets where the other six Vcell,Enzyme are drawn from a random distribution conditioned to produce uniform sampling (over the interval [0,6]) of the logarithmic variation of the parameter set, k:

(10)

Clearly, this implies that the individual Vcell,Enzyme parameters are not sampled uniformly (see the inset of Figure 4A). Effectively, our sensitivity analysis considers variations of up to almost two orders of magnitude in each of the six Vcell,Enzyme, and an overall uniform variation of six orders of magnitude for the complete parameter set. The fixed point for each of the 106 parameter sets is obtained and the corresponding is calculated.

Display large version of this figure
Display high quality version of this figure
Figure 4
(A) Histogram of the distribution PMGS() of the steady-state for a sampling of 106 random parameter sets, where Vcell,MGS is fixed and the other six Vcell,Enzyme values are varied. k is defined in Eq. (10) and measures the logarithmic variation of the parameter set. The individual Vcell,Enzyme distribution used (solid, top inset) ensures that the logarithmic variation k is sampled uniformly (solid, bottom inset). If uniform individual Vcell,Enzyme distributions were used (dashed, top inset), then k would approach a Gaussian distribution (dashed, bottom inset). (B) Marginal probability distributions PEnzyme() of all seven enzymes, where one Vcell,Enzyme is fixed, whereas the other six Vcell,Enzyme are varied. The dashed vertical line indicates the critical value below which bilayers do not form. (C) Cumulative probability that the steady state will be viable () against the logarithmic variation k of the state.

Figure 4A shows the results for the enzyme MGS as a two-dimensional histogram of , the distribution of the 106 random parameter sets. As expected, the distribution is centered on the reference value of and becomes broader as the variation of the parameter set, k, grows. Since biophysical experiments show that membranes with do not form bilayers, we can consider such compositions as biologically nonviable 17. This is marked as a dashed line in Figure 4A. Therefore, the distribution , or its marginal , quantifies how likely it is for the system to evolve to a nonbilayer state when a particular enzyme is kept fixed at its reference value and all other enzymes have uncertain Vcell,Enzyme values. Essentially, this is also a measure of the relevance of the particular enzyme for the controllability of the system, as it quantifies the variability of the output when a given parameter is kept tightly controlled.

The complete results for the system are summarized in Figure 4B, where we plot the marginal distributions for the seven enzymes. All the distributions are centered on the reference value of −1/7.9nm−1; however, the variance and the tails of the distributions differ. Specifically, DGS and MGS (marked with symbols in Figure 4B) have much smaller variances and sharper decay tails than the other five enzymes. The left tail of the distribution is relevant as it gives the proportion of biologically inviable stationary states that do not form bilayers. Figure 4C shows the fraction of viable steady-state lipid compositions as a function of the variability k. The results clearly show that keeping the rate of DGS fixed to its reference value is the most efficient way of guaranteeing a viable system.

The reason for this sensitivity is clear if we examine the network in Figure 2B. The enzyme DGS catalyzes the synthesis of DGlcDAG, with a small negative Js, from MGlcDAG, with a large negative Js. Since the membrane binding of amphipathic α-helices is increased by lipids with large negative Js, this provides a direct link with the biophysical feedback mechanism described above. The implication is that enzymes that catalyze reactions that result in a large positive change in (i.e., DGS and to a lesser extent MGS) are strong candidates to exert the biophysical feedback control of membrane curvature, in accordance with kinetic data. In the next section, the effect of the biophysical feedback provided by the amphipathic motifs is investigated in detail.


Effect of the biophysical feedback on

We now study the effect of the biophysical feedback, mediated by amphipathic α-helices, on the control of the steady-state Our model encodes this mechanism through the association constants Ka,MGS and Ka,DGS collected in the matrix

A sensitivity analysis similar to that performed in the preceding section is carried out to measure how much the variability of the system is reduced in the presence of feedback. We draw 106 parameter sets from a random distribution of Vcell,Enzyme such that the logarithmic variation k, defined in Eq. (10), is uniform over the interval [0,7]. Figure 5A shows four marginal distributions of the steady-state of the system without feedback and with different combinations of feedback on MGS and DGS. In particular, we model MGS to have a 19-residue amphipathic α-helix (see Figure 3B) and we hypothesize that a similar α-helix is responsible for the curvature-sensitive activity of DGS, as suggested by secondary-structure predictions 30. The numerics show that the effect of feedback is noticeable in the reduction of the left tail of the distribution. This means that inviable, nonbilayer steady states are less probable when feedback is present. This is especially prominent for DGS, although MGS also contributes to the control of as shown in Figure 5B. The combined feedback of MGS and DGS reduces the fraction of inviable oleoyl acyl lipid compositions by 19%.

Display large version of this figure
Display high quality version of this figure
Figure 5
(A) Effect of the biophysical feedback on the control of the steady-state The four marginal distributions P(), obtained from 106 random parameter sets in which all seven enzyme rates Vcell,Enzyme are varied, show a reduction of the probability of inviable states with <−1/6nm−1. In the absence of feedback, we fix Ka,MGS at the reference value of Ka,MGS () shown in Figure 3B (×). MGS is modeled to have a 19-residue α-helix (□), DGS is modeled with a 19-residue α-helix (○), both MGS and DGS are modeled with 19-residue α-helices (△). (B) The cumulative probability that the steady state forms a bilayer () as a function of the logarithmic variation k.

The robustness of the lipid compositions

A central feature of the biophysical feedback mechanism is the fact that the enzymatic network controls the physical property and not the steady-state lipid concentrations L*. However, is a function of the lipid concentrations, and it is important to study the underlying variability of L*, with respect to the reference lipid composition Lexp, when the parameters of the model are uncertain. This point can be illustrated with the data obtained in the preceding section through our sensitivity analysis. Fig. 6 shows the probability distribution of steady-state compositions L* as a function of and ||L*−Lexp||1, the distance to the experimental lipid composition, in the absence and in the presence of feedback.

Display large version of this figure
Display high quality version of this figure
Figure 6
Histogram of the distribution of steady-states L*, showing the probability distribution of and the distance of L* to the experimental concentration Lexp. The data are obtained by sampling 106 random parameter sets as in Figure 5A Distribution in the absence of feedback. The vertical dashed line at separates the nonbilayer and bilayer fractions. The dashed lines correspond to the steady-state that results from adding a particular lipid to Lexp until one reaches a monocomponent state, marked by crosses. (B) Distribution when both MGS and DGS exert biophysical feedback through a 19-residue amphipathic α-helix.

In the absence of feedback (Figure 6A), the data show that a majority of L* are close to Lexp, but a range of lipid mixtures is allowed by the system. Note that the system does not by default evolve toward pure, monocomponent compositions. In addition, the L* in the biologically viable region () consist mostly of mixtures of PG, DGlcDAG, and GPDGlcDAG, and not of the intermediates PA, CDP-DAG, or PGP. On the other hand, the inviable L*, with appear through the increase in MGlcDAG, and not of DAG. These trends stem from the constraints that the network structure imposes on the control mechanism and highlight how the steady-state is bounded by the simplex of the Js values of the individual lipids.

Figure 6B shows that the most notable effect of the feedback is to reduce the likelihood of inviable states with values of specifically those with high concentrations of MGlcDAG. In addition, our numerical results clearly indicate that the control is not achieved by targeting specific (possibly monocomponent) compositions, as the overall shape of the probability distribution remains broadly unchanged and mixed states are the norm. In fact, the proportion of monocomponent states is reduced when the feedback is on. Another effect of the feedback is the increase of the probability of states with close to the boundary between bilayer and nonbilayer states. This follows from the functional form of the feedback that is centered on the value cbound/2=−1/6.3nm−1. This is clearly observable in Fig. 6. The implication is that lipid compositions with large negative are more improbable, but at the same time there is an increase in lipid compositions with close to the bilayer to nonbilayer transition.


The effect of the length of the amphipathic α-helices

As seen in Figure 3B, the proposed biophysical feedback exerted by a 19-residue amphipathic α-helix gives rise to a fourfold regulation, which results in a modest reduction of the likelihood of inviable states. The functional form of our feedback implies that the magnitude of the feedback depends strongly on the length of the amphipathic α-helix. We now investigate how the length and number of amphipathic α-helices on MGS and DGS affect regulation.

Biological evidence indicates that the size of amphipathic motifs should be studied parametrically as it is an important feature in protein-membrane interactions. Secondary-structure predictions suggest that both MGS and DGS have multiple α-helices that may insert into the membrane 20,30. However, it is difficult to identify amphipathic α-helices and determine their length from genetic sequences. It is also important to mention that some enzymes with amphipathic α-helices can act in concert. For instance, there is evidence that CCT acts as a dimer that is only active when the 58-residue amphipathic α-helices of both monomers are bound to the membrane 58. In a simplified picture, this dimer would be viewed as having an effective amphipathic α-helix of 116 residues.

Fig. 7 plots the dependence on the length of the amphipathic α-helices of the marginal distribution P() in Figure 5A, where both MGS and DGS are curvature sensitive. Clearly, as the length of the α-helices is increased, the distribution becomes sharper and the likelihood of getting inviable states is reduced. As the inset of Fig. 7 shows, the reduction in the proportion of nonbilayer states increases to more than 50% when both MGS and DGS have 60-residue α-helices.

Display large version of this figure
Display high quality version of this figure
Figure 7
Effect of amphipathic α-helix length on steady-state Shown are three marginal distributions P() of the steady-state under uncertain parameters where both MGS and DGS are curvature sensitive with amphipathic α-helices of different lengths: 19-residue (blue, same as Figure 5A), 57-residue (green), and 114-residue (red). Vertical dashed line at indicates the nonbilayer region. (Inset) Probability that the steady state does not form a bilayer against α-helix length, when MGS has an α-helix (×); DGS has an α-helix (+); and MGS and DGS both have α-helices (□). Colored squares correspond to the three plotted P() in the main figure.


Discussion

This article outlines a bottom-up modeling framework that couples a biochemical network with an intrinsic biophysical feedback mechanism. The central idea behind the biophysical regulation is that the activity of certain enzymes involved in lipid biosynthesis is dependent on the spontaneous curvature of the lipids, which is itself a function of the lipid composition. This introduces a feedback loop that can regulate a property that must be kept within a narrow window to allow for cellular activity and survival 3,7.

The numerical results of our model show that the system evolves toward biologically plausible mixtures of lipids. When the biophysical regulation is present, it decreases the likelihood of inviable steady-state lipid compositions that would not be expected to form a lamellar bilayer. Moreover, our sensitivity analysis indicates that two enzymes in the network (DGS and MGS) have the largest effect on the steady-state in agreement with kinetic data. This fact follows from both the intrinsic properties of the enzymatic reactions (the corres