This section documents the comparative accuracy evaluating our raw, gently filtered, and stringently filtered variant calling datasets against the ground truth set using the evaluate_snps_v4.py.

In genomics, the ground truth represents a high-confidence, reference catalog of true markers physically present within the target Saccharomyces cerevisiae genomes. By cross-referencing my predictions against vs verified baseline, I can quantify the accuracy and errors during this processing chain.

It conducts a position-by-position comparison between our VCF calls and the reference file. By calculating genomic site into standard statistical validation matrices:

• True Positives (TP): Variant positions present in the organism that our pipeline successfully identified.
• False Positives (FP): Artifacts or sequencing errors that our pipeline misidentified as variants, acting as false alarms.
• False Negatives (FN): Real biological mutations present in the organism that our pipeline missed due to insufficunet depth or quality thresholds.

#!/bin/bash
# scripts/06_evaluation.sh
set -euo pipefail

# Create the folder to save outputs
mkdir -p results/06_evaluation

# Referene paths
TRUTH_VCF="/lustre/isaac24/proj/UTK0505/test2/check/truth_snps.vcf"
EVAL_SCRIPT="data/reference/evaluate_snps_v4.py"

# --- Phase 1: Raw set Cross-Comparison ---
echo "--- Phase 1: Raw set Cross-Comparison ---"
python3 $EVAL_SCRIPT \
    --truth "$TRUTH_VCF" \
    --calls results/04_variant_calling/yeast_data.raw.vcf.gz \
    --all-samples \
    --tag raw \
    --outdir results/06_evaluation/raw_eval

# --- Phase 2: Gentle set Cross-Comparison ---
echo "--- Phase 2: Gentle set Cross-Comparison ---"
python3 $EVAL_SCRIPT \
    --truth "$TRUTH_VCF" \
    --calls results/05_filtering/yeast_data.gentle.vcf.gz \
    --all-samples \
    --tag gentle \
    --outdir results/06_evaluation/gentle_eval

# --- Phase 3: Strict Set Cross-Comparison ---
echo "--- Phase 3: Strict Set Cross-Comparison ---"
python3 $EVAL_SCRIPT \
    --truth "$TRUTH_VCF" \
    --calls results/05_filtering/yeast_data.stringent.vcf.gz \
    --all-samples \
    --tag strict \
    --outdir results/06_evaluation/strict_eval

# --- Phase 4: Directory final structure ---
echo "--- Phase 4: Directory final structure ---"
tree -LF 2 . > results/06_evaluation/directory_structure.txt

echo "Step 6 evaluation commands completed successfully."

# Launch the master python evaluation process in the background
nohup ./scripts/06_evaluation.sh > results/06_evaluation.log 2>&1 &

# Follow cross-comparison status outputs in real time
tail -f results/06_evaluation.log

├── data/
│   ├── raw/
│   └── reference/
├── README.md
├── results/
│   ├── 01_fastqc_raw/
│   ├── 01_multiqc_raw_report_data/
│   ├── 01_multiqc_raw_report.html
│   ├── 02_fastp_new.log
│   ├── 02_fastqc_trimmed/
│   ├── 02_multiqc_trimmed_report_1.html
│   ├── 02_multiqc_trimmed_report_data/
│   ├── 02_multiqc_trimmed_report_data_1/
│   ├── 02_multiqc_trimmed_report.html
│   ├── 02_trimmed/
│   ├── 03_alignment/
│   ├── 03_alignment.log
│   ├── 03_multiqc_alignment_report_data/
│   ├── 03_multiqc_alignment_report.html
│   ├── 04_variant_calling/
│   ├── 05_filtering/
│   ├── 06_evaluation/
│   ├── 06_evaluation.log
│   ├── environment_lock.yml
│   ├── evaluate_snps_v4.py -> /lustre/isaac24/proj/UTK0505/test2/check/evaluate_snps_v4.py*
│   └── truth_snps.vcf -> /lustre/isaac24/proj/UTK0505/test2/check/truth_snps.vcf*
└── scripts/
    ├── 01_fastqc_multiqc.sh*
    ├── 02_read_trimming.sh*
    ├── 03_alignment.sh*
    ├── 04_variant_calling.sh*
    ├── 05_filtering.sh*
    └── 06_evaluation.sh*

16 directories, 17 files

The accuracy scorecard table below outlines the TP, FP, FN, precision, recall, and F1-Scores calculated across all sample genotypes:

To locate the evaluation files, isolate suspicious variants, and identify target file dependencies for JBrowse deployment, the following terminal commands were executed on the cluster:

# 1. Audit output log matrices to find non-empty evaluation summaries
ls -lh results/06_evaluation/raw_eval/

# 2. Extract real genomic coordinates from Sample B False Positive records
head -n 3 results/06_evaluation/raw_eval/Sample_B_raw_false_positives.tsv

# 3. Extract real genomic coordinates from Sample B False Negative records
head -n 3 results/06_evaluation/raw_eval/Sample_B_raw_false_negatives.tsv

# 4. Uncover exact filename for the genomic reference file assembly
ls -l data/reference/*fast*

How did samples differ from each other in read count versus in read quality? Why does each matter for SNP calling?

The biometric profiles extracted across my samples indicated a difference between absolute read depth volume and single-nucleotide base calling accuracy. Samples A, C, and D maintained optimal, robust sequence count volumes, and presenting a sufficient raw molecular depth to cover the S. cerevisiae genome. Conversely, Sample B suffered from a low-coverage constraint, despite a good quality markers. Sample D presented the inverse defect: deep sequence volume compromised by uniform quality capping. I think that both features are highly important for downstream variant calling but operate under different paradigms. Sufficient read count volume is needed to generate enough overlapping fragments to catch both ALT at heterozygous positions. Moreover, insufficient depth, the statistic threshold fails, leading to a drop in sensitivity. On the other hand, high read quality is require to observe true polymorphic base adjustments from automated sequencing errors or machine noise; low quality thresholds could lead to discard true variants or introduce false positive artifacts.

What are the relative merits and drawbacks of filtering? Interpret the above chart in terms of the original quality issues, and how those yielded different outcomes. When is filtering a good idea vs when it is it a bad idea? Think specifically about false positives vs false negatives.

Bioinformatic variant filtration works between precision (the called mutations) and recall (the sensitivity of detection). The main aim of stringency filtering is to be able of delete the FP, filtering out the noise, multi-mapping read confusion, and sequencing errors, which ensures that every called variant accurate. However, we also have to consider the FN, where genuine biological variants are accidentally discarded because their local depth or quality scores fall slightly below the hard numerical cutoffs. In my results, the precision scores=1 across Samples A, B, and C demonstrate that our bcftools joint calling configuration was already robust, leaving on zero false positives to clean up. In low-coverage scenarios like Sample B, a strict thresholds amplify data missingness and destroy biological discovery.

3. Visual verification of anomalies Using JBrowse

I inspected locally via the JBrowse 2 alignment track viewer using IndexedFastaAdapter streaming.

• Sample B - locus 1 (NC_001136.10 : 139214): This site represents a computational artifact called with a low Phred score of 3.4 and a shallow local mapping depth of 6 reads. Visual inspection in JBrowse 2 we can see that 3 reads randomly showed an "A" instead of the reference "T". Since the raw data didn't have any filters, bcftools confused this machine noise with a real SNP.

Figure 1. Sample B raw false positive at coordinate NC_001136.10:139214 showing machine base-calling noise.

• Sample C - locus 2 (NC_001137.3 : 42): This false positive is an artifact where the software called a mutation based on a single read (1X depth) with a poor quality score of 3.2. JBrowse 2 confirms that this position contains machine noise/error or an unclipped adapter remnant. By running our filtration process we removed the low-confidence site of this sample.

*Figure 2. Sample C raw false positive at coordinate NC_001137.3:42 caused by unclipped adapter artifacts.

• Sample B - locus 1 (NC_001133.9 : 8027): This is a real mutation that the pipeline completely missed because sample B suffered from a massive shortage of data. When checking this coordinate in JBrowse 2, the alignment track is empty, showing 0X coverage. The sequencer simply never passed through this spot, making the real mutation not visible.

Figure 3. Sample B filtered false negative at coordinate NC_001133.9:8027 showing complete data starvation and 0x coverage.

• Sample C - locus 2 (NC_001135.5 : 221): This real mutation became a false negative because caused by over-filtering. Our strict thresholds were so aggressive that the positioin was discarded.

figure 4. sample c filtered false negative at coordinate NC_001135.5:221 dropped due to strict quality thresholds.