Hi everyone,

I encounter an unexpected behavior of seqkit amplicon: seqkit amplicon -p returns a different result than seqkit amplicon -F -R.

seq.fasta

>seq1
ATGCGCTATATATATTTT
>seq2
GGGGAGTGTGTGTGTTTT

These commands produce the expected output:

seqkit amplicon -F "ATGC" -R "AAAA" seq.fasta --bed
# seq1    0       18      .       0       +       ATGCGCTATATATATTTT

seqkit amplicon -F "GGGG" -R "AAAA" seq.fasta --bed
# seq2    0       18      .       0       +       GGGGAGTGTGTGTGTTTT

However, if we introduce a file to store the primers, we do not get the same results:
arr.txt

arr1    ATGC    AAAA
arr2    GGGG    AAAA

seqkit amplicon -p arr.txt seq.fasta --bed
# seq2    0       18      arr1    0       +       GGGGAGTGTGTGTGTTTT
# seq2    0       18      arr2    0       +       GGGGAGTGTGTGTGTTTT

It looks like the sequence seq2 is misreported as a valid amplicon of arrangement arr1; whereas I would expect seq1 to be reported, as in the -F -R example. Tested on 2.7.0. The same behavior may be observed without the --bed option.

Am I missing something?

HTH