Entering edit mode
Dear Anthoula,
You are getting an error message from mroast() because the vector of
gene
weights must be of the same length as the number of the genes in the
set,
not of length equal to all the genes in your dataset. This is what
the
error message is telling you. Hence you need
gene.weights=fitA.genes$Amean[inds]
rather than
gene.weights=fitA.genes$Amean
There is no need to change the code of the function itself.
Best wishes
Gordon
> Date: Thu, 22 Dec 2011 14:18:11 +0100
> From: Anthoula Gaigneaux <anthoula.gaigneaux at="" lbmcc.lu="">
> To: bioconductor at r-project.org
> Subject: [BioC] limma: adding weights to mroast
>
> Dear List,
>
> First thank you for this very nice package.
>
> I'm currently using mroast in limma for differential geneset testing
on
> the basis of a MaList. Adding expression-related weights like
explained
> in the paper (Wu 2010) generated errors.
>
> When adding weights as a vector, dimensions did not fit (see below).
> Then I tried to build a weight list similar to iset object, but the
same
> error prompted. I've searched in several documentations, including
> mailing list and changelog, but have seen nothing related.
>
> Looking at the code, I've fixed the issue by adding a subscript in
the
> mroast function (gene.weights[[i]]). Is this fix is ok, or if
there's
> another way to define weights for mroast.
>
> Here is the code:
> ##### test
> identical(ma.genes$genes$SYMBOL, fitA.genes$genes$SYMBOL )
> gmt <- getGmt("c2.cp.v3.0.symbols.gmt" )
> universe = ma.genes$genes$SYMBOL
> gmt.list <- geneIds(gmt)
> inds <- symbols2indices(gmt.list, universe)
> design <- fitA.genes$design
> contrast <- fitA.genes$contrasts
> test <- mroast(iset=inds, y= ma.genes, design, contrast=
> contrast[,4],gene.weights=fitA.genes$Amean, nrot=10000)
> Error in roast(iset = iset[[i]], y = y, design = design,
contrast =
> contrast, :
> length of gene.weights disagrees with size of set
>
> weights <- lapply(inds, function(x) w <- fitA.genes$Amean[x] )
> test <- mroast(iset=inds, y= ma.genes, design, contrast=
contrast[,4],
> gene.weights= weights, nrot=10000)
> Error in roast(iset = iset[[i]], y = y, design = design,
contrast =
> contrast, :
> length of gene.weights disagrees with size of set
>
> ##### adapt function
> mroast2 <- function (iset = NULL, y, design, contrast =
ncol(design),
> set.statistic = "mean",
> gene.weights = NULL, array.weights = NULL, block = NULL,
> correlation, var.prior = NULL, df.prior = NULL, trend.var =
FALSE,
> nrot = 999, adjust.method = "BH")
> {
> if (is.null(iset))
> iset <- rep(TRUE, nrow(y))
> if (!is.list(iset))
> iset <- list(set = iset)
> nsets <- length(iset)
> if (is.null(names(iset)))
> names(iset) <- paste("set", 1:nsets, sep = "")
> fit <- lmFit(y, design = design, weights = array.weights,
> block = block, correlation = correlation)
> covariate <- NULL
> if (trend.var) {
> covariate <- fit$Amean
> if (is.null(covariate))
> covariate <- rowMeans(as.matrix(y))
> }
> sv <- squeezeVar(fit$sigma^2, df = fit$df.residual, covariate =
> covariate)
> var.prior <- sv$var.prior
> df.prior <- sv$df.prior
> pv <- adjpv <- active <- array(0, c(nsets, 3), dimnames =
> list(names(iset),
> c("Mixed", "Up", "Down")))
> if (nsets < 1)
> return(pv)
> for (i in 1:nsets) {
> out <- roast(iset = iset[[i]], y = y, design = design,
> contrast = contrast, set.statistic = set.statistic,
> gene.weights = gene.weights[[i]], array.weights =
> array.weights,
> block = block, correlation = correlation, var.prior =
> var.prior,
> df.prior = df.prior, nrot = nrot)[[1]]
> pv[i, ] <- out$P.Value
> active[i, ] <- out$Active.Prop
> }
> adjpv[, "Mixed"] <- p.adjust(pv[, "Mixed"], method =
adjust.method)
> adjpv[, "Up"] <- p.adjust(pv[, "Up"], method = adjust.method)
> adjpv[, "Down"] <- p.adjust(pv[, "Down"], method =
adjust.method)
> list(P.Value = pv, Adj.P.Value = adjpv, Active.Proportion =
active)
> }
>
> test <- mroast2(iset=inds, y= ma.genes, design, contrast=
contrast[,4],
> gene.weights= weights, nrot=10000) #ok
>
> ##### sessionInfo()
> R version 2.14.0 (2011-10-31)
> Platform: x86_64-apple-darwin9.8.0/x86_64 (64-bit)
>
> locale:
> [1] en_GB.UTF-8/en_GB.UTF-8/en_GB.UTF-8/C/en_GB.UTF-8/en_GB.UTF-8
>
> attached base packages:
> [1] stats graphics grDevices utils datasets methods base
>
> other attached packages:
> [1] GSEABase_1.16.0 graph_1.32.0 annotate_1.32.0
> [4] AnnotationDbi_1.16.0 Biobase_2.14.0 limma_3.10.0
>
> loaded via a namespace (and not attached):
> [1] DBI_0.2-5 IRanges_1.12.1 RSQLite_0.10.0 tools_2.14.0
> [5] XML_3.4-3 xtable_1.6-0
>
> Thanks for helping and merry Christmas,
>
> Anthoula Gaigneaux, PhD
> Bioinformatician
> Laboratoire de Biologie Mol?culaire et Cellulaire du Cancer (LBMCC)
> H?pital Kirchberg
> 9, rue Steichen
> L-2540 Luxembourg
> anthoula.gaigneaux at lbmcc.lu
______________________________________________________________________
The information in this email is confidential and
intend...{{dropped:5}}