Filtering & Frequency Calculations

Step 7 takes the joint VCF that step 6 produced and turns it into the frequency tables you actually read. Between the two sit the filters that decide which of the called sites survive.

The order is fixed — stages 5 to 9 of The Filter Chain: normalize to the major allele, drop alleles too few samples support, drop sites too thin or too poorly called, split SNPs from indels, then convert allelic depths to frequencies. Only two of those stages have parameters, and both are analysis-affecting.

Parameter	What it does
poolSize, ploidy	Set the smallest allele frequency worth believing, which is what the false-positive filter tests against. Both can be set per sample or per run
filterFalsePositives.sampleThreshold	Removes alternate alleles without support across enough samples. The filter that makes the pipeline pool-aware, and the one most worth understanding before changing anything
vcffilter.minDP, vcffilter.minQUAL	Remove whole sites that are too shallow in any one sample, or too poorly called
vcf.fileName	Removes nothing — it names the files this step writes

filterFalsePositives.sensitivity is computed from poolSize and ploidy rather than set, and can be overridden per sample in metadata.csv — see Metadata for the per-sample form.

For the stage-by-stage account with counts, read The Filter Chain.

poolSize and ploidy

poolSize = 50    // individuals in one pool
ploidy   = 2     // copies of the genome per individual

These two are the pipeline's most consequential settings. They are not used to model anything — they set the minimum credible allele frequency:

\[S = \frac{1}{2 \times \text{ploidy} \times \text{poolSize}}\]

An alternate allele must reach S in a sample for that sample to count as supporting it.

poolSize is per pool, not per run. With eight pools of 50 individuals each, it is 50, not 400.

The factor of two is deliberate. A pool of 50 diploids holds 100 chromosomes, so one chromosome is frequency 0.01. S comes out at 0.005 — half of that — so a genuine singleton clears the threshold with margin rather than sitting exactly on it. The observed fraction of a singleton is itself noisy, and a threshold placed exactly at the expected value would reject half of them.

poolSize	ploidy	One chromosome is	S
10	2	0.0500	0.0250
25	2	0.0200	0.0100
50	2	0.0100	0.0050
100	2	0.0050	0.0025
200	2	0.0025	0.00125
50	1	0.0200	0.0100
50	4	0.0050	0.0025

If your pools differ in size, give each one its own param_poolSize in metadata.csv and every column is judged at its own threshold. The global poolSize here stays as the default for pools that do not state one. See Pool size belongs to the pool.

ploidy applies to a whole run. It describes the sequence rather than the sample, so a diploid animal's haploid mitochondrial genome belongs in a separate run with its own value — see Sequences with a different ploidy.

sampleThreshold

filterFalsePositives {
    sampleThreshold = 0.2
}

The fraction of samples that must independently support an allele at frequency S or above for it to be kept. This is what makes the filter a corroboration test rather than a frequency cutoff, and it is what allows S to be set so low without drowning in sequencing error: errors do not recur at the same position across independent libraries.

The comparison is count >= n_samples × sampleThreshold, so the effective requirement is the next whole number up:

Samples	M at 0.2	Must support
1	0.2	1
4	0.8	1
5	1.0	1
6	1.2	2
8	1.6	2
10	2.0	2
12	2.4	3
20	4.0	4

The default removes population-private alleles

With eight samples, an allele found in one pool is discarded regardless of how frequent it is there — one does not reach two. If private or population-specific variation is the subject of your study, this default is wrong for you.

Value	Behavior	Suits
0.2 (default)	Needs roughly a fifth of samples	Shared, evolving variation across comparable pools
Low enough to require 1 sample	Any single pool can carry an allele	Private variants, small sample counts, discovery runs
Higher, e.g. 0.5	Needs half the samples	Conservative core-variant sets; high-noise data

Remember that the denominator is the number of VCF columns, which is the number of distinct RG_Sample values — not the number of FASTQ pairs. Merging replicates by sharing an RG_Sample changes this threshold as a side effect (see Metadata).

Depth and quality

vcffilter {
    minDP   = 20
    minQUAL = 30
}

Both are site-level filters, applied as two commands in sequence:

bcftools view -e "FMT/DP<20" -Ov -o <name>_dp.vcf <input>
vcftools --vcf <name>_dp.vcf --minQ 30 --recode --recode-INFO-all --out <name>_dq

minDP → bcftools view -e "FMT/DP<N"

Removes a site if any sample falls below the depth. Read the negation carefully: the site survives only when every sample meets the floor.

minQUAL → vcftools --minQ

Removes sites whose QUAL falls below the value.

The weakest library sets the threshold for every site

Because the test is "any sample below minDP", one under-sequenced pool removes sites for all of them. Three pools at depths 50/15/55, 60/12/28 and 10/12/20, with minDP = 20:

minDP = 20  ->  0 of 3 sites kept   (poolB is 15/12/12, so it fails everywhere)
minDP = 12  ->  2 of 3 sites kept

Check Output/Reports/Coverage/ for your weakest sample before choosing a value, and after the first run compare the site count in <name>.vcf against the distinct positions that reached the frequency tables. A near-total wipeout is this filter, not a broken pipeline.

If "every sample" is too strict for your design, the alternatives are one-line swaps in 7_vcf2freq.nf. Tested against the depths above at minDP = 20:

Expression	Semantics	Sites kept
-e "FMT/DP<20" (current)	Every sample must pass	0 of 3
-i "COUNT(FMT/DP>=20)>=2"	At least 2 samples pass	2 of 3
-i "MEAN(FMT/DP)>=20"	Mean depth across samples	2 of 3
-i "INFO/DP>=20"	Cohort total depth	3 of 3
-i "FMT/DP>=20"	Any one sample passes — not a depth floor	3 of 3

The last row is worth noting as a trap: -i "FMT/DP>=20" reads like the obvious inverse of the current expression and is not, because bcftools evaluates a FORMAT condition per sample and keeps the site if it holds for any of them.

Genotype-level filtering is not available here

A per-sample depth floor — blanking one pool's frequency while keeping the row — cannot be done in the VCF, because vcftools and bcftools both express a genotype-level verdict by rewriting GT, and GT is set to ./. throughout by major-allele normalization (why). Frequency conversion reads AD. If you need per-sample blanking rather than whole-site removal, it has to happen in bin/depth2freq.awk, which already computes each sample's depth as the denominator.

Output naming

vcf {
    fileName = 'Test'
}

Base name for every VCF and frequency table. With the default, a finished run leaves Test.vcf, Test_snp_freq.tsv, Test_indel_freq.tsv, and Test_annotated.vcf if annotation ran. Worth setting to something descriptive — it is the name your results carry from here on.

Changing it after a run does not rename anything; it makes the resume checks look for files that do not exist, and the whole VCF and frequency branch runs again alongside the old files.