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)[52], 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[50].
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)[52] 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)[51], 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.[55]. Energies were calculated using the DREIDING forcefield[54]. NEIMO dynamics calculations were performed using software written at JPL and Caltech, and interfaced to the BIOGRAF program.