Dinh/Dinh 2013/NOTES/2013-1-7: Difference between revisions

From ZhangLabWiki
Jump to navigation Jump to search
>Dinh
mNo edit summary
>Dinh
mNo edit summary
Line 4: Line 4:
* First map reads using a shell script:
* First map reads using a shell script:
* Edit '''Go.mapBisulfite.sh''' as follows:
* Edit '''Go.mapBisulfite.sh''' as follows:
  ### EDIT below ####
  ### EDIT below ####
  # current directory: change the following path to your working directory.
  # current directory: change the following path to your working directory.
  '''cur_dir="/oasis/triton/scratch/ddiep/Working/WGBS_Noi_hg19/HELLO"'''
  '''cur_dir="/oasis/triton/scratch/ddiep/Working/WGBS_Noi_hg19/HELLO"'''
  # reads directory: change the following path to the reads directory
  # reads directory: change the following path to the reads directory
  # reads_dir and cur_dir doesn't have to be the same
  # reads_dir and cur_dir doesn't have to be the same
  '''reads_dir="/oasis/triton/scratch/ddiep/Working/WGBS_Noi_hg19/HELLO"'''
  '''reads_dir="/oasis/triton/scratch/ddiep/Working/WGBS_Noi_hg19/HELLO"'''
  # List the index names (or unique names to each read)
  # List the index names (or unique names to each read)
  '''INDX="Indx1 Indx2 Indx3"'''
  '''INDX="Indx1 Indx2 Indx3"'''
Line 32: Line 29:
  soap="/home/ddiep/softwares/soap2.21release/soap"
  soap="/home/ddiep/softwares/soap2.21release/soap"
  ############################################################
  ############################################################
  cd $cur_dir
  cd $cur_dir
  for n in ${INDX}
  for n in ${INDX}
Line 51: Line 47:
   echo "#PBS -A zhang-lab" >> $n.job
   echo "#PBS -A zhang-lab" >> $n.job
   echo "cd $cur_dir" >> $n.job  
   echo "cd $cur_dir" >> $n.job  
   echo "$bisReadMapper -r $reads_dir/$f,$reads_dir/$g -m $mm -W $template_fwd -C $template_rev -g $template_fa -a $soap -b $qual -p $p -n $n -q $qtrim -5 $fivetrim -3 $threetrim > $n.status" >> $n.job
   echo "$bisReadMapper -r $reads_dir/$f,$reads_dir/$g -m $mm -W $template_fwd -C $template_rev -g $template_fa -a $soap -b $qual -p $p -n $n -q $qtrim -5 $fivetrim -3 $threetrim > $n.status" >> $n.job
   echo "rm *encoded" >> $n.job
   echo "rm *encoded" >> $n.job

Revision as of 02:53, 11 January 2013

bisReadMapper pipeline

Triton

  • First map reads using a shell script:
  • Edit Go.mapBisulfite.sh as follows:
### EDIT below ####
# current directory: change the following path to your working directory.
cur_dir="/oasis/triton/scratch/ddiep/Working/WGBS_Noi_hg19/HELLO"
# reads directory: change the following path to the reads directory
# reads_dir and cur_dir doesn't have to be the same
reads_dir="/oasis/triton/scratch/ddiep/Working/WGBS_Noi_hg19/HELLO"
# List the index names (or unique names to each read)
INDX="Indx1 Indx2 Indx3"

email="ddiep@ucsd.edu"
mm=2 # max number of allowable mismatches
qual=64 # base quality offset value
qtrim=20 # Phred quality value for soft-trimming
threetrim=0 # number of bases to trim from 3' end (rightmost)
fivetrim=0 # number of bases to trim from 5' end (leftmost)
p=8

### DO NOT EDIT, unless you know where these files are ####
# the following paths should stay the same on triton. 
bisReadMapper="/home/ddiep/scripts/MethylationPipeline/smartBisReadMapper.pl"
template_fwd="/projects/zhang-lab/ddiep/LatestGenome/BisRef/bisHg19/hg19.fa.bis.fwd.index"
template_rev="/projects/zhang-lab/ddiep/LatestGenome/BisRef/bisHg19/hg19.fa.bis.rev.index"
template_fa="/projects/zhang-lab/ddiep/LatestGenome/BisRef/bisHg19/hg19.fa"
soap="/home/ddiep/softwares/soap2.21release/soap"
############################################################
cd $cur_dir
for n in ${INDX}
do
 ###
 f="s_1_1_$n.txt" # name format for read 1, make sure it matches the file names
 g="s_1_2_$n.txt" # name format for read 2, make sure it matches the file names
 ###
 echo "#!/bin/csh" > $n.job
 echo "#PBS -q small" >> $n.job
 echo "#PBS -l nodes=1:ppn=8" >> $n.job
 echo "#PBS -l walltime=36:00:00" >> $n.job
 echo "#PBS -o $n.log" >> $n.job
 echo "#PBS -e $n.err" >> $n.job
 echo "#PBS -V" >> Idx$n.job
 echo "#PBS -M diep.hue.dinh@gmail.com" >> $n.job
 echo "#PBS -m abe" >> $n.job
 echo "#PBS -A zhang-lab" >> $n.job
 echo "cd $cur_dir" >> $n.job 
 echo "$bisReadMapper -r $reads_dir/$f,$reads_dir/$g -m $mm -W $template_fwd -C $template_rev -g $template_fa -a $soap -b $qual -p $p -n $n -q $qtrim -5 $fivetrim -3 $threetrim > $n.status" >> $n.job
 echo "rm *encoded" >> $n.job
 qsub $n.job
done
  • Run the script:
sh Go.mapBisulfite.sh > job_ids
  • Now, there will be a different job for each sequence file that you have.
  • If you have multiple sequencing lanes, just change the f and g variables to s_2, s_3, etc, and re-run the shell script with the command above.