# Perform a 2 group comparison with SAM method
DeSam<-function(mtrx, grps, paired=FALSE, logged=TRUE, nperm=100, ...) {
# mtrx A numeric matrix of gene expression data. Rows are unique genes and columns include 2 groups of samples to be compared
# grps A 2-vector list, each vector has the column indexes (numeric vectors) or column names (character vectors) of a group; vectors are named by group names
# paired Paired test
# logged If the data has been log-transformed
# nperm SAM parameter
require(DEGandMore);
require(samr);
prepared <- PrepareDe(mtrx, grps, paired);
mtrx <- prepared[[1]];
grps <- prepared[[2]];
paired <- prepared[[3]];
e1 <- mtrx[, grps[[1]], drop=FALSE];
e2 <- mtrx[, grps[[2]], drop=FALSE];
if (paired) {
y <- c(-(1:ncol(e1)), 1:ncol(e2));
type <- "Two class paired";
} else {
y <- rep(c(1, 2), sapply(grps, length));
type <- "Two class unpaired";
}
data <- list(x=cbind(e1, e2), y=y, geneid=rownames(mtrx), genenames=rownames(mtrx), logged2=logged);
capture.output(sam<-samr(data, resp.type=type, nperms=nperm, ...))->x;
p <- samr.pvalues.from.perms(sam$tt, sam$ttstar);
capture.output(delta<-samr.compute.delta.table(sam))->x;
x <- samr.compute.siggenes.table(sam, 0, data, delta, all.genes=T);
x <- rbind(x[[1]], x[[2]]);
rownames(x) <- as.vector(x[,2]);
x <- x[rownames(mtrx), ];
fc <- as.numeric(x[,7]);
q <- as.numeric(x[,8])/100;
m1 <- rowMeans(e1);
m2 <- rowMeans(e2);
nm <- paste(names(grps), collapse='-');
if (logged) lgfc<-m2-m1 else lgfc<-log2(fc);
lgfc[is.na(lgfc)]<-0;
s <- cbind(m1, m2, m2-m1, lgfc, pmin(1, pmax(0, p)), pmin(1, pmax(0, q)));
colnames(s) <- c(paste('Mean', names(grps), sep='_'), 'Mean_Change', 'LogFC', 'Pvalue', 'FDR');
list(stat=s[rownames(mtrx), ], group=grps, sam=sam, sig=x)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.