Problems with methRead function when reading a .txt report

[Uxue Lazkano]

Hello everyone, 

I'm having problems to read a repor in .txt format (dowloaded from a published article).

Looks like this:

chr1 10469 + 0 0 CG CGC
chr1 10470 - 0 0 CG CGA
chr1 10471 + 0 0 CG CGG
chr1 10472 - 0 0 CG CGC
chr1 10484 + 5 0 CG CGG
chr1 10485 - 0 0 CG CGG
chr1 10489 + 4 1 CG CGC
chr1 10490 - 0 0 CG CGG
chr1 10493 + 2 3 CG CGC
chr1 10494 - 0 0 CG CGG

When using the readMeth() function insted of creating a methylRawList object it stores the information in an array.

GSE5879988 = methRead("2.Data/scRRBS/GSE196387/GSM5879988_Res_3_T.CpG_report.txt", sample.id = "test", treatment = 0,
           assembly = "hg38", header = FALSE, context = "CpG", resolution = "base",
           pipeline = list(fraction = FALSE, chr.col = 1, start.col = 2, end.col = 2,
                           coverage.col = 4, strand.col = 3, freqC.col = 5))

I don't know if I'm missing some parameter.

Thanks in advance for your kind help, I'm new analyzing RRBS data and perhaps I'm missing basic things.

Uxue

[Alexander Blume]

Hi Uxue,

this file looks like a bismark cytosine report, which is a format supported by methRead. Please use "bismarkCytosineReport" as the pipeline argument. 

Best,
Alex

[Uxue Lazkano]

I'm so sorry to bother you again but I have tried changing the pipeline argument to "bismarkCytosineReport" and the problem remains the same.

GSE5879988 = methRead("2.Data/scRRBS/GSE196387/GSM5879988_Res_3_T.CpG_report.txt.gz", sample.id = "test", treatment = 0,
           assembly = "hg38", context = "CpG", resolution = "base",
           pipeline = "bismarkCytosineReport")

The output instead of being a methylRawList is a simple table.

Thanks for your help!

Uxue

--
You received this message because you are subscribed to a topic in the Google Groups "methylkit_discussion" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/methylkit_discussion/vhitFWYGvos/unsubscribe.
To unsubscribe from this group and all its topics, send an email to methylkit_discus...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/methylkit_discussion/35462953-1b98-4a75-9a84-e4ff33182bd8n%40googlegroups.com.

[Alexander Blume]

Hi Uxue,

Since you are only providing a single file, the output should be a methylRaw object. 

Could you please post the output of

head(GSE5879988) 

Best,
Alex

You received this message because you are subscribed to the Google Groups "methylkit_discussion" group.
To unsubscribe from this group and stop receiving emails from it, send an email to methylkit_discus...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/methylkit_discussion/CAK0x_ZE94jUY6NTs-uA5Tn%2BDxbu83treE3beWOOHrxqyAf6e6w%40mail.gmail.com.

[Alexander Blume]

Hi Uxue,

You are right, the content of the objects is correct, only the formatting is messed up. It is possible that some other packages you are using are interfering with the representation. 

The error you get for the calculateDiffMeth is a bit cryptic, but it basically tells you that the methylRaw object is not the correct input object. Also, for running the differential analysis you need more than one sample, and they have to be merged beforehand, as described in the vignette (https://bioconductor.org/packages/devel/bioc/vignettes/methylKit/inst/doc/methylKit.html#3_Comparative_analysis). 

Best,
Alex

Am Mo., 11. Sept. 2023 um 15:34 Uhr schrieb Uxue Lazkano <uxu...@gmail.com>:

Hello Alex, 

> head(GSE5879988) methylRaw object with 6 rows -------------- chr start end strand coverage numCs numTs 1 chr1 10497 10497 + 106 95 11 2 chr1 10525 10525 + 106 103 3 3 chr1 10542 10542 + 106 102 4 4 chr1 10563 10563 + 106 100 6 5 chr1 10571 10571 + 104 96 8 6 chr1 10577 10577 + 88 72 16 -------------- sample.id: test assembly: hg38 context: CpG resolution: base

It seems that even in the environment panel does not appear as a methylRaw object, in reallity is stored correctly.

Anyway, when I continue to follow the tutorial and try to filterByCoverage(GSE5879988) it works and filters but the next steps as calculateDiffMeth(GSE5879988.filt) do not work and arise the following error:

GSE5879988.filt <- filterByCoverage(GSE5879988,
                               lo.count=10,
                               lo.perc=NULL,
                               hi.count=NULL,
                               hi.perc=99.9)
getCorrelation(GSE5879988.filt,plot=TRUE)

myDiff <- calculateDiffMeth(GSE5879988.filt,
                            overdispersion = "MN",
                            adjust="BH")
myDiff

Error in (function (classes, fdef, mtable) : unable to find an inherited method for function ‘calculateDiffMeth’ for signature ‘"methylRaw"’

Thanks for your infinite patience, I really appreaciate your help.

Uxue

To view this discussion on the web visit https://groups.google.com/d/msgid/methylkit_discussion/ddddde30-deae-4eef-99ac-5f13ae134048n%40googlegroups.com.

[Uxue Lazkano]

Thanks Alex, all solved!

I am very grateful for the help you have given me.

All the best, 

Uxue