User Tools

Site Tools


Writing /users/jbklauda/pub/wiki/data/cache/6/6193a558dc90ed4c190004bf67ecff70.captchaip failed
tilt_angle

Differences

This shows you the differences between two versions of the page.

Link to this comparison view

Both sides previous revisionPrevious revision
Next revision
Previous revision
tilt_angle [2016/11/16 15:42] – admintilt_angle [2025/01/17 09:46] (current) – admin
Line 2: Line 2:
  
 **Brief instruction to obtain tilt angle distribution** \\ **Brief instruction to obtain tilt angle distribution** \\
-Authors: Xiaohong Zhuang and Viviana Monje+Authors: Xiaohong Zhuang and Viviana Monje\\ 
 +Updated: Dr. Jeff Klauda (2021)
  
  
-This brief instruction is used to calculate tilt angle distribution of sterols/cholestrols. The scripts that you need are in DT1 path: **/lustre/xzhuang/simulation/script_folders/lipids_analyses/tilt** +This brief instruction is used to calculate tilt angle distribution of sterols/cholestrols. The scripts that you need are in ZT1 path: **/afs/shell.umd.edu/project/energybio/shared/jbklauda/scripts/tilt
-Or DT2 path: **/lustre/xzhuang/scripts/lipids_analyses/tilt+
 ** **
 +or using zip file: {{ :tilt.gz |}}
  
-The example of the outputs is in the DT2 path: **/lustre/xzhuang/soybean/1_analyses/hybo-b1/tilt** +In each sbatch file, the ZT1 charmm are provided. 
- +
- +
-In each sbatch file, the DT1 and DT2 charmm are provided. So if you work in DT1, please comment out the line with DT2 CHARMM by # sign, and uncomment the DT1 CHARMM (more detail is shown in the step 5 below). +
  
 Copy the scripts folder “tilt” from path highlighted above to your directory that have all the dyn files. Copy the scripts folder “tilt” from path highlighted above to your directory that have all the dyn files.
  
 +1. Update the //rtfprm.str//, if needed and copy toppar files.
  
-The following scripts are the example setup used to calculate the title angle of beta-sitosterol (SITO) of the hypocotyl-b1 system. +Copy the toppar.str and the toppar subdirectory from the CHARMM-GUI build directory (for the CHARMM files not NAMD). Then, update the PSF and CRD file in the rtfprm.str associated with your files, if needed.
- +
- +
-====== 1. ====== +
-1.1 Update the //rtfprm.str// +
-Update the path for the topology stream file and the psf and pdb files +
- +
- +
-e.g. +
-1.1  In //rftprm.str//, specify the path of the //toppar.str// file (it is toppar_less_1.str for hypo system), psf and crd file  +
- +
- +
-1.2 In the //toppar_less_1.str// file, +
- Replace the short path: //toppar///  +
- by the full path to your CHARMM toppar folder  +
-e.g. **/lustre/xzhuang/soybean/hypocotyl-b1/toppar/** +
- +
- +
-Note: +
-“Replace” command in WinSCP would save time to replace the short paths by full path. Since the toppar.str is streamed in the current path, the system will try to find the toppar folder from the current path. So in order to let system to use the correct path of the *rtf and *prm files, the full path is provided in toppar.str file.+
  
  
 2. Update the //set-range.str //  2. Update the //set-range.str // 
 <code commend> <code commend>
-set top1 = 145 ! first resid of lipid of interest on the top leaflet +set rnc = CHL1   ! Update the resname of sterol 
-set top2 = 165 ! last resid of lipid of interest on the top leaflet +set ti = 10      ! Update the dyn number first trajectory 
-set bot1 = 166 ! first resid of lipid of interest on the bottom leaflet +set tf = 11      ! Update the dyn number of last trajectory 
-set bot2 = 186 ! last resid of lipid of interest on the bottom leaflet +set tra = 18     ! Update the value of (first trajectory plus 8) to define unit# to be read
-set lip = 21   ! total number of lipids per leaflet +
-set ti = 41     ! first trajectory of interest +
-set tf = 75    ! last trajectory of interest +
-set tra = 49    ! first trajectory ti plus 8 to define the first unit# +
 </code> </code>
  
-Note: +Note:\\ 
-a. The resid range can be found in the //step5_assembly.crd// file, usually the first half resid is on the top leaflet, and rest half on bottom leaflet in hypo system. You may confirm the top and bottom leaflet by VMD +a. Choose ti and tf based on the equilibrated data. So if SA/lipid shows that you have equilbration at dyn10 and higher then use ti=10 and tf is the number of DCD files.\\ 
- +b. In this example the sterol of interest is CHL1. Other sterols can be used and rnc is changed based on their resname.\\ 
- +c. We use @tra=@ti+8 because in tilt.inp, the first unit number is calculated by calc u = 8 + @i \\
-b. The lipid of interest for soybean system is beta-sitosterol (SITO) and stigmasterol (STIG) for hypo system. The above setup is for SITO, which has 21 SITO in each leaflet. Assume system reach equilibrium by dyn 41, so dyn range 41 to 75 are used to calculate average angle distribution.\\ +
-c. We use @tra=@ti+8 because in tilt.inp, the first unit number is cacluated by calc u = 8 + @i \\+
 … \\ … \\
 ''traj firstu @tra nunit @N''\\ ''traj firstu @tra nunit @N''\\
Line 64: Line 37:
 3. Update the  //tilt.inp//  3. Update the  //tilt.inp// 
  
-1.1 Update the data points per DYN/dcd file\\+3.1 Update the data points per DYN/dcd file\\
 <code commend> <code commend>
 calc mxs = 10 + 16*100 ! TIME SERIES PER MOLECULE +10   calc mxs = 10 + 16*100 ! TIME SERIES PER MOLECULE +10  
-calc mxt = 1000 * @N  ! update first number to data points per DYN file+calc mxt = 500 * @N  ! update first number to data points per DYN file
 </code> </code>
  
  
-Note: +Note:\\ 
-For mxs, it is just used to save space, it depends on the lines of calculation performed (in 1.2 below). It does not need to be exact. Can increase the number if the error occurs because the space (time seriers) limitation.+For mxs, it depends on the lines of calculation performed (in 1.2 below) and for most cases does not need to change. You can increase the value if the error occurs because the space (time series) limitation.
    
  
-For mxt, 1000 is # of data points/frames per DYN file (calculated by # of total time steps divided by dcd frequency, both can be obtained from dyn2-inp, e.g. 1000000/1000=1000)+For mxt, 500 is # of data points/frames per DYN file (calculated by # of total time steps divided by dcd frequency, both can be obtained from dyn2-inp, e.g. numsteps/dcdfreq=500)
  
  
-1.2 Update the reference atoms to calculate the tilt angles +3.2 Update the reference atoms to calculate the tilt angles 
  
  
Line 87: Line 60:
  
  
-Note:+Note:\\
 MEMB is segid @t is resid (which is based on the  resid range of  top1,top2,bot1 and bot 2 set up in set_range.str), C3 (first C close to the OH group near head end) or C17 (carbon at the "base" of the last rings near tail end) are atom type names selected (which is more clear when you look at the topology file). The reference vectors are created based on connecting the reference atoms C3 and C17. Then, the angle between reference vector and the bilayer normal are calculated.  MEMB is segid @t is resid (which is based on the  resid range of  top1,top2,bot1 and bot 2 set up in set_range.str), C3 (first C close to the OH group near head end) or C17 (carbon at the "base" of the last rings near tail end) are atom type names selected (which is more clear when you look at the topology file). The reference vectors are created based on connecting the reference atoms C3 and C17. Then, the angle between reference vector and the bilayer normal are calculated. 
  
 +3.3  Update location of DCD file if needed.
 +
 +<code commend>
 +calc u = 8 + @i
 +open unit @u file read name ../dyn@i.dcd
 +incr i by 1
 +</code>
 +
 +Update the path for the DCD file if needed.
  
 4. Update the //top.inp// and //bot.inp// 4. Update the //top.inp// and //bot.inp//
Line 100: Line 82:
 </code> </code>
  
-The variable values in top.inp are the following (which can be found in dist.c file):\\ +The variable values in top.inp are the following (which can be found in dist.c file): 
-35000        ! # of total data points/frames (= # frames per dcd times # dyn, e.g. 1000 X 35)\\ +<code commend> 
-90.10         ! Maximum (value, which is the one step larger than the real maximum (90+0.10)\\ +1000        ! # of total data points/frames (= # frames per dcd times # dyn, e.g. 500 X 2) 
-0                    ! Minimum (angle) value\\ +90.10         ! Maximum (value, which is the one step larger than the real maximum (90+0.10) 
-t-all-nc1.dat  ! input file name,\\ +0                    ! Minimum (angle) value 
-t-tilt.dat        ! output file name with the distribution\\+t-all-nc1.dat  ! input file name 
 +t-tilt.dat        ! output file name with the distribution 
 +</code>
  
  
 Note:\\ Note:\\
-//t-all-nc1.dat// is same as //t-all.dat// but with the first time step columns removed. Since //dist.c// calculate distribution (with bin size 0.1) from the first column, the input file is //t-all-nc1.dat// (you won’t see it as it is deleted after the calculation) is used instead of t-all.dat which having the first time steps columns.+The only change you likely need is the fist line based on the number of frames per DCD and number of DCD file. The other lines are described here. //t-all-nc1.dat// is same as //t-all.dat// but with the first time step columns removed. Since //dist.c// calculate distribution (with bin size 0.1) from the first column, the input file is //t-all-nc1.dat// (you won’t see it as it is deleted after the calculation) is used instead of t-all.dat which having the first time steps columns.
  
  
Line 119: Line 103:
    
  
-5.1 Update the charm location depends on you run the analyses on DT1 or DT2. Use # to comment out the one you don’t use\\ +5.1 Update the charm location if needed but the current working code is typically in this file:\\ 
-# In DT1\\ + 
-# set chm = **/lustre/home/jbklauda/charmm/c39a1/exec/gnu/gfort-xxlg** \\+<code commend>
 # In DT2\\ # In DT2\\
-set chm = **/lustre/xzhuang/software/charmm/gnu/gfort-xxlg** \\+set chm = **/homes/jbklauda/charmm/c43b1/exec/gnu/charmm** \\ 
 +</code>
  
 +5.2 Update the path to where you want to save your final tilt distributions file, it can be a common folder than you want to save all your tilt angle data, or it can be current folder path.
  
-5.2Update the path to where you want to save your final tilt distributions file, it can be a common folder than you want to save all your tilt angle data, or it can be current folder path. 
  
- +''set all = "/lustre/jbklauda/lipid/analysis-test/gl-chl1/tilt"''
-''set all = "/lustre/xzhuang/soybean/all_analyses_summary/tilt"''+
  
  
Line 135: Line 119:
  
  
-Note: \\ +Note:\\ 
-In the following commands, \\ +In the following commands: 
-# INVOKE CHARMM\\ + 
-$''chm N=$1 < tilt.inp > tilt.out'' \\ +<code commend> 
-# GET DISTRIBUTION \\ +# INVOKE CHARMM 
-# especify last and first tilt data file \\+$''chm N=$1 < tilt.inp > tilt.out''  
 +# GET DISTRIBUTION  
 +# especify last and first tilt data file 
 ''./dist.csh $2 $3'' ''./dist.csh $2 $3''
 +</code>
  
  
Line 157: Line 144:
    
  
-''sbatch -A energybio-hi gettilt.csh 35 21 1 hypo-b1'' \\+''sbatch -A energybio-hi gettilt.csh 2 20 1 chl1'' \\
  
  
Line 163: Line 150:
  
  
-As shown in the get.scr file, the four variables are \\ +As shown in the 1_get.scr file, the four variables are \\ 
-75: Total number of dcd files, which will be value of N in //gettilt.csh// \\ +2: Total number of dcd files to be read, which will be value of N in //gettilt.csh// \\ 
-21: last tilt file # which is the # of lipid of interest, \\ +20: last tilt file number, which is the number of lipids of interest per leaflet (not the total), \\ 
-1: first tilt file # \\ +1: first tilt file number \\ 
-hypo-b1: name to copy final dat file to common directory \\+chl1: name to copy final dat file to common directory \\
  
  
-After you run '' get.scr command $chm N=$1 < tilt.inp > tilt.out'', the files// t1.dat …t#.dat, b1.dat…b#.dat// will be created. Each t#.dat is a tilt file, the # of tilt files depends on the number of lipid of interest. Therefore, last tilt file # is the # of lipid of interest.+7.  Submit the job by command: ''1_get.scr'' 
  
  
-7.  Submit the job by command: ''1_get.scr'' \\ 
 //T1.dat …t#.dat and b1.dat … b#.dat// will be obtained, and the final top and bottom leaflet averaged tilt angle distribution is provided in //avg-tilt-nh.dat//. The first column is angle, and the second column is normalized probability in each angle. //T1.dat …t#.dat and b1.dat … b#.dat// will be obtained, and the final top and bottom leaflet averaged tilt angle distribution is provided in //avg-tilt-nh.dat//. The first column is angle, and the second column is normalized probability in each angle.
  
Line 180: Line 166:
  
  
-8.Plot the //avg-tilt-nh.dat// to visualize the distribution+8.Plot the //avg-tilt-nh.dat// to visualize the distribution in a program of your choice.  
 + 
 + 
 +//**Some additional details of the codes here:**//
  
  
-Additonal note: +This section is to aid in the understanding the codes for this tilt calculation. 
- (This section only written for helps understanding the code, so it is fine if you don’t understand since no modification is required)+
  
  
Line 194: Line 182:
 </code> </code>
  
-2. The numbers of 21 and 1 in //get.scr// are defined as the variables+2. The numbers of 20 and 1 in //get.scr// are defined as the variables
  $2 $3 in the command line in the //gettilt.csh// file,   $2 $3 in the command line in the //gettilt.csh// file, 
  ''./dist.csh $2 $3 ''\\  ''./dist.csh $2 $3 ''\\
Line 206: Line 194:
  
 because the variable number are defined based on the order of specification when a script is executed. Since 4 variables are specified when gettilt.csh is excecuted by command  line in //1_get.scr// \\ because the variable number are defined based on the order of specification when a script is executed. Since 4 variables are specified when gettilt.csh is excecuted by command  line in //1_get.scr// \\
-''sbatch -A energybio-hi gettilt.csh 35 21 1 hypo-b1''+''sbatch -A energybio-hi gettilt.csh 2 20 1 chl1''
  
  
tilt_angle.1479328925.txt.gz · Last modified: 2016/11/16 15:42 (external edit)