In an attempt to simulate the expansion effect, two different models of the TBSV asymmetric unit were developed. The first contains the protein atoms and calcium ions as they appear in the protein database coordinate file (2TBV), with hydrogen atoms added to nitrogen, oxygen and sulfur atoms to allow for hydrogen bonding. In addition, Na and Cl ions were added to balance the charges of unpaired acidic and basic residues, respectively. This structure is termed the ``PH7'' model. The second representation is the ``NOCA'' model, in which the six Ca ions were removed and the free aspartic acid residues were allowed to form salt bridges with basic residues, or were given Na counterions. In this model, the 15 Asp residues are no longer held together by interactions with Ca, but are free to move independently. This is believed to be the major factor in the expansion of the virus particle.
Inclusion of explicit waters in a calculation of this type can improve its accuracy but can also greatly expand the computational cost. We do not, therefore, have waters included in the system. However, we are able to compensate for this partially in two ways: 1) by using a distance-dependent dipole, and 2) by using counterions to balance lone charges. The first mimics the charge-shielding capacity of water by using an effective dielectric constant between charges and proportional to the distance between the two charges,
Inclusion of counterions in the simulation can provide charge stabilization for lone charged groups. In nature, such stabilization is provided both by water dipoles and free ions. The 893 residues of the simulated triad include 69 acidic (48 aspartic acid and 21 glutamic acid) and 58 basic (38 arginine and 20 lysine) residues. Of the 48 aspartic acid residues, 15 are involved in the binding of the 12 Ca ions. At each protein interface, two calcium ions and five aspartic acids are present; at the A/B interface Asps 181, 183, and 186 from A are present along with Asps 153 and 225 from B. A basic residue, Lys 232 from A, is also present, approximately 3.1 Å from Asp 183. Excluding these six residues per protein, there were a total of 109 free charged residues.
Not every charged residue requires a counterion, since many are involved in salt bridges. We eliminated (+,-) paired sidechains by looking at all pairwise distances between the central atoms of oppositely charge sidechains. These atoms are C in aspartic acid, C in glutamic acid, C in arginine, and N in lysine. We checked not only for (+,-) pairs within the three proteins of the asymmetric unit but expanded the viral coat into its entire 180 proteins and checked for pairs between residues in different triads. Pairs less than 10 Å apart were considered to stabilize one another, and were not given countercharges. There were 30 such pairs in addition to the Lys 232-Asp 183 pair mentioned in the preceding paragraph. Of these, 23 occurred within the asymmetric unit and 10 occurred between residues of different triads. As can be seen in Figure , each P domain is closely paired with a P domain from a different triad, and several salt bridges occur in this region. After eliminating charge-paired residues, 25 Cl were still needed to balance 10 arginines and 15 lysines while 24 Na were needed for 18 aspartic acids and 6 glutamic acids. These counterions were placed in idealized locations, determined by previous calculations on individual sidechains. In the NOCA model, the Ca ions were removed, and the need for counterions was recalculated. In this model, the Asp 53 sidechains from B formed weak (+,-) pairs with Lys 182 of B, so counterions were not needed for these residues. This was the case at all three interfaces. Therefore, the NOCA model required a total of 33 Na counterions and 22 Cl and there was no net change in the total number of atoms in the system.
The asymmetric unit of the TBSV viral coat contains three copies of the coat protein in slightly different conformations. The S (residues 102-274) and P (275-387) domains are resolved crystallographically in all three conformations. In each case, the two domains were built independently, so mismatches exist in the hinge region (residues 273-275). The crystal structure (2TBV) lists alternate S and P coordinates for the residues in the overlap region. For these calculations, the two alternates were averaged and re-optimized by energy minimization. The R domains are not resolved, except for residues 67-101 in the C conformation. Only these are included in the calculations. The RNA is also not included in the current calculations, since no structure is available for it. The simulated PH7 system contains 893 protein residues having a total of 8083 atoms in addition to six Ca, 24 Na and 25 Cl for a total of 8138 atoms. As explained above, the NOCA model has zero Ca but 33 Na and 22 Cl, so the total number of atoms remains the same.
In order to accurately model the capsid environment, the nonbonded forces acting on the asymmetric unit included interactions with the other 177 proteins. This was made possible by the Cell-Multipole Method (CMM), an extremely fast and accurate method for calculating nonbonds in large systems. CMM divides the simulation space into a hierarchy of cubic cells, the smallest of which contains, ideally, 4 or 5 atoms and the largest of which contains the entire system. Each cell at level contains eight cells from level . For the triad alone, four levels are used. There are 4096 () cells at this level, measuring 6.397 Å on a side. Since the system is much flatter than it is cubic, 81.4%of these cells are empty, leaving 762 populated cells with an average of 10.7 atoms per cell. When the triad is expanded into the full 180 protein capsid, the cell multipole method uses six levels for the 488,820 atom system. At level 6, there are 262,144 cells measuring 5.340 Å on a side. 87.5%of these are empty, leaving 8.0 atoms per populated cell. Both the dimension and the average population in this case are better than those in simulations of the triad, alone.
In CMM calculations, the total charge, dipole, and quadrupole are calculated for each cell at each level. Exact pairwise interactions, like those in Equation () are only calculated for atoms in adjacent cells. Interactions with distant cells are calculated as interactions with the cumulative charge, dipole and multipole of those cells. For nearby cells, just beyond the nearest neighbors, interactions are calculated with the lowest-level cells, e.g., level six cells in the capsid simulation. Interactions at increasing distance are calculated with larger, higher-level cells, up to level 1. This hierarchical approach makes the nonbond calculation proportional to , the number of atoms, rather than , as would be the case for traditional methods, yet is highly accurate because the errors introduced are proportional to the strength of the interaction.
All calculations were performed on one processor of a Silicon Graphics 4D/380 workstation using the BIOGRAF software from Molecular Simulations, Inc.. Energies were calculated using the DREIDING forcefield. NEIMO dynamics calculations were performed using software written at JPL and Caltech, and interfaced to the BIOGRAF program.