Systematic protein-protein docking and molecular dynamics studies of HIV-1 gp120 and CD4: insights for new drug development.

BACKGROUND AND THE PURPOSE OF THE STUDY: The interactions between HIV-1 gp120 and mutated CD4 proteins were investigated in order to identify a lead structure for therapy based on competitive blocking of the HIV binding receptor for human T-cells. Crystal structures of HIV gp120-CD4 complexes reveal a close interaction of the virus receptor with CD4 Phe43, which is embedded in a pocket of the virus protein.
METHODS: This study applies computer simulations to determine the best binding of amino acid 43 CD4 mutants to HIV gp120. Besides natural CD4, three mutants carrying alternate aromatic residues His, Trp and Tyr at position 43 were investigated. Several docking programs were applied on isolated proteins based on selected crystal structures of gp120-CD4 complexes, as well as a 5 ns molecular dynamics study on the protein complexes. The initial structures were minimized in Gromacs to avoid crystal packing effects, and then subjected to docking experiments using AutoDock4, FireDock, ClusPro and ZDock. In molecular dynamics, the Gibbs free binding energy was calculated for the gp120-CD4 complexes. The docking outputs were analyzed on energy within the respective docking software. RESULTS AND
CONCLUSION: Visualization and hydrogen bonding analysis were performed using the Swiss-PdbViewer. Strong binding to HIV gp120 can be achieved with an extended aromatic group (Trp). However, the sterical demand of the interaction affects the binding kinetics. In conclusion, a ligand for an efficient blocking of HIV gp120 should involve an extended but conformational flexible aromatic group, i.e. a biphenyl. A docking study on biphenylalanine-43 confirms this expectation.

INTRODUCTION

HIV-1 infection starts with the binding of the viral surface-protein gp120 to the human T-cell receptor CD4 via a hydrophobic pocket on gp120, which binds effectively to Phe43 in CD4 (1). This initial association is followed by membrane fusion and transfer of the viral genetic material into the cell, thus enabling the reproduction of the virus whilst destroying the host. The crucial role of the gp120-CD4 interaction has been demonstrated by an experiment leading to HIV-1 vulnerable rats after insertion of human CD4-CDR2 protein sequences (2), while the impact of Phe43 was confirmed by significant loss of binding affinity of CD4 for gp120 after mutation around the binding site (3, 4). Therefore the amino acid sequence in gp120 is highly conserved and less prone to mutations (1, 2), thus providing a promising target for the development of new drugs (5).

Blocking the cellular entry of the virus is more practical and exhibit less side effects compared to a drug operating at intracellular level, e.g. a protease inhibitor (6). BMS-378806 inhibits the gp120-CD4 binding, as shown in an enzyme-linked immunosorbent assay (ELISA), without inhibitory activity against HIV-1 reverse transcriptase, protease and integrase and it has been proposed as a good candidate against HIV-1 infection (7, 8). A clinical study of gp120 inhibition has not been conducted so far, however an oligopeptide CD4 mimic has been used to study the inhibition of the HIV-1 entry using a cell-based fusion assay (9).

Several computational studies have been performed on docking and molecular dynamics of small ligands with either gp120 (5, 10) or CD4 (11) to identify suitable inhibitor candidates. Another investigation applied molecular dynamics on the gp120-CD4 complex targeting to predict a mimic for the natural Phe43 conformation in the complex (1). This is the first study that conducts computational protein-protein interactions to propose the design of a therapeutic strategy.

MATERIAL AND METHODS

Generic approach

PDB files for gp120-CD4 protein complexes of HIV-1 (12) from up to 2010 were downloaded and checked for the presence of the indicated hydrophobic interaction via Phe43. As none of the structures had additional molecules at the binding site, they were stripped from all antibodies, water and ligands to lower the computational time. In order to avoid crystal packing and antibody induced conformational effects, the protein complexes were minimized in Gromacs (13) using steepest-descent, followed by conjugate gradient approach to an energy convergence of 0.01 kJ/mol. Out of 22 crystal structures, only 1rzk (14) and 1g9n (15) remained intact and were used for the study. The root mean square deviations (RMSD) of the structures before and after minimization were computed using the Swiss-PdbViewer (16). Finally, the complex was separated into individual proteins.

The Phe43 cap was mutated in the original protein complex using the Swiss-PdbViewer to His43 (neutral and protonated-generated by ADT (17)), Trp43 and Tyr43. Docking studies applied a keylock protocol with receptor (gp120) and ligand (CD4) in ten replicates on AutoDock4 (18), PatchDock & FireDock (19, 20), ZDock (21) and ClusPro (22). The latter three applications were run on the respective websites. Protein binding energies were determined within the software packages for AutoDock and PatchDok/FireDock, while ZDock and ClusPro output files were minimized in Gromacs and the energy was compared with the accumulated energy of the individual proteins.

Rigid & flexible receptors

Docking studies were performed with both rigid and flexible receptors and conformational flexible ligand at position 43. The amino acid selection for flexible domains of the receptor was based on closest proximity, including Thr108, Asp190, Glu192, Ile193, Asn239, Met240, Trp241 and Gly287 for 1g9n and Asp190, Glu 192, Asn239 and Trp241 for 1rzk (Figs. 1a and 1b). Exceptionally, ClusPro uses a blind docking approach without selection of flexibility.

Figure 1

Stick representation shows the binding sites containing flexible amino acid residues in the hydrophobic pocket of gp120 and the CD4 Phe43 cap. For 1g9n eight amino acid residues are in close proximity to the Phe43 cap (a), while only four amino acid residues fulfill the criteria in 1rzk (b).

Analysis and visualization of docked conformations

Final docking output files were analyzed for hydrogen bonds using Swiss-PdbViewer with default settings. FireDock provides an additional analysis tool to quantify long and short range electrostatic interactions, van der Waals bindings and p-p stacking energies. The RMSDs of docked conformations with respect to the minimized crystal structure of the corresponding protein complex were calculated using Swiss-PdbViewer for the amino acid around the binding site, which were used as criteria for the gap-cap interaction between CD4 and gp120 to complement the qualitative visual analysis of the docking structure.

Molecular dynamics simulations (MD) of gp120-CD4 complex

Molecular dynamics were applied on both individual proteins as well as on the protein complexes. The initially minimized structures were soaked in SPC explicit water solvent (23) and minimized in Gromacs by steepest-descent and conjugate gradient approach to energy convergence of 0.01 kJ/mol. MD was performed for 5ns under constant pressure of one atmosphere with a pressure coupling constant of 1.0 ps and at 310 K with a temperature coupling constant of 0.1 ps at 1 femtosecond time step. The g_lie module in Gromacs was used to calculate the Gibbs free energy of binding (DGbind) of Aqvist (24) as shown in Eq. 1:

The approach is based on the assumption that a receptor-ligand protein-protein interaction can be estimated by determining the difference in solvent (water) interaction energy of the free ligand (Efree) and the complex (Ebound). The remaining variable, i.e. the receptor and its water interaction, can be ignored, as it remains constant in the comparison, thus only affects the total energy, but not the difference between the ligands. The hydrogen bonding was calculated using the g_hbond module in Gromacs. The RMSDs and potential energies were calculated by Gromacs modules g_rms and g_energy, respectively. All analyses were computed over the full 5ns MD trajectories.

RESULTS AND DISCUSSION

The selection of PDB crystal structures was based on a continuous peptide structure after minimization. 1rzk and 1g9n are the only structures that do not expose unnatural binding distances due to missing amino acids. Minimization of the protein complexes only led to minor changes in RMSDs, i.e. below 2 Å, as shown in table 1. The results for the docking binding energy and compliance with the gap-cap interaction are summarized in table 2. The latter is demonstrated in figures 2a and 2b. The binding energy corresponds to best out of ten runs that complies with the gap-cap interaction. Only if such structure was not found, the global energy minimum was selected.

Table 1

Root Mean Square Deviations (RMSDs) of the minimized structures with respect to the initial crystal structure for 1rzk and 1g9n protein complexes with their respective amino acid mutants as calculated by Swiss-PdbViewer.

Amino acid mutants		RMSD [Å]
		whole structure		binding site
		Ca	all atoms	Ca	all atoms
1rzk	Phe	1.18	1.56	0.26	0.51
	His	1.15	1.52	0.31	0.57
	Trp	1.24	1.63	0.39	0.79
	Tyr	1.27	1.67	0.32	0.63
1g9n	Phe	1.16	1.51	0.48	0.60
	His	1.26	1.62	0.47	0.61
	Trp	1.25	1.61	0.57	0.81
	Tyr	1.26	1.62	0.52	0.66

Table 2

Docking energies [kJ/mol] and gap-cap conformation fractions for 1rzk and 1g9n using AutoDock 4, ClusPro, ZDock and FireDock. Values in parenthesis indicate the rate of proper pocket gap-cap conformations over a total of ten dockings. Best docking conformation refers to lowest energy with highest gap-cap conformation fraction, when possible; lowest docking energies and highest gap-cap conformation fractions are highlighted (in asterisks). As FireDock provided best results, the data are specially highlighted.

Amino Acid-43			AutoDock4		ClusPro	ZDock	FireDock
		Rigid	Flexible	Blind	Directed	Directed	Blind
	Phe	-2.0	-3.1	-5798*	303164	-88.5	-87.3
1rzk		(0/10)	(0/10)	(1/10)	(3/10)	(5/10)*	(2/10)
	His(neutral)	-0.9	-5.7	-4851	291689	-85.0	-79.1
		(0/10)	(0/10)	(0/10)	(3/10)	(2/10)	(1/10)
	His(+H+)	-	-	-	-	-87.3	-84.4
							(3/10)	(1/10)
	Trp	-1.7	-4.3	-5057	281363	-95.0*	-101.3*
		(0/10)	(0/10)	(1/10)	(5/10)	(3/10)	(1/10)
	Tyr	-1.0	0.7	-5258	288568	-77.0	-93.0
		(0/10)	(0/10)	(1/10)	(3/10)	(3/10)	(2/10)
1g9n	Phe	-3.9	-0.2	-4859	289117	-84.2	-84.8
		(0/10)	(0/10)	(1/10)	(3/10)	(2/9)	(2/10)
	His(neutral)	-4.2	9.7	-5063	281750	-85.7	-84.6
		(0/10)	(0/10)	(1/10)	(3/10)	(6/10)*	(2/10)
	His(+H+)	-	-	-	-	-86.6	-86.6
						(6 /10)*	(1/10)
	Trp	-5.1	-6.8	-5220*	285966	-96.9*	-101.7*
		(0/10)	(0/10)	(1/10)	(2/10)	(1/10)	(1/10)
	Tyr	-4.7	-3.6	-4422	262119	-6.5	-90.3
		(0/10)	(0/10)	(0/10)	(3/10)	(0/6)	(3/10)

Figure 2

Conformational analysis of gp120-CD4. Docking conformation (a) reflects the proper hydrophobic gap-cap interaction, while (b) refers to a non-gap-cap binding, as listed in table 2.

FireDock and ClusPro provided both, proper conformations and energies for the dockings, while ZDock appeared to be disfavored by endothermic energies. However, the external energy calculation for ZDock and ClusPro may be misleading. Therefore, FireDock proved to be the best docking tool. The RMSD of the docked conformation with respect to the minimized crystal structure was used to supplement the visual conformation inspection. Table 3 lists the RMSD results for all atoms at the binding site. For qualified gap-cap conformations the limit was set to 1 Å. FireDock results indicate easiest bindings for small aromatic amino acids (phenylalanine and histidine) but stronger binding for the larger aromatic tryptophan.

Table 3

Root mean square deviations (RMSDs) of the docked conformations (1-10 for FireDock) for 1rzk and 1g9n protein complexes (directed and blind approaches) with reference to the minimized crystal structures as calculated by Swiss-PdbViewer. All atoms of the binding site were considered.

Docked conformations	1 rzk [Å]							1g9n [Å]
	directed					blind		directed					blind
	His		Phe	Trp	Tyr	His		His		Phe	Trp	Tyr	His
	±0	+H+				±0	+H+	±0	+H+				±0	+H+
1	0.36	0.41	0.43	0.42	0.39	0.41	0.39	0.37	0.32	0.40	0.33	3.68	0.39	0.32
2	0.38	0.40	0.34	0.43	0.41	x	x	0.37	0.38	0.36	x	x	0.37	x
3	1.61	0.36	0.41	0.40	x	x	x	0.39	0.34	x	x	x	x	x
4	x	2.86	0.41	x	x	x	x	0.32	0.38	x	x	x	x	x
5	x	x	0.31	x	x	x	x	0.18	0.36	x	x	x	x	x
6	x	2.47	1.36	x	x	x	x	3.41	0.40	x	x	x	x	x
7	x	x	x	x	x	x	x	0.49	x	x	x	-	x	x
8	x	x	x	x	x	x	x	x	x	x	x	-	x	x
9	9.69	2.92	x	x	x	x	x	x	x	x	1.79	-	x	x
10	x	x	x	x	1.02	x	x	2.04	2.00	-	x	-	x	x

x : not analyzed due to non-gap-cap conformation when inspected by eyes;

-: unresolved.

Italic numbers refer to RMSDs for structures that are likely to fail the visual inspection.

Bold numbers indicate proper gap-cap interaction; criteria: RMSD ≤ 1.0 Å

Table 4 shows a breakdown of the total binding energy for FireDock. The sums of individual energies exceed the total binding energy in table 2 indicating protein folding. The strong binding of gp120 and Trp43 mutated CD4 corresponds with large van der Waals interaction and p-p-stacking. This finding matches previous reports on interactions of Phe-43 with its neighboring Trp and Tyr (25, 10). Long range electrostatic interactions are irrelevant. No intermolecular hydrogen bond in the gap-cap region could be found, except a hydrogen bond between the phenolic hydroxyl group of CD4 Tyr43 and the Asncarbonyl in gp120.

Table 4

Contribution of different interaction types (electrostatic, hydrogen bonding, vdW and π-π stacking) [kJ/mol] to the total binding energy for protein-protein complexes of gp120 and (mutated) CD4 in FireDock. Lowest energies are highlighted (in asterisks).

Amino Acid-43			Electrostatic		Hydrogen bonding	vdW	π-π stacking
			Short-range	Long-range
Directed Bound & Unbound	1rzk	Phe	−106.1*	7.9	−14.3	−57.6	0.0
		His	−80.6	−10.9	−13.7	−42.3	0.0
		His (+H+)	−99.3	4.1	−14.5	−54.0	0.0
		Trp	−82.4	8.3	−15.0*	−62.6*	−0.5*
		Tyr	−78.7	2.6	−14.4	−57.2	0.0
	1g9n	Phe	−83.3	4.2	−15.7*	−55.4*	0.0
		His	−88.6	−3.6	−13.7	−54.0	0.0
		His (+H+)	−118.2*	4.1	−12.8	−57.6*	0.0
		Trp	−98.4	−4.6	−12.1	−55.2*	−1.5*
		‡Tyr	−8.4	0.4	−3.8	−11.7	−3.0
Blind Bound & Unbound	1rzk	Phe	−98.1*	0.7	−11.3	−49.8	0.0
		His	−95.1	2.7	−11.6	−51.2	0.0
		His (+H+)	−97.2	5.1	−12.1	−56.8	0.0
		Trp	−73.0	−2.0	−13.4	−64.1*	−0.5*
		Tyr	−92.4	−1.2	−13.6	62.2	0.0
	1g9n	Phe	−96.5	5.1	−11.0	−57.7	0.0
		His	−122.5*	6.9	−15.1	−53.1	0.0
		His (+H+)	−118.2*	4.1	−12.8	−57.6	0.0
		Trp	−98.2	−7.2	−13.7	−62.3*	−1.5*
		Tyr	−104.4	−4.0	−11.6	−59.7	0.0

The data for the directed binding of the Tyr43 mutated CD4 (boxed) were not considered due to the non-gap-cap based protein-protein binding.

Molecular dynamics investigations of the minimized complexes converged energetically but the RMSD indicates partial equilibrium. Figures 3a and 3b show one example. The Gibbs free binding energies, DGbind, for 1rzk indicate slightly higher binding for tryptophan and tyrosine (Table 5), while no significant differences are found for 1g9n. This may be considered as a confirmation of FireDock results (Table 2) and matches with the significant decrease of affinity for gp120 by Tyr43(F43Y) mutated CD4 compared to the Trp43(F43W) analog (26). No intermolecular hydrogen bonding between ligand and receptor was found around the binding site.

Figure 3

Autocorrelation functions of molecular dynamics – energy profile (a) and root mean square deviations (RMSDs of protein backbone) (b) for 1rzk-Trp.

Table 5

Gibbs free binding energy in solvent water of gp120-CD4 complex and CD4 after 5ns of molecular dynamics simulation. ∆Gbind was obtained from the g_lie module applying factors a = 0.181 and = 0.5. All energies are averaged, lowest binding energies are highlighted (in asterisks).

Amino acid-43 mutants	Energy for proteins in water [kJ mol-1]				∆Gbind [kJ mol-1] (g_lie)
	gp120-CD4		CD4		Absolute	Normalized
	EVdW,bound	ECoul, bound	EVdW,free	ECoul, free
1rzk	Phe	750.3	51959.2	202.0	19806.1	-16176 ± 20	90 ± 38
	His	737.6	51951.5	226.6	19834.8	-16151 ± 28	115 ± 46
	Trp	704.9	52133.2	214.7	19779.3	-16266 ± 18*	0*
	Tyr	754.9	52093.8	220.9	19769.5	-16259 ± 26	7 ± 44

In order to evaluate the prediction of a flexible aromatic inhibitor, an additional docking of gp120 by biphenylalanine (BiPhe)-mutated CD4 was performed. The biphenylalanine was modeled with Chem3D and the aromatic moiety was attached to the backbone of the corresponding CD4-Phe43. FireDock produced only 1/10 gap-cap docked conformations for 1rzk and 1g9n with docking energies of -103.9 kJ mol-1 and -92.3 kJ mol-1, respectively. This refers to best binding for 1rzk, and second best for 1g9n after Trp. These results agree well with the experimental result of enhanced binding affinity of gp120 for biphenylalanine Bip43-mutated CD4 (9). For a drug design, the recommendation is a larger aromatic system with conformational flexibility, to avoid inefficient binding kinetics.

CONCLUSIONS

Hydrophobic interactions and p-stacking are the dominating features in the binding of gp120 and CD4. While extended aromatic systems on the mutated CD4 can enhance the binding strength, the docking kinetics was reduced, leading to less efficient binding. However, conformational flexibility in the aromatic residue may overcome that obstacle. Additional hydrogen bonding based interactions involving a hydrogen donor on the CD4 may further improve the efficiency of inhibition, provided that the hydrogen donor does not affect the hydrophobic interaction.

Out of the investigated four docking tools, only the network FireDock and ZDock applications provided good results with respect to the docking conformation, while AutoDock failed to provide any suitable structure. The missing feature of an integrated energy evaluation in ZDock makes FireDock the software of choice for protein-protein docking studies. Molecular dynamics can only support docking studies qualitatively, as the results depend on the input structure.