-
Notifications
You must be signed in to change notification settings - Fork 3
/
Copy pathChIPSeqFPro.SR.hg19.sh
74 lines (52 loc) · 1.63 KB
/
ChIPSeqFPro.SR.hg19.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
#!/bin/bash
SECONDS=0
###fastqc quality control - requires fastqc installed and placed in PATH
ls -1 *fastq.gz > commands.1
sed -i 's/^/.\/FastQC\/fastqc /g' commands.1
source commands.1
mkdir FastQC_OUTPUT
mv *zip FastQC_OUTPUT
mv *html FastQC_OUTPUT
###mapping with BWA
files=(*fastq.gz)
for (( i=0; i<${#files[@]} ; i+=1 )) ; do
mkdir "${files[i]}.BWA"
done
GenomeDir='~/reference_genomes/hg19/'
GenomeFasta='~/reference_genomes/hg19/hg19.fa'
mkdir STATS
files=(*fastq.gz)
for (( i=0; i<${#files[@]} ; i+=1 )) ; do
echo $(pwd)/${files[i]}
Reads="$(pwd)/"${files[i]}""
Cores=$1
echo $Reads
cat >> commands.2.${files[i]}.tmp <<EOL
#!/bin/bash
echo Proccessing `pwd`: ${files[i]}
# run bwa mem
bwa mem -M -t $Cores $GenomeFasta $Reads > ${files[i]}.sam
# run samtools to convert sam to bam
samtools view -Sb ${files[i]}.sam > ${files[i]}.bam
# get stats on bam files
samtools flagstat ${files[i]}.bam > ${files[i]}.flagstats
mv ${files[i]}.flagstats ../STATS/
#convert bam to bigwig
cd ..
./bam2bigwig.sh ${files[i]}.BWA/${files[i]}.bam
#running MACS2 on bam
macs2 -t ${files[i]}.BWA/${files[i]}.bam -g hs -n ${files[i]}.MACS2
EOL
done
files=(*fastq.gz)
for (( i=0; i<${#files[@]} ; i+=1 )) ; do
sed -i "3i\\\tcd ${files[i]}.BWA" commands.2.${files[i]}.tmp
sed -i "4i\\\t# enter the correct folder" commands.2.${files[i]}.tmp
done
for (( i=0; i<${#files[@]} ; i+=1 )) ; do
sed -i "2i\ Reads=\"`pwd`/${files[i]} \"" commands.2.${files[i]}.tmp
done
for (( i=0; i<${#files[@]} ; i+=1 )) ; do
source commands.2.${files[i]}.tmp
done
rm *tmp