tests/testMatrixOps.R

library(oompaBase)
# a key point is to test the matrixPairedT code
# we create a simple example
suppressWarnings( RNGversion("3.5.3") )
set.seed(372284)
nPairs <- 3
nGenes <- 7
m <- matrix(rnorm(2*nPairs*nGenes), ncol=2*nPairs)
v <- factor(rep(c("A","B"), each=nPairs))
pf <- rep(1:nPairs, 2)
colnames(m) <- paste(as.character(v), pf, sep='')
round(m, 2)
# Now we run the package code
mpt <- matrixPairedT(m, v, pf)
# We also loop over the standard t.test function
# for each gene
realt <- sapply(1:nGenes, function(i) {
  x <- m[i, v=="A"]
  y <- m[i, v=="B"]
  t1 <- t.test(y, x, paired=TRUE)
  t1$statistic
})
# The differences should all be less than 1e-15
all(round(mpt-realt, 15) == 0)

# for completeness, we also check the unpaired result

mtt <- matrixT(m, v)
realt <- sapply(1:nGenes, function(i) {
  x <- m[i, v=="A"]
  y <- m[i, v=="B"]
  t1 <- t.test(x, y, paired=FALSE, var.equal=TRUE)
  t1$statistic
})
# The differences should all be less than 1e-15
all(round(mtt-realt, 15) == 0)

# and the unequal variance case
mut <- matrixUnequalT(m,v)
realt <- sapply(1:nGenes, function(i) {
  x <- m[i, v=="A"]
  y <- m[i, v=="B"]
  t1 <- t.test(x, y, paired=FALSE, var.equal=FALSE)
  t1$statistic
})
# The differences should all be less than 1e-15
all(round(mut$tt-realt, 15) == 0)

Try the oompaBase package in your browser

Any scripts or data that you put into this service are public.

oompaBase documentation built on Aug. 24, 2019, 5:02 p.m.