5. Filtering (bcftools ‐ two scenarios) - AreliDGS/EPP622-2026-SNPCalling GitHub Wiki

5. Filtering (bcftools - two scenarios)

This section documents the comparative evaluation between a gentle and a highly stringent variant filtration layout applied to the joint Saccharomyces cerevisiae raw dataset using bcftools filter.


Workspace and coding

# 1. Make sure you are in the project root directory
cd /lustre/isaac24/proj/UTK0505/agonza84/EPP622-2026-SNPCalling

# 2. Create and open the filtering script using vi
vi scripts/05_filtering.sh

# 3. Give execute permissions to the script
chmod +x scripts/05_filtering.sh

Execution script code

#!/bin/bash

# Exit immediately if any command or pipe fails
set -euo pipefail

# Create the folder to save the clean variant files
mkdir -p results/05_filtering

# --- First filtering (gentle) ---
# -e means "exclude": we throw away sites with low quality or low read support.
# QUAL<20 filter variants with a quality score below 20 
# INFO/DP<5 deletes sites with fewer than 5 overlapping reads across samples.
# -Oz saves the output as a compressed VCF file (.vcf.gz).
bcftools filter \
  -e 'QUAL<20 | INFO/DP<5' \
  -Oz \
  -o results/05_filtering/yeast_data.gentle.vcf.gz \
  results/04_variant_calling/yeast_data.raw.vcf.gz

# Run statistics 
bcftools stats results/05_filtering/yeast_data.gentle.vcf.gz > results/05_filtering/yeast_data.gentle.stats.txt


# --- Second filter (strict) ---
# We try keep only the highest-confidence mutations.
# QUAL<30 drops variants with quality under 30.
# INFO/DP<10 drops positions supported by fewer than 10 total reads.
# INFO/MQ<40 drops sites with a mapping quality score under 40 (from poor alignment).
bcftools filter \
  -e 'QUAL<30 | INFO/DP<10 | INFO/MQ<40' \
  -Oz \
  -o results/05_filtering/yeast_data.stringent.vcf.gz \
  results/04_variant_calling/yeast_data.raw.vcf.gz

# Run statistics
bcftools stats results/05_filtering/yeast_data.stringent.vcf.gz > results/05_filtering/yeast_data.stringent.stats.txt

Usage: 05_filtering.sh

# Launch the script in the background 
nohup ./scripts/05_filtering.sh > results/05_filtering.log 2>&1 &

# See the execution progress in real time
tail -f results/05_filtering.log

Results

  • bcftools stats report path and name for gentle filtering: results/05_filtering/yeast_data.gentle.stats.txt
  • bcftools stats report path and name for stringent filtering: results/05_filtering/yeast_data.stringent.stats.txt

Data extraction commands

To verify the exact number of SNPs surviving both filtering approaches without parsing bias, the following explicit diagnostic queries were executed directly in the terminal:

# Extract SNP count from the gentle filtration summary report
grep "number of SNPs:" results/05_filtering/yeast_data.gentle.stats.txt

# Extract SNP count from the stringent filtration summary report
grep "number of SNPs:" results/05_filtering/yeast_data.stringent.stats.txt

VCF variant progression summary table

VCF Dataset Source File Tracker Number of SNPs
raw yeast_data.raw.vcf.gz 4,920
gentle yeast_data.gentle.vcf.gz 4,869
stringent yeast_data.stringent.vcf.gz 4,774