Search
Question: featureAlignedSignal() error: deparse(sig): object 'sig' not found.
0
gravatar for michael.alexander.finger
14 months ago by

Hello all.
I am attempting to use "featureAlignedDistribution" in the "ChIPpeakAnno" package to get the binding curves for some subset of annotated peaks I have stored as GRanges object.
My experimental test function is as follows:

bindingPattern<-function(gr, bwpath, outprefix, ups=2000, downs=2000) {
    wid <- width(gr)
    gr_recentered <- gr_center <- gr
    start(gr_center) <- start(gr) + floor(wid/2)
    width(gr_center) <- 1
    start(gr_recentered) <- start(gr_center) - ups
    end(gr_recentered) <- end(gr_center) + downs
    sig <- import(bwpath, format="BigWig", which=gr_recentered, as="RleList")
    aligned_signal <- featureAlignedSignal(sig, gr_center, upstream=ups, downstream=downs)
    distout <- paste(outprefix, "distribution.png", sep="")
    checkMakeFile(distout)
    png(filename=distout)
    featureAlignedDistribution(aligned_signal, gr_center, upstream=ups, downstream=downs, type="l")
    dev.off()
}

(Note that bwpath is the filepath to a bigwig file I have made from a bedgraph file which in turn was generated from MACS2 peakcalling. Both show up fine in IGV)

However, featureAlignedSignal() fails with the error: 

Error in deparse(sig) : object 'sig' not found

I thought perhaps sig was null or the wrong type, so after the sig<-import... line I checked it out:

1.class(sig) shows that sig is of type RleSimpleList, which is one of the accepted object types for featureAlignedSignal

2.Printing it out shows:
RleList of length 17
$chr1
numeric-Rle of length 258084 with 80751 runs
  Lengths:            39199                1                1                1                1 ...                1                2                3                1            34199
  Values :                0 67.9800796508789 69.5986633300781 67.9800796508789 66.3615036010742 ...   6.474289894104 4.85572004318237 3.23714995384216 1.61856997013092                0

$chr10
numeric-Rle of length 766161 with 310066 runs
  Lengths:            32592                1                1                2                2 ...               29                4                3                7             2282
  Values :                0 51.7943496704102 50.1757698059082  48.557201385498 46.9386291503906 ... 9.71144008636475 8.09286975860596   6.474289894104 4.85572004318237                0

...etc...

I tried looking at the source code for featureAlignedSignal() but I did not see any obvious reason why this error should occur. Intestingly, I found that if I convert sig to a list:

sig <- list(sig)

then featureAlignedSignal() runs and sets aligned_signal as a list object, which prints as:
[1] "Signal_aligned class after featureAlignedSignal: list"
[[1]]
              [,1]        [,2]        [,3]       [,4]        [,5]       [,6]       [,7]        [,8]       [,9]        [,10]       [,11]       [,12]      [,13]       [,14]       [,15]
   [1,]  90.721035 126.2082518 150.6082455 133.289509 114.9996338 103.952873  90.478250  98.6115808  94.241433  92.01589451  82.9518829  78.5412695  93.391681 118.7628157 124.5087505
   [2,] 213.489816 227.0858353 241.6529968 241.491138 249.9886482 254.601581 261.601912 261.1163364 248.774718 227.69279861 217.9408966 208.3103859 225.507725 229.5541584 225.4672607

...etc...

but that leads to a problem in featureAlignedDistribution():

Error in legend("topright", legend = colnames(density), col = col, lty = lty,  :
  'legend' is of length 0
In addition: Warning message:
In featureAlignedDistribution(aligned_signal, gr_center, upstream = ups,  :
  cvglists contain NA values. NA value will be omit.


I imagine it's because featureAlignedSignal wasn't meant to take a 'list' object in the first place, as its source code shows that it only wants RleList, SimpleRleList, or CompressedRleList.

I'm stumped. Can anyone help me out? Your time and help is much appreciated!

(Just a little extra info: Running on Ubuntu. R version is 3.3.1-1xenial, chippeakanno version is 3.6.5. Using the latest version of bioconductor)

ADD COMMENTlink modified 23 days ago by jun.yu0 • written 14 months ago by michael.alexander.finger0
0
gravatar for jun.yu
23 days ago by
jun.yu0
jun.yu0 wrote:

Alex, I ran your code with some of my test file. It ran successfully. The ChIPpeakAnno has been updated many times since your question. I hope the new version solved your problem.

 

> gr_center
GRanges object with 51385 ranges and 2 metadata columns:
                seqnames           ranges strand |     tx_id     tx_name
                   <Rle>        <IRanges>  <Rle> | <integer> <character>
      [1]           chr1 [ 11874,  11874]      + |         1  uc001aaa.3
      [2]           chr1 [ 69091,  69091]      + |         4  uc001aal.1
      [3]           chr1 [321084, 321084]      + |         5  uc001aaq.2
      [4]           chr1 [321146, 321146]      + |         6  uc001aar.2
      [5]           chr1 [322037, 322037]      + |         7  uc009vjk.2
      ...            ...              ...    ... .       ...         ...
  [51381] chrUn_gl000237   [ 2687,  2687]      - |     82956  uc011mgu.1
  [51382] chrUn_gl000241   [36876, 36876]      - |     82957  uc011mgv.2
  [51383] chrUn_gl000243   [11501, 11501]      + |     82958  uc011mgw.1
  [51384] chrUn_gl000243   [13608, 13608]      + |     82959  uc022brq.1
  [51385] chrUn_gl000247   [ 5817,  5817]      - |     82960  uc022brr.1
  -------
  seqinfo: 93 sequences (1 circular) from hg19 genome

> bwpath=bwfiles[1]
> sig <- import(bwpath, format="BigWig", which=gr_recentered, as="RleList")
> aligned_signal <- featureAlignedSignal(sig, gr_center, upstream=ups, downstream=downs)
> class(sig)
[1] "SimpleRleList"
attr(,"package")
[1] "IRanges"

ADD COMMENTlink modified 23 days ago • written 23 days ago by jun.yu0
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.2.0
Traffic: 127 users visited in the last hour