Some Feedback and Questions

[tjprins]

Hello, really looking forward to giving this a try with our lab's dataset. Detecting m6A modifications has been quite a challenge to get up and running, but m6Anet looks like a great alternative to other already established methods.

After getting things setup, I just have a few questions. I'm pretty sure I know the answer to these, but I just want to verify:

• For input data, you must always use transcriptome aligned data, correct? It would not be valid to attempt to basecall and align nanopore raw data to the genome and then attempt to put that data through m6Anet?
• For basecalling requirements, previous m6a detection methods processing single samples required the utilized basecalling algorithms to match those that were used for the trained models (for instance, EpiNano-SVM requires data basecalled with Guppy 3.1.5, which is very outdated). Since .fast5 files generated from more recent, accurate, and faster versions of Guppy are not backwards compatible with older Guppy versions (see the nanopore thread here), I'm curious if there are any 'basecalling algorithm' requirements for using m6anet? Can your input data be basecalled with any Guppy version?

Also, a little feedback:

• The command found in the documentation for m6anet-run_inference --input_dir demo_data --out_dir demo_data --infer_mod-rate --n_processes 4 contains a typo. --infer_mod-rate should actually be --infer_mod_rate, otherwise it throws an error.
• I had segmentation fault issues when doing an installation of m6anet using the setup.py method, which was due to running the installation within an already existing conda environment. I instead created a new empty conda environment (conda create --name m6anet) prior to the install which fixed the issue.
• (This is not related to m6anet, but useful for people here) Guppy caused a few headaches for me when attempting to run nanopolish index and eventalign, because of folder paths and file names in the sequencing_summary.txt file. Nanopolish originally thought I had corrupted or missing data, but it was simply looking for it in the wrong place. Although slower, if you have trouble with the first pre-processing step, try again without including the sequencing_summary.txt file in your pipeline and see if that gets you through to the next step.

Thank you so much!

[chrishendra93]

hi @tjprins, sorry for the delay in my response. Thanks a lot for your feedback, especially the --infer_mod_rate command which I was sure I had already fixed but must be somehow overwritten. I've just fixed it after reading this, thanks a lot for bringing this to my attention. As for your questions:

• We must always use transcriptome-aligned data. Genome mapping might result in a single position with multiple different 5-mers, mentioned in here and this will throw off an assertion error. If you still want to map to the genome, filtering the splice junctions might help
• There is no explicit requirement for the version of the basecalling algorithm used for m6Anet. Nevertheless, Nanopolish might not work effectively with some basecalling algorithms with many errors or RNA reads that are too short .

Anyway thanks for including m6Anet in your m6A detection protocol. Let me know if you have more questions / feedbacks and I will try to respond to them as soon as I can

Best Regards

Christopher Hendra

[tjprins]

Awesome, thank you for your response! I have a few more general questions then:

• Newer versions of Guppy separate out reads based on qscore into "pass" and "fail" folders. With other methods, such as EpiNano-SVM, I know that this is a problem since sometimes modifications can be mistakenly classified as errors, and so therefore you should generate your fast5 files with --disable_qscore_filtering option enabled. For m6ANet, would you recommend disabling qscore filtering of reads? Or just run m6ANet on "pass" reads?
• In your quick start documentation, you specify that "We also recommend a threshold of 0.9 for selecting m6A sites based on the probability_modified column, which can be relaxed at the expense of having lower model precision." I'm curious how that threshold was determined?
• Finally, to what extent do you consider mod_ratio vs. probability_modified? Based upon your paper and the documentation, I presume that the mod_ratio is the percentage of reads that were determined to have a modification with the probability specified by probability_modified, but is it possible to have a high mod probability and low mod ratio? Do you think it should be normalized in some fashion?

Thanks for all your help!!

[chrishendra93]

Hi @tjprins,

Thank you for the questions, I am converting this issue into discussion since it feels more like a discussion at this point. But to answer your questions:

• During training and validation, I have always included the failed fast5 files into nanopolish eventalign so as to obtain as many candidate sites as possible, again, also because some of the modified reads might have low qscore
• If you refer to the manuscript Figure 3e, we have tried several different thresholds, benchmarked against whether the predicted modified sites appear in m6Ace, miCLIP, or predicted by xPore to show statistically significant different in terms of signal-level between the wild-type HEK293T samples and the METTL3-KO samples. At p=0.9, we found an acceptable level of precision (as well as strong 3'UTR enrichment in our samples), but of course this can be adjusted depending on the research aim
• I haven't actually looked at sites with high mod probability and low mod ratio, although in my opinion this is unlikely to happen unless by some chance we keep on sampling modified reads during inference. This can be minimized by setting num_iteration parameters to a huge number (100 or 1000 for example), which is slow with the current m6anet implementation but we have fixed this issue in the development branch which should be released soon. Also, can you explain what do you mean by normalizing the probability here?

[chrishendra93]

issue converted as discussion @tjprins. Please reply in this discussion thread instead. Thanks!

[tjprins]

Thanks for moving it, this is a much better home for this thread. Thank you for your reply and comments, they are really helpful.

Also, can you explain what do you mean by normalizing the probability here?

I suppose this is mostly due to my ignorance about what m6ANet does exactly under the hood (I am still trying to wrap my head around everything included in your paper). I am mostly just curious about how differences in total read number between samples are taken into account with m6ANet? It looks like the minimum requirement is 20 reads (based upon my samples), but what if you are comparing that mod probability to the same transcript in another sample that has 500 reads? Is your confidence in the mod probability for both the same, or does m6ANet take the number of reads into account?

Basically I want to compare modifications between two sample conditions and I'm trying to figure out the best way to do that.

[chrishendra93]

The actual probability modified is derived by sampling 20 reads without replacement, m6Anet then assigns a probability score to each of these reads and calculates the probability that at least one of them is modified. So with enough sampling iterations (controlled by --num_iterations during inference), this should somewhat equalize the probability for sites with a high number of reads and sites with a low number of reads. Also in my experience, the probability-modified score is not correlated with the number of reads

I'm still not sure about the best way to compare modifications between two sample conditions. [xPore (https://github.com/GoekeLab/xpore/tree/master/xpore) provides an easy alternative to this since it was designed with sample comparison in mind (you can further filter the sites using m6Anet probability modified score before checking the results from xPore for example), another way till be to use the mod_ratio number and performs some kind of statistical test I suppose