Topological dynamics of an intrinsically disordered N‐terminal domain of the human androgen receptor

2.1. Modeling of three‐dimensional structure of AR‐NTD

Due to its disordered conformation, there is no resolved structure of the NTD domain deposited on the Protein Data Bank (PDB). To start from a computationally efficient initial structure, I‐TASSER server ³⁸ was used for modeling the three‐dimensional structure of the AR‐NTD (residues 1–538). I‐TASSER is a hierarchical approach to protein structure prediction and structure‐based function annotation. It is ranked the best protein structure prediction method by the Critical Assessment of Protein Structure Prediction (CASP) community. The Z‐score and normalized Z‐score >1 indicates a template with a good alignment quality. The Confidence‐score (C‐score, values −5 to 2) reflects the threading alignment and convergence of structure refinement simulations with higher score corresponding to a high confidence model. A Template Modeling (TM) score >0.5 indicates correct conformation of the model, while scores <0.17 corresponds to a model of random topology. Therefore, models with the highest Z‐score, C‐score and TM‐score were selected for our MD studies.

In order to further optimize and validate the initial structure, an energy minimization step using steepest descent method was performed followed by a conjugate gradient with a ff99SB all atom force field to perform a total of 100,000 steps using the GROMACS software package. ³⁹ For validation of energy refined structures, WHATIF ⁴⁰ was used for anomalous bond angle and bond length detection. RAMPAGE ⁴¹ was used for Ramachandran plot analysis and detection of incorrectly modeled protein regions was done using ERRAT. ³⁰ The compatibility of the three‐dimensional protein model to amino acid sequence was verified using VERIFY3D ⁴² and the global quality of the models was calculated using QMEAN. ⁴³ Structural disorder was analyzed using the PONDR web server and plots directly obtained from the server. ⁴⁴

2.2. Molecular dynamics setup

In this study, 5 μs of MD simulations was carried out using the SIRAH coarse‐grained force field for proteins in combination with the WT4 explicit coarse‐grained water model. ⁴⁵ The advantage of using this model is to reproduce the roughly tetrahedral ordering of bulk water due to the existence of hydrogen bonds between atomistic water molecules. Three MD simulation with different random seeds started from energy minimized predicted structure obtained from structure prediction pipeline. The entire AR‐NTD (residues 1–538) was mapped to a coarse‐grained representation according to the standard SIRAH mapping. Rhombic dodecahedron box was used to dissolve NTD structure by adding WT4 water molecules. Electroneutrality and physiological concentration of salt were achieved by replacing corresponding amount of water molecules with NaW and ClW (coarse‐grained representations of Na+ and Cl− ions, respectively). All systems were minimized using the steepest descent algorithm and then went through a 5 ns NVT equilibration, a 5 ns NPT equilibration, and a NPT production run. The leapfrog integrator with a 20 fs time step was used throughout. Protein beads were constrained with the LINCS algorithm ⁴⁶ during the equilibration, and no constraints were used during the minimization and production. The temperature was kept at 310 K with a velocity rescale thermostat, ⁴⁷ and the pressure at 1 bar with the Parrinello−Rahman barostat. ⁴⁸ τ_T for the thermostat was set to 1.0 ps during the equilibration phases and to 2.0 ps during the production. τ_P for the barostat was set to 10.0 ps during both the NPT equilibration and the production. Both nonbonded cut‐offs (van der Waals and shortrange electrostatics) were set to 1.2 nm. Long‐range electrostatics were treated with the PME method with a 0.2 nm grid spacing during the equilibration and 0.25 nm during the production. Nonbonded interactions were calculated using a 1.2 nm cut‐off neighbor list, updated every 25 steps (in the production and the NPT equilibration) or 10 steps (in the NVT equilibration). Both energy and pressure dispersion corrections were applied. Periodic boundary conditions and the minimum image convention were used. Snapshots were collected every 1,000 steps (20 ps). All simulations and analyses were carried out with the GROMACS 2020 software. ³⁹

2.3. Molecular docking

The binding interactions of NTD with LBD and DBD was studied using the ClusPro web server ⁴⁹ with the aim of identifying interaction sites. This docking suite has been consistently rated among the best global docking procedures in the CAPRI challenge (Critical Assessment of Predicted Interaction). ⁴⁹ ClusPro first performs rigid body docking by sampling a wide range of conformations and subsequently clusters the thousand lowest‐energy structures based on root‐mean‐square deviation (RMSD) to identify the largest clusters that represent the most likely models. It furthermore refines the represented structures by energy minimizing the CHARMM potential. ⁴⁹ We supplied ClusPro with previously reported crystal structures for DBD (PDB code: 1R4I), LBD (PDB code: 1T7T), and RAP74‐CTD (PDB code: 1I27) to dock with our energy minimized structure of the NTD. The top 10 best ranked models in ClusPro were selected for further analysis.

In order to find the possible binding site of the EPI‐001 compound on the CR region of the NTD, drug–protein interaction analysis was performed. The mol file of the EPI‐001 was created in online molecular drawing module and then transferred to PRODRG server ⁵⁰ to generate the required topologies and structure coordinates. The representative structures of three clusters obtained from clustering the CR region of NTD were prepared for docking by adding Kollman charges and hydrogen atoms to the polar groups of the protein using AutoDock tools. ⁵¹ AutoDock 4.2 package ⁵¹ was used for docking the EPI‐001 to the target molecules. Blind docking was performed for the drug using the default parameters on Autodock.

2.4. Preparing the structures for circuit topology analysis

After the trajectory of the system was obtained, atomic positions of amino acids were generated from the location of CG beads. Backmapping was done using sirah_vmdtk.tcl plugin followed by 100 steps of steepest‐descent and 50 steps of Conjugated Gradient minimization in vacuum using the sander module of AmberTools. ⁴² This process made the procedure robust and independent from the fine details of the backmapping library. The obtained atomistic coordinates were used for circuit topology analysis.

2.5. Circuit topology analysis

Circuit topology parameters were retrieved using custom‐made Python code, modified for energy, length filtering, and circuit decomposition. For doing this, 1,000 frames were extracted from the trajectory containing only the protein molecule. Contacts between residues were calculated based on a distance cut‐off of 4.5 Å. Residue pairs where more than five atoms were found to be at a respective distance below the distance cut‐off were computed as contacts. The three closest neighbors of each residue were excluded from analysis.

In order to perform dynamic CT analysis, frames from the AR‐NTD simulations were taken 5 ns apart and transformed into binary contact maps. The stack of contact maps was processed to create the contact life‐time kymographs: all combinations of residues which never participate in contacts were removed. In order to create the contact lifetime distribution, the maximum lifetime was chosen for each particular combination of residues.

Time filtering for topological study was performed by creating life‐time filters matrices. These matrices were then multiplied to regular contact matrices as a customization of the CT Python tool. All dynamic CT analysis was performed using a custom‐made Jupyter lab data analysis file.

Accession numbers:

PDB: 1R4I, 1T7T, 1I27.

3.1. AR‐NTD folds into two disjoint highly dynamic regions

First, we performed molecular dynamics simulations of the AR‐NTD in an aqueous solution with physiological salt concentration. The initial structure, shown in Figure 1a, was built using I‐TASSER server. ³⁸ The model was chosen as the best ranked model obtained from the server. The model was superior to conformations predicted by the AlphaFold and RosettaFold severs, based on confidence measures and comparative analysis with experimental data (Figure S2). Comparison between α‐helix content of the models predicted by different algorithms with Circular Dichroism analysis of the NTD ¹⁸ revealed that among the predicted models, the best ranked model from I‐TASSER and the AlphaFold model showed closer values to CD data. However, given the very low confidence level of the AlphaFold model, the best ranked model obtained from I‐TASSER was chosen for performing MD simulations.

After minimization and relaxation, we performed MD simulations of the NTD structure. Visual examination of the trajectory and RMSD plots show that the initial conformation has undergone an extensive structural change (Figures S3 a and 1b). Figure 1c provides a kymograph of all the secondary structures within the NTD, which is mostly composed of random coils, turns, and loops despite the apparent changes in the global 3D shape of the system (Figure S3 b and c). We repeated the MD simulation three times using different initial velocities to ensure that our observations are not dependent on the initial conditions. Importantly, all three independent runs consistently led to the emergence of two disjoint regions in 3D within 2 μs of simulations: an extended N‐terminal sub‐region (NR, residues 1–224), and a C‐terminal sub‐region (CR, residues 225–538). Both sub‐regions showed high level disorder and high solvent accessibility, which reached a plateau after 2 μs of simulations (Figure 1d,e). We monitored the dynamics for an additional 3 μs and computed root‐mean‐square fluctuations (RMSF) to quantify the fluctuations of the chain. We observed an order of magnitude larger RMSF for AR‐NTD in comparison to the folded AR‐LBD (Figure 1f, dashed line at 0.15 nm). Considering the average fluctuations of the chain we observed lower values of RMSF in the vicinity of residue 400. Quantification of solvent accessibility revealed that both CR and NR are highly solvent accessible. The CR however has a central segment (residues 355–469) which is relatively less solvent accessible (Figures 1e and S4), corresponding to the reduction seen in the RMSF profile. Disorder prediction data produced by PONDR analysis agrees with the solvent accessibility and RMSF profiles, with the central region within CR having less disorder than the rest of the CR (Figure S5). We denote this segment as CR core, to differentiate it from the rest of the CR, CR shell. Interestingly the core includes the well‐known transcriptional activation unit 5 (Tau‐5) motif of the AR‐NTD. Comparing the radii of gyration of CR and NR, one can clearly see that the CR is significantly more compact than NR (Figure 1g). This is consistent with our model that CR, despite its disorder, carries residual structures in its core.

To further demonstrate that the formation of disjoint regions does not depend on the choice of initial configurations, we melted the compact NTD at 350 K and then ran the simulation for 1 μs at physiological temperature. Melting opened up the compact structures and led to disappearance of separation between regions. The subsequent simulations, however, recovered the NR/CR in all three experiments and the low solvent accessibility of the CR core. Moreover, we find that our observation of NTD folding to two disjoint regions does not depend on the F‐parameter of the force field, which encodes the interactions between the polypeptide and the surrounding solvent (Figure S6). We note that compaction of the chain to form the CR/NR segregation and the emergence of the CR core are AR‐NTD specific, as synthetic models of AR‐NTD using large poly‐glutamine (polyQ) or poly‐alanine (polyA) stretches, generated by mutating all residues of the structure shown in Figure 1a to Q or A respectively, did not produce these features (Figure S7, S8).

It is noteworthy that we have used these two models as control analysis to better understand the complex behavior of the AR‐NTD (particularly during the first 1 μs). The synthetic polyQ model used to demonstrate the sequence specificity is not to be confused with a natural short polyQ motif that is embedded within the AR‐NTD. Q and A residues are known as disorder promoting and order promoting amino acids and were thus chosen to build these synthetic models with extreme behaviours. Although it has recently been shown that short lengths of the polyQ tract (up to 25 residues) can form a stable alpha helical structure ¹⁴ there is evidence that extended lengths of the polyQ chain show a largely disorderd structure, ⁵² while polyA shows a high propensity to form ordered structures. ⁵³

3.2. Mapping CR folding dynamics to the topology landscape showed a directional conformational evolution

A natural way to express how the NTD peptide chain folds dynamically is to describe its topological folding patterns over time. In contrast to the conventional geometric approaches, topological approaches are specifically designed to quantify and categorize residual structural (shape) properties of constantly deforming polypeptides. Thus, we decided to map the chain dynamics onto the chain topology space, to identify patterns or motifs of the folding process. ⁵⁴ Contacts in the chain are defined based on the spatial proximity of residues. The arrangement of residue–residue contacts provides a first‐order topological representation of a folded chain, which can be highly informative on folding processes. ²⁵ Here, we [46] exploited the CT framework to provide a dynamic topological description of a disordered protein for the first time. In order to map the topological dynamics of the NTD over time, we calculated the percentages of P, S, and X in chain conformations at 5 ns time points with a total of 1,000 frames over 5 μs. In this way, we followed the topological dynamics of the chain in time (Figure S9).

The topological analysis revealed critical differences in the folding process of NR and CR. While a clear topological trajectory is present for CR, the NR topology map is not showing any readily visible pathway. As seen in Figure S9 c, the starting conformation of CR was rich in series (S) but low in parallel (P) and cross (X) configurations. However, over time, while the percentage of P was almost constant, major changes were observed for S and X content, over all three runs. This topological profile is consistent with the observed emergence of the CR core. Ideally, an unfolded or semi‐unfolded chain will present local contacts at first, ⁵⁵ which have little interaction with far away residues along the chain. Non‐interacting contacts belong by definition to the S topological configuration. However, as the folding process develops, the chain will fold onto itself, creating long‐range interactions between residues, with a more complex topological profile. This process is true for any folding process, as exemplified by the topological evolution of our polyQ model of AR‐NTD, see Methods section for details about the polyQ model.

In contrast to the CR topological trajectory, the NR region explores the topological space in what appears to be a random fashion. As expected from such a highly dynamic region, the opportunities of forming and breaking contacts with high frequency leads to a very wide configurational space. This results in a globular pattern in the topological space. The lack of a preferred direction for topological evolution is nicely visualized by the width and overall shape of the pattern. This is in contrast with the narrow topological evolution of CR, which indicates a number of transient states with significantly different configurations: once a new configuration is reached, the set of possible contacts changes. This phenomenon can be visualized on the topological space as a “leap” from one cluster to another, resulting in a narrow trajectory with high directionality. The lower percentage of series relations relative to the CR indicates that the NR evolution is truly dynamic in nature, with long‐range contacts continuously forming and breaking and the absence of sub‐structures. Indeed, contacts internal to a sub‐structure are most likely in series with contacts embedded into another sub‐structure. This motivated us to look more closely into the topological dynamics of NR and search for patterns using our circuit topology approach.

3.3. Dynamic CT analysis reveals temporal regimes for topological development

Dynamic topological analysis of the AR‐NTD shows the presence of three different time regimes for topological evolution, characterized by the time scales of residue–residue contact lifetimes. We envisioned a dynamic CT‐based topological analysis to identify the time scale of topological evolution of the two previously identified NTD sub‐regions, NR and CR. This, it is possible to follow the temporal evolution of every residue–residue contact formed at any moment of the MD simulation by plotting a kymograph of all possible residue‐residue contact combinations (Figure 2a). By looking at the kymograph (and zooming in on regions of interest), we see how some contacts are fairly short‐lived, often breaking and reforming with high frequency. On the other hand, other contacts live for a few microseconds. Therefore, we might have one (or more) lifetimes for each residue–residue contact: by picking the maximum lifetime of each particular residue combination, we can build a lifetime distribution of all contacts (Figure 2b).

As expected by such a dynamic and heterogeneous system, contact lifetimes span across the whole range of the simulation, with the vast majority of them being short lived (lifetime distribution maxima are 35, 18, and 15 ns, respectively for the three runs). The distribution in Figure 2b respects to some extent the well‐known scale free law ($P\left(\tau \right)=A{\tau }^{-\gamma }$), indicating there might not be a preferential duration for contact lifetimes. ⁵⁶ We see, however, that the degree with which the lifetime data match the power law is not the same for all time scales: we identify a “short life” range ($\tau <0.5\mathrm{\mu s}$), where matching is maximum. After 0.5 μs the distribution starts deviating from the scale‐free behavior ($0.5\mathrm{\mu s}<\tau <2\mathrm{\mu s}$, “middle life”). After 2 $\mathrm{\mu s}$, instability prevails, and deviations from the law appear with wide fluctuations, probably due to poor fitting of the distribution caused by the lack of statistics. The scale free distribution might represent a useful ideal model to quantify the behavior of IDPs. The adherence to this law is to a somewhat higher degree present for our polyQ case study (Figure S10). This regime behavior is highly conserved over all three runs, for both NR and CR (Figure S11). We can calculate the C‐factor for each residue, that is, the average contact lifetime among all contacts in which the residue participates (Figures 2c, S12,13). Here, this heterogeneity becomes apparent. If we exclude residue 23–27, where the only alpha helix present within the NTD resides, higher lifetimes are way more prevalent in the CR region, indicating higher stability and compaction with respect to the NR region. This stark difference in behavior between the two regions cannot be replicated in the synthetic polyQ (Figure S10), indicating this asymmetry is specific to AR‐NTD.

The lifetime analysis allows the identification of connectivity properties within the NTD under three different time regimes. This information can be exploited to disentangle the role of fast forming, short‐lived contacts from that of longer‐lived contacts during topological evolution. Next, we applied a lifetime filter to NTD contact maps and repeated the topological analysis. Results for short and middle life NR contacts can be seen in Figure 2d. Short life contacts seem overall to occupy only one topological state, over all the simulation, which is represented as a globule in the triangular S/P/X topological space. While these contacts oscillate around what looks like a wide topological equilibrium state, with no preferred direction for topological evolution, middle life contacts hop from one state to the other, creating narrow trajectories in the topological space. It is interesting to point out how this narrowing of the topological evolution pattern for middle lifetime contacts is not as evident in the synthetic polyQ model. This indicates there is some AR‐specific coding within its secondary structures that leads to a well‐defined configurational change within the middle lifetime. This is contrary to the view of the NTD as having a completely random conformation. The three runs show a clear pattern towards higher fractions of cross relation, indicating an evolution towards higher degrees of entanglement. The topology of long‐life contacts is shown in Figure S14. As the scale τ > 2 $\mathrm{\mu s}$ is close to the duration of the simulation, this scale is subject to heterogeneous behavior caused by kinetic traps, which manifest in the topological space as globules similar to those we have encountered for short lifetime. The same analysis carried out on CR (Figure S15) reveals a similar topological evolution, although the topological make‐up of the two regions is quite different. While NR and CR have similar cross fractions, NR consistently scores higher in parallel, and the CR in series fractions.

In the context of lifetime analysis, one can consider following the evolution of a specific set of residues, and selectively visualize their contact formation dynamics. For example, cysteine residues play an important role in the folding and conformational dynamics of the proteins by forming disulfide bonds and covalently connecting two distinct parts of the chain. ⁵⁷ There are 11 cysteine residues in AR‐NTD, which make up for a total of 55 possible cys‐cys contacts. Interestingly, only five of such possible contacts actually occur across our simulation runs. These contacts are by nature very dynamic and often form and break multiple times during a run; however, it is possible to calculate the cumulative time a certain contact is detected across frames in the simulation, in order to see which configurations are most likely to occur. Among the five outstanding contacts, we found contacts ²³⁹Cys‐Cys²⁶⁶, ²³⁹Cys‐Cys³²⁶ in which Cys239 stand out as a main contributor. This agrees with previous finding that Cys239 plays an important role in oxidative hemostasis of the AR and aggregation behaviour of the chain. ⁵⁸

One might also consider including spatial information encoded in sequence separation, to shed a light on these topological differences. For this purpose, we calculate the average range of contacts for each residue. The distribution P(d) of such spatial range (Figure S16 a) peaks at about 50 residue–residue distances, for both regions and all MD runs, indicating once again the high dynamicity of this domain. The most probable range of interaction in NTD is already larger than what is commonly identified as interaction range threshold defined for long‐range contacts in proteins (typically 24 residues). Thus, there is a need for defining a relevant threshold for contact range. Based on these distributions, we identified a threshold for long‐range contacts ($\lambda \ge <d>+{\sigma }_d$ ~ 82 residues). The NR region consistently shows a higher density of such long‐range contacts (Figure S16 b). This might explain the different topological profiles observed for the two regions: short range contacts tend to be in series with each other, as the contact sites are too close together to form entangled loops. Moreover, NR long‐range contacts do not display series fraction at all (Figure S16 c), which might indicate the propensity of the region to consistently interact with itself along its full length (i.e., without creating self‐interacting sub‐domains). Furthermore, middle life long‐range contacts do not consistently show the trajectories in the topological space, which we have identified in Figure 2d. This excludes long‐range contacts as main drivers for conformational change, which therefore happens for smaller spatial scales.

Finally, we note that topological arrangements of contacts have implications for overall shape of the NTD domain. We sampled P‐rich, S‐rich, and X‐rich topologies from topology space of middle‐life contacts from last 500 ns of each simulation. Aligning the structures, we found highly dynamic region composed of residues 50–140 with an obvious structural shift in X‐rich conformations (Figure S17). This also was approved by our RMSD analysis which showed a larger deviation for X‐rich conformations in comparison to P‐ and S‐rich ones. The most obvious finding to emerge from the analysis is that high content of X configuration had a key role in this large structural deviation. In the following section, we take a complementary approach and first identify the NTD shape prototypes and will subsequently study their differences in topology and function.

3.4. AR‐NTD adopts several prototypic shapes capable of targeting AR‐LBD

Our spatial and temporal analysis study shows that AR‐NTD is dynamic, yet there are regional differences in adopting certain topologies. As a next step, we asked whether these topological transitions contribute globally to the overall shape of the NTD, and whether this would have any impacts on the so‐called N/C interaction which is mediated by ²³FQNLF²⁷ motif and a hydrophobic pocket on the surface of the AR‐LBD. To address this question, we modeled the dynamics of the protein after regional compaction was achieved and sampled the protein conformations. From the sampled conformations, we identified six different classes by performing structural clustering. Figure 3a represents the clusters together with their corresponding occupancy. Looking closer at the circuit topology of these structures, we could categorize them in three groups for whole NTD and two groups for NR region (Figure S18). Structures from the same circuit topology groups were similar in global shape and distinct from other groups. One can readily map the well‐known NTD motifs onto the surface of these structures. The motifs are surrounded by residues that may change during shape transitions (Figure S19). For example, a motif, which is surrounded by positive charges in one NTD conformation, may be surrounded by negative charges in a different shape state. These transitions could be the regulatory mechanism for the binding interactions of these motifs with their partners and be important for their function.

Strikingly, a cleft forms between the CR and NR regions to adopt orientations that resemble opening and closing movements within the AR‐NTD structure (Figure 3a). We performed docking of the identified NTD structures (Figure S20) with LBD (PDB code: 1T7T) and found that the shape of the complex shows striking resemblance to the recently published cryo‐EM image of the full‐length AR. ⁴ Particularly, the highest score docking pose from cluster 2 was nearly identical to the reported experimental structure (Figure S21). Interestingly, in the highest score docking pose from cluster 4 (Figure 3b), the ²³FQNLF²⁷ motif is located at a close distance to the AF‐2 pocket, which agrees with the well‐established models of N/C interactions. ¹ Studying best docking poses from all six clusters (a total of 600 poses) reveals that LBD interacts with NTD primarily via the NR region. This interaction dominantly formed either with only NR region or at the cleft. In the latter cases, we do find poses in which the ⁴³³WHTLF⁴³⁷ motif is in close proximity with the AF‐2 site of the LBD. Nonetheless, we do find a small number of clusters where the LBD also interacts with NTD‐CR part only.

3.5. AR NTD‐DBD interaction suggests a regulatory role via two distinct binding modes

The AR DNA‐binding domain (DBD) bind to AR response elements, and through this mediates the critical role in the AR transcriptional activity. ⁵⁹ We investigated whether and how the identified AR‐NTD conformations could interact with DBD. We used the previously reported crystal structure of the DBD (PDB code: 1R4I) and docked it against the AR‐NTD conformations we have identified, as shown in Figure 4a,b. For each AR‐NTD conformation, highest ranked poses (ranked by cluster size) were selected, and their interaction interfaces were analyzed (Figures 4 and S22). Our data revealed three interaction modes namely (1) interaction via NR, (2) interaction via CR, and (3) interactions that involves the cleft and neighboring NR/CR regions.

Strikingly, CR‐mediated binding typically covers DNA‐binding surface of DBD (residues His553, Tyr554, Ser561, Val564, Arg568, Tyr576, Arg591, Lys592, and Pro595), while NR‐mediated interactions leave this surface exposed. When DBD interacts with NR region both D and P boxes surfaces are accessible. DBD interactions with CR region limited the access to the P box, which is responsible for DNA binding, however, the D box remained exposed and accessible for dimerization. This suggests a model in which CR‐mediated NTD‐DBD interactions could occlude DNA binding, while NR‐mediated interaction would allow for DNA binding. Figure 4c summarizes NTD–DBD interaction modes and their overlap with DNA‐binding site. Our analysis also showed that no cysteine residues were involved in these interactions.

3.6. Compact CR binds to an EPI drug

Due to its medical importance, developing drugs that target AR NTD is an area of intense research; however, the conformational disorder significantly hinders drug development. Only recently, EPI‐001 has been developed as the first inhibitor of the AR‐NTD, ⁶¹ yet we do not know how the drug molecule interacts with the receptor. Interestingly, biochemical studies revealed that the binding happens in the CR part of the AR‐NTD. ⁵² To study the drug binding, we clustered the CR region after the regional compaction (last 3 μs) and identified three representative structures. Note that the clusters found in Figure 3a are the result of clustering of the whole NTD and the NR contributes to the heterogeneity due to its dynamics. It is thus expected that CR clustering would lead to a smaller number of representative structures. Modeling the three CR conformations show how the core is enclosed by the shell and how the dynamic shell can mask or expose the core. We performed docking study on the identified CR structures and the EPI‐001 drug and identified the most probable binding modes (Figure 5). The most probable drug‐CR complex is established by three hydrogen bonds. No π‐π or T‐shaped π interactions were observed between the EPI‐001 and the CR region of the NTD. In this binding mode, more than 90% of the interacting residues are located in the Tau‐5 region, which suggests a high affinity of the drug to this key functional motif, in agreement with previous reports. ⁶² , ⁶³ NMR analysis suggested that residues 354–448 are essential for drug binding, ⁶² which coincides with core of CR. This is consistent with our most probable model (interacting with residues 372–415). Furthermore, we see that residues from CR shell are also involved in drug binding (residues 527 and 529). This has not been seen in NMR study as the NTD construct used in those studies (residues 142–448) was shorter than the wild‐type and so the shell interaction could not be probed. Moreover, we studied other likely modes and, interestingly, the third CR clusters also showed binding configurations that engage both the CR core (interacting with residues 374–381) and the CR shell (residues 348 and 450) which has a complete overlap with Tau‐5 that is consistent with the NMR data. However second CR cluster showed interactions only with CR shell (residues 492–535) which involves NTD parts not evaluated in the NMR study (and thus cannot be compared).

As a case study to demonstrate the drug effect, we asked whether our identified drug‐NTD interaction can interfere with AR‐NTD for RAP74 interactions. Previous research has found Tau‐5 to be an important interaction site on the AR‐NTD for RAP74 interactions. ⁷ Tau‐5 is found to bind helices 2, 2.5, and 3 on the C‐terminal domain of RAP74 (RAP74‐CTD), which form a groove where AR can bind. ² , ⁷ Importantly, our docking analysis shows that RAP74‐CTD binds the highest ranked CR structure and the binding site overlaps both with the Tau‐5 and the aforementioned EPI‐001 drug‐binding sites (Figure S23). We see a consistent picture with other CR models as well. The data clearly suggest that the drug binding can interfere with RAP74‐CTD binding, in agreement with previous studies.