limma: adding weights to mroast
0
0
Entering edit mode
@anthoula-gaigneaux-5021
Last seen 9.6 years ago
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 -- Anthoula Gaigneaux, PhD Laboratoire de Biologie Moléculaire et Cellulaire du Cancer (LBMCC) Hôpital Kirchberg 9, rue Steichen L-2540 Luxembourg tel: 00352/ 2468 4046 anthoula.gaigneaux@lbmcc.lu [[alternative HTML version deleted]]
Cancer limma Cancer limma • 879 views
ADD COMMENT

Login before adding your answer.

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