diffHic - incorporating NB-loess, CNV, and ICE offset into DGEList object
Entering edit mode
wliao • 0
Last seen 3.9 years ago

I'm trying to incorporate NB-loess, CNV, and ICE offset in my DGEList object before running estimateDisp(). 

ice.corrected <- correctedContact(data,winsor.high=0.02,ignore.low=0.02,exclude.local=1,average=FALSE)
anchor1.bias <- ice.corrected$bias[anchors(data,type="first",id=TRUE),]
anchor2.bias <- ice.corrected$bias[anchors(data,type="second",id=TRUE),]
ice.offset <- log(anchor1.bias*anchor2.bias)

nb.offset <- normOffsets(data,type="loess")

cnv.offset <- normalizeCNV(cnv.data,cnv.margin)

design <- model.matrix(~factor(groups))
colnames(design) <- c("Intercept", "Cell-type")

y <- asDGEList(data)
y$offset <- ice.offset + nb.offset + cnv.offset

y <- estimateDisp(y,design)

This returns the following error:

Error in y[, j] <- fitted(locfit::locfit(y[, j] ~ x, weights = weights,  :
  number of items to replace is not a multiple of replacement length

It works fine without including the ICE offsets. Here is the sessionInfo():

R version 3.3.0 (2016-05-03)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: CentOS Linux 7 (Core)

 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C
 [9] LC_ADDRESS=C               LC_TELEPHONE=C

attached base packages:
[1] stats4    parallel  stats     graphics  grDevices utils     datasets
[8] methods   base

other attached packages:
 [1] reshape2_1.4.1             gridExtra_2.2.1
 [3] ggplot2_2.1.0              depmixS4_1.3-3
 [5] Rsolnp_1.16                MASS_7.3-45
 [7] nnet_7.3-12                rtracklayer_1.32.0
 [9] csaw_1.6.0                 edgeR_3.14.0
[11] limma_3.28.2               diffHic_1.4.0
[13] InteractionSet_1.0.2       SummarizedExperiment_1.2.0
[15] Biobase_2.32.0             GenomicRanges_1.24.0
[17] GenomeInfoDb_1.8.0         IRanges_2.6.0
[19] S4Vectors_0.10.0           BiocGenerics_0.18.0

loaded via a namespace (and not attached):
 [1] Rcpp_0.12.4.5           plyr_1.8.3              XVector_0.12.0
 [4] GenomicFeatures_1.24.0  bitops_1.0-6            tools_3.3.0
 [7] zlibbioc_1.18.0         biomaRt_2.28.0          statmod_1.4.24
[10] gtable_0.2.0            RSQLite_1.0.0           rhdf5_2.16.0
[13] lattice_0.20-33         BSgenome_1.40.0         Matrix_1.2-6
[16] DBI_0.4-1               stringr_1.0.0           Biostrings_2.40.0
[19] locfit_1.5-9.1          grid_3.3.0              AnnotationDbi_1.34.0
[22] XML_3.98-1.4            BiocParallel_1.6.0      magrittr_1.5
[25] Rhtslib_1.4.2           scales_0.4.0            Rsamtools_1.24.0
[28] GenomicAlignments_1.8.0 splines_3.3.0           colorspace_1.2-6
[31] KernSmooth_2.23-15      stringi_1.0-1           munsell_0.4.3
[34] RCurl_1.95-4.8          truncnorm_1.0-7
diffhic • 580 views
Entering edit mode
Aaron Lun ★ 27k
Last seen 20 hours ago
The city by the bay

That's probably because the ICE offsets have some NA's in them (because low-coverage regions get excluded). There's no reliable information to estimate the genomic biases for such regions, which is why they become NA. If that's the case in your ice.offset, you'll need to remove the offending rows prior to further analyses. 

Also, as I mentioned to you off-list, I don't think that adding offsets together is necessary as genomic biases from ICE should (largely) cancel out or be removed by the other normalisation strategies. Moreover, it's not clear that adding offsets together makes sense - especially if the genomic biases corrected for by ICE are somehow correlated with the abundances, in which case you'll basically be over-normalising by adding the ICE offsets to the offsets from everything else. Similarly, it's wrong to add the loess offsets to the normalizeCNV offsets. By default, normalizeCNV will smooth on the abundances as well, so any abundance-dependent trend will already be removed. Adding in the offsets from loess normalisation will result in over-normalisation.


Login before adding your answer.

Traffic: 250 users visited in the last hour
Help About
Access RSS

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

Powered by the version 2.3.6