0
|
1 asca_w4m <- function(datamatrix, samplemetadata, factors, variablemetadata, threshold, scaling="none", nPerm)
|
|
2 {
|
|
3 ## Transpose
|
|
4 # datamatrix <- t(datamatrix)
|
|
5
|
|
6 # Check sample ID's
|
|
7 rownames(datamatrix) <- make.names(rownames(datamatrix), unique = TRUE)
|
|
8 colnames(datamatrix) <- make.names(colnames(datamatrix), unique = TRUE)
|
|
9 rownames(samplemetadata) <- make.names(rownames(samplemetadata), unique = TRUE)
|
|
10 rownames(variablemetadata) <- make.names(rownames(variablemetadata), unique = TRUE)
|
|
11
|
|
12 if(!identical(rownames(datamatrix), rownames(samplemetadata)))
|
|
13 {
|
|
14 if(identical(sort(rownames(datamatrix)), sort(rownames(samplemetadata))))
|
|
15 {
|
|
16 cat("\n\nMessage: Re-ordering dataMatrix sample names to match sampleMetadata\n")
|
|
17 datamatrix <- datamatrix[rownames(samplemetadata), , drop = FALSE]
|
|
18 stopifnot(identical(sort(rownames(datamatrix)), sort(rownames(samplemetadata))))
|
|
19 }else {
|
|
20
|
|
21 cat("\n\nStop: The sample names of dataMatrix and sampleMetadata do not match:\n")
|
|
22 print(cbind.data.frame(indice = 1:nrow(datamatrix),
|
|
23 dataMatrix=rownames(datamatrix),
|
|
24 sampleMetadata=rownames(samplemetadata))[rownames(datamatrix) != rownames(samplemetadata), , drop = FALSE])
|
|
25 }
|
|
26 }
|
|
27
|
|
28 # Check feature ID's
|
|
29 if(!identical(colnames(datamatrix), rownames(variablemetadata)))
|
|
30 {
|
|
31 if(identical(sort(colnames(datamatrix)), sort(rownames(variablemetadata))))
|
|
32 {
|
|
33 cat("\n\nMessage: Re-ordering dataMatrix variable names to match variableMetadata\n")
|
|
34 datamatrix <- datamatrix[, rownames(variablemetadata), drop = FALSE]
|
|
35 stopifnot(identical(sort(colnames(datamatrix)), sort(rownames(variablemetadata))))
|
|
36 }else {
|
|
37 cat("\n\nStop: The variable names of dataMatrix and variableMetadata do not match:\n")
|
|
38 print(cbind.data.frame(indice = 1:ncol(datamatrix),
|
|
39 dataMatrix=colnames(datamatrix),
|
|
40 variableMetadata=rownames(variablemetadata))[colnames(datamatrix) != rownames(variablemetadata), , drop = FALSE])
|
|
41 }
|
|
42 }
|
|
43
|
|
44 # Design
|
|
45 design <- data.matrix(samplemetadata[, colnames(samplemetadata) %in% factors])
|
|
46
|
|
47 # Scaling if scaling!=none
|
|
48 datamatrix <- prep(datamatrix, scaling)
|
|
49
|
|
50 # Computation of the A-SCA model
|
|
51 data.asca <- ASCA.Calculate_w4m(datamatrix, design, scaling=scaling)
|
|
52
|
|
53 # Permutation test
|
|
54 data.asca.permutation <- ASCA.DoPermutationTest(data.asca[[1]], perm=nPerm)
|
|
55 p <- c(data.asca.permutation, 0)
|
|
56
|
|
57
|
|
58 # % of explained variance
|
|
59 ssq <- (data.asca[[2]]$summary.ssq)
|
|
60 ssq <- cbind(round(rbind(ssq[2], ssq[3],ssq[4],ssq[5])*100, 2), p)
|
|
61 rownames(ssq) <- c(factors[1], factors[2], "Interaction", "Residuals")
|
|
62 colnames(ssq) <- c("% of explained variance", "Permutation p-value")
|
|
63
|
|
64 # Add Scores and loadings at the end of meatadata files
|
|
65 noms <- colnames(samplemetadata)
|
|
66 samplemetadata <- cbind(samplemetadata, (data.asca[[1]]$'1'$means.matrix + data.asca[[1]]$remainder) %*% data.asca[[1]]$'1'$svd$v[, 1:2],
|
|
67 (data.asca[[1]]$'2'$means.matrix + data.asca[[1]]$remainder) %*% data.asca[[1]]$'2'$svd$v[, 1:2],
|
|
68 (data.asca[[1]]$'12'$means.matrix + data.asca[[1]]$remainder) %*% data.asca[[1]]$'12'$svd$v[, 1:2])
|
|
69 colnames(samplemetadata) <- c(noms, paste(factors[1],"XSCOR-p1", sep="_"), paste(factors[1],"XSCOR-p2", sep="_"),
|
|
70 paste(factors[2],"XSCOR-p1", sep="_"), paste(factors[2],"XSCOR-p2", sep="_"),
|
|
71 "Interact_XSCOR-p1", "Interact_XSCOR-p2")
|
|
72
|
|
73 noms <- colnames(variablemetadata)
|
|
74 variablemetadata <- cbind(variablemetadata, data.asca[[1]]$'1'$svd$v[, 1:2], data.asca[[1]]$'2'$svd$v[, 1:2], data.asca[[1]]$'12'$svd$v[, 1:2])
|
|
75 colnames(variablemetadata) <- c(noms, paste(factors[1],"XLOAD-p1", sep="_"), paste(factors[1],"XLOAD-p2", sep="_"),
|
|
76 paste(factors[2],"XLOAD-p1", sep="_"), paste(factors[2],"XLOAD-p2", sep="_"),
|
|
77 "Interact_XLOAD-p1", "Interact_XLOAD-p2")
|
|
78
|
|
79 l <- list(data.asca[[1]], data.asca.permutation, ssq, samplemetadata, variablemetadata)
|
|
80 names(l) <- c("ASCA","p-values", "ssq", "samplemetadata", "variablemetadata")
|
|
81 return(l)
|
|
82 }
|