Rational coupled dynamics network manipulation rescues disease-relevant mutant cystic fibrosis transmembrane conductance regulator

A novel approach identifying networks of residues involved in trans-protein dynamic coupling is applied to rescue mutant CFTR.


Introduction
Protein allostery is a ubiquitous mechanism central to the regulation of many cellular processes, such as enzyme catalysis and signal transduction. 1,2 An allosteric change involves the protein population ensemble redistributing itself among the available conformations. [2][3][4] NMR studies of dynamic coupling between residues in proteins support the idea that allostery is a common intrinsic property of many proteins. 3,5,6 Exerting control over such dynamic coupling interactions could open the door for novel therapeutic strategies in diseases such as cystic brosis. However, predicting and designing the effect of coupling interactions on protein structure and function remains a challenge. Determining interaction networks that couple distant sites of the protein and identifying hub ("hot spot") residues that control the dynamics of such coupling interactions are essential steps to uncovering the molecular mechanisms of protein allostery. 7 Previous efforts toward these ends have focused mainly on searching for new allosteric sites and mechanisms; X-ray crystallographic studies of bound and unbound structures have offered particularly important structural insights into allosteric regulation. 7 However, analysis of the static structures of allosteric end states cannot provide a complete picture of the inter-residue dynamics and interactions involved in allosteric conformational change or the coupling of dynamics at distal sites. NMR studies of protein dynamics have been pivotal in identifying "hidden" networks of residues with strong dynamic coupling. 6,8 Thermodynamic mutation cycles, 9 which measure the coupling between two mutation sites by their mutual contribution to protein stability, provide a direct method to systematically probe such relations between protein sites. However, due to experimental limitations and practical considerations, a large-scale study using these methods would be prohibitively laborious and time consuming.
Computational methods can be used to probe the coupling of amino acids and to identify networks of residues controlling protein conformational changes. For example, sequence-based approaches can reveal co-evolving residues likely to be energetically or functionally coupled. [10][11][12][13][14] However, the application of sequence-based approaches, which rely on evolutionary information, is limited by the availability of homologous sequences and complicated by the fact that evolutionary conservation is oen driven by factors other than function (e.g. stability, folding kinetics). 15,16 Recently, we have proposed that coupled dynamics plays a critical role in the pathophysiology of cystic brosis (CF). 17,18 Over 1500 mutations in the CFTR gene have been identied in patients with CF. The protein product of this gene, CFTR, plays a fundamental role in epithelial ion transport, providing a ratelimiting step in the regulation of salt secretion and reabsorption. In approximately 90% of CF patients, the deletion of a single phenylalanine (DF508) in its rst nucleotide-binding domain (NBD1) results in misfolding and misassembly of the protein (DF508-CFTR). 17 Cheng et al. determined that the DF508 mutation prevents maturation and trafficking of CFTR to the cell membrane, 19 although maturation and partial function of DF508-CFTR have been observed at sub-physiological temperature (<30 C) in various cell types. [20][21][22][23] Second site mutations, like the commonly used I539T substitution, can improve maturation, 24,25 but DF508-CFTR rescued in this manner exhibits poor function at physiological temperature, indicating that the fundamental defect in DF508-CFTR maturation and function is a consequence of reduced thermal stability. 26 This hypothesis is further supported by the fact that F508 plays a critical role in interfacing NBD1 with the fourth cytoplasmic loop of CFTR, thus contributing signicantly to the structural integrity and stability of the protein. 27,28 Therefore, effective therapeutic rescue of mutant CFTR requires restoration of thermal stability, oen mediated by protein uctuations and dynamics, making mutant CFTR an ideal system for development of methods for dynamic control.
Here, we develop a widely applicable computational methodology that utilizes concepts from graph theory to identify specic residues that propagate dynamic coupling effects between structurally distant sites in proteins. We apply our approach to the rescue of DF508-CFTR, using data from discrete molecular dynamics simulations, performable on a personal computer, to identify hot spot sites. We then utilize mutagenesis to demonstrate control over dynamic coupling between two distant regions of NBD1. We show that rational mutagenesis of the identied bottleneck sites dramatically rescues aberrant dynamics, as well as dramatically improving maturation and function of DF508-CFTR in mammalian cells. We further demonstrate the ability of our method to identify sites important to wild-type protein function and stability by showing that our designed mutation can also rescue the adjacent but distinct and more severe disease mutant, DI507-CFTR.

Results
Heightened thermal uctuations in the regulatory insertion propagate to the SDR and ATP-binding sub-domain We and others have reported that the NBD1 regulatory insertion (RI, residues 404-435) inuences the dynamics of NBD1 in CFTR, 18,26,29 but the nature and mechanism of this inuence are still unknown. To elucidate the effects of RI within NBD1, we perform discrete molecular dynamics (DMD) simulations 30,31 of wild type (WT), DF508-, and DI507-NBD1. We nd that thermal uctuations increase and/or are shied signicantly in several key regions of NBD1 upon deletion of F508 or I507: the RI, structurally diverse region (SDR), residues 532-552, 29 RI-SDR bridge (residues 492-502), F508-loop (residues 507-514), and part of the ATP-binding sub-domain (residues 570-600) (regions indicated by color-coded arrows, Fig. 1a). These heightened Fig. 1 Deletion of F508 or I507 results in altered coupling between regions of NBD1. (a) The root mean square fluctuations (RMSF) over the simulation of each residue in WT NBD1 compared with those in DF508and DI507-NBD1 suggests increased flexibility in several regions upon deletion mutation. (b) Difference map derived from the correlation maps of WT and DF508 NBD1. Blue denotes lost correlation upon deletion mutation, while red denotes gained correlation. Tube representation of NBD1 protein dynamics highlights changes in fluctuations in NBD1 upon deletion of F508. The thickness and color of the tube represent the extent of change in dynamics of the corresponding region during the simulation. Warmer colors indicate greater increase in flexibility, while colder colors indicate no change. (c) Same as (b), for DI507-NBD1.
uctuations suggest the possibility of dynamic coupling between the affected regions, which would transfer a perturbation in one region to any of the others. For example, deletion of F508 or I507 in the F508 loop affects the dynamics of RI, the RI-SDR bridge, the SDR, and the ATP-binding sub-domain.
F508 deletion results in non-native internal coupling of the a-subdomain and ATP-binding subdomain In order to conrm our hypothesis of coupled dynamics between key structural regions of NBD1, we perform covariation analysis of dynamic uctuations observed in simulations. 32 In DF508-NBD1, we nd a marked decoupling of the RI from the entire a-subdomain and parts of the ATP-binding subdomain as compared to the wild type ( Fig. 1b-i). In addition, we observe an increase in coupling between the various key regions of the a-subdomain and the ATP-binding subdomain ( Fig. 1b-ii), and an increase in internal coupling in the ATP-binding subdomain ( Fig. 1b-iii). From these ndings, we conclude that the F508 loop, the SDR, the RI-SDR bridge, the RI, and the ATP-binding subdomain are inter-connected in a network of coupled dynamics within NBD1. These changes suggest that the deletion of F508 from the F508 loop could cause dynamic instability and increased uctuations in this region, which lead to decoupling from the RI. The RI decouples from the entire a-subdomain, leading to the gain of non-native coupling within that region. Ultimately, the observed increase in uctuations and changes to dynamic coupling within NBD1 may result in domain instability and consequent misprocessing of DF508-CFTR. Deletion of the RI eliminates this transfer of dynamics through the domain, rescuing the protein as observed. 26 I507 deletion results in non-native coupling between the SDR and the ATP-binding sub-domain Instead of the complete decoupling of the RI from the a-subdomain seen with the F508 deletion (discussed above), DI507-NBD1 features a shi in the regions of RI that undergo coupling, losing coupling in one area only to gain it in adjacent residues ( Fig. 1c-i). We detect a similar shi of coupling interactions in the SDR, with its coupling to the RI ( Fig. 1c-i), and in addition note that the SDR loses coupling within itself ( Fig. 1cii), likely due to the increase in uctuations in this area (Fig. 1a). The perturbation to the SDR also includes a gained non-native coupling with part of the ATP-binding subdomain ( Fig. 1c-ii), which decouples from the a-subdomain (Fig. 1c-iii). Lastly, we nd a strong non-native coupling of the SDR with the F508 loop ( Fig. 1c-iv), likely the perturbation caused by the deletion mutation. Combining these ndings, we conclude that the deletion of I507 causes a perturbation to coupled dynamics that is transferred by coupled dynamics through the SDR to the RI and the ATP-binding subdomain, shiing coupling interactions all throughout NBD1, which severely affects the dynamic stability of the domain. Likewise, the non-native coupling of the ATP-binding subdomain to the SDR and the loss of internal coupling in the ATP-binding subdomain could affect the binding of ATP, which is stabilizing to NBD1 under normal conditions. In the absence of ATP-binding, NBD1 is signicantly destabilized, potentially compounding the effects of the DI507 mutation.

Network formalism of NBD1 dynamics reveals bottleneck residues that control dynamic coupling
In order to deduce the pathways through which these changes in dynamic coupling occur, for each system of wild type, DF508, and DI507 NBD1, we represent the pairwise correlation coefficients between residues as a complete graph, with the correlation between every possible pair of residues represented as a connection between nodes (an edge). We then weight this graph (G(N,E), where N is the set of nodes and E is the set of edges) such that each edge is enumerated by the correlation coefficient between the corresponding pair of residues. Weighting enables us to isolate those edges having the most signicant impact on the transduction of uctuations across the protein. To identify the specic residues mediating dynamic coupling within NBD1, we impose a correlation cutoff to our graph, eliminating edges with a weight below our cutoff until the largest component of the resulting disconnected sub-graph comprises approximately 50% of the total number of nodes in the graph (Fig. S1 †). This cutoff is at the critical threshold, where the network exhibits critical properties 33 and transitions from connected to disconnected. 34 Our rationale with this construction is to create an algorithm for network mapping that is without free parameters. From the resulting network, we determine whether a node is critical for the connectivity of the largest component in the sub-graph by monitoring the topological changes in the sub-graph as we iteratively remove each node from the sub-graph. If removing node a splits the largest subgraph into two or more components, then node a is a critical node, or "bottleneck," in the network. A bottleneck residue is a residue with the ability to affect the entire network because of its positioning as the lone connection between at least two regions of a network (i.e., an unavoidable step in the coupling of dynamics between one region and the other). Mutation or modication of a bottleneck residue can affect the dynamics, stability, and allosteric or binding behavior of a protein or domain. We nd several such bottleneck residues in the dynamic networks of DF508, DI507, and wild type NBD1, which we perturb to ameliorate the pathological dynamics of the mutant NBD1s.
Network bottlenecks control the transfer of heightened thermal uctuations from the RI of DF508-NBD1 We have previously reported that thermal uctuations are increased in the SDR and dynamic coupling is lost between the F508 loop and the ATP-binding sub-domain upon deletion of F508, both of which are ameliorated upon deletion of the RI. 26 We hypothesize that manipulating the coupling between the RI and the 508-loop and between RI and the SDR may ameliorate uctuations and improve dynamic stability in these regions. In order to identify bottleneck residues that mediate the transfer of uctuations in the RI across the domain, we conduct a network analysis of the correlated dynamics between residues in NBD1 (Methods). We conclude from this analysis that the residues K464, T465, L468, S492, I601, and V603 are critical nodes involved in the transfer of the thermal uctuations from the RI to other regions of NBD1 via dynamic coupling (Fig. 2a). K464 and V603 mediate coupling with the b-strands forming the ATP-binding core sub-domain (Fig. S2, † red residues), L468 and L601 with the regulatory extension (RE) of NBD1 (Fig. S2, † blue residues), and T465 with both of these regions ( Fig. 2b and S2 †). However, S492 in b5 is the only residue mediating dynamic coupling between the RI and the 508-loop and SDR ( Fig. 2b: yellow nodes; Fig. S2: † green residues). In addition, we note that S492 has been identied as the site of a rescue mutation for CFTR in an independent study. 18 We conclude that S492 is a network bottleneck that dynamically couples the RI with these regions in NBD1.
DI507-NBD1 network bottlenecks reside in regions essential for CFTR maturation and function Unlike DF508-CFTR, removal of the RI in NBD1 does not rescue DI507-CFTR, so we instead adopt the strategy of ameliorating the non-native coupling of the SDR and ATP-binding subdomain. Using network analysis of correlated dynamics between DI507-NBD1 residues (Methods), we nd that residues T547, L548, K564, D565, A566, D567, and L568 are critical nodes in the cross-domain coupling of dynamics in DI507-NBD1, and also that each of these critical nodes is central to coupling between the SDR and the ATP-binding sub-domain ( Fig. 3a and  b). However, residues T547 and L548 reside directly adjacent to or within the signature motif responsible for ATP binding, 35 and residues D565, A566, and D567 comprise the di-acidic exit code recognized by coat complex II (COPII) for transport of CFTR from the endoplasmic reticulum to the membrane, 36,37 with K564 and L568 directly adjacent (Fig. 3c). Although several reversion mutations exist in the signature motif, 24,38,39 T547 and L548 appear to be in position to participate in the binding of ATP (Fig. 3c), one of the major contributing factors to CFTR stability. 40,41 Using computational methods, we nd potential stabilizing mutations at some of these sites (although only marginally stabilizing, DDG $1-2 kcal mol À1 ), but experiments conrm that these mutations result in lack of mature CFTR ( Fig. S3 †), reecting the crucial nature of these conserved residues to full-length CFTR in the cell. While we therefore cannot utilize mutagenesis at these positions to rescue DI507-CFTR, we demonstrate that these important sites are hubs of coupled dynamics networks capable of affecting domain dynamics.

S492P substitution restores dynamic stability and function to DF508-CFTR
To conrm that DF508-NBD1 dynamics can be improved by restoring coupling between the F508 loop and the ATP-binding sub-domain, we perform computational mutagenesis of S492 in the DF508 background. In an alignment of NBD1 with NBD2 and nucleotide-binding domains from many other ABC proteins, we nd that the position corresponding to S492 contains exclusively either serine or proline (Fig. S4 †). 18 We therefore test the effect of a proline substitution at this position in NBD1. We perform DMD simulations and RMSF analysis of the NBD1 variant DF508-S492P in the context of I539T, a previously identied maturation reversion mutation 24,25 (DF508-PT, generated using Medusa 15,42,43 ). We nd that DF508-PT NBD1 exhibits a thermal uctuation prole comparable to that of wild type NBD1, suggesting potential thermal rescue of DF508-CFTR (Fig. 4a, in comparison to Fig. 1a). We demonstrate that this rescue of dynamics is due to the S492P substitution, as I539T mutation alone does not return wild type dynamics to DF508-NBD1 (Fig. S5 †). Notably, the network of correlated dynamics in DF508-PT does not contain any bottlenecks; no single residue is so crucial to the dynamics that without it the network splits into non-communicating pieces (Fig. 4b). Such radical topological rearrangement could be a consequence of decreased uctuations in the RI due to strengthened dynamic coupling. Alternately, in DF508-PT NBD1 the uctuations disconnected sub-graph after disconnecting the corresponding residue (x-axis). We consider residues that partition the sub-graph into connected components of ten or more residues as critical nodes. Different colors represent different components. (b) Network representation of the disconnected subgraph. S492 forms a critical node that partitions the sub-graph into two connected components of nearly equal size (yellow and gray circles). Yellow nodes with black outline represent residues in the RI. Smaller components are represented as white circles.
originating from the RI might dissipate throughout the domain such that they are not channeled towards the 508-loop or the SDR.
To determine whether DF508-PT CFTR is functionally rescued at physiologically relevant temperatures, we perform experimental characterization of the DF508-PT variant of CFTR. As predicted from our simulations, Western blot analysis of lysates from HEK293 cells transiently expressing wild type and mutant (DF508, DF508-PT) CFTR reveals that maturation of DF508-PT is improved to a level similar to that of wild type CFTR (Fig. 4c). Using single channel recordings, we found in a previously published study 18 that DF508-PT CFTR is functional at 35 C while DF508-I539T is not, suggesting that the S492P substitution not only improves maturation but also restores function to DF508-CFTR (Fig. 4d). Furthermore, this result indicates that it is the dynamic stability of the protein, not evasion of the cell's quality control machinery such as attained by the I539T substitution, that guarantees CFTR function. These single channel tracings reveal that the fast, ickering gating mode characteristic of DF508-CFTR is interrupted only rarely by normal full conductance gating transitions in the I539T variant, whereas normal gating transitions are the dominant mode in DF508-PT CFTR. When assessed near human physiological temperature, DF508-PT features a level of channel activity that provides conductance similar to wild type CFTR (Fig. 4d). In parallel, we monitor iodide efflux from BHK cells stably expressing wild type, DF508, or DF508-PT CFTR (Fig. 4e). Together, the relative single channel properties of DF508-PT and wild type CFTR (Fig. 4d) and iodide efflux behavior (Fig. 4e) provide substantial evidence that DF508-PT CFTR trafficked to the plasma membrane is functional, and that the S492P mutation identied with our graph theoretical approach successfully rescues DF508-CFTR.
S492P substitution restores partial wild-type function in DI507-CFTR Because of the similarities in increased uctuations and nonnative coupling of the a-subdomain between DI507-NBD1 and DF508-NBD1, we hypothesize that the S492P mutation may have a similar rescue effect in DI507-NBD1. In simulations of DI507-PT NBD1, we nd that uctuations are ameliorated in comparison with DI507-NBD1, returning to a wild type-like prole (Fig. 4a). Notably, this dynamic stabilization results in wild type-like maturation (Fig. 4c). These results suggest that stiffening the RI-SDR bridge with a proline mutation at position 492 prevents the transfer of correlated dynamics from the SDR to both the RI and the ATP-binding subdomain, alleviating the changes in coupled motions and rescuing the mutant.

Discussion
Allostery in proteins is viewed as an effect that propagates tens of angstroms across the structure, an action caused by a perturbation at a distant site. Described 50 years ago, this phenomenon is ubiquitously observed in various biological processes, from gene regulation to protein biogenesis. However, the mechanism of coupled dynamics transfer (caused by ligand binding, for example) over large distances (>10Å, spanning multiple bond lengths) remains unknown, and determining disconnected sub-graph after disconnecting the corresponding residue (x-axis). We consider residues that partition the sub-graph into connected components of ten or more residues as critical nodes. Different colors represent different components. (b) Network representation of the disconnected sub-graph. Several residues form critical nodes that partition the sub-graph into two connected components of nearly equal size (yellow and blue circles; grey circles represent nodes that switch components depending on which critical node is removed). (c) Structure of NBD1 with critical nodes highlighted in red stick representation. Color-coded by sub-structure, the RI is pink, ATP-binding subdomain is green, SDR is blue, and F508 loop is orange. I539 (blue) and S492 (grey) are highlighted in stick representation. this mechanism is a fundamental pursuit of structural biology. 4 We utilize graph theoretical methods to determine the coupling between regions that mediates dynamics transfer and apply our methodology to rescue cystic brosis-causing mutants of the protein CFTR. We present two major advances: (i) development of a generally-applicable method to determine inter-residue interactions mediating trans-domain coupled dynamics in proteins (Fig. S6 †) and (ii) rescue of CF-causing mutant forms of the CFTR chloride ion channel. Using our methodology, we not only detect the networks of inter-residue interactions that mediate dynamic coupling of distal protein regions, but also redesign these networks, demonstrating the ability to modulate protein dynamics to achieve specic goals. The applicability of our methodology is not limited to CFTR, but can be extended to any macromolecule of interest. Since the method involves representation of a protein using a generic mathematical construct, a graph, the concept is applicable to all macromolecules that are known to feature coupled dynamics, including RNA.
During the preparation of this work, Roy and Post 44 have employed a similar method to examine the long-distance effects of drug binding in viral capsids. However, that work studies only local networks, dening edges only between those residues featuring separation less than 14Å, regardless of the correlation between those residues. Our method of including all pairwise edges and culling edges featuring a correlation less than a cutoff or DF508-CFTR (solid triangles). Stimulation cocktail is added at time 0 to activate iodide efflux through CFTR channels. The values represent the mean AE standard deviation of the amount of iodide released from the cells during a one minute interval (n ¼ 3). Efflux buffer containing 0.1% NP40 is added at the end of each assay (arrow) to release remaining iodide. We find that DF508-PT CFTR partially recovers WT activity as compared to DF508-CFTR.
value gives the dual advantage of increased network coverage with the inclusion of distal correlations and faster computation of shortest path, since weak connections are eliminated and are not included in calculations of shortest path.
In the present study, we apply our computational methodology to study dynamic coupling within NBD1 of CFTR. We identify hubs in NBD1 that mediate dynamic coupling between various sites across the protein (Fig. S2 †), and nd that the residue S492 alone mediates dynamic coupling between the dynamic RI region and the F508 loop where disease-relevant deletion occurs, making this residue a hot spot for rational design via mutagenesis. We further nd that mutation of S492 to proline (an evolutionarily-conserved amino acid in this position among NBDs) successfully rescues DF508and DI507-CFTR. To evaluate whether the rescue of CFTR upon S492P mutation is due simply to thermodynamic stabilization, we calculate the DDG of mutation at this position using Eris (Methods). 42,43 Intriguingly, we nd that the thermodynamic stability of DF508-NBD1 is not improved by proline substitution, nor by any amino acid substitution at this position, emphasizing the role of dynamic coupling within NBD1 in CFTR structural stability (Fig. S7 †). The conservation of S492 in CFTR among various species further indicates an evolutionary basis for an dynamic coupling mechanism of regulation of NBD1 dynamics (Fig. S4 †). Furthermore, we nd both computationally and experimentally that mutagenesis of S492 in the context of I539T (referred to as DXXXX-PT), a previously iden-tied maturation reverting mutation, 24 signicantly improves the dynamic stability, maturation, and function of CFTR deletion mutants. The biological signicance of the combined I539T and S492P substitutions is also apparent in nature. For instance, some non-mammalian (e.g. Xenopus laevis and Squalus acanthias) CFTR orthologs naturally possess T539 and P492, and the DF508 mutation does not compromise their folding and maturation ( Fig. S4 and S8 †). Notably, mammalian species (e.g. Equus caballus and Ovis aries) do not contain these variations, and their folding and maturation is affected by F508 deletion (Fig. S4 and S8 †).
In order to demonstrate the importance of dynamic behavior in the trans-domain "communication" between various regions, for contrast we perform a similar analysis of the static structure of human NBD1, using the number of contacts made between each pair of residues as the weight of each edge in the graph. In the graph of the static system, S492 does not appear in any of the paths between the RI and the F508 loop, further supporting the idea that an understanding of protein dynamics is critical to rational manipulation of proteins.
In the case of CF, the main cause of disease (approximately 90% of cases) is the deletion of F508 from NBD1 of CFTR. 17 The deletion of I507, the adjacent residue, also results in disease. 45 We nd that deletion of these residues from the F508 loop substructure leads to non-native dynamic coupling interactions within NBD1. In DF508-NBD1, the deletion mutation causes a perturbation that decouples the RI from the rest of NBD1, and the loss of coupling with the RI results in non-native coupling within the a-subdomain and the ATP-binding domain. In DI507-CFTR, the deletion of I507 changes the interactions of the mutation site with the SDR, from which perturbations affect coupling throughout the domain. This difference explains why the removal of RI rescues DF508-CFTR, but does not affect DI507-CFTR (Fig. S9 †). However, in both CFTR constructs, the perturbed SDR plays a major role in mediating transfer of aberrant dynamics from the mutation site. The substitution of a proline for a serine at position 492 in NBD1 stiffens the protein backbone and inhibits uctuations, inhibiting the coupling between the SDR and both the RI and the ATP-binding subdomain, interrupting the transfer of non-native dynamics and rescuing both deletion mutants. Notably, the phenomenon of lost native coupling in the ATP-binding sub-domain of DI507-NBD1, but not in DF508-NBD1, can potentially explain why the latter, but not the former, mutant can be rescued by decreasing the temperature from 37 C to 27 C (Fig. S9 †); because the binding of ATP to the NBDs is known to be a signicant stabilizing factor in CFTR, 40,41 the loss of coupling between the different areas of the ATP-binding subdomain could interfere with the binding of this important ligand.
Given these ndings and conclusions, we are faced with the question: if the S492P mutation rescues both DF508and DI507-CFTR, why do we not identify S492 in our computational analysis of the DI507-NBD1 dynamic network? We note that coupled motions are markedly lower in DI507-NBD1 as opposed to DF508, as evidenced by the much lower correlation coefficient cutoff in our networks (0.50 for DF508, 0.32 for DI507). We hypothesize that the low level of correlation in the DI507-NBD1 domain and sharp transition in node inclusion (Fig. S1 †) make the choice of network cutoff much more critical to network properties. We nd that S492 is indeed located in the bottleneck region, in fact correlated (correlation coefficient 0.30, as compared to the cutoff 0.32) with one of the identied bottleneck residues, L568. We note that in using our method, an expert researcher could manually ne-tune the correlation coefficient cutoff based on prior knowledge or on regions or residues of interest. Manual intervention may be particularly useful in systems such as DI507-NBD1, where the sharp transition in network composition with correlation coefficient cutoff makes choice of this parameter particularly critical.
The sensitivity of the DI507-NBD1 system leads us to emphasize that our method serves as a guide to narrow the wide eld of possible rescue mutations, highlighting residues with potential inuence over protein dynamics, and that not all identied residues will be the site of viable rescue mutations, nor that all rescue mutation sites will be identied. The computational methods that we employ necessarily are molecular in scale, and cannot integrate all cellular processes that inuence the expression, maturation, and interactions of particular mutants. Also due to molecular detail, our computational method is only applicable in cases where the structure of the protein or protein domain of interest has been solved to high resolution (<4Å). Lower resolution structures would require additional computational methods to model atomic resolution before our approach could be applied. Finally, the structure must be of a size that is practical for extensive dynamic sampling during molecular simulations.

NBD1 models and discrete molecular dynamics simulations
Crystal structures exist of wild type NBD1 (PDB ID: 2BBO) and DF508-NBD1 (PDB ID: 2BBT), but because of the inherent exibility of the RI this region is not resolved in either structure. We reconstruct the missing residues in the RI ab initio using discrete molecular dynamics (DMD) simulations 30 with the Medusa force eld. 15 We perform multiple iterations of replica exchange simulation and energetic minimization using Chiron 46 to obtain nal minimized models of full-length wild type and DF508-NBD1. No structure of DI507-NBD1 is currently available, therefore starting from our minimized model of fulllength wild type NBD1, we manually delete the I507 residue and reseal the protein backbone using DMD simulations. We reseal the protein backbone by imposing peptide-bonding constraints on the cleaved ends until the backbone is once again intact. We hold the majority of the protein static during this process, allowing only three residues on either side of the deletion to move freely and reseal the gap caused by deletion. From these three initial full-length models (wild type, DF508, and DI507), we utilize the Eris suite 42,43 to create S492P/I539T mutants of each construct. We perform energy minimization of all constructs using Chiron.
We perform long timescale single-temperature DMD simulations using the nal minimized models of full-length NBD1. We perform 10 randomized simulations for each construct at a temperature of 0.4 kcal mol À1 k B À1 , with each individual simulation having a length of 10 6 DMD time steps (approximately 50 ns), summing to a total of approximately 500 ns for each construct. The results that we present here are averaged over all simulations for each construct, unless otherwise specied.

Construction of dynamic networks and determination of optimal paths
We utilize the resulting simulation trajectories to determine dynamic coupling between the various regions of NBD1 by computing correlation coefficients of motion between all residues. 32 We represent the pairwise correlation map as a complete weighted graph G (N,E), with the C a atoms in NBD1 representing the nodes (N) of the graph and the edges (E) being the connections between these nodes. The edge weight E ij between any two nodes i and j is the correlation coefficient C ij between the corresponding pair of residues.
To determine optimal paths of dynamics transfer through coupled residues, we use the same network with edge weights of E ij ¼ 1 À |C ij | and apply Dijkstra's algorithm. 47 Since dynamic coupling is mediated by physical interactions, we reduce the complete graph by removing edges between nodes representing residues that are not in contact. We consider two residues to be in contact if their C b atoms (C a for Glycine) are within 7.5Å of one another. Our aim is to determine the transfer of correlated motion across the protein. Hence, we emphasize the effect of local interactions by retaining only those edges in which the participating nodes feature a contact frequency (u c ) of at least 0.5 over the simulation.

Flexible backbone redesign
We perform iterative exible backbone redesign and structural relaxation using the Medusa suite 15 in order to ensure optimum backbone conguration for computational mutagenesis of human NBD. We perform 15 randomly-seeded iterations of backbone redesign Monte Carlo simulations with 20 Monte Carlo iterations per replicate. We select the lowest energy structure from each iteration to use as the starting point for the next iteration. Finally, we choose the structure with lowest energy and smallest Kabsch root mean square deviation (KRMSD) from the initial structure for xed-backbone mutation analysis using the Eris suite. 42,43 CFTR construction and expression We express human CFTR cDNAs encoding wild type or mutant proteins transiently in HEK 293 cells or stably in BHK-21 cells, with pcDNA3 or pNUT vectors, respectively, as previously described. 28 We use the QuickExchange protocol (Stratagene) to generate mutant CFTR constructs in pcDNA3 and pNUT vectors from human WT CFTR cDNA, and conrm sequences by automated DNA sequencing (UNC-CH Genome Analysis Facility). We carry out transfection using jetPEI transfection reagent (Fermentas, Glen Burnie, MD) according to the manufacturer's instructions. For stable cell line establishment, we select and maintain BHK cells expressing CFTR in methotrexate-containing media as previously described. 48

Membrane isolation
We harvest BHK or HEK 293 cells expressing CFTR or its variants by scraping, and then homogenize the cells on ice in 10 mM Hepes, pH 7.2, 1 mM EDTA (ethylenediaminetetraacetic acid) containing a protease inhibitor cocktail (benzamidine at 120 mg ml À1 , E64 at 3.5 mg ml À1 , aprotinin at 2 mg ml À1 , leupeptin at 1 mg ml À1 and Pefablock at 50 mg ml À1 ). We centrifuge the resulting samples at 600 g for 15 minutes to remove nuclei and undisrupted cells, followed by centrifugation at 100 000g for 60 minutes to pellet the membranes, which we then resuspend in phosphorylation buffer (10 mM Hepes, pH 7.2, containing 0.5 mM EGTA (ethylene glycol bis(b-aminoethyl ether), N,N 0 -tetraacetic acid), 2 mM MgCl 2 , and 250 mM sucrose). We utilize brief (3 Â 20 s) bath sonications to generate vesicles of uniform size. For single-channel recordings, we phosphorylate membrane vesicles by incubating with 50 nM PKA catalytic subunit (Promega) and 2 mM Na 2 ATP (Sigma) in phosphorylation buffer for 20 min at +4 C. We then aliquot the membranes and store at À80 C for later use.

Single-channel measurements
To prepare planar lipid bilayers, we drill a 0.2 mm hole in a Teon cup and paint the hole with a phospholipid solution containing a 3 : 1 mixture of 1-palmitoyl-2-oleoyl-sn-glycero-3phosphoethanolamine and 1-palmitoyl-2-oleoyl-sn-glycero-3phosphoserine (Avanti Polar Lipids) in n-decane. The lipid bilayer separates 1 ml of solution in the Teon cup (cis side) from 5 ml of solution in the outer glass chamber (trans side). Both chambers are magnetically stirred and thermally insulated. We utilize a temperature control system (TC2BIP, Cell Micro Controls).
We transfer CFTR ion channels into the pre-formed lipid bilayer by spontaneous fusion of membrane vesicles containing the CFTR variants. To maintain uniform orientation and functional activity of the CFTR channels transferred into the bilayer, we add 2 mM ATP, 50 nM PKA, and membrane vesicles into the cis compartment only. We perform all measurements in symmetrical salt solution (300 mM Tris-HCl, pH 7.2, 3 mM MgCl 2 and 1 mM EGTA) under voltage-clamp conditions using an Axopatch 200B amplier. We maintain a membrane voltage potential of À75 mV, the difference between cis and trans (ground) compartments. We analyze the resulting data as previously described. 50 Iodide efflux assay We grow BHK cells stably expressing wild type and mutant CFTR to $100% conuence in six-well plates and incubate in an iodide loading buffer (136 mM NaI, 3 mM KNO 3 , 2 mM Ca(NO 3 ) 2 , 11 mM glucose, and 20 mM Hepes, pH 7.4) for one hour at room temperature. We rinse the cells with iodide-free efflux buffer (which is the same as the loading buffer except that NaI is replaced by NaNO 3 ) to remove extracellular iodide. We collect samples by completely replacing the efflux buffer (1 ml volume) with fresh solution at one-minute intervals. We use results from the rst four samples to establish a baseline. We measure iodide efflux upon stimulation with PKA agonists (10 mM forskolin, 100 mM dibutyl-cAMP and 1 mM 3-isobutyl-1methylxanthine) using an iodide-selective electrode LIS-1461CM (Lazar Res. Lab., Inc.), as previously described. 48

Conclusions
We demonstrate here that NBD1 of CFTR, for which transdomain coupled dynamics was not previously described, features patterns of correlated motions that form a network throughout the domain and allow structural uctuations to be transferred to distal sites through dynamic coupling, supporting the general notion that all folded proteins are to some extent allosteric in nature. 3 Moreover, we show that these networks can be rationally redesigned to affect a desired outcome. Identifying the residues involved in trans-domain dynamics can provide a map to greatly improve the stability and/or modulate the function of enzymes for numerous biotechnological applications, as well as contribute to our understanding of the many human diseases caused by protein dysfunction. Determining dynamic coupling networks and hot spots is also of practical importance in rational drug design. 3 Conventionally, structurebased virtual screening utilizes experimentally-validated binding sites, such as an enzyme catalytic site or a ligand-bound pocket, but the determination of pre-existing dynamic communication pathways has great potential for identifying novel "druggable" sites that can be targeted to modulate protein function in human diseases. Thus, beyond the signicant advance in understanding the fundamental defect in cystic brosis, our approach has broad application in the elucidation of mechanisms of protein function and dysfunction in disease.