This repository has been archived by the owner on Aug 31, 2020. It is now read-only.
forked from rpolicastro/ChIPseq
-
Notifications
You must be signed in to change notification settings - Fork 1
/
main.sh
106 lines (81 loc) · 2.2 KB
/
main.sh
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
#!/bin/bash
cd $PBS_O_WORKDIR
source settings.conf
#########################################
## ChIP-seq Analysis of Multiple Samples
#########################################
## loading conda environment
## -------------------------
source activate chipseq-automation
## download files from SRA if necessary
## ------------------------------------
# creating sequence directory
mkdir -p $SEQDIR
# downloading from SRA if necessary
Rscript ${BASEDIR}/bin/getSRA.R \
--download $DOWNLOAD \
--samplesheet $SAMPLE_SHEET \
--seqdir $SEQDIR
## fastqc of reads
## ---------------
# creating directory to output fastqc results
mkdir -p ${OUTDIR}/fastqc
# saving fastq file names to variable
FASTQ_FILES=$(find ${SEQDIR} -name "*fastq")
# running fastqc
fastqc \
-t $CORES \
-o ${OUTDIR}/fastqc \
$FASTQ_FILES
## bowtie 2 read alignment
## -----------------------
# create directory for genomic index
mkdir -p ${OUTDIR}/genome/index
# create bowtie2 index
bowtie2-build \
-f --threads $CORES \
$GENOME_FASTA \
${OUTDIR}/genome/index/genome
# create directory to output alignments
mkdir -p ${OUTDIR}/aligned
# read alignment
Rscript ${BASEDIR}/bin/alignReads.R \
--outdir $OUTDIR \
--seqdir $SEQDIR \
--threads $CORES \
--samplesheet $SAMPLE_SHEET
# converting to coordinate sorted bam with index
for SAM in ${OUTDIR}/aligned/*sam; do
samtools sort \
-O BAM -@ $CORES \
-o ${OUTDIR}/aligned/$(basename $SAM .sam).bam \
$SAM
done
for BAM in ${OUTDIR}/aligned/*bam; do samtools index $BAM; done
## peak calling and annotation
## ---------------------------
# creating directory to output peaks
mkdir -p ${OUTDIR}/peaks
# calling peaks
Rscript ${BASEDIR}/bin/callPeaks.R \
--outdir $OUTDIR \
--threads $CORES \
--samplesheet $SAMPLE_SHEET \
--genomesize $GENOME_SIZE
# creating directory to output annotated peak files
mkdir -p ${OUTDIR}/annotated_peaks
# annotating peaks
Rscript ${BASEDIR}/bin/annoPeaks.R \
--outdir $OUTDIR \
--upstream $UPSTREAM \
--downstream $DOWNSTREAM \
--genomegtf $GENOME_GTF
## bam to bigwig
## -------------
# creating directory to output bigwigs
mkdir -p ${OUTDIR}/bigwigs
# converting bams to bigwigs
Rscript ${BASEDIR}/bin/bamToBigwig.R \
--outdir $OUTDIR \
--threads $CORES \
--samplesheet $SAMPLE_SHEET