Prediction of protein complex structure using surface-induced dissociation and cryo-EM

Dataset

The dataset used for protein-protein docking with SID was described previously.¹³ In short, the dataset contained triose phosphate isomerase (homodimer, 8TIM), streptavidin (homotetramer, 1SWB), hemoglobin (heterotetramer, 1GZX), cholera toxin B (homopentamer, 1FGB), C-reactive protein (homopentamer, 1GNH), and serum amyloid P (homopentamer, 1SAC). For the complexes in the dataset, experimental SID AE values for specific interfaces were extracted from ERMS plots and normalized by the number of inter-subunit protein-protein interfaces as described previously.^{12, 13} Simulation of the 14 Å resolution density maps with noise is described in the SI.

Input structures for docking

Rather than solely redocking crystal structures for full complexes, here, we obtained docking input structures three different ways, which are visually depicted in Figure 1A and described in further detail in the SI. In short, when possible, we built homology models, used unbound crystal structures, and extracted bound crystal structures.

Rather than performing completely rigid docking simulations with single input structures, we performed extra backbone sampling for all inputs prior to docking. For each input structure, 10 models were built using Rosetta relax, normal mode analysis (NMA), and backrub, resulting in a total of 30 models with slightly different backbones. These structures were input into docking as described in the following sections.

Ensemble docking

Using a somewhat similar approach to our previous rigid docking method with SID,¹³ we performed docking in order to build (sub)complexes corresponding to each interface where AE was extracted. However, here we performed ensemble docking, where backbones were allowed to swap during the simulation from a pool of inputs. As described in the SI, these docking simulations were performed for all of the different types of input structures: 5 for HM, 2 for unbound, and 8 for bound+perturbed. A list of the specific interfaces and partners built for each input structure type is provided in Table S1. Additionally, for each input structure, the 30 extra backbone structures were pooled prior to docking. The orientation of the mobile chain was randomized using the -randomize2 flag. Ensemble docking in RosettaDock⁵² was performed where the prepacked monomers were allowed to swap during the simulation, allowing for backbone flexibility based on the sampled tertiary structures. For each input structure, 10,000 docked structures were built. The workflow for ensemble docking is shown in Figure 1B (left).

Symmetric docking

In addition to ensemble docking, symmetric docking was also performed to build full complexes in all cases based on the symmetry derived from stoichiometry provided by SID data. De novo symmetry was used when possible, however due to poor observed sampling when using D2 symmetry, complexes for streptavidin (1SWB and 5N8T) were docked using noncrystallographic point symmetry. A list of applied symmetries for each input structure is provided in Table S2. Note that for hemoglobin (1GZX, dimer of heterodimers) in both HM and bound+perturbed, the input heterodimer was the predicted structure from the ensemble docking (using the combined score described in the following section) and was subsequently docked with C2 symmetry to form the heterotetramer. A definition of C2 symmetry was also used in all other cases where dimers were used to build tetramers (2HBC, 1GZX_dimers, and 1SWB_dimers). For symmetric docking, backbones of input structures were not allowed to swap during the simulations, therefore separate trajectories of 500 structures for each of the 30 inputs were performed for a total of 15,000 docked structures. The workflow for symmetric docking is shown in Figure 1B (right).

Scoring strategy

Structures from both the ensemble and symmetric docking were scored using the Rosetta interface score (Isc). Then, cryo-EM score and SID scores were subsequently included in the overall scoring. The components corresponding to the score terms are visually depicted in Figure 1C. Rosetta Isc represents the Rosetta energy of the interface and the model with the lowest Isc represents the predicted model using Rosetta without experimental data. The Rosetta cryo-EM score (elec_dens_fast) was calculated based on the generated low-resolution density maps with noise.^{36, 53, 54} Calculation of SID scores was described in depth previously.¹³ In short, for each docked structure, the AE was predicted using the model (with rigidity factor [RF]¹² specifically calculated and averaged for the ensemble of input structures) and compared to the experimental AE. Structures were scored based on this difference using a fade function, with no penalty assigned to structures for which the predicted AE was within 100 eV of the experimental lab frame energy. If multiple AEs were available for a complex where all interfaces corresponding to the aforementioned AEs were built, the SID score used was the average SID score for those AEs. For the combined scores (Isc+cryo-EM and Isc+cryo-EM+SID), each individual score was first normalized based on the maximum observed absolute value for the score of the docked structures for that protein. The total score was a weighted sum of the three normalized terms, as shown in Equation 1. The weights (w_Isc, w_cryo-EM, and w_SID) used for ensemble and symmetric docking are provided in Table S3. We always chose the lowest scoring model as the predicted structure.

Confidence metric

A confidence metric was calculated in order to classify the prediction for each protein as either high or low confidence without the knowledge of the native structure. To derive the confidence metric, we first calculated RMSD (all atom, refined, and aligned) of each model with respect to the top overall scoring model for each protein. The confidence metric was defined as the P_near (lambda=10.0 and k_BT=2.0, using the RMSD values with respect to the top scoring model). We hypothesized that when accurate predictions (RMSD₁₀₀ of top scoring model less than 4 Å) were made, a larger number of favorably scoring models would be more similar to the top scoring model and superior funneling would be observed, thus P_near (with respect to the top scoring model) would be high. Proteins for which P_near was greater than 0.1 were identified as high confidence.

Previous work showed SID AE could improve model selection

While SID has typically been used to determine stoichiometry and connectivity for protein complexes, we established previously that a structurally-dependent measure, appearance energy (AE), for protein-protein interfaces could be extracted from SID.¹² We observed correlations between AE and interface features such as size and other measures relating to interface energy in Rosetta. From these correlations, methods to predict AE from structure were developed. AE prediction methods were then used to score sets of docked structures based on the agreement between predicted and observed AE for each structure.¹³ We showed that the structural information contained in the prediction model was sufficient to facilitate the discrimination between nativelike and non-nativelike poses. When the SID scoring function was used, for 3/9 cases, the RMSD of the selected model improved by more than 18 Å with respect to when no experimental data were included, and the inclusion of SID was never detrimental. We observed that the SID scoring function generally scored a variety of structures equally well (notably including many of the top models) but penalized the majority of bad structures. While this study laid the foundation of SID usage in protein complex structure prediction, showing for the first time that SID data could facilitate model selection for a diverse set of docked poses, there were a few shortcomings in our initial implementation that we addressed here. Originally, tertiary structure obtained from crystal structures of the complexes was used as input into the rigid docking. Additionally, for the majority of our docking simulations (6/9), we predicted structures of subcomplexes rather than entire complexes. For these reasons, we sought to perform a set of docking simulations that was much more realistic, by modifying the input structures to decrease the amount of knowledge required and also building structures of the entire complexes in more cases.

Different strategies for obtaining input structures for docking require less prior knowledge

To perform more realistic docking simulations and thus reduce the necessary prior knowledge for docking by not relying on the bound structures, we performed docking using alternative input structures. First, when possible, we built homology models of monomers. The monomer building blocks were accurately built (RMSD < 2 Å) using RosettaCM multi-template modeling for the following five proteins: 1GNH (RMSD of 1.36 Å), 1GZX (chain A: 1.70 Å, chain B: 1.74 Å), 1SAC (1.10 Å), 1SWB (0.59 Å), and 8TIM (1.12 Å). Next, we used unbound structures for the following two proteins: hemoglobin (heterodimer, 2HBC, corresponding to substructure of 1GZX) and streptavidin (monomer, 5N8T, corresponding to subunit of 1SWB). Unbound input structures were not available for the remaining proteins in the benchmark set. Unsurprisingly, the RMSD of both unbound input structures (with respect to bound) was small (1.48 Å for 2HBC to 1GZX and 0.51 Å for 5N8T to 1SWB). Finally, for all proteins in the dataset, we perturbed the backbones of the bound structures using several backbone-altering methods so as to not run the simulation with the tertiary structure that exactly matched the bound crystal structures. In addition to adjusting the backbones of bound structures, the same approach was used for all HM and unbound input structures with the purpose of adding flexibility into all docking as well. For each input structure, different backbones were sampled by generating 10 models each using Rosetta relax, backrub, and normal mode analysis (NMA), for a total of 30 models. Overall, when using these sampling methods, the structures did not change drastically, with RMSDs comparing inputs and outputs of typically ~0.5 Å and rarely above 1.5 Å, as shown in Figure S1. These ensembles of structures with slightly different backbone conformations were input into docking as described in the following sections.

Scoring strategy combines SID and cryo-EM

In addition to using data from SID to score structures, here, simulated low-resolution cryo-EM density maps with noise were also included in scoring. To evaluate the complex structure prediction when including experimental data, scoring was performed using three different scoring functions: Rosetta interface score (Isc, i.e., score with no experimental data included), Isc with Rosetta cryo-EM score, and a combined Isc, cryo-EM, and SID score. The combined score was a linear combination of the three normalized terms and represented the score when cryo-EM and SID data were both included. The weights used for Isc, cryo-EM, and SID were the same for each protein, however, different sets of weights were used for ensemble and symmetric docking (HM and crystal). The predicted structure for each scoring function was identified as the top scoring model. The ensemble and symmetric docking are described in detail in the following sections, but we also compared the score vs. RMSD plots for the SID scores individually to previous work. Example SID score vs. RMSD plots for 1SWB are shown in Figure S2 for bound rigid redocking (previous work),¹³ HM ensemble docking, unbound (5N8T) ensemble docking, perturbed+bound ensemble docking, HM symmetric docking, unbound (5N8T) symmetric docking, and perturbed+bound symmetric docking. We note that while SID scores were not normalized in previous work and RMSD values were calculated differently (they were essentially ligand RMSD [LRMSD] values), the results are generally comparable. In all cases, the general trend that most accurate (low RMSD) structures were scored favorably and most inaccurate structures (high RMSD) were scored unfavorably continued to be observed. Although the SID scores are not enough to predict structures, they can be beneficial when combined with other scores (such as Isc and cryo-EM scores described here), when structures are scored slightly incorrectly. This example demonstrates that the relationship between SID AE and interface composition are similar here to previous work.

Ensemble docking resulted in accurate modeling of subcomplexes, strengthened by SID data

To increase the conformational space sampled during docking, we performed flexible ensemble docking (as opposed to rigid docking) using RosettaDock. This docking was performed for all input structures (5 for HM, 2 for unbound, and 8 for bound+perturbed; see Table S1) and was used to build complexes across all interfaces with measured experimental AE. Figure 2A shows the RMSD₁₀₀ of the top scoring model using each scoring function. When no experimental data were included (Isc, gray), the RMSD₁₀₀ of the predicted model was less than 4 Å for only a few cases (2/5 for HM, 1/2 for unbound, and 4/8 for bound+perturbed). When cryo-EM data were included in scoring (Isc+cryo-EM), the RMSD₁₀₀ of the predicted structure was less than 4 Å for 12/15 cases. After the inclusion of SID data (Isc+cryo-EM+SID), all (15/15) cases resulted in accurate prediction. For predictions with this level of accuracy, the method was able to identify the correct complex topology. Remarkably, for 5 cases, the lowest RMSD model was also the top scoring model for the combined Rosetta, cryo-EM, and SID score. Predicted structures (and RMSD₁₀₀) are shown in Figure 2C (aligned to native shown in green) for 1GZX (HM), 5N8T (unbound), and 1SAC (bound+perturbed). Agreement with native crystal structure improved significantly when all experimental data were included.

The average model quality of the predicted structure generally improved with the inclusion of additional experimental data (both cryo-EM and then SID) when considering the CAPRI (Critical Assessment of PRediction of Interactions) metrics, as shown in Table S4 (left) and Figure 2B. For example, the average DockQ score,⁵⁵ a quantitative measure of model quality inspired by CAPRI categories, improved from 0.30 when no experimental data were included to 0.48 when cryo-EM was included, and further improved to 0.58 when SID data were also included. The same trend was observed for average RMSD (all atom, aligned), RMSD₁₀₀, F_nat, F_non-nat, iRMSD, and LRMSD. Using CAPRI model quality designations, 9/15 predicted structures using the combined score were identified as medium or high quality, with no cases identified as incorrect (i.e., all were at least acceptable quality), compared to 9/15 incorrect when scoring without experimental data. Improvement of funneling for the score vs. RMSD distributions was also observed when experimental data were included, with the P_near improving by an average of 0.27 (on a scale of 0 to 1). The RMSD for the top 50 scoring models improved by an average of 2.8 Å when SID and cryo-EM data were included as well. While much of the improvement when experimental data were included was due to the inclusion of the cryo-EM score, the SID score did play a significant role in a few cases, namely 1GNH (HM), 1GZX (HM), and 5N8T (unbound). Furthermore, accurate structures were predicted in all cases.

Symmetric docking resulted in accurate modeling overall

In order to predict structures of entire complexes, rather than building only across interfaces identified with AE, we performed symmetric docking using Rosetta SymDock. In comparison to our previous docking studies, the symmetric docking was new and extensively increased the breadth of SID-guided modeling. In addition to each input structure (5 for HM, 2 for unbound, and 8 for perturbed bound; see Table S2), as the dimer of heterodimers could not be directly built with the structures of the monomers and symmetric docking, input structures for two cases were obtained from the output of the ensemble docking (1GZX HM and bound+perturbed, with RMSD₁₀₀s of 1.36 Å and 0.33 Å respectively). Figure 3A shows RMSD₁₀₀ of the top scoring model using each scoring function. When no experimental data were included (Isc, gray), the RMSD₁₀₀ of the predicted model was less than 4 Å for only 5/15 cases (0/5 for HM, 1/2 for unbound, and 4/8 for bound+perturbed). When data from both SID and cryo-EM were included, the RMSD₁₀₀ of the predicted structure was less than 4 Å for 14/15 cases (4/5 for HM, 2/2 for unbound, and 8/8 for perturbed bound). Predicted structures are shown in Figure 3C (aligned to native shown in green) for 1GNH (HM), 2HBC (unbound), and 1FGB (bound+perturbed). Agreement with native crystal structure improved significantly when all experimental data were included for these cases.

When considering the CAPRI metrics, as shown in Table S4 (right) and Figure 3B, the results were generally best when cryo-EM and SID data were included. The averages for all RMSD-based metrics improved when more data were included. The average results appeared slightly worse for f_nat, f_non-nat, and DockQ score upon the inclusion of SID (comparing Isc+cryo-EM to Isc+cryo-EM+SID), however this was due to a single case, 8TIM (average DockQ score improves slightly from 0.54 to 0.55 for the other proteins). To illustrate this minor, non-deleterious difference in structure for 8TIM, the structures predicted from Isc+cryo-EM and Isc+cryo-EM+SID are shown in Figure S3. Noticeably, there are very few differences in the accuracy of the predicted structure.. When assessing model quality using DockQ score, 10/15 models were identified as inaccurate without experimental data. However, when the combined score including SID was used, 7/15 were identified as medium or high quality (only two inaccurate). For symmetric docking, improvement in funneling was observed when SID and cryo-EM data were included. P_near, goodness of funneling in score vs. RMSD distributions, improved over Isc by an average of 0.37 (on a scale of 0 to 1). Improvement of the average RMSD of 7.3 Å for the top 50 scoring models was also observed when SID and cryo-EM data were included. In summary, the symmetric docking results were slightly improved by SID, but the inclusion of both cryo-EM and SID produced the best results. Importantly, for the first time, all full complexes were built using SID data and input structures required much less knowledge than previous work.

Dependence of prediction accuracy on input structure followed expected trends

Due to the difference in accuracy of the input structures (with respect to their structures in the complex), the overall docking results were expected to depend on the method used to obtain them. While Table S4 shows average metrics over all proteins for ensemble and symmetric docking, we also examined averages of each type of input structure (HM, unbound, and bound+perturbed) separately. The category-dependent results of F_nat, iRMSD, LRMSD, and DockQ score are shown in Figures 2B and 3B for ensemble and symmetric docking respectively. First, similar to the total average results, the accuracy of the predicted structures was generally best when using all experimental data (cryo-EM and SID). When considering the difference between input structure categories, we expected the structures built from bound+perturbed to be the best followed by unbound and HM based on accuracy of input structures prior to docking. The results typically followed the expected trend for both ensemble and symmetric docking. The results were best for bound+perturbed, followed by unbound and HM, respectively. We also note that there were no cases where the inclusion of SID improved the prediction results for the bound+perturbed docking simulations (Isc+cryo-EM+SID compared to Isc+cryo-EM), which is indicated by the very similar results when comparing Isc+cryo-EM and Isc+cryo-EM+SID for bound+perturbed. However, this did not indicate poor performance from the SID score, rather this occurred because the scoring was already excellent. On the other hand, SID improvements were much more common for the docking simulations using homology models since the scoring was slightly off without SID. To summarize, whether SID can improve prediction results depends on the scoring before including SID (in that the SID score is most beneficial when the scoring was close to being accurate).

Confidence metric effectively separated nativelike and non-nativelike structures

As not all predictions resulted in accurate structures (RMSD₁₀₀ less than 4 Å), we sought to quantify a prediction confidence metric that could be calculated without knowledge of native structure. For the confidence metric (applied to both ensemble and symmetric docking), P_near of the score vs. RMSD distribution with respect to the top scoring model was calculated. We hypothesized that when accurate predictions were made, a larger quantity of similar structures would score well and thus the funneling of the score vs. RMSD distribution would be stronger (as quantified by a higher P_near).²² Figure 4 shows the RMSD₁₀₀ of the top scoring model as a function of the confidence metric. Proteins for which the P_near was greater than 0.1 were identified as high confidence, as indicated by the vertical solid line. The prediction was accurate for all proteins identified as high confidence (26/26), as indicated by the dotted horizontal line. All inaccurate predictions were flagged as low confidence (1/1), with only three accurate predictions flagged as low confidence. All ensemble docking predictions (black, all accurate) were identified as high confidence (15/15). Overall, this confidence metric served as an effective indicator of prediction accuracy.