Training an R10.4.1 methylation model

[hasindu2008]

This is closely related to #1059. As that issue was for the nucleotide models, I am opening another issue to discuss about the cpg models.

to train a cpg methylation model, I have done what you suggested at #825 some time ago. I made a hybrid reference genome - one with all CGs changed into MG, and the other with the CGs left unchanged. Then reads were mapped to the corresponding reference and the BAMs were merged together.

Now the doubt is about what I should be using as an input cpg model. I created an input model like below:

#model_name r10_450bps.cpg.9mer.template.model
#kit    r10_450bps
#strand template
#k      9
#alphabet       cpg
kmer    level_mean      level_stdv      sd_mean sd_stdv weight
AAAAAAAAA       54.204906       3.0     1       1       1       1
AAAAAAAAC       58.590797       3.0     1       1       1       1
AAAAAAAAG       55.952068       3.0     1       1       1       1
AAAAAAAAM       58.590797       3.0     1       1       1       1
AAAAAAAAT       58.433520       3.0     1       1       1       1
AAAAAAACA       63.659959       3.0     1       1       1       1
AAAAAAACC       65.805750       3.0     1       1       1       1
AAAAAAACG       64.072506       3.0     1       1       1       1
AAAAAAACM       65.805750       3.0     1       1       1       1
...

Here what I did was for aby methylated k-mers, I put the same values as for the corresponding nucleotide k-mers from the ONT's base model. Then I did the following command:

nanopolish  train -r positive_and_negative_pass.fastq -g hg38noAlt_hybrid.fa -b merged.bam  -t 40 --train-kmers=all --input-model-filename r10.4.1_400bps.cpg.9mer.model -d ../meth/

Unfortunately, none of the k-mers was trained. Note that now the nucleotide training works after the modification you made in [#1059], but still not for cpg. Is there a problem with the approach I am following or wonder if there is something else in the code that prevents the training?

[jts]

There's a few different reasons this could fail. Were the nucleotide (ACGT only) k-mers trained, or nothing at all trained? Were any status messages printed to the terminal?

[hasindu2008]

Oooohhh. I added fprintfs here andthere and finally found the potential problem - as the base model is using the inbuilt old r10 model. I manually replaced the .inl model and seems like itis now actually doing the training. Will keep you posted when the training is done.

[jts]

Ah yes, good find!

[hasindu2008]

well, will take a while it seems and the server is crying :D