voltage - dependent anion 
 channels ( vdacs ) are integral membrane 
 proteins in the mitochondrial outer membrane .
vdacs are related to 
 ancient porin proteins commonly seen in bacteria , and the eukaryotic 
 channels regulate ion and metabolite flux between the cytosol and 
 the mitochondrial intermembrane space .
crystallographic and nmr structures of mammalian vdac revealed a 
 -barrel with 19 strands connected by loops and an n - terminal 
 helix inside the channel . at low membrane potentials
, vdac assumes an open 
 conformation with modest anion selectivity , whereas applied voltages 
 greater than 3040 mv cause vdac to gate into a low conducting , 
 cation - selective state .
the vdac - mediated 
 confluence of the cytosol and mitochondrial intermembrane 
 space suggests low membrane potentials across the mitochondrial outer 
 membrane and predominantly open channels in vivo .
however , 
 in addition to voltage , vdac function is modulated by proteins that 
 include bcl - xl , tbid , hexokinase , and tubulin , and the interactions of these proteins with 
 vdac critically regulate apoptosis and cellular respiration .
lipids can also affect channel 
 activity by altering the properties of the surrounding membrane ( e.g. , 
 lateral pressure ) or through direct interactions with the protein .
in addition to endogenous regulators , vdac is a target of drugs 
 that include erastin and general anesthetics , although the functional influence of these ligands is not entirely 
 clear . here , we used atomistic molecular dynamics simulations 
 to investigate 
 vdac interactions with cholesterol .
cholesterol is known to bind mammalian 
 vdac in vivo and in vitro , and specific cholesterol sites on the protein have been 
 identified with nmr .
functionally , cholesterol 
 enhances the structural integrity of isolated vdac and aids channel 
 insertion into bilayers , and cholesterol promotes uniformity of open channel conductances .
together , these 
 effects suggest vdac function might be modulated when mitochondrial 
 outer membrane cholesterol content increases , such as in cancers .
the majority of proteins that bind cholesterol are integral 
 membrane 
 proteins ; however , cholesterol interactions 
 with -barrels , which in mammals consist of vdac and other mitochondrial 
 proteins , have not been extensively characterized . 

this is in contrast to extensive efforts investigating cholesterol 
 binding to other ion channels that reside on the plasma membrane .
therefore , this work represents an early step toward understanding 
 how cholesterol interacts with vdac , the effects of cholesterol binding 
 on channel structure and dynamics , and the potential interplay associated 
 with ligand and protein binding .
four residues in the mouse crystal structure 
 of vdac1 ( pdb code 3emn ) , which was solved at 2.3  resolution , 
 were mutated in pymol to humanize the 
 protein ( asp55  thr55 , val129  met129 , ser160  ala160 , val227  ile227 ) .
we specifically used 3emn because the conformation of the 
 n - terminal helix , which differs among other vdac1 structures , was confirmed with solid - state nmr .
the n - terminal methionine , water , and 
 detergent were removed before loading the structure into autodocktools . to the protein ,
hydrogens were added , nonpolar hydrogens were merged , 
 and kollman charges were added . molecular coordinates for cholesterol 

were downloaded from the charmm small molecule library ; in autodocktools , geisteiger charges were added , nonpolar hydrogens 
 were merged , and six torsions were allowed ( i.e. , cholesterol was 
 fully flexible ) . for docking ,
grid boxes targeted the specific cholesterol 
 sites and biased the exterior , membrane side of the -barrel . 

residue side chains in the site of interest and projecting outside 
 the barrel were flexible during docking runs with autodock vina .
the similarity of some docking 
 poses within each site often led to sampling of multiple poses during 
 preliminary simulations used to build the cholesterol - bound vdac model 
 ( see below ) . for this reason ,
multiple docking calculations with different 
 grid centers were performed on some sites to increase heterogeneity 
 in the starting poses .
after 
 docking , pdb files 
 of vdac , the docked cholesterol molecules , and corresponding flexed 
 residues were created .
the protein pore was oriented along the z - axis , normal to the membrane , using the ppm server of 
 orientations of proteins in the membrane , and the output models were loaded into the charmm - gui membrane builder . in all simulations , the n - terminal residue ( ala2 )
ala2 was acetylated , and residues were protonated according to their 
 standard states at ph 7.4 .
the water thickness ( minimum height on 
 the top and bottom of the system ) was 15.0  . for preliminary 
 simulations used to examine different orientations of cholesterol 
 binding to vdac , and for production simulations of the final cholesterol - bound 
 model ,
5 cholesterols were docked to the protein , and 160 1,2-dioleoyl - sn - glycero-3-phosphocholine ( dopc ) molecules ( 80 in each 
 leaflet ) and 11 randomly distributed cholesterols comprised the membrane . 

as a reference , phosphatidylcholines are the most 
 abundant mitochondrial outer membrane phospholipids at 4350 
 mol percent , with cholesterol 811 mol percent .
a total of 23 k and 25 cl ions were 
 randomly placed by the charmm - gui to neutralize the system and provide 
 a salt concentration of 0.15 m kcl .
kcl is generally favored over 
 nacl in vdac experiments to minimize the effects of ionic size and 
 diffusivity on measurements .
systems were separately generated with 
 the charmm - gui server , hence system details ( number of atoms , starting 
 coordinates , etc . ) differed between simulations .
systems were approximately 
 87  , 87  , and 74  in x , y , and z dimensions , respectively , with 
 about 56  000 total atoms that included about 9100 tip3p waters .
the charmm36 
 model was used for protein and phospholipid parameters , whereas the modified 
 ( charmm ) 36c model was used for cholesterol ; parameters for tip3p waters and ions are well - established .
. 
 interactions between nonbonded atoms were cut off at 12  , and 
 bonds involving hydrogen were constrained using the shake / rattle algorithms . 
 a langevin thermostat and barostat were used to maintain a temperature 
 and pressure of 303.15 k and 1 atm , respectively , and no surface tension 
 was imposed .
prior to production , 
 20  000 minimization steps , and a 0.725 ns equilibration protocol 
 detailed elsewhere were performed on 
 each system to gradually release restraints on the protein . in total , 
 22 separate simulations totaling over 1.4 s contributed to 
 this work , including 10 simulations of 100 ns .
generally , 
 the trajectories were aligned to and centered around the backbone 
 of the crystal structure before performing the following analyses : rmsf : root - mean - square fluctuation was measured for residue 
  carbons and describes the fluctuation of the atom relative 
 to its average position throughout each simulation . ion diffusion : 
 cl / k permeability ratios
were determined 
 after normalizing the number of ions that diffused through the pore 
 to the total number of that type of ion in the system ( 23 for k and 25 for cl ) .
diffusion was defined 
 as traversing the membrane through the channel pore from 20 
  to 20  , or vice versa , along z. the 
 protein spanned these coordinates in every frame of every simulation , 
 as detected with the hole software.channel radius : the minimum channel radius for each frame 
 was measured with the hole software , which fit consecutive spheres 
 inside the channel with the centers spaced every 0.5  along 
 the z axis .
hydrogen bonds : hydrogen 
 bonds were calculated with the hydrogen bond plugin of vmd with a 
 3.3  and 20 donor  acceptor cutoff .
electrostatics : the protein , membrane , and water contribution to the electrostatic 
 potential was generated with the pme electrostatics 
 plugin ( pmepot ) using 88  88  80 grid position counts 
 ( 1  spacing ) .
the generated 
 maps contained an identical number of data points from each simulation 
 ( the dimensions of the systems were also essentially identical ) .
potential of mean force : averaged multi - ion potential of 
 mean force plots were made for each trajectory . for each frame , we 
 measured both the number of ions and the number of water molecules 
 in a cylinder of radius 15  that ran parallel to the z - axis through the center of the channel .
this was used 
 to calculate the concentration of ions , c(z ) , in 1  bins , and the relative free energy was estimated 
 with the equation:1where r is the gas 
 constant , t is the temperature , and cbulk is the concentration of ions in the bulk water.ionic density : ionic density maps for entire trajectories 
 were generated with the volmap tool plugin of vmd .
structural figures 
 were generated with vmd and pymol , and graphical figures were generated 
 with graphpad prism ( graphpad software inc .
, la jolla , ca ) and kaleidagraph 
 ( synergy software , reading , pa ) ; movies were also generated with vmd . 

statistical analyses were performed as described in the text within 
 the graphpad prism software . for each analysis ,
the mean from each 
 simulation was considered a single experimental unit ; therefore , where 
 applicable , values are reported as the mean measurement from multiple 
 simulations with standard error .
published 
 nmr chemical shift 
 mapping experiments revealed resonances of nine vdac backbone amides 
 that had a significant chemical shift induced by the presence of cholesterol , 
 with the structure of the protein essentially unchanged ( figure s1 ) .
this 
 indicated a change in the chemical environment surrounding the backbone 
 amides and suggested specific binding by cholesterol . in that experiment , 

cholesterol was dissolved in a hydrophobic micellar phase and interacted 
 with vdac from outside the barrel , analogous to membrane cholesterol 
 approaching the protein . with binding demonstrated experimentally
, 
 we aimed to identify the effects of cholesterol site occupancy on 
 vdac . to initiate this
, we began a combined docking and molecular 
 dynamics approach targeting cholesterol to the specific vdac backbone 
 atoms implicated in binding .
the topology of the protein and initial 
 docking calculations suggested a model in which five cholesterol molecules 
 with unique orientations were required to simultaneously occupy all 
 residues that bind cholesterol .
a 5 to 1 mol ratio of cholesterol 
 to protein also reflected conditions under which chemical shift mapping 
 was performed , and is consistent with the number of cholesterols that 
 bind vdac in detergent .
for simplicity , the targeted 
 cholesterol sites are referred to numerically with ascending residue 
 composition ( site 1 : lys96 , thr116 , gly117 , asp128 ; site 2 : phe169 , 
 ala170 ; site 3 : gly172 ; site 4 : lys236 ; site 5 : ser260 ) ( figure s1 ) .
we docked cholesterol to each 
 site using autodock , which gave multiple 
 poses and reasonable docking scores for each site ( g  6.0 to 8.3 kcal / mol ) .
fifteen 
 atomistic molecular dynamics simulations , averaging 37.3 ns each , 
 were computed in an iterative approach .
we sampled multiple starting 
 poses in each site with the goal of optimizing a cholesterol - bound 
 vdac model ( figure s2 ) .
the primary criterion 
 for selecting favorable poses was sustained residency of cholesterol 
 in the sites .
unbinding of cholesterol oriented unfavorably and subsequent 
 diffusion into the membrane was observable on these time scales . the final cholesterol - bound vdac model used for subsequent production 
 simulations is shown in figure 1a .
the consistency 
 and fidelity of the model were measured in five 100 ns simulations . 
 as a correlate to chemical shift mapping , we measured the distance 
 between cholesterol and the backbone amide atoms in the sites throughout 
 the trajectories ( figure 1b and figure s3 ) .
in each 100 ns simulation , cholesterol 
 maintained a 36  intermolecular distance to 
 an amide atom from at least three sites .
additionally , each site was 
 essentially occupied for the full 100 ns in at least three simulations 
 ( figure 1b ) .
the unsustained binding we observed 
 in some simulations , however , suggested low affinity interactions , 
 and simultaneous occupation of these five sites can not be determined 
 at the nmr time scale at which the sites were initially identified .
( a ) optimized cholesterol - bound vdac model , with cholesterol 
 colored 
 magenta .
the n - terminal helix is colored dark blue , and the residues 
 within each site are colored yellow .
( b ) minimum distance between 
 cholesterol and an amide atom ( n or h ) of a residue within that site . 

multiple intermolecular contacts 
 contributed to stable binding 
 in each site ( table 1 and movies s1s3 ) .
the cholesterol molecules primarily 
 oriented in grooves defined by ridges of hydrophobic and sometimes 
 aromatic side chains that alternately project into the membrane from 
 adjacent -strands ; sites 2 and 3 are defined by a common ridge , 
 as are sites 4 and 5 .
while we docked cholesterol to the protein in vacuo , we anticipated favorable binding in simulations 
 would be enhanced by the sterol hydroxyl forming hydrogen bonds .
the 
 hydroxyls of sites 4 and 5 cholesterols extended into the membrane , 
 contacting water and lipid head groups , and the cholesterol in site 
 1 formed significant hydrogen bonds with tyr153 .
cholesterol in sites 
 1 and 3 were also accessible to water funneling from the edge of the 
 protein .
specifically , water was available to site 1 through tyr153 
 and the adjacent ser167 , which constitutively contacted water in the 
 channel , and could reach site 3 between the side chains of loop residues 
 his122 and tyr146 , which also contacted bulk water .
in contrast 
 to the above , cholesterol bound the protein parallel 
 to the membrane with the hydroxyl among lipid tails in site 2 , stably 
 wedged between side chains , and inaccessible to bulk water and lipid 
 head groups .
the initial coordinate of this cholesterol hydroxyl was 
 1.3  from the bilayer midplane .
membrane cholesterol commonly 
 assumes an upright orientation in membranes of saturated phospholipids , 
 with the hydroxyl at the hydrophilic interface 16  from 
 the bilayer midplane .
however , in the 
 presence of polyunsaturated fatty acids , membrane cholesterol can 
 assume a flat orientation in the middle of the bilayer , not unlike cholesterol in site 2 ; therefore , cholesterol could assume 
 this orientation and bind protein without initial placement , especially 
 in more complex membranes .
site 1 cholesterol hydrogen bonded 
 with tyr153 for 4050 ns in separate simulations , and site 
 5 cholesterol hydrogen bonded with asn269 for 12 ns in one simulation .
our initial simulations were 
 designed to identify favorable orientations 
 of cholesterol in the experimentally identified sites . validating 
 our approach , we observed replacement of a docked cholesterol by a 
 randomly placed membrane cholesterol , which assumed an identical orientation 
 in site 4 ( movie s4 ) .
this should theoretically 
 be observed at all sites ; however , measuring multiple unbinding and 
 rebinding events at each site is beyond currently feasible time scales 
 and comes at an uncertain computational cost . while cholesterol potentially 
 binds in other orientations or sites with physiologically relevant 
 affinities , the present model identifies several binding modes that 
 reasonably and reproducibly satisfy experimental binding data .
we used the cholesterol - bound 
 vdac model to investigate the effects of the sterol on the protein . 

five 100 ns simulations with an identical mole fraction of cholesterol 
 in the membrane , but not docked to the protein , were computed ( the 
 protein in these simulations is referred to as  apo - vdac  
 as opposed to  cholesterol - bound vdac  using our model ) . 
 cholesterol did not achieve binding poses in apo - vdac simulations 
 that were equivalent to cholesterol - bound vdac .
our experimental approach 
 therefore provided comparable systems that decoupled the effects of 
 cholesterol as a ligand , which was investigated here , from its influence 
 on the membrane .
five 100 ns replicates of vdac with and without bound 
 cholesterol allowed for a measure of reproducibility and biased our 
 sampling toward the intended apo- or cholesterol - bound states .
the protein backbone equilibrated rapidly and remained stable in 
 the presence or absence of bound cholesterol ( figure s4 ) , and the rmsd of backbone  carbons from 
 averaged structures of apo- or cholesterol - bound vdac was 0.6  . 

we calculated the average rmsf of the residue  carbons throughout 
 the trajectories ( figure 2a , b ) . with cholesterol 
 bound , 64  7% of residues per simulation
residue stabilization 
 by cholesterol was generally global but not absolute ; for example , 
 residues 130138 ( between -strands 8 and 9 ) were consistently 
 more dynamic ( figure 2c ) . for any individual 
 residue ,
however , the absolute difference in averaged rmsf ( i.e. , 
 the average rmsf from cholesterol - bound simulations minus that from 
 apo - vdac simulations ) was quantitatively small and between 0.63 
  and + 0.55  .
generally , the loops that connect -strands 
 were most dynamic 
 ( figure 2a , b ) , as previously observed . in an aqueous environment
, these loops lack the structural support 
 provided by the membrane and inter - residue electrostatic interactions . 
 even without direct cholesterol binding ,
the ability of distant loops to be 
 affected by changes in global protein dynamics has been well - characterized 
 in soluble proteins . here , vdac residues that bind cholesterol might form 
 networks with adjacent amino acids that relay the presence of cholesterol 
 to the loops . across these loops ,
the largest percent increases in 
 averaged rmsf were seen in residues ala134 and lys274 , and the largest 
 percent decreases were seen in residues pro105 , gly187 , and asn269 
 ( figure 2d ) .
( a ) the average rmsf value for each residue 
 from five simulations 
 of apo - vdac and ( b ) cholesterol - bound vdac systems .
standard error 
 ( n = 5 ) is indicated , and when not visible , the error 
 was smaller than the size of the point .
( c ) the number of simulations 
 ( out of 5 ) that each residue from cholesterol - bound vdac simulations 
 had increased  carbon rmsf relative to the average from apo - vdac 
 simulations .
( d ) percent change in rmsf averaged across five cholesterol - bound 
 and five apo - vdac simulations .
points below the dotted lines indicate 
 that the residue was less dynamic in cholesterol - bound vdac simulations 
 relative to apo - vdac simulations .
we investigated further the effects of cholesterol on the 
 n - terminus , 
 as the helix and cholesterol bind the same wall of the -barrel , 
 but inside and outside of the pore , respectively ( figure s1 ) .
bonding pairs between helix residues 225 
 and the barrel wall were variable both within and between the two 
 groups ( with or without bound cholesterol ) ( table 
 s1 ) , and this variability extended to residues that hydrogen 
 bond in the crystal structure .
however , 
 the helix maintained 3.1  0.1 and 2.9  0.1 hydrogen bonds 
 per frame with the barrel wall throughout apo - vdac and cholesterol - bound 
 vdac simulations , suggesting a dynamic network maintained helical 
 stability .
separate from hydrogen bonding , channel stability is also 
 both dependent on and indicated by a hydrophobic contact between helix 
 residue leu10 and barrel residue val143 . the average distance ( 3  ) and maximum distance ( 5 
  ) between these residues throughout simulations were further 
 consistent with a constitutively open channel for both apo- and cholesterol - bound 
 vdac .
interestingly , we detected a cholesterol - induced increase in 
 the average minimum channel radius from 5.5  0.1  ( apo - vdac ) 
 to 5.8  0.1  ( cholesterol - bound vdac ) using an algorithm 
 that fits spheres inside the channel along the z axis ( comparison of minimum channel radii : p < 0.05 , two - tailed t test , n = 5 simulations each group ) . the increased channel radius 
 with cholesterol bound does not necessarily imply greater area for 
 water and ion accessibility ; however , it does suggest larger circular 
 or less elliptical constriction sites .
we generated 
 averaged electrostatic potential maps , which revealed entirely positive 
 net potentials inside the anion channel .
cholesterol - induced protein 
 rigidity resulted in larger regions of strong positive potential , 
 with the strongest potentials in the channel arising from residues 
 with decreased positional fluctuation in cholesterol - bound vdac ( figure 3 ) .
mutagenesis of multiple residues in these positions 
 modulates ion permeability under low voltage in silico and in vitro .
similarly , negatively 
 charged residues that were also stabilized by cholesterol are located 
 along the barrel wall opposite the n - terminal helix , providing a path 
 for k ions through cholesterol - bound vdac .
this was seen 
 adjacent to glu84 ( figure s5 ) , previously 
 shown to influence k permeability through vdac in silico .
despite changes in charge distribution , there was no significant 
 change in the average number of cl or k ions that diffused through the channel in each simulation ( number 
 of ions that translocate per simulation  cl : 26      3 ( cholesterol - bound ) and 23  
 1 ( apo - vdac ) ; k : 10  2 ( cholesterol - bound ) and 7 
  1 ( apo - vdac ) ) .
the cl / k permeability 
 ratios from apo- and cholesterol - bound vdac were 3.3  0.4 and 
 2.8  0.7 , in reasonable agreement with the experimental value 
 of 2 for the open channel .
we 
 estimated the potential of mean force for cl and 
 k as a function of the ion height z in the pore ( figure 4 ) . as described previously , free energy wells for cl were observed inside 
 the channel on either side of the n - terminal helix , which is the source 
 of the strongest positive potentials ( figure 3 ) . the free energy well at 10  < z < 2.5  coincided with the location of the minimum 
 pore radius , which was also the location of the free energy maxima 
 for k. the averaged magnitude and shape of the k profiles are similar to previous computational work despite variability between simulations at 0  < z < 10  , possibly due to under - sampling of k inside the anion channel .
electrostatic potential maps showing positive 
 potential inside 
 apo - vdac and cholesterol - bound vdac .
blue electrostatic potential 
 maps correspond to ( a ) + 64 kt / e and ( b ) + 38 kt / e . the n - terminal helix 
 is colored green , and the protein backbones correspond to the averaged 
  carbon structures from the simulations with or without bound 
 cholesterol . the n- and c - termini are at the left entrance ( in lateral 
 structures ) or on the near side of the channel ( when looking into 
 the barrel ) . in lateral structures ,
the orange spheres correspond 
 to the  carbons of leu69 and ser101 which , for reference , are 
 located at approximately z  8 and z  8 , respectively . 
 in ( a ) , lys12 , lys20 , and lys236
are shown in blue sticks , and they 
 are the charged residues that predominantly contribute to the electrostatic 
 field .
arg15 , lys32 , lys119 , lys174 , and lys224 are also shown in 
 ( b ) .
all of the mentioned residues are stabilized by cholesterol ( i.e. , 
 have a decreased average rmsf ) .
( b ) k ions along the z axis , with and 
 without bound cholesterol , estimated using eq 1 .
by 
 examining the interactions of vdac with cholesterol , which is 
 both a membrane component and endogenous ligand , we demonstrate stable binding and alteration of the dynamics 
 of the mitochondrial anion channel . guided by experimental data , we 
 identified reasonable binding orientations for cholesterol in specific 
 vdac sites .
although we observed replicable cholesterol binding , it 
 is unclear from our simulations whether the five investigated sites 
 would be simultaneously occupied in vivo .
however , occupation of all 
 five sites does not appear necessary to affect vdac dynamics .
experimentally , 
 in addition to modifying the starting positions of molecules , the 
 stoichiometry of site occupancy can be modulated by the replacement 
 rate of vacated sites by cholesterol from the membrane , which in turn 
 is affected by membrane content .
mitochondrial outer membrane cholesterol 
 content varies in disease states to potentially modulate vdac function , 
 with 210-fold increases in cholesterol per mg of protein observed , suggesting higher site occupancy under these conditions .
however , 
 varying membrane cholesterol content might also affect vdac through 
 membrane - mediated mechanisms in addition to direct effects of specific 
 protein binding , which was investigated 
 here .
regulation of membrane protein function by lipid composition 
 can likely occur by exerting the physical properties of the membrane 
 and its components through specific protein sites .
studies of -barrels , 
 whether computational or experimental , should therefore occur in the 
 presence of at least physiologically relevant cholesterol concentrations
. 
 this will be particularly critical as investigation proceeds toward 
 clinically relevant and structurally similar sterols that also bind 
 these proteins , potentially at similar 
 sites .
 the 
 mitochondrial voltage - dependent anion channel ( vdac ) allows 
 passage of ions and metabolites across the mitochondrial outer membrane . 
 
 cholesterol binds mammalian vdac , and we investigated the effects 
 of binding to human vdac1 with atomistic molecular dynamics simulations 
 that totaled 1.4 s . 
 we docked cholesterol to specific sites 
 on vdac that were previously identified with nmr , and we tested the 
 reliability of multiple docking results in each site with simulations . 
 
 the most favorable binding modes were used to build a vdac model with 
 cholesterol occupying five unique sites , and during multiple 100 ns 
 simulations , cholesterol stably and reproducibly remained bound to 
 the protein . for comparison 
 , vdac was simulated in systems with identical 
 components but with cholesterol initially unbound . 
 the dynamics of 
 loops that connect adjacent -strands were most affected by 
 bound cholesterol , with the averaged root - mean - square fluctuation 
 ( rmsf ) of multiple residues altered by 2030% . 
 cholesterol 
 binding also stabilized charged residues inside the channel and localized 
 the surrounding electrostatic potentials . despite this , 
 ion diffusion 
 through the channel was not significantly affected by bound cholesterol , 
 as evidenced by multi - ion potential of mean force measurements . 
 although 
 we observed modest effects of cholesterol on the open channel , our 
 model will be particularly useful in experiments that investigate 
 how cholesterol affects vdac function under applied electrochemical 
 forces and also how other ligands and proteins interact with the channel . 
Introduction
Methods
Results
and Discussion
Conclusions