Molecular dynamics study on the effects of charged amino acid distribution under low pH condition to the unfolding of hen egg white lysozyme and formation of beta strands.

Aggregation of unfolded or misfolded proteins into amyloid fibrils can cause various diseases in humans. However, the fibrils synthesized in vitro can be developed toward useful biomaterials under some physicochemical conditions. In this study, atomistic molecular dynamics simulations were performed to address the mechanism of beta-sheet formation of the unfolded hen egg-white lysozyme (HEWL) under a high temperature and low pH. Simulations of the protonated HEWL at pH 2 and the non-protonated HEWL at pH 7 were performed at the highly elevated temperature of 450 K to accelerate the unfolding, followed by the 333 K temperature to emulate some previous in vitro studies. The simulations showed that HEWL unfolded faster, and higher beta-strand contents were observed at pH 2. In addition, one of the simulation replicas at pH 2 showed that the beta-strand forming sequence was consistent with the 'K-peptide', proposed as the core region for amyloidosis in previous experimental studies. Beta-strand formation mechanisms at the earlier stage of amyloidosis were explained in terms of the radial distribution of the amino acids. The separation between groups of positively charged sidechains from the hydrophobic core corresponded to the clustering of the hydrophobic residues and beta-strand formation.

1. Introduction

Lysozymes are a family of globular enzymes in the immune systems of animals. A lysozyme molecule is a single polypeptide of around 130 residues that can partially hydrolyze the peptidoglycans of gram-positive bacterial cell walls [1, 2]. Lysozyme is one of the protein types associated with the formation of amyloid fibrils, as the failure of specific peptides or proteins to fold or to remain correctly folded triggered non-functional protein aggregation [3]. Amyloid aggregation is a hallmark of several degenerative diseases, such as Alzheimer’s disease, Parkinson’s disease, type II diabetes, Creutzfeldt Jakob, Huntington’s, amyotrophic lateral sclerosis (ALS) [4, 5]. However, non-toxic forms of amyloid fibrils could be useful in some areas of applications, such as bioengineering, biosensor, drug delivery, regenerative medicine, cell-encapsulating materials, tissue engineering, molecular and electronic devices, etc. [6-17]. Amyloid fibril formation in vitro can be controlled under several factors that either stimulate or inhibit aggregation [18]. Hen egg-white lysozyme (HEWL) is one the most commonly used proteins in protein aggregation research due to its well-characterized structure and low cost. A high temperature and an acidic condition can accelerate the unfolding and aggregation of HEWL, as a previous study showed that amyloid fibrils were formed when incubating HEWL under the temperature of 65°C and pH 2 for 196 hours [1].

Conformational changes of the unfolded and aggregated proteins can be monitored by using several methods, such as atomic force microscopy (AFM), Raman spectroscopy, and electrochemical impedance spectroscopy (EIS) [19-21]. Under an elevated temperature and acidic condition, AFM could provide images of lysozyme aggregation at various incubation times by depicting the time-dependent sizes of self‐assembled spheroidal oligomers. Changes in both secondary and tertiary structures of aqueous lysozymes could also be observed by shifting Raman spectra in the range of 650–1875 cm-1 from its native state [1]. Signals from EIS spectroscopy are highly sensitive to structural changes of proteins due to their unique charge transfer resistances (Rct) [21]. These methods, along with dynamic light scattering, size exclusion chromatography, transmission electron microscopy (TEM), surface plasmon resonance (SPR) [1, 22, 23], have their advantages but still lack understanding of the molecular details of protein structural changes.

Computer simulations have become an alternative tool to provide more insight into the structures and dynamics of proteins at different states under different conditions in atomistic details. The accuracy of atomistic molecular dynamics (MD) simulations to predict the molecular behavior of proteins under extreme conditions [24, 25] has been improved with the continuing development of molecular mechanics forcefield parameters [26, 27]. Therefore, some computational studies were conducted to observe the effects of temperature, solvents, and external perturbation on the unfolding and aggregation of human lysozymes and HEWL [25, 28–34]. From MD simulations performed by Moraitakis et al., replacing an aspartic acid residue with a histidine caused human lysozyme to unfold significantly faster under high temperature and could reproduce an experimental result. The mutation was proposed as a possible seed for amyloidosis [25]. Jafari et al. demonstrated in molecular detail that the lysozyme unfolds better in high concentrations of the sodium dodecyl sulfate (SDS) surfactant at 370 K, higher than the thermal denaturation midpoint temperature (Tm) [28]. Jiang et al. reported that high electric fields could enhance the possibility of protein unfolding due to the heterogeneous nature of charge distribution within proteins [29].

In our study, the effects of pH on the propensity of beta-strand formation that might lead to amyloidosis were addressed in atomistic details as the low pH condition was previously reported to facilitate amyloidosis. In addition, the altered electrostatic properties by the protonation of some amino acids should affect protein denaturation and beta-strand formation. Firstly, a series of MD simulations were performed under a high temperature at both pH conditions to accelerate unfolding processes, followed by simulations under an optimum temperature for amyloidosis reported in previous studies. Then, the conformational analysis was performed to characterize the simulated protein structures and provide the detailed mechanisms of beta-sheet formation under low pH and molecular insight of accelerating amyloidosis in vitro.

2. Methodology

GROMACS 2019.6 simulation package was used to carry out the Molecular Dynamics simulations. The initial atomistic structure of hen-egg white lysozyme was obtained from Protein Data Bank (PDB ID: 1AKI). In order to emulate pH conditions, protein coordinates were input to the PropKa software [35, 36] to estimate pKa values for aspartic acid, glutamic acid, and histidine residues. If the specified pH was lower than the pKa of a residue, the protonation state was assigned to that residue. At the pH 2 condition, all histidines, glutamic acids, and aspartic acids were fully protonated. Meanwhile, at pH 7, all the aforementioned amino acids were deprotonated. All protein structures were parameterized by the GROMOS54a7 forcefield [37-39], suitable for the prediction of protein unfolding under high temperatures [37]. Then, the lysozyme structures at both pH conditions were solvated by using the SPC water model [38, 39] within simulation boxes of size 6.9 x 6.9 x 6.9 nm3, which was large enough to cover a whole protein molecule with 1 nm buffer distance. As HEWL in both protonated and deprotonated forms were positively charged, all systems were neutralized by adding Cl- counter-ions. For each simulation, the energy minimization was performed with the steepest descent algorithm for the maximum number of steps of 50000. The short-range electrostatic cutoff distance and the short-range Van der Waals cutoff distance were 1.0 nm. PME treatment was used to calculate long-range electrostatic interactions. Then, an NPT equilibration stage was performed at T = 300 K and P = 1 bar for 200 ns, which was long enough to accommodate the conformation changes from protonation states introduced at pH 2. After that, three replicas of MD simulations at 450 K were performed for 100 ns to accelerate unfolding processes of both protonated (pH 2) and deprotonated (pH 7) structures. Then, all six simulations were continued at 333 K for 200 ns to emulate the in vitro condition that beta-strand formation was observed.

All MD trajectory replicas at both pH conditions were then analysed. Root mean square deviation (RMSD) was calculated to quantify the level of global conformational changes compared with the reference PDB structure. The time-dependent information of the secondary structure content of lysozymes at both pH 2 and pH 7 were analyzed by the DSSP algorithm [40], which identified the type of secondary structure for all regions within the protein. During the unfolding process at 450 K and the refolding processes at 333 K, the states of backbone torsions and alpha-beta transitions were monitored by Ramachandran plots for all simulation replicas at both pH conditions. The radius of gyration (Rg) was calculated as a function of time for different amino acid groups of the HEWL to monitor the distribution of charged amino acids and hydrophobic amino acids.

3. Results and discussions

A series of atomistic molecular dynamics (MD) simulations on HEWL under two different pH conditions and conformational analysis of all simulation replicas were performed to understand the effect of low pH on the unfolding of HEWL and the propensity of beta-strand formation. For each replica, HEWL was explicitly simulated at 450 K for 100 ns to accelerate the unfolding process within a feasible timescale for MD simulations, followed by a 200-ns run at 333 K temperature to emulate the condition where amyloidosis was observed in vitro. In Fig 1, the root mean square deviation (RMSD) of HEWL was calculated from all MD trajectories to monitor the global conformational changes of HEWL, along with the snapshots at the start, the middle (50 ns), and the end (100 ns) of the 450-K unfolding simulation. For the pH 2 condition at 450 K (Fig 1A), the RMSD values of lysozyme in all repeats increased rapidly within the first 30 ns. The replica R0 displayed the most abrupt change as the RMSD exceeded 1.2 nm after 10 ns, while the RMSD of replicas R1 and R2 increased beyond 1.2 nm after 30 ns and 20 ns, respectively. The conformational snapshots of HEWL replicas showed that most of the alpha-helical structures disappeared in R0 and R2 after 50 ns. Only one short alpha-helix was found for the R1 replica at 450 K and pH 2 after 100 ns. Interestingly, additional beta-strands were formed in all three replicas, especially in R0 and R2. The averaged RMSD from all replicas at 450 K and pH 2 after 50 ns was 1.35 ± 0.15 nm. Less fluctuation was then observed for the RMSD values when the systems at pH 2 were cooled down to 333 K, and the averaged RMSD was found at 1.30 ± 0.10 nm. For simulations at pH 7 and 450 K (Fig 1B), RMSD of the replicas R0, R1, and R2 exceeded 1.2 nm after 35 ns, 20 ns, and 30 ns, respectively. Similar to the pH 2 simulations, most of the alpha-helices of HEWL became lost during 450 K simulations. However, more helical structures were still left at pH 7 as two short helices were observed for the R1 replica, while one helical structure was observed for the R0 and R2 replicas. No beta-strand was observed after 100 ns for all pH 7 simulations. Averaged RMSD at pH 7 and 450 K after 50 ns was found around 1.30 ± 0.15 nm and became less fluctuating when the temperature was decreased to 333 K. Transition of alpha helices to either beta strands or random coils could also be visualized through the Ramachandran plots (see S1 and S2 Figs in S1 File).

Fig 1

RMSD calculated along all three replicas of MD simulations at 450 K and 333 K–R0 (black), R1 (red), and R2 (green)—at (a) pH 2 and (b) pH 7. Vertical dashed lines represent the time t = 100 ns where the temperature was switched from 450 K to 333 K. Conformational snapshots after 0 ns, 50 ns, and 100 ns of 450 K simulations were also presented for each replica, in which alpha helices were represented in purple and beta strands were represented in yellow.

RMSD results from 450 K simulations at pH 2 and 7 suggested that HEWL unfolded faster at pH 2, incorporating additional beta-strand formation. The DSSP (Define Secondary Structure of Protein) algorithm could quantitatively measure the secondary structure content for all simulations by identifying the types of structures for all the amino acid residues from patterns of hydrogen bonding network. Time evolution of the secondary structure content within the simulated HEWL at pH 2 (Fig 2) showed that most of the alpha-helical structures were lost during the 100-ns simulations at 450 K. In the case of the R0 replica at pH 2, three of the four alpha-helices were denatured after only 10 ns, corresponding to the most rapidly increasing RMSD. The R2 replica displayed similar behavior at the early 450 K stage but with a slightly longer timescale than the R0 replica. For the R1 replica, an alpha-helix of the sequence TASVNCAKKIVS at residues 89–100 could withstand the 450 K temperature up to 90 ns. The same helical region was also the last alpha-helix to denature for the R0 and R2 replicas. Beta strand formation and deformation occurred at different times and locations under 450 K for all replicas at pH 2. After the temperature was switched from 450 K to 333 K at 100 ns, the regions identified as beta-strands or at the beginning of 333 K simulations were extended, corresponding to the beta-strands observed in the final snapshots of all three replicas for pH 2 simulations. Only short beta-strands were found for the replica R1, while a long beta-strand of sequence VCAAKFE at residues 29–35 was found for the replica R2. Meanwhile, the more extended SLGNWVCAAKFES beta-strand at residues 24–36 was found for the replica R0 about the same region, along with the YGILQINSRWW beta-strand at residues 53–63. The time evolution of the secondary structures of simulated HEWL at pH 7 (Fig 3) showed significantly less amount of beta-strand formed after the temperature switching at 100 ns. On the other hand, alpha-helices were formed at the 333 K simulation stages of the R1 and R2 replicas around the region containing residues 89–100, where the alpha helix of sequence TASVNCAKKIVS displayed high thermostability at both pH conditions. Another alpha-helix near the N-terminus was formed during the 200-ns simulation at 333 K and pH 7 for the replica R1.

Fig 2

Time evolution of the secondary structures contents determined by the DSSP algorithm along the 450-K and 333-K simulation replicas (a) R0, (b) R1, and (c) R2 at pH 2. Snapshots of the final structure after the 333 K simulation was also given for each replica.

Fig 3

Time evolution of the secondary structures contents determined by the DSSP algorithm along the 450-K and 333-K simulation replicas (a) R0, (b) R1, and (c) R2 at pH 7. Snapshots of the final structure after the 333 K simulation was also given for each replica.

Tables 1 and 2 summarized the amount of alpha-helix and beta-strand secondary structures for all simulation replicas at both pH conditions. The DSSP analysis on the HEWL structures after a 200-ns equilibration at 300 K showed that the alpha-helix content within the native HEWL was 32.5% of the whole structure at pH 2, and 34.1% at pH 7. Meanwhile, the beta-sheet or beta-strand content within the HEWL was 6.2% of the whole structure at both pH conditions. After all the 300-ns simulations finished, the final percentage of alpha-helix content was found between 0.0% - 4.4% at pH2 and 0.1% - 13.0% at pH 7, suggesting that the full protonation on HEWL at pH 2 corresponded to more alpha helix loss. However, the percentage of beta-sheet content was found between 9.1% - 19.1% at pH 2 and 0.6% - 3.8% at pH 7, confirming that alpha-beta transition was more likely to occur at pH 2.

Table 1

The percentage of alpha-helix and beta-strand secondary structures for all simulation replicas at pH 2.

pH2 replicas	Structure	0 ns (%)	50 ns (%)	100 ns (%)	300 ns (%)
R0	α	32.5	4.1	0.1	0
	β	6.2	5.5	11.9	16
R1	α	32.5	20.9	3.3	3.5
	β	6.2	4.8	7.1	9.1
R2	α	32.5	0.8	0.2	4.4
	β	6.2	3.6	16.1	19.1

Table 2

The percentage of alpha-helix and beta-strand secondary structures for all simulation replicas at pH 7.

pH7 replicas	Structure	0 ns (%)	50 ns (%)	100 ns (%)	300 ns (%)
R0	α	34.1	11.9	5.9	0.9
	β	6.2	0.5	1.4	0.6
R1	α	34.1	2.6	8.3	13.2
	β	6.2	6.5	0.4	3.8
R2	α	34.1	2.2	4.3	0.1
	β	6.2	4.5	8.7	3.1

To further consider the effects of pH on the beta-strand formation of the simulated HEWL, the radius of gyration (Rg) was determined from different groups of amino acids. The distribution of positively-charged (Arg, Lys, and His) and negatively-charged (Glu and Asp) amino acids was affected by the protonation of all eight Glu and Asp residues at pH 2, which neutralized all the negatively-charged sidechains. Meanwhile, the Rg was calculated from all hydrophobic amino acids to monitor the propensity to form clusters. Fig 4 displayed the radius of gyration calculated along with all the 333 K simulations at both pH conditions for groups of hydrophobic amino acids, positively-charged amino acids (ARG-LYS-HIS), and negatively-charged amino acids (GLU-ASP). For all simulations at pH 2 (Fig 4A and 4C), Rg of the positively-charged ARG-LYS-HIS group was found above the Rg of the whole HEWL, while Rg of the groups of hydrophobic residue was found below the Rg of the whole HEWL. Rg of the protonated GLU-ASP group was located above the Rg of the whole HEWL for the R0 replica and was found below the Rg of the whole HEWL for the R1 and R2 replicas. At pH 7 (Fig 4D and 4F), Rg of the positively-charged ARG-LYS-HIS group was found slightly above the Rg of the whole HEWL, while the only slight difference in Rg was observed for the negatively-charged GLU-ASP group, hydrophobic group, and the whole HEWL.

Fig 4

Radius of gyration (Rg) calculated for the groups of whole HEWL structure (black), positively-charged residues (ARG-LYS-HIS; blue), negatively-charged residues (GLU-ASP; green), and hydrophobic residues (red) along all three 333-K simulation replicas of HEWL at (a-c) pH 2 and (d-f) pH 7. The negatively-charged residues were protonated at pH 2. Final snapshots highlighting C-alpha atoms of hydrophobic residues (red sphere) were also shown. Dash circles indicated the cluster of hydrophobic residues around the beta-strand regions at pH 2.

The significant difference between the Rg of the group of hydrophobic residues and the group of positively-charged residues (ARG-LYS-HIS) at pH 2 signified that hydrophobic residues tended to form a hydrophobic cluster within HEWL cores. In contrast, the positively charged residues tended to be more exposed to the solution. Final snapshots highlighting C-alpha atoms of hydrophobic residues in Fig 4A and 4C showed that the clusters of hydrophobic residues were presented along with the formed beta strands at pH 2. At pH 7, as Glu and Asp residues were deprotonated and the total charge of HEWL was reduced from +16e to +8e, the difference between the Rg of hydrophobic and charged groups became less significant. As a result, hydrophobic became more exposed to the solution and less likely to form a cluster at pH 7. The relationship between hydrophobicity and beta-strand formation could be discussed regarding the beta-strand sequences observed from simulations. Table 3 summarized the regions where beta-strands were found in final structures at pH 2, along with their amino acid sequences. Beta strands from final structures of all simulations at pH 2 involved 57 amino acids, of which 28 of them (49.12%) contained hydrophobic sidechains. The hydrophobic content within the beta-strands was clearly higher than the 35.66% of hydrophobic content within the whole HEWL. The YGILQINSRWW beta-strand found in the R0 replica at pH 2 contained 54.55% of hydrophobic content, even higher than the average hydrophobic content for beta-strands. This YGILQINSRWW sequence at residues 53–63 in one of our simulations covered the K-peptide region (residues 54–62) previously proposed as the core for amyloidosis [41-44]. An experimental study by Tokunaga et al. later showed that K-peptides were most likely to form fibrils under acidic conditions. Adding the STDY sequence from the residues 50–53 of HEWL further enhanced the amyloidosis under low pH [45].

Table 3

Residue numbers and amino acid sequences of beta-strand found in the final structure of HEWL from all simulation replicas at pH 2.

pH 2 Replica	Residue number of Beta
R0	24–36 (SLGNWVCAAKFES), 53–63 (YGILQINSRWW), 65–66 (ND), 69–70 (TP), 103–104 (NG), 112–113 (RN)
R1	29–31 (VCA), 83–85 (LLS), 99–101 (VSD), 107–109 (AWV)
R2	2–3 (VF), 6–10 (ELAA), 18–20 (DNY), 29–35 (VCAAKFE), 38–39 (FN), 75–77 (LCN), 109–110 (VA)

4. Conclusions

In this study, the beta-strand formation mechanism at the early stage of amyloid fibrilization has been proposed. Lysozyme is one of the amyloid proteins that can misfold into fibrils and cause some diseases. However, the controlled synthetic amyloid fibrils can become useful biomaterials for bioengineering and biosensing applications. Our molecular dynamics simulations showed that beta-strands were more likely to form when HEWL was unfolded at pH 2. Mechanisms of beta-strand formation were explained by the radial distribution of charged and the hydrophobic amino acids. The beta-strand forming sequences from the HEWL simulations contained a significantly higher hydrophobic content than the whole HEWL. This result was consistent with the relationship between hydrophobic clustering and amyloidosis previously reported by Mishima et al. [46]. Furthermore, at pH 2, a significant positive charge added by the protonation of all glutamic acid and aspartic acid residues corresponded with the separation of positively-charged amino acids from the group of hydrophobic residues, which promoted both hydrophobic clustering and beta-strand formation. However, further validation of this mechanistic scheme is still needed for the beta-strand formation of other amyloid proteins towards the control of fibril production and the prediction of nucleation sites.

Supplementary information file contains Ramachandran plots at pH 2 and pH 7.

(DOCX)

Click here for additional data file.

14 Jun 2021

PONE-D-21-08882

Molecular Dynamics Study on the Effects of Charged Amino Acid Distribution Under low pH Condition to the Unfolding of Hen Egg White Lysozyme and Formation of Beta Strands.

PLOS ONE

Dear Dr. Sutthibutpong,

Thank you for submitting your manuscript to PLOS ONE. After careful consideration, we feel that it has merit but does not fully meet PLOS ONE’s publication criteria as it currently stands. Therefore, we invite you to submit a revised version of the manuscript that addresses the points raised during the review process.

Please submit your revised manuscript by Jul 29 2021 11:59PM. If you will need more time than this to complete your revisions, please reply to this message or contact the journal office at plosone@plos.org. When you're ready to submit your revision, log on to https://www.editorialmanager.com/pone/ and select the 'Submissions Needing Revision' folder to locate your manuscript file.

Please include the following items when submitting your revised manuscript:

A rebuttal letter that responds to each point raised by the academic editor and reviewer(s). You should upload this letter as a separate file labeled 'Response to Reviewers'.

A marked-up copy of your manuscript that highlights changes made to the original version. You should upload this as a separate file labeled 'Revised Manuscript with Track Changes'.

An unmarked version of your revised paper without tracked changes. You should upload this as a separate file labeled 'Manuscript'.

If you would like to make changes to your financial disclosure, please include your updated statement in your cover letter. Guidelines for resubmitting your figure files are available below the reviewer comments at the end of this letter.

If applicable, we recommend that you deposit your laboratory protocols in protocols.io to enhance the reproducibility of your results. Protocols.io assigns your protocol its own identifier (DOI) so that it can be cited independently in the future. For instructions see:  http://journals.plos.org/plosone/s/submission-guidelines#loc-laboratory-protocols . Additionally, PLOS ONE offers an option for publishing peer-reviewed Lab Protocol articles, which describe protocols hosted on protocols.io. Read more information on sharing protocols at  https://plos.org/protocols?utm_medium=editorial-email&utm_source=authorletters&utm_campaign=protocols .

We look forward to receiving your revised manuscript.

Kind regards,

Human Rezaei

Academic Editor

PLOS ONE

Journal Requirements:

When submitting your revision, we need you to address these additional requirements.

1. Please ensure that your manuscript meets PLOS ONE's style requirements, including those for file naming. The PLOS ONE style templates can be found at

https://journals.plos.org/plosone/s/file?id=wjVg/PLOSOne_formatting_sample_main_body.pdf and

https://journals.plos.org/plosone/s/file?id=ba62/PLOSOne_formatting_sample_title_authors_affiliations.pdf

[Note: HTML markup is below. Please do not edit.]

Reviewers' comments:

Reviewer's Responses to Questions

Comments to the Author

1. Is the manuscript technically sound, and do the data support the conclusions?

The manuscript must describe a technically sound piece of scientific research with data that supports the conclusions. Experiments must have been conducted rigorously, with appropriate controls, replication, and sample sizes. The conclusions must be drawn appropriately based on the data presented.

Reviewer #1: Partly

**********

2. Has the statistical analysis been performed appropriately and rigorously?

Reviewer #1: No

**********

3. Have the authors made all data underlying the findings in their manuscript fully available?

The PLOS Data policy requires authors to make all data underlying the findings described in their manuscript fully available without restriction, with rare exception (please refer to the Data Availability Statement in the manuscript PDF file). The data should be provided as part of the manuscript or its supporting information, or deposited to a public repository. For example, in addition to summary statistics, the data points behind means, medians and variance measures should be available. If there are restrictions on publicly sharing data—e.g. participant privacy or use of data from a third party—those must be specified.

Reviewer #1: Yes

**********

4. Is the manuscript presented in an intelligible fashion and written in standard English?

PLOS ONE does not copyedit accepted manuscripts, so the language in submitted articles must be clear, correct, and unambiguous. Any typographical or grammatical errors should be corrected at revision, so please note any specific errors here.

Reviewer #1: Yes

**********

5. Review Comments to the Author

Please use the space provided to explain your answers to the questions above. You may also include additional comments for the author, including concerns about dual publication, research ethics, or publication ethics. (Please upload your review as an attachment if it exceeds 20,000 characters)

Reviewer #1: In this paper, Zein et al provide molecular dynamics simulations of Hen egg white lysozyme partial denaturation and refolding at high temperature. The protocol they present is a 100-ns simulation at temperature 450 K followed by 200 ns at 333 K. The authors use two protonation states, one where all titrable side chains (including glutamate, aspartate, histidine) are protonated ('pH 2'), and one where glutamates, aspartates and histidines are deprotonated ('pH 7'). They make three replicas of each state and analyze several parameters along the simulations : RMSD, secondary structures and various radii of gyration. The actual content of the paper is small and quite incremental with respect to existing data, both experimental and computational, on lysozyme unfolding and refolding and on amyloidosis. Still the authors make one convincing and interesting observation: That under their protocol protonated lysozyme consistently transitions from alpha-helix- to beta-sheet-containing structures, while there is no beta-sheet formation in the deprotonated form. A shortened paper focusing on this observation would be more suitable for publication.

Pages are not numbered so the following refers to the pages of the consolidated pdf, with the introduction starting on page 9.

1) p. 10 'The accuracy of atomistic molecular dynamics (MD) simulations to predict the

molecular behavior of proteins under extreme conditions [24,25] has been improved with the

continuing development of molecular mechanics forcefield parameters [26,27].'

I agree with this, but the authors should at least indicate in the methods which forcefield they used, and if there was a reason (suitability for high temperature ?) for the choice.

Similarly p.11 'the SPC216 water model': Was there a reason for this water model ?

2) 'Table 1: All simulations in this study' is not useful. As stated there are just two nearly identical systems simulated in the same way with three replicas each. Please remove.

3) p. 11 'At the pH 2 condition, all histidines, glutamic

acids, and aspartic acids were fully protonated. Meanwhile, at pH 7, all the aforementioned amino acids were deprotonated.'

Was the C-terminus considered ?

4) p. 14

pH 2 'At this stage, the slight drop of RMSD values denoted the

refolding of proteins.'

pH7 : 'No significant drop of RMSD was seen when switching the temperature to 333 K, and hence

no protein refolding.'

a) This is not obvious from the RMSD and analyses provided. The parts about protein refolding are not substantiated and should be removed.

b) Generally speaking the results are full of discussion elements and speculations such as these. The manuscript thus can and should be drastically shortened (see below points 8 and 9).

5) p. 15 'The transitions from alpha-helixes, which were the major part of native lysozymes, into

beta-strands and random coils were quantified through the Ramachandran plot'

There is no quantification given from figures 2 and 3 (Ramachandran plots). Quantifications are given from figures 4 and 5 (DSSP plots) only. All elements given in the Ramachandran part are qualitative, such as :

'The alpha-beta transitions occurred less frequently for the non-protonated proteins at

pH7 (Figure 3), in which a large number of amino acid'

This part should be shortened, removing particularly interpretations such as 'signified the incomplete refolding process as observed from the conformational snapshots'

6) DSSP plots

p. 16 'the percentage of beta-sheet content was

found between 9.1% - 19.1% at pH 2 and 0.6% - 3.8% at pH 7,'

This is the interesting and substantial result the authors bring and should be the focus of the results section. As it is, it is clear from the DSSP plots that the same residues are not found in the same conformation in different replicas. This should at least be stated in the text. Preferably residues that tend to be incorporated in beta-strands should de identified, in relation to the 'radii of gyration' plots (see below).

7) Radii of gyration

That part is the least clear and convincing of the paper. A major problem is that it is not clear which selections statements such as 'negatively-charged

amino acids affected by the protonations.' 'positively charged, negatively-charged, and hydrophobic compositions' refer to.

a) I take it from Figure 6 and 7 labels ('Negative-protonated' instead of 'Negative') that the same selections (namely asp+glu) are considered at pH 2 and pH 7 ? Please then refer explicitly to these as 'asp+glu', both in the text and the figures.

b) What about 'positive' residues ? Are histidine in the 'positive' selection at pH 2 and pH 7 ? Only at pH 2 ? 'Rg of the positively charged group was the highest Rg due to a large number of positive charge residues.' suggests that it is the latter. But then, no meaningful comparison can be made between 'positive' radii of gyration at pH 7 vs pH 2 (different selections). At any rate, please also make explicit descriptions of the 'positive' selection (whether arg+lys or arg+lys+his).

c) p. 17 'However, at pH 7

(Figure 6d, 6e and 6f), Rg values of different amino acid groups tended to converge as the

simulation progressed.'

This is simply false. Only for R1 is it verified, for R0 and R2 'positive' radii stay above others as for pH 2.

d) p. 17

'the formation of beta-strands

or beta-sheets from backbone parts of hydrophobic amino acids was facilitated by the increased

compactness (low Rg) of the hydrophobic clusters'

No results are given as to the composition of beta strands (see point 6 above) and whether they tend to be made of hydrophobic residues. Figure 8 is illustrative at best. The second part of the sentence is interpretation.

8) The whole section '4. Discussions' is an extended repetition of results. For instance

p. 20 'positively charged residues formed a hydrophilic shell with larger Rg than the averaged Rg' (and as stated this particular point is moot anyway since it is also the case at pH 7 for 2 replicas out of 3). Furthermore, the 'results' section itself is riddled with interpretations. The easiest way to amend this is to remove section '4. Discussions' completely, and change the section '3. Results' into '3. Results and discussion'

9) In this new 'Results and discussion', please avoid speculations and interpretations that are not substantiated by the results themselves. For instance

a) p. 20 'The higher amount of Coulombic repulsion at

lower pH had driven most of the positively charged sidechain further from the backbone,

leaving the backbone to stay at the middle between hydrophobic and hydrophilic shells. The

beta-strands were finally formed by nucleation of the ordered backbone part.'

No analyses are provided to support this, not even the justification of a visual inspection of the backbone.

b) p. 18

'to refold into betastrands

at pH 2, while misfolded into alpha-helices and random coils at pH 7'

Why speak of 'refolding' into beta-strands and 'misfolding' into alpha-helices for an alpha-helical protein ?

10)Minor: There are many typos throughout the manuscript such as'to observed the effects', 'analyzedanalyzed by the DSSP algorithm'. Please proofread carefully.

**********

6. PLOS authors have the option to publish the peer review history of their article (what does this mean?). If published, this will include your full peer review and any attached files.

If you choose “no”, your identity will remain anonymous but your review may still be made public.

Do you want your identity to be public for this peer review? For information about this choice, including consent withdrawal, please see our Privacy Policy.

Reviewer #1: No

[NOTE: If reviewer comments were submitted as an attachment file, they will be attached to this email and accessible via the submission site. Please log into your account, locate the manuscript record, and check for the action link "View Attachments". If this link does not appear, there are no attachment files.]

While revising your submission, please upload your figure files to the Preflight Analysis and Conversion Engine (PACE) digital diagnostic tool, https://pacev2.apexcovantage.com/. PACE helps ensure that figures meet PLOS requirements. To use PACE, you must first register as a user. Registration is free. Then, login and navigate to the UPLOAD tab, where you will find detailed instructions on how to use the tool. If you encounter any issues or have any questions when using PACE, please email PLOS at figures@plos.org. Please note that Supporting Information files do not need this step.

30 Jul 2021

We thank the Editor-in-chief for giving us an opportunity to revise the manuscript and the reviewer for his/her comments and suggestions. The paper is revised according to comments and suggestions of the reviewer and we have highlighted the changes made in the revised manuscript with the red color.

**************************************************************************

Reviewer 1’s: In this paper, Zein et al provide molecular dynamics simulations of Hen egg white lysozyme partial denaturation and refolding at high temperature. The protocol they present is a 100-ns simulation at temperature 450 K followed by 200 ns at 333 K. The authors use two protonation states, one where all titrable side chains (including glutamate, aspartate, histidine) are protonated ('pH 2'), and one where glutamates, aspartates and histidines are deprotonated ('pH 7'). They make three replicas of each state and analyze several parameters along the simulations : RMSD, secondary structures and various radii of gyration. The actual content of the paper is small and quite incremental with respect to existing data, both experimental and computational, on lysozyme unfolding and refolding and on amyloidosis. Still the authors make one convincing and interesting observation: That under their protocol protonated lysozyme consistently transitions from alpha-helix- to beta-sheet-containing structures, while there is no beta-sheet formation in the deprotonated form.

A shortened paper focusing on this observation would be more suitable for publication. Pages are not numbered so the following refers to the pages of the consolidated pdf, with the introduction starting on page 9.

**************************************************************************

Comment 1:

p. 10 (p.5 in the new version) 'The accuracy of atomistic molecular dynamics (MD) simulations to predict the molecular behavior of proteins under extreme conditions [24,25] has been improved with the continuing development of molecular mechanics forcefield parameters [26,27].' I agree with this, but the authors should at least indicate in the methods which forcefield they used, and if there was a reason (suitability for high temperature ?) for the choice.

Similarly p.11 'the SPC216 water model': Was there a reason for this water model ?

Response: Our choice of forcefield is GROMOS54a7, which has been added into the methodology section (p.5 in the new version). The reason behind choosing this forcefield comes from a benchmarking study (also cited in the paper) of Kamenik et al. 2020, where proteins simulated by GROMOS54a7 forcefield showed a good agreement with experimental state population and unfolding time. Then, the SPC water model was used along with the GROMOS forcefield as the GROMOS forcefield itself was parameterized to reproduce the free enthalpy of hydration in SPC water (‘SPC216’ was the name of the file so we changed it back to the actual name of the model ‘SPC’). Even if SPC is known to incorrectly reproduce the nature of lipid bilayer, there should be no problem as the lipid bilayers do not exist in our models.

**************************************************************************

Comment 2:

'Table 1: All simulations in this study' is not useful. As stated there are just two nearly identical systems simulated in the same way with three replicas each. Please remove.

Response: The table has been removed.

**************************************************************************

Comment 3:

p. 11 'At the pH 2 condition, all histidines, glutamic acids, and aspartic acids were fully protonated. Meanwhile, at pH 7, all the aforementioned amino acids were deprotonated.' Was the C-terminus considered ?

Response: We missed out considering the terminus but we believe that it will not significantly affect the protein that already contains eight protonable sites.

**************************************************************************

Comment 4:

p. 14 pH 2 'At this stage, the slight drop of RMSD values denoted the refolding of proteins.' pH7 : 'No significant drop of RMSD was seen when switching the temperature to 333 K, and hence no protein refolding.' a) This is not obvious from the RMSD and analyses provided. The parts about protein refolding are not substantiated and should be removed. b) Generally speaking the results are full of discussion elements and speculations such as these. The manuscript thus can and should be drastically shortened (see below points 8 and 9).

Response: a) We removed all explanations about protein refolding from the manuscript, but chose to present the events as the formation of beta-strands. RMSD results were discussed along with the snapshot of the HEWL at the start, 50 ns, and 100 ns of 450-K simulations to monitor the conformational change of the proteins at different pH (see Figure 1 and p.7). b) In this version, we already tried to avoid speculative assumptions, e.g. protein refolding after the temperature switching and Coulombic force that drive positively charged residues away. These changes are described more in the comments below.

**************************************************************************

Comment 5:

p. 15 'The transitions from alpha-helices, which were the major part of native lysozymes, into beta-strands and random coils were quantified through the Ramachandran plot' There is no quantification given from figures 2 and 3 (Ramachandran plots). Quantifications are given from figures 4 and 5 (DSSP plots) only. All elements given in the Ramachandran part are qualitative, such as : 'The alpha-beta transitions occurred less frequently for the non-protonated proteins at pH7 (Figure 3), in which a large number of amino acid' This part should be shortened, removing particularly interpretations such as 'signified the incomplete refolding process as observed from the conformational snapshots'

Response: We discussed more about quantitative measurement on the secondary structure from DSSP results (Table 1 and Table 2) including the percentage, and the amino acid residues that formed beta strand in the final structure (see Table 3). And we removed the discussion about the Ramachandran plot, as it was already covered by the DSSP result - which should shorten the manuscript.

**************************************************************************

Comment 6:

DSSP plots p. 16 'the percentage of beta-sheet content was found between 9.1% - 19.1% at pH 2 and 0.6% - 3.8% at pH 7,' This is the interesting and substantial result the authors bring and should be the focus of the results section. As it is, it is clear from the DSSP plots that the same residues are not found in the same conformation in different replicas. This should at least be stated in the text. Preferably residues that tend to be incorporated in beta-strands should be identified, in relation to the 'radii of gyration' plots (see below).

Response: We have provided the tables that summarized the percentage of all replicas and the amino acid residues that tend to form beta-strands. We also discussed what happened for each replica at both pH (see p.8 in the new version), the sequences of all beta-strands formed and their hydrophobicity (see p.10 in the new version). As the beta-strand forming sequences contained more hydrophobic residues than other regions. Hydrophobic residues tends to form clusters, resulting in the enhanced compactness but lower radius of gyration (Rg)

**************************************************************************

Comment 7:

Radii of gyration That part is the least clear and convincing of the paper. A major problem is that it is not clear which selections statements such as 'negatively-charged amino acids affected by the protonations.' 'positively charged, negatively-charged, and hydrophobic compositions' refer to.

a) I take it from Figure 6 and 7 labels ('Negative-protonated' instead of 'Negative') that the same selections (namely asp+glu) are considered at pH 2 and pH 7 ? Please then refer explicitly to these as 'asp+glu', both in the text and the figures.

Response: Negatively charged amino acids refer to GLU and ASP amino acids at both pH. At pH2 GLU and ASP amino acid were affected by the protonation of all eight Glu and Asp residues.

b) What about 'positive' residues ? Are histidine in the 'positive' selection at pH 2 and pH 7 ? Only at pH 2 ? 'Rg of the positively charged group was the highest Rg due to a large number of positive charge residues.' suggests that it is the latter. But then, no meaningful comparison can be made between 'positive' radii of gyration at pH 7 vs pH 2 (different selections). At any rate, please also make explicit descriptions of the 'positive' selection (whether arg+lys or arg+lys+his).

Response: Positive residues refer to ARG, LYS, and HIS and were mentioned as ARG-LYS-HIS in the Figure 4 and the whole manuscript.

c) p. 17 'However, at pH 7 (Figure 6d, 6e and 6f), Rg values of different amino acid groups tended to converge as the simulation progressed.' This is simply false. Only for R1 is it verified, for R0 and R2 'positive' radii stay above others as for pH 2.

Response: We have now removed this discussions. In Figure 4 of the revised version, we discussed only the Rg from 333 K simulations. At pH 7, Rg of every groups was not significantly different. Meanwhile, the significant difference in Rg between the positively charged group and the hydrophobic group suggested that the positively charged group was separated from the hydrophobic core. This was discussed in p.10 of the revised version

d) p. 17 'the formation of beta-strands or beta-sheets from backbone parts of hydrophobic amino acids was facilitated by the increased compactness (low Rg) of the hydrophobic clusters' No results are given as to the composition of beta strands (see point 6 above) and whether they tend to be made of hydrophobic residues. Figure 8 is illustrative at best. The second part of the sentence is interpretation.

Response: In the revised version, %content of hydrophobic residues within the beta-strands was analysed in comparison with the whole HEWL structure (p.10 of the new version). Beta-strands contained more hydrophobic content than the whole HEWL structure.

**************************************************************************

Comment 8:

The whole section '4. Discussions' is an extended repetition of results. For instance p. 20 'positively charged residues formed a hydrophilic shell with larger Rg than the averaged Rg' (and as stated this particular point is moot anyway since it is also the case at pH 7 for 2 replicas out of 3). Furthermore, the 'results' section itself is riddled with interpretations. The easiest way to amend this is to remove section '4. Discussions' completely, and change the section '3. Results' into '3. Results and discussion'

Response: We have now merged section ‘3. Result’ and ‘4. Discussion’ into section ‘3. Result and Discussion’.

**************************************************************************

Comment 9:

In this new 'Results and discussion', please avoid speculations and interpretations that are not substantiated by the results themselves. For instance

a) p. 20 'The higher amount of Coulombic repulsion at lower pH had driven most of the positively charged sidechain further from the backbone, leaving the backbone to stay at the middle between hydrophobic and hydrophilic shells. The beta-strands were finally formed by nucleation of the ordered backbone part.' No analyses are provided to support this, not even the justification of a visual inspection of the backbone.

Response: We left out the discussion on the Coulombic interactions but only mentioning the separation between positively-charged amino acids from the hydrophobic core (judging from significant difference in Rg at pH 2; p.10 of the revised version). We also tried to visualize the tendency for hydrophobic clustering in Figure 4.

b) p. 18 'to refold into betastrands at pH 2, while misfolded into alpha-helices and random coils at pH 7' Why speak of 'refolding' into beta-strands and 'misfolding' into alpha-helices for an alpha-helical protein ?

Response: We already left out this statement, only reporting time-dependent secondary structure content from DSSP.

**************************************************************************

Comment 10:

Minor: There are many typos throughout the manuscript such as 'to observed the effects', 'analyzedanalyzed by the DSSP algorithm'. Please proofread carefully.

Response: The draft has been thoroughly proofread. Thank you very much.

**************************************************************************

Submitted filename: Reviewers-1.docx

Click here for additional data file.

24 Feb 2022

Molecular Dynamics Study on the Effects of Charged Amino Acid Distribution Under low pH Condition to the Unfolding of Hen Egg White Lysozyme and Formation of Beta Strands.

PONE-D-21-08882R1

Dear Dr. Sutthibutpong,

We’re pleased to inform you that your manuscript has been judged scientifically suitable for publication and will be formally accepted for publication once it meets all outstanding technical requirements.

Within one week, you’ll receive an e-mail detailing the required amendments. When these have been addressed, you’ll receive a formal acceptance letter and your manuscript will be scheduled for publication.

An invoice for payment will follow shortly after the formal acceptance. To ensure an efficient process, please log into Editorial Manager at http://www.editorialmanager.com/pone/, click the 'Update My Information' link at the top of the page, and double check that your user information is up-to-date. If you have any billing related questions, please contact our Author Billing department directly at authorbilling@plos.org.

If your institution or institutions have a press office, please notify them about your upcoming paper to help maximize its impact. If they’ll be preparing press materials, please inform our press team as soon as possible -- no later than 48 hours after receiving the formal acceptance. Your manuscript will remain under strict press embargo until 2 pm Eastern Time on the date of publication. For more information, please contact onepress@plos.org.

Kind regards,

Hannes C Schniepp, Dr. sc. nat.

Academic Editor

PLOS ONE

Additional Editor Comments (optional):

Reviewers' comments:

Reviewer's Responses to Questions

Comments to the Author

1. If the authors have adequately addressed your comments raised in a previous round of review and you feel that this manuscript is now acceptable for publication, you may indicate that here to bypass the “Comments to the Author” section, enter your conflict of interest statement in the “Confidential to Editor” section, and submit your "Accept" recommendation.

Reviewer #1: All comments have been addressed

**********

2. Is the manuscript technically sound, and do the data support the conclusions?

The manuscript must describe a technically sound piece of scientific research with data that supports the conclusions. Experiments must have been conducted rigorously, with appropriate controls, replication, and sample sizes. The conclusions must be drawn appropriately based on the data presented.

Reviewer #1: Yes

**********

3. Has the statistical analysis been performed appropriately and rigorously?

Reviewer #1: N/A

**********

4. Have the authors made all data underlying the findings in their manuscript fully available?

The PLOS Data policy requires authors to make all data underlying the findings described in their manuscript fully available without restriction, with rare exception (please refer to the Data Availability Statement in the manuscript PDF file). The data should be provided as part of the manuscript or its supporting information, or deposited to a public repository. For example, in addition to summary statistics, the data points behind means, medians and variance measures should be available. If there are restrictions on publicly sharing data—e.g. participant privacy or use of data from a third party—those must be specified.

Reviewer #1: Yes

**********

5. Is the manuscript presented in an intelligible fashion and written in standard English?

PLOS ONE does not copyedit accepted manuscripts, so the language in submitted articles must be clear, correct, and unambiguous. Any typographical or grammatical errors should be corrected at revision, so please note any specific errors here.

Reviewer #1: Yes

**********

6. Review Comments to the Author

Please use the space provided to explain your answers to the questions above. You may also include additional comments for the author, including concerns about dual publication, research ethics, or publication ethics. (Please upload your review as an attachment if it exceeds 20,000 characters)

Reviewer #1: (No Response)

**********

7. PLOS authors have the option to publish the peer review history of their article (what does this mean?). If published, this will include your full peer review and any attached files.

If you choose “no”, your identity will remain anonymous but your review may still be made public.

Do you want your identity to be public for this peer review? For information about this choice, including consent withdrawal, please see our Privacy Policy.

Reviewer #1: No

15 Mar 2022

PONE-D-21-08882R1

Molecular Dynamics Study on the Effects of Charged Amino Acid Distribution Under low pH Condition to the Unfolding of Hen Egg White Lysozyme and Formation of Beta Strands.

Dear Dr. Sutthibutpong:

I'm pleased to inform you that your manuscript has been deemed suitable for publication in PLOS ONE. Congratulations! Your manuscript is now with our production department.

If your institution or institutions have a press office, please let them know about your upcoming paper now to help maximize its impact. If they'll be preparing press materials, please inform our press team within the next 48 hours. Your manuscript will remain under strict press embargo until 2 pm Eastern Time on the date of publication. For more information please contact onepress@plos.org.

If we can help with anything else, please email us at plosone@plos.org.

Thank you for submitting your work to PLOS ONE and supporting open access.

Kind regards,

PLOS ONE Editorial Office Staff

on behalf of

Dr. Hannes C Schniepp

Academic Editor

PLOS ONE