Single molecule m6A probability predictions

[Stakaitis]

Hi @chrishendra93,

Is there a script for single molecule m6A probability prediction? It would be very useful to know the percentage of modified/unmodified sites in a transcript.

Thanks,
Rytis Stakaitis

[chrishendra93]

hi @Stakaitis, apology for the late reply. I have yet to provide the script for this functionality but thinking of providing this in some format in the future.

Currently I'd say that the single molecule prediction is a bit noisy and not well calibrated (i.e the probability 0.5 does not reflect 50% of the times it is modified), and maybe the single molecule representation (the representation extracted from the last layer before pooling and classification in the network) might be more useful if you can perhaps do some clustering on it, then again I have yet to explore this

May I know what specific use cases you are interested to use this for so that I can take this into account for the next release?

Regards

Christopher Hendra

[Stakaitis]

Thank you for the reply,

I understand that the m6A modification calling in an individual read can have a lot of false positives/negatives, I think the same is with the Epinano 1.2. However, some sort of estimation on the proportion of modified/unmodified transcripts in a given sample, would be useful in treatment efficiency, disease progression, or tumor malignancy research.

In my particular case, I want to determine whether I can detect changes in m6A modification levels in specific transcripts during tumor progression. In other words, I want to see if transcript/gene X has more m6A modifications in high malignancy cancer vs low malignancy cancer.

If I understand correctly, now the final result from m6Anet shows a probability that a specific site of a transcript is modified. In this case, the information on m6A modification percentage is lost. Then, I assume, that in the scenario where there are 2 samples, and sample A has 5/100 m6A modifications at a site Y of a transcript X, and sample B has 35/100 m6A modifications at the same location, both of these samples will show the probability that site Y of a transcript X is unmodified. But these additional 30 modifications could be the cause of initiation of a certain biological processes - which could help detect early stages of various diseases, like cancer.

Rytis Stakaitis

[chrishendra93]

Hi @Stakaitis , sorry for the slow response, I was only back recently

I am planning to add this feature where you can output the individual read probability along with the site probability, hopefully in the near future (1 or 2 weeks time if everything goes well). I am thinking of outputting individual read probability and let the user decides the threshold to determine whether a read is modified or not (we will be providing a recommended value of course). Hopefully it will be useful

Regards

Christopher Hendra

[chrishendra93]

hi @Stakaitis, we have just released a new version of m6anet that can perform single molecule quantification. Do try it and let us know if it's working for you

Regards

Christopher Hendra

[Stakaitis]

Thank you! I will try to run it on my data soon and let you know if I notice any issues.

[Stakaitis]

Also, is this update also available on nf-core/nanoseq -r dev or is it more complicated to update that pipeline?

[chrishendra93]

hi @Stakaitis , the update is currently not available through nanoseq. You have to install this separately from github / pypi for now, but we're working to make the update available through nanoseq, hopefully within this month or next month

[MaestSi]

Hi @chrishendra93 and @Stakaitis ,
I just installed m6anet following installation instructions and ran m6anet-dataprep and m6anet-run_inference on the demo dataset. However, I am not able to find read-level modification results. Is there another script I should run? I was not able to find additional information in the docs.
Thanks in advance,
Simone

[WeipengMO]

hi @MaestSi, if you're looking for a feature that can label each read as modified, you can try the branch named read_id. It is not fully tested yet, mind you, but this will give you an additional file with the read index of each read (as per nanopolish read index, which if you have summary.txt from nanopolish eventalign you can map it back to the actual read ids) and the probability modified for each one. You might need to set a threshold to call the modification on the read level (our threshold is 0.033379376)

Hi @chrishendra93 , thank you for developing such an excellent tool. I have tried the "read_id" branch. And got the "data.indiv_proba.csv.gz" results file. I want to know how to map the read index to the corresponding read_id?

This is the eventalign.txt file from nanopolish eventalign:

head eventalign.txt | column -t
contig       position  reference_kmer  read_index  strand  event_index  event_level_mean  event_stdv  event_length  model_kmer  model_mean  model_stdv  standardized_level  start_idx  end_idx
AT1G01020.2  81        GTGAA           7           t       1            92.95             4.060       0.00332       GTGAA       101.50      4.49        -1.79               48992      49002
AT1G01020.2  81        GTGAA           7           t       2            101.64            3.559       0.00664       GTGAA       101.50      4.49        0.03                48972      48992
AT1G01020.2  81        GTGAA           7           t       3            108.44            3.056       0.00299       GTGAA       101.50      4.49        1.46                48963      48972
AT1G01020.2  82        TGAAC           7           t       4            131.81            3.539       0.00398       TGAAC       118.66      5.56        2.23                48951      48963
AT1G01020.2  82        TGAAC           7           t       5            124.34            4.994       0.00232       TGAAC       118.66      5.56        0.96                48944      48951
AT1G01020.2  82        TGAAC           7           t       6            129.49            5.299       0.00531       TGAAC       118.66      5.56        1.83                48928      48944
AT1G01020.2  82        TGAAC           7           t       7            109.05            8.522       0.00232       TGAAC       118.66      5.56        -1.63               48921      48928
AT1G01020.2  82        TGAAC           7           t       8            126.58            8.079       0.01129       TGAAC       118.66      5.56        1.34                48887      48921
AT1G01020.2  83        GAACA           7           t       9            96.93             2.319       0.00531       GAACA       95.52       3.45        0.39                48871      48887

Thank you very much!

[WeipengMO]

I think I got the solution. Add --summary summary.txt parameter to nanopolish eventalign command.

[chrishendra93]

hi @WeipengMO, you are absolutely right. The summary.txt provides a mapping between read_id and read index, probably it is a good idea to provide this mapping in the next release.

For the two of you that need this feature @WeipengMO and @MaestSi , can I just check what are you planning to use the individual read probability for? Just curious

[WeipengMO]

Thank you for your kind response @chrishendra93 . I'm trying some experimental test. I want to use the m6A methyltransferase to label RNA and check whether I can check this signal in single-molecule level.

[MaestSi]

Hi @chrishendra93 , no concrete plan at the moment yet, but it will be nice to find a biological application for this new exciting feature!

[MaestSi]

Dear @chrishendra93 ,
I just saw that read_id branch of m6anet was removed. Does the main branch provide single-read results now?
Thanks,
Simone

[chrishendra93]

hi @MaestSi

Yeah I have just merged the read_id branch into the development branch. It will support single-read results by default now (no need additional flag) and the latest implementation also improves the speed for high --num_iterations options so if you want to filter the read modification by selecting sites with high probability_modified number, you can set higher num_iterations to make it more reliable

[MaestSi]

Thanks for the information!
SM