molecular dynamics simulations 
 have been used as a powerful tool 
 to study lipid membranes but are limited by the length and time scales 
 of the probing systems . in particular , the lateral diffusion rate 
 of a lipid in a membrane , although varying among different experimental 
 techniques , is in the range 1010 cm s. this means that it will take a lipid about 
 hundreds of nanoseconds to microseconds time scale to cover a 1 nm area .
furthermore , lipids are observed to move collectively 
 with their neighbors , and therefore , the rate at which a lipid swaps 
 position with its neighbor is even slower .
this 
 poses a serious problem when we simulate lipid bilayers with multiple 
 components , since the lateral organization of different components 
 requires extensive simulation to attain equilibrium and common accessible simulation time scales may provide 
 configurations ensemble biased toward the initial conditions .
temperature replica exchange molecular dynamics ( t - remd ) is one 
 of the methods that has been widely used to accelerate equilibration 
 in simulations and achieved numerous success in protein folding .
however , its application to lipid bilayers is rare since the number 
 of replicas scales with the degrees of freedom ( dofs ) of the whole 
 system , and lipid membranes usually have many more dofs than the protein 
 systems , making it computationally prohibitive to use t - remd to study 
 lipid bilayers .
one promising method to get around the poor 
 scalability of t - remd 
 with system size is replica exchange with solute tempering ( rest ) , 
 which was initially proposed by berne and co - workers .
rest is a specific variation of a generalized hamiltonian 
 replica exchange method . by changing 
 the solute  solute and solute  solvent interactions in 
 the system
although rest has been demonstrated to sample the 
 conformational ensemble of the alanine dipeptide successfully , its 
 efficiency in folding larger proteins remains unclear , which impedes 
 a wide adoption of the method . in this work , 

we show that rest can be used as an efficient way 
 to accelerate lateral equilibration in a mixed lipid bilayer .
we applied 
 constant pressure rest to a 1,2-dipalmitoyl - sn - glycero-3-phosphocholine 
 ( dppc ) bilayer with 50 mol % cholesterol ( chol ) .
constant pressure 
 in rest simulation is used because volume expansion at high temperature 
 increases lipids lateral diffusion . by carefully choosing the tempering 
 solute
, we managed to simulate the system with only 12 rest replicas , 
 which compares favorably to the requirement for 100 replicas 
 with t - remd at the same replica exchange rate .
the relative diffusion 
 rate between molecules in rest is , on average , an order of magnitude 
 faster than the rate in standard md simulation .
we also show that 
 the lateral radial distribution function between all molecular pair 
 types ( chol  chol , chol - dppc , dppc 
finally , we use 
 rest to obtain the gibbs free energy profiles between different molecular 
 pair types from the corresponding lateral radial distribution functions .
let us start with a brief review of replica exchange 
 with solute 
 tempering ( rest ) .
rest can be derived 
 from a more general form of replica exchange called hamiltonian replica 
 exchange . in hamiltonian replica exchange , 

replica m is simulated with potential energy em , at temperature tm ( the corresponding inverse temperature 
 will be referred as m , where m = 1/kbtm and kb is the boltzmann constant ) and constant pressure pm . in an isothermal 
isobaric 
 ensemble , the probability of configuration xm with volume vm in replica m is1where zm is the corresponding partition function .
the exchange between 
 replica m and n can be treated as 
 the change from state i to state f in the generalized ensemble,2and we use t(i  f ) to denote the transition probability 
 for i  f. applying the detailed 
 balance condition3gives the ratio of the transition 
 probabilities4where5 if the replicas are simulated at the 
 same temperature t0 and pressure p0 but a different potential energy , eq 5 can be further reduced to6 the 
 detailed balance condition guarantees the boltzmann sampling 
 of all replicas when sufficiently long simulations are performed . in rest , each replica is simulated at the same temperature ( t0 ) and pressure ( p0 ) but a different potential energy function .
we can divide the system 
 potential energy into three terms:7where each term , in order , represents solute  solute , 
 solute  solvent , solvent  solvent interactions . for replica m ,
we scale its potential energy according to8 following from eq 1 , the boltzmann distribution 
 for replica m becomes9equation 9 shows that , 
 thermodynamically , we can interpret the solute of the system as if 
 it is simulated at the scaled potential ( m/0)ess at temperature t0 , or at the original potential ess at an  effective temperature  tm . in each rest setup , we should 
 always ensure that there is one replica simulated with tm = t0
. we 
 will refer t0 also as the  target 
 temperature  because it is the temperature of interested for 
 the studied system .
replicas simulated with tm  t0 are 
 used for the sole purpose of enhancing sampling at t0 and configurations obtained from these additional replicas 
 do not represent any thermodynamic ensembles that have experimental 
 counterparts .
first , the 
 definition of 
 solute and solvent in a system is not absolute .
the solute can generally 
 present a part of the system whose sampling we want to accelerate , 
 while the solvent is defined as the rest of the system .
second , the 
 potential function of replica m can be of a form 
 different from eq 8 .
rest is just a specific 
 form of hamiltonian replica exchange , and there is no restriction 
 in the form of the potential energy function used in each replica 
 in hamiltonian replica exchange .
the advantage of using eq 8 is that there is a physical interpretation associated 
 with it , but this is not required .
therefore , the choice of the prefactor 
 in front of esw is not unique .
we choose 
 ( m/0 ) for the ease of implementation , as suggested by terakawa et al .
following eq 6 , 
 the exchange ratio between 
 replica m and n is determined by10it is clear from the above equation 
 that the 
 exchange ratio is independent of solvent  solvent interactions 
 ( eww ) .
therefore , one can enhance the 
 solute dynamics and simultaneously reduce the number of tempered degrees 
 of freedom , which results in a reduction of the number of required 
 replicas .
we want to point out that in our implementation of rest 
 we use eq 6 to calculate exchange rate rather 
 than eq 10 for practical reasons . specifically , 
 the total potential energy is easily available from the simulation 
 code . in the following text
, we describe how we scale the potential 
 energy 
 function according to eq 8 , using separate force 
 field parameters for each replica . in common molecular dynamics force 
 fields , the potential function consists of11where the first two terms are coulomb and 
 lennard - jones interactions and the third and the following terms are 
 bonded interactions which define the bond length potential , the bond 
 angle potential , the torsion potential , and so forth . the scaling 
 of the potential in replica m descried in eq 8 can be done as follows : ( 1 ) for the bonded interactions 
 in the solute , scale the spring constants by ( m/0 ) , ( 2 ) scale the charges in solute 
 by ( m/0 ) , and ( 3 ) scale ij by 
 m/0 if both i and j are in the solute , and by ( m/0 ) if only i or j is in the solute .
we scale ij directly because it can be applied to any 
 combination rule . by these
the system we studied 
 consisted of 144 dppcs , 144 cholesterols ( chol ) , and 14k water 
 molecules .
each monolayer in the bilayer was built independently by 
 randomly placing 72 dppcs and 72 chols on 12  12 planar grids . 

we define the z axis as the bilayer normal and refer 
 the values of z > 0 nm and z < 
 0 nm as the  upper  and  lower  monolayers . 
 then , we equilibrated the system for 20 ns at 323 k and 1 atm .
the 
 resulting configuration was used as the initial structure for all 
 replicas in rest and the standard md simulations .
we note here that 
 every replica in rest has the same starting configuration as in the 
 md simulation . in rest , we chose dppc molecules as the solute 
 and treated cholesterols and waters together as solvent .
the explanation 
 for this choice is provided in the results section . in total , 12 replicas were used and a 25% exchange rate 
 was achieved between the neighboring replicas .
each replica was simulated 
 at 323 k ( t0 ) , while the  effective 
 temperatures  of dppc were set at 323 , 341 , 360 , 380 , 400 , 
 421 , 445 , 471 , 500 , 531 , 564 , and 600 k ( tm ) .
each replica in rest was simulated for 60 ns , 
 while the md simulation was conducted for 400 ns . we implemented 
 hamiltonian replica exchange in gromacs 4.5.7 software package to conduct rest .
the default hamiltonian replica exchange in gromacs 
 is done through thermodynamic integration , and it suffers from a great 
 performance loss if the potential function involving a large portion 
 of the system is altered .
systems were 
 simulated under periodic boundary conditions , at constant temperature 
 and pressure . for temperature coupling , dppc and chol molecules
we note here that in rest , each replica is simulated at 323 k and 
 the  heating  of solute is done by reducing the solute 
rahman 
 barostat at 1 atm was used and the pressure 
 in the plane of the bilayer was coupled separately from the pressure 
 normal to the bilayer .
the temperature and pressure time constants 
 of coupling were 0.1 and 2.0 ps , respectively .
the spc / e model was used 
 for water and the 43a1-s3 force field was used for dppc and chol .
settle was 
 used to constrain water molecules and lincs was used to constrain all other bond lengths in the system .
the 
 sixth - order particle mesh ewald ( pme ) method was used for electrostatic interactions with a fourier spacing of 
 0.15 nm .
the real space coulomb interactions and pair - list calculations 
 were set to 1.0 nm .
a 1.0/1.6 nm twin - range cutoff scheme was used 
 for vdw interactions and the pair - list was updated every 10 steps . to study 
 how molecules diffuse relative to one another , we calculated the relative 
 diffusion coefficient dij between each pair of molecules i , j. the relative diffusion removes the contribution from collective 
 motions of molecules to the diffusion and should give a better estimate 
 of how fast the system samples various lateral configurations .
we 
 define the mean squared relative displacement between molecule i , j in a time interval t as1213where ri(t ) is the position of molecule i at time t and the bracket means an average 
 over different starting times t. then we linearly 
 fit rij2(t ) as a function 
 of t for larger t .
dij was assigned as /4 ( two - dimensions ) the slope of the curve .
we note here 
 that in rest , the relative diffusion rate is an average over temperatures 
 ( due to the exchange among replicas , replicas jump in the temperature 
 ladder ) ; however , it should still provide a meaningful description 
 of how fast the simulations sample bilayer lateral configurations 
 in general .
as mentioned in 
 the methods session , the choice of solute 
 for rest simulations is not absolute .
the solute can generally be 
 the part of the system whose dynamics we want to accelerate , while 
 the solvent is the rest of the system . since in this study
we want 
 to accelerate the dynamics of lipid bilayers , it is natural to choose 
 both the dppc and chol as the solute . before we conduct rest , it is 
 always a good practice to test the system at the highest temperature 
 that we want to simulate in rest .
we found that when we simulated 
 the dppc and chol at 600 k ( tm ) , the chol moved out of the monolayers and formed a third 
 layer sandwiched by the dppc bilayer .
a snapshot of the system is 
 shown in figure 1 . a possible explanation of 
 this can be obtained by carefully examining chol molecular structures . 

chol has a small hydrophilic alcohol headgroup and a bulky hydrophobic 
 body . at low temperature , the hydrophilic interactions between the 
 chol alcohol group and water favors aligning chol along the bilayer 
 normal .
when the effective temperature of chol increases , the entropic 
 effect becomes dominant and chol will gain entropically by placing 
 itself in the middle of the bilayer . on the other hand
, dppc has a 
 larger hydrophilic headgroup than chol ; thus , even at 600 k , it still 
 can anchor itself upright to form the bilayer .
snapshot of the system 
 where both dopc and chol molecules are tempered 
 at 600 k. water is colored blue . the oxygen atom in chol and the phosphate 
 atom in dppc
one way is to lower the highest temperature in rest to keep chol aligned 
 with the bilayer normal .
another way is to use only dppc as solute 
 and treat both chol and water as solvent .
we want to note that we 
 should be able to just conduct rest with both dppc and chol set to 
 600 k , even though chol moved out of the monolayers at this temperature . 
 however , this would be very inefficient .
configurations similar to 
 figure 1 obtained from high solute temperature 
 replicas will have vanishingly small probabilities in the target temperature 
 replica ensemble ( t0 = 323 k ) due to the 
 detailed balance condition ( eq 6 ) . as a result
neale et al . has observed such phenomena in another hamiltonian 
 replica exchange system . in this case , 
 high temperature replicas will not enhance the sampling in the target 
 temperature replica but consume computing time . in this work ,
we take 
 the second approach in which only dppc is chosen as the tempered solute . 

since this reduces the dofs in the solute , it further decreases the 
 number of replicas required to span our temperature range . by 
 choosing dppc as 
 solute
, we managed to heat dppc from 323 to 600 k with 12 replicas . 

the exchange rates between neighboring replicas are 25% , 22% , 24% , 
 26% , 27% , 22% , 23% , 22% , 22% , 24% , and 22% . in order to compare 
 the efficiency between rest and t - remd
, we estimate how many replicas 
 we need in our system to conduct t - remd to maintain a 25% exchange 
 rate between neighboring replicas .
figure 2 shows that if 12 replicas are used in t - remd , we 
 can only cover the temperature ranging from 323 to 340 k. this would 
 hardly accelerate the simulation as t - remd usually gains simulation 
 efficiency by increasing the enthalpy barrier crossing rate in high 
 temperature replicas .
comparatively , rest has a significant 
 advantage over t - remd as we can heat the dppc to a much higher temperature 
 with the same number of replicas .
number of replicas required to obtain 
 a  25% exchange rate 
 between neighboring replicas in rest and t - remd ( temperature replica 
 exchange ) .
below , we compare the sampling 
 efficiency between rest and a single 
 long md simulation . in our work , we ran rest for 60 ns and standard 
 md for 400 ns .
the 
 equilibrium sampling of lipid bilayers , especially of bilayers with 
 different components , depends on the ability to sample different lateral 
 organizations .
the faster the system can explore various lateral configurations , 
 the quicker equilibrium will be reached .
therefore , the lateral diffusion 
 coefficient plays a key role in determining the equilibrium rate . 

however , lipids usually move collectively in bilayers . this collective motion generally does not facilitate the 
 sampling of various lateral configurations but contributes significantly 
 to the lateral diffusion coefficient of each individual molecule . 
 therefore , we calculated the relative diffusion coefficient between 
 every pairs of molecule as it removes the contribution from collective 
 motion among lipids .
dppc ) of relative diffusion coefficients 
 are calculated from the md and rest simulations .
it clearly indicates that in all molecular pair types , molecules in 
 rest diffuse an order of magnitude faster than in standard md .
also , 
 we observe that the diffusion between chol and chol is the fastest 
 while the diffusion between dppc and dppc is the slowest .
dppc , on the other hand , has two 
 acyl tails that are usually entangled with other dppc , which reduces 
 the diffusion .
it is also worth noting that even though we only increase 
 the effective temperature of dppc , the diffusion of chol increases 
 as well . by omitting the chol from the tempered solute , we further 
 reduced dofs in the solute and therefore reduced the number of replicas 
 required for the system .
this also suggests that the efficiency of 
 rest may be further optimized by carefully choosing the tempering 
 solute , which is also pointed out in the original paper .
we note that the total simulation time in rest 
 ( 60  12 = 720 ns ) is almost twice as much as the time in md 
 ( 400 ns ) .
however , considering an order of magnitude increase in the 
 lipid diffusion , rest is quite efficient .
probability density of 
 all pairwise relative diffusion coefficients 
 between chol  chol ( c  c ) , dppc  chol ( d  c ) , 
 and dppc  dppc ( d  d ) in md and rest simulations .
the radial 
 distribution function ( rdf ) quantitatively describes the lateral organization 
 of molecules in a bilayer .
we define the lateral distance between 
 two molecules as the center of mass ( com ) distance between the molecules 
 projected in the x  y plane . 

as the coupling between separate monolayers is weak , the rdf calculated 
 from the upper and lower monolayers should converge to the same distribution . 

therefore , we can judge the convergence of a simulation by the rdf 
 difference calculated from separate monolayers .
figure 4 and supporting information 
 figures s1 and s2 show the rdfs among chol  chol , chol 
the rdfs 
 are calculated from separate monolayers ( blue , upper monolayer ; red , 
 lower monolayer ) from the rest and md simulations using different 
 block sizes . in all cases
chol rdfs ( figure 4 ) for 
 example , with a 10 ns block size , the rdfs in rest show reasonable 
 convergence between separate monolayers ( figure 4a  c ) . however , this is not the case for the md simulation 
 ( figure 4d  f ) .
the rdf difference in 
 md simulation in the last 10 ns block ( figure 4f ) is even larger than the difference in the first 10 ns block ( figure 4d ) .
we reason that as molecules diffuse slowly in 
 the md simulation , the sampled lateral configurations are highly correlated 
 in time ; thus , a larger block size is required to obtain independent 
 configurations .
therefore , we increased the time block size to 120 
 ns to calculate the rdf from the md simulation
. figure 4g  i show that rdf converges better when the block size 
 is increased .
the converged rdf in md ( figure 4i ) has a very similar shape with the rdf obtained from rest ( figure 4c ) .
the blue / red line in each subplot represents the rdf calculated from 
 the upper / lower monolayer , respectively . supporting information figures s1 and 

s2 show the rdfs between chol - dppc and dppc  dppc pairs , 
 respectively . 

we observe that the rdf of dppc  dppc does not converge as well 
 as the rdf of chol  chol or chol - dppc .
dppc pairs diffuse the slowest in both md 
 and rest ; therefore , we can expect that the rdf of dppc 
dppc 
 takes the longest time to converge . to quantitative analyze 
 the convergence rate
, we define the rdf 
 difference between separate monolayers ( rdfdiff ) as14where rdfu(i ) and rdfl(i ) are the value of the ith bin of rdfs 
 from the 
 upper and lower layer , respectively .
figure 5 shows the 
 convergence of rdf as a function of time block size in rest and md 
 simulations .
figure 5 shows that the 
 rdf in rest converges an order of magnitude faster than in md .
difference 
 between rdfs calculated from the upper and lower monolayers 
 as a function of block size .
error bars are estimated from consecutive 
 blocks with the same block size . in this session
, 
 we compare several bilayer structural properties calculated from the 
 rest and md simulations . for rest
, only the replica simulated with 
 solute temperature at tm = 323 k is used for analysis , as it represents the system with the 
 original hamiltonian .
normally , the area per lipid serves a 
 good indicator on bilayer structural properties .
however , since the 
 bilayer we studied has both chol and dppc , we calculated the average 
 area per molecule ( aapm ) instead .
aapm is defined as the projected 
 area of the bilayer in the x  y plane divided by the number of molecules in a monolayer .
the aapm 
 are 43.1  0.6  in rest and 42.8  0.4 
  in md .
another important bilayer 
 property is the deuterium order parameter 
 ( scd ) of the lipid acyl tails .
the order 
 parameter of a methylene at position i is defined 
 as15where 
 i is the angle between a c 
d 
 vector of the ith methylene in an acyl chain and 
 the normal of the bilayer ( z axis ) .
our calculation 
 suggests that rest and md simulations have similar bilayer properties , 
 which is a good validation of the rest method .
deuterium order parameters 
 |scd| of 
 dppc palmitoyl chains at 323 k calculated from the rest and md simulations .
chol is well - known for its condensing 
 effect on lipid bilayers 
 composed of lipids with saturated acyl tails .
it smooths the lipid liquid / gel phase transition to the lipid disordered / ordered 
 phase transition .
it is reported that 
 the average value of |scd| for dppc at 
 liquid ordered phase and disordered phase are 0.36 and 0.21 , respectively.supporting information figure 
 s3 shows the |scd| of dppc at different 
 solute temperatures in rest . at high solute temperatures ,
the dppc 
 acyl tails are more flexible and the |scd| indicates that the bilayer is in the liquid disordered phase .
the 
 range of sampled |scd| in rest suggests 
 that rest works well even when the system experiences the liquid disordered / ordered 
 phase transition .
based on 
 the high resolution rdf obtained from rest , we calculate the total 
 gibbs free energy g(r ) and excess gibbs free energy g(r ) profiles between chol  chol , 
 chol  dppc , and dppc  dppc .
g(r ) and g(r ) are defined as16and17g(r ) is 
 the lateral radial distribution function ( rdf ) , and r0 is the reference distance where we set g(r0 ) = 0 .
g(r ) is the contribution to the g(r ) due to the jacobian 
 or area effect in the two - dimensional space .
this means that neither chol or dppc 
 tends to aggregate at the 50% chol concentration .
chol gibbs free energy profile in dilute conditions 
 and found that g(r ) drops below zero in the range 1.0 < r < 
 1.5 nm .
dppc g(r ) is almost flat when r > 1.0 nm , suggesting a random distribution of dppc 
 at 
 large distance .
chol 
 g(r ) , indicating 
 preferential interacting locations for chol  chol pairs .
this 
 supports some phenomenological models , such as the supperlattice model and umbrella model , which suggest long - range ordering for chol . however , 
 the barriers between the free energy minimums are of the order of 
 kt scale , suggesting that the ordering of cholesterols is sensitive 
 to the temperature . as the derivative of g(r )
is the mean force , the g(r ) we obtained can be used 
 as a reference for various coarse grained models for this system .
total 
 free energy g(r ) ( solid line ) and the excess free energy g(r ) ( dashed line ) profiles 
 between different molecular types as a function of the lateral molecular 
 center of mass distance .
a general guideline is to choose the part of the system for which 
 we want to accelerate the dynamics as the solute .
for example , if we study the lipid  protein 
 interactions and are interested in the affinity of different lipid 
 components to the protein , we can temper the lipids to accelerate 
 their diffusion . on the other hand ,
if we are interested in how the 
 protein adapt its conformation to the bilayer environment , we can 
 temper the protein instead . with a good choice of solute
, rest 
 can accelerate the dynamics 
 of the system with fewer replicas ( compare to t - remd ) .
the sole purpose of the replicas with the scaled potentials ( tm  t0 ) is to sample configurations that are likely to occur 
 in the target temperature ensemble ( tm = t0 ) .
therefore , if 
 we scale the potential in such a way that the system samples configurations 
 with low probability to populate at the target temperature ensemble 
 ( configurations such as figure 1 , in our case ) , 
 rest will be less efficient .
however , no matter what the choice of 
 the solute is , we always have one replica in rest that is simulated 
 at the original unscaled potential . choosing the solute that makes 
 rest efficient may require an intuitive trial and error approach . in the past
usually , coarse grained systems have fewer degrees 
 of freedom than their atomic counterparts , which results in smoother 
 free energy landscape and faster lipid diffusion .
another method , 
 developed by tajkhorshid and co - workers , is a membrane mimetic model , 
 which separates the lipid headgroup from its hydrophobic tails .
this 
 method facilitates headgroup diffusion while maintaining a hydrophobic 
 core and has been used to study the insertion of peripheral proteins .
wang et al . also developed a method based on the accelerated molecular 
 dynamics ( amd ) method , which accelerates lipid diffusion by adding 
 a boost potential to the original system .
rest provides an efficient way of accelerating the equilibrium of 
 lipid bilayer systems while simulating at least one copy of an unperturbed 
 potential and maintaining atomistic details .
in this work , we applied replica exchange with solute tempering 
 ( rest ) to a cholesterol  dppc bilayer system . in rest ,
part 
 of the system is chosen as solute and the solute  solute and 
 solute  solvent interactions are scaled such that thermodynamically , 
 the solute is effectively sampling at a different temperature .
we 
 found that choosing both cholesterols and dppcs as solute is not efficient 
 because most of the cholesterols moved out of the monolayers to form 
 a third layer at  high temperature  .
therefore , we chose 
 to temper the dppc molecules only . since the number of replicas for 
 rest only scales with the degrees of freedom in the solute , we managed 
 to use 12 replicas to sample dppc at  temperature  ranging 
 from 323 to 600 k , which , otherwise , would require 100 replicas 
 in the traditional temperature replica exchange molecular dynamics 
 ( t - remd ) .
the relative diffusion coefficients between all molecular 
 pair types ( chol  chol , chol 
dppc ) 
 in rest are , on average , an order of magnitude larger than in standard 
 md simulation , indicating a better sampling of lateral structures 
 in rest .
dppc , and 
 dppc  dppc radial distribution function ( rdf ) between separate 
 monolayers in the rest and md simulations .
since the coupling between 
 monolayers is weak , the rdf should converge to the same distribution 
 from different monolayers .
our results show that the rdf converges 
 much faster in the rest than in the md simulation .
bilayer structural 
 properties such as average area per molecules and deuterium order 
 parameters are similar between the rest and md simulations .
finally , 
 we obtained the lateral free energy profile between different molecular 
 types from rest , which could be used as a reference to coarse - grained 
 models of the system .
chol lateral free energy has 
 several local minima and shows a long - range ordering , but the free 
 energy barriers between minima are on the kt scale , indicating the 
 ordering may be sensitive to the temperature . while we see a significant 
 advantage of using rest to accelerate lateral equilibrium in mixed 
 lipid bilayers , we believe that rest , or more generally , hamiltonian 
 replica exchange will have broader applications .
for example , the 
 preferential interactions between different lipids and membrane proteins 
 can be studied by tempering the lipids while leaving the protein and 
 solvent at the target temperature .
also , combined with umbrella sampling , 
 rest can be used to accelerate relaxation on degrees of freedom orthogonal 
 to the reaction coordinate , which is reported as a hurdle in free 
 energy calculations for lipid membranes .
 the 
 lateral heterogeneity of cellular membranes plays an important 
 role in many biological functions such as signaling and regulating 
 membrane proteins . 
 this heterogeneity can result from preferential 
 interactions between membrane components or interactions with membrane 
 proteins . 
 one major difficulty in molecular dynamics simulations aimed 
 at studying the membrane heterogeneity is that lipids diffuse slowly 
 and collectively in bilayers , and therefore , it is difficult to reach 
 equilibrium in lateral organization in bilayer mixtures . here 
 , we 
 propose the use of the replica exchange with solute tempering ( rest ) 
 approach to accelerate lateral relaxation in heterogeneous bilayers . 
 
 rest is based on the replica exchange method but tempers only the 
 solute , leaving the temperature of the solvent fixed . 
 since the number 
 of replicas in rest scales approximately only with the degrees of 
 freedom in the solute , rest enables us to enhance the configuration 
 sampling of lipid bilayers with fewer replicas , in comparison with 
 the temperature replica exchange molecular dynamics simulation ( t - remd ) 
 where the number of replicas scales with the degrees of freedom of 
 the entire system . 
 we apply the rest method to a cholesterol and 1,2-dipalmitoyl - sn - glycero-3-phosphocholine ( dppc ) bilayer mixture and find 
 that the lateral distribution functions of all molecular pair types 
 converge much faster than in the standard md simulation . 
 the relative 
 diffusion rate between molecules in rest is , on average , an order 
 of magnitude faster than in the standard md simulation . 
 although rest 
 was initially proposed to study protein folding and its efficiency 
 in protein folding is still under debate , we find a unique application 
 of rest to accelerate lateral equilibration in mixed lipid membranes 
 and suggest a promising way to probe membrane lateral heterogeneity 
 through molecular dynamics simulation . 