Frustration in the protein-protein interface plays a central role in the cooperativity of PROTAC ternary complexes

Abstract

Targeted protein degradation using proteolysis-targeting chimeras (PROTACs) offers a promising strategy to eliminate previously undruggable proteins. PROTACs are bifunctional molecules that link a target protein with an E3 ubiquitin ligase, enabling the formation of a ternary complex that promotes ubiquitination and subsequent proteasomal degradation. Although many ternary complex structures are available, understanding how structural features relate to PROTAC function remains challenging due to the dynamic nature of these complexes. Here we show that the interface between the target protein SMARCA2 and the E3 ligase VHL is conformationally flexible and stabilized by interactions involving disordered loops. Using molecular dynamics simulations and X-ray crystallography of SMARCA2–VHL complexes bound to five different PROTACs, we find that interfacial residues often adopt energetically suboptimal, or ‘frustrated,’ configurations. We further show that the degree of frustration correlates with experimentally measured cooperativity for a set of 11 PROTACs. These findings suggest that quantifying interface frustration provides a rational, structure-based approach to guiding PROTAC design.

Introduction

Proteolysis-targeting chimeras (PROTACs) are an emerging class of therapeutics that have opened opportunities for targeting a wider range of proteins, as they enable an alternative mode of action when compared to conventional small molecule drugs, e.g., inhibitors¹. PROTACs are heterobifunctional molecules consisting of two binding motifs each of which binds to the POI and an E3 ubiquitin ligase (E3), and are joined by a linker. Upon formation of a ternary complex (TC) composed of the POI, the PROTA,C and the E3 ligase (POI::PROTAC::E3), the E3 catalyzes the transfer of ubiquitin from the E2 ubiquitin conjugating ligase to the POI. The ubiquitinylated POI is then subjected to the ubiquitin–proteasomal degradation machinery, leading to the removal of the POI from the cell². Unlike conventional protein inhibitors that rely on an occupancy-driven mode of action, PROTACs operate through a different mechanism where only a brief interaction is needed to trigger the target protein’s degradation, essentially acting as a “catalytic” event rather than relying on constant occupancy. Very often, the ternary complex formation by a PROTAC is mediated through proximity-induced formation of de novo interactions within the complex³. The tendency of a PROTAC to stabilize such a ternary complex can be quantified using the cooperativity value α, i.e., the binding affinity ratio of a PROTAC in the binary complex (PROTAC::POI, PROTAC::E3) versus the ternary complex. The cooperativity is positive (α > 1) when the binding affinity to the ternary complex is greater than the affinities to the binary complexes, thus stabilizing the ternary complex. The cooperativity is negative (α < 1) when the binding affinity of the PROTAC to the ternary complex is lower than to the binary complexes, hence destabilizing the ternary complex. A non-cooperative cooperativity (α = 1) means no change in the binding affinity for the ternary complex when compared to the binary complexes^4,5.

Although the three-dimensional structures of TCs inform us about the nature of the POI-E3 interface^6,7,8,9, structure-based design of PROTACs still remains challenging. This is in part because some ternary complexes are highly dynamic^10,11,12,13, and the crystal structure of an individual ternary complex represents only one of, in some cases, many energy minima in the energy landscape of the TC conformational ensemble. Multiple thermodynamic properties are needed to characterize the measured cooperativity and degradation properties of PROTACs prior to synthesis due to the highly dynamic energy landscape of the TCs.

Computational methods play a critical role in understanding the dynamic behavior of the PROTAC-mediated TCs. Multiple methods for the modeling or prediction of ternary complex structures have been published^14,15,16, many of which are based on protein-protein docking. A recent study used weighted ensemble molecular dynamics (MD) simulations to model the formation of the ternary complex ensembles¹⁷. The calculated cooperativity of the PROTACs alone is not applicable to all PROTAC systems. There is a need for more descriptors derived from the PROTAC ternary complex ensembles to correlate with the cooperativity, even when the 3D-structures of the ternary complexes are known. Often, the protein-protein interface (PPI) residues of the POI and the E3 ligase are in unstructured or disordered regions, and hence are highly dynamic^10,11,18. Furthermore, PROTAC-mediated PPIs are not driven through coevolution of interface residues, which, in naturally occurring PPIs, is often the basis of the PPI stabilization. Given the lack of correlation between calculated binding affinities of PROTACS and cooperativity¹⁷, we looked for structural descriptors that characterize the PPIs in unstructured regions and their role in regulating the protein activity. Previous studies have shown that distinct residue positions that play a critical functional role are often in a suboptimal energetic state before being involved in intramolecular or intermolecular protein interaction¹². The presence of such locally perturbed thermodynamic states in proteins is referred to as ‘frustration’^9,15.

Frustration was developed as a concept to quantify the conformational variations and heterogeneity often seen in dynamic PPIs formed through disordered regions^19,20. Since the TCs are highly dynamic, we hypothesized that the frustration of the interactions at the POI-E3 ligase interface could be a relevant factor in determining the cooperativity of the ternary complex.

We investigated a series of PROTACs that mediate the degradation of SMARCA2 as the POI through recruitment of the E3-ligase VHL, to assess whether frustration could predict ternary complex stability. The helicases SMARCA2 and SMARCA4 are close paralogs and crucial components of the BAF chromatin remodeling complex. Based on a synthetic lethality concept, SMARCA2 has been proposed as a drug target for certain tumors exhibiting SMARCA4-deficiency (e.g., non-small cell lung cancer (NSCLC))²¹. The von Hippel-Lindau tumor suppressor (VHL) is an E3 ubiquitin ligase in which small molecule binders with suitable exit vectors for PROTAC design have been described before²². Farnaby, Koegl et al. utilized a SMARCA bromodomain (BD) binder reported by Genentech (GEN-1)²³ to synthesize the first VHL-based SMARCA2/4-degrading PROTACs and solved three ternary complex crystal structures to guide PROTAC design⁸. Later, Kofink, Trainor, Mair et al. discovered alternative SMARCA^BD binders with improved physicochemical properties. Guided by a further three TC crystal structures, SMARCA2-selective PROTACs with oral bioavailability could be designed⁹.

Here, we show that protein-protein interface frustration correlates with the cooperativity of PROTAC-induced ternary complex formation. We determine X-ray crystal structures of SMARCA2–VHL complexes bound to four PROTACs at resolutions ranging from 2.2 Å to 3.74 Å. Using these structures, we perform all-atom molecular dynamics simulations to characterize the conformational dynamics of the complexes and calculate residue-level frustration at the SMARCA2–VHL interface. We find that PROTACs with higher cooperativity exhibit a greater number of frustrated residue pairs at the interface. Extending this analysis to modeled ternary complexes of 11 GEN-1–based PROTACs, we observe the same trend. These results suggest that interfacial frustration may serve as a structure-based metric for prioritizing PROTAC candidates when cooperativity is an important driver of degradation. However, this approach may not be applicable to systems where degradation occurs independently of cooperativity, such as PROTACs that engage HBF-type E3 ligases²⁴.

Results

Synthesis of PROTACs for SMARCA2-VHL complexes and measurement of cooperativity

We investigated 16 PROTACs (shown in Fig. 1A and Supplementary Fig. 1) in this work and determined their cooperativity (α). Eleven of these PROTACs, P6 to P20, are based on the original GEN-1 SMARCA^BD binder²³, while P1 to P5 incorporate alternative SMARCA binder versions⁹. All PROTACs contain the same VHL binder, VH101 (compound 10)²², except for P3 in which the VH101 fluorine was replaced by a dimethyl amino function⁹. A phenolic hydroxyl group on VH101 serves as an exit vector to the linker, except for P2, P3, and P18, which use the benzylic position instead. P19 and P20 contain the cis-diastereomer of the hydroxy-pyrrolidine moiety in VH101, rendering the binder inactive to VHL⁸. Seven PROTACs, P1, P3, P4, P9, P10, P11, and P18, are disclosed, and their synthesis is described in the Supplementary Information.

Fig. 1: Ternary complex crystal structures of SMARCA2^BD, VCB, and PROTACs P1 to P5.

A The VHL-binding part of the PROTACs is oriented in similar positions. SMARCA2^BD and VHL are shown in green, PROTACs in white. The electron density map of the refine 2f_of_c map is depicted at a contour level of 1.0 σ for the PROTACs. Omit maps for ligands P1 to P5 (f_of_c maps at a contour level of 3.0 σ) are provided in Supplementary Fig. 8. B Protein residues that interact with PROTAC in the five crystal structures include VHL residues with indices below 1000 and SMARCA2 residues with indices above 1000.

A previously developed time-resolved fluorescence resonance energy transfer (TR-FRET) competition assay was used to assess the cooperativity of PROTAC binding to SMARCA2 and VHL⁸. Briefly, a biotinylated SMARCA2 probe²³ was incubated with the His₆-tagged SMARCA2^BD protein, and the proximity of both components was measured from a FRET donor/acceptor pair tagged to streptavidin or an anti-His antibody, respectively (experimental details are provided in the Supplementary Methods). The loss of the proximity signal, after addition of a concentration series of each PROTAC in the presence or absence of saturating concentrations of VCB (a pre-formed complex of VHL, Elongin-C, and Elongin-B), yielded the ternary and binary IC₅₀, respectively. The cooperativity (α) was defined as the ratio: IC₅₀(binary) / IC₅₀(ternary), so that the affinity gains of a PROTAC to SMARCA2 after addition of VCB, i.e., a lower ternary IC₅₀, lead to higher α values. Binary and ternary IC₅₀ values, as well as α values for all 16 PROTACs, are reported in Supplementary Table 1.

Crystal structures of five PROTAC ternary complexes

For two of the 16 PROTACs, P2 and P13, ternary complex crystal structures have previously been published (PDB codes 7Z77⁹ and 6HAY⁸, respectively). We solved the crystal structures of the PROTAC ternary complexes: P1, P3, P4, and P5, to enable comparison of their binding mode with existing crystal structures (Fig. 1A and Supplementary Tables 2 to 6). In all cases, the SMARCA binders bound in the acetyl-lysine binding site as previously described, i.e., with the key interactions of the quinazolinone (P4: benzoxazinone) core to Leu1456 and Asn1464 of SMARCA2^BD9. The VHL binder was consistently found to bind as described previously²². The interactions between SMARCA2^BD and VHL varied substantially across the PROTACs, as detailed later in the manuscript.

Comparison of the crystal structures of SMARCA2-VHL PROTAC ternary complexes shows that the ternary complexes are dynamic

The cooperativity of the five PROTACs, P1 to P5, ranges from 0.2 to 93 (see Fig. 2A, Supplementary Table 1). Structural alignment of the five TC crystal structures onto the VHL domain shows the high flexibility of the ternary complexes: SMARCA2 and VHL adopt very different conformations relative to each other with different PROTACs (Fig. 2B). These differences can stem both from the differences in the SMARCA2 binding motif and/or the linker lengths of the PROTACs. The root mean square deviation (RMSD) in the coordinates of the backbone atoms of the SMARCA2 domain ranges from 2.8 Å to 15.6 Å when aligned by the VHL domain of the P1 ternary complex crystal structure. The SMARCA2 domain orientations vary even between the two low cooperativity PROTACs (P4, α = 1.3 and P5, α = 0.2).

Fig. 2: Structural comparison of five ternary complexes.

A 2D Chemical structures of the five PROTACs. B Crystal structures of the five PROTAC-bound ternary complexes aligned onto their VHL domains.

MD simulations show that each PROTAC ternary complex samples multiple conformational states

Starting from the respective crystal structures, we performed all-atom MD simulations with the ternary complexes of the five PROTACs shown in Fig. 2A, for 3 μs each using the AMBER18 MD simulation package²⁵ (details in the Methods section). To analyze and visualize the dominant motions, Principal Component Analysis (PCA) was performed on the combined trajectories of all five complexes (see Methods for details). The top three weighted principal components capture nearly 80% of the variance (as shown in Supplementary Fig. 2). Principal component 1 (PC1) represents a rotation of the SMARCA2 domain with respect to the VHL domain, while PC2 represents a swinging motion of the SMARCA2 with respect to VHL (Fig. 3A, Supplementary Movies 1 and 2). Projection of the snapshots from MD simulations of all the ternary complexes (Fig. 3B and Supplementary Fig. 3) in the PC1 and PC2 coordinates shows that the high cooperativity PROTAC P1 complex samples a relatively compact conformational ensemble compared to the broader ensembles of the low cooperativity PROTAC complexes, for example, P4 and P5. The respective crystal structures thus only represent single low-energy snapshots in this wide conformational space, as seen in Fig. 3B.

Fig. 3: Principal component analysis of ternary complex dynamics.

A Domain motions represented by the Principal Component 1 and 2 (PC1 and PC2) resulting from PC analysis of the MD trajectories of all the five PROTAC ternary complexes. B Individual contour maps for all five PROTACs, projected onto the common PC landscape. C representative conformations from the highest populated cluster of each from the PROTAC bound TC MD simulations overlaid, after aligning their VHL domains. D Scatter plot of Shannon entropy (x-axis) versus the experimentally determined cooperativity values (y-axis) is shown as dots. PROTACs with α > 10 are colored in green, those with 2 <α < 10 in orange, and α < 2 cases are colored in red.

We clustered the conformations from the MD simulation trajectories of each ternary complex by PC1 and PC2, and Supplementary Table 2 shows the number of conformation clusters and their relative population for all five complexes. Each TC shows multiple conformation clusters indicative of their dynamic nature. The PROTACs with high and medium cooperativity: P1, P2, and P3 show fewer conformational clusters when compared to the PROTACs with low and negative cooperativity: P4 and P5, where the ensemble population is distributed over a larger number of clusters. This indicates that PROTACs with low cooperativity impart more flexibility to the ternary complex than PROTACs with high cooperativity. Figure 3C shows the alignment of representative structures extracted from the most populated conformational cluster (see Methods) for all five PROTACs. The positional spread of the SMARCA2 domain in relation to VHL is wider than observed when comparing only crystal structures (Fig. 2B). We then calculated the Shannon entropy using the conformation cluster population distributions (see Methods) for the five PROTAC complexes as a measure of their conformational diversity. As seen in Fig. 3D, the Shannon entropy shows weak but significant (p value: 0.036)²⁶ inverse correlation with cooperativity. The two PROTACs with low cooperativity, P4 and P5, show significantly higher Shannon entropy than the PROTACs, P1-P3, with high and medium cooperativity. We calculated the theoretical cooperativity using MMGBSA calculations with the MD simulation trajectories that include desolvation penalties (see Methods for details). The results in Supplementary Fig. 4 show no correlation to the measured cooperativity.

The Spatiotemporal residue contact heat map from MD simulations reflects the cooperativity of the PROTACs

In the previous section, we analyzed the conformational diversity of TCs with five different PROTACs using MD simulations. To understand how PROTAC flexibility could have contributed to these differences, we performed backbone RMSD clustering on the SMARCA2-VHL complex conformations from the MD simulations and extracted one representative structure from each cluster (Fig. 4A and B). As seen in Fig. 4A, PROTAC P1 with high cooperativity showed closely related PROTAC conformations, in contrast to PROTAC P5 with low cooperativity (Fig. 4B), which showed multiple disparate conformations. The higher flexibility of the PROTACs with low cooperativity could be due to the linker structure (Supplementary Fig. 5), which in turn could lead to increased flexibility in the PPI. To characterize the PROTAC flexibility further, we identified the VHL and SMARCA2 residues that are in contact with the PROTACs in the crystal structures. As seen in Fig. 4C, the number of PROTAC-protein contacts gleaned from the crystal structures are similar for all five PROTACs and do not distinguish between low and high cooperativity PROTACs. On the other hand, the contact frequency heatmap between the PROTAC and VHL or SMARCA2, as derived from MD simulation trajectories (shown in Fig. 4C and D), shows that the high and medium cooperativity PROTACs (P1 to P3) have a higher number of persistent PROTAC:protein contacts (P1:21, P2:23, P3:24) when compared to the low cooperativity PROTACs (P4:15, P5:19). In addition, the PROTAC contacts with SMARCA2 are varied for the different PROTACs (Fig. 4D). These frequency heatmaps reflect the dynamic PPI mediated by the PROTACs. The interaction fingerprints derived from MD simulations are therefore more indicative of PROTAC cooperativity when compared to fingerprints derived from crystal structures (Fig. 4D).

Fig. 4: Conformational ensembles and interface contacts of ternary complexes.

A The representative structures of the top three most populated TC conformation clusters of P1, superimposed onto the VHL domain. The three conformations are colored from dark to light color, SMARCA2 colored in blue, VHL colored in grey, PROTAC colored in green. The unstructured loop regions are colored in pink. B The representative structures of the three most populated clusters of P5. C Total number of residue contacts in PROTAC-VHL and PROTAC- SMARCA2 interfaces, as observed in crystal structure and MD simulations (contact formed for >50% of simulation time). In crystal structures and MD simulation ensembles, the number of contacts between the PROTAC and VHL are labeled by a yellow line on the bar plot. D The contact heatmap between PROTAC and proteins (POI and E3) obtained from the MD ensembles. Red “+” symbols indicate contacts formed in crystal structures. Residue numbers for VHL are below 1000, and for SMARCA2 are above 1000. Residues located in the unstructured loop regions are shown in magenta labels.

Residue contacts in the SMARCA2-VHL interface in PROTAC ternary complexes are dynamic and located in unstructured regions

Previous structural and biophysical studies have shown that the contact surface area between the POI and E3 ligase is critical to the cooperativity and can impact the degradation properties of PROTACs^27,28. We analyzed the pairwise residue interactions in the VHL-SMARCA2 interface, both in the crystal structures and in the MD trajectories of the five PROTACs (Fig. 5A). Comparison of the contact heatmaps in Figs. 4D and 5A, shows that most of the protein-protein interface contacts are mediated by the PROTAC, but also a few direct protein-protein contacts exist. The number of such direct contacts in the crystal structures does not allow discrimination between the high cooperativity PROTAC P1 from the low cooperativity PROTAC P4, although the negatively cooperative PROTAC P5 shows a negligible number of SMARCA2-VHL contacts. Interestingly, during the MD simulations, additional contacts were formed and many of the contacts observed in the crystal structures were weakened. The frequency heatmap calculated from the MD simulation trajectories shows that the high and medium cooperativity PROTACs exhibit more direct PPI contacts than the low cooperativity PROTACs, P4 and P5 (Fig. 5B). This loss of initial residue contacts in the SMARCA2-VHL interface can be attributed to the high flexibility of the ternary complexes. In contrast, the PROTAC-mediated SMARCA2-VHL contacts are maintained in the MD simulations, demonstrating the critical role of the PROTAC in stabilizing the ternary complex. Other MD simulation studies on PROTAC ternary complexes have also shown this behavior^17,29. We also observed that the SMARCA2-VHL interface contacts are assembled from unstructured loop regions in both VHL and SMARCA2 (POI), as shown in Fig. 5C and D. In summary, the number of PPI residue contacts in the crystal structures does not correlate with the cooperativity. Next, we investigated PPI properties with a focus on disordered regions, which are often in suboptimal, or ‘frustrated’, energy states^13,30. The SMARCA2-VHL interface is dynamic and is formed by residue contacts built from the unstructured loop regions, as shown in Fig. 5A, C and D. We calculated the frustration of these interface residues to examine whether or not frustration can distinguish between high and low cooperativity PROTACs.

Fig. 5: Interprotein interactions between SMARCA2 and VHL in ternary complexes.

A Heatmap of residue contacts made between SMARCA2 and VHL domains, from both crystal structures and MD simulations (contact present if formed for >50% simulation time). Red “+” signs indicate residue contacts observed in the crystal structures. Residue numbers for VHL are below 1000, and for SMARCA2 are above 1000. Residues located in the unstructured loop regions are labeled in magenta. B The total number of contacts formed between residues in SMARCA2 and VHL in the crystal structure and MD simulations. The residue pairs with more than 50% contact frequency in P1 (C) and P5 (D) are shown as yellow dashed lines. SMARCA2 is colored in blue, VHL is colored in grey, and the loop regions of the SMARCA2-VHL interface are colored in magenta.

High cooperativity PROTACs show a greater number of frustrated residue contacts in the protein-protein interface of the ternary complex

As described in the previous section, part of the SMARCA2-VHL interface interactions involves disordered loops in SMARCA2, ranging from residues N1396 to F1431 and from L1465 to S1468. Residues in unstructured loops contribute to the SMARCA2-VHL interface contacts during MD simulations (Fig. 5A). As discussed by Marton et al.³¹, disordered regions sample a wider range of conformations and thus play an important role in protein function. Chen et al.¹² have shown that functionally critical structural regions in a protein often have patches of ‘frustrated’ residues that adopt suboptimal energetic states. Based on this concept, Gonzalo-Parra et al. developed the ‘Frustratometer’ toolkit³² to calculate the frustration of residue contacts within proteins or between interacting protein domains. The Frustratometer toolkit is based on the coarse-grained forcefield AWESEM³³. Briefly, in the Frustratometer procedure, each amino acid pair is mutated in silico to all other amino acids and the energy distribution for all the mutated pairs (decoy pairs) is calculated with respect to the energy of the wild type pair. The wildtype residue pair energy is then converted into a z-score (${{zscore}}_{{wildtype}}=\frac{{E}_{{wildtype}}-\mu }{\sigma }$, where μ and σ are the mean and standard deviation of the decoy energy distribution, respectively). The z-score is a measure of the wildtype contact energy relative to all other amino acid substitutions at these residue positions, and hence an estimation as to whether a residue is in an energetic optimal or suboptimal state. Each residue pair was classified as highly frustrated, neutral, or minimally frustrated based on the criteria employed in Frustratometer V2 (a z-score of mutational frustration > 0.78 is minimally frustrated, < −1 is highly frustrated and between 0.78 and −1 is neutral). More details on the calculation of frustration are given in the Methods section (see also Supplementary Data 1).

For each PROTAC ternary complex shown in Fig. 2, we used the Frustratometer³² on the MD trajectories and calculated the number of frustrated residue pairs in the SMARCA2-VHL interface for each MD snapshot. As shown in Fig. 6A, the ratio of the number of highly frustrated pairs over the total number of residue pairs in the interface for each ternary complex correlates with the measured cooperativity. Since frustrated residue pairs are often located at active sites and protein interfaces¹², we speculate that the frustrated interface residues could show a high-energy transient state that is indicative of catalytically active conformations and hence PROTAC cooperativity.

Fig. 6: Frustration in ternary complexes.

A The ratio of the average number of highly frustrated residue pairs in the SMARCA2-VHL interface to the total number of residue pairs in the PPI, plotted against the measured cooperativity. PROTACs with cooperativity greater than 10 are colored in green, between 5 and 10 are colored as orange, and below 5 are colored red. B The highly frustrated residue pairs from the representative structures from the most populated cluster of PROTAC P1 and the representatives extracted from the top 2 populated conformational clusters of P5 are shown as cyan dashed lines. SMARCA2 is colored in blue, VHL in white. The unstructured loop region of the SMARCA2-VHL interface is colored in magenta. C Heatmap of frustration level for residue pairs located in the SMARCA2-VHL interface. The frustration level was calculated by averaging the frustration scores of all frames from the MD simulation. Residues located in the unstructured loop region are labeled in magenta.

Most frustrated residues are in the loop regions of either VHL or SMARCA2 as shown in Fig. 6B for high (P1) and low (P5) cooperativity PROTACs. This is also shown in the frustration heat map for the SMARCA2-VHL interface residue pairs (Fig. 6C). Supplementary Movies 3 and 4 showing the three-dimensional view of this region are given in the Supplementary Information. Using double mutant cycles, Jemth et al.³⁴ have shown that the interactions in a protein complex involving unstructured domains, especially in promiscuous protein complexes, are likely to be energetically suboptimal and exhibit some strain due to unfavorable interactions. This is what we observe here in the PROTAC ternary complexes. In most MD simulation runs, the number of frustrated residues in the SMARCA2 protein, is higher than in VHL (Supplementary Fig. 6). Interestingly, Proline is present in 23% of highly frustrated pairs in all five simulations. Both glutamine and glutamic acid are present in 21% of highly frustrated pairs, and asparagine in 11% of the frustrated residue pairs in the interface. These frequently occurring residues imply that the chemical nature of amino acids may also play a role in regulating the frustration level of PPIs, and it has been shown in literature that certain residues may be more critical in forming a catalytically effective interface^31,35. In summary, we show that, for the SMARCA2/VHL system under study, the frustrated residue pairs in the PPI of PROTAC ternary complexes can successfully distinguish high from low cooperativity PROTACs, and we will use the same high/low terminology in future discussions. It is the focus of future research to establish whether this can be utilized for all POI-E3 pairings in which cooperativity is a relevant optimization parameter.

Frustration involving the POI-E3 interface residues can be used to assess the cooperativity of PROTACs prior to synthesis

One goal of this work was to assess the utility of the frustration in the PPI for predicting the cooperativity of PROTACs, and ideally to assist synthesis prioritization. As a further validation of the frustration measure, we studied 11 GEN-1-based PROTACs (shown in Supplementary Fig. 1) that have no crystal structures except P13 (pdb ID: 6HAY). Using the methods detailed in the Methods section, we generated ternary complex structural models using a now publicly available PROTAC (P2) ternary complex template structure (pdb ID: 6HAX)⁸. We speculated that TC conformations, present in high-cooperativity PROTACs such as the one in 6HAX structure, would lead to more pronounced differences in the number of frustrated residue pairs. The 11 PROTACs feature the same VHL binder as the PROTACs P1, P2, P4, and P5, except for P19 and P20, which contain a non-binding distomer of VH101. The SMARCA2 binder in these PROTACs, GEN-1, binds into the same pocket and provides a similar exit vector²¹. Starting from the structural models of the 11 PROTAC TCs, we performed MD simulations using the same protocol as applied to the crystal structures (see Methods). Using the MD simulation trajectories for each TC model, we calculated the average number of frustrated residue pairs located at the interface between SMARCA2 and VHL. To strengthen the error estimation, we increased the MD simulation sample size by generating 10 different conformation ensembles using bootstrapping methods (see Methods for more details). As shown in Fig. 7, the average number of frustrated residue pairs shows the same trend as the cooperativity, with the higher cooperativity PROTACs being more likely to show higher frustration than weak cooperativity PROTACs.

Fig. 7: Frustration at the protein–protein interface correlates with cooperativity across modeled PROTAC ternary complexes.

Average number of highly frustrated residue pairs in the PPI of SMARCA2–VHL calculated from MD simulations of docked models for 11 GEN-1–based PROTACs, shown alongside their experimentally measured cooperativity (dots). Each system was simulated with n = 3 independent MD replicates, each 100 ns long, starting from randomized initial velocities. For each replicate, conformations with a SMARCA2–VHL interface RMSD < 1.5 Å from the starting structure were pooled. To estimate variability in the number of highly frustrated residue pairs, we applied a bootstrapping approach: 100 frames were randomly sampled 10 times from each replicate ensemble, and the average and standard deviation across these 10 trials were calculated. The standard deviations shown on the X-axis represent this bootstrap-based estimate. PROTACs with a cooperativity greater than 10 are colored blue, and those below 10 are colored orange (α refer to cooperativity).

To test the robustness of the frustration calculations, we used a bootstrapping approach. For each PROTAC, we generated 10 random conformation ensembles from the aggregated MD trajectories. The number of highly frustrated pairs was calculated from each of the 10 samples, and the average and standard deviation were estimated to assess the level of robustness in our sampling protocol. This analysis revealed that, while some variation in frustration is expected when selecting random frames, the overall trend in the number of frustrated pairs remains consistent and effectively differentiates between strong and weak PROTACs, as reflected in the error bars on the x-axis in Fig. 7. The trend from the docked models concurs with the frustration pattern observed among the high, medium and low cooperativity PROTACs in the MD simulations starting from the crystal structures. This indicates that frustration can be used to predict cooperativity and hence to prioritize PROTACs for synthesis and experimental testing. We calculated the cooperativity using MMGBSA calculations with the same ensemble as was used for frustration calculations for the 11 PROTAC ternary complexes (see Methods for details). The results in Supplementary Fig. 7 show no correlation to the measured cooperativity. In summary, frustration can be used as a measure to identify strong PROTACs.

Discussion

Despite some notable successes, the design of PROTACs can be challenging, since conventional structure-based drug design principles using single crystal structures do not represent the highly dynamic nature of some ternary complexes formed with PROTACs. PROTAC properties are governed by the conformational ensemble and the dynamic interacting interfaces in the ternary complexes. While the calculated cooperativity of PROTACs showed correlation with the measured cooperativity α in some systems, this hypothesis was not supported by our data in this study (Supplementary Figs. 4 and 7). Therefore, we used residue pair frustration to characterize the dynamic nature of the PPI in PROTAC ternary complexes. Our major findings based on studying the 5 crystal structures of VHL-SMARCA2 ^BD PROTAC system, four of which were solved in this study, and 11 homology modeled PROTAC systems are:

• 1. Pairwise residue frustration in the protein-protein interface of PROTAC ternary complexes (independent of the starting structure - homology model or crystal structure) can be used as a metric to characterize the PPI in TCs.

• 2. PROTACs with strong cooperativity show a greater number of highly frustrated residues in the interface compared to PROTACS with weak cooperativity. PPIs can therefore be characterized by frustration as an in silico descriptor and used to identify potent PROTACs prior to synthesis, given that the cooperativity calculated from MD simulations does not always correlate with the measured cooperativity.

• 3. The highly frustrated residue pairs are in the disordered regions in the VHL-BRD4 ternary complexes.

Comparison of the five SMARCA2::PROTAC::VHL TC crystal structures reveals large domain motions. We derived the free energy surfaces of these from MD simulations and identified multiple free energy minima, representing distinct TC conformations. We also observed that the TCs formed with low-cooperativity PROTACs are more flexible than those with high-cooperativity PROTACs. The spread of the relative orientations in the SMARCA2 domain is wider in the MD simulations than observed when comparing crystal structures only. The orientation of the E3 ligase domains relative to the POI has been shown to play a critical role in substrate selectivity by PROTACs³⁶. Such dynamic behavior of TCs with other POIs or E3 ligases, such as cereblon, has been shown by other MD simulation studies^17,37. We observed that the number of residue contacts formed during MD simulations in the protein-protein interacting interface is larger in high-cooperativity PROTACs when compared to low-cooperativity PROTACs (Supplementary Fig. 6 and Fig. 5B). Such differences have also been seen in cereblon-based PROTACs paired with other POIs²⁸.

Four new crystal structures, along with the MD simulations, also showed that, for different PROTACs, the SMARCA2-VHL interface is formed by residues that are mostly located in the unstructured regions. Frustration in protein-protein interfaces arises when residues at the interface cannot adopt optimal conformations or achieve dense packing due to conflicting interactions³⁰. While it is generally understood that frustration decreases the stability of a protein, the relationship between frustration and the functional activity of a protein remains vague. It is known that in some enzymes, frustration in the active site can promote catalysis by destabilizing the reactants, products, and facilitating the transition state³⁸. Here, we observe that the number of highly frustrated residue pairs in the SMARCA2-VHL interface of the ternary complex is higher for high-cooperativity than low-cooperativity PROTACs. The frustrated residue pairs have a higher occurrence of prolines, glutamine and glutamic acid residues compared to other amino acid residues that are also more frequently found in intrinsically disordered regions of proteins³⁵. Therefore, we hypothesize that frustration in the POI-E3 interface can be utilized for characterizing the nature of the TCs and for describing cooperativity. We have provided initial evidence that the number of frustrated residue pairs can be used as a measure to distinguish high from low cooperativity PROTACs, based on the TC models that were generated by comparative modeling in combination with MD simulations. In the future, we aim to use the specific frustrated residue pairs in the PPI of TCs for the rational design of PROTACs.

Recently, we published a study on frustration in BRD4-PROTAC-cereblon ternary complexes, where we observed similar trends³⁹. Five PROTAC-mediated ternary complexes with known degradation efficiencies were simulated, and the results revealed that higher degradation efficiency correlated with reduced POI-E3 domain motion and increased frustration in the PPI. A key difference between the cereblon and the VHL systems is that the BRD4-PROTAC-cereblon ternary complex features a distinct structured hydrophobic patch in the PPI interface. The frustration of this hydrophobic patch is regulated by the PROTAC linker, ultimately influencing degradation efficiency. In contrast, no such hydrophobic patch is observed in the SMARCA2-PROTAC-VHL complex. Despite the differences in the E3 ligases (VHL vs. cereblon) and the focus on cooperativity versus degradation efficiency, our findings highlight frustration as a valuable computational tool for evaluating PROTAC function across different systems³⁹. How is frustration related to degradation? The ubiquitination process involves multiple steps after the formation of the PROTAC ternary complexes and depends on multiple complex factors. Previous studies using frustration have shown that highly frustrated residues are often located in the active site of proteins^12,40. The limited speculation we can make from this study is that the high number of frustrated residue pairs in the strong PROTAC interface could be indicative of efficient catalytic activity of ubiquitin transfer, leading to degradation. This would be the subject of a future study.

Methods

MD simulation protocol

The VHL::P2::SMARCA2^BD ternary complex crystal structure has previously been published and is available in the protein databank (PDB ID: 7Z77)⁹. Each PROTAC ternary complex was initially prepared using the Maestro Protein Preparation Wizard module⁴¹ (Schrödinger Release 2023-1: Maestro, Schrödinger, LLC, New York, NY, 2021). Histidine protonation states were determined under neutral pH conditions. PROTAC forcefield parameters were derived using the Antechamber module from AmberTools²⁵, while the partial charges were calculated using the Jaguar module in Maestro⁴². The protein atoms were parameterized using the FF99SB-ILDN force field⁴³. The PROTAC TC complex was solvated in a truncated octahedron SPC⁴⁴ water box with at least 14 Å margin from the protein surface using tleap from the AmberTools 18 package²⁵. The simulation box was neutralized using an adequate number of Na⁺/Cl^- ions. The simulation box went through a two-step minimization procedure. In the first step, a 2000-step steepest descent minimization was conducted, with a 10 kcal/mol constraint applied to all atoms, excluding water. During the second minimization phase, all atoms were permitted to move, and an additional 2000-step steepest descent minimization was carried out. Following minimization, an equilibration simulation was executed. The initial equilibration step involved a brief 0.5 ns NVT heating phase, raising the temperature from 0 K to 300 K with a 10 kcal/mol restraint applied to the protein and PROTAC’s heavy atoms. The system underwent a seven-step NPT equilibration to gradually release restraints at 1 atmosphere pressure and 300 K temperature. Restraints were gradually reduced from 10 kcal/mol to 5 kcal/mol, then 4, 3, 2, 1, with each step lasting for 5 ns. Finally, a 10 ns equilibration without restraints was conducted at the end. The final snapshot of the equilibration phase was used as the initial conformation for production runs. Six random velocities were employed to generate six independent production runs. Each production run lasted for 500 ns, and the complete production trajectories (30,000 snapshots in total) were used for subsequent analysis. Throughout the simulations, the temperature and pressure were controlled by the Nosé-Hoover method^45,46. The particle mesh Ewald method was utilized to compute the electrostatic interactions⁴⁷, while the non-bound interactions were determined using a 12 Å cutoff.

Principal component analysis and Shannon entropy calculation

We conducted principal component analysis (PCA) on the combined production trajectories, considering all Cα atom positions for the RMSD calculation in the PCA, and projected the resulting principal components (PCs). The first two principal components constituted approximately 80% of the total eigenvector variance (Supplementary Fig. 2). The gmx anaeig module⁴⁸ was employed to examine the protein motions represented by the first two principal components. The two TC conformations with maximal PC1 and PC2 values were aligned onto the VHL domain to assess the conformational changes in the SMARCA2 domain (Fig. 3A). The combined sampling points were projected onto the first two PC axes and displayed as a contour map (Supplementary Fig. 3). The individual contour map for each PROTAC is displayed in Fig. 3B. The combined structural ensemble of all five PROTACs was clustered in the top two PC dimensions using the agglomerative clustering algorithm, generating seven clusters (yielded the best Silhouette Coefficient). The conformations corresponding to the centroids of all seven clusters were aligned by their VHL domains to compare the orientations of the SMARCA2 domains (Fig. 3C). For each PROTAC, we calculated the Shannon entropy as -$\sum {P}_{i}\mathrm{ln}({P}_{i})$, where P_i denotes the fraction of conformations that fall into cluster i (Fig. 3D).

Calculation of PROTAC cooperativity from MD simulation trajectories

We calculated the cooperativity using MMGBSA calculations in AMBER (MMGBSA.py) using the entire trajectories of the 5 crystal structures and the homology models for the 11 PROTAC ternary complexes. These calculations include desolvation penalties. The cooperativity was calculated using the following formula:

$$ {Theoretical\; Cooperativity} \\ =\exp \left(\frac{\left(-{G}_{{ternary}}+{G}_{{target}-{PROTAC}}+{G}_{E3-{PROTAC}}-{G}_{{PROTAC}}\right)}{{RT}}\right)$$

(1)

Where, G_ternary is the energy of the ternary complex, G_{target-PROTAC} is the energy of the binary complex of BRD4 + PROTAC, G_E3-PROTAC is the energy of the binary complex, VHL + PROTAC. Supplementary Figs. 4 and 7 show the results of the calculated cooperativity for the 5 crystal structures and the 11 homology models in the Supplementary Information.

RMSD based clustering and contact analysis

Each TC ensemble was clustered based on the RMSD among the Cα atom positions using the program CPPTRAJ⁴⁹. We selected a cutoff between 1.5 to 2.5 Å such that the top three clusters accounted for approximately 60% of the overall population. The representative structures of the top three populated clusters were extracted and aligned onto their VHL domain to display the conformational differences in the SMARCA2 domain (Fig. 4A and B). The RMSD values among the top three representative structures were calculated as follows: P1: 0.5~0.8 Å, P2: 0.5~0.8 Å, P3: 0.5~0.8 Å, P4: 1.4~4.5 Å, P5: 0.4~2.7 Å. We computed the contact frequencies between the proteins and the PROTAC, as well as between SMARCA2 and VHL, for both the crystal structures and the simulation trajectories. The contact analysis was carried out using the program get_contact (https://getcontacts.github.io/), considering every atom and all types of contacts in the interface^16,34,35.

Frustration calculation using the MD conformations (Fig. 6A)

We conducted frustration calculations on the trajectory of every TC, using Frustratometer³². We initially identified the SMARCA2-VHL interface residue pairs from every snapshot by locating all pairs with a minimum distance of less than 4.5 Å between heavy atoms. We then calculated the frustration levels for these identified pairs in each MD snapshot. Each pair was classified as highly frustrated, neutral, or minimally frustrated based on the criteria employed in Frustratometer V2 (a z-score of mutational frustration > 0.78 is minimally frustrated, <−1 is highly frustrated). Then we counted the number of residue pairs that were highly/neutral/minimally frustrated in over 40% of the cluster population. To compare the level of frustration among different PROTAC simulations, we calculated the ratio of highly frustrated pairs over the total number of frustrated pairs (including highly/neutral/minimally).

Comparative modeling, simulation, and analysis of PROTAC ternary complexes without crystal structures (Fig. 7)

To investigate the 11 GEN-1-based PROTACs, we utilized the publicly available SMARCA2^BD: PROTAC 2: VHL crystal structure (PDB ID: 6HAX)⁸ as a template to construct ternary complex models, as 6HAX displays a similar SMARCA2-VHL conformation to that found in the P1 and P2 TCs. This was one of the two crystal structures of SMARCA2-VHL complexes available at the start of this project. The purpose of generating these ternary structural models was not to predict the complex structure that each of the 11 PROTACs would form. Instead, our goal was to place each of the 11 PROTACs in a known high cooperativity PROTAC ternary model conformation and perform short MD simulations to calculate the frustration in the protein-protein interface.

We constructed the 3D coordinates of each of the 11 PROTACs based on the P2 in the 6HAX conformation as reference using Maestro, and these were subjected to geometry optimization using Hartree-Fock (HF) quantum mechanical calculations and the STO-3G basis set as implemented in the Jaguar module of Maestro⁴². Then these structures were aligned with PROTAC 2 in the 6HAX ternary complex formed by SMARCA2 and VHL using the pair-fit utility in PyMOL 3.0^50,51,52. Subsequently, the local clashes between the protein side chain atoms and the newly inserted PROTACs were resolved using the side chain replacement feature in Maestro’s PRIME module^53,54. Partial charges calculated from the HF calculations were used in the molecular topology of PROTACs created by Antechamber for MD simulations without any modifications.

Each TC was prepared and simulated using the same protocol as described in the previous sections. For each system, we conducted three 100 ns-long simulations starting with random velocities. After the production runs, we performed RMSD clustering and extracted frames with a SMARCA2-VHL interface RMSD below 1.5 Å with the initial structure, for a frustration calculation using Frustratometer. This step is to ensure that the frustration metric can be compared across all 11 PROTACs, sharing a consistent ensemble of PPIs.

To account for sampling variability in the frustration calculations, we implemented a bootstrapping approach. For each PROTAC, from the subset of conformations with RMSD < 1.5 Å, we randomly selected 100 frames, repeating this process 10 times. The number of highly frustrated pairs was calculated from each of the 10 trials, and the average and standard deviation were estimated to assess the level of robustness in our sampling protocol.

Chemical synthesis

A list of final compounds is compiled in Supplementary Table 2 and Supplementary Fig. 4. Full details of synthetic procedures for P1, P3, P4, P9, P10, P11, and P18 (including Supplementary Figs. 9 to 34) with NMR and HPLC spectra of final compounds are provided as Supplementary Methods. The syntheses of P2 and P5 have been described in ref. ⁹. (compounds 6 and 23 therein respectively), the syntheses of P6, P7, and P16 have been described in ref. ⁵⁵. the syntheses of P8, P13, P19, and P20 have been described in ref. ⁸ (ACBI1, PROTAC 1, cis-ACBI1 and cis-PROTAC 2 therein respectively).

Protein production and crystallography

Protein production for SMARCA2^BD and the VCB complex was done as previously described^8,9 and experimental details are provided in the Supplementary Methods. For P1 and P3 ternary complex crystallization, SMARCA2^BD and VCB complex were incubated with 1 mM of the ligand in a 1:1:1 ratio overnight and crystallized in 96-well sitting drop plates. P1 was crystallized with a reservoir solution of 0.1 M K₂HPO₄, 16% PEG 8000, and 0.2 M sodium chloride. For P3, the reservoir solution was 0.1 M Bis/Tris propan,e pH 7.5, 32% PEG 335,0 and 0.1 M tri-sodium citrate. Datasets were measured at SLS beamline X10SA and processed with autoPROC⁵⁶. Resolution cutoffs were determined using the default parameters of autoPROC. The structures were solved by molecular replacement using Phaser⁵⁷ with PDB: 7Z76⁹ as a search model. The final model was prepared in iterative cycles of model building with COOT⁵⁸ and refinement with autoBUSTER (Global Phasing Ltd). Data collection and refinement statistics are provided in Supplementary Table 3 and Supplementary Table 4, respectively. For the P4 ternary complex: SMARCA2^BD, VCB complex and P4 were mixed in a 1:1:1 molar ratio, concentrated to 10 mg/mL in a buffer containing 0.1 M HEPES pH 7.5, 100 mM sodium chloride, 1 mM TCEP, 2% DMSO, and incubated at room temperature. Drops were prepared by mixing 1 μL of the ternary complex with 1 μL of the well solution and crystallized at 18 °C using the hanging-drop vapor diffusion method. Crystals were obtained in 0.1 M HEPES pH 7.0, 8% ethylene glycol, 14% PEG 8 K. Harvested crystals were flash-cooled in liquid nitrogen following gradual equilibration into a cryoprotectant solution consisting of 20% PEG 8 K, 0.1 M HEPES pH 7.0, 20% ethylene glycol. Diffraction data were collected at Diamond Light Source beamline I24 (λ = 0.9686 Ã…) using a Pilatus3 6 M detector and processed using XDS. The crystals belonged to space group P21 with unit cell parameters a = 48.7, b = 89.8, c = 64.2 Ã… and α = 90.0°, β = 96.9°, γ = 90.0° and contained one copy of the ternary complex per asymmetric unit. The structure was solved by molecular replacement using PHASER with VCB coordinates derived from the VCB: MZ1: Brd4^BD2 complex (PDB: 5T35) and SMARCA2^BD (PDB: 4QY4) as search models. Subsequent iterative model building and refinement were done according to standard protocols using CCP4, COOT, and autoBUSTER (Global Phasing Ltd). For the VCB:P5:SMARCA2^BD ternary complex: VCB, P5, and SMARCA2^BD were mixed in a 1:1:1 stoichiometric ratio in 20 mM HEPES, pH 7.5, 100 mM sodium chloride, 1 mM TCEP, 2% DMSO, incubated at room temperature and concentrated to a final concentration of approximately 16 mg/ml. Drops were prepared by mixing 1 μL of the ternary complex with 1 μL of well solution and crystallized at 18 °C using the hanging-drop vapor diffusion method. Crystals were obtained in 28% PEG 3350, 0.2 M potassium thiocyanate, 0.1 M Bis-Tris Propane, pH 7.5. Harvested crystals were flash-cooled in liquid nitrogen following gradual equilibration into cryoprotectant solution consisting of 40% PEG 3350, 0.2 M potassium thiocyanate, 0.1 M Bis-Tris Propane, pH 7.5. Diffraction data were collected at Diamond Light Source beamline I24 (λ = 0.9686 Ã…) using a Pilatus3 6 M detector and processed using XDS. The crystals belonged to space group P21 with unit cell parameters a = 46.6, b = 81.6, c = 58.1 Ã… and α = 90°, β = 97.8°, γ = 90.0° and contained one copy of the ternary complex per asymmetric unit. The structure was solved by molecular replacement using PHASER with VCB coordinates derived from the VCB: MZ1: Brd4^BD2 complex (PDB: 5T35) and SMARCA2^BD (PDB: 4QY4) as search models. Subsequent iterative model building and refinement were done according to standard protocols using CCP4, COOT, and autoBUSTER (Global Phasing Ltd).

Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.