calculate_ks or the Earth Mover's Distance algorithm using calculate_emd.The algorithm is used to compare genomics data between any number of groups. Usually the data will be gene expression values from array-based or sequence-based experiments, but data from other types of experiments can also be analyzed (e.g. copy number variation).
Traditional methods like Significance Analysis of Microarrays (SAM) and Linear Models for Microarray Data (LIMMA) use significance tests based on summary statistics (mean and standard deviation) of the two distributions. This approach tends to give non-significant results if the two distributions are highly heterogeneous, which can be the case in many biological circumstances (e.g sensitive vs. resistant tumor samples).
The Cramer von Mises (CVM) algorithm generates a test statistic that is the sum of the squared values of the differences between two cumulative distribution functions (CDFs). As a result, the test statistic tends to overestimate the similarity between two distributions and cannot effectively handle partial matching like EMD does. However, it is one of the most commonly referenced nonparametric two-class distribution comparison tests in non-genomic contexts.
The CVM-based algorithm implemented in EMDomics has two main steps. First, a matrix (e.g. of expression data) is divided into data for each of the groups. Every possible pairwise CVM score is then computed and stored in a table. The CVM score for a single gene is calculated by averaging all of the pairwise CVM scores. Next, the labels for each of the groups are randomly permuted a specified number of times, and an CVM score for each permutation is calculated. The median of the permuted scores for each gene is used as the null distribution, and the False Discovery Rate (FDR) is computed for a range of permissive to restrictive significance thresholds. The threshold that minimizes the FDR is defined as the q-value, and is used to interpret the significance of the CVM score analogously to a p-value (e.g. q-value < 0.05 is significant.)
calculate_cvm(data, outcomes, nperm = 100, pairwise.p = FALSE, seq = FALSE, quantile.norm = FALSE, verbose = TRUE, parallel = TRUE)data matrix. The names should be the sample identifiers provided in data.FALSE.TRUE, if passing transcripts per million (TPM) data or raw
data that is not scaled. If TRUE, data will be normalized by first multiplying by 1E6, then adding
1, then taking the log base 2. If FALSE, the data will be handled as is (unless
quantile.norm is TRUE). Note that as a distribution comparison function, K-S will
compute faster with scaled data. Defaults to FALSE.TRUE, then the normalize.quantiles function is used.
Defaults to FALSE.TRUE.CVMomics object.
CVMomics CramerVonMisesTwoSamples
# 100 genes, 100 samples
dat <- matrix(rnorm(10000), nrow=100, ncol=100)
rownames(dat) <- paste("gene", 1:100, sep="")
colnames(dat) <- paste("sample", 1:100, sep="")
# "A": first 50 samples; "B": next 30 samples; "C": final 20 samples
outcomes <- c(rep("A",50), rep("B",30), rep("C",20))
names(outcomes) <- colnames(dat)
results <- calculate_cvm(dat, outcomes, nperm=10, parallel=FALSE)
head(results$cvm)
Run the code above in your browser using DataLab