Shians / NanoMethViz

Apache License 2.0
21 stars 2 forks source link

Finding genes with highest and lowest probability of methylation #12

Open 0Bernhard opened 2 years ago

0Bernhard commented 2 years ago

Hello

Great tool! It's the easiest tool I could find for analysing methylation data from megalodon. I am working with bacterial species. Is there any way to rank genes by methylation probability? Sorry if this is a really trivial question.

Thank you for your help, Bernhard

Shians commented 2 years ago

Hi Bernhard,

Good question, I don't provide a way to do this that easily since I've been focused on differential methylation (which I outsource to the base package). It would make sense to look at methylation levels within a sample and be able to rank them. I'll have a think about this and work it into the next release. For now the prototype of this looks like this:

library(dplyr)
library(purrr)
library(NanoMethViz)
# REPLACE WITH YOUR OWN DATA
nmr <- load_example_nanomethresult()

get_region_methy_stats <- function(index, regions) {
    query_methy(nmr, regions$chr[index], regions$start[index], regions$end[index]) %>%
        summarise(
            mean_methy_prob = mean(e1071::sigmoid(statistic)), # sigmoid transforms LLH to probability
            prevalence = mean(e1071::sigmoid(statistic) > 0.5)
        )
}

region_methy_stats <- function(nmr, regions) {
    bind_cols(
        regions,
        map_df(1:nrow(regions), get_region_methy_stats, regions = regions)
    )
}

# REPLACE WITH YOUR OWN ANNOTATION
gene_anno <- exons_to_genes(NanoMethViz::exons(nmr))

region_methy_stats(nmr, gene_anno)

Cheers, Shian

0Bernhard commented 2 years ago

Hello Shian

Thank you very much for this - it works like a dream with your trail data but i am having trouble adapting it for my data which is derived from megalodon. I am not sure if I am doing something wrong though or if I need to adapt it somehow. Sorry for the very naive question.

Thank you, Bernhard

Shians commented 2 years ago

It shouldn't matter what your source data is from if you've imported it into the NanoMethViz format using create_tabix_file() and have a NanoMethResult object. Are you getting any specific error?

BaylorLyu commented 2 years ago

Sir, i am confuse with the ruselt produced by code. Using NanoMethViz to visualing 6ma detected by the megelodon, the col prevalence have some probability is 1, is this probability trustworthy?

Shians commented 2 years ago

Sir, i am confuse with the ruselt produced by code. Using NanoMethViz to visualing 6ma detected by the megelodon, the col prevalence have some probability is 1, is this probability trustworthy?

Please open a new issue for this. I have not tested this package extensively outside of 5mC methylation. I have done a single experiment with 6mA using Megalodon and my impression was that the prevalence is much higher than expected (though only in the 5-10% range, not 100%). I can't say what the issue is without more detail about what kind of experiment it is, and how often you have 100% prevalence.

BaylorLyu commented 2 years ago
font{
    line-height: 1.6;
}
ul,ol{
    padding-left: 20px;
    list-style-position: inside;
}

Thanks for you reply. I open a new issue on git  already.

                ***@***.***

On 2/24/2022 07:55,Shian ***@***.***> wrote: 

Sir, i am confuse with the ruselt produced by code. Using NanoMethViz to visualing 6ma detected by the megelodon, the col prevalence have some probability is 1, is this probability trustworthy?

Please open a new issue for this. I have not tested this package extensively outside of 5mC methylation. I have done a single experiment with 6mA using Megalodon and my impression was that the prevalence is much higher than expected (though only in the 5-10% range, not 100%). I can't say what the issue is in your data without more detail about what kind of experiment it is, and how often you have 100% prevalence.

—Reply to this email directly, view it on GitHub, or unsubscribe.Triage notifications on the go with GitHub Mobile for iOS or Android.

You are receiving this because you commented.Message ID: @.***>