Asymmetric Ligand Binding Facilitates Conformational Transitions in Pentameric Ligand-Gated Ion Channels
The anesthetic propofol inhibits the currents of the homo-pentameric ligand-gated ion channel GLIC, yet the crystal structure of GLIC with five propofol molecules bound symmetrically shows an open-channel conformation. To address this dilemma and determine if symmetry of propofol binding sites affects the channel conformational transition, we performed a total of 1.5 As of molecular dynamics simulations for different GLIC systems with propofol occupancies of 0, 1, 2, 3, and 5. GLIC without propofol binding or with five propofol molecules bound symmetrically showed similar channel conformation and hydration status over multiple replicates of 100-ns simulations. In contrast, asymmetric binding to one, two or three equivalent sites in different subunits accelerated the channel dehydration, which was accompanied by increased conformational heterogeneity of the pore and shifted the lateral and radial tilting angles of the pore-lining TM2 towards a closed-channel conformation. The results differentiate two groups of systems based on the propofol binding symmetry. The difference between symmetric and asymmetric groups is correlated with the variance in the propofol-binding cavity adjacent to the hydrophobic gate and the force imposed by the bound propofol. Asymmetrically bound propofol produced greater variance in the cavity size that could further elevate the conformation heterogeneity. The force trajectory generated by propofol in each subunit over the course of a simulation exhibits an ellipsoidal shape, which has the larger component tangential to the pore. Asymmetric propofol binding creates an unbalanced force that expedites the channel conformation transitions. The findings from this study not only suggest that asymmetric binding underlies the propofol functional inhibition of GLIC, but also advocate for the role of symmetry breaking in facilitating channel conformational transitions.
INTRODUCTION<!>METHODS<!>Asymmetric propofol binding increased the probability of channel dehydration<!>Asymmetric propofol binding facilitated the pore-lining TM2 toward a closed-channel conformation<!>Propofol motion and imposing forces affect the channel hydration status<!>Additional remarks<!>CONCLUSIONS
<p>Structural symmetry of a protein assembly results from its functional evolution.1–4 The symmetry feature is often required to retain global structural stability and cooperative functionality.1,3,5,6 Cys-loop receptors mediate fast synaptic signal transmission. Each Cys-loop receptor is comprised of five homologous subunits that form a pentameric ligand gated ion channel (pLGIC). For a homo-pLGIC, a fivefold symmetry around the central pore is assumed. Each subunit contains an extracellular (EC) domain, a transmembrane (TM) domain of four TM helices (TM1 to TM4) with TM2 lining the pore, and an intracellular domain that links TM3 and TM4. Agonist binding to the EC domain induces channel opening and allows ions to move through the pore. Despite the existence of five identical agonist-binding sites in a homo-pLGIC, occupancy of three nonconsecutive sites by agonists was found to induce maximal mean channel open time.7 The maximum channel-gating efficacy could be reached when only three potential binding sites were occupied in the homomeric α1 glycine receptors (GlyRs)8 and ρ1 GABAA receptors.9 Skepticism remains as to whether asymmetric agonist binding induces spontaneous asymmetric conformational changes. In addition, it remains unclear whether ligands other than agonists also bind asymmetrically to these proteins to produce functional impact. Interestingly, an asymmetric intermediate conformation of a homo-trimeric transporter was captured recently in the crystal structure,10 indicating the involvement of asymmetric conformational change in biological function.</p><p>General anesthetics modulate the functions of pLGICs. At pharmacologically relevant concentrations, general anesthetics potentiate anion-selective GlyRs and GABAA receptors, but inhibit cation-selective nicotinic acetylcholine receptors (nAChRs) and serotonin receptors.11,12 Occupancy of a single binding site in the homo-pLGIC α1GlyR by volatile anesthetics and alcohols was found to be sufficient to potentiate channel currents.13 General anesthetics14 and alcohols15 can also modulate the function of the bacterial Gloeobacter violaceus pLGIC (GLIC). Crystal structures of GLIC in complex with propofol or desflurane16 reveal intra-subunit binding sites in the TM domain of all five subunits (Fig. 1). The structures of the anesthetic-GLIC complexes are virtually identical to the apo GLIC structure.16,17 GLIC is a proton-gated cationic channel with a half maximal effective concentration (EC50) at ~pH 5.18 It shares many pharmacological properties with eukaryotic cationic pLGICs, including inhibition by the intravenous anesthetic propofol.19–22 Propofol inhibits GLIC currents at concentrations used clinically.14,16 A higher propofol concentration could completely close the GLIC channel and inhibit ion conductance.14,16 Therefore, an apparently open channel structure of GLIC under symmetric anesthetic occupancy in the crystal structure seems incongruent with the potent inhibition observed in the functional measurements. The relevance of the identified anesthetic site in the structure, however, is supported by functional studies on various mutants.16 The question is whether inhibition of GLIC requires anesthetics to occupy all five subunits simultaneously or only a few subunits similar to channel activations by asymmetric agonist binding.7–9</p><p>To address this apparent disagreement between structure and function, we performed multiple sets of molecular dynamics (MD) simulations on the crystal structures of GLIC and the propofol-GLIC complex. The number of propofol molecules bound to GLIC was varied in different simulation systems by either keeping all five propofol molecules or deleting propofol from some of the subunits before the simulations. Two groups with distinct channel hydration states and conformations emerged over the course of MD simulations. GLIC with symmetric propofol occupancy in all five sites and the apo GLIC acted as one group, while GLIC with asymmetric propofol binding belonged to another group. The study suggests that symmetry of ligand binding has a profound effect on conformational transitions. Symmetry breaking by ligand binding facilitates conformational transitions. In general, symmetry breaking is a prevalent process in biology and symmetry breaking along well-defined axes is often linked to functional diversification on every scale.23 The current simulation study along with previous knowledge of asymmetric agonist binding in activation of Cys-loop receptors7–9,13 demand thorough characterizations of symmetry breaking by ligand binding in the functions of pLGICs, just as those characterized for many other biological systems.23</p><!><p>Crystal structures of the open-channel apo GLIC (PDB code: 3EAM) and the open-channel propofol-GLIC complex (PDB code: 3P50) were used for MD simulations. Crystal structures of the locally closed GLIC (PDB codes: 3TLS and 3TLW) were used as references for the closed-channel conformations observed from the simulations. Five simulation systems were generated by varying the propofol occupancy in GLIC: (i) no propofol bound to GLIC (0PFL); (ii) five propofol molecules bound to GLIC (5PFL) as shown in the X-ray structure of the propofol-GLIC complex;16 (iii) three propofol molecules bound to non-consecutive subunits (3PFL); (iv) two propofol molecules bound to consecutive sites (2PFL); and (v) one propofol molecule bound to GLIC (1PFL). For each system, three parallel runs were performed using different seed numbers and each run lasted for 100 ns. Fig. 1 shows all five systems and the propofol binding sites.</p><p>Protonation states of titratable residues in GLIC at pH = 4.6 were assigned based on the results reported by Bocquet et al.17 Some modifications, including deprotonation of five H235 residues and two E222 residues, were made based on our recent calculations.24 The TM domain of GLIC was inserted into a pre-equilibrated and solvated POPE/POPG (3:1) binary lipid mixture. Each simulation system has a hexagonal boundary condition of 104.6 Å × 104.6 Å × 129.8 Å, one GLIC, 167 POPE, 54 POPG, and approximately 23,700 TIP3 water molecules.</p><p>MD simulations were performed using NAMD 2.7b1.25 CHARMM27 force field with CMAP corrections (version 31) was used for protein, water, and lipids.26,27 Propofol parameterization was done following the protocol of the CHARMM General Force Field (CGenFF) for drug-like molecules.28 Details of propofol parameters are provided in the online supporting material (Fig. S1, Table S1–S4). All simulation systems followed the same simulation protocol. Each system was energy minimized for 20,000 steps before equilibration for 2 ns, during which the backbone constraint on GLIC was gradually reduced from 10 kcal·mol−1·Å−2 to zero. Each system underwent three runs up to 3×100 ns of unconstrained simulations under constant pressure (P = 1 Bar) and temperature (T = 310 K).29,30</p><p>Periodic boundary conditions, water wrapping, hydrogen atoms constrained via SHAKE, and evaluation of long-range electrostatic forces via the Particle Mesh Ewald (PME) algorithm31 were used in the simulations. Bonded interactions and short-range, non-bonded interactions were calculated every time step (2 fs). Electrostatic interactions were calculated every two time steps (4 fs). The cutoff distance for non-bonded interactions was 12 Å. A smoothing function was employed for the van der Waals interactions at a distance of 10 Å. The pair-list of the nonbonded interaction was calculated every 20 time-steps with a pair-list distance of 13.5Å.</p><p>VMD32 was used for visualization and most parts of data analysis. Unless otherwise specified, snapshots every 20 ps of the simulation trajectories (a total of 5000 snapshots) were used for data analyses for each 100-ns simulation. The number of water molecules (Nwater) inside the hydrophobic gate region was obtained by counting water inside the pore between I233 (I9′) to I240 (I16′). The channel hydration state was defined based on Nwater: fully hydrated if Nwater ≥10, partially dehydrated if 0<Nwater <10 and fully dehydrated if Nwater = 0. Ten water molecules inside the hydrophobic gate region are equivalent to 65% of the bulk water density and were used previously as a threshold in the evaluation of the channel hydration status.33 A histogram of Nwater was calculated based on the data from three 100-ns parallel runs for each system. Each run was named according to the duration of time that the channel remained hydrated. For example, the channel hydration time in the 5PFL system follows the order: 5PFL-1 > 5PFL-2 > 5PFL-3.</p><p>Orientation of the pore-lining TM2 helix was characterized by the radial (θ) and lateral (ϕ) tilting angles of the TM2 helices relative to the membrane normal, as defined in previous publications.34,35 The same method as detailed previously35 was used to calculate the radial and lateral tilting angles. For each system, distributions and histograms of the radial and lateral tilting angles were calculated for each channel hydration state and averaged over all five subunits in three replica simulations.</p><p>A normalized histogram for the joint events of (θ, ϕ) was used to estimate heterogeneity of the TM2 tilting angles using MATLAB®. The joint Shannon entropy S(i,j) was calculated by:36</p><p>
(1)S(i,j)=-∑j=1N∑i=1Npijlnpij where pij is the joint probability of the event (θi, ϕj) obtained from the normalized histogram of the tilting angles, and N is the number of bins. For each system, the radial and lateral tilting angles were collected over all five subunits in three replica simulations. To generate the histogram, a bin size of 0.1° was used to sample angles ranging from −12° to 15° for both θ and ϕ. Standard deviations of the entropies were estimated with the bootstrap method,37 using 100 sets of randomly sampled data points from the original 2600 × 3 snapshots for each system. The analysis of variance with post-hoc Tukey HSD multiple comparison test between groups was performed using SPSS v20.</p><p>The anisotropic network model (ANM)38 was used for structure-based analysis of GLIC dynamics. The Hessian matrix was built using all Cα atoms and a pairwise interaction cutoff of 15 Å. The ten lowest frequency modes of ANM were calculated and visualized using ProDy.39</p><p>The cavity size of the intra-subunit propofol-binding pocket was calculated using the POVME algorithm.40 Frames taken at 100 ps intervals from the simulation trajectories were used in the calculation. A grid encompassing the entire binding pocket was generated with 0.5-Å spacing. The grid points not occluded by protein atoms and attached through a series of adjacent grid points to the center of the binding site were used for the cavity volume calculation. If no grid points remained in the center of the site, the volume was reported as 0 Å3.</p><p>The force, resulting from VDW and electrostatic interaction of propofol with GLIC, was calculated in x, y, z directions using the pair Interaction module implemented in NAMD 2.7b1.25 Residues within 5 Å of propofol were selected in the calculations. In order to quantify the primary component of the force, we performed principal component analyses (PCA) on the force trajectory (5,000 frames over 100 ns) using MATLAB. For visualization, the principal component of the force is presented within the context of the GLIC structure. Each principal component was scaled by its eigenvalue and multiplied by the inverse of the eigenvector matrix. Each resultant vector was then centered on the propofol position and plotted.</p><p>For statistical analyses of the channel opening probability, we categorized a channel as "open" if the number of water molecules inside the channel gate is equal to or greater than 10; otherwise, the state of the channel is labeled as "closed". We pooled all 3 replicated simulations for each system (snapshots every 0.5 ns of the simulation trajectories were used) and represented the number (X) of simulation snapshots that assumed an open state using a binomial distribution: X ~ Bin(N, p), where N is the total number of simulation snapshots and p is the channel opening probability. The estimate of the channel open probability (p̂), confidence interval for p̂, the adjusted confidence interval for comparing two estimated channel opening probabilities, (p̃1 − p̃2), and therein the p-value were derived using an R language script and the methods described by Agresti and Caffo.41 The open-channel probability was estimated by p̂ = X/N, the confidence interval for the estimate of channel opening probability is 
p^±zα/2p^(1-p^)/N, and the adjusted confidence interval for the estimated difference between two channel opening probabilities, (p̃1 − p̃2), is</p><!><p>The channel hydration status is strongly related to the ion conductance. Ions cannot pass through a dry channel even if the pore is not yet geometrically closed.42 An evolution from a fully hydrated to a completely dehydrated channel in simulations signifies a transition of the channel functional state. The GLIC crystal structures used for our simulations, in the absence and presence of propofol, have the same open-channel conformation.16,17,43 At the beginning of all simulations, the GLIC channel was fully hydrated. Over the course of the simulations, GLIC exhibited only a small deviation from the crystal structure. The Cα root mean square deviations (RMSD) of the overall protein and the TM domain reached a plateau after 5 to 8 ns. The RMSDs remained under 2Å for the TM domain thereafter (Fig. S2). The duration of the fully hydrated state, however, varied by runs, particularly by systems. As shown in Fig. 2a, 5PFL-1 remained full hydration almost for the entire 100 ns, while 5PFL-2 and 5PFL-3 were fully hydrated for 65 and 50 ns, respectively. The system with no propofol, 0PFL, showed almost the same results (Figs. S3 and S4). In contrast, 3PFL, 2PFL, and 1PFL experienced much more rapid and extended channel dehydration (Figs. S3 and S4).</p><p>To objectively evaluate and compare the channel hydration status among different systems, we made histograms for each system based on the data from all three replicate runs. Figs. 2b–2f show that the distribution of the number of water molecules inside each channel gate across different simulations is bimodal, with one peaking near 20 waters and the other peaking at zero water. We hypothesize that these two peaks correspond to the "open" and "closed" states of the channel, respectively. Based on this assumption, we used the value of 10 waters per channel—the groove between the two modes—as a threshold to categorize the state of the channel in each snapshot as "open" (≥ 10 waters per channel) or "closed" (otherwise). The channel open probabilities for the symmetric systems 0PFL and 5PFL are 68±6% and 68±6%, respectively. In contrast, the open probabilities for the asymmetric systems 3PFL, 2PFL, and 1PFL are only 42±7%, 28±6%, and 32±6%, respectively. The errors represent 99.9% confidence intervals of the estimates. Using the same level of confidence interval for the differences between the estimated probabilities,41 we found that the channel open probability is significantly higher in the symmetric systems than in the asymmetric systems with a p-value < 0.001. The results clearly differentiate two groups of systems that are divided based on the propofol binding symmetry. Asymmetric propofol binding facilitated the transition from a fully hydrated channel to a dehydrated channel.</p><p>GLIC is expected to be mostly open and hydrated at pH 4.6.18 The channel in 0PFL or 5PFL was indeed hydrated for most of the simulation time. The dehydration occurred and resembled the observations from several previous simulations on GLIC.35,44–46 Since GLIC is a proton-gated channel, imperfect imitation of pH conditions within the simulations could be one of the reasons to cause dehydration. In reality, we do not have a complete set of titratable residues that has been experimentally proven to be responsible for pH activation of GLIC. Accurately predicting pKas in membrane proteins remains challenging. Other imperfections in the simulation environment may also have compromised the timescale of channel hydration. However, for all of the simulated systems reported here, everything was set identical except the number of propofol molecules. Any consequences induced by system imperfections are systematic and common to all the systems. Thus, it is reasonable to conclude that the observed propensity of changes in channel hydration status in Fig. 2 result primarily from variation in propofol occupancy.</p><p>A significantly higher probability of dehydration promoted by asymmetric propofol binding is consistent with the inhibitory effect of propofol on GLIC.14,16 It is known that propofol also inhibits cationic currents of nAChRs19,21 and 5HT3 receptors.22 The results observed here on GLIC are likely relevant to the inhibitory effect of propofol on these receptors as well.</p><p>The finding that asymmetric propofol occupancy facilitated the channel transition more effectively than 5PFL does not contradict the concentration dependence of anesthetic inhibition of GLIC.14 The crystal structure presents the maximum binding sites. In order to make propofol observable in the crystal structure, propofol was added to GLIC in a saturating amount for crystallization. In functional measurements, propofol inhibited GLIC in concentrations (IC50: 0.5~10 μM) several orders of magnitude lower than that used for crystallization.14,16 In addition, the Hill coefficient of propofol inhibition was only 0.42,14 an indication of negative cooperative process, in which one propofol bound to GLIC reduces the ability of another propofol to inhibit channel. All of these suggest that propofol completely inhibits GLIC at a concentration well below the saturated concentration and before it occupies all five sites. Once the channel falls into a closed state, an additional amount of anesthetics cannot resume the GLIC current (our own unpublished data). The crystal structure of the open channel GLIC with five propofol molecules symmetrically bound may only reflect a preferred conformation of GLIC at the crystallization conditions. The discrepancy between the functional state and the channel conformation captured in the crystal structures has also been observed recently on the Erwinia chrysanthemi pLGIC (ELIC).47</p><!><p>It has long been proposed that TM2 helix tilting underlies the channel gating of pLGICs.48 The lateral (δ) and radial (θ) tilting angles of the pore lining TM2 helices give quantitative measurements of pore conformational changes.34,35 Compared to the closed-channel ELIC (δ ≈ −7.9°; θ ≈ −3.5°), the open-channel GLIC (δ ≈ 0.7°; θ ≈ 6.5°) has no more than 10° difference on both angles.43,46,49,50 A combined crystallographic and functional study revealed the locally closed conformations of GLIC, demonstrating that a few degree changes in the lateral and radial tilting angles of TM2 are sufficient to stop GLIC current.51</p><p>Fig. 3 shows the lateral and radial tilting angles of TM2 for each of the five systems simulated to 100 ns. The angles are colored in green, purple, and black to represent three channel states: fully hydrated, transitional, and fully dehydrated at the hydrophobic gate region. For comparisons, the lateral and radial tilting angles of TM2 in the crystal structures of the open-channel GLIC43,46 and the locally-closed-channel GLIC51 are also marked in Fig. 3. For the fully hydrated state, the lateral and radial tilting angles of TM2 are highly populated near 0.7° and 6.5°, respectively, more or less the same as those in the open channel GLIC. Quantitative information about populations of tilting angles for each system is shown by histograms in Fig. S5. A deviation by a few degrees in either the lateral or radial tilting angle could result in channel dehydration. The TM2 tilting angles associated with the dehydrated state are largely shifted toward the angles in the locally closed GLIC.51</p><p>To evaluate conformational changes in other regions of GLIC, we also calculated lateral and radial tilting angles for TM1, TM3, and TM4 in all of the simulations. The same as indicated in the crystal structures of the open- vs. locally closed-channel GLIC (Fig. 3g), TM2 correlates most sensitively to the channel state, while other TM helices show much smaller conformational differences for different channel states, especially TM3 and TM4 (Fig. S6). A higher sensitivity of TM2 conformations to channel state is also shown in structural RMSD clustering analysis (Fig. S7). Pair-wise RMSD matrices over backbone atoms in the TM domain show good separation between different channel states. Bundles of clustered structures from open and closed channels show distinct separation of side chains of residues I240 and A237 at the hydrophobic gate region as well as TM2 between the hydrated and dry channels. Although the conformation differences for different channel states are also visible in other TM helices, they are less distinct than that observed in TM2.</p><p>The symmetry of propofol binding influenced the TM2 tilting angles (Fig. 3 and Fig. S5). Compared to the symmetric 0PFL and 5PFL, the asymmetric 3PFL and 2PFL as well as 1PFL had higher populations of TM2 whose lateral- and radial-tilting angles shifted more towards the locally closed-channel conformation.51 In addition, 1PFL, 2PFL and 3PFL had broader distributions of the TM2 tilting angles than 0PFL and 5PFL, particularly in the lateral angles. The broadness of the distributions reflects the conformational heterogeneity, which can be quantified by the joint Shannon entropy.36 The calculation of the joint Shannon entropies of the TM2 tilting angles using Eq. 1 yielded the values of 7.00±0.01, 6.81±0.01, 6.82±0.01, 6.04±0.01, and 6.31±0.01 for 1PFL, 2PFL, 3PFL, 5PFL, and 0PFL, respectively. One-way ANOVA with respect to propofol occupancy shows that the Shannon entropies for the five systems are significantly different (p < 0.001). Post hoc Tukey HSD comparison tests indicate that 1PFL, 2PFL, and 3PFL are significantly different from 0PFL and 5PFL (p < 0.001). Clearly, asymmetric propofol binding increased the conformational heterogeneity of TM2.</p><p>Although the quaternary twist motion is thought to dominate channel conformation transitions in pLGICs,46,52–54 asymmetric motion has been observed to lead to channel opening55 and closing.52,53 Experimental data also support the role of asymmetric motion in channel functions.56,57 An asymmetric and independent contribution of the TM2 residues to gating of nAChR was observed in a single-channel study.56 Thus, it is not surprising to see spontaneous asymmetric motion in our simulations. After carefully examining individual trajectories, we found that TM2 helices in different subunits experienced different degrees of tilting at a given time point. Moreover, inward radial tilting (or contraction) and/or lateral tilting of the TM2 helix in one or two subunits was sufficient to alter the channel hydration state (Figs. S8–S9). The same phenomena were also observed in other studies on GLIC.35,44,46 We also performed ANM analysis on the crystal structure of GLIC. Among the 10 lowest frequency modes, which are coupled with the large-scale global domain motions, mode 3 is the only mode showing symmetric twisting motion (Fig. S10). Asymmetric motions in other modes, especially asymmetric inward/outward motion among subunits, are expected to contribute to functional changes in channel conformations.53 Thus, asymmetric anesthetic binding does not create new modes of motion, but rather shifts the population of asymmetric motion.</p><!><p>Propofol bound asymmetrically to GLIC facilitates a population shift of TM2 towards the closed-channel conformation and increases conformational entropy. To further understand why asymmetric propofol binding facilitates the conformational transition, we examined changes in the cavity size under different propofol occupancies and the forces imposed by propofol on each subunit.</p><p>The size of the cavity for propofol was quantified using the previously reported method.40 As shown in Fig. 4, for a cavity initially bound with propofol, it expanded to accommodate propofol penetration deeper over the course of simulations. In the absence of propofol, the cavity collapsed over simulations and became noncontiguous (~0 Å3) due to lateral pinching by residues of TM1 and TM3 (Y197, I202, T255, and I258). Histograms of cavity sizes in Fig. 5 show a predominate population at 0 for 0PFL or near 100 Å3 for 5PFL. In contrast, at least two significant populations (0 and ~100 Å3) are observed in 3PFL, 2PFL, and 1PFL. Larger cavity variations in the systems of asymmetric propofol binding are consistent with greater conformational heterogeneity of the TM2 tilting angles found in these systems. It is reasonable to speculate that a larger cavity variance due to asymmetric propofol binding is one of the primary causes for greater TM2 conformational heterogeneity and consequent changes in channel hydration states.</p><p>Propofol binding not only affects the cavity size, but also imposed a force on GLIC. We calculated the force between propofol and residues within 5 Å of propofol and examined the time trajectory of the vector sum of all forces imposed on each subunit (Fig. 6 and Fig. S11). Several characteristics about the force are noteworthy. The force trajectory generated by each propofol over simulation times assumes the shape of an ellipsoid. Based on principal component analyses, the primary component of the force is tangential to the pore and substantially larger than the second component that is mostly radial to the pore. Both the tangential and radial forces are well balanced in 5PFL (Fig. 6), but obviously uneven in the systems of asymmetrical propofol binding. In addition to the force trajectories on all residues within 5 Å of propofol, we also examined the forces on individual residues in the binding cavity, T255 of TM3 and Y197 of TM1. Each residue experienced a force imposed by propofol at any given time. The force fluctuated along a narrow range of directions over the course of simulations. The time averaged net force resulted from dividing an accumulated force by number of snapshots that were sampled evenly from each simulation trajectory. In most cases, the averaged net force became smaller when simulation time became longer. For example, the averaged net force on T255 of 3PFL was 1.7, 1.3, and 1.0 kcal/(mol·Å) at the 10-ns, 50-ns, and 100-ns simulation, respectively. Even though the time averaged net force is not large, the propofol force at any given time is substantial to prevent the binding cavity shrunk as that observed for the empty cavities in the simulations (Fig. 4). Intuitively, an unevenly distributed force creates an unstable condition that could facilitate transitions, either to a direction leading channel closure (such as the case of anesthetic binding in GLIC) or to the direction of channel opening (such as the case of agonist binding in pLGICs). Indeed, when one of the propofol molecules in 5PFL-3 migrated out of the cavity after ~45 ns simulation, the channel was dehydrated soon thereafter.</p><p>The force analysis suggests that propofol imposes a larger force along the tangential direction than the radial direction to the pore. Impact of the force to conformational transition can be substantiated once the force becomes unbalanced among five subunits. Although it seems uncommon to link asymmetrically distributed force with conformational and ultimately functional changes, there is engrained experimental support for the biological significance of symmetry breaking.23 For instance, an asymmetric protrusive force resulting from symmetry breaking in the actin assembly drives directional cell mobility.58,59 The fact that asymmetric agonist binding activates Cys-loop receptors7–9 also speaks for the involvement of a symmetry breaking in the process.</p><!><p>The structures of GLIC bound with anesthetics are so far determined only in the open channel conformation.16,60 Thus, the simulations reported in this study has focused on propofol binding only to the open channel. Can the drug bind to a closed channel? Does the binding to the closed channel induce a structure change? These questions will be addressed with new experiments. To crystallize a fully closed GLIC remains challenging. However, a locally closed channel conformation of GLIC, in the absence of anesthetics, was captured recently in crystal structures by manipulations of cysteine cross-linking.51 This locally closed channel conformation shares most of the features for the open form except a locally closed pore. Its functional relevance during gating transitions has been substantiated by experiments.51 The dehydrated channels observed in our simulations show similar pore conformation of the locally closed channel.</p><p>The propofol binding site in the current study is directly adapted from the reported crystal structure.16 Can anesthetics bind to other parts of GLIC? A more recently published crystal structure reveals that ketamine binds to an inter-subunit cavity in the EC domain of GLIC, which is distinctively distant from the intra-subunit cavity for propofol binding.60 Other experimental and computational studies also suggest a possibility of multiple anesthetic binding sites in GLIC.35,45,61,62 A recent computational study by LeBard and colleagues has not only predicted propofol inhibition by blocking the open pore of GLIC, but also offered quantitative assessment of propofol binding affinity in a low μM concentration.62 The pore-blocking mechanism is worth noting and will be verified by future experiments.</p><!><p>The most important conclusion from this study is that transient symmetry breaking by asymmetric ligand binding in pLGICs facilitates changes in channel conformation. The simulations for systems that are otherwise identical except the number of propofol occupancy provide compelling evidence to support the conclusion. Binding without perturbing symmetry in 5PFL conserved the open-channel conformation as observed in 0PFL. The result is consistent with the consensus that the symmetry feature retains global structural stability.1,3,5,6 In contrast, asymmetric propofol binding perturbed the symmetry and facilitated conformational changes. The distinct difference resulting from asymmetric ligand binding does not come as a total revelation. As Blundell and his colleague stated based on their examinations of several enzymes, mild perturbation from perfect symmetry may be essential in some systems for dynamic functions.1 It is also known that asymmetric agonist binding without occupying all five equivalent sites can produce the maximal opening of Cys-loop receptor channels.7–9 It is likely that, no matter channel activation or inhibition, asymmetric ligand binding works more effectively to induce transitions from one state to another.</p><p>Our multi-seeded, parallel simulations exceeded 1.5 μs in total. Although extending each set to the μs timescale46,63 is desirable, such extensions for multiple μs simulations to cover all possible anesthetic-binding scenarios demand much more computational power that has not been available to us. Fortunately, the current simulation timescale is able to cover the transitions between different channel hydration states. Our statistical approaches with multiple independent runs have sufficient power to unequivocally differentiate the functional propensities of GLIC under different scenarios of propofol binding. To generalize the functional role of asymmetric ligand binding in pLGICs requires further experimental investigations, which may be challenging but are not impossible. The most encouraging examples are the elegant experimental demonstrations that asymmetric agonist binding activates homo-pLGICs7–9 and that anesthetic or alcohol binding to a single subunit in the homomeric α1 GlyR is sufficient to alter channel functions.13</p>