Skip to content

Instantly share code, notes, and snippets.

@arq5x
Last active December 21, 2015 01:58
Show Gist options
  • Save arq5x/6231203 to your computer and use it in GitHub Desktop.
Save arq5x/6231203 to your computer and use it in GitHub Desktop.
GBM
export SAMPLES="BLV4 NCH411GBM_CD133high NCH411GBM_CD133low NCH537P54_CD133neg NCH537P54_CD133pos NCH620P55_CD133neg NCH620P55_CD133pos NCH644GBM_CD133high NCH644GBM_CD133low NCH7Md_P43_CD133neg NCH7Md_P43_CD133pos NPC-v"
######################################
# Sort the original BAM files by name:
######################################
export GBMHOME=/net/midtier18/vol79/cphg-quinlan2/projects/gbm-seq-abounader
export STEPNAME=gbm-nmsrt
for sample in `echo $SAMPLES`
do
export QSUB="qsub -q cphg -W group_list=CPHG -V -l select=1:mem=6000m:ncpus=1 -N $STEPNAME -m bea -M [email protected]"
echo "cd $GBMHOME; samtools sort -n -m 2000000000 bam/$sample.rmdup.bam bam/$sample.rmdup.namesorted" | $QSUB
done
############################################################
# Create new FASTQ files for every sample
############################################################
#export GBMHOME=/net/midtier18/vol79/cphg-quinlan2/projects/gbm-seq-abounader
#export STEPNAME=gbm-b2fq
#for sample in `echo $SAMPLES`;
#do
#export QSUB="qsub -q cphg -W group_list=CPHG -V -l select=1:mem=8000m:ncpus=1 -N $STEPNAME -m bea -M [email protected]";
# echo "cd $GBMHOME; java -Xmx2g -jar /home/arq5x/cphg-home/bin/SamToFastq.jar \
# INPUT=bam/$sample.rmdup.namesorted.bam \
# VALIDATION_STRINGENCY=LENIENT \
# F=fastq/$sample.1.fq \
# F2=fastq/$sample.2.fq" | $QSUB;
#done
############################################################
# Create new FASTQ files for every sample
############################################################
export GBMHOME=/net/midtier18/vol79/cphg-quinlan2/projects/gbm-seq-abounader
export STEPNAME=gbm-b2fq
for sample in `echo $SAMPLES`;
do
export QSUB="qsub -q cphg -W group_list=CPHG -V -l select=1:mem=8000m:ncpus=1 -N $STEPNAME -m bea -M [email protected]";
echo "cd $GBMHOME; ~/shared/bin/bamUtil_1.0.9/bamUtil/bin/bam bam2FastQ --in bam/$sample.rmdup.bam \
--outBase fastq/$sample" | $QSUB;
done
############################################################
# Novoalign
############################################################
export GENOME=/home/arq5x/cphg-home/shared/genomes/hg19/bwa/gatk/hg19_gatk.fa.novo.k14.s1.idx
export GBMHOME=/net/midtier18/vol79/cphg-quinlan2/projects/gbm-seq-abounader
export STEPNAME=gbmnovo
for sample in `echo $SAMPLES`;
do
export QSUB="qsub -q cphg -W group_list=CPHG -V -l select=1:mem=32000m:ncpus=16 -N $STEPNAME -m bea -M [email protected]";
echo "cd $GBMHOME; novoalign -d $GENOME -o SAM $'@RG\tID:$sample\tSM:$sample' -r Random \
-f fastq/${sample}_1.fastq fastq/${sample}_2.fastq \
-c 15 \
| samtools view -Sb - > bam/$sample.novoalign.rg.bam" | $QSUB
done
############################################################
# sort
############################################################
export GBMHOME=/net/midtier18/vol79/cphg-quinlan2/projects/gbm-seq-abounader
export STEPNAME=gbmsort
for sample in `echo $SAMPLES`;
do
export QSUB="qsub -q cphg -W group_list=CPHG -V -l select=1:mem=6000m:ncpus=8 -N $STEPNAME -m bea -M [email protected]";
echo "cd $GBMHOME; novosort -c 8 -t ./ -m 4G -f bam/$sample.novoalign.rg.bam > bam/$sample.novoalign.rg.sorted.bam" | $QSUB
done
############################################################
# index
############################################################
export GBMHOME=/net/midtier18/vol79/cphg-quinlan2/projects/gbm-seq-abounader
export STEPNAME=gbmidx
for sample in `echo $SAMPLES`;
do
export QSUB="qsub -q cphg -W group_list=CPHG -V -l select=1:mem=20000m:ncpus=1 -N $STEPNAME -m bea -M [email protected]";
echo "cd $GBMHOME; samtools index bam/$sample.novoalign.rg.sorted.bam" | $QSUB
done
############################################################
# Identify realignment targets.
############################################################
export GBMHOME=/net/midtier18/vol79/cphg-quinlan2/projects/gbm-seq-abounader
export STEPNAME=gbm-realtgts
export GENOME=/home/arq5x/cphg-home/shared/genomes/hg19/bwa/gatk/hg19_gatk.fa
export BIN=/home/arq5x/cphg-home/shared/bin/GenomeAnalysisTK-2.4-9-g532efad
for sample in `echo $SAMPLES`;
do
export QSUB="qsub -W group_list=cphg_arq5x -q arq5xlab -V -l select=1:mem=8000m:ncpus=3 -N $STEPNAME -m bea -M [email protected]";
echo "cd $GBMHOME; java -Xmx2g -jar $BIN/GenomeAnalysisTK.jar \
-T RealignerTargetCreator \
-R $GENOME \
-I bam/$sample.novoalign.rg.sorted.bam \
-o bam/$sample.novoalign.rg.sorted.bam.intervals" | $QSUB;
done
############################################################
# Use GATK to realign the BAMs at the suspicious targets.
############################################################
export GBMHOME=/net/midtier18/vol79/cphg-quinlan2/projects/gbm-seq-abounader
export STEPNAME=gbm-real
export GENOME=/home/arq5x/cphg-home/shared/genomes/hg19/bwa/gatk/hg19_gatk.fa
export BIN=/home/arq5x/cphg-home/shared/bin/GenomeAnalysisTK-2.4-9-g532efad
for sample in `echo $SAMPLES`;
do
export QSUB="qsub -W group_list=cphg_arq5x -q arq5xlab -V -l select=1:mem=16000m:ncpus=2 -N $STEPNAME -m bea -M [email protected]";
echo "cd $GBMHOME; java -Xmx2g -Djava.io.tmpdir=$GBMHOME/bam -jar $BIN/GenomeAnalysisTK.jar \
-T IndelRealigner \
--num_threads 1 \
-R $GENOME \
-targetIntervals bam/$sample.novoalign.rg.sorted.bam.intervals \
-I bam/$sample.novoalign.rg.sorted.bam \
-o bam/$sample.novoalign.rg.sorted.realigned.bam" | $QSUB;
done
############################################################
# Merge BAM files
############################################################
export GBMHOME=/net/midtier18/vol79/cphg-quinlan2/projects/gbm-seq-abounader
export STEPNAME=gbm-merge
export BIN=/home/arq5x/cphg-home/shared/bin/picard-tools-1.60
export QSUB="qsub -q cphg -W group_list=CPHG -V -l select=1:mem=8000m:ncpus=4 -N $STEPNAME -m bea -M [email protected]";
echo "cd $GBMHOME; java -Xmx2g -jar $BIN/MergeSamFiles.jar \
O=bam/all.novoalign.rg.sorted.bam \
I=bam/BLV4.novoalign.rg.sorted.bam \
I=bam/NCH411GBM_CD133high.novoalign.rg.sorted.bam \
I=bam/NCH411GBM_CD133low.novoalign.rg.sorted.bam \
I=bam/NCH537P54_CD133neg.novoalign.rg.sorted.bam \
I=bam/NCH537P54_CD133pos.novoalign.rg.sorted.bam \
I=bam/NCH620P55_CD133neg.novoalign.rg.sorted.bam \
I=bam/NCH620P55_CD133pos.novoalign.rg.sorted.bam \
I=bam/NCH644GBM_CD133high.novoalign.rg.sorted.bam \
I=bam/NCH644GBM_CD133low.novoalign.rg.sorted.bam \
I=bam/NCH7Md_P43_CD133neg.novoalign.rg.sorted.bam \
I=bam/NCH7Md_P43_CD133pos.novoalign.rg.sorted.bam \
I=bam/NPC-v.novoalign.rg.sorted.bam \
TMP_DIR=$RSHOME/bam \
VALIDATION_STRINGENCY=SILENT \
USE_THREADING=true" | $QSUB;
############################################################
# Counts.
############################################################
export GBMHOME=/net/midtier18/vol79/cphg-quinlan2/projects/gbm-seq-abounader
export STEPNAME=gbmsort
for sample in `echo $SAMPLES`;
do
export QSUB="qsub -W group_list=cphg_arq5x -q arq5xlab -V -l select=1:mem=10000m:ncpus=1 -N $STEPNAME -m bea -M [email protected]";
echo "cd $GBMHOME; bedtools-v2.18 intersect -c -b bam/$sample.novoalign.rg.sorted.realigned.bam -a annotations/hg19.1Kb.windows.bed -sorted -g annotations/hg19.genome > coverage/$sample.1Kb.coverage.bedg" | $QSUB
done
############################################################
# Call SNPs and INDELs with GATK W/O using BAQ.
############################################################
export GBMHOME=/net/midtier18/vol79/cphg-quinlan2/projects/gbm-seq-abounader
export STEPNAME=gbm-varbaq
export GENOME=/home/arq5x/cphg-home/shared/genomes/hg19/bwa/gatk/hg19_gatk.fa
export BIN=/home/arq5x/cphg-home/shared/bin/GenomeAnalysisTK-2.4-9-g532efad/
export QSUB="qsub -q cphg -W group_list=CPHG -V -l select=1:mem=16000m:ncpus=16 -N $STEPNAME -m bea -M [email protected]";
echo "cd $GBMHOME; java -Xmx2g -jar $BIN/GenomeAnalysisTK.jar \
-T UnifiedGenotyper \
-glm BOTH \
--num_threads 14 \
-baq OFF \
-R $GENOME \
-I bam/BLV4.novoalign.rg.sorted.realigned.bam \
-I bam/NCH411GBM_CD133high.novoalign.rg.sorted.realigned.bam \
-I bam/NCH411GBM_CD133low.novoalign.rg.sorted.realigned.bam \
-I bam/NCH537P54_CD133neg.novoalign.rg.sorted.realigned.bam \
-I bam/NCH537P54_CD133pos.novoalign.rg.sorted.realigned.bam \
-I bam/NCH620P55_CD133neg.novoalign.rg.sorted.realigned.bam \
-I bam/NCH620P55_CD133pos.novoalign.rg.sorted.realigned.bam \
-I bam/NCH644GBM_CD133high.novoalign.rg.sorted.realigned.bam \
-I bam/NCH644GBM_CD133low.novoalign.rg.sorted.realigned.bam \
-I bam/NCH7Md_P43_CD133neg.novoalign.rg.sorted.realigned.bam \
-I bam/NCH7Md_P43_CD133pos.novoalign.rg.sorted.realigned.bam \
-I bam/NPC-v.novoalign.rg.sorted.realigned.bam \
-o varCalling/all.raw.nobaq.vcf" | $QSUB;
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment