The investigation of the effects of counterions in protein dynamics simulations

P. Drabik1,2,3, A. Liwo1, C. Czaplewski1 and J. Ciarkowski1

1 Faculty of Chemistry, University of Gdansk, ul. Sobieskiego 18, 80-952, Gdansk and 2 Department of Public Health, University School of Physical Education, ul. Wiejska 1, 80-336, Gdansk, Poland


    Abstract
 Top
 Abstract
 Introduction
 Computational methods
 Results and discussion
 References
 
Molecular simulations able to exactly represent solvated charged proteins are helpful in understanding protein dynamics, structure and function. In the present study we have used two different starting structures of papain (a typical, stable, globular protein of intermediate net charge) and different modeling procedures to evaluate some effects of counterions in simulations. A number of configurations have been generated and relaxed for each system by various combinations of constrained simulated annealing and molecular dynamics procedures, using the AMBER force field. The analysis of trajectories shows that the simulations of solvated proteins are moderately sensitive to the presence of counterions. However, this sensitivity is highly dependent on the starting model and different procedures of equilibration used. The neutralized systems tend to evince smaller root mean square deviations regardless of the system investigated and the simulation procedure used. The results of parameterized fitting of the simulated structures to the crystallographic data, giving quantitative measure of the total charge influence on the stability of various elements of the secondary structure, revealed a clear scatter of different reactions of various systems' secondary structures to counterions addition: some systems apparently were stabilized when neutralized, while the others were not. Thus, one cannot unequivocally state, despite consideration of specific simulation conditions, whether protein secondary structures are more stable when they have neutralized charges. This suggests that caution should be taken when claiming the stabilizing effect of counterions in simulations other than those involving small, unstable polypeptides or highly charged proteins.

Keywords: constrained simulated annealing/counterions/molecular dynamics/papain/protein stability


    Introduction
 Top
 Abstract
 Introduction
 Computational methods
 Results and discussion
 References
 
The dynamic behavior of proteins in solution is of fundamental interest due to the implicated relationships between protein dynamics and the structure and function of proteins (Karplus and McCammon, 1983Go; Petsko and Ringe, 1984Go). The conformational stability of a native protein structure results from a spectrum of interactions which include hydrogen bonding, van der Waals and electrostatic interactions, hydrophobic interactions and the internal freedom of motion of amino acids residues (Wada et al., 1985Go). Their individual contributions to protein stability have been extensively studied both experimentally and theoretically (Bodkin and Goodfellow, 1995Go; Honig and Yang, 1995Go; Kohn et al., 1997Go; Nolting et al., 1997Go). The ionic charges and partial atomic charges (which are defined by their Mulliken-type electron population) in proteins produce a unique electrostatic field pattern corresponding to the defined secondary structure and higher-order folding of the polypeptide chains (Wada et al., 1985Go). Electrostatic effects not only play an important role in stabilizing proteins' compact structures (Wada and Nakamura, 1981Go) but the specificity of structure recognition is also determined to a significant extent by electrostatics (Nakamura and Wada, 1985Go; Jug and Gerwens, 1998Go).

Thus, the electrostatic fields in and around biological macromolecules appear to be of central importance for both structure and function. The development of computational methods for the accurate determination of protein structure and dynamics is of both fundamental and applied importance (Zheng et al., 1993Go). However, the modeling and calculation of these fields presents a serious challenge (Gilson et al., 1987Go). While molecular simulations have proved to be very successful, there are still many questions to be answered concerning the various approximations that are commonly used in order to make the calculations computationally practical (Smith and Pettitt, 1991Go). Among molecular dynamics (MD) simulations, one of the challenges is the incorporation of the counterions into protein simulations (Ibragimova and Wade, 1998Go; Marti-Renom et al., 1998Go). While experimental knowledge is essential to understand the effects of ions on the structure and dynamic properties of proteins in solution, theoretical studies and computer simulations are necessary approaches to supplement the experimental data (Marti-Renom et al., 1998Go).

The problem of computationally representing charged biological macromolecules can be clearly exemplified in the treatment given to DNA whose conformational stability is modulated by the surrounding environment and is dependent upon the neutralization of its charged groups. Although it is still not clear if counterions allow a much better reproduction of the properties of DNA, it seems that the addition of counterions improves behavior of the system in MD simulations (Cheatham and Kollman, 1996Go; Tapia and Velazquez, 1997Go; Flatters and Lavery, 1998Go; Pastor et al., 1999Go). Also, studies were conducted to investigate counterions in simulations of not only DNA alone, but in complexes with polypeptide in protein–nucleic acid recognition research (Bishop et al., 1997Go; Kosztin et al., 1997Go; Roxstrom et al., 1998aGo,bGo). However, in the case of proteins, the value of counterion simulations has not been studied to such a great extent as for DNA. For proteins, charge-balancing counterions are positioned close to charged groups in the regions with the most favorable electrostatic potential for binding on the basis of either geometric or energetic criteria (Ibragimova and Wade, 1998Go). The addition of ions can strongly alter the thermodynamic and physical properties of proteins in solution. This is experimentally evident in the process of denaturation and salting proteins in and out of solution. The conformation of proteins is sensitive to the composition of the solution environment (Marlow et al., 1993Go).

Here, we investigate sensitivity to the placement and modeling of counterions in protein simulations. It is hoped that the latter would help in analyzing the factors which allow us to simulate a more physically correct and biologically realistic behavior. For the simulation in the present study, papain (EC 3.4.22.2) from Carica papaya was chosen. It is a globular 212-residue protein, with a net charge of +10e, three disulfide bridges, and has been investigated several times for charge distribution and related studies (Wada and Nakamura, 1981Go; Nakamura and Wada, 1985Go; Wada et al., 1985Go).


    Computational methods
 Top
 Abstract
 Introduction
 Computational methods
 Results and discussion
 References
 
All simulations were carried out using the AMBER 5.0 program (Case et al., 1997). The all-atom AMBER force field was used (Cornell et al., 1995Go).

Modeling procedures

Starting co-ordinates for all heavy atoms of two papain structures were obtained from the X-ray data [Protein Data Bank (PDB) files 1pe6 and 1pad for systems A and B, respectively]. The ligands were removed and the C25 sulfhydryl was left unmodified. Both systems were subjected to two different procedures denoted as Procedure I and II, respectively (Table IGo). The details of Procedure I are described in our previous article (Czaplewski et al., 1999Go). Two variants have been generated for this procedure via choosing a different seed for the random number generator, giving eight systems (two different structures, two different seeds and the presence of the counterions). In Procedure II, the crude systems were initially minimized in vacuo to remove close van der Waals contacts. This (and all the following) minimization(s) consisted of 20 000 steps (after 1000 steps the minimizer was switched from the steepest descent to the conjugate gradient mode). Afterwards, the system was heated in constrained simulated annealing (CSA) in vacuo from 0 to 1000 K for 1 ps, then, during thermal stabilization of 10 ps duration, snapshots were generated every 1 ps. From the two sets of 10 snapshots, three snapshots of each set (Nos 5, 8, 9 and 1, 6, 10 for systems A and B, respectively) were randomly selected for further simulations. The snapshots were subsequently cooled down for 17 ps (Table IIGo; the temperature relaxation times for particular snapshots were assigned individually to assure the slowest cooling possible). During the CSA general constraints (positional: all C{alpha} atoms in papain; all peptide bonds; and all improper dihedrals to retain proper chirality) were used. Subsequently, the systems were diluted in TIP3P water (Jorgensen et al., 1983Go), minimized, heated and equilibrated (5 ps of heating, 12 ps for water equilibration, and finally 10 ps for equilibration of the whole system; 27 ps in total) as described in Table IIIGo. The 12 ps time for water equilibration was selected, because the solvent equilibration phase should be sufficiently extensive to allow the solvent to completely readjust to the potential field of the solute. For MD this implies that the length of this phase should be longer than the relaxation time of the solvent—the time taken for a molecule to lose any ‘memory’ of its original orientation, which for water is about 10 ps (Leach, 1996Go). Thus, 12 ps should be enough for water surrounding a protein, which seems to be confirmed in our study by monitoring the physical parameters of a box (data not shown). All simulations were performed with and without adding counterions. Before surrounding the protein by water molecules, 10 Cl- charge-balancing counterions were added in positions with favorable ion–residue interactions to neutralize charges on the protein surface using the AMBER program LEAP, which adds counterions in a shell around a molecule using a Coulombic potential on a grid. Grid resolution was 1 Å. Thus, in total, 20 systems were generated for simulations: eight for Procedure I (A1-1, A1-1-I, A1-2, A1-2-I, B1-1, B1-1-I, B1-2, B1-2-I) and 12 for Procedure II (A2-5, A2-5-I, A2-8, A2-8-I, A2-9, A2-9-I, B2-1, B2-1-I, B2-6, B2-6-I, B2-10, B2-10-I). The first number denotes Procedure number, the second number denotes two different seeds in Procedure I and snapshot number in Procedure II, and I denotes inclusion of the counterions in the simulation.


View this table:
[in this window]
[in a new window]
 
Table I. The general protocols of Procedures I and II
 

View this table:
[in this window]
[in a new window]
 
Table II. The protocol for cooling (generating snapshots via simulated annealing in vacuo)
 

View this table:
[in this window]
[in a new window]
 
Table III. The protocol for heating
 
Molecular dynamics productive runs

Potentially equilibrated systems were subjected to MD productive runs. Simulations of complexes in solution were performed under periodic boundary conditions in a closed, isothermal, isobaric (NTP) ensemble. Throughout the simulation the solute and solvent were coupled to a constant-temperature (T = 300 K) heat bath and a constant-pressure (P = 1 bar) bath (Berendsen et al., 1984Go). All hydrogen-containing bonds were constrained using the SHAKE algorithm (Ryckaert et al., 1977Go) for holonomic constraints with a relative geometric tolerance for co-ordinate resetting of 0.0005 Å, allowing a time step of 1 fs. The Leapfrog version (Hockney, 1970Go) of the Verlet algorithm (Verlet, 1967Go) was employed to integrate the equations of motion. A double residue-based cutoff distance of 10/14 Å was used for non-bonded interactions. The TIP3P model was used for water molecules. A typical rectangular box size was 78x66x60 Å. Approximately 8200 TIP3P water molecules were in the box, i.e. 28 000 atoms in total. Co-ordinates were saved every 1000 steps. A set of 230 and 300 ps MD runs was carried out during Procedures I and II, respectively.

Parametric mean square approximation

To quantitatively compare time-averaged C{alpha} atoms mobilities as a function of sequence in dependence on simulation conditions, the minimum of {varepsilon}, where {varepsilon} is a function of one variable (a), was sought:

(1)
where {varepsilon} is a measure of similarity between x and y2), {x1, . . . , xn} are the reference data: crystallographic temperature factors, {y1, . . . , yn} are the simulation data: time-averaged residue-based mobilities, n is the number of residues and a is a parameter.

Parameter a was introduced because average C{alpha} mobilities differ from X-ray temperature factors by approximately one order of magnitude. Parameter a was calculated to minimize {varepsilon}:

(2)

Statistical analysis

Statistical calculations were performed using the STATISTICA 5.0 program (StatSoft, Inc., 1995). One-way analysis of variance (ANOVA) of repeated measurements (two levels of within-subjects factor: ions versus non-ions) was used to determine the differences in MD simulation behavior with and without counterions. Normality of the distribution of the data was assessed by a Kolmogorow–Smirnov test with Lilliefors probabilities (Lilliefors, 1967Go) and Shapiro–Wilks analysis (Shapiro et al., 1968Go) with Royston's algorithm (Royston, 1982Go). Homogeneity of variance was evaluated by means of Levene's test. Results of normality and homogeneity of variance tests allowed the usage of ANOVA statistics (data not shown). Statistical significance was defined at the p < 0.05 level.


    Results and discussion
 Top
 Abstract
 Introduction
 Computational methods
 Results and discussion
 References
 
Computer simulations are an approach, under continual development, for studying protein structure, dynamics and stability. In the present study, the conformational space was explored using MD with CSA. In CSA very high temperatures are used in the initial stages to permit the system to range widely over the energy surface. The temperature is then gradually reduced as the structure settles into a conformation which has low energy. The production of long and stable trajectories using MD is necessary to obtain interpretable results for such studies.

Experimentally, it has been shown elsewhere that the increase in the concentration of added salts leads to significant variations in the stability of proteins, usually to stabilization of non-chaotropic salts, probably owing to changes in the apparent hydrophobic effect and charge screening (Kohn et al., 1997Go). In the theoretical studies, improved protein stability has been observed several times when adding counterions (York et al., 1993Go; Yelle et al., 1995Go) which screen charged side chains, preventing unnatural electrostatic interactions between neighboring parts of the macromolecule (Ibragimova et al., 1998). A similar conclusion was reached from the studies of the highly charged activation domain of procarboxypeptidase B, which suggests that the use of counterions significantly contribute to the maintenance of native secondary structure (Marti-Renom et al., 1998Go). Ibragimova and Wade (Ibragimova and Wade, 1998Go) investigated counterion effects in simulations of the WW domain of YAP (Yes kinase-associated protein). The small size and weak stability of this domain make it very sensitive to conditions in MD simulations, and thus, in the opinion of the authors, they constitute a good system on which different protocols and models could be tested. The usual procedure of only placing charge-balancing counterions near charged groups on the protein surface did not result in a stable protein fold, nor did simulation without any counterions. For maintaining the stability of the WW domain, explicit modeling of full ionic strength (at 0.2 M) was necessary.

As it may seem quite natural that adding ions would help to stabilize highly charged polypeptides or proteins with marginal stability, the aim of our study was to investigate if such a procedure influences a typical protein of the type used widely in MD research, it not being a sophisticated example. Thus, we have used a different approach in choosing the model to simulate. Papain is a typical 212-residue globular protein and, when well equilibrated, is very stable during long MD simulations giving relaxed stable trajectories (data not shown). Four systems (two different starting structures and two different modeling procedures) in two variants each (with and without the counterions) were generated and investigated to maximize different simulation conditions and, in that way, to eliminate (at least partly) the modeling setup.

As MD simulations are chaotic and critically dependent on the starting conditions (Braxenthaler et al., 1997Go), to strengthen the statistical inference 10 sets of different runs (with a total time of about 5.5 ns) have been performed with differences in seeds for the random number generator, investigated systems and equilibration procedures. Figure 1Go shows the root mean square deviations (RMSD) for all simulated systems. Subjected to Procedure I non-neutralized systems failed to equilibrate, the RMSD values tended to increase until the very end of the MD run reaching 5 Å in the B1-1 simulation. The influence of added counterions on these systems is somewhat ambiguous. In system A ions efficiently helped the protein to equilibrate, drastically reducing the overall RMSD (Table IVGo). Surprisingly, the system B reacts on the same procedure in a different fashion: regardless of net charge, the system fluctuations rapidly increased throughout the simulation, although in the case of neutralized simulations they were much smaller in size (Table IVGo). The equilibration protocol used in Procedure I apparently gives no steady-state phase (in terms of RMSD). However, it was used successfully in many MD simulations performed on similar systems (see for example: Kazmierkiewicz et al., 1997Go; Czaplewski et al., 1999Go). It should be stressed that although there is an insufficiency in the applied simulation protocol, it still gives some insight into the problem of the addition of counterions.




View larger version (30K):
[in this window]
[in a new window]
 
Fig. 1. The RMSDs (Å) versus simulation time (ps) of C{alpha} atoms along the MD trajectories. Symbols of MD runs are as defined in the Computational methods section.

 

View this table:
[in this window]
[in a new window]
 
Table IV. The one-way ANOVA results for individual simulations
 
Analyzing Procedure II one can see that both systems look almost perfectly relaxed from the very beginning of the MD productive runs. The addition of counterions does not change this situation significantly. Although not influencing the curvature of the RMSD trajectory, the overall RMSD is lower (giving statistically significant difference in most of the simulations) in neutralized systems (Table IVGo). The one-way ANOVA calculations have been performed also on populations of runs (denoted X1 and X2 for collective sets of simulations in Procedures I and II, respectively). The differences between families of MD runs with and without the counterions turned out to be statistically significant, both in Procedure I and in Procedure II (Table IVGo). This application of ANOVA to compare not single trajectories, but families of them, gives more reliable information about the role of counterions of the structure.

It seems reasonable that simulations with counterions should behave better; they render more accurately biophysical conditions (solvent composition) in living organisms. When performing simulations without counterions, strong electrostatic interactions between charged side chains could disturb the structure of the protein during the simulation.

The secondary structure was examined to characterize conformational changes of the papain during simulations. The results of parameter fitting (see Computational methods section) of the calculated structures to the crystallographic data, giving quantitative measure of the total charge influence on the stability of various elements of the secondary structure, are shown in Table VGo. The scatter of different responses of various systems' secondary structures to counterions addition is clearly seen. It is hardly possible to find any order in these responses: some (6 of 10) systems apparently are stabilized by ions in the regions of both regular ({alpha}-helices, ß-sheets and turns) and irregular (mostly loops) secondary structures (A1-1, A1-2, B1-1, B1-2, A2-5 and B2-6 simulations), the A2-8 system seems to be stabilized by counterions in all residues except those forming {alpha}-helices, whereas neutralization driven improvement of A2-9 and B2-1 systems is clear only in mobile (irregular) residues. The B2-10 system is stabilized only in terms of turn regions. It seems that there is no significant superiority of counterions usage on the basis of the analysis of the four pairs of MD runs presented. It is rather hard to evaluate counterions effects on the local fluctuations of molecule mobility because this strongly depends on the details of specific simulations that were performed.


View this table:
[in this window]
[in a new window]
 
Table V. The measures of similarities (as defined in Computational methods section) between secondary structures in individual simulations versus experimental temperature factors for X-ray structure (PDB entry: 1pe6)
 
In the present study we have used two different starting structures of papain (a typical, stable, globular protein with three disulfide bridges) and two different modeling procedures. Our results show that counterion effects highly depended on the model variant and procedure used. It seems that one cannot unequivocally state, despite consideration of specific simulation conditions, whether protein secondary structures are more stable when they have neutralized charges. This suggests that caution should be taken when claiming the stabilizing effect of counterions in simulations other than those involving small, unstable polypeptides (Ibragimova and Wade, 1998Go) or highly charged proteins lacking disulfide bonds (Marti-Renom et al., 1998Go). Thus, the problem of neutralization of charged proteins needs further investigation that will enhance our understanding of protein structure and function.


    Notes
 
3 To whom correspondence should be addressed. E-mail: drabol{at}chemik.chem.univ.gda.pl Back


    Acknowledgments
 
Supported by grant DS 8245-4-0124-1 from the Polish State Committee for Scientific Research (KBN). The computations were carried out at the Academic Computer Center in Gdansk (TASK), Poland, and the Interdisciplinary Center for Mathematical Modeling (ICM) in Warsaw, Poland.


    References
 Top
 Abstract
 Introduction
 Computational methods
 Results and discussion
 References
 
Berendsen,H.J.C., Postma,J.P.M., van Gunsteren,W.F., DiNola,A. and Haak,J.R. (1984) J. Chem. Phys., 81, 3684–3690.[ISI]

Bishop,T.C., Kosztin,D. and Schulten,K. (1997) Biophys. J., 72, 2056–2067.[Abstract]

Bodkin,M.J. and Goodfellow,J.M. (1995) Protein Sci., 4, 603–612.[Abstract/Free Full Text]

Braxenthaler,M., Unger,R., Auerbach,D., Given,J.A. and Moult,J. (1997) Proteins, 29, 417–425.[ISI][Medline]

Case,D.A. et al. (1997) Amber 5.0. University of California, San Francisco.

Cheatham,T.E. and Kollman,P.A. (1996) J. Mol. Biol., 259, 434–444.[ISI][Medline]

Cornell,W.D., Cieplak,P., Bayly,C.I., Gould,I.R., Merz,K.M.,Jr, Ferguson,D.M., Spellmeyer,D.C., Fox,T., Caldwell,J.W. and Kollman,P.A. (1995) J. Am. Chem. Soc., 117, 5179–5197.[ISI]

Czaplewski,C., Grzonka,Z., Jaskólski,M., Kasprzykowski,F., Kozak,M., Politowska,E. and Ciarkowski,J. (1999) Biochim. Biophys. Acta, 1431, 290–305.[ISI][Medline]

Flatters,D. and Lavery,R. (1998) Biophys. J., 75, 372–381.[Abstract/Free Full Text]

Gilson,M.K., Sharp,K.A. and Honig,B.H. (1987) J. Comput. Chem., 9, 327–335.[ISI]

Hockney,R.W. (1970) Methods Comput. Phys., 9, 136–211.

Honig,B. and Yang,A.-S. (1995) Adv. Protein Chem., 46, 27–58.[ISI][Medline]

Ibragimova,G.T. and Wade,R.C. (1998) Biophys. J., 74, 2906–2911.[Abstract/Free Full Text]

Jorgensen,W.L., Chandreskhar,J., Madura,J.D., Imprey,R.W. and Klein,M.L. (1983) J. Chem. Phys., 79, 926–935.[ISI]

Jug,K. and Gerwens,H. (1998) J. Phys. Chem. B, 102, 5217–5227.[ISI]

Karplus,M. and McCammon,J.A. (1983) Annu. Rev. Biochem., 52, 263–300.[ISI][Medline]

Kazmierkiewicz,R., Czaplewski,C., Ciarkowski,J. (1997) Acta Biochim. Pol., 44, 453–466.[ISI][Medline]

Kohn,W.D., Kay,C.M. and Hodges,R.S. (1997) J. Mol. Biol., 267, 1039–1052.[ISI][Medline]

Kosztin,D., Bishop,T.C. and Schulten,K. (1997) Biophys. J., 73, 557–570.[Abstract]

Leach,A.R. (1996) Molecular Modelling. Principles and Applications. Longman Singapore Publishers (Pte) Ltd, Singapore, p. 326.

Lilliefors,H. (1967) J. Am. Stat. Assoc., 64, 399–402.

Marlow,G.E., Perkyns,J.S. and Pettitt,B.M. (1993) Chem. Rev., 93, 2503–2521.[ISI]

Marti-Renom,M.A., Mas,J.M., Oliva,B., Querol,E. and Aviles,F.X. (1998) Protein Eng., 11, 881–890.[Abstract]

Nakamura,H. and Wada,A. (1985) J. Phys. Soc. Jpn, 54, 4047–4052.[ISI]

Nolting,B., Golbik,R., Neira,J.L., Soler-Gonzalez,A.S., Schreiberg,G. and Ferscht,A. (1997) Proc. Natl Acad. Sci. USA, 94, 826–830.[Abstract/Free Full Text]

Pastor,N., MacKerell,A.D. and Weinstein,H. (1999) J. Biomol. Struct. Dynam., 16, 787–810.[ISI][Medline]

Petsko,G.A. and Ringe,D. (1984) Annu. Rev. Biophys. Bioeng., 13, 331–371.[ISI][Medline]

Roxstrom,G., Velazquez,I., Paulino,M. and Tapia,O. (1998a) J. Phys. Chem. B., 102, 1828–1832.[ISI]

Roxstrom,G., Velazquez,I., Paulino,M. and Tapia,O. (1998b) J. Biomol. Struct. Dynam., 16, 301–312.[ISI][Medline]

Royston,J.P. (1982) Appl. Stat., 31, 115–124.[ISI]

Ryckaert,J.P., Ciccotti,G. and Berendsen,H.J.C. (1977) J. Comput. Phys., 23, 327–341.[ISI]

Shapiro,S.S., Wilk,M.B. and Chen,H.J. (1968) J. Am. Stat. Assoc., 63, 1343–1372.[ISI]

Smith,P.E. and Pettitt,B.M. (1991) J. Chem. Phys., 95, 8430–8441.[ISI]

StatSoft, Inc. (1995) STATISTICA for Windows (Computer program manual). StatSoft, Inc., Tusla, OK.

Tapia,O. and Velazquez,I. (1997) J. Am. Chem. Soc., 119, 5934–5938.[ISI]

Verlet,L. (1967) Phys. Rev., 159, 98–103.[ISI]

Wada,A. and Nakamura,H. (1981) Nature, 293, 757–758.[ISI][Medline]

Wada,A., Nakamura,H. and Sakamoto,T. (1985) J. Phys. Soc. Jpn, 54, 4042–4046.[ISI]

Yelle,R.B., Park,N.-S. and Ichiye,T. (1995) Proteins, 22, 154–167.[ISI][Medline]

York,D.M., Darden,T.A., Pedersen,L.G. and Anderson,M.W. (1993) Biochemistry, 32, 1443–1453.[ISI][Medline]

Zheng,Q., Rosenfeld,R., Vajda S. and DeLisi,C. (1993) Protein Sci., 2, 1242–1248.[Abstract/Free Full Text]

Received June 20, 2001; revised July 6, 2001; accepted July 10, 2001.





This Article
Abstract
FREE Full Text (PDF)
Alert me when this article is cited
Alert me if a correction is posted
Services
Email this article to a friend
Similar articles in this journal
Similar articles in ISI Web of Science
Similar articles in PubMed
Alert me to new issues of the journal
Add to My Personal Archive
Download to citation manager
Search for citing articles in:
ISI Web of Science (1)
Request Permissions
Google Scholar
Articles by Drabik, P.
Articles by Ciarkowski, J.
PubMed
PubMed Citation
Articles by Drabik, P.
Articles by Ciarkowski, J.