Research Article  Open Access
Toward Structure Prediction for Short Peptides Using the Improved SAAP Force Field Parameters
Abstract
Based on the observation that Ramachandrantype potential energy surfaces of single amino acid units in water are in good agreement with statistical structures of the corresponding amino acid residues in proteins, we recently developed a new allatom force field called SAAP, in which the total energy function for a polypeptide is expressed basically as a sum of single amino acid potentials and electrostatic and LennardJones potentials between the amino acid units. In this study, the SAAP force field (SAAPFF) parameters were improved, and classical canonical Monte Carlo (MC) simulation was carried out for short peptide models, that is, Metenkephalin and chignolin, at 300 K in an implicit water model. Diverse structures were reasonably obtained for Metenkephalin, while three folded structures, one of which corresponds to a nativelike structure with three native hydrogen bonds, were obtained for chignolin. The results suggested that the SAAPMC method is useful for conformational sampling for the short peptides. A protocol of SAAPMC simulation followed by structural clustering and examination of the obtained structures by ab initio calculation or simply by the number of the hydrogen bonds (or the hardness) was demonstrated to be an effective strategy toward structure prediction for short peptide molecules.
1. Introduction
In conformational analysis of short peptides, Monte Carlo (MC) and molecular dynamics (MD) simulation techniques have been widely applied [1, 2], in which a set of potential energy functions (a socalled force field), such as ECEPP [3, 4], AMBER [5–8], CHARMM [9], OPLS [10, 11], and GROMOS [12], is employed to define the relation between the structure and the potential energy. However, it is still not practical to search for the possible conformers comprehensively within a limited computing time without using advanced sampling methods, such as multicanonical [13, 14] and replicaexchange [15–17] methods, which have been developed to overcome the sampling problem. In recent years, various types of coarsegrained force fields, such as AMBERUA [18], UNRES [19–23], and MARTINI [24, 25], have also been developed in order to reduce time to compute the potential energy. For instance, the UNRES force field successfully simulated the aggregation process and stable structure of βamyloid protein oligomers that are regarded as pathologic factors of Alzheimer’s disease [26].
On the other hand, we recently discovered the interesting feature of protein structures that Ramachandrantype potential energy surfaces of single amino acid units in water obtained by ab initio molecular orbital calculation are almost identical with statistical distributions of the Ramachandran plots obtained from protein databank (PDB) [27, 28]. In other words, the amino acid residues in proteins seem to statistically follow Boltzmann distributions on the single amino acid potential (SAAP) surfaces in water. Although physicochemical implications of this unexpected similarity are not yet clear, the finding strongly suggests prominent importance of SAAP as a determinant of protein structures. This point would also be supported by the previous experimental result by Dobson and coworkers [29, 30] that and dihedral angle distributions of individual amino acid residues of a polypeptide in a randomcoil state can be expressed by use of the Ramachandran plots, as well as by the theoretical approaches by Sakae and Okamoto [31] and Kamiya et al. [32] to the optimization of conventional allatom force field parameters by using the Ramachandran plots.
A discovery of the similarity between the SAAP in water and the statistical structure of the amino acid residues in folded proteins has prompted us to develop a new force field called SAAP for polypeptide molecules [33, 34]. The SAAP force field (SAAPFF) is entirely different from conventional force fields in that a polypeptide is divided to the amino acid units, not to the atomic units, solvent effects are implicitly included in the parameters, and the atomic charges are not constant but variable depending on the conformation of the amino acid units. In the SAAPFF, the total potential energy () for a polypeptide is expressed by (1) as the following: where is a sum of potential energies (SAAP) for the individual amino acid units, and are electrostatic and LennardJones potentials between the amino acid units, and is the other correlation term, which is ignored in a current version of SAAPFF. Detailed description of these terms and the simulation program was given previously (see also Supplementary Material available online at http://dx.doi.org/10.1155/2013/407862) [33, 34]. Canonical MC simulation using the SAAPFF reasonably reproduced randomly fluctuating conformation of Metenkephalin. However, for chignolin the SAAPMC simulation did not afford the native βhairpin structure when the simulation was started from the extended structure [34]. Thus, the previous SAAPFF parameters were a little too rough to be applied for structure prediction of the short peptide. In this study, the SAAPFF parameters have been improved in several points, and accuracy of the improved SAAP simulation suite to search for possible conformers of short peptides has been evaluated by performing the molecular simulation for Metenkephalin and chignolin. The results were then compared with those obtained by the conventional AMBERMD simulation method using the amber99sb force field [1] combined with the generalized Born solvent model [35–37]. A protocol of SAAPMC simulation followed by structural clustering is suggested to be an effective strategy toward structure prediction for short peptide molecules.
2. Method of Calculations
2.1. Model Peptides
Metenkephalin and chignolin without N and Cterminal protecting groups were employed as model short peptides. Metenkephalin is a short peptide, which consists of five residues (Tyr–Gly–Gly–Phe–Met) as shown in Figure 1. According to previous experimental [38] and theoretical [39–50] studies, the structure of this short peptide randomly fluctuates in water, while the presence of some characteristic conformations with intramolecular hydrogen bonds was indicated in vacuo [51–53]. Metenkephalin is a peptide widely used for testing performance of developed molecular simulation methodologies.
Chignolin has a unique folded βhairpin structure in water, although it is a short peptide consisting of ten amino acid residues (Gly–Tyr–Asp–Pro–Glu–Thr–Gly–Thr–Trp–Gly) [54]. Molecular simulation for this miniprotein has been performed actively by using AMBER and OPLS force fields [55–64]. As shown in Figure 2, chignolin folds around Pro4–Glu5–Thr6 to form a βhairpin structure, which is stabilized by three intramolecular hydrogen bonds, Gly7(N)Asp3(O) (Hbond I), Thr8(N)Asp3(O) (Hbond II), and Asp3(N)Thr8(O) (Hbond III), and a hydrophobic interaction between the aromatic rings of Tyr2 and Trp9.
2.2. Improvement of the SAAPFF Parameters
The SAAPFF parameters are comprised of a potential energy, atomic coordinates, atomic charges, and LennardJones potential parameters for all possible conformations of the single amino acid unit (HCOXaaNH_{2}). The conformations are defined by dihedral angles, , , , , and so forth, each with an interval of 15 degrees. In this study, the SAAPFF parameters previously developed [33, 34] were improved in the following points.
First, the parameters for mainchain, Pro, and Val in water, which were obtained in the previous version [34] by singlepoint calculation with the IEFPCM [65–67] solvent model using the structures optimized in vacuo, were replaced with those obtained by geometry optimization in water at the HF/631+G(d) level using the IEFPCM model in a Gaussian03 program (rev. E.01) [68]. The similar calculation level has been widely applied for small biomolecules [69, 70]. Relaxation of the potential surfaces in the implicit water model should improve accuracy of the SAAPFF.
Second, atomic charges of the SAAPFF parameters, which are variables as a function of the structure of an amino acid unit, were switched for all amino acids from Mulliken charges to electrostatic potential (ESP) charges [71] that were obtained at the recommended higher MP2/631G(d,p) level [72] using the IEFPCM solvent model. The ESP charges are defined to reproduce the electrostatic potential surface and are used in conventional force fields as reliable atomic charges. Therefore, the use of ESP charges would reproduce electrostatic interactions, that is, the term, more appropriately than the Mulliken charges.
Third, the parameters of potential energies for the mainchain unit were further improved as follows. The SAAPFF parameters for most amino acids, except for Gly, Ala, Pro, and Val, are divided into the mainchain and sidechain units according to the sidechain separation approximation method [34]. In this approximation, the mainchain unit was capped with a tertbutyl group, which was fixed eclipsed to the mainchain part during geometry optimization to obtain the SAAPFF parameters. However, this method brought a statistical error that the SAAP energies thus obtained are slightly larger than those obtained without the approximation, due probably to distortion of the bond angles derived from the fixed geometry. Indeed, such a statistical error was evident for the case of Val [34]. Therefore, the energies were replaced with those calculated for the relaxed structures, which were obtained by geometry optimization with the tertbutyl group relaxed in water, while for structural parameters (i.e., atomic coordinates and atomic charges) the geometry with the tertbutyl group fixed eclipsed was still employed in order to avoid abnormally short atomic contacts between the mainchain and sidechain units after the connection.
2.3. Simulation Conditions
Conformational properties of Metenkephalin and chignolin in water were studied by SAAPMC and AMBERMD [73] simulations. The SAAPMC simulation was carried out at 300 K in water by using the improved SAAPFF parameters, in which solvation effects were implemented implicitly. The conventional Metropolis method [74, 75] and the Mersenne Twister random number generator [76] were employed for the MC simulation. At each MC step, one dihedral angle, randomly chosen from all dihedral angles including mainchain and sidechain dihedral angles, was changed randomly with a maximum displacement angle of ±32 degree. The new structure was constructed by the interpolation of the SAAPFF parameters and subsequent connection of the single amino acid units, and the potential energy was calculated by (1). The potential energy thus obtained was then compared with the previous one. Ten trajectories with total MC steps of 100 million were obtained for Metenkephalin, while twenty trajectories with total MC steps of 2 billion were obtained for chignolin. As an initial structure of Metenkephalin for the SAAPMC simulation, the extended structure from PDB (PDB ID: 1PLW) [77] was employed. On the other hand, extended (defined by all mainchain and dihedral angles of −180 degree) and native folded (PDB ID: 1UAO) [54] structures were selected for chignolin as the initial structures. An Intel Xeon W3550 3.06 GHz processor with 12 GB memory was employed as a platform for calculation. The computing time of SAAPMC simulation for Metenkephalin was about 6 hours and for chignolin was about 11 days.
For comparison, ten trajectories of AMBERMD [73] simulation were also calculated for Metenkephalin and chignolin, respectively, in water at 300 K. The simulation time was 1 s with a time step of 2 fs, and a cutoff distance of nonbonded interactions was 12 Å. Solvation effects were implicitly considered by using the generalized Born model () [35–37] with a dielectric constant of 78.5, following the recent study by Götz et al. [78]. The AMBER10 program [73] with the amber99sb force field [1] was used. The computing time of the AMBERMD simulation to obtain one trajectory was about 165 hours for Metenkephalin, and that for chignolin was about 18 days, on the same calculation platform used for the SAAPMC simulation.
2.4. Data Analysis
The 10,000 structures of Metenkephalin that were extracted in every 10,000 step from the SAAPMC simulation trajectory were classified to twenty structural clusters by using a clustering algorithm called the means method [79] based on the RMSD for the mainchain C, O, and N atoms or all heavy atoms except for the hydrogen atoms. Similarly, ten structural clusters were obtained for chignolin from the 20,000 structures extracted in every 100,000 step from the SAAPMC simulation trajectories. On the other hand, 10,000 structures were extracted in every 100 ps from the AMBERMD trajectories for both Metenkephalin and chignolin. The obtained structures were classified to ten or twenty clusters by application of the means method under the same conditions to those applied for the SAAPMC simulation results.
For chignolin, the RMSD values with respect to the native structure were also calculated for the structures extracted from the SAAPMC and AMBERMD trajectories by using the amberptraj program [7, 8] based on the mainchain atoms or all heavy atoms. The freeenergy surfaces of chignolin were further analyzed by using all trajectories obtained from the SAAPMC simulation according to (2): where is a probability of the structures with the corresponding structural parameters, is the probability of the structure that had the lowest energy, is a difference in the free energy, is a gas constant, and is an absolute temperature. Ramachandrantype freeenergy surfaces were obtained for each amino acid residue by dividing the surface into 10 × 10 degree units, while mainchain RMSD versus hydrogen bond or hydrogen bond versus hydrogen bond freeenergy surfaces were obtained by dividing the surface into 0.1 × 0.1 Å units.
2.5. Ab Initio Calculation
Relative energies of the representative structures obtained for chignolin by the SAAPMC simulation were calculated in water by the ab initio molecular orbital method. The calculation was performed at the HF/IEFPCM/631+G(d,p) level by using a Gaussian03 program [68]. The polarizable continuum model using the integral equation formalism variant (IEFPCM), which is equipped as the default selfconsistent reaction field method, was employed for the calculation in water.
3. Results
3.1. MetEnkephalin by SAAPMC and AMBERMD
Energetic trajectories of SAAPMC and AMBERMD simulations for Metenkephalin were obtained in an implicit water model at 300 K (see Figure S1). In SAAPFF, the total energy is basically expressed as a sum of , , and . The values of these energetic terms were maintained stable during the SAAPMC simulation with large fluctuation for and , while the value of was almost zero, as observed previously [34]: the large dielectric constant of water would make the electrostatic interaction between the amino acid residues ignorable. Similarly, the total energy was maintained stable for the AMBERMD simulation. However, the distribution of the distance between the terminal Cα atoms () was significantly different for the two simulation methods. The value of dispersed in a range from 5 to 12 Å in the SAAPMC simulation, while those converged at around 6 Å in the AMBERMD simulation (Figure 3). The results suggested that the structures obtained from the SAAPMC simulation contain diverse conformations and are more extended than those obtained from the AMBERMD simulation on the average. Indeed, when the structures were classified into twenty clusters based on the RMSD values for the mainchain atoms, various types of structures were obtained from the SAAPMC trajectory. The representative examples are shown in Figure 4. The structures are different from each other with a diverse value, and there is no hydrogen bond in the three representative structures. In contrast, the structures of Metenkephalin obtained from the AMBERMD trajectory seemed to be trapped in local energy minimums during the simulation (see Figure S2).
(a)
(b)
(a)
(b)
(c)
(d)
3.2. Chignolin by SAAPMC and AMBERMD
Twenty trajectories were obtained for chignolin by SAAPMC simulation. The structures in each trajectory were classified to ten clusters based on the mainchain RMSD by using the means method. According to the mainchain folds, that is, combination of the dihedral angles for the amino acid residues, the obtained clusters were further classified to three representative structures, which are called here A, B, and C (Figure 5). Populations of the three structures are summarized in Table 1. Structure A has a nativelike mainchain fold stabilized by three hydrogen bonds, which are the same as those observed in the native structure (i.e., Hbonds I–III). However, the hydrophobic interaction between the aromatic side chains of Tyr2 and Trp9 is not present. Structures B and C have nonnative conformation with a common type II βturn around Thr6 and Gly7. A weak hydrogen bond exists between Glu5(O) and Thr8(N) with a distance about 4.0 Å. However, the conformations of the Cterminals are different from each other. Other minor structures have various local and global conformations, such as a type I βturn, a helical conformation, and other conformations.

In the trajectories from 1 to 12, the existence ratio of structure A was high (31.2 to 58.5%), whereas the ratio was much lower than that of the misfolded structure B in the remaining trajectories. The mean ratio for structure A averaged for all trajectories was 29.2%. On the other hand, when the clustering analysis was performed based on allatom RMSD, the nativelike structure (A′) was obtained from 0 to 19.5% (Table 1). The existence ratio of structure A′ averaged over all trajectories was 6.7%. Structure A′ maintains three native Hbonds as well as a hydrophobic interaction between Tyr2 and Trp9 (Figure 5). The allatom and mainchain RMSD values of structure A′ from the native structure were 2.6 and 2.0 Å, respectively, suggesting that structure A′ keeps the native fold. The trace of the mainchain RMSD with respect to the native structure is shown in Figure 6(a) for the case of trajectory 10, where structures A and B are involved in almost equal amounts. The RMSD values fluctuated between 1 and 6 Å, and repeats of the structural transitions were observed. Similarly, in the AMBERMD simulation, the mainchain RMSD value fluctuated between 1 and 6 Å (Figure 6(b)). The structural clustering analysis showed that the native structure with three native hydrogen bonds was involved about 30% (Figure S5).
(a)
(b)
The convergence of the SAAPMC simulation is not clear from Table 1. However, the similar relative populations were obtained for structures A (~30% including 6% of A′) and B (~60%) when the simulation was started from the native structure. Moreover, in our preliminary results of the SAAPMC simulation using ten trajectories, almost the same mean populations were obtained for the representative structures (data not shown). Based on these observations, we employed the twenty trajectories shown in Table 2 for further analysis.

The structure with the smallest allatom RMSD (1.3 Å) obtained from the twenty trajectories is superimposed on the native structure in Figure 7. The mainchain and Cα atom RMSD values of this structure with respect to the native structure were 0.9 and 0.5 Å, respectively, confirming that the structure is almost identical to the native structure.
The Ramachandrantype freeenergy surfaces obtained for Tyr2–Trp9 residues from all trajectories of the SAAPMC simulation are shown in Figure 8. The profiles for Tyr2, Pro4, Glu5, and Thr6 fit well to the native structures, while the freeenergy minimums for Asp3, Gly7, and Thr8 are slightly moved from the positions of the native structures. For Trp9, the stable areas corresponding to α and β regions tend to merge. This is roughly consistent with diversity in the Cterminal conformation of the native structure. The dihedral angles for structures A and A′ are located around the energy minimum points on the freeenergy surfaces. In contrast, structure B has different geometry from structures A and A′ for Glu5 and Thr6: the angles are ca. + degree. Structure C has the geometry with angles of Glu5, Thr6, Thr8, and Trp9 different from those of structures A and A′. During the SAAPMC simulation, mismatch of the mainchain dihedral angles for Glu5, Thr6, and Gly7 to the native structure occurred at high ratios.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
The freeenergy surface of chignolin projected on Hbonds I versus III plane and that projected on the mainchain RMSD versus Hbond II plane obtained from all trajectories are shown in Figure 9. In Figure 9(a), the most stable area was found at Å and Å. This energy minimum corresponds to the misfolded structure B. Another potential energy minimum was found at Å and Å, corresponding to structures A and A′. Similarly, two stable minimums were located in Figure 9(b); Å and Å, corresponding to structure B, and Å and Å, corresponding to structure A. The freeenergy surfaces showed that structure A maintains three native hydrogen bonds, but it has larger mainchain RMSD from the native structure than misfolded and the most stable structure B. This conflict is ascribed to the N and Cterminals of structure A, which are significantly apart due to the absence of a hydrophobic interaction between Tyr2 and Trp9 side chains. It should be noted that structure A′ is separated from structure A in Figure 9(b).
(a)
(b)
3.3. Evaluation of the Relative Energies of Structures A–C by Ab Initio Calculation
To evaluate relative stabilities of structures A, A′, B, and C more accurately, singlepoint ab initio calculation was carried out. Table 2 lists the SCF and relative energies calculated in the IEFPCM solvent model. Although structure B was most populated in the SAAPMC simulation, it is significantly unstable compared to the other structures, suggesting that structure B would not be populated in practical solutions. On the other hand, structure A′ with a native fold was found to be the most stable. According to the ab initio calculation results, it is likely that the structures having more hydrogen bonds and/or a hydrophobic interaction between the aromatic rings are more stable.
4. Discussion
4.1. Improvement of SAAPFF Parameters
In this study, we have modified the SAAPFF parameters in the following points. The mainchain, Pro, and Val parameters in water were replaced with those obtained by geometry optimization in water. The atomic charges were replaced from Mulliken to ESP charges. The energies of the mainchain unit were replaced with those obtained for the geometry with the tertbutyl cap being relaxed during the optimization. The effects of these modifications were found to be significant as the SAAPMC simulation using the new parameters reproduced the native structure of chignolin (structure A′), while it was not obtained when the simulation was performed by using the previous parameters [34].
The clustering of the structures obtained for chignolin by SAAPMC simulation using the improved parameters showed that structures A and B are dominant in the implicit water model as shown in Table 1. Structure A (29.2%) maintains three native hydrogen bonds, indicating the native mainchain fold, while nonnative structure B (43.2%) is stabilized by a type II βturn. In addition, the structural clustering based on allatom RMSD demonstrated that native structure A′, which maintains not only the three native hydrogen bonds but also a hydrophobic interaction between the aromatic side chains of Tyr2 and Trp9, is involved in 6.7% of the total structures. The freeenergy analysis clearly indicated that the mainchain dihedral angles of the each amino acid residue and the distances of Hbonds I–III for structures A and A′ are consistent with the native structure (see Figures 8 and 9), although their populations were less than that of the misfolded structure B (vide infra).
In the meantime, randomly fluctuating structure was obtained for Metenkephalin by using both the previous [34] and improved SAAPFF parameters (see Figures 3 and 4). This would be due to the two Gly residues, whose SAAP parameters have been modified only for the atomic charges in this study, being involved in the mid of the amino acid sequence.
There are a number of reports in the literature on the molecular simulation for chignolin by using conventional force fields to obtain the native structure from the unfolded state [55–64]. In these studies, various stateoftheart techniques of molecular simulation have been employed, including a multicanonical MD method [55, 56], a replica exchange MD method [57–60], a longterm conventional MD method by using a highspeed computer system [62], and originally developed MHOP and MSES methods [63, 64]. The values of allatom RMSD for the obtained native structure were in a range from 1.3 to 1.8 Å with respect to the experimentally determined structure, which is smaller than that of structure A obtained as a major nativelike structure by the SAAPMC simulation in this study. Moreover, the AMBERMD simulation carried out in this study showed that the native structure of chignolin is a little more frequently obtained by using the amber99sb force field than SAAPFF. Thus, it is obvious that further improvement of the SAAPFF parameters is required in future works.
4.2. Efficiency of SAAPMC Simulation for Conformational Sampling
Previous studies applying an explicit water molecule model [43] or a generalizedensemble method with RISM theory [39–42] have shown that Metenkephalin has extended conformations in water. This feature was reasonably reproduced in this study by application of SAAPFF as shown in Figures 3 and 4. Clustering of the structures obtained from the SAAPMC simulation indicated diverse conformations of Metenkephalin in water. On the other hand, the AMBERMD simulation for Metenkephalin carried out in this study produced rather uniform folded conformations with intramolecular hydrogen bonds (see Figure S2), probably due to insufficient computing time. Thus, the SAAPMC method would be useful for sampling possible conformations for short peptide molecules in an implicit water model, although the accuracy is not satisfactory in terms of the reproducibility of the relative energies.
The efficiency for conformational sampling arises probably from a less number of variables used in the SAAPFF than that in the conventional force fields. In SAAPFF, the structure of peptides is defined only by the dihedral angles of the each amino acid unit, not by the cartesian coordinates of each atom. Acceleration of conformational sampling by reducing the number of structural parameters is a common strategy of coarsegrained force fields, such as unitedatom force fields [18–26, 80, 81]. Efficiency of the SAAPMC simulation for conformational sampling can also be explained by the fact that the solvent effects are implicitly included in the parameters of SAAPFF. This enables the SAAPMC simulation in an implicit water model to be performed in the same computing time as that in vacuo.
4.3. Prediction of Native Structures
Although the accuracy of SAAPFF is not yet satisfied, the efficiency for conformational sampling would be advantageous for prediction of the stable structures of peptides in water. Therefore, we subsequently explored the application of the SAAPMC method to predict the native structure of chignolin.
Indeed, it was found by singlepoint ab initio calculation in water at the HF/IEFPCM/631+G(d,p) level that structure A′ is the most stable one in water among the structures shown in Figure 5 (see Table 2). The procedure should involve a large error in the energy because the structure is not optimized. However, the tendency of the relative energies would be reliable to some extent as previously demonstrated for HCOAla_{2}NH_{2} in vacuo [33]. Thus, a protocol of SAAPMC simulation followed by structure clustering and examination of the obtained structures by ab initio calculation would be a useful strategy toward the prediction of the native structure for short peptides. The relative energies obtained by ab initio calculation also suggested the possibility that structure B can be produced as an artifact of the force field because the relative energy was quite large (>20 kcal/mol).
As for prediction of the native structure, another simpler method would be possible based on the hardness (i.e., the number of interamino acid interactions or the depth of the potential hole on the freeenergy surfaces) of the structures. As seen in Figure 5, structure B of chignolin maintains only one hydrogen bond, while structures A and A′ hold three hydrogen bonds. In terms of the number of hydrogen bonds, it can be assumed that structures A and A′ are more nativelike than structure B, although the population of structure B was larger than those of structures A and A′. Indeed, the areas around the freeenergy minimums corresponding to structures A and A′ in Figures 9(a) and 9(b) are slightly narrower than those corresponding to structures B and C. Such robust structures on the conformational surface would be predicted to be more nativelike than the structures with more conformational flexibility. Thus, we suggest that a protocol of SAAPMC simulation followed by structural clustering and examination of the obtained structures by the number of hydrogen bonds (or the hardness) may be another simple strategy toward structure prediction for short peptide molecules.
5. Conclusions
The parameters of SAAPFF, which was previously developed to analyze the structures and folding of polypeptides, have been improved in several points in this study. The mainchain, Pro, and Val parameters in water were replaced with those obtained by geometry optimization in the IEFPCM solvent model. The atomic charges were replaced from Mulliken to ESP charges. (3) The energies of the mainchain unit were replaced with those obtained for the geometry with the tertbutyl cap being relaxed during the optimization. To investigate the accuracy of the improved SAAPFF parameters, the SAAPMC simulation was carried out for Metenkephalin and chignolin as model peptides. Analysis of the obtained structures revealed that the SAAPMC simulation reasonably reproduced randomly fluctuating conformation for Metenkephalin, while three folded structures, one of which corresponds to a nativelike structure with three native hydrogen bonds, were obtained for chignolin. The latter result supported that the accuracy of the SAAPFF parameters has been remarkably improved from the previous ones, although the native structure of chignolin could not be assigned to the most stable structure in the conformational space defined by SAAPFF.
In the meantime, efficiency of SAAPMC simulation for conformational sampling was demonstrated for Metenkephalin as the SAAPMC simulation afforded diverse structures. The feature, combined with the structural clustering analysis, was subsequently applied to the structure prediction of chignolin. Among the representative structures obtained by the clustering, structure A′ with a native fold was assigned to the most stable structure according to ab initio calculation. Thus, a protocol of SAAPMC simulation followed by structural clustering and examination of the obtained structures by ab initio calculation would be a useful strategy toward prediction of the native structure for short peptides. Alternatively, the structure with a maximum number of hydrogen bonds or with the hardest conformation on the potential surface may be considered to be the most nativelike among the clustered structures. Applications of this method to structure prediction of other peptide molecules as well as further improvement of the SAAPFF parameters are being conducted by our group.
Acknowledgment
This work was supported by GrantinAid for Scientific Research on Innovative Areas (no. 2120005) from the Ministry of Education, Culture, Sports, Science and Technology. The SAAP force field parameters are available for download at http://saap.sc.utokai.ac.jp/.
Supplementary Materials
The SAAP force field and the simulation algorithm, the results of SAAPMC and AMBERMD simulations for Metenkephalin and chignolin, the structural parameters of structures A, A’, B, and C as well as the smallest allatom RMSD and the lowest total potential energy structures obtained for chignolin, and the results of singlepoint ab initio calculation for the representative structures.
References
 V. Hornak, R. Abel, A. Okur, B. Strockbine, A. Roitberg, and C. Simmerling, “Comparison of multiple amber force fields and development of improved protein backbone parameters,” Proteins: Structure, Function and Genetics, vol. 65, no. 3, pp. 712–725, 2006. View at: Publisher Site  Google Scholar
 J. W. Ponder and D. A. Case, “Force fields for protein simulations,” Advances in Protein Chemistry, vol. 66, pp. 27–85, 2003. View at: Publisher Site  Google Scholar
 F. A. Momany, R. F. McGuire, A. W. Burgess, and H. A. Scheraga, “Energy parameters in polypepltides. VII. Geometric parameters, partial atomic charges, nonbonded interactions, hydrogen bond interactions, and intrinsic torsional potentials for the naturally occurring amino acids,” Journal of Physical Chemistry, vol. 79, no. 22, pp. 2361–2381, 1975. View at: Google Scholar
 G. Némethy, K. D. Gibson, K. A. Palmer et al., “Energy parameters in polypeptides. 10. Improved geometrical parameters and nonbonded interactions for use in the ECEPP/3 algorithm, with application to prolinecontaining peptides,” Journal of Physical Chemistry, vol. 96, no. 15, pp. 6472–6484, 1992. View at: Google Scholar
 S. J. Weiner, P. A. Kollman, D. A. Case et al., “A new force field for molecular mechanical simulation of nucleic acids and proteins,” Journal of the American Chemical Society, vol. 106, no. 3, pp. 765–784, 1984. View at: Google Scholar
 S. J. Weiner, P. A. Kollman, D. T. Nguyen, and D. A. Case, “An all atom force field for simulations of proteins and nucleic acids,” Journal of Computational Chemistry, vol. 7, no. 2, pp. 230–252, 1986. View at: Publisher Site  Google Scholar
 D. A. Pearlman, D. A. Case, J. W. Caldwell et al., “AMBER, a package of computer programs for applying molecular mechanics, normal mode analysis, molecular dynamics and free energy calculations to simulate the structural and energetic properties of molecules,” Computer Physics Communications, vol. 91, no. 1–3, pp. 1–41, 1995. View at: Google Scholar
 D. A. Case, T. E. Cheatham III, T. Darden et al., “The Amber biomolecular simulation programs,” Journal of Computational Chemistry, vol. 26, no. 16, pp. 1668–1688, 2005. View at: Publisher Site  Google Scholar
 B. R. Brooks, R. E. Bruccoleri, B. D. Olafson, D. J. States, S. Swaminathan, and M. Karplus, “CHARMM: a program for macromolecular energy, minimization, and dynamics calculations,” Journal of Computational Chemistry, vol. 4, no. 2, pp. 187–217, 1983. View at: Publisher Site  Google Scholar
 W. L. Jorgensen and J. TiradoRives, “The OPLS potential functions for proteins. Energy minimizations for crystals of cyclic peptides and crambin,” Journal of the American Chemical Society, vol. 110, no. 6, pp. 1657–1666, 1988. View at: Google Scholar
 W. L. Jorgensen and J. TiradoRives, “Potential energy functions for atomiclevel simulations of water and organic and biomolecular systems,” Proceedings of the National Academy of Sciences of the United States of America, vol. 102, no. 19, pp. 6665–6670, 2005. View at: Publisher Site  Google Scholar
 W. R. P. Scott, P. H. Hünenberger, I. G. Tironi et al., “The GROMOS biomolecular simulation program package,” Journal of Physical Chemistry A, vol. 103, no. 19, pp. 3596–3607, 1999. View at: Google Scholar
 U. H. E. Hansmann, Y. Okamoto, and F. Eisenmenger, “Molecular dynamics, langevin and hybrid Monte Carlo simulations in a multicanonical ensemble,” Chemical Physics Letters, vol. 259, no. 34, pp. 321–330, 1996. View at: Google Scholar
 N. Nakajima, H. Nakamura, and A. Kidera, “Multicanonical ensemble generated by molecular dynamics simulation for enhanced conformational sampling of peptides,” Journal of Physical Chemistry B, vol. 101, no. 5, pp. 817–824, 1997. View at: Google Scholar
 K. Hukushima and K. Nemoto, “Exchange monte carlo method and application to spin glass simulations,” Journal of the Physical Society of Japan, vol. 65, no. 6, pp. 1604–1608, 1996. View at: Google Scholar
 Y. Sugita and Y. Okamoto, “Replicaexchange molecular dynamics method for protein folding,” Chemical Physics Letters, vol. 314, no. 12, pp. 141–151, 1999. View at: Publisher Site  Google Scholar
 Y. Okamoto, “Generalizedensemble algorithms: enhanced sampling techniques for Monte Carlo and molecular dynamics simulations,” Journal of Molecular Graphics and Modelling, vol. 22, no. 5, pp. 425–439, 2004. View at: Publisher Site  Google Scholar
 L. Yang, C. Tan, M.J. Hsieh et al., “Newgeneration amber unitedatom force field,” The Journal of Physical Chemistry B, vol. 110, no. 26, pp. 13166–13176, 2006. View at: Publisher Site  Google Scholar
 A. Liwo, S. Ołdziej, M. R. Pincus, R. J. Wawak, S. Rackovsky, and H. A. Scheraga, “A unitedresidue force field for offlattice proteinstructure simulations. I. Functional forms and parameters of longrange sidechain interaction potentials from protein crystal data,” Journal of Computational Chemistry, vol. 18, no. 7, pp. 849–873, 1997. View at: Google Scholar
 A. Liwo, M. R. Pincus, R. J. Wawak, S. Rackovsky, S. Ołdziej, and H. A. Scheraga, “A unitedresidue force field for offlattice proteinstructure simulations. II. Parameterization of shortrange interactions and determination of weights of energy terms by Zscore optimization,” Journal of Computational Chemistry, vol. 18, no. 7, pp. 874–887, 1997. View at: Google Scholar
 A. Liwo, R. Kaźmierkiewicz, C. Czaplewski et al., “Unitedresidue force field for offlattice proteinstructure simulations: III. Origin of backbone hydrogenbonding cooperativity in unitedresidue potentials,” Journal of Computational Chemistry, vol. 19, no. 3, pp. 259–276, 1998. View at: Google Scholar
 A. Liwo, C. Czaplewski, J. Pillardy, and H. A. Scheraga, “Cumulantbased expressions for the multibody terms for the correlation between local and electrostatic interactions in the unitedresidue force field,” Journal of Chemical Physics, vol. 115, no. 5, pp. 2323–2347, 2001. View at: Publisher Site  Google Scholar
 A. V. Rojas, A. Liwo, and H. A. Scheraga, “Molecular dynamics with the unitedresidue force field: Ab initio folding simulations of multichain proteins,” Journal of Physical Chemistry B, vol. 111, no. 1, pp. 293–309, 2007. View at: Publisher Site  Google Scholar
 S. J. Marrink, H. J. Risselada, S. Yefimov, D. P. Tieleman, and A. H. De Vries, “The MARTINI force field: coarse grained model for biomolecular simulations,” Journal of Physical Chemistry B, vol. 111, no. 27, pp. 7812–7824, 2007. View at: Publisher Site  Google Scholar
 L. Monticelli, S. K. Kandasamy, X. Periole, R. G. Larson, D. P. Tieleman, and S. Marrink, “The MARTINI coarsegrained force field: extension to proteins,” Journal of Chemical Theory and Computation, vol. 4, no. 5, pp. 819–834, 2008. View at: Publisher Site  Google Scholar
 A. Rojas, A. Liwo, D. Browne, and H. A. Scheraga, “Mechanism of fiber assembly: treatment of Aβ peptide aggregation with a coarsegrained unitedresidue force field,” Journal of Molecular Biology, vol. 404, no. 3, pp. 537–552, 2010. View at: Publisher Site  Google Scholar
 M. Iwaoka, D. Yosida, and N. Kimura, “Importance of the single amino acid potential in water for secondary and tertiary structures of proteins,” Journal of Physical Chemistry B, vol. 110, no. 29, pp. 14475–14482, 2006. View at: Publisher Site  Google Scholar
 M. Iwaoka, M. Okada, and S. Tomoda, “Solvent effects on the φψ potential surfaces of glycine and alanine dipeptides studied by PCM and IPCM methods,” Journal of Molecular Structure: THEOCHEM, vol. 586, pp. 111–124, 2002. View at: Publisher Site  Google Scholar
 L. J. Smith, K. M. Fiebig, H. Schwalbe, and C. M. Dobson, “The concept of a random coil : residual structure in peptides and denatured proteins,” Folding and Design, vol. 1, no. 5, pp. R95–R106, 1996. View at: Google Scholar
 K. M. Fiebig, H. Schwalbe, M. Buck, L. J. Smith, and C. M. Dobson, “Toward a description of the conformations of denatured states of proteins. Comparison of a random coil model with NMR measurements,” Journal of Physical Chemistry, vol. 100, no. 7, pp. 2661–2666, 1996. View at: Google Scholar
 Y. Sakae and Y. Okamoto, “Optimization of protein forcefield parameters with the Protein Data Bank,” Chemical Physics Letters, vol. 382, no. 56, pp. 626–636, 2003. View at: Publisher Site  Google Scholar
 N. Kamiya, Y. S. Watanabe, S. Ono, and J. Higo, “AMBERbased hybrid force field for conformational sampling of polypeptides,” Chemical Physics Letters, vol. 401, no. 1–3, pp. 312–317, 2005. View at: Publisher Site  Google Scholar
 M. Iwaoka and S. Tomoda, “The SAAP force field. A simple approach to a new allatom protein force field by using Single Amino Acid Potential (SAAP) functions in various solvents,” Journal of Computational Chemistry, vol. 24, no. 10, pp. 1192–1200, 2003. View at: Publisher Site  Google Scholar
 M. Iwaoka, N. Kimura, D. Yosida, and T. Minezaki, “The SAAP force field: development of the single amino acid potentials for 20 proteinogenic amino acids and monte carlo molecular simulation for short peptides,” Journal of Computational Chemistry, vol. 30, no. 13, pp. 2039–2055, 2009. View at: Publisher Site  Google Scholar
 G. D. Hawkins, C. J. Cramer, and D. G. Truhlar, “Parametrized models of aqueous free energies of solvation based on pairwise descreening of solute atomic charges from a dielectric medium,” Journal of Physical Chemistry, vol. 100, no. 51, pp. 19824–19839, 1996. View at: Google Scholar
 G. D. Hawkins, C. J. Cramer, and D. G. Truhlar, “Pairwise solute descreening of solute charges from a dielectric medium,” Chemical Physics Letters, vol. 246, no. 12, pp. 122–129, 1995. View at: Google Scholar
 V. Tsui and D. A. Case, “Theory and applications of the generalized Born solvation model in macromolecular simulations,” Biopolymers, vol. 56, no. 4, pp. 275–291, 2001. View at: Google Scholar
 W. H. Graham, E. S. Carter II., and R. P. Hicks, “Conformational analysis of Metenkephalin in both aqueous solution and in the presence of sodium dodecyl sulfate micelles using multidimensional NMR and molecular modeling,” Biopolymers, vol. 32, no. 12, pp. 1755–1764, 1992. View at: Publisher Site  Google Scholar
 A. Mitsutake, M. Kinoshita, Y. Okamoto, and F. Hirata, “Combination of the replicaexchange Monte Carlo method and the reference interaction site model theory for simulating a peptide molecule in aqueous solution,” Journal of Physical Chemistry B, vol. 108, no. 49, pp. 19002–19012, 2004. View at: Publisher Site  Google Scholar
 A. Mitsutake, M. Kinoshita, Y. Okamoto, and F. Hirata, “Multicanonical algorithm combined with the RISM theory for simulating peptides in aqueous solution,” Chemical Physics Letters, vol. 329, no. 34, pp. 295–303, 2000. View at: Publisher Site  Google Scholar
 M. Kinoshita, Y. Okamoto, and F. Hirata, “Firstprinciple determination of peptide conformations in solvents: combination of Monte Carlo simulated annealing and RISM theory,” Journal of the American Chemical Society, vol. 120, no. 8, pp. 1855–1863, 1998. View at: Publisher Site  Google Scholar
 M. Kinoshita, Y. Okamoto, and F. Hirata, “Solvation structure and stability of peptides in aqueous solutions analyzed by the reference interaction site model theory,” Journal of Chemical Physics, vol. 107, no. 5, pp. 1586–1599, 1997. View at: Google Scholar
 L. Ramya and N. Gautham, “Effects of hydration on the conformational energy landscape of the pentapeptide Metenkephalin,” Journal of Chemical Theory and Computation, vol. 5, no. 8, pp. 2180–2190, 2009. View at: Publisher Site  Google Scholar
 K. Y. Sanbonmatsu and A. E. García, “Structure of Metenkephalin in explicit aqueous solution using replica exchange molecular dynamics,” Proteins: Structure, Function and Genetics, vol. 46, no. 2, pp. 225–234, 2002. View at: Publisher Site  Google Scholar
 D. Van Der Spoel and H. J. C. Berendsen, “Molecular dynamics simulations of leuenkephalin in water and DMSO,” Biophysical Journal, vol. 72, no. 5, pp. 2032–2041, 1997. View at: Google Scholar
 L. Zhan, J. Z. Y. Chen, and W. Liu, “Comparison of predicted native structures of metenkephalin based on various accessiblesurfacearea solvent models,” Journal of Computational Chemistry, vol. 30, no. 7, pp. 1051–1058, 2009. View at: Publisher Site  Google Scholar
 L. Su and R. I. Cukier, “Hamiltonian and distance replica exchange method studies of metenkephalin,” Journal of Physical Chemistry B, vol. 111, no. 42, pp. 12310–12321, 2007. View at: Publisher Site  Google Scholar
 D. A. Evans and D. J. Wales, “The free energy landscape and dynamics of metenkephalin,” Journal of Chemical Physics, vol. 119, no. 18, pp. 9947–9955, 2003. View at: Publisher Site  Google Scholar
 M. M. Palian, Boguslavsky, D. F. O'Brien, and R. Polt, “Glycopeptidemembrane interactions: glycosyl enkephalin analogues adopt turn conformations by NMR and CD in amphipathic media,” Journal of the American Chemical Society, vol. 125, no. 19, pp. 5823–5831, 2003. View at: Publisher Site  Google Scholar
 B. A. Berg and H. Hsu, “Metropolis simulations of MetEnkephalin with solventaccessible area parametrizations,” Physical Review E, vol. 69, no. 2, Article ID 026703, 9 pages, 2004. View at: Publisher Site  Google Scholar
 K. Vengadesan and N. Gautham, “Energy landscape of Metenkephalin and Leuenkephalin drawn using mutually orthogonal Latin squares sampling,” Journal of Physical Chemistry B, vol. 108, no. 30, pp. 11196–11205, 2004. View at: Publisher Site  Google Scholar
 S. G. Itoh and Y. Okamoto, “Effective sampling in the configurational space of a small peptide by the multicanonicalmultioverlap algorithm,” Physical Review E, vol. 76, no. 2, Article ID 026705, 2007. View at: Publisher Site  Google Scholar
 A. Mitsutake, U. H. E. Hansmann, and Y. Okamoto, “Temperature dependence of distributions of conformations of a small peptide,” Journal of Molecular Graphics and Modelling, vol. 16, no. 4–6, pp. 226–238, 1998. View at: Publisher Site  Google Scholar
 S. Honda, K. Yamasaki, Y. Sawada, and H. Morii, “10 Residue folded peptide designed by segment statistics,” Structure, vol. 12, no. 8, pp. 1507–1518, 2004. View at: Publisher Site  Google Scholar
 D. Satoh, K. Shimizu, S. Nakamura, and T. Terada, “Folding freeenergy landscape of a 10residue miniprotein, chignolin,” FEBS Letters, vol. 580, no. 14, pp. 3422–3426, 2006. View at: Publisher Site  Google Scholar
 T. Terada, D. Satoh, T. Mikawa, Y. Ito, and K. Shimizu, “Understanding the roles of amino acid residues in tertiary structure formation of chignolin by using molecular dynamics simulation,” Proteins: Structure, Function and Genetics, vol. 73, no. 3, pp. 621–631, 2008. View at: Publisher Site  Google Scholar
 D. Van Der Spoel and M. M. Seibert, “Protein folding kinetics and thermodynamics from atomistic simulations,” Physical Review Letters, vol. 96, no. 23, Article ID 238102, 2006. View at: Publisher Site  Google Scholar
 M. M. Seibert, A. Patriksson, B. Hess, and D. Van Der Spoel, “Reproducible polypeptide folding and structure prediction using molecular dynamics simulations,” Journal of Molecular Biology, vol. 354, no. 1, pp. 173–183, 2005. View at: Publisher Site  Google Scholar
 W. Xu, T. Lai, Y. Yang, and Y. Mu, “Reversible folding simulation by hybrid Hamiltonian replica exchange,” Journal of Chemical Physics, vol. 128, no. 17, Article ID 175105, 2008. View at: Publisher Site  Google Scholar
 S. Kannan and M. Zacharias, “Enhanced sampling of peptide and protein conformations using replica exchange simulations with a peptide backbone biasingpotential,” Proteins: Structure, Function and Genetics, vol. 66, no. 3, pp. 697–706, 2007. View at: Publisher Site  Google Scholar
 X. Dou and J. Wang, “Folding free energy landscape of the decapeptide chignolin,” Modern Physics Letters B, vol. 22, no. 31, pp. 3087–3098, 2008. View at: Publisher Site  Google Scholar
 A. Suenaga, T. Narumi, N. Futatsugi et al., “Folding dynamics of 10residue betahairpin peptide chignolin,” Chemistry, vol. 2, no. 5, pp. 591–598, 2007. View at: Publisher Site  Google Scholar
 S. Roy, S. Goedecker, M. J. Field, and E. Penev, “A minima hopping study of allatom protein folding and structure prediction,” Journal of Physical Chemistry B, vol. 113, no. 20, pp. 7315–7321, 2009. View at: Publisher Site  Google Scholar
 K. Moritsugu, T. Terada, and A. Kidera, “Scalable free energy calculation of proteins via multiscale essential sampling,” Journal of Chemical Physics, vol. 133, no. 22, Article ID 224105, 2010. View at: Publisher Site  Google Scholar
 B. Mennucci, E. Cancès, and J. Tomasi, “Evaluation of solvent effects in isotropic and anisotropic dielectrics and in ionic solutions with a unified integral equation method: theoretical bases, computational implementation, and numerical applications,” Journal of Physical Chemistry B, vol. 101, no. 49, pp. 10506–10517, 1997. View at: Google Scholar
 E. Cancès, B. Mennucci, and J. Tomasi, “A new integral equation formalism for the polarizable continuum model: theoretical background and applications to Isotropic and anisotropic dielectrics,” Journal of Chemical Physics, vol. 107, no. 8, pp. 3032–3041, 1997. View at: Google Scholar
 G. Scalmani and M. J. Frisch, “Continuous surface charge polarizable continuum models of solvation. I. General formalism,” Journal of Chemical Physics, vol. 132, no. 11, Article ID 114110, 2010. View at: Publisher Site  Google Scholar
 M. J. Frisch, G. W. Trucks, H. B. Schlegel et al., Gaussian 03, Revision C.02, Gaussian, Wallingford, UK, 2004.
 Y. K. Kang, “Ab initio and DFT conformational study of proline dipeptide,” Journal of Molecular Structure: THEOCHEM, vol. 675, no. 1–3, pp. 37–45, 2004. View at: Publisher Site  Google Scholar
 M. Makowski, J. Makowska, and L. Chmurzyński, “Ab initio study of the energetics of protonation, deprotonation and homocomplexed cations and anions formation in systems modeling side chains of biomolecules,” Journal of Molecular Structure: THEOCHEM, vol. 674, no. 1–3, pp. 61–67, 2004. View at: Publisher Site  Google Scholar
 B. H. Besler, K. M. Merz Jr., and P. A. Kollman, “Atomic charges derived from semiempirical methods,” Journal of Computational Chemistry, vol. 11, no. 4, pp. 431–439, 1990. View at: Publisher Site  Google Scholar
 V. V. Kostjukov, N. M. Khomytova, A. A. H. Santiago, R. L. Ibarra, D. B. Davies, and M. P. Evstigneev, “Calculation of the electrostatic charges and energies for intercalation of aromatic drug molecules with DNA,” International Journal of Quantum Chemistry, vol. 111, no. 3, pp. 711–721, 2011. View at: Publisher Site  Google Scholar
 D. A. Case, T. A. Darden, T. E. Cheatham III et al., AMBER 10, University of California, San Francisco, Calif, USA, 2008.
 N. Metropolis, A. W. Rosenbluth, M. N. Rosenbluth, A. H. Teller, and E. Teller, “Equation of state calculations by fast computing machines,” The Journal of Chemical Physics, vol. 21, no. 6, pp. 1087–1092, 1953. View at: Google Scholar
 K. Binder, “Applications of Monte Carlo methods to statistical physics,” Reports on Progress in Physics, vol. 60, no. 5, pp. 487–559, 1997. View at: Publisher Site  Google Scholar
 M. Matsumoto and T. Nishimura, “Mersenne Twister: a 623dimensionally equidistributed uniform pseudorandom number generator,” ACM Transactions on Modeling and Computer Simulation, vol. 8, no. 1, pp. 3–30, 1998. View at: Google Scholar
 I. Marcotte, F. Separovic, M. Auger, and S. M. Gagné, “A multidimensional ^{1}H NMR investigation of the conformation of methionineenkephalin in fasttumbling bicelles,” Biophysical Journal, vol. 86, no. 3, pp. 1587–1600, 2004. View at: Google Scholar
 A. W. Götz, M. J. Williamson, D. Xu, D. Poole, S. La Grand, and R. C. Walker, “Routine microsecond molecular dynamics simulations with AMBER on GPUs. 1. Generalized Born,” Journal of Chemical Theory and Computation, vol. 8, no. 5, pp. 1542–1555, 2012. View at: Publisher Site  Google Scholar
 J. Shao, S. W. Tanner, N. Thompson, and T. E. Cheatham III., “Clustering molecular dynamics trajectories: I. Characterizing the performance of different clustering algorithms,” Journal of Chemical Theory and Computation, vol. 3, no. 6, pp. 2312–2334, 2007. View at: Publisher Site  Google Scholar
 H. A. Scheraga, M. Khalili, and A. Liwo, “Proteinfolding dynamics: overview of molecular simulation techniques,” Annual Review of Physical Chemistry, vol. 58, pp. 57–83, 2007. View at: Publisher Site  Google Scholar
 Y. Sakae and Y. Okamoto, “Optimisation of OPLSUA forcefield parameters for protein systems using protein data bank,” Molecular Simulation, vol. 36, no. 14, pp. 1148–1156, 2010. View at: Publisher Site  Google Scholar
Copyright
Copyright © 2013 Kenichi Dedachi et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.