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

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.