comparison ASCA.Calculate_w4m.R @ 0:93312041f1d5 draft default tip

planemo upload for repository https://github.com/workflow4metabolomics/ascaw4m commit 7ea9b0f8abc5a60c2c04fd2098788497f14766b6
author marie-tremblay-metatoul
date Fri, 21 Sep 2018 05:51:14 -0400
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:93312041f1d5
1 ASCA.Calculate_w4m <- function (data, levels, equation.elements = "", scaling, only.means.matrix = FALSE, use.previous.asca = NULL)
2 {
3 ASCA.GetEquationElement <- function(asca, evaluation, previous.asca) {
4 s <- list()
5 s$factors.evaluated <- evaluation
6 if (!is.null(previous.asca)) {
7 s$level.combinations <- previous.asca[[paste(evaluation,
8 collapse = "")]]$level.combinations
9 }
10 else {
11 s$level.combinations <- ASCA.GetRowRepeats(asca$levels[,
12 s$factors.evaluated, drop = FALSE])
13 }
14 s$means.matrix <- matrix(nrow = dim(asca$data)[1], ncol = dim(asca$data)[2])
15 for (p in 1:dim(s$level.combinations$row.patterns)[1]) {
16 mean.for.this.level.combination <- colMeans(asca$data[s$level.combinations$indices.per.pattern[[p]],
17 , drop = FALSE])
18 for (i in s$level.combinations$indices.per.pattern[[p]]) {
19 s$means.matrix[i, ] <- mean.for.this.level.combination
20 }
21 }
22 s
23 }
24 s <- list()
25 dataAdjusted <- MetStaT.ScalePip(data, center = FALSE, scale = FALSE,
26 quietly = TRUE)
27 s$ssq.mean <- sum(rep(dataAdjusted$center.vector/dataAdjusted$scale.vector,
28 nrow(data))^2)
29 s$ssq <- sum(data^2)
30 s$data <- dataAdjusted$data
31 if (!is.numeric(levels)) {
32 stop("The supplied levels are not numeric.")
33 }
34 s$levels <- levels
35 if (!only.means.matrix) {
36 s$svd <- PCA.Calculate(s$data)
37 }
38 s$ee.names <- c()
39 if (identical(equation.elements, "")) {
40 equation.elements <- ASCA.GetPowerSet(c(1:dim(s$levels)[2]),
41 exclude.empty.set = TRUE)
42 }
43 if (is.character(equation.elements))
44 equation.elements <- lapply(strsplit(strsplit(equation.elements,
45 split = ",")[[1]], split = ""), as.numeric)
46 for (ee in equation.elements) {
47 for (f in ee) if (f > dim(levels)[2] || f < 1) {
48 stop(paste("Factor ", f, " is beyond scope of study-design",
49 sep = ""))
50 }
51 }
52 if (dim(data)[1] != dim(levels)[1]) {
53 stop(paste("Number of rows in data (", dim(data)[1],
54 ") and study design (", dim(levels)[1], ") do not match",
55 sep = ""))
56 }
57 order.to.evaluate.ee <- sort(as.numeric(unlist(lapply(equation.elements,
58 paste, collapse = ""))), index.return = TRUE)$ix
59 s$remainder <- s$data
60 for (ee in order.to.evaluate.ee) {
61 new.equation.element <- ASCA.GetEquationElement(s, equation.elements[[ee]],
62 use.previous.asca)
63 reductions <- ASCA.GetPowerSet(equation.elements[[ee]],
64 exclude.empty.set = TRUE, exclude.complete.set = TRUE)
65 for (r in reductions) {
66 new.equation.element$means.matrix <- new.equation.element$means.matrix -
67 s[[c(paste(r, collapse = ""))]]$means.matrix
68 }
69 new.equation.element$ssq <- sum(new.equation.element$means.matrix^2)
70 if (!only.means.matrix) {
71 s$remainder <- s$remainder - new.equation.element$means.matrix
72 new.equation.element$reduced.matrix <- s$remainder
73 new.equation.element$svd <- PCA.Calculate(new.equation.element$means.matrix)
74 }
75 ee.name <- paste(equation.elements[[ee]], collapse = "")
76 s$ee.names <- c(s$ee.names, ee.name)
77 s[[ee.name]] <- new.equation.element
78 }
79 s$ssq.remainder <- sum(s$remainder^2)
80 if (!only.means.matrix)
81 asca.summary <- ASCA.GetSummary(s, quietly = TRUE)
82 return(list(s, asca.summary))
83 }
84