methylation coverage in the gene region
0
0
Entering edit mode
4.5 years ago
rthapa ▴ 90

Hi,

I want to estimate the coverage of methylation in the gene-body region. I have coverage output from bismark.

Chr05   241 241 100 9   0
Chr05   258 258 100 1   0
Chr05   260 260 88.8888888888889    8   1
Chr05   261 261 100 2   0
Chr05   263 263 88.8888888888889    8   1
Chr05   419 419 100 15  0
Chr05   421 421 100 7   0
Chr05   438 438 100 15  0
Chr05   440 440 100 10  0
Chr05   441 441 100 15  0
Chr05   443 443 90  9   1
Chr05   597 597 95.2380952380952    20  1
Chr05   599 599 100 6   0
Chr05   616 616 95.2380952380952    20  1
Chr05   618 618 100 5   0

I also have bed file of gene region.

chr       start  end       geneid
Chr01   0   1951    Sb.001G000100
Chr01   9180    11180   Sb.001G000200
Chr01   21399   23399   Sb.001G000300
Chr01   20391   22391   Sb.001G000400
Chr01   50891   52891   Sb.001G000501
Chr01   51781   53781   Sb.001G000700
Chr01   60892   62892   Sb.001G000800
Chr01   77159   79159   Sb.001G000900
Chr01   79932   81932   Sb.001G001000
Chr01   101029  103029  Sb.001G001066
Chr01   104847  106847  Sb.001G001132
Chr01   107968  109968  Sb.001G001200
Chr01   110766  112766  Sb.001G001300
Chr01   114116  116116    Sb.001G001400
Chr01   150558  152558  Sb.001G001700

I tried using my python script to estimate the coverage in the gene region but since the coverage output file is very big, no success. I would appreciate any suggestion.

Thanks

methylation • 800 views
ADD COMMENT
0
Entering edit mode

methylkit (R/Bioconductor) has a convenience function to read Bismark outputs and aggregate calls over user-specified regions.

ADD REPLY
0
Entering edit mode

Thank you for your reply. I went through methylkit package. I tried to read my file in but got an error. I am not sure what I am doing wrong. Could you please suggest?

myobj = readBismarkCoverage( location = "/work/WGBS/methylkit",sample.id = "CHG_context_SRR5748817_1_bismark_bt2_pe.bed.gz.bismark.cov",assembly="unknown",treatment,context="CHG",min.cov=1)


Error in fread(filepath, ...) : 
  File '/work/WGBS/methylkit' is a directory. Not yet implemented.
ADD REPLY
0
Entering edit mode

From the manual: location => a list or vector of file paths to coverage files

You entered a directory (as the error describes) but need to give the path to the file.

ADD REPLY

Login before adding your answer.

Traffic: 2711 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6