lower base quality and more indels than actual data

[npdungca]

moved from Psy-Fer/buttery-eel#76 (comment)

Hi, Hasindu! Here you go:

I generated synthetic data using squigulator then basecalled with buttery-eel. There seems to be much more indels and lower base quality scores with the generated synthetic data than the actual data. Below is an IGV screenshot (upper panel: squiggulator+buttery-eel data; bottom panel:actual amplicon data sequenced in a R10 flowcell, lib prep kit NBD114, basecalled with SUP)

Here are the commands used to generate the synthetic data:

config=dna_r10.4.1_e8.2_400bps_sup.cfg

#create artificial datasets
time $squigulator -x dna-r10-prom -f ${d} -t 8 -r ${r} -q $outdir/${i}_${d}"x"_${r}ideal_${n}.fasta \
--bps 400 --ont-friendly=yes $ref/${i}.fasta -o $datadir/${i}_${d}"x"_${r}_${n}.blow5

#basecall
time buttery-eel -g $basecaller --config $config --device cuda:1 -i $datadir/${i}_${d}"x"_${r}_${n}.blow5 -o $outdir/${i}_${d}"x"_${r}_buttery-eel_${n}.fastq \
--port auto --use_tcp --dorado_download_path $dorado_download_path

The mean baseQs are 13.3 for the synthetic data and 34.8 for the actual data.

Thank you in advance for your help.

[hasindu2008]

Hi @npdungca

This is likely to be because the r10 pore model in Squigulator (which is derived from ONT's R10 model) is not of the greatest quality.

How did you generate $outdir/${i}_${d}"x"_${r}ideal_${n}.fasta that you are using? Is that generated from the CMVstrains.fasta? If you can share this CMVstrains.fasta (here or via email), I can see if adjusting any squigulator parameters would have an impact.

[npdungca]

Hi, Hasindu!

It's generated by squigulator. It has the reads without the errors. Here's the link for the ideal.fasta and the reference sequence from which it was derived.
https://udemontreal-my.sharepoint.com/:f:/g/personal/nina_theresa_dungca_umontreal_ca/Erj-ekYDPrZFoSRG4FjDMb4B4KyM5PRj9joUUaGj0KKtaw?e=2wSdcX

[hasindu2008]

Thank you. Just did a quick check.

If you add --dwell-std 0.0, the mean qscore improves a tiny bit to around 13.7.

Then, checking how --amp-noise factor can affect using the following snippet

for i in 0.9 0.95 1.0 1.1; 
do 
squigulator -x dna-r10-prom -f 20 --ont-friendly=yes -o $i.blow5 --dwell-std 0.0 --amp-noise $i CMV.fasta; 
/install/buttery-eel-0.4.2+dorado7.2.13/scripts/eel --config dna_r10.4.1_e8.2_400bps_sup.cfg  --device 'cuda:0' -i $i.blow5 -o  $i.fastq
done

--amp-noise 0.95 improves another tiny bit to 13.9 qscore.

Interestingly, using the v5.0 sup models in slow5-dorado gives an average qscore bump to 19.1. But your real data also might see a similar bump in qscore with these models.

Command used:

/install/slow5-dorado-0.8.3/bin/slow5-dorado basecaller /install/slow5-dorado-0.8.3/models/dna_r10.4.1_e8.2_400bps_sup@v5.0.0 a.blow5 --emit-fastq -x cuda:0 > a.fastq

But I had to manually patch the slow5 file to make dorado happy like #27.

Nevertheless, the next step for improving the model would be to polish the pore models. Do you happen to have modification-free real data from this reference?

[npdungca]

Hi! I uploaded the real data to the same folder.

Yes, we do get better mean baseq with v5. I had a dataset that was basecalled with 4.3 and the mean baseq was 26.7 then it improved to 30 after re-basecalling with v5.

Sorry, haven't had the chance yet to train a pore model with uncalled4.

[hasindu2008]

Is this real data free of any DNA modifications?

[npdungca]

that's amplicon data, sorry. but I added another dataset that wasn't PCR'd (Merlin_GFP_Merlinmapped_hd.fasta). DNA was extracted from cultures then sequenced. I don't have a lot of that kind of data, though. The usual workflow always involves PCR.

[hasindu2008]

So to make sure I understood correctly:
Data without modification - CMV.fasta
Data with modification - Merlin_GFP_Merlinmapped_hd.fasta
Is that right?

[npdungca]

Hi, Hasindu.

Yes, you can take CMV.fasta as modification-free, it is the universal reference sequence for CMV.

Then Merlin_GFP_Merlinmapped is also modification-free. As it was directly sequenced viral DNA using Nanopore.

Merlin_SF+ONT might have some modification as it was amplified with PCR, so there might be some artifacts/errors.

[hasindu2008]

If you can share some BLOW5/POD5 for this sequenced CMV.fasta dataset, we can try generating a pore model using poregen. You can share privately. When I say modification-free, I refer to base modifications.

[npdungca]

Unfortunately, I don't have the raw reads for CMV.fasta since I only downloaded it from Genbank. I do have some pod5s but since the samples were sequenced in May 2024, all the data is saved in a single pod5. I'm still trying to find a way on how to share a 240GB pod5.

[npdungca]

@hasindu2008 Hi! I shared some pod5s via globus. Let me know if you'll encounter any issues in accessing the directory. Thank you!

[hasindu2008]

Can you send me an email from yours? I can then share the globus ID with you. The email alias you currently added seem not to give me access to the dataset.