Mercurial > repos > devteam > dwt_ivc_all

annotate execute_dwt_IvC_all.R @ 1:506ae7b0d85d draft default tip

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
"planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
author devteam
date Mon, 06 Jul 2020 20:31:56 -0400
parents
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
1
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
1 ###########################################################################################
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
2 ## code to do wavelet Indel vs. Control
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
3 ## signal is the difference I-C; function is second moment i.e. variance from zero not mean
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
4 ## to perform wavelet transf. of signal, scale-by-scale analysis of the function
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
5 ## create null bands by permuting the original data series
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
6 ## generate plots and table matrix of correlation coefficients including p-values
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
7 ############################################################################################
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
8 library("wavethresh");
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
9 library("waveslim");
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
10
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
11 options(echo = FALSE)
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
12
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
13 ## normalize data
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
14 norm <- function(data) {
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
15 v <- (data - mean(data)) / sd(data);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
16 if (sum(is.na(v)) >= 1) {
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
17 v <- data;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
18 }
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
19 return(v);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
20 }
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
21
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
22 dwt_cor <- function(data_short, names_short, data_long, names_long, test, pdf, table, filter = 4, bc = "symmetric", wf = "haar", boundary = "reflection") {
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
23 print(test);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
24 print(pdf);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
25 print(table);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
26
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
27 pdf(file = pdf);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
28 final_pvalue <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
29 title <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
30
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
31 short_levels <- wavethresh::wd(data_short[, 1], filter.number = filter, bc = bc)$nlevels;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
32 title <- c("motif");
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
33 for (i in 1:short_levels) {
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
34 title <- c(title, paste(i, "moment2", sep = "_"), paste(i, "pval", sep = "_"), paste(i, "test", sep = "_"));
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
35 }
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
36 print(title);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
37
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
38 ## loop to compare a vs a
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
39 for (i in seq_len(length(names_short))) {
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
40 wave1_dwt <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
41 m2_dwt <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
42 diff <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
43 var_dwt <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
44 out <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
45 out <- vector(length = length(title));
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
46
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
47 print(names_short[i]);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
48 print(names_long[i]);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
49
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
50 ## need exit if not comparing motif(a) vs motif(a)
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
51 if (names_short[i] != names_long[i]) {
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
52 stop(paste("motif", names_short[i], "is not the same as", names_long[i], sep = " "));
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
53 }
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
54 else {
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
55 ## signal is the difference I-C data sets
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
56 diff <- data_short[, i] - data_long[, i];
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
57
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
58 ## normalize the signal
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
59 diff <- norm(diff);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
60
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
61 ## function is 2nd moment
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
62 ## 2nd moment m_j = 1/N[sum_N(W_j + V_J)^2] = 1/N sum_N(W_j)^2 + (X_bar)^2
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
63 wave1_dwt <- waveslim::dwt(diff, wf = wf, short_levels, boundary = boundary);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
64 var_dwt <- waveslim::wave.variance(wave1_dwt);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
65 m2_dwt <- vector(length = short_levels)
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
66 for (level in 1:short_levels) {
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
67 m2_dwt[level] <- var_dwt[level, 1] + (mean(diff)^2);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
68 }
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
69
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
70 ## CI bands by permutation of time series
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
71 feature1 <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
72 feature2 <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
73 feature1 <- data_short[, i];
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
74 feature2 <- data_long[, i];
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
75 null <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
76 results <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
77 med <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
78 m2_25 <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
79 m2_975 <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
80
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
81 for (k in 1:1000) {
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
82 nk_1 <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
83 nk_2 <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
84 m2_null <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
85 var_null <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
86 null_levels <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
87 null_wave1 <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
88 null_diff <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
89 nk_1 <- sample(feature1, length(feature1), replace = FALSE);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
90 nk_2 <- sample(feature2, length(feature2), replace = FALSE);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
91 null_levels <- wavethresh::wd(nk_1, filter.number = filter, bc = bc)$nlevels;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
92 null_diff <- nk_1 - nk_2;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
93 null_diff <- norm(null_diff);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
94 null_wave1 <- waveslim::dwt(null_diff, wf = wf, short_levels, boundary = boundary);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
95 var_null <- waveslim::wave.variance(null_wave1);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
96 m2_null <- vector(length = null_levels);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
97 for (level in 1:null_levels) {
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
98 m2_null[level] <- var_null[level, 1] + (mean(null_diff)^2);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
99 }
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
100 null <- rbind(null, m2_null);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
101 }
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
102
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
103 null <- apply(null, 2, sort, na.last = TRUE);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
104 m2_25 <- null[25, ];
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
105 m2_975 <- null[975, ];
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
106 med <- apply(null, 2, median, na.rm = TRUE);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
107
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
108 ## plot
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
109 results <- cbind(m2_dwt, m2_25, m2_975);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
110 matplot(results, type = "b", pch = "*", lty = 1, col = c(1, 2, 2), xlab = "Wavelet Scale", ylab = c("Wavelet 2nd Moment", test), main = (names_short[i]), cex.main = 0.75);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
111 abline(h = 1);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
112
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
113 ## get pvalues by comparison to null distribution
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
114 out <- c(names_short[i]);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
115 for (m in seq_len(length(m2_dwt))) {
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
116 print(paste("scale", m, sep = " "));
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
117 print(paste("m2", m2_dwt[m], sep = " "));
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
118 print(paste("median", med[m], sep = " "));
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
119 out <- c(out, format(m2_dwt[m], digits = 4));
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
120 pv <- NULL;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
121 if (is.na(m2_dwt[m])) {
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
122 pv <- "NA";
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
123 }
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
124 else {
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
125 if (m2_dwt[m] >= med[m]) {
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
126 ## R tail test
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
127 tail <- "R";
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
128 pv <- (length(which(null[, m] >= m2_dwt[m]))) / (length(na.exclude(null[, m])));
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
129 }
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
130 else{
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
131 if (m2_dwt[m] < med[m]) {
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
132 ## L tail test
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
133 tail <- "L";
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
134 pv <- (length(which(null[, m] <= m2_dwt[m]))) / (length(na.exclude(null[, m])));
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
135 }
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
136 }
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
137 }
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
138 out <- c(out, pv);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
139 print(pv);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
140 out <- c(out, tail);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
141 }
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
142 final_pvalue <- rbind(final_pvalue, out);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
143 print(out);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
144 }
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
145 }
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
146
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
147 colnames(final_pvalue) <- title;
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
148 write.table(final_pvalue, file = table, sep = "\t", quote = FALSE, row.names = FALSE);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
149 dev.off();
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
150 }
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
151 ## execute
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
152 ## read in data
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
153 args <- commandArgs(trailingOnly = TRUE)
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
154
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
155 input_data <- read.delim(args[1]);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
156 input_data_names <- colnames(input_data);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
157
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
158 control_data <- read.delim(args[2]);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
159 control_data_names <- colnames(control_data);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
160
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
161 ## call the test function to implement IvC test
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
162 dwt_cor(input_data, input_data_names, control_data, control_data_names, test = "IvC", pdf = args[3], table = args[4]);
506ae7b0d85d "planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/dwt_ivc_all commit f929353ffb0623f2218d7dec459c7da62f3b0d24"
devteam
parents:
diff changeset
163 print("done with the correlation test");