J. Phys. Chem. 1996, 100, 2475-2479
2475
ARTICLES Homology Modeling and Molecular Dynamics Simulations of the Gla Domains of Human Coagulation Factor IX and Its G[12]A Mutant Leping Li Department of Chemistry, UniVersity of North Carolina, Chapel Hill, North Carolina 27599-3290
Tom Darden National Institute of EnVironmental Health Sciences, Research Triangle Park, North Carolina 27709
Richard Hiskey and Lee Pedersen*,† Department of Chemistry, UniVersity of North Carolina, Chapel Hill, North Carolina 27599-3290 ReceiVed: August 1, 1995; In Final Form: September 20, 1995X
We have tested molecular dynamics as a predictive tool for stability of protein structure. We first used the X-ray crystallographic coordinates of bovine prothrombin fragment 1 to model the structure of the Gla domain of human factor IX, an essential coagulation protein. In addition, a mutant protein that differs in only a single amino acid (G[12]A) but is associated with one form of Hemophilia B was also modeled. Simulations performed in an identical fashion for both proteins, with considerable care to accommodate long-range forces in these highly ionic systems, indicated that substantial distortions occur in the mutant structure relative to the wild type factor IX. Surprisingly, the Gla domain (residues 1-34) of the wild type and mutant remain similar. Instead, it is the interaction of the added Ala at position 12 in the mutant that repels the post-Gla region (residues >33, the C-terminal helix) region.
1. Introduction Molecular dynamics (MD) simulation techniques have been widely used to obtain time-dependent atomic-level structural information of liquids and macromolecules. Because of the limitations of experimental techniques, MD simulations provide an alternative tool for biological problems that are complimentary to experimental techniques. The successful applications of long MD simulations have become increasingly possible with the awareness of the importance of accommodating long-range ionic forces.1-11 In our MD simulation of bovine prothrombin fragment 1 bound with Ca2+ (bf1/Ca),12 we found that unfolding of the calcium bound region of the molecule occurred unless the long-range ionic forces were properly accommodated. In this work, we focus on the Gla domain of human coagulation factor IX and that of a single substitution mutant protein known to have devastating consequences. Human coagulation factor IX is a vitamin K-dependent plasma glycoprotein which plays an important role in both the extrinsic and intrinsic pathways of blood coagulation.13,14 The zymogen form of factor IX is a single-chain peptide of 416 amino acids.15 Upon activation by either factor XIa or the VIIatissue factor (VIIa-TF) complex, factor IX is converted into a two-chain active serine protease, factor IXaβ. The resultant active factor IXaβ proteolytically activates factor X in the presence of Ca2+, phospholipids (PL), and cofactor VIIIa. While the heavy chain of factor IXaβ contains the serine protease domain, the light chain consists of a Gla domain (residues 1-34) † Also at the National Institute of Environmental Health Sciences, Research Triangle Park, NC 27709. X Abstract published in AdVance ACS Abstracts, January 15, 1996.
0022-3654/96/20100-2475$12.00/0
and two EGF-like domains. As with other vitamin K-dependent proteins, the Gla domain of factor IX is rich in γ-carboxyglutamic acids (a total of 12 Gla residues). These Gla residues provide ligands for Ca2+ and are essential for the Ca2+dependent properties of the Gla domains such as PL binding. Recently, a two-dimensional NMR study of the Gla domain of human factor IX showed that the metal-free Gla domain is flexible with the majority of the peptide solvent accessible.16 Similarly, in a NMR study of the conformation of the GlaEGF domains of bovine factor X, Sunnerhagen et al. found that the Gla domain of the Ca-free form of bovine factor X is mainly helical with the helical regions connected by flexible hinges.17 Comparison of the structures of the Ca-free Gla domain (NMRderived) with the Ca-bound Gla domain (modeled) of bovine factor X showed that there is a striking difference between the two structures.17 In the Gla domain of the Ca-free form, the hydrophobic residues Phe4, Leu5, and Val8 in the Ω-loop are buried inside the domain and the hydrophilic residues Gla7 and -8 are exposed to solvent. In the Ca-bound form of factor X, these relationships are switched. It was thus concluded that Ca-induced structural changes were essential for exposing the PL binding epitope of the Ω-loop region.17 It was also suggested that these structural changes in factor X induced by calcium ions were general for all vitamin K-dependent coagulation proteins. Although the Gla domains of vitamin Kdependent coagulation proteins share high homology with evidence for a common metal ion-induced conformational change in the Gla domains of several vitamin K-dependent proteins,18,19 differences in PL binding20-23 are apparent. Although it would be most useful to have three-dimensional structures of the Gla domains of all of the vitamin K-dependent © 1996 American Chemical Society
2476 J. Phys. Chem., Vol. 100, No. 7, 1996
Figure 1. Sequences of bovine prothrombin fragment 1 (residues 1-47) and human factor IX (residues 1-47). Identical residues are in bold. X, γ-carboxyglutamic acid.
coagulation proteins, only one well-defined structure is currently known.24 A previous model structure for the Gla domain of factor IX with seven calcium ions was built that employed distance geometry constraints.22 The model structure was refined using MD simulation with a small shell of water. However, no analysis was provided. It has been suggested that the hydrophobic helical stack (HS) (Thr39-Asp47) adjacent to the Gla domain (residues 1-34) is important for the Gla-dependent properties of vitamin Kdependent coagulation proteins.25-29 Thus, we modeled the three-dimensional structure of the Gla domain plus the HS of human factor IX (residues 1-47) using the X-ray crystallographic coordinates bf1/Ca.24 In addition, we also computationally mutated Gly12 to Ala in a separate simulation. Hemophilia B is a classification for bleeding disorders that result from mutations in factor IX; G[12]A is one of these naturally occurring mutations.30 Since Gly12 is largely conserved in human vitamin K-dependent proteins, it is of considerable interest to examine the effect of the mutation on the global structure. We then performed MD simulations on both the wild type and G[12]A mutant of hfix/Ca (1-47) in order to understand the structural basis for Hemophilia B that results from the mutation. 2. Methodology 2.1. Model Construction. The X-ray crystallographic coordinates of bf1/Ca (1-47) were extracted as a template for the construction of the Gla domain of human factor IX (residues 1-47) [hfix/Ca (1-47) ] in a manner similar to that reported in the modeling of human prothrombin fragment 1 with calcium ions (hf1/Ca).31 The sequence identity between bf1 (1-47) and hfix (1-47) is 53% (Figure 1). With such high identity, homology modeling can be used to construct a reliable model with confidence.32 Briefly, all identical residue pairs were assumed to initially have the same conformations as those in the template and were kept fixed throughout the initial model building stages. Side chains that differed from the template were mutated to the corresponding amino acids in the target in such a way that the conformation of the common part of the pair was retained when possible. Conformational searches for nonidentical residues were performed when necessary to avoid repulsive van der Waals contacts. Crystal water molecules from bf1/Ca structure were preserved as much as possible. Two additional calcium ions (Ca-8 and -9) were placed at the side chain(s) of Gla33, and Gla36, and -40, respectively. These two calcium ions are assumed to form a malonate complex with the side chain oxygen atoms of these Gla residues as is the case for most of calcium ions in bf1/Ca. The introduction of these two calcium ions to hfix/Ca (1-47) is based on the following observations: (1) the pair Gla36 and -40 is unique to factor IX, and the side chain oxygen atoms of these residues are near enough to form a metal ion binding site; (2) mutation of Gla33 to Asp in factor IX results in Hemophilia B,33 suggesting that Gla33 is important, reasonably through its Ca2+ interaction; (3) the net charge of the system is zero after the addition of these
Li et al. calcium ions. The potential for the side chain of Gla40 to bridge to an EGF domain calcium has also been suggested.34 The resulting model structure was energy-optimized using the AMBER 3.0A all-atom force field with a distance-dependent dielectric function in a manner similar to that reported previously.31 2.2. Molecular Dynamics Procedure. The AMBER version 3.0A program was used for the molecular dynamics simulations. Several modifications were made to improve throughput.12 Briefly, a modified Newton method suggested by Ryckaert et al.35 was used to replace SHAKE. Computation of the nonbonded list was modified by implementation of a grid technique. At every time step, the box was recentered to the center of the protein coordinates. These modifications increase the computational speed by a factor of 2. The partial charges and parameters for Gla and calcium were taken from Maynard et al.36 and Charifson et al.37 The simulations were performed according to the protocol used in our simulations of bf1/Ca12 and hf1/Ca.31 Briefly, the protein and crystal water were placed in a large box of Monte Carlo TIP3P water with each side being at least 14.0 Å away from the nearest protein atom. Water molecules with oxygen closer than 2.8 Å or hydrogens closer than 2.0 Å to any protein atoms were removed. Approximately 7000 water molecules were necessary to solvate the system. The water was then energy-minimized using the distancedependent dielectric function at constant volume. A 40 ps amount of dynamics were performed on water only while holding the protein fixed. The water was then reminimized followed by full energy minimization of the whole system. The simulations utilized the AMBER all-atom force field (dielectric constant ) 1) with a step size of 1 fs and were performed at 300 K after an initial stepwise heating over 3 ps. Nonbonded interactions were updated every 10 steps. A twin range cutoff (10 and 21 Å) was applied to calculate van der Waals and electrostatic interactions. Inclusion of the long-range interactions was found to be essential for stabilizing the protein structure in our previous simulations.6,12 A Cartesian position constraint of 5 (kcal/mol)/Å on the CR’s of the last three residues (Tyr45-Asp47) was applied. All covalent bonds involving hydrogen atoms were constrained using the Ryckaert procedure.35 The calculations were performed on a CRAY Y-MP supercomputer at the National Cancer Institute and in parallel on a Silicon Graphics Iris 4D/380-VGX workstation. 3. Results and Discussions 3.1. Comparison of hfix/Ca (1-47) and bf1/Ca (1-47). The root mean square deviations (rmsd) of the backbone atoms of the starting model (energy-minimized) structure of hfix/Ca (1-47) with respect to that of the X-ray homologous bf1/Ca (1-47) is 0.44 Å. The N-terminal networks for both hfix/Ca (1-47) and bf1/Ca (1-47) are similar (Table 1), despite the fact that position 1 is Tyr for hfix while it is Ala for bf1. The bulky side chain of Tyr1 is not buried in the interior of the Gla domain. Rather, it is partially exposed to solvent and interacts with the hydrophobic part of side chains Thr3 and Lys22. The hydroxy group is not hydrogen-bonded to any protein atoms. Interestingly, replacement of Tyr1 by Ala in factor IX showed no effect on PL binding and coagulant activity of factor IX.22,38 As was the case for bf1/Ca, Ca-4 and -5 are in close contact with the N-terminus (Table 1). The overall shapes of the Ω-loops and the entire Gla domains for both proteins are similar. The secondary structure elements found in bf1/Ca (1-47) are also present in hfix/Ca (1-47) when analyzed using the program DSSP.39 In both structures, the C-terminal long helix (Thr35-Val46) is preceded by two short helices (Phe25-Phe32 and Leu14-Gla17).
Gla Domains of Human Coagulation Factor IX and Its Mutant TABLE 1: Comparison of the Tyr1-N Network for hfix/Ca (1-47), G[12]A hfix/Ca (1-47), av-hfix/Ca (1-47), and G[12]A av-hfix/Ca (1-47)a
hfix/Ca (1-47)
av-hfix-Ca (1-47)
G[12]A hfix/ Ca (1-47)
G[12]A av-hfix/ Ca (1-47)
D-A
D‚‚‚A distance (Å)
D-H‚‚‚A angle (deg)
Tyr1N-Glal7-OE3 Tyr1N-Lys22-O Tyr1N-Gla27-OE3 Tyr1N-Gla27-OE4 Tyr1N-Ca-4 Tyr1N-Ca-5 Tyr1-N-Gla17-OE3 Tyr1-N-Gla21-OE4 Tyr1-N-Gla27-OE4 Tyr1-N-Ca-4 Tyr1-N-Ca-5 Tyr1-N-Gla17-OE3 Tyr1-N-Lys22-O Tyr1-N-Gla27-OE3 Tyr1-N-Gla27-OE4 Tyr1-N-Ca4 Tyr1-N-Ca5 Tyr1-N-Gla21-OE4 Tyr1-N-Gla27-OE4 Tyr1-N-Ca-4 Tyr1-N-Ca-5
2.7 2.9 3.0 2.8 4.5 4.3 2.5 2.6 2.7 4.6 3.8 2.7 2.8 2.9 2.7 4.4 4.4 2.5 2.9 4.2 4.1
144 140 136 148
J. Phys. Chem., Vol. 100, No. 7, 1996 2477 TABLE 2: Comparison of Calcium Ion Coordinations for hfix/Ca (1-47), G[12]A hfix/Ca (1-47), av-hfix/Ca (1-47), and G[12]A av-hfix/Ca (1-47)a hfix/Ca Ca-1
Ca-2 125 154 138 168 143 152 131 162 158
OE2, Gla26 OE3, Gla26 OE1, Gla30 OE3, Gla30 OE2, Gla8 OE4, Gla8 OE2, Gla27 OE3, Gla30 OE4, Gla30
G[12]A hfix/Ca
av-hfix/Ca
G[12]A av-hfixCa
OE2, Gla26 OE3, Gla26 OE1, Gla30 OE3, Gla30 OE2, Gla8 OE4, Gla8 OE2, Gla27 OE3, Gla30
OE2, Gla26 OE3, Gla26 OE1, Gla30 OE3, Gla30 OE2, Gla8 OE4, Gla8 OE2, Gla27 OE3, Gla30
OE2, Gla26 OE3, Gla26 OE1, Gla30 OE3, Gla30 OE2, Gla8 OE4, Gla8 OE2, Gla27 OE3, Gla30
OE3, Gla26 OE4, Gla26 OE1, Gla27 Ca-3
Ca-4
D ) donor, A ) acceptor, and H ) hydrogen. Hydrogen bonding criteria: D-H‚‚‚A angle g 125°, D‚‚‚A distance e 3.0 Å. a
OE3, Gla8 OE4, Gla8 OE1, Gla17 OE2, Gla17 OE1, Gla27 OE4, Gla30 OD1, Asn2 OE1, Gla7 OE4, Gla7 OE4, Gla7 OE1, Gla17 OE1, Gla27 OE4, Gla27
OE3, Gla8 OE4, Gla8 OE1, Gla17 OE2, Gla17 OE1, Gla27 OE4, Gla30 OD1, Asn2 OE1, Gla7 OE4, Gla7 OE1, Gla17 OE1, Gla27 OE4, Gla27
OE4, Gla8 OE1, Gla17 OE2, Gla17 OE1, Gla27 OE4, Gla30 OD1, Asn2 OE1, Gla7 OE4, Gla7 OE4, Gla7 OE1, Gla17 OE4, Gla27 OE2, Gla27
OE3, Gla8 OE4, Gla8 OE1, Gla17 OE2, Gla17 OE1, Gla27 OD1, Asn2 OE1, Gla7 OE4, Gla7 OE1, Gla17 OE1, Gla27 OE4, Gla27 OE2, Gla7 OE3, Gla17
Ca-5
Ca-6
Ca-7
Ca-8 Figure 2. Root mean square deviations of the backbone atoms of avhfix/Ca (1-47) and G[12]A av-hfix/Ca (1-47) compared to hfix/Ca (1-47) (0 ps) and G[12]A hfix/Ca (1-47) (0 ps), respectively: (solid light curve) G[12]A hfix/Ca (1-47); (solid dark curve) hfix/Ca (147). Corresponding comparisons: (dotted light curve) G[12]A hfix/Ca (1-34), (dotted dark curve hfix/Ca (1-34) are shown for residues 1-34.
3.2. Comparison of Average Structure [av-hfix/Ca (147)] with hfix/Ca (1-47). Both fragments 1-47 and the Gla domain (1-34) reached equilibration after 80 ps on the basis of the rmsd of the backbone atoms (Figure 2). An average structure [av-hfix/Ca (1-47)] was generated from the atomic coordinates of the period 88-115 ps. The rmsds for the backbone atoms for both av-hfix/Ca (1-47) and the Gla domain compared to those of the 0 ps structure are 2.1 and 1.3 Å, respectively. The Gla domain reached equilibration quickly and showed small rms fluctuations which suggest that the model structure is stable upon solvation. On the other hand, the relatively large rmsd for the backbone atoms of the whole fragment (1-47) compared to that of the Gla domain (1-34) may reflect movement between the Gla domain and the C-terminal long helix (Thr35-Val46). During the course of the simulation, Tyr1-N remained tightly associated with the interior of the Gla domain through hydrogen
Ca-9
a
OE2, Gla7 OE3, Gla17 OE4, Gla17 OE3, Gla21 OE4, Gla21 OE1, Gla20 OE2, Gla20 OE2, Gla21 OE3, Gla21 OE3, Gla15 OE4, Gla15 OE3, Gla20 OE4, Gla20
OE2, Gla7 OE3, Gla17 OE4, Gla17 OE3, Gla21 OE4, Gla21 OE1, Gla20 OE2, Gla20 OE2, Gla21 OE3, Gla21 OE3, Gla15 OE4, Gla15 OE3, Gla20 OE4, Gla20
OE2, Gla33 OE4, Gla33 OE2, Gla36 OE3, Gla36 OE2, Gla40 OE3, Gla40
OE2, Gla33 OE4, Gla33 OE2, Gla36 OE3, Gla36 OE2, Gla40 OE3, Gla40
OE1, Gla7 OE2, Gla7 OE3, Gla17 OE4, Gla17 OE3, Gla21 OE4, Gla21 OE1, Gla20 OE2, Gla20 OE2, Gla21 OE3, Gla15 OE4, Gla15 OE4, Gla20 O, Met19 OE2, Gla33 OE4, Gla33 OE2, Gla36 OE3, Gla36 OE2, Gla40 OE3, Gla40
OE2, Gla7 OE3, Gla17 OE4, Gla17 OE3, Gla21 OE4, Gla21 OE1, Gla20 OE2, Gla20 OE2, Gla21 OE3, Gla15 OE4, Gla15 OE3, Gla20 OE4, Gla20 OE2, Gla33 OE4, Gla33 OE2, Gla36 OE3, Gla36 OE2, Gla40 OE3, Gla40
Ca-O distance