medaka tools resolve_model --auto_model not able to get correct model

[replikation]

Hi

Medaka resolve model is not working on fastq files. Even though the model is listed in the read and it is a valid model in your medaka code.

Running medaka tools resolve_model in ontresearch/medaka:mr324_shae731e46af6f89b80b1590a9696e6f8630512e358-amd64

medaka tools resolve_model --auto_model consensus one_read.fastq

causes error:

Cannot import pyabpoa, some features may not be available.
Traceback (most recent call last):
  File "/home/epi2melabs/conda/bin/medaka", line 8, in <module>
    sys.exit(main())
  File "/home/epi2melabs/conda/lib/python3.8/site-packages/medaka/medaka.py", line 801, in main
    args = parser.parse_args()
  File "/home/epi2melabs/conda/lib/python3.8/argparse.py", line 1768, in parse_args
    args, argv = self.parse_known_args(args, namespace)
  File "/home/epi2melabs/conda/lib/python3.8/argparse.py", line 1800, in parse_known_args
    namespace, args = self._parse_known_args(args, namespace)
  File "/home/epi2melabs/conda/lib/python3.8/argparse.py", line 1988, in _parse_known_args
    positionals_end_index = consume_positionals(start_index)
  File "/home/epi2melabs/conda/lib/python3.8/argparse.py", line 1965, in consume_positionals
    take_action(action, args)
  File "/home/epi2melabs/conda/lib/python3.8/argparse.py", line 1874, in take_action
    action(self, namespace, argument_values, option_string)
  File "/home/epi2melabs/conda/lib/python3.8/argparse.py", line 1159, in __call__
    subnamespace, arg_strings = parser.parse_known_args(arg_strings, None)
  File "/home/epi2melabs/conda/lib/python3.8/argparse.py", line 1800, in parse_known_args
    namespace, args = self._parse_known_args(args, namespace)
  File "/home/epi2melabs/conda/lib/python3.8/argparse.py", line 1988, in _parse_known_args
    positionals_end_index = consume_positionals(start_index)
  File "/home/epi2melabs/conda/lib/python3.8/argparse.py", line 1965, in consume_positionals
    take_action(action, args)
  File "/home/epi2melabs/conda/lib/python3.8/argparse.py", line 1874, in take_action
    action(self, namespace, argument_values, option_string)
  File "/home/epi2melabs/conda/lib/python3.8/argparse.py", line 1159, in __call__
    subnamespace, arg_strings = parser.parse_known_args(arg_strings, None)
  File "/home/epi2melabs/conda/lib/python3.8/argparse.py", line 1800, in parse_known_args
    namespace, args = self._parse_known_args(args, namespace)
  File "/home/epi2melabs/conda/lib/python3.8/argparse.py", line 2006, in _parse_known_args
    start_index = consume_optional(start_index)
  File "/home/epi2melabs/conda/lib/python3.8/argparse.py", line 1946, in consume_optional
    take_action(action, args, option_string)
  File "/home/epi2melabs/conda/lib/python3.8/argparse.py", line 1874, in take_action
    action(self, namespace, argument_values, option_string)
  File "/home/epi2melabs/conda/lib/python3.8/site-packages/medaka/medaka.py", line 50, in __call__
    model = medaka.models.model_from_basecaller(input_file, variant=variant)
  File "/home/epi2melabs/conda/lib/python3.8/site-packages/medaka/models.py", line 140, in model_from_basecaller
    raise ValueError(
ValueError: Input file did not contain precisely 1 basecaller model reference.

Read file cat

@22179a98-884c-4c5a-a3fb-ffac2c1b50f1 runid=9042ff1a5dd45db8143188f898c539b3b2152c37 read=26 ch=311 start_time=2023-11-21T15:08:27.358316+01:00 flow_cell_id=FAX59214 protocol_group_id=20231121_PCR_random_isolates sample_id=20231121_PCR_random_isolates_15-22 barcode=barcode15 barcode_alias=barcode15 parent_read_id=22179a98-884c-4c5a-a3fb-ffac2c1b50f1 basecall_model_version_id=dna_r10.4.1_e8.2_400bps_sup@v4.2.0
ACTGGTGGAACGATTTTATTGCTGGCCTGGTATTGTTCGGCATACTGGTATTCGATGGGCGGCTGCGTTGCGCATTACAACGCAATCTGCACCACCAGAATATGCCCGTTTTATATCACCACCCACTCCACTACAGGCGGAAGCAAAAACGCACGCACAACAGAATAAAAACAAAGAGGTGGCATGATGAATCCATGGCGACGCTATAGCTGGGAAATTGCGCTGGCAGCCTTATTGATCTTTGAAATTCTGGCTTTCGGTCTGATTAACCCACGTTTATTAGATATTAATGTCTTACTTTTTAGCACTAGCGATTTTATTTGTATCGGTATTGTCGCTTTGCCGCTAACAATGGTTATTGTCAGCGGTGGTATGGATATTTCATTTGGTTCTACAATCGGGCTATGCGCGATTACCCTGGGTGTGTATTTAACTCATATGCGCTACCTTTAGCGATTATTCTTTACCTATACTTCGGCGCAATATGTAGGTGGTTCAATTGCCGGACTTATTATTTATACCGGCGTAAACCCGTTGGTGATTACCCTGGGAACCATGTATTTATTTGGTGGTAGCGCATTATTGTTATCTGGTATCGCTCGCGCCACGGGTTATGAAGGTATTGGTGGATTTG
+
67?77755955668==7577667CEEDFHLFFEFEGF221101+***,321,))*66DEGSGNJJFISIFGRSMPFHLMSJIIECB7776))'%%%%%*1119899><=BCH@KFENMGHSJOGIKMOSOQLSJJMGPHDFGFSEDCDKSLSLLPHMFISGIHMLIILMJLH???>633335>=====88888;CBGIIGGHFHSEE<//1<777IHJHIFKEHKA:...-02,+/01633443<EJFMKSSJOLJSLMGSSSHNSJIKIGIFSFLMKKOSHNFHSNHJSMSKISIEDE@@@?>AA@AACECB@FFGL@@A@@D>>9:::;BB98899ADIEDECCCHSKLKIKGSKKJKJOISKJIJIHJSSMGEB@<===9;==?FLEOQG>=<<=@DDEEE==<;332224/0150054055,*))'&$%%%$%&%%'('()//1..8..:7888795)(&&&&(,+***+7444:>CDCB3,)'&%$$$%%%%%&)))))88GFJIHOMJHKHSKSHFEDHJGFCC<60--,)+688>86633355AEIMNOSHFHEHMLMLKJJKMSJFEEHFNISKSNNGSQLFD;;;;;4,+,)'('()'(-334<GFBBA?ADCBB;:;977430'

The basecaller model is written down in the read file (basecall_model_version_id=dna_r10.4.1_e8.2_400bps_sup@v4.2.0). based on a quick glance in your medaka code, you search for some RG comment tags instead of the basecaller_model field?

[cjw85]

Hi @replikation,

Thank you for reporting this. The parsing code is written to handle SAM outputs that have been converted to fastq with samtools fastq. I realised after implementing that outputs from MinKNOW format the comment string in a different fashion but haven't yet had time to implement the logic to parse such data.

[replikation]

Hi @cjw85 Thanks for the quick reply. Yes, it would be great to add this, as it's the default output of MinKNOW. It would also be great to let medaka Exit if it can't find a model instead of using the default one (when you run e.g., medaka_consensus).

Edit: models are currently always listed like this in each read header which is "space" separated: basecall_model_version_id=dna_r10.4.1_e8.2_400bps_sup@v4.2.0

[cjw85]

I've done a quick implementation of parsing the Guppy/MinKNOW-style headers. This will appear shortly has version 1.11.2.