New Insights into the Binding Mechanism of Co-regulator BUD31 to AR AF2 Site: Structural Determination and Analysis of the Mutation Effect.

INTRODUCTION: Androgen Receptor (AR) plays a pivotal role in the development of male sex and contributes to prostate cancer growth. Different from other nuclear receptors that bind to the co-regulator LxxLL motif in coregulator peptide interaction, the AR Ligand Binding Domain (LBD) prefers to bind to the FxxLF motif. BUD31, a novel co-regulator with FxxLF motif, has been demonstrated to suppress wild-type and mutated AR-mediated prostate cancer growth.
METHODS: To find out the interaction mechanisms of BUD31 with WT/T877A/W741L AR complex, molecular dynamics simulations were employed to study the complex BUD31 and WT/mutant ARs. The molecular mechanics Poisson-Boltzmann surface area (MM-PBSA) results demonstrated that T877A and W741L point mutations can reduce the binding affinity between BUD31 and AR. The RMSF and dynamic cross-correlation analysis indicated that amino acid point mutations can affect the motions of loop residues in the AR structure.
RESULTS: These results indicated that AR co-regulator binding site AF2 can serve as a target for drug discovery to solve the resistance problem. Copyright© Bentham Science Publishers; For any queries, please email at epub@benthamscience.net.

INTRODUCTION

Androgen Receptor (AR) plays a crucial role in the development and progress of male reproductive tissues [1, 2]. As a member of ligand-inducible transcription factor of the Nuclear Receptor (NR) superfamily, AR shares a common structural and functional organization with other members, including an N-Terminal Domain (NTD), a DNA-binding domain (DBD) and a C-terminal Ligand-Binding Domain (LBD) containing Activation Function 2 (AF-2) [3-5]. Upon binding of the endogenous ligands, testosterone or 5α-dihydrotestosterone (DHT), AR LBD undergoes conformational changes leading to the formation of a hydrophobic co-regulator binding groove in the AF-2 of LBD [6]. This groove acts as a high-affinity binding site for a short α-helix FxxLF motif, Presenting in the AR N-terminal domain and in AR co-regulators, which can be used as a harbor to accommodate the bulky peptide side chains. AR regulates the expression of target genes through binding to the specific Androgen Response Elements (ARE) after translating into the nucleus. Meanwhile, co-regulators bind to the AF-2 to enhance AR transcription function by histone modifications, chromatin remodeling, proteasomal degradation, DNA repair and bridging of the receptor to other components of the transcription initiation process, including general transcription factors and RNA polymerase II [7-9].

It has become clear that the transcriptional activity of AR is modulated by a co-regulatory protein. Co-regulators interact with the AR to either enhance transcription or reduce the transcription of target genes but do not alter the essential transcription rate. Previous research demonstrated that AR prefers to bind aromatic-rich FxxLF motif [10-12]. Studies of the first (+1) and the fifth (+5) amino acid positions of the FxxLF motif indicate that these residues are involved in mediating high affinity and specific interactions with the AR ligand-binding domain [13]. In fact, the FxxLF-like motif exists in many AR co-regulators and plays a significant role in modulating co-regulator-mediated AR functions [14-17]. However, peptides containing a tyrosine in the +5 position of the motif were frequently identified in the screened peptides [17] by the combined phage-display screening and X-ray crystal structure analysis of the wild-type AR-DBD-LBD (WT AR-DBD-LBD) and T877A mutant AR-DBD-LBD (T877A AR-DBD-LBD) cofactor-binding grooves. An anti-AR BUD31 containing an FxxFY motif was also identified, which was found to interact tightly with AR (Fig. ). The corresponding experimental results indicate that BUD31 can block AR transcriptional activity, the AR N-C interaction and AR-mediated cancer cell growth [18-21].

Here we report the structural basis of FxxFY motif binding to both WT and mutant ARs. Molecular dynamics simulations were employed as a promising tool to understand the influence of point mutations on the binding of peptide BUD31 (Sequence: FDLFY) with AR. Key residues interacting with BUD31 were obtain and the important H-bonds were analyzed. Furthermore, the binding affinity was also calculated.

MATERIALS AND METHODS

System Preparations

The initial AR structure together with peptide BUD31 was obtained from the Protein Data Bank (PDB code: 4OED, 4OH6 and 4OKB) [22]. Crystal water molecules and heteroatom were deleted, and missing residues were complemented by preparing a ligand module in Discover studio 2.5. All of the hydrogen atoms were added using LEaP module in the AMBER 12 software [23]. The AMBER ff99 force field [24] was used as the parameter for all of the amino acid residues. The appropriate number of chloride counterions were added to maintain the electric neutrality of all the systems and then each system was immersed into a hexahedral box using TIP3P water molecules with at least 10Å distance around the complex [25].

Molecular Dynamics Simulations

To carry out MD simulations, the topology and coordinate files of the three complexes were built with LEaP module of the AMBER 12 package. Subsequently, the first minimization was performed using the steepest descent method that switched to a conjugate gradient every 2500 steps for a total of 5000 steps with 0.1 kcal/mol Å-2 restraints on the atoms of the complexes. The other two energy minimization steps were performed by only restraining the protein and further releasing all the restraints for 5000 steps of each round. After that, the system was gradually heated from 0 to 310 K over 500 ps using a Langevin thermostat with a coupling coefficient of 2.0 ps. Later on, MD equilibration of 50 ps each was performed with restraint weight 0.5 kcal/mol 
Å-2. Eventually, 50 ns production MD simulations were performed without any restraint on these three systems in the isothermal-isobaric (NPT) ensemble at a temperature of 310 K and a pressure of 1 atm. All MD phases were carried out in the NPT ensemble using a Berendsen barostat [26] with a target pressure of 1 bar and a pressure coupling constant of 2.0 ps. Long-range Coulombic interactions were analysed using the Particle Mesh Ewald (PME) [27, 28] summation. For the equilibrium and production phase of simulation, the SHAKE algorithm [29] was employed on all atoms covalently bonded to a hydrogen atom, allowing for an integration time step of 2 fs.

Binding Free Energy Calculations

The binding free energy and corresponding components of WT and its mutations were calculated by MM-PBSA methods [30]. At first, 5000 snapshots were taken from the last 5ns trajectory with an interval of 10 ps and the counterions and water molecules were stripped. Subsequently, MM-PBSA methods were employed to compute the difference in the binding free energy between WT and its mutations. In general, the complex approximate free energy was calculated by using the following equations:

ΔGbind = ΔEMM + ΔGsol – TΔS (1)

ΔEMM = ΔEval + ΔEele + ΔEvdw (2)

ΔGsol = ΔGPB + ΔGnp (3)

ΔGnp = γSASA + β (4)

Here, ΔGbind was considered as the sum of the change in the molecular mechanical (MM) gas-phase binding energy (ΔEMM), solvation energy (ΔGsol), and entropy term (– TΔS) (Eq. 1). The ΔEMM is calculated by (Eq. 2), where ΔEval, ΔEele and ΔEvdw represent the internal energy contribution from bonds, angles and torsions, electrostatic and van der Waals interactions, respectively, computed by using the same parameter set as used in the MD simulation. The solvation free energy

consist of two parts, ΔGPB and ΔGnp (Eq. 3). The electrostatic contribution to the solvation free energy is calculated by the Poisson-Boltzmann (PB) method using the Delphi program. The grid spacing was set to 0.5 Å. The nonpolar solvation contribution is described in (Eq. 4), where SASA is the solvent-accessible surface area and the solvation parameters, γ and β, which were set to 0.00542 kcal/mol Å-2 and 0.92 kcal/mol, respectively [31].

Simulation Analysis

In order to evaluate the correlation motion of residues caused by T877A and W741L mutations, cross-correlation analysis was employed to explore the fluctuation in the coordinates of the Cα atoms relative to their average positions from the last 5 ns of the simulations. Moreover, hydrogen bond analysis was also used to investigate the binding affinity, including bond, angle and probability. The hydrogen bond criteria include an acceptor-donor distance < 3.0 Å, and acceptor…H-donor angle > 135°. In addition, the average structure and hydrophobic surface were used to explore the binding mechanism of the protein and the peptide.

RESULTS

The Stability of the Studied Systems

To obtain the equilibrated and stable systems, fifty nanoseconds molecular dynamics simulations for three systems (WT-BUD31, T877A-BUD31 and W741L-BUD31) were performed. The root mean square deviation (RMSD) of backbone Cα atoms of the protein reference the structure of production dynamics were calculated and plotted in Fig. (). In Fig. (), the RMSD of the backbones of WT-BUD31, T877A-BUD31 and W741L-BUD31 show small fluctuations in the whole production process, indicating that the studied system reached stabilization.

To further explore the conformation change of the protein induced by the two mutations, Root Mean Square Fluctuations (RMSF) analysis was also employed to study the stability of the systems. The structure of the AR ligand binding domain consists of 250 residues: the loop region consist of 681-695, 722-723, 759-770, 779, 798-799, 814-822, 842-851, 889-891, and 909-918, while the binding pocket includes L712, V713, V716, K717, A719, K720. F725, H729, V730, D731, Q733, M734, I737, Q738, W741, E893, 894, E897, I898, V901, and Q902. It can be seen from Fig. () that the flexibility of the residues 690-696, 756-761, 766-800 and 841-865, significantly increased, which was mainly concentrated in the loop region. However, the flexibility of residues for T877A-BUD31 and W741L-BUD31 decreased compared to wild type, indicating that T877A and W741L mutations reduced the flexibility of the residues, especially W741L mutation. For W741L-BUD31, it maintains the lowest flexibility among these three systems. However, at the binding site, the RMSF values remained low in all these systems, such as residues 700-750 and 890-905, which demonstrates point mutation, thereby, changing the flexibility of the Loop region, which did not change the flexibility of the 
active site.

Dynamic Cross-correlation Analysis

To further investigate the residue motions caused by point mutants, cross-correlation matrices of the Cα atoms fluctuations among last 10 ns were calculated and plotted in Fig. (). Herein, highly positive regions (red and yellow) are associated with strongly correlated motions of specific residues, whereas negative regions (dark blue) are indicative of a strong anticorrelation in the specific residue movements. In Fig. (, there are a very few highly correlated motions except for the y=x diagonal, which inflect the correlation of a certain residue with itself.

Fig. (4)

Cross-correlation matrices of the fluctuations of coordinates for Cα atoms around their mean positions during the last 10 ns of MD simulation. The extent of correlated motions and anticorrelated motions are color-coded. a: WT-BUD31. b: T877A-BUD31. c: W741L-BUD31.

Compared with WT, T877A and W741L point mutations increased the anticorrelated motions near T877 and W741 residues. This result suggests that the conformation near T877 and W741 residues changed by point mutations, while the co-regulator binding site was not affected.

MM-PBSA Calculation

In order to deeply analyze the effects of point mutations between AR and BUD31, the binding free energies and corresponding components were calculated by the MM-PBSA method. Table shows that the ΔG of all three complexes are almost similar, which do not contradict with the previous report [22]. It can be seen (Table ) that ΔEele andΔGpol are the dominant factors for regulating the binding affinity. To obtain the detailed interaction profile between AR and co-regulator BUD31, MM-PBSA method was further used to decompose the interaction energies to each residue. The contribution profile of the residues in the active site (residues within 5Å around BUD31) for ΔEele, ΔGpol and ΔG are plotted in Fig. (. As mentioned above, electrostatic interaction, as one of the key components, on one hand, decreases the contribution of residues K717 and K720, on the other hand, it significantly increases the contribution of residues D731 and E897. For solvation free energy, point mutations decrease the contribution of residues D731, E893 and E897, but increase the contribution of residues K717 and K720. However, as a whole, the contribution values of the key residues were observed to be slightly different in the active site.

Fig. (5)

The residues contribution of AR to BUD31 binding. a: per residue contribution profile for ΔGele. b: per residue contribution profile for ΔGpol. c: per residue contribution profile for ΔG.

Hydrogen Bond Analysis

In the process of molecular dynamics simulations, BUD-31 can form hydrogen bonds with certain key residues, such as K720, Q733 and E897. All hydrogen bonds with a probability higher than 10% are listed in Table . From Table , it can be seen that in these three systems, N-H in the amino group of K720, as well as O-H in the hydroxyl group of Q733 and E897, can form a hydrogen bond with polypeptide BUD-31. Compared with WT AR systems, the probabilities of hydrogen bonds after T877A and W741L mutations decreased, but the distance increased. These results indicate that T877A and W741L point mutations reduce hydrogen bond stability between polypeptide with Q733. Nevertheless, hydrogen bond probability as a benchmark only reflects electrostatic interactions between the residues and polypeptide BUD-31, which cannot fully dominate the binding affinity between AR and polypeptide BUD-31.

The Structure Comparison of WT and Mutant Complexes

To reveal the contribution of the vital residues from a structural perspective, three structures with peptide and key residues in the binding pocket taken from the average MD-simulated structure are shown in Fig. (. Besides, in order to fully understand the actual distribution of key peptide residues in AR active sites, hydrophobic surfaces between AR and polypeptide BUD-31 are plotted in Fig. (. It can be seen that the peptide, composed of Phe+1, Leu+4, and Phe+5, binds to the hydrophobic groove formed by the helix3, 4, 5 and 12 of AR. As expected, BUD-31 is coordinated by residues K720, M734, M894 and E897, suggesting that AR can accommodate short peptides through these conserved residues.

Fig. (6)

The average structures of WT and its mutants take from the last 10 ns of the molecular dynamics simulations with the key residues of the binding pocket of the complexes. a: WT-BUD31. b: T877A-BUD31. c: W741L-BUD31.

Fig. (7)

The hydrophobic surfaces of the complexes of WT and its mutants: a: WT-BUD31. b: T877A-BUD31. c: W741L-BUD31. Herein, orange represents hydrophobic regions, whereas blue represents hydrophilic regions.

Point mutations lead to slight differences in the interactions between the peptide and AF2 site, such as binding affinity, and orientation of key residues in +4, +5 positions of the peptides. In the +4 position, although hydrogen bond contract was not formed, the orientation of this side chain was influenced by hydrophobic residues, such as L712, V713 and V716. As previously reported, +4 position of the groove is more suitable for accommodating bulky side chains [10]. It can be seen in Fig. ( that the side chain of

The hydrogen bonds are determined by the acceptor-donor atom distance of 3.0Å and acceptor-donor angle of 135.

+4 residues were exposed to the broad groove, indicating that the surrounding hydrophobic residues that form an appropriate orientation of the peptide side chain. At the same time, T877A and W741L mutations make the phenyl group to be distant from V713, and close to M894, which may explain to some extent that the point mutation can decrease the binding affinity.

Compared with +4 position of BUD31, +5 residues are deeply embedded into the hydrophobic groove. By comparing the hydrophobic surface (Fig. 7), we could see that both the hydrophobicity and the size of the hydrophobic pocket were obviously different in these three complexes. Although all the complexes were stabilized by two hydrogen bond contracts between the phenolic hydroxyl group of +5 peptide, hydroxyl group of Q733 and the amino group of K720, the length of the bonds changed after point mutations. This result indicates that point mutations cannot change the formation of hydrogen bond interaction, which also contributes to the difference of the bond energy between a protein and a peptide.

DISCUSSION

BUD31 is a novel co-regulator that regulates the transcriptional activity of the AR. However, the mechanism of BUD31 interaction with AR is not clear. The present work reveals how point mutations affect the efficacy of a peptide associated with a protein. The action mechanism was explored using molecular dynamics simulation, which provids structural insights into the binding mechanism of the co-regulator BUD31 to AR AF2. Though the result of binding free energy analysis suggests that T877A and W741L point mutations can decrease the binding affinity between co-regulator BUD31 and protein, the difference in the interaction energies between BUD31 and protein is not significant. Combined with RMSF analysis, one possible explanation is that the T877A and W741L mutations were unable to disturb the scaffolds of the AR, and just increased the fluctuation of the loop regions. This view was also proved by the cross-correlation analysis. The hydrogen bond analysis suggests that T877A point mutation reduced hydrogen bond stability by influencing the formation of the hydrogen bond between polypeptide and Q733.

CONCLUSION

In conclusion, when binding to the co-regulator AR AF2 site, no matter WT or mutant protein, the binding poses of peptide BUD31 have no significant difference, the possible reason is that AF2 is a highly reserved binding site. The structural insights obtained from this study are important for us to design novel peptide drugs to overcome the resistance caused by T877A and W741L mutations.

Table 1

The calculated binding free energies and corresponding components between AR and BUD31 using MM-PBSA method (unit: kcal/mol).

Terms	WT-BUD31	T877A -BUD31	W741L -BUD31
				ΔEvdwa	-50.1±1.69	-47.85±2.18	-46.86±2.26
ΔEeleb	-217.83±5.18	-238.84±8.40	-242.7±10.03
ΔGpolc	221.46±7.70	241.95±7.27	245.64±10.98
ΔGnon-pold	25.86±0.38	12.37±0.71	24.75±0.70
ΔGbinde	-20.61±5.91	-19.99±2.15	-19.17±1.14

avan der Waals interaction energies between AR and BUD31.

bElectrostatic interaction energies between AR and BUD31.

cPolar contributions to the solvation free energy.

dNonpolar contributions to the solvation free energy.

eΔGbind = ΔEele + ΔEvdw + ΔGpol + ΔGnon-pol.

Table 2

The hydrogen bonds between AR and polypeptide BUD31.

Complex	Acceptor	Donor	Distance (Å)	Angle (。)	Probability (%)
WT	+5Y@O	K720@HZ1	2.80	156.81	26.00
-	+5Y@O	K720@HZ2	2.80	156.66	19.40
-	+5Y@O	K720@HZ3	2.80	155.84	20.63
-	Q733@OE1	+5Y@OH	2.75	164.89	86.80
-	E897@OE2	+1F @N	2.85	157.79	43.45
-	E897@OE1	-2Y@N	2.86	153.26	40.00
T877A	+5Y@O	K720@HZ1	2.83	155.64	14.70
-	+5Y@O	K720@HZ2	2.81	155.58	15.87
-	+5Y@O	K720@HZ3	2.82	155.47	14.33
-	Q733@OE1	+5Y@OH	2.76	162.21	77.80
-	E897@OE1	+1F @N	2.89	165.59	30.47
-	E897@OE2	-1I@H	2.83	154.04	68.00
-	E897@OE2	-2Y@H1	2.78	154.63	25.00
-	E897@OE2	-2Y@H2	2.79	154.25	23.00
-	E897@OE2	-2Y@H3	2.79	154.06	25.00
W741L	+5Y@O	K720@HZ1	2.81	155.63	18.93
-	+5Y@O	K720@HZ2	2.81	157.00	18.83
-	+5Y@O	K720@HZ3	2.80	155.87	18.07
-	Q733@OE1	+5Y@OH	2.75	164.62	87.57
-	I242@O	+4F @N	2.89	155.59	25.33
-	Y241@O	+3L@H	2.90	156.44	25.00
-	E897@OE2	+1F@H	2.87	164.27	33.00