code in supplementary 2 of DESeq
1
0
Entering edit mode
Yuan Jian ▴ 250
@yuan-jian-2603
Last seen 8.3 years ago

I am trying to make the code (in supplementary 2 of DESeq paper) to work in current DESeq (not DESeq2).

I got error in function  "baseVarFunc":

Error in baseVarFunc(cdsFly, "A")(10^xg) : could not find function "rvf"

 

##here is the code (file fly_RNA_counts.tsv was downloaded from supplementary 3).
baseVarFunc <- function( cds, cond ) {
rvf <- rawVarFunc( cds, cond )
sf <- sizeFactors(cds)[ conditions(cds) == cond ]
xim <- sum(1/sf) / length(sf)
function( q ) rvf( q ) + xim * q
}

rawVarFunc <- function( cds, condOrName ) {
stopifnot( is( cds, "CountDataSet" ) )
res <- cds@fitInfo[[ as.character(condOrName) ]]
if( is.null(res) ) {
res <- cds@fitInfo[[ cds@dispTable[ as.character(condOrName) ] ]]
if( is.null(res) )
stop( sprintf( "No base variance function found for condition or with name '%s'.", condOrName ) )
}
res
}


countsTableFly <- read.delim( "fly_RNA_counts.tsv" )
condsFly <- c( "A", "A", "B", "B" )
rownames( countsTableFly ) <- paste( "Gene", 1:nrow(countsTableFly), sep="_" )
cdsFly <- newCountDataSet( countsTableFly, condsFly )
cdsFly <- estimateSizeFactors( cdsFly )
cdsFly <- estimateDispersions( cdsFly )
gg<-log10( baseVarFunc(cdsFly,"A")( 10^xg ) ) ## this sentence got error


> sessionInfo()
R version 3.2.3 (2015-12-10)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 7 x64 (build 7601) Service Pack 1

locale:
[1] LC_COLLATE=English_United States.1252 LC_CTYPE=English_United States.1252
[3] LC_MONETARY=English_United States.1252 LC_NUMERIC=C
[5] LC_TIME=English_United States.1252

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

other attached packages:
[1] DESeq_1.22.1 lattice_0.20-33 locfit_1.5-9.1 Biobase_2.30.0 BiocGenerics_0.16.1

loaded via a namespace (and not attached):
[1] IRanges_2.4.7 XML_3.98-1.3 grid_3.2.3 xtable_1.8-2 DBI_0.3.1
[6] stats4_3.2.3 RSQLite_1.0.0 genefilter_1.52.1 annotate_1.48.0 S4Vectors_0.8.11
[11] splines_3.2.3 RColorBrewer_1.1-2 geneplotter_1.48.0 survival_2.38-3 AnnotationDbi_1.32.3

deseq • 1.2k views
ADD COMMENT
0
Entering edit mode
@wolfgang-huber-3550
Last seen 7 days ago
EMBL European Molecular Biology Laborat…

Dear Yuan Jian

Thank you for posting the question, and in such a reproducible way. The implementation of the fitInfo slot of the CountDataSet class in the DESeq package is different in version 1.22.1 than in version 1.1.12, which the paper refers to (as indicated in the sessionInfo of Supplement S2). (It is no longer a function, but a list one of whose elements is a function.) There are two options for you:

  • revert back to R 2.12 and DESeq 1.1.12
  • understand the DESeq code well enough to adapt what you want to do to the more recent versions of DESeq.

 

ADD COMMENT

Login before adding your answer.

Traffic: 716 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