Command used:
modkit extract calls pure_duplex_mapped.bam duplex_read_mods.tsv.gz
--reference /path/to/GRCm39.primary_assembly.genome.fa
--mapped-only
--filter-threshold C:0.7
--mod-threshold m:0.8
--mod-threshold h:0.9
--bgzf
-t 18
--io-threads 4
--out-threads 4
--log modkit_extract_stats.log
According to the documented threshold logic (see "Threshold examples" in the modkit book, which uses this exact command as its worked example), a call with call_code = - (canonical) and call_prob = 0.36 should be flagged as failed, since 0.36 is below the specified --filter-threshold C:0.7.
However, in my output I'm seeing rows like this where fail = false despite call_prob being well below the threshold:
call_code = -
call_prob = 0.36
canonical_base = G
modified_primary_base = C
fail = false
The key detail is that this row comes from a read aligned to the negative strand. Because the BAM SEQ field is reverse-complemented for reverse-strand alignments, canonical_base (which reflects the literal SEQ-orientation base) shows G, while modified_primary_base (the semantically correct primary base carrying the modification, accounting for the MM tag's strand indicator) correctly shows C.
Command used:
modkit extract calls pure_duplex_mapped.bam duplex_read_mods.tsv.gz
--reference /path/to/GRCm39.primary_assembly.genome.fa
--mapped-only
--filter-threshold C:0.7
--mod-threshold m:0.8
--mod-threshold h:0.9
--bgzf
-t 18
--io-threads 4
--out-threads 4
--log modkit_extract_stats.log
According to the documented threshold logic (see "Threshold examples" in the modkit book, which uses this exact command as its worked example), a call with call_code = - (canonical) and call_prob = 0.36 should be flagged as failed, since 0.36 is below the specified --filter-threshold C:0.7.
However, in my output I'm seeing rows like this where fail = false despite call_prob being well below the threshold:
call_code = -
call_prob = 0.36
canonical_base = G
modified_primary_base = C
fail = false
The key detail is that this row comes from a read aligned to the negative strand. Because the BAM SEQ field is reverse-complemented for reverse-strand alignments, canonical_base (which reflects the literal SEQ-orientation base) shows G, while modified_primary_base (the semantically correct primary base carrying the modification, accounting for the MM tag's strand indicator) correctly shows C.