-
Notifications
You must be signed in to change notification settings - Fork 0
/
bcftools_filtration.sh
executable file
·54 lines (42 loc) · 1.68 KB
/
bcftools_filtration.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
#!/bin/bash
vcf=''
print_usage() {
echo "Creating files for genome vizualization."
echo "Usage: '-i' your VCF file"
}
while getopts 'i:' flag; do
case "${flag}" in
i) vcf="${OPTARG}" ;;
*) print_usage
exit 1 ;;
esac
done
# paths
# export PATH=${PATH}:/mnt/tank/scratch/atomarovsky/tools/bcftools-1.12/bin:/mnt/tank/scratch/atomarovsky/tools/parallel-20210622/bin:/mnt/tank/scratch/atomarovsky/tools/vcftools-0.1.16/bin
# vcf prefix
pref=${vcf%.*.*}
# script
mkdir bcftools_filtration;
cd bcftools_filtration/;
echo "bcftool filter..."
bcftools filter -S . -o ${pref}.filt.vcf --exclude 'QUAL < 20.0 || (FORMAT/SP > 60.0 | FORMAT/DP < 5.0 | FORMAT/GQ < 20.0)' ../${vcf}
fvcf=${pref}.filt.vcf
pref=${fvcf%.*}
touch tmp.list
echo "bcftools query..."
for SAMPLE in `bcftools query -l ${fvcf}`; do
echo "Handling ${SAMPLE}..." ;
vcf-subset --exclude-ref -c ${SAMPLE} ${fvcf} > ${pref}.${SAMPLE}.vcf ;
echo ${pref}.${SAMPLE}.vcf >> tmp.list ;
done
for file in $(cat tmp.list); do
echo "bcftools filter...";
bcftools filter -i 'TYPE="indel"' -O v ${file} > ${file%.*}.indels.vcf && bcftools filter -i 'TYPE="snp"' -O v ${file} > ${file%.*}.SNPs.vcf;
echo "get hyterozygous variants..."
bcftools filter -i 'FMT/GT = "het"' -O v ${file%.*}.indels.vcf > ${file%.*}.indels.hetero.vcf && bcftools filter -i 'FMT/GT = "het"' -O v ${file%.*}.SNPs.vcf > ${file%.*}.SNPs.hetero.vcf
echo "get homozygous variants..."
bcftools filter -i 'FMT/GT = "hom"' -O v ${file%.*}.indels.vcf > ${file%.*}.indels.homo.vcf && bcftools filter -i 'FMT/GT = "hom"' -O v ${file%.*}.SNPs.vcf > ${file%.*}.SNPs.homo.vcf
done
rm tmp.list
# Draw densities of heterozygous files
# $TOOLS/bashare/draw_densities_of_hetero.sh