in order 
 to characterize the radial undulations from a simulated vesicle , we 
 start by recasting the positions of the bilayer beads into an equi - angular 
 discrete surface representation using a  , -grid ( colatitude and longitude respectively ) , where   [ 0 ,  ] and   [ 0 , 
 2 ] , with  = 0 at the northern pole . 

the angular resolution is defined as1 where s is the arc - length 
 where the bilayer fluctuations transition between a continuum mode 
 to molecular ( s  4 nm ) , r0 is the average vesicle radius , and ns is the number of sampling points per s with ns  2 to satisfy nyquist 
 sampling theorem . to render the surface
, we define the origin 
 at the com of the vesicle ( lipids + membrane inclusions ) and then 
 transform the system into spherical coordinates ( xi , yi , zi  ri,i , i , where ri,i , i corresponds 
 to radius , colatitude , and longitude for each bead i , respectively )
. the direct use of an equi - angular  , -grid to bin bead positions introduces 
 discontinuities at the poles .
we mitigate these discontinuities by 
 implementing an arc - length low - pass filter with filter cutoff qarc = 2/s .
our previous 
 work with undulation analysis on flat patch systems identified a cutoff 
 wavenumber , q0 = 1.5 nm , at which the long - wavelength undulations transition into molecular 
 structure fluctuations .
we use this q0 as a guide in defining an appropriate qarc for the vesicle systems .
figure 1a shows a snapshot from an equilibrated , 
 34 nm dmpc vesicle , and figure 1b presents 
 a cutaway view , emphasizing the vesicle s diameter relative 
 to its thickness .
we parse all lipids into the inner and outer monolayer 
 using the average vector between the acyl - chain terminal beads and 
 the phosphate bead , where the radial component defines which monolayer 
 the lipid belongs to ( i.e. , r < 0 , inner monolayer ; r > 0 , outer monolayer ) . using an arc - length filter with qarc = 2.5 nm ,
we define 
 the local  , -membrane surface for both 
 inner and outer monolayers , rin( ,  ) and rout( ,  ) , respectively .
figure 1c schematizes the first stages of our method , with the inner 
 ( red ) and outer ( blue ) monolayer beads shown for numerous  , -positons . on the left hemisphere of figure 1c , regions of rin( ,  ) and rout( ,  ) are presented for illustration .
the average 
 of the two surfaces define the undulating radial surface as2 a. 
 snapshot of a 34 nm dmpc - lipid vesicle .
c. an arc - length filter is used to define a surface for both inner 
 and outer monolayers ( blue and red selection beads respectively ) . 

d. rund( ,  ) is determined as the average of the inner and outer monolayers 
 ( color - map indicates fluctuations about the average radius , r0 , ( units in nm ) .
figure 1d shows rund(, ) calculated from 
 one frame of the 
 dmpc trajectory .
the color - map displays the deviations about the average 
 vesicle radius , r0 ( units in nm ) , which is typically smaller 
 than the ideal sphere radius r0 . from rund( ,  )
, we then 
 define the normalized radial fluctuations as3 next , spectral 
 decomposition of rund( , 

spherical harmonics are standing waves on a sphere . given 
 that rund( ,  ) is a discrete surface defined on a spherical manifold
, we can represent 
 it as a linear combination of spherical harmonics with degree ( l ) and order ( m ) , corresponding to the 
 number of waves in the  , -dimensions , 
 respectively .
spha decomposes rund( , 
  ) into a series of standing waves on a sphere with 
 degree l and order m. helfrich s formulism for spha expands the 
 normalized radial 
 fluctuations in spherical harmonics as4where alm are the spherical harmonic coefficients 
 of degree l and order m ( l  
 0,1,lmax and m  l, ... ,l ) with lmax defined by the number of 
 colatitude grid points.ylm are the spherical harmonic basis functions5where plm( ) are the fully normalized 
 associated legendre polynomials with the normalization factor , nlm , defined as6 because f( ,  ) is a real - valued function alm = al , m * , the complex conjugate of alm .
we thereby redefine alm such that7allowing us to write8where now m  0, ... ,l .
we generate the 
 matrix p with dimensions 
 ( 2n * lmax )  n for a given  , -distribution , where n is the number of colatitude parallels , such that for each li  0, ... ,lmax , the corresponding 2 m + 1 rows are 
 defined as9from which the transformation 
 matrix y follows as10 from the matrix y we can write 
 the 
 spherical harmonic forward transform as11where alm1 is a recasting 
 of the spherical harmonic coefficients alm with dimension n  
 1 corresponding to the row construction of y , and f,1 is a matrix of positions with dimension 
 ( 2n * lmax )  1 . whereas the forward transform is exact , the inverse transform ( i.e. , 
 spha ) is overdetermined . using the factorize package in matlab ,
we compute the pseudoinverse of y using qr - factorization and apply a least - squares 
 approximation to determine alm1 as12 implementation of the spha 
 algorithm results in error propagation 
 from numerous sources ( i.e. , truncation error , round - off error , and 
 least - squares approximation ) .
we evaluated the magnitude of these 
 effects by calculating the root - mean - squared - difference ( rmsd ) between rund( ,  ) and the 
 transformed rund( , 
  ) , where13 rmsd 
 calculated for the dmpc system is 3.5  10 nm and for the dmpc + cholesterol system is 1.8  10 nm .
these errors are 2 orders of magnitude lower than the magnitude 
 of the radial fluctuations , o ( 0.5  10 nm ) , 
 allaying concerns regarding the inherent propagating errors .
membrane structure 
 ( (z ) and al ) for 
 flat - patch membrane systems was determined using our previously published 
 method . briefly
the urs was determined 
 using the phosphate atoms and subsequently filtered with a 1.5 nm cutoff .
the filtered surface was then used to reference 
 every atom to the local bilayer midplane .
these distances were subsequently 
 binned and volume normalized to extract the number density profile . the standard method for determining 
 the radial membrane number density profile , rbin(r ) , references every bead / atom in the system relative 
 to the com of the vesicle in spherical coordinates , averaging across 
 all  , -angles .
we have previously shown 
 that the membrane structure profile is smoothed out in systems where 
 long - wavelength undulations develop . by 
 isolating the long wavelength undulations and referencing every bead / atom 
 to the local undulating reference surface
following the same principles for our method developed 
 for flat - patch systems , we apply a low - pass ideal filter to the spherical 
 harmonic coefficients with an order cutoff lcut = qarcr0  
 0.5 such that14where lm1 are the filtered coefficients . for the dmpc 
 vesicle lcut = 26 , whereas for the slightly 
 larger dmpc + cholesterol 
 vesicle , lcut = 29 .
we next apply an the 
 inverse spherical harmonic transform to resolve a filtered radial - undulation 
 surface rund( ,  ) 
 as15 with 
 the filtered surface we resolve the vesicle s local 
 membrane structure , uc(r ) , by referencing 
 every bead / atom i relative to it is position on the 
 surface , such that16where rund( , )nearest corresponds 
 to the closest equal - angular grid point on rund( ,  ) .
the grid definition samples the corresponding 
 filter wavelength at least 2 times ( satisfying nyquist sampling theorem ) . 
 as a default
we employ 4-samples per wavelength ( at qarc = 2.5 nm the arc - length resolution 
 is 0.4 nm ) , providing a close approximation to the unique rund ,  i( ,  ) .
the resulting ri( ,  ) 
 is then binned with a bin - width of dr = 0.1 nm and 
 subsequently normalized by the differential volume shell17with r0the 
 radius of the ideal sphere with equal volume as rund( ,  ) defined as18and adjusted for every bin index b by the corresponding distance from the surface .
in addition to number density profiles , 
 we can identify the vesicle s 
 area per lipid , al , both as a whole and 
 independently for both monolayers . instead of defining r0 with 
 all vesicle beads / atoms
, we parse them independently as inner and 
 outer monolayers . then using eq 18
, we obtain 
 two additional ideal sphere radii ( three in total ) , r0,vesicle , r0,inner , and r0,outer .
the corresponding al for 
 each ideal radius is simply192021where nl ,  inner and nl ,  outer refer to the number 
 of lipids in the respective monolayer . from alm1 , we obtain the undulation power 
 spectrum by binning the modulus of the spherical harmonic coefficients 
 across degree l. the resulting profile can be interpreted 
 according to the helfrich continuum model for undulations on a sphere 
 with vanishing spontaneous curvature as22where t is temperature and l  2, ... ,lmax . in the flat - patch spectral method , the 
 undulations power spectrum ,
|u(q)| , determined by direct fourier transformation 
 of the bilayer selection atoms can be modeled as23where n is the number of 
 selection atoms , a is the projected area per lipid , 
 and s(q ) is the 
 in - plane structure factor .
s(q ) is subtracted from n|u(q)| to provide a broader range 
 of modes within the q regime .
kc is then determined by fitting the low - q , long - wavelength modes . using a filter to 
 define the urs attenuates the spectral intensity of the undulations . 

this attenuation can span a broad frequency bandwidth , bleeding through 
 into the desired frequencies of the signal .
we explored a range of 
 arc - length filter cutoffs to characterize the filter s frequency 
 response , with the goal of identifying a cutoff where signal attenuation 
 is limited to frequencies above the crossover wavenumber , q0 = 1.5 nm .
figure 3 presents the power spectra for the dmpc vesicle 
 system using four different arc - length filters with wavenumber cutoff , qcut , ranging from 0.5 nm to 2.5 nm .
undulation spectra for a range of arc - length 
 cutoff wavenumbers , 
 ( qcut = 0.5 nm : blue ; 
 1.5 nm : green ; 2.0 nm : red ; 
 and 2.5 nm : cyan ) for the dmpc vesicle system . 

the frequency response of the arc - length filter is far from 
 ideal , 
 with significant bleed - through extending below the desired cutoff 
 wavenumber .
this is most noticeable when comparing the crossover wavenumber q0 = 1.5 nm ( black dashed - line)the 
 wavenumber where continuum undulations transition into molecular fluctuations  to the 1.5 nm filter 
 cutoff ( green triangles ) . there is significant loss of undulation 
 intensity for degrees 1225 , all before the 1.5 nm crossover wavenumber .
this loss of intensity skews the kc - fit , resulting in a larger kc .
furthermore , increasing the filter cutoff increases the number 
 of degrees that comprise the linear region below q0 , thereby 
 improving the fit to eq 22 . bending
with qcut = 2.5 nm , all 26 degrees ( spanning the full range , q0  1.5 nm ) can be used to determine kc .
all 
 simulations were 
 run using the martini force field and using the gromacs v4.5.3 program .
isobaric 
 ( npt ) ensemble at constant pressure and temperature ( 1 bar and 303 
 k , respectively ) using either a 25 fs time step for pure lipid or 
 10 fs time step for cholesterol systems .
pressure coupling was isotropic 
 for the vesicle systems and semi - isotropic for flat - patch , with independent xy- and z - barostats , resulting in a tensionless 
 bilayer .
initial equilibration for vesicle systems included 100  000 
 steps of steepest descent minimization followed by 500 ns of dynamics 
 using the velocity - rescaling thermostat and a berendson barostat ( flat - patch 
 systems underwent 10  000 steps steepest descent minimization 
 and 100 ns velocity - rescaling dynamics ) .
hoover 
 thermostat and the parrinello  rahman barostat with a time constant 
 of 2.5 and 250 ps , respectively .
pressure coupling 
 was applied isotropically for vesicle simulations and semi - isotropically 
 for flat - patch systems .
vesicle production simulations were 2.5 s 
 ( 10 s scaled simulation time ) , flat - patch simulations were 
 5 s ( 20 s scaled time ) sampled every 1 ns .
the 
 starting configurations for the dmpc vesicle ( 11  126-dmpc , 
 1  123  315-coarse grained water + antifreeze ,  high - water  , 
 1  113  742-cg water + antifreeze  low - water  ) 
 and dmpc + cholesterol ( 12  771-dmpc + 5473-cholesterol , 1  313  698-cg 
 water + antifreeze ) were constructed by randomly seeding two opposing 
 monolayers with the appropriate number of lipids ( dmpc and/or cholesterol ) 
 based on the area per lipid of each species , the surface area of the 
 monolayer shell , and the mole fraction of the mixture . a second dmpc 
 vesicle was simulated using the 3 s frame as a starting configuration 
 with 10% of the internal water beads removed to explore effects of 
 system construction on structure and bending rigidity .
flat - patch 
 bilayer systems were constructed with 3200-dmpc lipids and 2240-dmpc 
 + 960-cholesterol , randomly seeded in flat monolayers with 70,400 
 coarse grained water + antifreeze particles . during the early 
 stages of equilibration , rapid pore formation 
 occurred throughout the vesicles .
these pores coalesced and closed 
 quickly  within 70 ns  providing an opportunity for both 
 lipid flip - flop and water exchange across the vesicle to equilibrate 
 the system .
for the low - water system , the vesicle collapsed into an 
 ellipsoidal shape within 100 ns , and maintained an ellipsoidal character 
 throughout the 3 s production run .
trajectories were manipulated and processed 
 using both the gromacs v.4.5.3 simulation package and the mdanalysis python library .
further 
 data analysis and figure rendering was performed using matlab ( v.r2012a ) 
 with use of the factorize library .
having established the framework 
 to define rund( ,  ) for vesicle simulations , 
 it now becomes possible to extract the underlying membrane structure 
 from the fluctuating vesicle s trajectory .
figure 4 presents the number density profiles for the dmpc 
 vesicle , rbin(r ) and uc(r ) , and the dmpc flat - patch , uc(z ) , systems .
vesicle profiles are illustrated 
 as dashed lines , flat - patch profiles as filled distributions , with 
 the component groups color - coded as described in the caption .
negative 
 radial positions correspond to the inner monolayer , positive values 
 correspond to the outer monolayer .
comparisons of component number densities 
 for dmpc vesicle and 
 flat - patch systems .
flat - patch 
 profiles are filled - distributions with headgroup ( black / gray ) , carbonyl - glycerol 
 ( red / pink ) , and acyl - chain ( blue / cyan ) .
( a ) , ( b ) 
 rbin for the high - water and low - water dmpc system , 
 respectively .
( c ) , ( d ) uc for high - water and low - water 
 dmpc system . the same flat - patch profile has been included in all 
 panels .
figure 4a compares the high - water vesicle 
 profile , rbin(r ) , to the flat - patch 
 profile , uc(z ) .
the uc(z ) profile was translated by the average 
 vesicle radius , r0 , to allow direct comparison to the 
 vesicle profiles . as expected ,
using a global reference frame with 
 a fluctuating bilayer results in a broadening of the component distributions 
 ( rmsd between rbin(r ) and uc(z ) is 3.08 nm ) .
this 
 broadening results in a loss of structural resolution in the vesicle 
 profile , most noticeably in the acyl - chain distributions , where the 
 typically pronounced terminal methyl trough is absent .
figure 4b compares rbin(r ) from the ellipsoidal , low - water dmpc system to the flat - patch profile . 
 as fluctuations increase ,
the low - water rbin(r ) displays 
 significant distortion ( rmsd is 7.12 nm ) across 
 all lipid - component distributions , highlighting the problem with using 
 a global com reference frame to define membrane structure .
figure 4c , d present comparisons of our new 
 method , uc(r ) , to the flat - patch 
 profile for both the high- and low - water dmpc systems , respectively . 
 in both cases ,
the resulting number density more closely resembles 
 the flat - patch result ( rmsd between uc(r ) and uc(z ) is 1.66 nm for low - water dmpc and 2.34 nm for high - water ) . 

each component distribution is tightened up relative to rbin(r ) , and key structural features ( e.g. , 
 the terminal methyl trough ) are properly characterized . in the high - water 
 system ( figure 4c ) , the headgroup distributions 
 are centered at the appropriate position ( i.e. , the bilayer thickness 
 matches that of the flat - patch ) . however , the amplitude and width 
 of the headgroup distributions do not fully agree with the flat - patch 
 result : the amplitudes are too high and the widths to narrow .
this 
 small discrepancy suggests that the lipids are under tension , either 
 due to water - density or lipid density asymmetry across the bilayer .
comparison of the low - water profiles ( figure 4b , d ) demonstrates the dramatic improvement in calculated membrane 
 structure when undulations are isolated and characterized during structure 
 determination .
even with a dramatic ellipsoidal geometry , the uc(r ) method extracts the underlying membrane 
 structure .
the low - water profile amplitudes agree with the flat - patch 
 profile , more so than the high - water result .
however , the membrane 
 appears thicker and the component distributions are slightly broader 
 than the flat - patch result . comparing figure 4c with 4d ( uc(r ) for the spherical 
 and ellipsoidal dmpc vesicle systems )
as expected , varying the number 
 of waters inside the vesicle can lead to different membrane tensions 
 and corresponding structures . with our new algorithm ,
this type of 
 comparison can serve as a guide in the initial construction of vesicle 
 simulations , evaluating the membrane structure for imbalances in water 
 density or lipid density across the bilayer .
analysis of the 
 dmpc + cholesterol system highlights additional 
 complexities that can develop due to the different system geometries .
supplemental figure 1 presents rbin(r ) and uc(r ) 
 for the dmpc + cholesterol system .
the undulation correction still 
 sharpens the component distributions , albeit to a lesser extent than 
 the pure dmpc system .
cholesterol s facile ability to flip - flop 
 during the time - scale of the simulation results in an asymmetric cholesterol 
 distribution across the bilayer in the vesicle geometry .
the development 
 of an asymmetric cholesterol distribution is not surprising for the 
 vesicle as the spontaneous curvature difference between inner and 
 outer monolayers will induce an asymmetric partitioning of lipids 
 with non - neutral spontaneous curvature .
close inspection of all vesicle number density profiles highlight 
 an asymmetry ( particularly in the headgroup distributions ) where the 
 inner monolayer has a higher density relative to the outer monolayer 
 ( figure 4 and supplemental 
 figure 1 ) .
we can not know for certain whether this asymmetry 
 is due to a redistribution of lipid components or an artifact from 
 a poorly equilibrated vesicle .
vesicle structure does not change 
 across the last 1.5 s 
 for both the dmpc and dmpc + cholesterol systems .
however , the differences 
 in structure between high - water and low - water dmpc vesicle systems 
 highlight the effect of water density imbalances across the bilayer . 

furthermore , in all simulations , lipid flip - flop was only observed 
 for cholesterol and never for dmpc .
significantly longer simulation 
 times ( > 100s of s ) of standard md would be necessary to 
 allow 
 for both water equilibration and sufficient lipid flip - flop events 
 to fully equilibrate the vesicle .
our focus here has been on 
 the development of a method to calculate 
 structure from vesicle simulations and not to perfectly refine these 
 particular simulations .
nevertheless , the problem of vesicle construction 
 and equilibration is an important challenge that must be answered 
 to facilitate further study of these more complex systems .
risselada 
 et al . employed an artificial pore allow for both water and lipid 
 exchange across the vesicle prior to the production md simulation . for simple lipid compositions ,
this approach 
 is preferable over an iterative refinement of the lipid and water 
 distributions ( altering them by hand on the basis of structural profile 
 imbalance ) and far more computationally efficient than extending the 
 simulation until sufficient lipid flip - flop is observed .
figure 5 presents 
 the area per lipid ( al ) time - series for 
 the dmpc vesicle systems .
we compare 
 values determined along the urs for the full vesicle ( all , blue ) and 
 both inner and outer monolayers ( green and red , respectively ) with 
 those of the flat - patch al ( cyan ) .
a summary 
 of the average al is presented in table 1 . in the high - water dmpc system ( figure 5a )
there are three distinct differences between 
 the vesicle and flat - patch result : ( 1 ) the inner and outer monolayer 
 have different al , with inner monolayer 
 lipids occupying more lateral area than those in the outer leaflet ; 
 ( 2 ) all measures of al are greater than 
 the flat - patch al ; and ( 3 ) the fluctuations 
 in al are significantly decreased in the 
 high - water vesicle system versus the flat - patch .
the imbalance between 
 inner and outer monolayers correlates with the structural asymmetry 
 in the headgroup uc(r ) .
( a ) al trajectory for the high - water 
 dmpc vesicle system with total vesicle ( blue ) , inner monolayer ( green ) , 
 outer monolayer ( red ) , and flat - patch system ( cyan ) .
we presume that this reflects that 
 the inner monolayer is under 
 greater tension than the outer , constraining the inner headgroups 
 in the radial dimension and resulting in a tighter position distribution 
 ( recall that the thickness is the same , figure 3b ) .
there are two potential causes of this tension , either an imbalance 
 in the water density or the lipid density across the vesicle .
the 
 area per unit cell , auc = 2 asurface / ndmpc , for the dmpc 
 + cholesterol system shows a similar trend as the dmpc vesicle ( see supplemental figure 2 ) .
figure 5b illustrates the same al analysis for the ellipsoidal , low - water dmpc system . 

it is immediately evident that the water content dramatically affects 
 the vesicle s structure and dynamics . the ideal sphere al is increased for all three al metrics , even though the average vesicle radius is smaller . 

there is a partial recovery in the al symmetry 
 across the monolayers ( i.e. , a reduction in al = al , out 
the residual discrepancy between al , out and al , in likely stems from 
 a lipid density imbalance across monolayers , which is still evident 
 in figure 4d .
although the structure from the 
 low - water system more closely matches the flat - patch profile , it is 
 the al of the high - water system that is 
 in better agreement .
the number density profile calculation isolates 
 the increased fluctuations in the low - water system , whereas the al does not .
figure 6 presents 
 the power spectra and subsequent model fits for dmpc high - water ( black ) 
 and dmpc + cholesterol ( red ) in both vesicle ( figure 6a - c ) and flat - patch ( figure 6d ) geometries . 
 in each case
the spectra was determined over the final 1.5 s 
 of the trajectory . as expected , in both
the flat - patch and vesicle 
 simulations cholesterol reduces the magnitude of the undulation intensity , 
 reflecting the sterol s well - know rigidifying effect.kc values were determined by fitting 
 the linear regime across low degrees corresponding to q0  1.5 nm to the appropriate helfrich 
 continuum model for vesicle or flat - patch geometries ( l  26 for dmpc and l  29 for the 
 slightly larger dmpc + cholesterol vesicle ) . as shown in table 2 , we observe a similar cholesterol - induced increase 
 in kc in the vesicle- and flat - patch geometry 
 ( the latter of which was obtained via our the previously published 
 method ) .
power spectra for both vesicle systems 
 ( a ) with corresponding kc fit for dmpc 
 ( b ) and dmpc + cholesterol ( c ) . 
 undulation power spectrum and kc - fits 
 for flat - patch systems . for all panels , dmpc ( black ) and dmpc + cholesterol 
 ( red ) .
the * denotes values from brandt 
 et al . ; experimental values ( * * ) from 
 pan et al .. the high - water ( hw ) dmpc vesicle kc is in better agreement with our previously published kc from a large flat patch , where a 30  000 
 lipid 
 dmpc coarse - grained flat - patch system had kc = 1.5  10 j. to ensure the use of an equi - angular grid did not introduce any 
 artifacts into the vesicle analysis , we repeated the spha on the high - water 
 dmpc vesicle in a 90 rotated reference frame .
the results were 
 identical to the nonrotated system ( see supplemental 
 figure 3 ) .
the increased kc from our smaller 3200 
 lipid dmpc flat - patch highlights potential system - size effects with 
 spectral methods to extract kc .
the low - water 
 ( lw ) system shows a reduced kc , relative 
 to the high - water kc , agreeing with the 
 increased al fluctuations ( see figure 5b , table 2 , and supplemental figure 4 ) .
as we discuss below , 
 other factors may also contribute to the differences in kc across the two system geometries , including hydration 
 levels .
for example , in the low - water system , we observed a dramatic 
 increase in fluctuations and concomitant decrease in kc .
the increase in rigidity for the dmpc + cholesterol 
 system , determined 
 as the ratio of kc between dmpc + cholesterol 
 and dmpc , shows excellent agreement for the high - water ( hw ) case , 
 1.8 for flat - patch versus 1.9 for vesicles .
for both cases ( pure lipid 
 and + cholesterol ) , the vesicle s kc is smaller than the corresponding flat - patch ( 1.3  10 j versus 2.1  10 j for 
 dmpc and 2.4  10 j versus 3.8  10 j for dmpc + cholesterol ) . in both the flat - patch 
 and vesicle systems ,
the calculated kc value for dmpc is an order of magnitude greater 
 than that determined experimentally ( 0.58  10 j , determined from diffuse x - ray scattering on oriented bilayer 
 stacks ) .
however , given that our goal is 
 to establish a reliable and robust ( i.e. , force - field independent ) 
 algorithm , we consider the agreement between the established flat - patch 
 method and our new vesicle method as the relevant measure of success , 
 rather than close agreement between our result and experiment .
that 
 said , we do note that there is much better agreement between the kc calculated from simulation and the experimental 
 values for both dmpc low - water ( 6.1  
 10 j versus 5.8  10 j ) and the dmpc + cholesterol systems ( 2.4  10 j versus 2.73  10 j , respectively ) . achieving 
 better agreement with experiment in these values remains a challenge 
 for force - field development , though it remains to be seen whether 
 any force field ( even fully atomistic ones ) will be able to accurately 
 capture kc across all lipid types .
because of the sensitivity of the spectrum to the arc - length filter 
 cutoff ( see figure 3 ) , it was important to 
 evaluate the sensitivity of the kc fit 
 to the loss of undulation intensity .
we did so by varying the number 
 of degrees used in the data - fitting , testing both the 1.5 nm and 2.5 nm filters .
supplemental 
 figure 5 illustrates the sensitivity of kc fits and time course of kc for 
 both dmpc and dmpc + cholesterol vesicle systems .
the increased linear 
 region for the 2.5 nm spectra results in a much 
 broader span of converged kc fits for 
 both vesicle systems ( supplemental figure 5b ) .
this added robustness confirms our choice of filter parameters , 
 specifically the 2.5 nm arc - length filter ( for 
 rendering the initial surface ) and the 1.5 nm ideal 
 filter ( to define the urs and truncate the sph coefficients for the kc - fit at the transition from continuum to molecular 
 length scales ) .
the evolution of kc shows 
 very small changes across the duration of the simulations for both 
 dmpc and dmpc + cholesterol systems ( supplemental 
 figure 5c ) .
we have developed 
 a method to determine membrane structure ( number 
 density and area per lipid ) and bending rigidity from md simulations 
 of lipid vesicles .
comparisons of number density profiles determined 
 from vesicle simulations using either a global or a local reference 
 frame highlight the broadening effect introduced by local bilayer 
 fluctuations , similar to what was previously observed in large flat - patch 
 simulations . using a urs
removes 
 the fluctuation - induced broadening effect and 
 recovers the underlying structure profile , even under extreme vesicle 
 deformations ( e.g. , elliptical vesicles ) .
our new spha method can 
 be applied to any md simulation where the topography is such that 
 radial lines intersect the topography only once ( e.g. , vesicles , bubbles , 
 micelles ) . in systems where the topography is more tortuous ( e.g. , 
 where radial lines intersect the topography multiple times ) the spha 
 method breaks down .
nevertheless , the local structure profile can 
 still be resolved via the urs , as long as the method used to define 
 the urs is modified to accommodate the more complex topography .
our choice of using an equi - angular grid in order to define the 
 urs as well as in implementing the spha was done to simplify the calculation 
 of the sph coefficients via matrix transformation ( eq 12 ) .
alternative surface definitions ( e.g. , geodesic gridding ) 
 are equally viable ; however , they may introduce increased complexity 
 and computational cost in determining the transformation matrix .
any 
 predefined  ,  surface reference frame 
 requires either a real - space filter ( e.g. , arc - length filter ) or a 
 direct spherical harmonic transformation that requires treating every 
 membrane bead as a radial delta function to obtain the urs .
although 
 the latter is the direct corollary to our flat - patch direct fourier 
 method , the increased computational cost 
 of the explicit sph transform makes that approach prohibitive .
the spherical harmonics analysis of the urs allows us to determine 
 the al for the equivalent ideal sphere 
 that samples the same area as the undulating bilayer .
we define a 
 unique al for the vesicle as well as both 
 monolayers .
these three al metrics describe 
 structural imbalances that exist due to initial vesicle construction 
 ( e.g. , water or lipid imbalances across the bilayer ) .
the convergence 
 of these three metrics may serve as a simple readout for vesicle equilibration .
this paper lays the groundwork for an iterative process to improve 
 vesicle construction . with flat - patch membranes ,
numerous studies 
 have been successful in matching experimental structure profiles by 
 tuning the al through altering the periodic 
 box dimensions .
the vesicle membrane geometry raises a new and more difficult challenge : 
 obtaining al , structure , and bending rigidity 
 correct by varying the lipid distribution in the two monolayers and 
 the water density inside the vesicle .
our algorithm can extract 
 these structural and mechanical properties 
 to guide the construction and equilibration of these complex vesicle 
 systems .
comparison of these structural and mechanical results with 
 experimental measurables is the ultimate goal .
although we are currently 
 limited by the accuracy of the simulation force fields , changes in 
 structure that correlate with changes in rigidity may provide insight 
 into the interpretation of experimental results .
 we have developed an algorithm to 
 determine membrane structure , 
 area per lipid , and bending rigidity from molecular dynamics simulations 
 of lipid vesicles . 
 current methods to extract structure from vesicle 
 simulations define densities relative to the global center of mass 
 of the vesicle . 
 this approach ignores the long - wavelength fluctuations 
 ( undulations ) that develop across the sphere and broaden the underlying 
 structure . 
 our method establishes a local reference frame by defining 
 a radially undulating reference surface ( urs ) and thereby removes 
 the broadening effect of the undulations . using an arc - length low - pass 
 filter 
 , we render the urs by defining the bilayer midplane on an equi - angular  , -grid ( colatitude , longitude ) . 
 this surface 
 is then expanded onto a truncated series of spherical harmonics . 
 the 
 spherical harmonic coefficients characterize the long - wavelength fluctuations 
 that define both the local reference frame  used to determine 
 the bilayer s structure  and the area per lipid ( al ) along the undulating surface . 
 additionally , 
 the resulting power spectrum of spherical harmonic coefficients can 
 be fit to a helfrich continuum model for membrane bending in spherical 
 geometry to extract bending rigidity ( kc ) . 
 kc values determined for both dmpc 
 and dmpc + cholesterol ( 30 mol % ) vesicles are consistent with values 
 from corresponding flat - patch systems determined using an independent , 
 previously published spectral method . 
 these new tools to accurately 
 extract structure , al , and kc should prove invaluable in evaluating the construction 
 and equilibration of lipid vesicle simulations . 