variant calling and filtering by family for melissa crosses

Post date: Apr 15, 2016 6:58:49 PM

Note that for some families multiple fathers are possible, my notes from before suggest that these are the actual fathers:

F 200, M = 202 (this one is not certain)

F 26, M = 28

F 257, M = 258

F 272, M = 273

samtools (1.2) and bcftools (1.3) were used for variant calling for each family. Here is the shell script (all families included):

#!/bin/sh

#SBATCH --time=48:00:00

#SBATCH --nodes=1

#SBATCH --ntasks=8

#SBATCH --account=gompert

#SBATCH --partition=kingspeak

#SBATCH --job-name=samtools

#SBATCH --mail-type=FAIL

#SBATCH --mail-user=zach.gompert@usu.edu

echo ------------------------------------------------------

echo SLURM: job identifier is $SLURM_JOBID

echo SLURM: job name is $SLURM_JOB_NAME

echo ------------------------------------------------------

module load samtools

cd /uufs/chpc.utah.edu/common/home/u6000989/data/lycaeides/melissa_mappingfams/alignments/

samtools mpileup -b fam26bams -C 50 -d 250 -f /uufs/chpc.utah.edu/common/home/u6000989/data/lycaeides/melissa_genome/final.assembly.fasta -q 20 -Q 20 -I -u -g -t DP,DPR -o fam26.bcf

bcftools call -v -c -p 0.01 -P 0.001 -O v -o fam26.vcf fam26.bcf

samtools mpileup -b fam123bams -C 50 -d 250 -f /uufs/chpc.utah.edu/common/home/u6000989/data/lycaeides/melissa_genome/final.assembly.fasta -q 20 -Q 20 -I -u -g -t DP,DPR -o fam123.bcf

bcftools call -v -c -p 0.01 -P 0.001 -O v -o fam123.vcf fam123.bcf

samtools mpileup -b fam143bams -C 50 -d 250 -f /uufs/chpc.utah.edu/common/home/u6000989/data/lycaeides/melissa_genome/final.assembly.fasta -q 20 -Q 20 -I -u -g -t DP,DPR -o fam143.bcf

bcftools call -v -c -p 0.01 -P 0.001 -O v -o fam143.vcf fam143.bcf

samtools mpileup -b fam181bams -C 50 -d 250 -f /uufs/chpc.utah.edu/common/home/u6000989/data/lycaeides/melissa_genome/final.assembly.fasta -q 20 -Q 20 -I -u -g -t DP,DPR -o fam181.bcf

bcftools call -v -c -p 0.01 -P 0.001 -O v -o fam181.vcf fam181.bcf

samtools mpileup -b fam200bams -C 50 -d 250 -f /uufs/chpc.utah.edu/common/home/u6000989/data/lycaeides/melissa_genome/final.assembly.fasta -q 20 -Q 20 -I -u -g -t DP,DPR -o fam200.bcf

bcftools call -v -c -p 0.01 -P 0.001 -O v -o fam200.vcf fam200.bcf

samtools mpileup -b fam257bams -C 50 -d 250 -f /uufs/chpc.utah.edu/common/home/u6000989/data/lycaeides/melissa_genome/final.assembly.fasta -q 20 -Q 20 -I -u -g -t DP,DPR -o fam257.bcf

bcftools call -v -c -p 0.01 -P 0.001 -O v -o fam257.vcf fam257.bcf

samtools mpileup -b fam272bams -C 50 -d 250 -f /uufs/chpc.utah.edu/common/home/u6000989/data/lycaeides/melissa_genome/final.assembly.fasta -q 20 -Q 20 -I -u -g -t DP,DPR -o fam272.bcf

bcftools call -v -c -p 0.01 -P 0.001 -O v -o fam272.vcf fam272.bcf

samtools mpileup -b fam297bams -C 50 -d 250 -f /uufs/chpc.utah.edu/common/home/u6000989/data/lycaeides/melissa_genome/final.assembly.fasta -q 20 -Q 20 -I -u -g -t DP,DPR -o fam297.bcf

bcftools call -v -c -p 0.01 -P 0.001 -O v -o fam297.vcf fam297.bcf

Number of variants before filtering:

fam123.vcf:36460

fam143.vcf:66579

fam181.vcf:42990

fam200.vcf:87879

fam257.vcf:36261

fam26.vcf:46174

fam272.vcf:77487

fam297.vcf:36250

Filtering

Variants were first filtered using my vcfFilter.pl script, with the following options: no more than one alternative allele (relative to ref.), minimum coverage of 2X (actual number of reads adjusted by family), 10+ alternate reads, no fixed for alt. allele, minimum mapping quality of 30, no more than 4 individuals with missing data, and a minimum MAF of 0.05 (5%).

Here are sample sizes and numbers of filtered reads by family: