Mercurial > repos > marie-tremblay-metatoul > asca_w4m
diff galaxy/ASCA.Calculate_w4m.R @ 0:c5f11e6f8f99 draft
Uploaded
author | marie-tremblay-metatoul |
---|---|
date | Mon, 30 Jul 2018 07:29:40 -0400 |
parents | |
children |
line wrap: on
line diff
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/galaxy/ASCA.Calculate_w4m.R Mon Jul 30 07:29:40 2018 -0400 @@ -0,0 +1,84 @@ +ASCA.Calculate_w4m <- function (data, levels, equation.elements = "", scaling, only.means.matrix = FALSE, use.previous.asca = NULL) +{ + ASCA.GetEquationElement <- function(asca, evaluation, previous.asca) { + s <- list() + s$factors.evaluated <- evaluation + if (!is.null(previous.asca)) { + s$level.combinations <- previous.asca[[paste(evaluation, + collapse = "")]]$level.combinations + } + else { + s$level.combinations <- ASCA.GetRowRepeats(asca$levels[, + s$factors.evaluated, drop = FALSE]) + } + s$means.matrix <- matrix(nrow = dim(asca$data)[1], ncol = dim(asca$data)[2]) + for (p in 1:dim(s$level.combinations$row.patterns)[1]) { + mean.for.this.level.combination <- colMeans(asca$data[s$level.combinations$indices.per.pattern[[p]], + , drop = FALSE]) + for (i in s$level.combinations$indices.per.pattern[[p]]) { + s$means.matrix[i, ] <- mean.for.this.level.combination + } + } + s + } + s <- list() + dataAdjusted <- MetStaT.ScalePip(data, center = FALSE, scale = FALSE, + quietly = TRUE) + s$ssq.mean <- sum(rep(dataAdjusted$center.vector/dataAdjusted$scale.vector, + nrow(data))^2) + s$ssq <- sum(data^2) + s$data <- dataAdjusted$data + if (!is.numeric(levels)) { + stop("The supplied levels are not numeric.") + } + s$levels <- levels + if (!only.means.matrix) { + s$svd <- PCA.Calculate(s$data) + } + s$ee.names <- c() + if (identical(equation.elements, "")) { + equation.elements <- ASCA.GetPowerSet(c(1:dim(s$levels)[2]), + exclude.empty.set = TRUE) + } + if (is.character(equation.elements)) + equation.elements <- lapply(strsplit(strsplit(equation.elements, + split = ",")[[1]], split = ""), as.numeric) + for (ee in equation.elements) { + for (f in ee) if (f > dim(levels)[2] || f < 1) { + stop(paste("Factor ", f, " is beyond scope of study-design", + sep = "")) + } + } + if (dim(data)[1] != dim(levels)[1]) { + stop(paste("Number of rows in data (", dim(data)[1], + ") and study design (", dim(levels)[1], ") do not match", + sep = "")) + } + order.to.evaluate.ee <- sort(as.numeric(unlist(lapply(equation.elements, + paste, collapse = ""))), index.return = TRUE)$ix + s$remainder <- s$data + for (ee in order.to.evaluate.ee) { + new.equation.element <- ASCA.GetEquationElement(s, equation.elements[[ee]], + use.previous.asca) + reductions <- ASCA.GetPowerSet(equation.elements[[ee]], + exclude.empty.set = TRUE, exclude.complete.set = TRUE) + for (r in reductions) { + new.equation.element$means.matrix <- new.equation.element$means.matrix - + s[[c(paste(r, collapse = ""))]]$means.matrix + } + new.equation.element$ssq <- sum(new.equation.element$means.matrix^2) + if (!only.means.matrix) { + s$remainder <- s$remainder - new.equation.element$means.matrix + new.equation.element$reduced.matrix <- s$remainder + new.equation.element$svd <- PCA.Calculate(new.equation.element$means.matrix) + } + ee.name <- paste(equation.elements[[ee]], collapse = "") + s$ee.names <- c(s$ee.names, ee.name) + s[[ee.name]] <- new.equation.element + } + s$ssq.remainder <- sum(s$remainder^2) + if (!only.means.matrix) + asca.summary <- ASCA.GetSummary(s, quietly = TRUE) + return(list(s, asca.summary)) +} + \ No newline at end of file