Transmembrane exchange of materials is a fundamental process in biology. Molecular dynamics provides a powerful method to investigate in great detail various aspects of the phenomenon, particularly the permeation of small uncharged molecules, which continues to pose a challenge to experimental studies. We will discuss some of the recent simulation studies investigating the role of lipid-mediated and protein-mediated mechanisms in permeation of water and gas molecules across the membrane.
The primary function of cellular membranes is to separate the interior of a living cell from its immediate surroundings. At the same time, the cell relies heavily on continuous exchange of materials and information with its environment to sustain its activity and growth and to communicate and coordinate with other cells within a cellular community, e.g., within a tissue or in a whole organism. Transmembrane trafficking of materials, therefore, is one of the most fundamental and highly regulated processes in the biology of all living organisms. As a matter of fact, a cell dedicates a large fraction of its energy and its genetic material to devising and implementing diverse mechanisms facilitating the translocation of materials across the membrane, most prominently through specialized proteins that mediate selective translocation of various molecular species between the two sides of the membrane, namely membrane channels and transporters.
The cellular membrane is primarily composed of lipids, which, owing to their hydrophobic nature, create an energetic barrier against the diffusion of polar and charged species. Excluding exocytosis and endocytosis, which are entirely nonselective mechanisms accompanied by large-scale reorganization of the membrane structure, transmembrane exchange of molecules mainly occurs through one of the three mechanisms: 1) passive diffusion through the lipid bilayer; 2) passive diffusion through membrane channels; and 3) active (energy dependent) translocation mediated by membrane transporters. These three mechanisms along with some examples of the substrates and proteins involved in these processes are graphically presented in FIGURE 1. Translocation of most molecular species, including ions, amino acids and other nutrients, neurotransmitters, as well as larger molecules such as drugs, peptides, and even proteins, across the membrane relies on the function of membrane channels and transporters. In contrast, the transport of very small molecules, e.g., O2, NO, and water, was long deemed to be mediated primarily through nonselective “pores” within the lipid molecules and thus independent of membrane proteins. The discovery of aquaporins (AQPs) and the demonstration of their role in homeostasis of water in various cells and organs suggested that lipid-mediated transport of very small molecules, while providing an important means, might not be sufficient for the desired level of transport that is required, e.g., for physiologically optimal activity of certain cells (10, 77).
A detailed characterization of structural and dynamic properties of lipid bilayers and membrane proteins is necessary for studying transport processes in biological membranes at high resolution. Recent advances in experimental structural biology complemented by biophysical measurements have provided invaluable information on structural and dynamic elements involved in the function of membrane proteins and in characterizing membrane processes. However, a significant number of functional aspects remain unresolved. The challenge is to achieve simultaneously the high temporal and spatial resolutions required for a detailed description of the molecular events and processes involved in membrane transport. Processes such as the mechanism of permeation and selectivity of membrane channels, co-transport of ions along with a neurotransmitter by a secondary transporter, and passive diffusion of hydrophobic molecules through the space (gaps) between the lipid molecules in a membrane can only be “visualized” by methods that offer a dynamic view of the system at hand at an atomic resolution. Molecular dynamics (MD) simulations offer an unparalleled combination of these resolutions and, thus, provide information that cannot be obtained currently by experimental methodologies.
In this article, we will review the application of MD to the study of transport of very small molecules across the membrane. Both lipid-mediated and protein-mediated modes of transport of these molecular species have been studied with MD and will be discussed here. Permeation of water or O2 through membrane channels represent fast phenomena (on the order of nanoseconds) and, thus, are ideal applications for MD simulations (21, 102, 113, 114, 116, 117). Given the simplicity of their primary function, i.e., water transport, and the wealth of structural information (36), AQPs have been extensively studied with MD simulations (21, 102, 113, 114, 116, 117). These studies have significantly contributed to our understanding of the mechanisms of permeation, substrate selectivity, and gating in these channels. AQP research has particularly benefited from close collaborative studies between experiment and computation where either new experiments were designed based on the results of the simulations or hypotheses developed by the simulations were corroborated by experiments (24, 26, 32, 75, 105, 121). Keeping the focus of the review on the application of MD to study permeation processes, we have refrained from a comprehensive discussion on various functional aspects of AQPs. Rather, we use AQPs as an example to demonstrate the application of MD and to showcase the various types of information, both qualitative and quantitative, that can be derived from the method.
Basic Concepts of MD
MD and its application to molecular systems are based on basic concepts of physics and chemistry, in particular classical mechanics and statistical mechanics. From an algorithmic perspective, the method has evolved along with the developments in mathematics and computer science. MD is a numerical method that integrates the Newtonian equations of motion, e.g., for individual atoms in a given molecular system (3, 34). During a “classical” MD simulation, the interactions between the atoms are calculated using a pre-determined set of parameters that define a potential energy function (U), also known as a force field. The potential energy functions of various force fields available for biomolecular simulations could be different, but they all generally include terms that describe major bonded (bonds, angles, and dihedral angles), and nonbonded (van der Waals and electrostatic) modes of interaction: The parameters used for these energy terms are derived from a combination of experimental data and quantum mechanical calculations (41) and are tuned to optimally reproduce the structure and vibrational modes of the molecular systems of interest, as well as their thermodynamic properties. At each time step during an MD simulation, the total force on each atom is calculated as the negative derivative of the potential energy with respect to its coordinates. Using these forces, the Newtonian equations of motion are then integrated to produce the trajectory of individual atoms. The result of an MD simulation is, therefore, a trajectory (collection of snapshots or configurations) of the system over a certain period of time, usually tens to hundreds of nanoseconds. These snapshots can be used to describe the dynamics of the system and to calculate macroscopic properties using principles of statistical mechanics.
FIGURE 2 shows a typical membrane protein simulation system, featuring the water channel protein aquaporin-1 (AQP1) embedded in a lipid bilayer composed of 187 POPE lipids (∼80,000 atoms). Explicit representation of lipid and water in these simulations provides a natural environment for the protein and also allows one to perform the simulations under various pressures when desired. To mimic the experimental conditions, MD simulations are often carried out under constant-temperature and constant-pressure conditions.
Strengths and Limitations of MD
In most all-atom MD simulations, the integration has to be performed every 1–2 femtoseconds (10−15 s) to describe properly the motion of atoms and to avoid numerical instability. As a result, current MD simulations of macromolecular systems are limited to the nanosecond-microsecond time scales, although consistent development in software and advances in computer hardware continue to tackle this limitation, as evidenced by very recent studies reporting simulations on the order of hundreds of microseconds (94). A second shortcoming associated with “classical” MD simulations is their inadequate treatment of electronic polarization effects (74) and simplified potential energy functions, precluding a complete description of interatomic interactions and accurate dynamics of the system. With regard to specific application of the method to membrane systems, such shortcomings have resulted, e.g., in underestimation of area per lipid (30, 31) and, thus, inaccuracies in the density and even the phase of the lipid bilayer under investigation. Moreover, the slow diffusion of lipid molecules continues to keep us from satisfactory sampling of mixed lipid bilayers.
Despite these limitations, MD has been successfully employed in studying a wide range of biomolecular systems and phenomena (66). In particular, and in regard to the subject of this review, membranes and membrane proteins have been widely studied by MD simulations (68, 73, 93, 118).
MD has the advantage of providing dynamic information of a system at very high spatial and temporal resolutions. Starting from an experimentally determined structure, e.g., those solved by X-ray crystallography or NMR measurements, MD allows the researcher to “watch” the natural motion of the molecular system and to monitor its dynamic behavior in real time. Therefore, MD can be considered as the in silico version of a conventional microscope, with the advantage of providing a live view of the molecule at an atomic resolution. Time-averaged properties computed from an MD trajectory can be compared with macroscopic quantities that are measured experimentally (5, 96, 125).
Furthermore, MD offers researchers the opportunity to explore their systems under an unlimited number of artificial conditions that are often inaccessible experimentally. For instance, not only can one modify a residue to mimic the effect of an experimentally studied point mutation, but one can also selectively screen, or even completely neutralize, the electrostatic forces generated by any given group of atoms (11, 102).
Enhanced Sampling Techniques
In the past few decades, MD has been used in the study of a wide range of biological phenomena, e.g., protein folding (25), ion conduction (90), and muscle elasticity (97). The rapid increase of the computational power and recent developments in simulation software have enabled MD studies of significantly larger systems, e.g., the ribosome in complex with a protein-conducting channel (2.7 million atoms) (40), and of much longer processes (hundreds of microseconds) (94). Despite the striking progress, a substantial gap remains between the time scales currently accessible to MD simulations (on the order of nanosecond-microsecond) and the time scale of “slow” biological processes, e.g., those that involve large-scale domain motions of macromolecular systems (microseconds to milliseconds or even longer). To bridge this gap, special simulation and sampling techniques have been developed that allow an enhanced sampling of the configuration space, e.g., by enforcing certain molecular events in the system (5, 14, 38). Here, we briefly describe a number of enhanced sampling techniques that are used by studies reviewed in this article. For a more comprehensive discussion, we refer the reader to a number of excellent reviews on the subject (5, 14, 38).
In steered MD (SMD) (57), external forces are applied to one or a group of atoms to induce certain molecular events during the simulation or to drive the system from one state to another. Some of the first applications of the method to biomolecular processes include mechanical characterization of the muscle protein titin to interpret the results of AFM experiments (37) and describing the unbinding pathway of a ligand from its receptor (59). More relevant examples to this review include SMD simulation of permeation of glycerol through AQPs (62), as well as induction of hydrostatic pressure difference for either water (103, 124–126) or gas molecules (113, 114, 116, 117) across the membrane by adding forces that drive these molecules toward the membrane on one side and away from it on the other side, thus generating a chemical potential gradient across the membrane. Apart from probing the system qualitatively, SMD trajectories have also been used to calculate the free energy difference (ΔG) between the states of the system (62, 85).
Other enhanced sampling techniques that are developed particularly with the purpose of calculating ΔG include free energy perturbation (71, 127), thermodynamic integration (70, 98), and umbrella sampling (89), as well as the more recently developed methods of adaptive biasing force (ABF) (18, 47) and implicit ligand sampling (ILS) (15). The ILS method, which calculates a three-dimensional free energy map from simulations without explicit presence of the ligand (e.g., gas molecules), has proven particularly useful in characterizing gas migration pathways within proteins (15) and lipid bilayers (113, 117) that might not be sampled during a standard MD simulation.
Simulation Studies of Water Conduction by Aquaporins
AQPs are a family of membrane proteins best characterized as water channels, since they increase the membrane's permeability to water and other small hydrophilic substrates. They form tetramers in the membrane, with each monomer acting as an independent water pore (1, 7, 65). An AQP monomer consists of six transmembrane helices and two half-membrane spanning repeats related by a quasi-twofold symmetry. These two repeats, which line a large portion of the water pore, each contain a highly conserved asparagine-proline-alanine (NPA) motif (51) and meet approximately at the center of the pore (FIGURE 3A). In the periplasmic half of the water pore, a region lined by aromatic residues and a conserved arginine forms the narrowest part of the pore, often referred to as the selectivity filter (SF) (91). To date, over 14 structures have been determined for different members of the AQP family, including E. coli GlpF (35, 102) and AqpZ (91), human AQP1 (79), AQP2 (92), AQP4 (53), and AQP5 (54), bovine AQP1 (99) and AQP0 (39), rat AQP4 (52), sheep AQP0 (44), spinach SoPIP2;1 (105), archaeal AqpM (72), yeast AQY1 (32), as well as PfAQP from the malarial parasite Plasmodium falciparum (83). As a result, AQPs have become the richest family of membrane channels with regard to the abundance of atomic-resolution structures.
The first set of MD simulations of AQPs (19, 102) were prompted by the high-resolution structures of GlpF (35) and AQP1 (79) and aimed at understanding the dynamics and selectivity of water permeation. These studies, which employed equilibrium simulation conditions, i.e., in the absence of any biasing force, found that water molecules form a continuous single file and move in a highly correlated manner within the water pores of individual AQP monomers (19, 64, 102). Every water molecule forms on average two hydrogen bonds with neighboring water molecules in addition to a hydrogen bond with the protein (FIGURE 3A). This fine-tuned hydrogen bond network, where water-water interaction is elegantly supplemented by protein-water interaction, allows AQPs to achieve a high efficiency in water conduction. The importance of the natural motions of pore-lining residues in water conduction by AQPs was also clearly demonstrated by simulations in which artificially confining the motion of these residues resulted in a complete block of water conduction (95, 115).
Apart from elucidating the dynamics of water permeation, MD simulations have made a significant contribution to our understanding of the selectivity of AQPs. While conducting water at a rate close to its self-diffusion limit, AQPs strictly exclude protons, a property that is critical to the maintenance of the electrochemical potential across the membrane (122). However, given the high mobility of protons in bulk water, the mechanism of proton exclusion by AQPs remained a long-standing puzzle. A significant number of simulation studies (9, 11, 13, 20, 21, 67, 102) have provided a detailed insight into the fine-tuned electrostatic forces that form the basis for this selectivity mechanism. These studies have shown that water molecules within each water pore form a characteristic “bipolar” orientation (FIGURE 3A) (102). This unique orientation of water reflects a positive electrostatic potential at the center of protein, which is mainly generated by the dipole moments of the two half-membrane spanning helices (11, 20, 102). Together with the desolvation effect (9, 67), the positive electrostatic potential effectively blocks the passage of protons, while leaving the rapid water conduction unaffected (22).
Calculation of the Permeability Coefficient from MD Trajectories
The motion of individual water molecules in MD simulations can be closely monitored. Therefore, all the required information to calculate various kinetic properties of the process of water permeation through AQPs is fully accessible from an MD trajectory. Experimentally, the rate of water conduction is characterized by the osmotic permeability coefficient (pf), a quantity that is measured in the presence of an osmotic (or hydrostatic) pressure gradient across the membrane. Several MD studies reported the calculation of pf from the simulated trajectories (22, 23, 45, 46, 61, 124, 126). Initially, a pressure gradient across the membrane was generated by applying small forces to individual water molecules in the simulation system to bias their translocation and, hence, establish a net flow through an AQP channel (124, 126). Later studies (125) developed theories that related pf to the diffusive behavior of water inside the channel during equilibrium simulations, thereby allowing for the calculation of pf from equilibrium MD simulations. The calculated value of pf for AQP1 from MD simulations is 7.1–7.5 × 10−14·cm3·s−1 (22, 126), which is in good agreement with the experimental value of 5.43 × 10−14·cm3·s−1 (112). Experimentally, pf values on the order of 10−14·cm3·s−1 have been reported for different members of the AQP family (28, 120). A notable exception with this regard is AQP0, in which the water flow appears to be about an order of magnitude slower than in other AQPs (12, 120). Recent MD studies have suggested that the “phenolic barriers” formed by two luminal tyrosines are responsible for the reduced rate of water transport by AQP0 (42, 60).
Energetics of Substrate Permeation Through the Water Pores
Given adequate sampling, the energetics associated with the molecular phenomenon at hand can be readily calculated from MD simulations. Since water permeation through the water pores of AQPs happens on a nanosecond time scale, all reported equilibrium MD simulations (21, 45, 63, 84, 101, 102, 105, 113, 117, 123) have been able to collect ample sampling on the dynamics of water inside the water pores. The relatively fast translocation of water in the water pores has resulted in adequate sampling of all points along the pore axis in all reported equilibrium MD simulations (21, 45, 63, 84, 101, 102, 105, 113, 117, 123) from which a free energy profile for water permeation can be reconstructed based on the probability distribution of water along the pore. Calculations based on this method have shown that water permeation through the water pores of AQPs requires crossing barriers of ∼3 kcal/mol (see Table 1) (55, 117).
In addition to water, some AQPs are also permeable to other hydrophilic compounds. In particular, a subfamily of AQPs known as aquaglyceroporins are able to conduct glycerol and other small linear sugar molecules. The energetics associated with glycerol permeation through the E. coli aquaglyceroporin GlpF has also been calculated from MD simulations (48, 62). In this case, however, the significantly slower diffusion of glycerol along the pore did not allow the entire process, i.e., translocation of a glycerol molecule from one end of the channel to the other, to be captured in equilibrium MD simulations (63). As such, enhanced sampling methods had to be adopted to ensure glycerol occupancy of all points along the pore axis to reconstruct the energetics associated with its permeation (48, 62). Jensen et al. performed SMD simulations to induce the permeation of glycerol with an external force (62). The free energy profile calculated from these simulations revealed a barrier of 7.3 kcal/mol (62), which is comparable with the 9.6 ± 1.5 kcal/mol Arrhenius activation energy estimated from experiments (8). Using the ABF method (18, 47), Henin et al. calculated a comparable energy barrier of 8.7 kcal/mol (48). We note that a notably smaller energy barrier (3.2 kcal/mol) was reported in a recent umbrella sampling study (55). Both Jensen et al. and Henin et al. also examined the stereoselectivity of the channel and showed that the gauche-gauche conformation of glycerol is favored at the SF region, with the non-polar backbone of glycerol facing the hydrophobic side of the channel and its three hydroxyl groups forming hydrogen bonds with the polar residues (FIGURE 3C). Accordingly, it was suggested that the search for an optimal conformation at the SF region constitutes the rate-limiting step of glycerol permeation through GlpF (48).
Gating of the Water Pores
Although largely known as always-open water channels, gating of the water pores has been reported for a number of AQPs, most prominently for plant AQPs (39, 75, 82, 106). Gating of AQPs is believed to provide the cell with a mechanism to modulate their permeability under conditions such as drought and flooding. It has been shown, for instance, that plant AQPs can be regulated by phosphorylation of highly conserved serine residues or by variations of the intracellular pH (39, 82, 106). A recent combined experimental and MD study on yeast AQY1 suggests that AQPs might also be gated by changes in membrane curvature or surface tension, a mechanism reminiscent of the gating of the mechanosensitive ion channels (32). Regardless of the gating signal, a rather consistent gating mechanism has emerged from the structural and MD studies of these AQPs (FIGURE 4). In spinach AQP, SoPIP2;1, the phosphorylation-mediated gating has been linked to the coupled movement of a pore plug (a conserved leucine residue that lines the pore in the cytoplasmic half) and one of the cytoplasmic loops that detects the phosphorylation signal or the intracellular pH change (69, 84, 105). In the yeast AQY1, a conserved tyrosine acts as a plug, whereas the extended N-terminus has been suggested to act as the sensor of the membrane tension (32). As discussed in recent structural and MD studies (33, 69, 84, 109), gating of AQPs may be mediated by a diverse set of mechanisms that couple the motion of different signal-detecting domains (sensors) to the pore plug.
Lipid- and Protein-Mediated Transmembrane Gas Conduction
Since the first report on the CO2 permeability of human AQP1 (86), the involvement of AQPs in gas conduction across cellular membranes has attracted marked attention from both experimentalists (6, 16, 17, 27, 29, 43, 49, 50, 58, 78, 80, 81, 104, 107, 108, 110, 119) and theoreticians (55, 56, 113, 117). Several studies have suggested that AQPs can conduct gas molecules across biological membranes and that this conduction can be of physiological relevance. For instance, the effect of tobacco AQP, NtAQP1, on photosynthesis and stomatal opening has been attributed to its role in the facilitating transmembrane CO2 permeation (107). More recently, it was demonstrated that human AQP1 may conduct NO when reconstituted in liposomes (50). Later, through the measurement of acetylcholine-induced vasorelaxation in AQP1 knockout mice, the direct involvement of AQP1 in NO-dependent vasorelaxation was demonstrated in vivo (49). Apart from participating in gas permeation, relative CO2 and NH3 permeability measurements have suggested that AQPs might even exhibit selectivity for different gas molecules (80). However, contradicting bodies of evidence also exist (29, 78, 110, 119). For example, it has been shown that AQP1 deletion has no effect on CO2 transport in erythrocytes, lungs, or kidneys (29, 119).
At the heart of the debate over the role of AQPs in gas conduction are two questions: first, whether AQPs can conduct gas molecules, and second, whether AQP-mediated gas conduction is of physiological relevance. These questions turn out to be very difficult to answer experimentally, primarily due to the fact that gas permeation through AQPs is often masked by the high gas permeability of the surrounding artificial lipid bilayer used to reconstruct the protein. Although the answer to the second question (physiological relevance) remains to be determined experimentally, most likely in in vivo settings MD simulations provide a powerful technique to address the first question, namely the availability of gas conduction pathways in AQPs, or in other membrane proteins for that matter. A number of reported MD studies (55, 56, 113, 117) have investigated this problem, and it is now well established that AQPs indeed are permeable to gas molecules. Interestingly, these studies have shown that the central pore of AQPs (the pore formed between the four monomers; see FIGURE 2) provides a more favorable pathway for hydrophobic gas species than the water pores. The calculated energy barrier against O2 permeation through the water pores is found to be 5–6 kcal/mol, whereas the barrier for the same process through the central pore is only 3–4 kcal/mol. A summary of the energy barriers against gas permeation through AQPs calculated by MD simulations is given in Table 1. Below, after summarizing the results of the MD simulations of gas permeation through lipid bilayers and AQPs, we will discuss their implications with regard to the physiological relevance of AQP-mediated gas conduction.
Gas Permeability of the Water Pores
Using umbrella sampling, Hub and de Groot (56) showed that the highest energy barrier against CO2 permeation through the water pores of AQP1 is 5.5 kcal/mol and is located at the SF region (see FIGURE 3). Using a combination of complementary sampling techniques, i.e., explicit ligand sampling, implicit ligand sampling (113), and umbrella sampling (117), comparable energy barriers (∼5 kcal/mol) were calculated for the permeation of O2 and NO permeation through the water pores of AQP1 (113) and AQP4 (117) (FIGURE 5C). Comparison of the energy barriers against gas permeation through the water pores calculated by various MD simulations (Table 1) shows that, despite the difference in methodology, these studies have yielded a rather consistent picture, namely an energy barrier of 5–6 kcal/mol for the permeation of hydrophobic species (CO2, O2 and NO) through the water pores of AQPs and a lower barrier (4 kcal/mol) for the more hydrophilic NH3. These energy barriers are all located at the SF region (FIGURE 5C) and appear to arise from the loss of favorable water-water and water-protein interactions in the SF (55) on the insertion of the gas molecule into the single file of water within the water pores (FIGURE 5F). The more hydrophilic NH3 partially compensates for this loss by NH3-protein interactions and, hence, faces a smaller energy barrier during its permeation (55).
It is noteworthy that the water pore of the aquaglyceroporin GlpF is found to have a significantly lower energy barrier against gas permeation (55). As shown in Table 1, a barrier of only ∼3 kcal/mol is obtained for the permeation of different gas species through the water pores in GlpF (55). This low and non-discriminant barrier is explained by the wider and more hydrophobic SF region of GlpF compared with pure water channels such as AQP1 and AQP4. Indeed, it has been shown that introduction of more hydrophobic residues in the SF of AQP1 can dramatically reduce the energy barrier against gas permeation (55). In humans, AQP3, AQP7, and AQP9 are aquaglyceroporins (2). Although no atomic structures are currently available for these AQPs, they share a high level of sequence identity with GlpF, especially with regard to the SF region. Therefore, comparably high gas conduction properties can be expected for these human AQPs.
Gas Permeability of the Central Pore
Besides the four monomeric pores, a fifth pore, named the central pore, is formed in the middle of the four monomers in an AQP tetramer (FIGURE 2). Although the function of the water pores has been well characterized, no functional implication has been experimentally established for the central pore. The central pore is largely lined by hydrophobic residues, which, along with the small radius of the pore, result in an effective exclusion of water molecules. During the MD simulations, the central pore stays completely dry, effectively blocking the entrance and permeation of water molecules. Therefore, it is unlikely for the central pore to contribute to the experimentally measured overall water conduction rate of AQPs. The MD simulations, however, have revealed a novel property of the central pore, namely its high affinity for gas molecules; the pore acts as an ideal gas reservoir that is rapidly filled with gas molecules during the simulations. Furthermore, it serves as the most favorable pathway for permeation of hydrophobic gas molecules through the protein, specifically when compared with the water pores (56, 113). This property has been demonstrated by equilibrium MD simulations of AQP1 and AQP4, where the gas molecules initially placed in solution were shown to rapidly accumulate in the central pore (FIGURE 3, D–F) during the simulations (113, 117). Furthermore, free energy calculations have revealed an energy barrier of only 3–4 kcal/mol (56, 113, 117) against the permeation of CO2, O2, and NO through the central pore, supporting the view that the pore indeed serves as a favorable gas conduction pathway.
Physiological Relevance of AQP-Mediated Gas Conduction
The physiological significance of gas permeation through AQPs will depend on their relative gas permeability compared with the surrounding lipid bilayers. Recent experimental studies offer different views with this regard. Endeward et al. suggested that AQP1 is a major pathway for CO2 transport across the human erythrocyte membrane (27), whereas Missner et al. reported that the resistance to transmembrane CO2 permeation is mainly caused by the unstirred layer, instead of the membrane itself, thereby suggesting that facilitation of CO2 transport by AQPs is not physiologically important (78). MD simulations performed on artificial lipid bilayers containing, e.g., 100% POPE or POPC show that these pure phospholipid bilayers produce a very small energy barrier against the permeation of hydrophobic gas molecules such as O2 and CO2 (55, 56, 113, 117). For instance, the energy barrier of a pure POPE bilayer to O2 permeation is only 0.4 kcal/mol (117), which is essentially negligible compared with the barrier of O2 permeation through AQPs. The only exception to this seems to be the hydrophilic NH3, in which case a barrier of 4.5 kcal/mol for its permeation through a pure lipid bilayer has been calculated (55). These results support the notion that AQPs are unlikely to be of significance in the overall gas exchange across a “model” membrane. Nevertheless, “real” cellular membranes can differ dramatically, both in lipid composition and structure, from the model membranes used in the simulation and experimental studies. Therefore, the physiological relevance of AQP-mediated gas conduction will remain unknown until the composition of the lipid bilayers embedding AQPs in a living cell has been determined and used accordingly in experimental or computational studies.
It should be also noted that, despite a lower permeability of AQPs for gas molecules compared with pure lipid bilayers, AQPs might occupy a significant portion of the membrane surface. Under these conditions, the AQP-mediated gas permeation might become physiologically significant. AQP4 presents a case in which this scenario might apply. AQP4 is the predominant water channel in the central nervous system (CNS) (76, 111), where it is most abundantly present in astroglial cells that surround capillaries and form glia limitans, as well as in ependymal cells lining ventricles (4, 52, 100). Freeze fracture and immunogold labeling experiments have shown that AQP4 forms extensive two-dimensional square arrays in the membrane of these cells (87, 88), where large sets of tightly packed AQP4 proteins occupy a significant portion of the membrane. The high density of AQP4 in these cells suggests that a non-trivial fraction of gas exchange in certain regions in the CNS might take place through pathways provided by the protein. Interestingly, recent experimental studies report that AQP1 can mediate NO permeation both in vitro and in vivo (49, 50). Consistent with the intriguing notion that AQP4 might play a role in the conduction of this important signaling molecule in the CNS, MD simulations have found a small energy barrier against NO permeation through AQP4 (117).
In this article, we have reviewed some of the recent MD studies of membranes and membrane proteins, specifically seeking to showcase the power of the method in describing molecular events involved in the translocation of small molecules across the cellular membrane. As exemplified, a wide range of phenomena can be captured at a high resolution by MD simulations. The calculated trajectories have been used to reveal, e.g., a novel and unique configuration of water molecules in AQPs and the basis of substrate selectivity and proton exclusion. MD simulations continue to provide a powerful platform to develop and test new hypotheses regarding the relevant functional aspects of the system at hand, for instance, the involvement of electrostatic forces in the selectivity of AQPs. On the more quantitative side, MD simulations have been used to calculate key thermodynamic and kinetic properties, e.g., the permeability coefficients of different AQPs and the free energy profiles associated with the permeation of various substrates along different conduction pathways.
With regard to the exchange of very small molecules (water, CO2, O2, and NO) across the membrane, there is currently no evidence suggesting the involvement of active (energy-dependent) transporters in the process. Therefore, the exchange of these molecules across the cellular membrane appears to be merely mediated by their passive diffusion through transient or permanent pathways in the membrane. Such pathways can be in principle provided by either the lipid bilayer itself, or by proteins that are embedded in the membrane. For charged molecules, as well as for the majority of polar molecules, the hydrophilicity of the pathway is critical for their efficient permeation, an aspect that explains why the conduction of such species is exclusively through membrane proteins (channels). Very small molecules, in particular hydrophobic gas molecules, on the other hand, seem to favor hydrophobic translocation pathways, e.g., transient pores that constantly form between the lipid molecules due to their structural fluctuation, or hydrophobic pores and pathways within membrane proteins that are expected to be more permanent in nature due to the relatively stable structure of proteins compared with lipids. Extensive MD simulations of AQPs have provided strong evidence that AQPs, and likely other membrane proteins, can indeed conduct gas molecules, although the physiological significance of the phenomenon remains to be established. It is interesting that the central (tetrameric) pore of AQPs, for which no function has been demonstrated experimentally, turns out to be the most favorable conduction pathway for gas molecules. The formation of an additional pore that is involved in the permeation of substrates other than the main substrate of the channel (water) might even be speculated to be one of the reasons for oligomerization of AQPs, and maybe for other oligomeric membrane proteins.
The authors acknowledge support from the National Institute of General Medical Sciences (grants R01-GM-086749 and R01-GM-067887), as well as computational resources provided by TeraGrid (grant no. MCA06N060).
- ©2010 Int. Union Physiol. Sci./Am. Physiol. Soc.