Skip to main content
Advertisement
  • Loading metrics

A lipophilicity-based energy function for membrane-protein modelling and design

  • Jonathan Yaacov Weinstein,

    Roles Conceptualization, Formal analysis, Investigation, Methodology, Software, Validation, Visualization, Writing – original draft, Writing – review & editing

    Affiliation Department of Biomolecular Sciences, Weizmann Institute of Science, Rehovot, Israel

  • Assaf Elazar,

    Roles Conceptualization

    Affiliation Department of Biomolecular Sciences, Weizmann Institute of Science, Rehovot, Israel

  • Sarel Jacob Fleishman

    Roles Conceptualization, Funding acquisition, Methodology, Project administration, Resources, Software, Supervision, Validation, Writing – original draft, Writing – review & editing

    sarel@weizmann.ac.il

    Affiliation Department of Biomolecular Sciences, Weizmann Institute of Science, Rehovot, Israel

Abstract

Membrane-protein design is an exciting and increasingly successful research area which has led to landmarks including the design of stable and accurate membrane-integral proteins based on coiled-coil motifs. Design of topologically more complex proteins, such as most receptors, channels, and transporters, however, demands an energy function that balances contributions from intra-protein contacts and protein-membrane interactions. Recent advances in water-soluble all-atom energy functions have increased the accuracy in structure-prediction benchmarks. The plasma membrane, however, imposes different physical constraints on protein solvation. To understand these constraints, we recently developed a high-throughput experimental screen, called dsTβL, and inferred apparent insertion energies for each amino acid at dozens of positions across the bacterial plasma membrane. Here, we express these profiles as lipophilicity energy terms in Rosetta and demonstrate that the new energy function outperforms previous ones in modelling and design benchmarks. Rosetta ab initio simulations starting from an extended chain recapitulate two-thirds of the experimentally determined structures of membrane-spanning homo-oligomers with <2.5Å root-mean-square deviation within the top-predicted five models (available online: http://tmhop.weizmann.ac.il). Furthermore, in two sequence-design benchmarks, the energy function improves discrimination of stabilizing point mutations and recapitulates natural membrane-protein sequences of known structure, thereby recommending this new energy function for membrane-protein modelling and design.

Author summary

Membrane proteins comprise a third of the genome and have essential roles as intermediaries between the cell and its environment. Despite exciting recent progress in membrane-protein modelling and design, however, these fields lag far behind advances in soluble proteins, chiefly because of inaccurate modelling of the membrane environment. Recently, our lab developed an assay, called dsTβL, that used high-throughput experimental screening to infer the energetics of each amino acid across the bacterial plasma membrane. Here, we encode the dsTβL energetics in the Rosetta software suite for biomolecular modelling and design and subject the energy function to three structure prediction and design benchmarks. The new energy function consistently outperforms the previous Rosetta membrane energy function. Additionally, ab initio structure prediction of homooligomeric membrane proteins results in accurate predictions in ⅔ of the examples in our benchmark. Therefore, we present a web server, called TMHOP, to compute the structures of single-pass homooligomeric membrane proteins directly from sequence. The results suggest that the automated design of large and complex membrane proteins is within reach.

Introduction

Membrane proteins have essential biological roles as receptors, channels, and transporters. Over the past decade, significant progress has been made in the design of membrane proteins, including the first design of membrane-integral inhibitors [1], a transporter [2,3], and a de novo designed structure based on coiled-coil motifs [4]. Despite this exciting progress, however, modelling, design, and engineering of membrane proteins lag far behind those of soluble proteins [5]. This lag is due, in part, to the relatively small number of high-resolution membrane-protein structures [6] and is exacerbated by these proteins’ typically large size. Clearly, however, the heterogeneity of a membrane protein’s environment, comprising water, lipid, and polar headgroups, is the most significant complication [7]. Therefore, modelling solvation is a fundamental problem that impacts all membrane-protein structure prediction and design.

Current energy functions used in modelling and design incorporate simplified solvation models [8]. For instance, RosettaMP uses information inferred from water-to-cyclohexane partitioning [9] as a proxy for amino acid solvation in the plasma membrane [10,11]. Due to these simplifications, expert analysis has been a prerequisite for accurate membrane-protein modelling and design [12,13]. Automating modelling and design processes and extending them to complex membrane proteins will likely require an accurate energy function that correctly balances intra-protein interactions, membrane solvation and water solvation [1416].

To understand the contributions to membrane-protein solvation, we recently established a high-throughput experimental screen, called deep sequencing TOXCAT-β-lactamase (dsTβL), which quantified apparent amino acid transfer energies from the cytosol to the E. coli plasma membrane [17]. From the resulting data, we inferred apparent position-specific insertion profiles for each amino acid relative to alanine, reconciling previously conflicting lines of evidence [18]. Foremost, the lipophilicity inferred for hydrophobic residues, such as Leu, Ile, and Phe, was greater than previously measured in some membrane mimics, including the water-to-cyclohexane transfer energies that are the basis for membrane solvation in Rosetta [9,11,16] (approximately 2 kcal/mol according to dsTβL compared to ½ kcal/mol), and in line with theoretical considerations [19,20]. Second, the profiles exhibited a strong 2 kcal/mol preference for Arg and Lys in the intracellular side of the plasma membrane compared to the extracellular side. While this preference, known as the “positive-inside” rule, was revealed based on sequence analysis 30 years ago [2123], the dsTβL assay was the first to indicate a large energy gap favouring positively charged residues in the intracellular relative to the extracellular membrane leaflet. The accuracy and generality of the dsTβL apparent transfer energies were partly verified by demonstrating that they correctly predicted the locations and orientations of membrane spans directly from sequence even in several large and complex eukaryotic transporters [24]. Taken together, these results provided reassurance that the dsTβL apparent insertion energies correctly balanced essential aspects of membrane-protein solvation.

As the next step towards accurate all-atom membrane-protein modelling and design, we develop a new lipophilicity-based energy term based on the dsTβL amino acid-specific insertion profiles and integrate this energy term in the Rosetta centroid-level and all-atom potentials. We furthermore develop a strategy to enhance conformational sampling of membrane-spanning helical segments and of helix-tilt angles observed in naturally occurring membrane proteins. Encouragingly, the new energy function outperforms previous ones in three benchmarks essential to modelling and design: atomistic ab initio structure prediction starting from completely extended chains of single-spanning membrane homo-oligomers of known structure, prediction of mutational effects on stability, and sequence recovery in combinatorial sequence design. An automated web server for structure prediction in transmembrane homooligomeric proteins (TMHOP) is available at (http://tmhop.weizmann.ac.il) and may enable modelling of the membrane-spanning domains of receptors. We conclude that the combination of lipophilicity and energetics developed for soluble proteins provides a basis for accurate structure prediction and design of membrane proteins.

Results

A lipophilicity-based membrane-protein energy function

The recent all-atom energy function in Rosetta, ref2015, is dominated by physics-based terms, including van der Waals packing, hydrogen bonding, electrostatics and water solvation [25]. This energy function was parameterized on a large set of crystallographic structures and experimental data of water-soluble proteins and was shown to outperform previous energy functions in several structure-prediction benchmarks. For membrane-protein modelling and design, however, the ref2015 solvation potential is relevant only to the water-embedded regions of the protein; a different potential is required to model the energetics of amino acids near and within different regions of the plasma membrane.

Accordingly, we sought to replace the ref2015 solvation model with one that encodes a gradual transition from the default water-solvation that evaluates regions distant from the plasma membrane and the dsTβL insertion profiles near and within the plasma membrane. The dsTβL profiles were inferred from an experimental mutation analysis of a monomeric membrane span into which each of the 20 amino acids were individually introduced at each position [17]; the profiles were then normalized to express the apparent transfer energy for each amino acid at each position relative to a theoretical poly-Ala membrane span, yielding apparent ΔΔGAla—>mut at each position across the plasma membrane (Fig 1). As a first step to encoding these energy profiles in Rosetta, we smoothed these profiles and symmetrised them with respect to the presumed membrane midplane, except the profiles for Arg, His, and Lys, for which the “positive-inside” rule applies (S1 Fig).

thumbnail
Fig 1. The lipophilicity-based ref2015_memb energy function.

Membrane-insertion profiles for six representative amino acids are shown. Raw dsTβL data (purple dots), ref2015 (dashed green line), the ref2015_memb potential (dashed blue line) and the dsTβL profiles (red line). Negative and positive membrane depths indicate the inner and outer membrane leaflets, respectively; the presumed membrane midplane is at 0.

https://doi.org/10.1371/journal.pcbi.1007318.g001

Next, we implemented an iterative strategy to encode the dsTβL energetics in a modified ref2015 all-atom energy function which we called ref2015_memb. To enable efficient conformational search as required in ab initio structure prediction and de novo design, we also encoded this energetics in the centroid-level energy function [26]. As a reference state in both all-atom and centroid-level modelling, we generated an ideal poly-Ala α helix and placed it perpendicular to the membrane plane. At each position along the helix (including the aqueous and membrane phases), we introduced each of the 19 point mutations, relaxed the models using the all-atom or centroid-level energy function, and computed the energy difference due to each single-point mutation ΔΔGAla—>mut. In the first iteration of these calculations, the unmodified ref2015 or centroid-level energy functions were used, resulting, as expected, in large deviations from the apparent energies observed in the dsTβL profiles (red lines in Fig 1). We then added a new context-dependent 1-body energy term, called MPResidueLipophilicity, which encoded the difference between the computed and dsTβL energies for each mutation at each position, ΔΔΔGAla—>mut. We iterated mutation, relaxation, energy calculations, and MPResidueLipophilicity updates for each of the mutations at each position up to ten times, noting that the computed energies converged with the trends observed in the experiment (blue and red lines in Fig 1, respectively). Scripts for calibrating the all-atom and centroid energy functions are available at github.com/Fleishman-Lab/membrane_protein_energy_function to enable adapting future improvements of the Rosetta energy functions to encode the dsTβL energetics.

The dsTβL apparent energy profiles were inferred from a monomeric segment [17]. Consequently, the profiles express the lipophilicity of each amino acid relative to Ala across the membrane when that amino acid is maximally solvent-exposed. To account for amino acid burial in multispan or oligomeric membrane proteins, we derived a continuous, differentiable and easily computable weighting term that expresses the extent of a residue’s burial in other protein segments. For any amino acid, this weighting term is based on the number of heavy-atom neighbours within 6 and 12 Å distance of the amino acid’s Cβ atom (Eq 1) resulting in a weight that expresses the extent to which a residue is buried in other protein segments or exposed to solvent (0 to 1, respectively). Water-embedded and completely buried positions are treated with the ref2015 solvation energy; fully membrane-exposed positions are treated with the MPResidueLipophilicity energy, and positions of intermediate exposure are treated with a linearly weighted sum of the two terms.

In summary, the actual contribution from solvation of an amino acid is a function of its exposure to the membrane and depends on the amino acid’s lipophilicity according to the dsTβL apparent energy and the position’s location relative to the membrane midplane. Note that this energy term averages lipophilicity contributions in the plasma membrane and does not express atomic contributions to solvation that are likely to be important in calculating membrane-protein energetics in different types of biological membranes [11,27], in non-helical membrane-exposed segments, or surrounding water-filled cavities [28].

The dsTβL assay reports on residue-specific insertion into the plasma membrane. Ab initio modelling and de novo design, however, also require a potential that addresses the protein backbone solvation. Although the low-dielectric environment in the core of the membrane enforces a strong tendency for forming canonical α helices [7], deviations from canonical α helicity can make important contributions to membrane-protein structure and function [29]. Therefore, we encoded an energy term, called MPHelicality, that allows sampling backbone dihedral angles and penalises deviations from α helicity (Eq 5). MPHelicality enforces strong constraints on the dihedral angles in the lipid-exposed surfaces at the core of the membrane and is attenuated in regions that are buried in other protein segments and in the extra-membrane environment (using the same weighting as for lipophilicity, Eq 1); this term thus allows significant deviations from α helicity only in buried or water-embedded regions.

In preliminary ab initio calculations starting from a fully extended chain, we noticed that conformational sampling significantly favoured large helical tilt angles relative to the membrane normal (θ in Fig 2). By contrast, 50% of naturally observed membrane spans exhibit small tilt angles in the range 15–30°. The skew in conformational sampling towards large tilt angles is expected from previous theoretical investigations according to which the distribution of helix-tilt angles in random sampling is proportional to sin(θ), substantially preferring large angles compared to the distribution observed in natural membrane proteins [30]. To eliminate this skew in conformational sampling, we introduced another energy term, called MPSpanAngle (Eq 4 and Fig 2), that strongly penalised large tilt angles, guiding ab initio sampling to tilt angles observed in natural proteins.

thumbnail
Fig 2. Observed versus expected tilt angles in membrane-spanning helices relative to the membrane normal.

The distribution of helix tilt angles (θ in the inset sphere) in natural membrane proteins shows a strong preference for small angles (red bars, left), whereas the distribution resulting from random conformational sampling is proportional to sin(θ) (blue bars) [30] significantly overrepresenting large tilt angles. The MPSpanAngle energy term (green line; Eq 4) penalises large tilt angles and focuses ab initio conformational sampling on tilt angles observed in membrane-protein structures. inset The expected distribution of helix-tilt angles is proportional to the circumference of a circle plotted by that helix around an axis perpendicular to the membrane-normal (panel adapted from ref. [30]). The membrane plane is depicted as a grey circle.

https://doi.org/10.1371/journal.pcbi.1007318.g002

In summary, ref2015_memb encodes three new energy terms relative to the soluble energy function ref2015: (1) a lipophilicity term based on amino acid type, membrane-depth, and burial; (2) a penalty on deviations from α helicity in backbone-dihedral angles; and (3) a penalty on the sampling of large tilt angles with respect to the membrane-normal (S1 Table). In the calculations reported below, the penalties on deviations from α helicity and helix-tilt angles are implemented in all centroid-level ab initio structure prediction simulations; all-atom calculations use the ref2015 energy modified with the lipophilicity term (ref2015_memb).

Ab initio structure prediction in membrane proteins

Previous structure-prediction benchmarks started from canonical α helices or from monomers obtained from experimental structures of homodimers and used the bound-structures in grid search or rigid-body docking [11,16,3134]. Additionally, structure-prediction studies used experimental constraints, conservation analysis or correlated-mutation analysis to predict residue contacts in order to constrain conformational sampling [12,13,3540]. Several automated predictors dedicated to single-span dimers used shape complementarity [41,42], sequence-packing motifs [43] or comparative modelling [44], but to the best of our knowledge, ab initio modelling calculations, starting from a fully extended chain, have not been described. Given that deviations from canonical α helicity make important contributions to membrane-protein structure and function [29], we decided to apply a more stringent test using ab initio modelling, sampling all symmetric backbone, sidechain, and rigid-body degrees of freedom.

To test ab initio modelling using the new energy function, we applied the fold and dock protocol [45], which has been successfully applied in a variety of soluble-protein structure prediction and design studies [4649]. Briefly, fold and dock starts from an extended chain and conducts several hundred iterations of symmetric centroid-level backbone-fragment insertion and relaxation moves. It then applies symmetric all-atom refinement including all dihedral sidechain and backbone degrees of freedom (S1 Movie). To generate an energy landscape, we ran 5,000 independent trajectories (50,000 for high-order oligomers) for every 19 and 21 residue subsequence of each homooligomer, filtered the resulting models according to energy and structure parameters (Methods), and isolated the lowest-energy 10% of the models. Models were then clustered according to their energies and conformations, and five cluster representatives were compared to the experimental structures (Figs 3 and 4, Table 1). For comparison, we applied the described methodology using ref2015_memb, ref2015 and the current membrane-protein energy function in Rosetta, RosettaMP [11].

thumbnail
Fig 3. Energy landscapes for the TMHOP ab initio structure prediction benchmark.

All models that passed the energy and structure-based filters are shown as semi-transparent grey dots. Each of the five lowest-energy clusters is indicated by coloured circles (lowest-to-highest energy: red, blue, green, purple and black). The PDB entry is indicated on each panel, and the oligomeric state is specified by grey circles for oligomeric states than homodimers. Y-axes report the ref2015_memb energy normalised by the monomeric sequence length of each model.

https://doi.org/10.1371/journal.pcbi.1007318.g003

thumbnail
Fig 4. Structural comparison of the top-predicted model (by RMSD) produced by TMHOP and the experimentally determined structure (gold and grey, respectively).

PDB entry, RMSD and the model’s ranking (in energy) among the top-5 predicted models are indicated. Only accurately predicted structures (< 2.5 Å) are shown.

https://doi.org/10.1371/journal.pcbi.1007318.g004

thumbnail
Table 1. TMHOP Structure prediction benchmark.

Grey cells indicate RMSD < 2.5Å or fraction of native contacts using ref2015_memb > 0.7.

https://doi.org/10.1371/journal.pcbi.1007318.t001

The Protein Data Bank (PDB) contains 17 nonredundant (sequence identity <80%) NMR and X-ray crystallographic structures (adopted from Lomize et al. [44]) of natural single-span homodimers, two tetramers and one pentameric structure. Of the 20 cases in the benchmark, fold-and-dock simulations using ref2015_memb predicted near-native (<2.5 Å root-mean-square deviation [RMSD]) low-energy models for 14 homooligomers compared to nine using RosettaMP; of the 14 oligomers accurately predicted by ref2015_memb, the soluble energy function ref2015 also resulted in nine correct predictions. Prediction success rate using ref2015_memb was somewhat higher for right-handed relative to left-handed homodimers (80 and 50%, respectively; S2 Table), reflecting the tendency of right-handed homodimers to be more tightly packed [31], and in 11 cases, a near-native prediction was found among the top 3 lowest-energy predicted models (Fig 3). Of the three high-order oligomers tested, ref2015_memb successfully recapitulated the structures of the M2 tetramer and phospholamban pentamer. The PREDDIMER [42] and TMDIM [43] structure-prediction web servers, which do not use ab initio modelling, found models at <2.5 Å RMSD for nine and eight of the 17 homodimers, respectively. Thus, ab initio calculations using ref2015_memb accurately predict structures in two-thirds of the homooligomers in our benchmark, including high-order oligomers that cannot be predicted by other automated methods. Given the high success rate of the ab initio calculations, we developed a web-accessible server for predicting the structures of membrane-spanning homo-oligomers such as are observed in receptor tyrosine kinases and other membrane proteins (http://tmhop.weizmann.ac.il).

The successfully predicted homooligomers exhibit different structural packing motifs. The majority of the homodimer interfaces are mediated by the ubiquitous Gly-xxx-Gly motif [50], in which two small amino acids separated by three positions on the primary sequence enable close packing between the helices. There is uncertainty whether these motifs additionally form stabilising Cα hydrogen bonds [51,52]. Our structure-prediction analysis cannot resolve this uncertainty; note, however, that the new energy function ref2015_memb does not encode terms for Cα hydrogen bonds and yet recapitulates a large fraction of the homodimer structures (Figs 3 and 4, and Table 1). The underlying reason for successful prediction is that the dsTβL energetics encodes a strong penalty on exposing Gly residues to the lipid bilayer (approximately 2 kcal/mol/Gly at the membrane mid-plane; Fig 1), driving the burial of Gly amino acids within the homodimer interface (i.e., “solvophobicity”). Thus, lipophilicity and interfacial residue packing are sufficient for accurate structure prediction in a large fraction of the targets we examined.

In several single-spanning membrane receptors, conformational change in the membrane domain is thought to underlie receptor activation. For instance, past modelling of the ErbB2 membrane domain suggested two non-overlapping interaction sites involving two small-xxx-small motifs within the membrane domain and a molecular switching mechanism that underlies receptor activation [54]. The only experimental structure for ErbB2 involves the N-terminal small-xxx-small motif [55], which is recapitulated by the second predicted cluster (Fig 5A), whereas in the fourth predicted cluster, dimerisation is mediated via the C-terminal motif (Fig 5B), suggesting that in some cases, TMHOP may provide structural hypotheses for alternative binding modes for receptor homooligomeric domains.

thumbnail
Fig 5. Low-energy predicted models of ErbB2 recapitulate two binding modes thought to underlie receptor activation.

(A) The ErbB2 experimental structure (PDB entry 2JWA, grey) and a representative of the second-best cluster (orange). (B) The fourth-best cluster (orange), RMSD 5.8 Å from 2JWA, exhibits a right-handed conformation in which Gly68 and Gly72 are buried, in qualitative agreement with experiments [53] and modelling [54].

https://doi.org/10.1371/journal.pcbi.1007318.g005

Using the dsTβL assay, we also examined the effects of dozens of point mutations in glycophorin A on apparent association energy (ΔΔGbinding) in the bacterial plasma membrane [17]. As a stringent test of the new energy function, we conducted fold-and-dock calculations using both ref2015_memb and RosettaMP starting from the sequences of each of the point mutants. To reduce uncertainty in interpreting the experimental results, we focused on 32 mutations that exhibited large apparent energy changes in the experiment (|ΔΔGbinding| ≥ 2 kcal/mol) and compared the median computed ΔΔGbinding of the lowest-energy models to the experimental observation (Fig 6, S3 Table). ref2015_memb outperformed RosettaMP, correctly assigning 81% of mutations as stabilizing or destabilizing compared to 66% for RosettaMP. The six false-positive predictions using ref2015_memb are due to mutations at position Gly86, which is exposed to the membrane, explaining why our simulations predict these mutations to be neutral or favourable. Note that as observed in studies of mutational effects on stability in soluble proteins, the correlation coefficient between computed and observed values is low (Pearson r2 = 0.21 and 0.02 for ref2015_memb and RosettaMP, respectively) [5659]. Such low correlation coefficients provide an impetus for improving the energy function; however, as we previously demonstrated, discriminating stabilizing from destabilizing mutations is sufficient to enable the design of accurate, stable, and functionally efficient proteins [5964].

thumbnail
Fig 6. Predicted versus experimental ΔΔGbinding values of single-point mutations in glycophorin A.

The structure of every point mutant was predicted ab initio, and the median ΔΔGbinding relative to the wild type sequence is reported. Only point mutations that exhibited |ΔΔGbinding| ≥ 2 kcal/mol in the experiment were analysed. TP, TN, FP, and FN—true positive, true negative, false positive, and false negative, respectively.

https://doi.org/10.1371/journal.pcbi.1007318.g006

We next tested sequence-recovery rates using combinatorial sequence optimisation based on ref2015, ref2015_memb, and RosettaMP in a benchmark of 20 non-redundant structures (<80% sequence identity) ranging in size from 124–765 amino acids [65]. ref2015_memb outperformed the other energy functions, exhibiting 83% sequence recovery, on average, when each design was compared to the target’s natural homologs (Table 2). To our surprise, the soluble energy function ref2015 outperformed RosettaMP in this test and was almost as successful as ref2015_memb (78% overall success), implying that the packing and electrostatic models of ref2015 [25] enabled at least some of the improvement observed in sequence recovery by ref2015_memb (see S1 Table for a comparison of the energy functions). High sequence recovery in both buried and exposed positions implies that ref2015_memb may be applied effectively to design large and complex membrane proteins.

thumbnail
Table 2. Sequence recovery rates in Rosetta combinatorial sequence optimisation.

https://doi.org/10.1371/journal.pcbi.1007318.t002

Discussion

An accurate energy function is a prerequisite for automated modelling and design, and solvation makes a critical contribution to protein structure and function. The recent dsTβL apparent energies of insertion into the plasma membrane [17] enabled us to derive an empirical lipophilicity-based energy function for Rosetta. The results demonstrate that ref2015_memb outperforms RosettaMP in three benchmarks that are important for structure prediction and design. As ref2015_memb is based on the current state-of-the-art water-soluble Rosetta energy function, prediction accuracy is high for ref2015_memb both in soluble regions and in the core of the membrane domain. Thus, the lipophilicity preferences inferred from the dsTβL energetics together with the residue packing calculations in Rosetta enable accurate modelling in several ab initio prediction cases. The current energy function and the fold and dock procedure accurately model homooligomeric interactions in the membrane and the effects of point mutations, suggesting that they may enable the accurate design of homooligomeric single-span receptor-like transmembrane domains. The high accuracy models generated by the TMHOP method also suggest that laborious and often failed experiments to determine the structures of homooligomeric receptor membrane domains may be circumvented through ab initio modelling.

Nevertheless, certain important attributes of membrane-protein energetics are not yet addressed by ref2015_memb; for instance, atomic-level solvation and the impact on electrostatic interactions due to changes in the dielectric constant in various parts of the membrane are currently not treated [16,28] and warrant further research. Furthermore, the dsTβL profiles are based on measurements conducted on α-helical proteins in E. coli inner membranes. Ref2015_memb may, therefore, not perform as well in outer-membrane proteins or in proteins residing in membranes with a substantially different lipid composition. The benchmark reported here provides a basis on which improvements in the energy function can be verified. Furthermore, structure prediction in heterooligomers is important for understanding receptor cross-activation and for the design of membrane inhibitors [1]. In preliminary calculations, however, we found that fold and dock simulations of heterooligomeric systems fail to converge due to the much larger conformation space open to a non-symmetric system. A potentially exciting extension of the current work is to use the information on preferred crossing angles between membrane helices to constrain conformational sampling in heterooligomers [66,67].

We recently showed that evolution-guided atomistic design calculations, which use phylogenetic analysis to guide the design processes [68], enabled the automated, accurate and effective design of large and topologically complex soluble proteins. Designed proteins exhibited atomic accuracy, high expression levels, stability [59,60,69], binding affinity, specificity [64], and catalytic efficiency [62,63]. Membrane proteins are typically large and challenging targets for conventional protein-engineering and design methods. Looking ahead, we anticipate that evolution-guided atomistic design using ref2015_memb may enable reliable design in this important but often formidable class of proteins.

Methods

Rosetta source code

All code is available in the Rosetta release at www.rosettacommons.org (git version: b210d6d5a0c21208f4f874f62b2909f926379c0f). Command lines and RosettaScripts [70] are available in the supplement.

Membrane-insertion profiles

The original dsTβL insertion profiles [17] were modified to generate smooth and symmetric functions [24]. The polar and charged residues Asp, Glu, Gln and Asn, which exhibited few counts in the deep sequencing analysis, were averaged such that the insertion energy at the membrane core (-10 to 10 Å; negative values correspond to the inner membrane leaflet and positive values to the outer leaflet) was applied uniformly to the entire membrane span. The profile for His was capped at the maximal value observed in the experiment (2.3 kcal/mol) between 0Å (membrane midplane) and 20 Å. The dsTβL profile for Cys is unusually asymmetric. Cys residues are rare in membrane proteins [71] and are likely to have similar polarity to Ser. We, therefore, applied the profile measured for Ser to Cys. To convert the values from the dsTβL insertion profiles to Rosetta energy units (R.e.u.) they were multiplied by 2.94 following the interpolation reported in ref. [25]. The dsTβL profiles spanned 27 positions, and we correspondingly translated them to span -20 to +20 Å relative to the membrane midplane.

Residue lipophilicity

The context-dependent, one-body energy term MPResidueLipophilicity was implemented to encode the dsTβL insertion profiles in ref2015. Starting from an ideal poly Ala α helix embedded perpendicular to a virtual membrane, every position was mutated to all 19 identities, relaxed, and the energy difference between the ref2015 energy and the dsTβL energy was implemented in MPResidueLipophilicity. This process was repeated ten times to reach convergence, and the resulting energy profiles were fitted by a cubic spline [72], generating continuous, differentiable functions for all 19 amino acids relative to Ala, which was assumed to be 0 throughout the membrane. The splines were recorded in the Rosetta database and are loaded at runtime. Insertion -profile adjustments were done using a python3 script available at github.com/Fleishman-Lab/membrane_protein_energy_function.

Residue burial

The number of protein atoms within 6 and 12 Å of each amino acid’s Cβ atom is computed and transformed to a burial score (Eq 1). We used sigmoid functions which range from 0 to 1, corresponding to completely buried and completely lipid-exposed, respectively.

(1)

Where N is the number of heavy atoms and S and O determine the slope and offset of the sigmoids and are different for all-atom and centroid calculations. Each parameter has different thresholds at 6 or 12 Å. For all-atom calculations, S = 0.15 and 0.5 and O = 20 and 475, for 6 and 12 Å radii, respectively. For centroid-level calculations, S = 0.15 and 5 and O = 20 and 220 for 6 and 12Å radii, respectively. For each amino acid, the product of the 6 and 12Å sigmoid functions is taken, producing a continuous, differentiable function that transitions from buried to exposed states. These parameters were determined by visualising the burial scores of all amino acids in several polytopic membrane proteins of known structure.

Tilt-angle (Θ; Fig 2) penalty

All membrane-spanning helices reported in the PDBTM [73] dataset (version 20170210) were analyzed for their tilt angles with respect to the membrane normal. A second-degree polynomial was fitted to this distribution using scikit-learn [74].

(2)

As Bowie noted, the expected distribution function of helix-tilt angles is sin(Θ) [30]. We, therefore, used a partition function to convert the expected distribution (sin(Θ)) and observed one (Eq 2) to energy functions, finally subtracting the expected energy from the observed one to derive the helix-tilt penalty function: (3)

Where θ is given in degrees. In order to simplify runtime calculations, we approximated Eq 3 using a third-degree polynomial (using scikit-learn) (Fig 2).

(4)

Penalizing deviations from ideal α helicity

The MPHelicality energy term penalizes the energy of every position that exhibits ϕ-ψ torsion angles significantly different from ideal α helices. A paraboloid function was manually calibrated to express a penalty for any given (ϕ, ψ). The paraboloid centre, for which the penalty is 0, was set to the centre of the helical region according to the Ramachandran plot (ϕ = 60°, ψ = 45°) [75]. The paraboloid curvature was set to 25, such that the penalty is low throughout the ϕ-ψ torsion angles space observed for α helices [75]. As segments buried against the protein should not be penalized to the same extent as those completely exposed to the membrane, the burial approximation of Eq 1 is used to weight MPHelicality. Moreover, as the protein extends outside of the membrane, the penalty is attenuated with a function that follows the trend observed for the hydrophobic residues, Leu, Ile, and Phe (see Fig 1A). In effect, the MPHelicality term favours α helicity in lipid-exposed surfaces in the core of the membrane, thereby enforcing some of the electrostatic and solvophobic effects that are essential for correctly modelling the backbone but are not expressed in the residue-specific dsTβL energy profiles.

(5)

Where ϕ and ψ are given in degrees, z is the distance from the membrane midplane of residue i, and burial is calculated as in Eq 1.

A benchmark for structure prediction of single-span homooligomers

17 structures of single-span homodimers, two homotetramers and one pentamer were selected from the PDB (S2 Table). For each structure, a 20–30 residue segment comprising the membrane-spanning domain was manually chosen. A sliding window then extracted all 19 or 21 residue subsequences. For each subsequence, three and nine residue backbone fragments were generated using the Rosetta fragment picker application [76]. The fold-and-dock protocol [45] was used to compute 5000 models (50,000 models for tetramers and the phospholamban pentamer), and the lowest-energy 10% of the models were subsequently filtered using structure and energy-based filters (solvent accessible surface area >500 Å; shape complementarity [77] Sc>0.5; ΔΔGbinding <-5 R.e.u.; rotameric binding strain [78] < 4 R.e.u.; helicality <0.1 R.e.u. (computed using Eq 5); and closest distance between the interacting helices < 9 Å, as calculated by the filter HelixHelixAngle). For each target, the filtered models from all subsequences were then pooled together and clustered using a score-wise clustering algorithm. This is an iterative process, where each iteration calculates the RMSD of all unclustered models to the best-energy model, and removes the ones closer than 4 Å. RMSD to NMR structures were calculated with respect to the first model in the PDB entry.

A benchmark for ΔΔGbinding predictions of single-spanning homodimers

Glycophorin A mutants that exhibited |ΔΔGbinding| > 2 kcal/mol according to the dsTβL study [17] were modelled using the same fold-and-dock protocol described for the structure prediction of homodimers. To reduce computational load, we used a single sequence (73-ITLIIFGVMAGVIGTILLI-91), and the median of computed ΔΔGbinding for the top models was reported (models were filtered using structure and energy based filters; solvent accessible surface area > 600 Å, packing > 0.4, shape complementarity > 0.5, ΔΔGbinding < -10 R.e.u., binding strain < 4 R.e.u. MPHelicality < 0.1, minimal atomic distance between helices < 4.5 Å and minimal distance between helix vectors < 8 Å. Of these models, only the top 10% scoring models were used).

Sequence-recapitulation benchmark

20 structures of polytopic membrane-spanning proteins were taken from ref. [65], 11 of which were symmetric complexes. All were refined (eliminating sidechain conformation information before refinement), and for each protein, 100 designs were computed using combinatorial sequence design followed by sidechain and backbone minimization, and the lowest-energy 10 designs were checked for the fraction of mutations relative to the target protein. For each target protein, a multiple-sequence alignment was prepared: homologous sequences were automatically collected using BLASTP [79] on the nonredundant sequence database [80] with a maximal number of targets set to 3,000 and an e-value ≤ 10−4. All sequences were clustered using CD-hit [81] with a 90% sequence identity threshold. Sequences were then aligned using MUSCLE [82] with default parameters. A position-specific scoring matrix (PSSM) was calculated using PSI-BLAST [83]. In the sequence-recovery benchmark, where homologous sequences are considered, the substitution of a given position to an identity with a PSSM score ≥ 0 is considered a match.

Supporting information

S1 Table. Comparison of energy term weights.

ref2015_memb is identical to the ref2015 energy function, except for the addition of the mp_ terms, whereas RosettaMP is based on score12[10]. For a detailed explanation on energy terms see ref [25,26,84]. 1 used only in centroid-level sampling. 2 used in centroid-level sampling and full-atom sampling of single spanning proteins.

https://doi.org/10.1371/journal.pcbi.1007318.s001

(PDF)

S2 Table. Structures in the prediction benchmark.

1small-xxx-small in sequence.

https://doi.org/10.1371/journal.pcbi.1007318.s002

(PDF)

S3 Table. Raw data on the ΔΔGbinding prediction data.

https://doi.org/10.1371/journal.pcbi.1007318.s003

(PDF)

S4 Table. Structures for the sequence recovery benchmark.

https://doi.org/10.1371/journal.pcbi.1007318.s004

(PDF)

S1 Movie. A representative fold and dock simulation for glycophorin A.

In the first 25 seconds of the movie, simulations are in centroid mode, followed by all-atom refinement.

https://doi.org/10.1371/journal.pcbi.1007318.s005

(MP4)

S1 Fig. Adjustment and calibration of insertion profiles.

Each panel shows membrane insertion profiles for a different amino acid. Raw dsTβL data (purple dots, right-hand Y-axis), dsTβL adjusted profiles (red line), ref2015 insertion profiles (green) and ref2015_memb profiles (blue dashed line). Different residues affect the α helix differently, and therefore have different baselines. Note that Pro has no profile under ref2015_memb due to its effect on the backbone.

https://doi.org/10.1371/journal.pcbi.1007318.s006

(TIF)

S1 File. RosettaScripts XMLs and command lines for all benchmarks.

https://doi.org/10.1371/journal.pcbi.1007318.s007

(PDF)

S2 File. Top models from structure prediction benchmark.

https://doi.org/10.1371/journal.pcbi.1007318.s008.tar

(GZ)

Acknowledgments

We thank Rebecca Alford for help with the membrane protein framework in Rosetta and for suggesting the use of splines for fitting the dsTβL profiles and Shlomo-Yakir Hoch, Rotem Barzilay and Assaf Glick for developing the TMHOP web server. We also thank Adi Goldenzweig, Olga Khersonsky and Saar Shoer for helpful comments.

References

  1. 1. Yin H, Slusky JS, Berger BW, Walters RS, Vilaire G, Litvinon RI, et al. Computational Design of Peptides That Target Transmembrane Helices. Science. 2007;315: 1817–1822. pmid:17395823
  2. 2. Joh NH, Grigoryan G, Wu Y, DeGrado WF. Design of self-assembling transmembrane helical bundles to elucidate principles required for membrane protein folding and ion transport. Philos Trans R Soc Lond B Biol Sci. 2017;372. pmid:28630154
  3. 3. Feng X, Ambia J, Chen K-YM, Young M, Barth P. Computational design of ligand-binding membrane receptors with high selectivity. Nat Chem Biol. Nature Publishing Group; 2017;13: 715–723.
  4. 4. Lu P, Min D, DiMaio F, Wei KY, Vahey MD, Boyken SE, et al. Accurate computational design of multipass transmembrane proteins. Science. 2018;359: 1042–1046. pmid:29496880
  5. 5. Barth P, Senes A. Toward high-resolution computational design of the structure and function of helical membrane proteins. Nat Struct Mol Biol. europepmc.org; 2016;23: 475–480.
  6. 6. Koehler Leman J, Ulmschneider MB, Gray JJ. Computational modeling of membrane proteins. Proteins. Wiley Online Library; 2015;83: 1–24.
  7. 7. White SH, Wimley WC. Membrane protein folding and stability: physical principles. Annu Rev Biophys Biomol Struct. 1999;28: 319–365. pmid:10410805
  8. 8. Lazaridis T, Karplus M. Effective energy function for proteins in solution. Proteins: Structure, Function and Genetics. 1999;35: 133–152.
  9. 9. Lazaridis T. Effective energy function for proteins in lipid membranes. Proteins: Structure, Function and Genetics. Wiley Subscription Services, Inc., A Wiley Company; 2003;52: 176–192.
  10. 10. Yarov‐Yarovoy V, Schonbrun J. Multipass membrane protein structure prediction using Rosetta. Proteins: Struct Funct Bioinf. Wiley Online Library; 2006; https://onlinelibrary.wiley.com/doi/abs/10.1002/prot.20817
  11. 11. Alford RF, Koehler Leman J, Weitzner BD, Duran AM, Tilley DC, Elazar A, et al. An Integrated Framework Advancing Membrane Protein Modeling and Design. PLoS Comput Biol. 2015;11: 1–23.
  12. 12. Wang Y, Barth P. Evolutionary-guided de novo structure prediction of self-associated transmembrane helical proteins with near-atomic accuracy. Nat Commun. nature.com; 2015;6: 7196.
  13. 13. Ovchinnikov S, Kinch L, Park H, Liao Y, Pei J, Kim DE, et al. Large-scale determination of previously unsolved protein structures using evolutionary information. Elife. 2015;4: 1–25.
  14. 14. Mravic M, Hu H, Lu Z, Bennett JS, Sanders CR, Orr AW, et al. De novo designed transmembrane peptides activating the α5β1 integrin. Protein Eng Des Sel. 2018; pmid:29992271
  15. 15. Mravic M, Thomaston JL, Tucker M, Solomon PE, Liu L, DeGrado WF. Packing of apolar side chains enables accurate design of highly stable membrane proteins. Science. 2019;363: 1418–1423. pmid:30923216
  16. 16. Barth P, Schonbrun J, Baker D. Toward high-resolution prediction and design of transmembrane helical protein structures. Proc Natl Acad Sci U S A. 2007;104: 15682–15687. pmid:17905872
  17. 17. Elazar A, Weinstein J, Biran I, Fridman Y, Bibi E, Fleishman SJ. Mutational scanning reveals the determinants of protein insertion and association energetics in the plasma membrane. Elife. 2016;5. pmid:26824389
  18. 18. Shental-Bechor D, Fleishman SJ, Ben-Tal N. Has the code for protein translocation been broken? Trends Biochem Sci. 2006;31: 192–196. pmid:16530414
  19. 19. Karplus PA. Hydrophobicity regained. Protein Sci. 1997;6: 1302–1307. pmid:9194190
  20. 20. Vajda S, Weng Z, DeLisi C. Extracting hydrophobicity parameters from solute partition and protein mutation/unfolding experiments. Protein Eng. academic.oup.com; 1995;8: 1081–1092.
  21. 21. Gavel Y, Steppuhn J, Herrmann R, von Heijne G. The “positive-inside rule”applies to thylakoid membrane proteins. FEBS Lett. Wiley Online Library; 1991;282: 41–46.
  22. 22. von Heijne G. Control of topology and mode of assembly of a polytopic membrane protein by positively charged residues. Nature. 1989;341: 456–458. pmid:2677744
  23. 23. von Heijne G. The distribution of positively charged residues in bacterial inner membrane proteins correlates with the trans-membrane topology. EMBO J. Wiley Online Library; 1986;5: 3021–3027.
  24. 24. Elazar A, Weinstein J, Prilusky J, Fleishman SJ. The interplay between hydrophobicity and the positive-inside rule in determining membrane-protein topology. Proceedings of the National Academy of Sciences. 2016;113: 10340–10345.
  25. 25. Park H, Bradley P, Greisen P, Liu Y, Mulligan VK, Kim DE, et al. Simultaneous Optimization of Biomolecular Energy Functions on Features from Small Molecules and Macromolecules. J Chem Theory Comput. American Chemical Society; 2016;12: 6201–6212.
  26. 26. Rohl CA, Strauss CEM, Misura KMS, Baker D. Protein Structure Prediction Using Rosetta. Methods Enzymol. 2004;383: 66–93. pmid:15063647
  27. 27. Nolde DE, Arseniev AS, Vergoten G, Efremov RG. Atomic Solvation Parameters for Proteins in a Membrane Environment. Application to Transmembrane α-Helices. J Biomol Struct Dyn. Taylor & Francis; 1997;15: 1–18.
  28. 28. Lai JK, Ambia J, Wang Y, Barth P. Enhancing Structure Prediction and Design of Soluble and Membrane Proteins with Explicit Solvent-Protein Interactions. Structure. Elsevier; 2017;25: 1758–1770.e8.
  29. 29. Yohannan S, Faham S, Yang D, Whitelegge JP, Bowie JU. The evolution of transmembrane helix kinks and the structural diversity of G protein-coupled receptors. Proc Natl Acad Sci U S A. 2004;101: 959–963. pmid:14732697
  30. 30. Bowie JU. Helix packing angle preferences. Nat Struct Biol. 1997;4: 915–917. pmid:9360607
  31. 31. Fleishman SJ, Ben-Tal N. A novel scoring function for predicting the conformations of tightly packed pairs of transmembrane alpha-helices. J Mol Biol. Elsevier; 2002;321: 363–378.
  32. 32. Fleishman SJ, Unger VM, Ben-Tal N. Transmembrane protein structures without X-rays [Internet]. Trends in Biochemical Sciences. 2006. pp. 106–113. pmid:16406532
  33. 33. Fleishman SJ, Ben-Tal N. Progress in structure prediction of alpha-helical membrane proteins. Curr Opin Struct Biol. 2006;16: 496–504. pmid:16822664
  34. 34. Weiner BE, Woetzel N, Karakaş M, Alexander N, Meiler J. BCL::MP-fold: folding membrane proteins through assembly of transmembrane helices. Structure. 2013;21: 1107–1117. pmid:23727232
  35. 35. Matthews EE, Thévenin D, Rogers JM, Gotow L, Lira PD, Reiter LA, et al. Thrombopoietin receptor activation: transmembrane helix dimerization, rotation, and allosteric modulation. The FASEB Journal. Federation of American Societies for Experimental Biology Bethesda, MD, USA; 2011;25: 2234–2244.
  36. 36. Nugent T, Jones DT. Accurate de novo structure prediction of large transmembrane protein domains using fragment-assembly and correlated mutation analysis. Proc Natl Acad Sci U S A. 2012;109: E1540–7. pmid:22645369
  37. 37. Hopf TA, Colwell LJ, Sheridan R, Rost B, Sander C, Marks DS. Three-Dimensional Structures of Membrane Proteins from Genomic Sequencing. Cell. 2012;149: 1607–1621. pmid:22579045
  38. 38. Fleishman SJ, Harrington S, Friesner RA, Honig B, Ben-Tal N. An automatic method for predicting transmembrane protein structures using cryo-EM and evolutionary data. Biophys J. 2004;87: 3448–3459. pmid:15339802
  39. 39. Fleishman SJ, Harrington SE, Enosh A, Halperin D, Tate CG, Ben-Tal N. Quasi-symmetry in the cryo-EM structure of EmrE provides the key to modeling its transmembrane domain. J Mol Biol. 2006;364: 54–67. pmid:17005200
  40. 40. Fleishman SJ, Unger VM, Yeager M, Ben-Tal N. A Calpha model for the transmembrane alpha helices of gap junction intercellular channels. Mol Cell. 2004;15: 879–888. pmid:15383278
  41. 41. Polyansky AA, Volynsky PE, Efremov RG. Multistate organization of transmembrane helical protein dimers governed by the host membrane. J Am Chem Soc. 2012;134: 14390–14400. pmid:22889089
  42. 42. Polyansky AA, Chugunov AO, Volynsky PE, Krylov NA, Nolde DE, Efremov RG. PREDDIMER: a web server for prediction of transmembrane helical dimers.
  43. 43. Cao H, Ng MCK, Jusoh SA, Tai HK, Siu SWI. TMDIM: an improved algorithm for the structure prediction of transmembrane domains of bitopic dimers. J Comput Aided Mol Des. Springer; 2017;31: 855–865.
  44. 44. Lomize AL, Pogozheva ID. TMDOCK: An Energy-Based Method for Modeling α-Helical Dimers in Membranes. J Mol Biol. Elsevier B.V.; 2016; pmid:27622289
  45. 45. Das R, André I, Shen Y, Wu Y, Lemak A, Bansal S, et al. Simultaneous prediction of protein folding and docking at high resolution. Proceedings of the National Academy of Sciences. 2009;106: 18978–18983.
  46. 46. Morag O, Sgourakis NG, Baker D, Goldbourt A. Capsid model of M13 bacteriophage virus from Magic-angle spinning NMR and Rosetta modeling [Internet]. 2015.
  47. 47. DiMaio F, Leaver-Fay A, Bradley P, Baker D, André I. Modeling symmetric macromolecular structures in Rosetta3. PLoS One. 2011;6: e20450. pmid:21731614
  48. 48. Spreter T, Yip CK, Sanowar S, André I, Kimbrough TG, Vuckovic M, et al. A conserved structural motif mediates formation of the periplasmic rings in the type III secretion system. Nat Struct Mol Biol. 2009;16: 468–476. pmid:19396170
  49. 49. Boyken SE, Chen Z, Groves B, Langan RA, Oberdorfer G, Ford A, et al. De novo design of protein homo-oligomers with modular hydrogen-bond network-mediated specificity. Science. 2016;352: 680–687. pmid:27151862
  50. 50. Russ WP, Engelman DM. The GxxxG motif: a framework for transmembrane helix-helix association. J Mol Biol. 2000;296: 911–919. pmid:10677291
  51. 51. Senes A, Ubarretxena-Belandia I, Engelman DM. The Cα—H⋅⋅⋅O hydrogen bond: A determinant of stability and specificity in transmembrane helix interactions. Proceedings of the National Academy of Sciences. 2001;98: 9056–9061.
  52. 52. Yohannan S, Faham S, Yang D, Grosfeld D, Chamberlain AK, Bowie JU. A Cα- H⊙⊙⊙ O Hydrogen Bond in a Membrane Protein Is Not Stabilizing. J Am Chem Soc. ACS Publications; 2004;126: 2284–2285.
  53. 53. Jura N, Endres NF, Engel K, Deindl S, Das R, Lamers MH, et al. Mechanism for activation of the EGF receptor catalytic domain by the juxtamembrane segment. Cell. 2009;137: 1293–1307. pmid:19563760
  54. 54. Fleishman SJ, Schlessinger J, Ben-Tal N. A putative molecular-activation switch in the transmembrane domain of erbB2. 2002;
  55. 55. Bocharov EV, Mineev KS, Volynsky PE, Ermolyuk YS, Tkach EN, Sobol AG, et al. Spatial structure of the dimeric transmembrane domain of the growth factor receptor ErbB2 presumably corresponding to the receptor active state. J Biol Chem. 2008;283: 6950–6956. pmid:18178548
  56. 56. Potapov V, Cohen M, Schreiber G. Assessing computational methods for predicting protein stability upon mutation: good on average but not in the details. Protein Eng Des Sel. 2009;22: 553–560. pmid:19561092
  57. 57. Yin S, Ding F, Dokholyan NV. Eris: an automated estimator of protein stability. Nat Methods. 2007;4: 466–467. pmid:17538626
  58. 58. Kellogg EH, Leaver-Fay A, Baker D. Role of conformational sampling in computing mutation-induced changes in protein structure and stability. Proteins. 2011;79: 830–838. pmid:21287615
  59. 59. Goldenzweig A, Goldsmith M, Hill SE, Gertman O, Laurino P, Ashani Y, et al. Automated Structure- and Sequence-Based Design of Proteins for High Bacterial Expression and Stability. Mol Cell. 2016;63: 337–346. pmid:27425410
  60. 60. Campeotto I, Goldenzweig A, Davey J, Barfod L, Marshall JM, Silk SE, et al. One-step design of a stable variant of the malaria invasion protein RH5 for use as a vaccine immunogen. Proc Natl Acad Sci U S A. National Academy of Sciences; 2017;114: 998–1002.
  61. 61. Goldenzweig A, Fleishman SJ. Principles of Protein Stability and Their Application in Computational Design. Annu Rev Biochem. 2018;87: 105–129. pmid:29401000
  62. 62. Lapidoth G, Khersonsky O, Lipsh R, Dym O, Albeck S, Rogotner S, et al. Highly active enzymes by automated combinatorial backbone assembly and sequence design. Nat Commun. Nature Publishing Group; 2018;9: 2780.
  63. 63. Khersonsky O, Lipsh R, Avizemer Z, Ashani Y, Goldsmith M, Leader H, et al. Automated Design of Efficient and Functionally Diverse Enzyme Repertoires. Mol Cell. 2018;72: 178–186.e5. pmid:30270109
  64. 64. Netzer R, Listov D, Lipsh R, Dym O, Albeck S, Knop O, et al. Ultrahigh specificity in a network of computationally designed protein-interaction pairs. Nat Commun. 2018;9: 5286. pmid:30538236
  65. 65. Duran AM, Meiler J. Computational design of membrane proteins using RosettaMembrane. Protein Sci. 2018;27: 341–355. pmid:29090504
  66. 66. Zhang S-Q, Kulp DW, Schramm CA, Mravic M, Samish I, DeGrado WF. The membrane- and soluble-protein helix-helix interactome: similar geometry via different interactions. Structure. 2015;23: 527–541. pmid:25703378
  67. 67. Walters RFS, DeGrado WF. Helix-packing motifs in membrane proteins. Proc Natl Acad Sci U S A. 2006;103: 13658–13663. pmid:16954199
  68. 68. Khersonsky O, Fleishman SJ. Why reinvent the wheel? Building new proteins based on ready-made parts. Protein Sci. 2016;25: 1179–1187. pmid:26821641
  69. 69. Warszawski S, Katz A, Lipsh R, Khmelnitsky L, Ben Nissan G, Javitt G, et al. Optimizing antibody affinity and stability by the automated design of the variable light-heavy chain interfaces.
  70. 70. Fleishman SJ, Leaver-Fay A, Corn JE, Strauch E-M, Khare SD, Koga N, et al. RosettaScripts: a scripting language interface to the Rosetta macromolecular modeling suite. PLoS One. 2011;6: e20161. pmid:21731610
  71. 71. Senes A, Chadi DC, Law PB, Walters RFS, Nanda V, DeGrado WF. Ez, a depth-dependent potential for assessing the energies of insertion of amino acid side-chains into membranes: derivation and applications to determining the orientation of transmembrane and interfacial helices. J Mol Biol. Elsevier; 2007;366: 436–448.
  72. 72. Press WH, Vetterling WT, Teukolsky SA, Flannery BP. Numerical Recipes in C++: The Art of Scientific Computing. 2nd ed. New York, NY, USA: Cambridge University Press; 2002.
  73. 73. Kozma D, Simon I, Tusnády GE. PDBTM: Protein Data Bank of transmembrane proteins after 8 years. Nucleic Acids Res. Oxford University Press; 2012;41: D524–D529.
  74. 74. Pedregosa F, Varoquaux G, Gramfort A, Michel V, Thirion B, Grisel O, et al. Scikit-learn: Machine Learning in Python. J Mach Learn Res. jmlr.org; 2011;12: 2825–2830.
  75. 75. Page RC, Kim S, Cross TA. Transmembrane Helix Uniformity Examined by Spectral Mapping of Torsion Angles. Structure. 2008;16: 787–797. pmid:18462683
  76. 76. Gront D, Kulp DW, Vernon RM, Strauss CEM, Baker D. Generalized fragment picking in rosetta: Design, protocols and applications. PLoS One. 2011;6. pmid:21887241
  77. 77. Lawrence MC, Colman PM. Shape complementarity at protein/protein interfaces. J Mol Biol. 1993;234: 946–950. pmid:8263940
  78. 78. Fleishman SJ, Khare SD, Koga N, Baker D. Restricted sidechain plasticity in the structures of native proteins and complexes. Protein Sci. 2011;20: 753–757. pmid:21432939
  79. 79. Altschul SF, Gish W, Miller W. Basic local alignment search tool. Journal of molecular …. 1990;215: 403–410.
  80. 80. Wheeler DL, Chappey C, Lash AE, Leipe DD, Madden TL, Schuler GD, et al. Database resources of the National Center for Biotechnology Information. Nucleic Acids Res. 2000;28: 10–14. pmid:10592169
  81. 81. Li W, Godzik A. Cd-hit: A fast program for clustering and comparing large sets of protein or nucleotide sequences. Bioinformatics. 2006;22: 1658–1659. pmid:16731699
  82. 82. Edgar RC. MUSCLE: Multiple sequence alignment with high accuracy and high throughput. Nucleic Acids Res. 2004;32: 1792–1797. pmid:15034147
  83. 83. Altschul SF, Gertz EM, Agarwala R, Schäffer AA, Yu Y-K. PSI-BLAST pseudocounts and the minimum description length principle. Nucleic Acids Res. 2009;37: 815–824. pmid:19088134
  84. 84. Alford RF, Leaver-Fay A, Jeliazkov JR, O’Meara MJ, DiMaio FP, Park H, et al. The Rosetta All-Atom Energy Function for Macromolecular Modeling and Design. J Chem Theory Comput. 2017;13: 3031–3048. pmid:28430426