\chapter{Effect of Temperature on Peptide Aggregation}
The phase state of a multicomponent solution at a given temperature is  the result of  the interplay of various interactions between the constituent particles. As the temperature is increased, the kinetic energy should become dominant over the potential energy, and at a certain temperature, the solution should become macroscopically homogenous.
\par   
Amyloidogenic peptides are highly insoluble in liquid water and their critical concentration is generally expected to be in the micromolar to nanomolar range.\cite{Harper1997} Under ambient conditions, their aqueous solutions are separated into a water-rich phase, containing the critical peptide concentration, and an organic-rich phase, which appears as amyloid fibrils. It can be expected that amyloid fibrils dissolve in liquid water upon heating at a particular temperature and pressure, depending on the peptide concentration.  In this chapter, the temperature-induced demixing phase transition between the aggregated and the disaggregated states of the peptide system was examined with the REMD simulation method, using an explicit water model.
\label{ch:replica}
\section{System setup}
The simulated system consists of 12 peptide fragments having the amino acid sequence FLVHS, corresponding to residues 15-19 of hIAPP, and 3079 water molecules, resulting in a peptide concentration of 12.4 w/w\% . The N- and C-termini of the peptide fragments were capped with acetyl and methylamide groups, respectively. A modified AMBER force field \cite{Simmerling2002} was used for the peptides along with the TIP3P water model.\cite{Jorgensen:83a} The Replica exchange molecular dynamics (REMD) simulation was performed using the RPMDRUN package.\cite{paschek} The temperature spacing was calculated by performing initial uncoupled simulations for 0.5 ns at 36 different temperatures ranging from 275  to 991 K. The average energies were fitted with a polynomial and Eq. \ref{eq:rep1} was solved iteratively for the temperature distribution for which $P(acc) = 0.20$. Sixty replicas distributed over a temperature range of 286.7 K to 645.5 K were used with a state exchange probability of 0.1 leading to a time of about 3 ps for each replica between two state exchanges. The SHAKE algorithm was applied to constrain all bonds in order to allow a time step of 2.0 fs. The electrostatic interactions were treated by Particle Mesh Ewald (PME) summation and the simulation was performed for 49 ns, corresponding to a total simulation time of 2.94 \textmu s.
\par
The intrinsic volume $V_p$ of the peptides was calculated as:
\begin{equation}
V_p = L_{box}^3 - V_h - (N_0 - N_w)/ \rho_b,
\label{eq:volp}
\end{equation}
where $V_h$ is the volume of the hydration layer, calculated using the slab approximation (see section \ref{bulksetup} of Chapter \ref{ch:bulkphase} for details), $N_0$ and $N_w$ are the numbers of water molecules in the cubic simulation box (with the edge $L_{box}$) and in the hydration shell, respectively, and $\rho_b$ is the numerical density of bulk water determined at the same thermodynamic states in constant-volume simulations of TIP3P water. The radius of gyration, hydrophobic contacts, hydrogen bonds, and secondary structure were calculated using the methods described in section \ref{bulksetup} of Chapter \ref{ch:bulkphase}.
%
\section {Results}
In order to study the equilibrium thermodynamic states of the system, the part of the MD trajectory where system properties vary monotonically with time should be excluded. The time dependence of various system properties was analyzed by dividing the simulation run into five non-overlapping segments and calculating the average properties of each segment. Of all the measured properties, the inter-peptide hydrophobic contacts and hydrogen bonds needed the longest equilibration time ($\approx$ 25 ns). Other parameters, such as the SASA and the radius of gyration, equilibrated much faster. Figure \ref{fig:repequi} shows the average number of inter-peptide hydrophobic contacts and SASA of all peptides for each segment of the simulation run. 
\begin{figure}[htb]
\begin{center}
\includegraphics[trim= 0.5cm 0.6cm 0.5cm 0.2cm,clip,width=15.00 cm]{repequi}
\end{center}
\caption{Average number of inter-peptide hydrophobic contacts and SASA as a function of temperature in four non-equilibrium segments of the total MD simulation run and in the equilibrium run between 26 to 49 ns. }
\label{fig:repequi}
\end{figure}
Based on the time evolution of all parameters studied, the equilibration period was estimated to be 25 ns. Only the equilibrium results obtained from the last 24 ns of the REMD simulation run are further discussed.
\par
The radius of gyration $R_g$ of all peptides and its mean-square fluctuations ($\sigma$ = $\left\langle R_g^2 \right\rangle$ - $\left\langle R_g \right\rangle$ $^2$) were shown to be useful parameters for characterizing the phase transition in aqueous solution of peptides.\cite{Cecchini2004,Baumketner2005} The $R_g$ shows a sigmoid-like increase with increasing temperature, indicating a transition from the state with compact aggregate with $R_g \approx$ 1.25 nm to a disaggregated state with $R_g \approx$ 2 nm (Figure \ref{fig:reprg}). The fit of the $R_g$(\textit{T}) to a sigmoidal function (dashed line) indicates a midpoint of the transition at $T \approx$ 450 K (vertical dashed line) and a smearing of the transition within a range $\pm$ 55 K. A transition temperature can be estimated using the mean-square fluctuations $\sigma$, which pass through a maximum at $T = $ 486 K (vertical solid line). The discrepancy between the two estimations is $\sim$ 35 K, which could be due to notable deviations of $R_g$(\textit{T}) from a sigmoidal shape. The behavior of the radius of gyration as a function of the reversed temperature 1/\textit{T} was also analyzed (Figure \ref{fig:reprginvt}, solid squares). The $R_g$(1/\textit{T}) dependence can be almost perfectly described by a sigmoidal function, with a better goodness of fit as compared to  that of the $R_g$(\textit{T}) dependence.  The inflection point was found to be at $T \approx$ 465 K ( Figure \ref{fig:reprginvt}, solid vertical line).
\begin{figure}[htb]
\centering
\begin{minipage}[t]{7.5cm}
\includegraphics[height=5.6 cm]{reprg}
\caption{Temperature dependence of the radius of gyration $R_g$ of all peptides (squares) and its mean-square fluctuations $\sigma$ (circles). The temperature where $\sigma$ passes through a maximum is indicated by a vertical solid line. The sigmoid fit of $R_g(T)$ and location of its inflection point are shown by dashed lines.  }
\label{fig:reprg}
\end{minipage}
\hfill
\begin{minipage}[t]{7.5cm}
\includegraphics[height=5.6 cm,width=5.6 cm]{reprginvt}
\caption{Radius of gyration $R_g$ and SASA of all peptides as functions of reversed temperature $1/T$ (symbols). The fits to a sigmoid function and locations of the inflection points are shown by solid and dashed lines, respectively }
\label{fig:reprginvt}
\end{minipage}
\end{figure}
Similarly, the dependence of SASA on the reversed temperature 1/\textit{T} is much closer to a sigmoid (Figure \ref{fig:reprginvt}) and the goodness of fit is about three times better as compared to the SASA(\textit{T}) dependence. The inflection of the sigmoidal fit of SASA(1/\textit{T}) is located at $T \cong$ 480 K (Figure \ref{fig:reprginvt}, vertical dashed line).  
\par
Peptide aggregation can be quantitatively characterized utilizing cluster analysis. The probability distributions \textit{P}($S_{max}$) to find the largest peptide cluster containing $S_{max}$ peptides, obtained with a connectivity criterion based on the distance between the centers of mass of the peptides, are shown in Figure \ref{fig:repsmax} for selected temperatures. At the lowest temperature studied (\textit{T} = 287 K), the largest cluster consists of all 12 peptides in $\sim$ 50$\%$ of the observed configurations and only rarely contains less than 5 peptides. Conversely, at the highest temperature studied (\textit{T} = 646 K), largest clusters consisting of just 2 or 3 peptides dominate. At \textit{T} = 440 K, the largest peptide cluster strongly fluctuates, consisting of 4 to 10 peptides with comparable probabilities. The temperature dependence of the probabilities to have an aggregate including more than 2/3 of the peptides, based on the connectivity criteria of either the distance between the centers of mass of two peptides ($r_{crit} \leq$ 1.1 nm) or the number of hydrophobic contacts between the two peptides ($n_c \geq$ 10), are shown in Figure \ref{fig:repdrop}. Fits to a sigmoidal function give inflection points at 415  and 435 K, respectively.
\begin{figure}[htb]
\centering
\begin{minipage}[t]{7.5cm}
\includegraphics[trim= 0.5cm 0.5cm 0cm 0.2cm,clip,width=9.3 cm]{repsmax}
\caption{Probability that the largest cluster consists of $S_{max}$ peptides, calculated at $T$ = 287, 348, 376, 401, 440, 501, 603 and 646 K.}
\label{fig:repsmax}
\end{minipage}
\hfill
\begin{minipage}[t]{7.5cm}
\includegraphics[width=7.5 cm]{repdrop}
\caption{Temperature dependence of the probabilitiy that the largest peptide cluster includes $>$ 2/3 of all molecules (more than 8 peptides). The criterion for the connectivity between two peptides was based on the distance between the centers of mass of the peptides (red circles) and on the number of hydrophobic contacts (blue circles). Fits to the sigmoidal function are shown by solid lines. The inflection points are indicated by vertical lines.}
\label{fig:repdrop}
\end{minipage}
\end{figure}
\par
The cluster size distribution $n_S$, which is the probability to find a cluster containing \textit{S} peptides, can also serve for characterizing aggregation. Examples of the $n_S$ distributions are shown in Figure \ref{fig:repclustd} using a double logarithmic scale. Such distributions are very sensitive to a phase transition. In the one-phase region, where the aqueous solution is homogeneous at the macroscopic scale, the distribution of $n_S$ deviates downwards from the power law $n_S \sim S^{-2.2}$.\cite{Fisher1967} Such a behavior is seen from  \textit{T} = 646 K up to 566 K. Phase separation becomes evident with the appearance of a peptide droplet, and is indicated by a sharp increase of $n_S$ at  large \textsl{S} values, e.g., the distributions from \textit{T} = 287 K to 340 K. Near the phase transition, $n_S$ shows a power-law behavior $\sim S^{-2.2}$ in the widest range of \textit{S}. However, such a scenario is only true for dilute macroscopic systems. Due to the finite size of the simulation box, the probability distribution of  large clusters gets distorted and a hump appears on the $n_S$ distribution at large \textit{S} values. This hump does not indicate the presence of a droplet, however, and can also be observed at  temperatures much higher than that of the phase transition. These humps can be seen for temperatures $T <$ 557 K, therefore, 557 K should be considered as an absolute upper temperature limit for the presence of an aggregate.
\begin{figure}[htb]
\centering
\includegraphics[width=6.0 cm]{repclustd}
\caption{Size distribution $n_S$ of the peptide clusters for selected temperatures (symbols) in a double-logarithmic plot. The power law $n_S \sim S^{-2.2}$ is shown by lines. }
\label{fig:repclustd}
\end{figure}
\par
The average number of peptide-water hydrogen bonds per peptide, $n_H^{wp}$, closely follows the temperature dependence of SASA in the temperature range 287 K $< T <$ 430 K (Figure \ref{fig:repsasa}). In particular, one additional water-peptide H-bond corresponds to an increase of SASA of about 1.1 nm$^2$. The temperature induced increase of SASA above 430 K is not accompanied by the creation of new water-peptide H-bonds. This means that the peptide is not able to form H-bonds with more than about 14 water molecules. The effective hydrophilicity of the surface can be characterized by the number of H-bonds with water per unit SASA. The increase of the ratio $n_H^{wp}$/SASA upon cooling indicates an increasing hydrophilicity of the exposed surface upon aggregation (Figure \ref{fig:repsasa}, lower panel). These results are in agreement with previous studies of aggregation of the same peptide at \textit{T} = 330 K, proving that the aggregation of the peptide studied occurs mainly via contacts between hydrophobic groups. The average number of intermolecular peptide-peptide H-bonds per peptide, $n_H^{pp}$, is $\sim$ 2.5 at low temperatures (Figure \ref{fig:rephbond}, upper panel), continuously decreases and finally vanishes upon heating, indicating that  most of the peptides are in monomeric form at high temperatures. The average number of intra-peptide H-bonds, $n_H^{intra}$, increases slightly up to 450 K, beyond which it remains almost unchanged. The population of various elements of secondary structure is shown as a function of temperature in the lower panel of Figure \ref{fig:rephbond}. In the aggregated state, which dominates at low temperatures, $\beta$-strands are the dominating structural elements. The fraction of the PP$_{\texttt{II}}$ structure is about 25$\%$ at low temperatures and decreases by a few percents upon heating. The population of random coil conformations increases from 18 to 26 $\%$ upon heating, whereas the $\alpha$-helix content shows a non-monotonic behavior, i.e.,  it increases from $\sim$ 19$\%$ at \textit{T} = 287 K, passes through a maximum ($\sim$ 23$\%$) at $T \approx$ 475 K,  and drops to $\sim$ 21$\%$ at the highest temperature studied.
\begin{figure}[htb]
\centering
\begin{minipage}[t]{7.5cm}
\includegraphics[height=7.8 cm]{repsasa}
\caption{Temperature dependences of the SASA of all peptides (open squares), the number of peptide-water H-bonds per peptide $n_H^{wp}$ (open circles), and their ratio (open triangles).}
\label{fig:repsasa}
\end{minipage}
\hfill
\begin{minipage}[t]{7.5cm}
\includegraphics[trim= 0.5cm 0.5cm 0.8cm 0.2cm,clip,height=7.8cm]{rephbond}
\caption{ Number of inter- ($n_H^{pp}$) and intra-peptide H-bonds ($n_H^{intra}$), and  populations of various elements of secondary structure, as a function of temperature are shown in upper and lower panels, respectively.}
\label{fig:rephbond}
\end{minipage}
\end{figure}
The similarity of the temperature dependences of the $\beta$-strand content and the number $n_H^{pp}$ of inter-peptide H-bonds (Figure \ref{fig:rephbond}) indicates their coupling, which is expected as the $\beta$-sheets are formed due to the creation of inter-peptide H-bonds. The degree of such  coupling is clearly seen, when the $\beta$-strand content is shown as a function of inter-peptide H-bonds per peptide, $n_H^{pp}$ (Figure \ref{fig:repcorr}, lower panel). The dependence is linear with a slope of $\sim$ 0.036 over the entire  temperature range studied. As the presence of intra-molecular H-bonds is an important condition for the existence of $\alpha$-helical structures, the $\alpha$-helical content as a function of $n_H^{intra}$ is shown in the upper panel of Figure \ref{fig:repcorr}. Three temperature intervals with an almost linear relation between the $\alpha$-helical content and $n_H^{intra}$ (dashed line) can be distinguished. Up to about 450 K, both properties increase with temperature due to the dissolution of the peptide aggregate, and a qualitative change of the correlation between the $\alpha$-helical content and $n_H^{intra}$ is observed at $\sim$ 375 K. Upon heating beyond 475 K, the $\alpha$-helical content starts to decrease, whereas $n_H^{intra}$ remains almost unchanged.
\begin{figure}[htb]
\centering
\includegraphics[width=7.00 cm]{repcorr}
\caption{Correlations between $\alpha$-helix and $\beta$-strand contents and the numbers of intra- ($n_H^{intra}$) and inter-peptide H-bonds ($n_H^{pp}$), respectively. }
\label{fig:repcorr}
\end{figure}

\par
Additional information about the peptide surface and the character of peptide aggregation can be obtained from the behavior of the volumetric properties of the hydration water and the intrinsic volumetric properties of the peptides.\cite{Brovchenko:08b}  The density of water in the hydration shells of the  peptides, $\rho_h$ = $N_w$/$V_h \cdot m_w$, where $m_w$ is the mass of a water molecule, is shown in Figure \ref{fig:reprhoshell} as a function of temperature.
\begin{figure}[htb]
\centering
\includegraphics[width=14.0 cm]{reprhoshell}
\caption{Dependence of the water density in the bulk $\rho_b$ (solid circles) and in the hydration shell $\rho_h$ (circles) on temperature along the simulated thermodynamic path. Pressure-temperature diagram of REMD simulations (right panel). }
\label{fig:reprhoshell}
\end{figure}
$\rho_h$ was found to be notably  lower ($\sim$ 25-30 $\%$) than the bulk density, $\rho_b$, in the same thermodynamic states, in agreement with the pronounced hydrophobic nature of the peptide surface. The densities of both the bulk and the hydration water strongly depend on the thermodynamic state of the system, which is characterized by temperature and pressure in the present constant-volume REMD simulations. The density of bulk water increases with temperature, due to a strong increase of pressure with temperature (Figure \ref{fig:reprhoshell}). In contrast, the density of hydration water decreases with temperature, and the thermal expansion coefficient of the hydration water, $\alpha_h$, is estimated to be -2$\cdot$10$^{-4}$ K$^{-1}$ up to $T \approx$ 450 K, and about one order of magnitude smaller at higher temperatures. The thermal expansivity of the hydration water observed in the present REMD simulations is notably smaller than that obtained near other peptides using constant temperature molecular dynamics simulations at \textit{p} = 0.1 MPa,\cite{Brovchenko:08b,Mitra:08a} which is probabily due to the very high pressures at high temperatures in the REMD simulations. However, the qualitative change in the temperature dependence of $\rho_h$ that takes place between 450 and 550 K cannot be attributed to the specific change of the thermodynamic state of system, as the pressure\textit{ p} increases almost linearly with temperature.
\par
The intrinsic volume $V_p$ of the peptides varies non-monotonically along the thermodynamic path studied: it slightly decreases upon heating up to $T \approx$ 350 K and decreases by about 30$\%$ upon further heating. The temperature dependence of ln$V_p$ is shown in Figure \ref{fig:repintvol}. The slope is equal to the peptide's intrinsic thermal expansion coefficient $\alpha_p$. The temperature dependence of  ln$V_p$ is almost linear between 375  and 575 K, with $\alpha_p$ = -2.3$\cdot$10$^{-3}$ K$^{-1}$ (solid line in Figure \ref{fig:repintvol}). The temperature dependence of the intrinsic density $\rho_p$ of the peptides is shown in the lower panel of Figure \ref{fig:repintvol}. The density of the peptides is about 1.75 g/cm$^3$ and increases with temperature and pressure upon the dissolution of the aggregate.
%
\begin{figure}[htb]
\centering
\includegraphics[width=7.0 cm]{repintvol}
\caption{Temperature dependencies of ln$V_p$ (upper panel) and the average density $\rho_p$ of the peptides (lower panel).  }
\label{fig:repintvol}
\end{figure}
%++++++++++++++++++++++++++
\section{Discussion}  
The simulation studies of an aqueous solution of the amyloidogenic peptides show that upon cooling, the solution undergoes a transition from a one-phase state, with peptides completely dissolved in water, to a two-phase state, where a solid-like peptide aggregate is in equilibrium with a saturated solution of the peptide. Various properties of the peptide system show a behavior similar to a conventional first order phase transition of fluids, which is strongly smeared out due to the finite size of the simulation box. The rounding of the transition is caused by the destabilization of the minor phase (peptide aggregate) in any simulated finite system due to thermodynamic reasons.\cite{Binder:03a} This is an intrinsic problem of any simulation study and it can be minimized by using larger systems. Additional rounding of the transition may originate from the exchange of configurations between the one-phase and the two-phase thermodynamic states in REMD simulation methods. 
\par
Various parameters were used to characterize the peptide aggregation and to locate the temperature of the demixing transition which separates these two distinct states of the system. One type of parameters reflects the properties of the whole system of peptides, e.g. the total radius of gyration or SASA, the total number of hydrophobic contacts or H-bonds, etc. Such an approach is valid, when the critical peptide concentration in the water-rich phase (the solubility limit) is so small that it is effectively reproduced by pure water, while all the peptides in the simulation box almost permanently belong to the organic-rich phase. This situation corresponds to the aqueous solution of amyloidogenic peptides at low temperatures. However, upon heating, the critical peptide concentration increases and the essential part of the peptides at this point belongs to the water-rich phase. Hence, the ability of parameters that utilize the properties of all  the peptides to describe the formation of a peptide aggregate degenerates as they become less sensitive to the formation of the peptide aggregate in the organic-rich phase. This results in the smearing out of the temperature dependences of these parameters in a wide temperature range, and a shift of the midpoint of the transition towards higher temperatures (see dependences in Figures \ref{fig:reprg} and \ref{fig:reprginvt}). Therefore, the inflection points of these dependences (450 to 485 K) give overestimated temperatures for the corresponding phase transition. In a similar analysis of the temperature-induced transition of amyloidogenic peptides in implicit water, the transition was found to be smeared out in an essentially narrower temperature interval, despite the use of a smaller number of peptides,\cite{Cecchini2004,Baumketner2005} which could be due to the lower peptide concentrations used in these studies.
\par
The appearance of the aggregated phase can be seen in the size distribution of the peptide clusters (Figure \ref{fig:repclustd}). The distortion of these distributions due to the small system size containing a few peptides only  prevents an accurate location of the transition temperature. The clustering analysis can be used to recognize the largest peptide aggregate that represents an organic-rich phase and offers the possibility to analyze its properties. The probability distribution of the size of the largest peptide cluster (Figure \ref{fig:repsmax}) shows  largest fluctuations at $\sim$ 440 K. The temperature dependence of the probability to have an aggregate including more than 2/3 of peptides, based on the connectivity criterion of  the distance between the centers of mass of two peptides ($r_{crit}\leq$ 1.1 nm), gives an inflection point at \textit{T} = 415 K. As this is one  of the most robust parameters for characterizing aggregation of this particular peptide, the midpoint of the demixing phase transition can be expected to be located at $\sim$ 415 K.
\par
Similar to the previously reported studies of the same peptide at \textit{T} = 330 K, the hydrophilicity of the total exposed surface increases upon aggregation (Figure \ref{fig:repsasa}). This conclusion is also supported by the increase in the density of hydration water upon aggregation (left panel in Figure \ref{fig:reprhoshell}). It is to be noted that this is not related to a decrease in temperature, as a corresponding decrease in pressure leads to an even more drastic decrease in the density of bulk water.
\par
The thermal expansivity of the peptides in an aggregated state is close to zero (Figure \ref{fig:repintvol}), indicating a solid-like structure of the aggregate. As the aggregate dissolves with increasing temperaure, the thermal expansivity $\alpha_p$ of the peptides becomes negative. This is expected for a random chain, and was also observed for the amyloidogenic peptide  A$\beta_{42}$ and an elastin-like peptide in water.\cite{Brovchenko:08b} In the latter case, a random chain behavior was confirmed by analyzing the structural properties of the peptide.\cite{Krukau2007} Therefore, a negative value of $\alpha_p$ obtained for the peptide system studied indicates a random chain conformational behavior of single peptides in the disaggregated state, which dominates at high temperatures.  
\par
The  density of the peptide aggregate (1.75 g/cm$^3$) largely exceeds the density (0.66 g/cm$^3$) estimated for the same peptide at 330 K in the previous Chapter. This could be due to  the fact that the latter value was obtained by estimating the $V_p$  from the radius of gyration $R_g$ of the largest peptide cluster and assuming its shape to be spherical. Since a sphere provides the maximal volume at a fixed radius of gyration, the  value estimated was the lower limit for the aggregate density. However, in this study, $V_p$ was calculated from Eq. \ref{eq:volp}, where the volume of the hydration shell was calculated using a slab approximation with a cutoff 0.35 nm, which could lead to systematic errors in absolute values of $V_p$. The dissolution of the peptide aggregate upon heating leads to the disappearance of the loosely packed structure with numerous voids, inaccessible to water molecules. Accordingly, the observed intrinsic density of the peptides, $\rho_p$, increases upon heating (Figure \ref{fig:repintvol}, lower panel). 
\par
$\beta$-sheet conformations are predominant in the peptide aggregate and the dissolution of the aggregate upon heating leads to the formation of $\alpha$-helical structures (upper panel in Figure \ref{fig:rephbond}). Some qualitative changes in this process, presumably related to the demixing transition, occur between 350 and 375 K. A drastic change of the correlation between the $\alpha$-helical content and the number of intra-peptide H-bonds occurs at about 475 K and 2.5 kbar. This indicates a change in the conformational behavior of individual peptides in liquid water, which occurs essentially above the demixing transition temperature. Approximately at the same state point, the dependence of the density of hydration water also changes qualitatively (see the left panel in Figure \ref{fig:reprhoshell}), whereas this is not the case for bulk liquid water. These changes point towards a possible relation between the state of the hydration water and the conformational behavior of the peptides.

