158

Hayward

frequency modes that result from projecting the Hessian into the space defined by intrablock motions. It has been shown that the method scales with N2, where N is the total number of degrees of freedom.

Finally, any symmetry that a molecule possesses can be exploited by the methods of group theory to reduce the Hessian to a number of independent submatrices [20]. If the symmetry group to which the molecule belongs has been determined, then the procedures of group theory show how one can construct a basis (the basis for the irreducible representations) in which the Hessian comprises a number of smaller independent submatrices lined along the diagonal. These independent submatrices can be individually diagonalized and the eigenvectors and eigenvalues of the whole system reconstructed. This process gives a considerable saving on the diagonalization of the whole Hessian. For example, if the molecule has a tenfold rotational symmetry (i.e., it belongs to the group C10), the order of the blocks is one tenth that of the whole Hessian. Group theory has been applied to normal mode analyses of the gramicidin A dimer [21] and one layer of the tobacco mosaic virus protein disk [22]. The whole problem of performing normal mode analysis on symmetrical protein assemblies in dihedral angle space has been tackled by Gibrat et al. [14]. This promises to make viable the normal mode analyses of very large molecular assemblies such as those of virus capsids and crystals.

B. Normal Mode Analysis in Dihedral Angle Space

For most large biological molecules, most of the flexibility arises from torsional or dihedral angle rotations. In comparison, bond lengths and bond angles are comparatively rigid. This automatically leads one to the idea of performing a normal mode analysis in dihedral angle space. In such an analysis the bond lengths and angles are kept fixed at their minimum energy values. At least for bond length fluctuations, the fact that they correspond to a permanent ground quantum oscillator state at physiological temperatures suggests that constraining them is not unphysical. For proteins the ratio of the number of rotatable dihedral angles to Cartesian coordinates is about 1:8 in proteins and about 1:11 for nucleic acids [23]. This provides a considerable saving, but there is the practical disadvantage that a normal mode analysis in dihedral angle space is technically more involved.

The dihedral angle space procedure is complicated by the fact that the kinetic energy term cannot be expressed as a simple function of the dihedral angle variables as it can in the case of the mass-weighted Cartesian coordinates. There are six degrees of freedom for the overall rigid-body motion of the molecule, and formulating the kinetic energy in a way that ensures that changes in dihedral angles do not result in the movement of the center of mass of the molecule or an overall rotation (by use of the Eckart conditions) means that the equations of motion can be solved for the internal motion alone. This is usually done by calculating the so-called K matrix, which is the linear term in a Taylor expansion,

Kij	∂qi	(10)
	∂φj

where qi is the ith mass-weighted Cartesian coordinate, and φj the jth dihedral angle. The derivative is calculated at the energy minimum, subject to the six constraints [24]. The kinetic energy is then expressed by way of the mass tensor, often denoted H, as

T	1		˙ t	˙	(11)
			φ	Hφ

2

Normal Mode Analysis of Biological Molecules

159

where φ˙ is the column vector whose elements are the time derivatives of the dihedral angles, and H Kt K. The solution of the equations of motion proceeds in the same way as for the Cartesian coordinates, except that one must solve a generalized eigenvalue problem whereby the eigenvector matrix V now simultaneously diagonalizes H to the identity matrix I and the Hessian to the eigenvalue matrix Ω. To convert between the normal mode coordinates, which now describe collective variations in dihedral angles, and the mass-weighted Cartesian coordinates, one can use the matrix

W KV

(12)

where WWt I.

As a dihedral angle space analysis is equivalent to a Cartesian coordinate space analysis where the bond angles and lengths are kept fixed, one can compare the results from a conventional Cartesian coordinate space analysis to those from a dihedral angle space analysis to access the effects of bond angle and bond length variations [25,26]. Although low frequency modes can be represented solely by dihedral angle variations, allowing bond length and bond angle variations has the indirect effect of making their amplitudes of fluctuation larger. A study of the protein bovine pancreatic trypsin inhibitor (BPTI) showed that the important subspaces (the subspaces within which most of the fluctuation occurs) overlap considerably, implying that dihedral angle space normal mode analysis is a viable alternative to Cartesian coordinate space analysis [26]. In fact, it has been argued that if one includes second-order terms in the conversion from dihedral angle fluctuations to Cartesian coordinate fluctuations, then the results from a dihedral angle space analysis are, in fact, better than those from a Cartesian coordinate space analysis because the harmonic approximation is valid over a wider range in dihedral angle space than in Cartesian coordinate space [27]. This approach has been used to calculate nuclear magnetic resonance (NMR) order parameters in proteins [28].

Dihedral angle space normal mode analyses of a number of proteins (not all of them large by current standards) have been made, including BPTI [1,2,26,27,29–31], lysozyme [30,32–35], G-actin [36], myoglobin [37–39], epidermal growth factor [40], bacteriophage 434 Cro and 434 repressor [41], and subtilisin–Eglin c complex [42]; nucleic acids, including tRNA [43]; and a double-stranded DNA dodecamer [23]. In the latter, flexibility arising from the pseudorotation of the furanose rings was included [44].

C. Approximate Methods

The methods described so far are, within the framework of the analysis, accurate, in the sense that they are implemented by using basically the same detailed force fields for the bonded and non-bonded interactions that are used in molecular dynamics simulations. Tirion [45] showed that a ‘‘single parameter model’’ can give results similar to those of these detailed analyses. The model replaces the Hessian by a matrix whose elements are zero for any pair of atoms separated by a distance greater than a cutoff distance (equal to the sum of the van der Waals radii plus a distance parameter Rc that models the decay of the interaction) and have values according to a simple Hookean pairwise potential with the same force constant C for all pairs of atoms within the cutoff distance of each other. By adjusting C and Rc, the predictions from conventional normal mode analyses and this simplified model are found to be astoundingly similar. Interestingly, for good fits to the conventional normal mode analysis data, the product CR2c has a constant value, indicating

160

Hayward

a universal ‘‘bond strength’’ of 3 kJ/mol. As distances between atoms are calculated directly from their X-ray-determined coordinates, costly energy minimization is avoided.

Bahar et al. [46] have used this kind of approach to predict the B-factors of 12 X- ray structures. Elements in the ‘‘Hessian’’ corresponding to atom pairs separated by a

˚

distance of less than 7 A are set to zero, and the remainder have the same value dependent on a single adjustable parameter. Generally B-factor predictions for the α-carbons compare very well with the B-factors measured by X-ray crystallography. Figure 1 shows the result for the subunit A of endodeoxyribonuclease I complexed with actin.

Given that one does not need to perform an energy minimization and that the ‘‘Hessian’’ is very sparse, it is not surprising that the computation time is reported to be at least one order of magnitude less than for a conventional normal mode analysis.

IV. NORMAL MODE REFINEMENT

The normal mode refinement method is based on the idea of the normal mode important subspace. That is, there exists a subspace of considerably lower dimension than 3N, within which most of the fluctuation of the molecule undergoing the experiment occurs, and a number of the low frequency normal mode eigenvectors span this same subspace. In its application to X-ray diffraction data, it was developed by Kidera et al. [33] and Kidera and Go [47,48] and independently by Diamond [49]. Brueschweiler and Case [50] applied it to NMR data.

Figure 1 Experimental B-factors of α-carbon atoms (thin curve) compared with those predicted using their single-parameter model for subunit A of endodeoxyribonuclease I complexed with actin (thick curve). Calculations were performed using both subunits A and D comprising 633 residues. Their parameter is adjusted such that the area under the predicted curve equals the area under the experimental curve. (From Ref. 46.)

Normal Mode Analysis of Biological Molecules

161

A. Normal Mode X-Ray Refinement

In the procedure of X-ray refinement, the positions of the atoms and their fluctuations appear as parameters in the structure factor. These parameters are varied to match the experimentally determined structure factor. The term pertaining to the fluctuations is the Debye–Waller factor in which the atomic fluctuations are represented by the atomic distribution tensor:

∆xj ∆xj	∆xj ∆yj
∆yj ∆xj	∆yj ∆yj
Uj ∆zj ∆xj	∆zj ∆yj

∆xj ∆zj

∆yj ∆zj (13)∆zj ∆zj

where j is the atomic index. In a common method of refinement it is often assumed that atomic fluctuations are isotropic and independent. In such a case the atomic distribution tensor is a scalar quantity, the B-factor. It is the B-factors that are varied in the refinement process. The normal mode refinement method requires two refinement procedures. The first uses the isotropic B-factors to refine a structure well enough for a normal mode analysis to be performed. Then the process of normal mode refinement begins by expressing the atomic distribution tensor in terms of the fluctuations of a selected number, M, of the lowest frequency normal modes through

Uj Wj Λref Wjt

(14)

This equation is a variant of Eq. (9b), but here Wj is the 3 M portion of W corresponding to atom j and the first M normal mode eigenvectors, and Λref is the M M variance– covariance matrix of the M lowest frequency normal mode coordinates. Note that Λref would equal the corresponding portion of the diagonal matrix kB TΩ 1 from Eq. (9b) if the dynamics of the real molecule in the experiment were perfectly described by the normal mode analysis. In the real situation this is not the case, so the elements of Λref are the parameters to be varied in the refinement process. If M is small enough that M M is less than N, the number of isotropic B-factors, there can be a significant saving in the number of parameters and, as has been reported, an improvement in the R-factor [47] and Rfree-factor [51]. If the six modes describing external motion are included, then simulation showed that it could distinguish between internal and external fluctuations [48]. Perhaps the most significant advantage is that, unlike the normal atomic B-factor refinement methods, correlations between atoms can be determined. Kidera and Go, together with Inaka and Matsushima, tested it on real diffraction data from human lysozyme [33]. Figure 2 shows the result of a comparison of the internal root-mean-square fluctuations (RMSFs) of human lysozyme residues determined from the normal mode refinement procedure in comparison with the results from a normal mode analysis on that molecule. The correspondence between the two is remarkable and unexpected. Diamond [49] presented a similar analysis on the protein BPTI.

B. Normal Mode NMR Refinement

The normal mode NMR refinement method of Brueschweiler and Case [50] can be applied to experimentally measurable quantities such as order parameters or nuclear Overhauser spectroscopy (NOSEY) intensities. Unlike the X-ray case, the expression of these quanti-