-
Notifications
You must be signed in to change notification settings - Fork 1
/
individual-mapping-freeb-snp-calling.sh
104 lines (63 loc) · 2.65 KB
/
individual-mapping-freeb-snp-calling.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
#!/bin/bash --login
#PBS -l walltime=60:00:00
#PBS -l mem=200gb
#PBS -m abe
#PBS -N elk-mapping
#PBS -A ged
cd ${PBS_O_WORKDIR}
source hpcc.modules
# index reference genome
bwa index final.contigs.fa
### ELK 2919 ###
## MAPPING ##
# map raw reads from 2919 to reference genome
bwa mem final.contigs.fa ../2919-INSERT-FILE-NAME-HERE.gz > 2919-elk-reads.sam
# index reference genome with samtools
samtools faidx final.contigs.fa
# convert sam to bam
samtools import contigs.fa.fai 2919-elk-reads.sam 2919-elk-reads.bam
# sort the bam file
samtools sort 2919-elk-reads.bam 2919-elk-reads.sorted
# index the sorted bam file
samtools index 2919-elk-reads.sorted.bam
## SNP CALLING ##
freebayes -f final.contigs.fa --ploidy 2 --min-coverage 10 --no-mnps --no-complex --min-alternate-count 3 2919-elk-reads.sorted.bam | vcffilter -f "QUAL > 20" > 2919-elk-snps.vcf
### ELK 2920 ###
## MAPPING ##
# map raw reads from 2920 to reference genome
bwa mem final.contigs.fa ../2920-INSERT-FILE-NAME-HERE.gz > 2920-elk-reads.sam
# convert sam to bam
samtools import contigs.fa.fai 2920-elk-reads.sam 2920-elk-reads.bam
# sort the bam file
samtools sort 2920-elk-reads.bam 2920-elk-reads.sorted
# index the sorted bam file
samtools index 2920-elk-reads.sorted.bam
## SNP CALLING ##
freebayes -f final.contigs.fa --ploidy 2 --min-coverage 10 --no-mnps --no-complex --min-alternate-count 3 2920-elk-reads.sorted.bam | vcffilter -f "QUAL > 20" > 2919-elk-snps.vcf
### ELK 1339 ###
## MAPPING ##
# map raw reads from 1339 to reference genome
bwa mem final.contigs.fa ../1339-INSERT-FILE-NAME-HERE.gz > 1339-elk-reads.sam
# convert sam to bam
samtools import contigs.fa.fai 1339-elk-reads.sam 1339-elk-reads.bam
# sort the bam file
samtools sort 1339-elk-reads.bam 1339-elk-reads.sorted
# index the sorted bam file
samtools index 1339-elk-reads.sorted.bam
## SNP CALLING ##
freebayes -f final.contigs.fa --ploidy 2 --min-coverage 10 --no-mnps --no-complex --min-alternate-count 3 1339-elk-reads.sorted.bam | vcffilter -f "QUAL > 20" > 1339-elk-snps.vcf
### ELK 1351 ###
## MAPPING ##
# map raw reads from 1351 to reference genome
bwa mem final.contigs.fa ../1351-INSERT-FILE-NAME-HERE.gz > 1351-elk-reads.sam
# convert sam to bam
samtools import contigs.fa.fai 1351-elk-reads.sam 1351-elk-reads.bam
# sort the bam file
samtools sort 1351-elk-reads.bam 1351-elk-reads.sorted
# index the sorted bam file
samtools index 1351-elk-reads.sorted.bam
## SNP CALLING ##
freebayes -f final.contigs.fa --ploidy 2 --min-coverage 10 --no-mnps --no-complex --min-alternate-count 3 1351-elk-reads.sorted.bam | vcffilter -f "QUAL > 20" > 1351-elk-snps.vcf
cat ${PBS_NODEFILE}
env | grep PBS
qstat -f ${PBS_JOBID}