Nothing
#' Direction Analysis for Pathways
#'
#' The main function of direction Analysis. This function takes in a matrix of test statistics with
#' two (2-dimensional space) or three (3-dimensional space) columns, the direction of interests, and
#' the annotation list such as pathway annotation, and test for enrichment of pathways on the specified
#' direction.
#'
#' @usage directPA(Tc, direction, annotation, minSize=5, gene.method="OSP",
#' path.method="Stouffer", visualize=TRUE, ...)
#'
#' @param Tc a numeric matrix. Rows are genes and columns are treatments vs control statistics.
#' @param direction the direction to be tested for enrichment. Either specified as a degree for
#' two-dimensional analysis or as contrast (in a triplet) for three-dimensional analysis.
#' @param annotation a list with names correspond to pathways and elements correspond to genes belong to
#' each pathway, respectively.
#' @param minSize the size of annotation groups to be considered for calculating enrichment. Groups
#' that are smaller than the minSize will be removed from the analysis.
#' @param gene.method the method to be used for integrating statistics across treatments for each gene.
#' Available methods are Stouffer, OSP, Fisher, and maxP.
#' Default method is OSP.
#' @param path.method the method to be used for integrating statistics of all genes that belongs to a
#' pathway. Available methods are Stouffer, OSP, Fisher, and maxP.
#' Default method is Stouffer.
#' @param visualize whether to visualize the plot.
#' @param ... other visualization parameters to pass on.
#' @return a list that contains directional p-values for each gene and directional enrichment for each pathway.
#' @export
#' @examples
#'
#' # load the proteomics dataset
#' data(PM)
#'
#' # load pathway annotations
#' data(Pathways)
#'
#' # display reactome pathways. Could be replaced by any other pathway databases
#' Pathways.reactome[1:5]
#'
#' # direction pathway analysis in 3-dimensional space. Implemnted as rotating by contrast
#' # (1) test combined effect of all 3 treatments (stimulation and inhibitions) vs control (basal)
#' # on the original direction.
#' dPA <- directPA(Tc=PM, direction=c(1,1,1), annotation=Pathways.reactome)
#' dPA$gst[order(unlist(dPA$gst[,1])),][1:20,]
#' # rank substrates on the direciton of interest
#' sort(dPA$gene.pvalues)[1:20]
#'
#' # (2) test combined effect of all 3 treatments vs controls on direction c(1,-1, 0)
#' # this rotates Ins by 0 degree, Wmn by 90 degree, and MK by 45 degree.
#' dPA <- directPA(Tc=PM, direction=c(1,-1,0), annotation=Pathways.reactome)
#' dPA$gst[order(unlist(dPA$gst[,1])),][1:20,]
#'
#' # (3) test combined effect of all 3 perturbations vs controls on direction c(1,-1, 1)
#' # this rotates Ins by 0 degree, Wmn by 90 degree, and MK by 0 degree.
#' dPA <- directPA(Tc=PM, direction=c(1,-1,1), annotation=Pathways.reactome)
#' dPA$gst[order(unlist(dPA$gst[,1])),][1:20,]
#'
directPA <- function(Tc, direction, annotation, minSize=5, gene.method="OSP", path.method="Stouffer", visualize=TRUE, ...){
## spherical coordinates for three-dimensional rotation
if (length(direction) == 3) {
# step 1. convert statistics into z-scores
Tc.zscores <- apply(Tc, 2, function(x){qnorm(rank(x)/(nrow(Tc)+1))})
# step 2. rotate z-scores
Tc.rotated <- rotate3d(Tc.zscores, direction)
# step 3. integrate statistics across treatments
gene.pvalues <- apply(Tc.rotated, 1, geneStats, gene.method)
if (visualize == TRUE) {
#HC = rainbow(length(gene.pvalues)*1.2)
#plot3d(Tc, col=HC[rank(gene.pvalues)], size=5, ...)
#abclines3d(x=0, y=0, z=0, a=diag(3), col="black", lwd=3)
#abclines3d(x=0, a=direction, col="pink", lwd=5)
df <- data.frame(Tc)
colnames(df) <- c("x", "y", "z")
df$pvalue <- gene.pvalues[rownames(df)]
my_col = colorRampPalette(rainbow(12))(100)
p <- plotly::plot_ly(df, x=~x, y=~y, z=~z, color=~pvalue, colors=my_col, size=5)
p <- plotly::add_markers(p)
p <- plotly::layout(p, scene = list(xaxis = list(title = colnames(Tc)[[1]]),
yaxis = list(title = colnames(Tc)[[2]]),
zaxis = list(title = colnames(Tc)[[3]])))
print(p)
}
# step 4. integrate statistics for pathways
gene.zscores <- qnorm(gene.pvalues, lower.tail = FALSE)
gst <- t(sapply(annotation, pathwayStats, gene.zscores, minSize=5, path.method))
result <- list()
result$gene.pvalues <- gene.pvalues
result$gst <- gst
return(result)
}
## polar coordinates for two-dimensional rotation
if (length(direction) == 1) {
# step 1. convert statistics into z-scores
Tc.zscores <- apply(Tc, 2, function(x){qnorm(rank(x)/(nrow(Tc)+1))})
# step 2. rotate z-scores
Tc.rotated <- rotate2d(Tc.zscores, direction)
# step 3. integrate statistics across treatments
gene.pvalues <- apply(Tc.rotated, 1, geneStats, gene.method)
if (visualize == TRUE) {
HC = rainbow(length(gene.pvalues)*1.2)
plot(Tc, col=HC[rank(gene.pvalues)], pch=16, ...)
abline(v = 0,h = 0, lty=2, col="gold")
abline(a=0, b=1, col="darkgreen", lty=2)
abline(a=0, b=-1, col="darkgreen", lty=2)
}
# step 4. integrate statistics for pathways
gene.zscores <- qnorm(gene.pvalues, lower.tail = FALSE)
gst <- t(sapply(annotation, pathwayStats, gene.zscores, minSize=5, path.method))
result <- list()
result$gene.pvalues <- sort(gene.pvalues)
result$pathways <- gst[order(as.numeric(gst[,1])),]
return(result)
}
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.