Couldn't anchor the Sanger trace in the selected reference genome. error in tracy but not in Indigo

[jacorvar]

Hi,

I'm experiencing a very similar issue to #34 .

We've sequenced a chunk of the human MBL2 gene of ~250 nt. However, the machine sequences almost 1000 nucleotides; therefore, most of the sequence in the ab1 file is just rubbish.
For that reason I've decided to set -q 50 -u 750, so that the first 50 low-quality bases and the last 750 (false) bases are excluded.
I tried first with the whole GATK GRCh38 genome and got the error Couldn't anchor the Sanger trace in the selected reference genome. when running both Indigo (setting left and right trim sizes to 50 and 750, respectively) and tracy in the command line as follows:

tracy decompose -o forward -r Homo_sapiens_assembly38.fasta.gz -q 50 -u 750 MF-102_MBL2.ab1
[2023-May-04 12:10:37] tracy decompose -o forward -a homo_sapiens -r Homo_sapiens_assembly38.fasta.gz -q 50 -u 750 MF-102_2MBL2.ab1 
[2023-May-04 12:10:37] Load ab1 file
[2023-May-04 12:10:37] Find Reference Match
[2023-May-04 12:10:37] Load FM-Index
Couldn't anchor the Sanger trace in the selected reference genome.

As you pointed out here, that issue could be circumvented using a shorter sequence as a reference file.

Then I downloaded and indexed the fasta file for the MBL2 gene and repeated the process with the same parameters. Although it works well now with Indigo, tracy still fails with the same error message in the command line.

I'm using tracy v0.7.5 singularity container in CentOS 7.9.

[tobiasrausch]

Indigo just runs tracy in the backend:

tracy decompose -v -r MBL2.gene.slice.fa -q 50 -u 750 -p 0.3 MF-102_MBL2.ab1

[jacorvar]

OK, that works now

BTW, given that I must provide the gene sequence instead of the whole genome/chromosome sequence in order to anchor the trace, which strategy would you suggest to annotate (adding RS identifiers) the bcf provided its positions are related to the gene and not to the genome?
I know I could do that using awk and some scripting, but I'd like to know if there's a more elegant and (possibly) better alternative.

Thanks a lot

[tobiasrausch]

I am afraid we don't have a better solution for the annotation with RS identifiers with a gene sequence as input.