### R code from vignette source 'ALDEx2_vignette.Rnw'

###################################################
### code chunk number 1: ALDEx2_vignette.Rnw:49-53
###################################################
library(ALDEx2)
data(selex)
#subset for efficiency
selex <- selex[1201:1600,]


###################################################
### code chunk number 2: ALDEx2_vignette.Rnw:59-60
###################################################
conds <- c(rep("NS", 7), rep("S", 7))


###################################################
### code chunk number 3: ALDEx2_vignette.Rnw:65-70
###################################################
x <- aldex(selex, conds, mc.samples=16, test="t", effect=TRUE,
     include.sample.summary=FALSE, denom="iqlr", verbose=FALSE)

aldex.plot(x, type="MA", test="welch")
aldex.plot(x, type="MW", test="welch")


###################################################
### code chunk number 4: ALDEx2_vignette.Rnw:79-80
###################################################
x <- aldex.clr(selex, conds, mc.samples=16, denom="iqlr", verbose=TRUE)


###################################################
### code chunk number 5: ALDEx2_vignette.Rnw:87-88
###################################################
x.tt <- aldex.ttest(x, conds, paired.test=TRUE)


###################################################
### code chunk number 6: ALDEx2_vignette.Rnw:92-93 (eval = FALSE)
###################################################
## x.glm <- aldex.glm(x, conds)


###################################################
### code chunk number 7: ALDEx2_vignette.Rnw:98-99
###################################################
x.effect <- aldex.effect(x, conds, include.sample.summary=FALSE, verbose=TRUE)


###################################################
### code chunk number 8: ALDEx2_vignette.Rnw:103-104
###################################################
x.all <- data.frame(x.tt,x.effect)


###################################################
### code chunk number 9: ALDEx2_vignette.Rnw:120-124
###################################################

par(mfrow=c(1,2))
aldex.plot(x.all, type="MA", test="welch")
aldex.plot(x.all, type="MW", test="welch")


###################################################
### code chunk number 10: ALDEx2_vignette.Rnw:181-197
###################################################

# identify which values are significant in both the t-test and glm tests
found.by.all <- which(x.all$we.eBH < 0.05 & x.all$wi.eBH < 0.05)

# identify which values are significant in fewer than all tests
found.by.one <- which(x.all$we.eBH < 0.05 | x.all$wi.eBH < 0.05)

# plot the within and between variation of the data
plot(x.all$diff.win, x.all$diff.btw, pch=19, cex=0.3, col=rgb(0,0,0,0.3),
 xlab="Difference within", ylab="Difference between")
points(x.all$diff.win[found.by.one], x.all$diff.btw[found.by.one], pch=19,
 cex=0.5, col=rgb(0,0,1,0.5))
points(x.all$diff.win[found.by.all], x.all$diff.btw[found.by.all], pch=19,
 cex=0.5, col=rgb(1,0,0,1))
abline(0,1,lty=2)
abline(0,-1,lty=2)


###################################################
### code chunk number 11: ALDEx2_vignette.Rnw:225-241
###################################################
# this is the function to retrieve the denominator
# running it as illustration only
iqlr <- aldex.set.mode(selex, conds, denom="iqlr")
# show the first 5 rows that will be used for the clr denominator
iqlr[[1]][1:5]

# use the iqlr features for the analysis
y <- aldex.clr(selex, conds, denom="iqlr")
y.e <- aldex.effect(y, conds)
y.t <- aldex.ttest(y, conds)
y.all <- data.frame(y.e, y.t)

# find significant features in each
sig.x <- x.all$wi.eBH < 0.05
sig.y <- y.all$wi.eBH < 0.05



###################################################
### code chunk number 12: ALDEx2_vignette.Rnw:247-268
###################################################

# plot it
par(mfrow=c(1,2))
plot(x.all$diff.win,x.all$diff.btw, pch=19, col=rgb(0,0,0,.3), cex=0.3, main="denom=all") 
points(x.all$diff.win[sig.x],x.all$diff.btw[sig.x], pch=19, col="red", cex=0.4) 
points(x.all$diff.win[iqlr[[1]]],x.all$diff.btw[iqlr[[1]]], col=rgb(0,0,1,.7), cex=0.5) 
abline(0,0,lty=3, col="grey")
abline(0,1,lty=2, col="grey")
abline(0,-1,lty=2, col="grey")
abline(0,2,lty=2)
abline(0,-2,lty=2)

plot(y.all$diff.win,y.all$diff.btw, pch=19, col=rgb(0,0,0,.3), cex=0.3, main="denom=iqlr") 
points(y.all$diff.win[sig.y],y.all$diff.btw[sig.y], pch=19, col="red", cex=0.4) 
points(y.all$diff.win[iqlr[[1]]],y.all$diff.btw[iqlr[[1]]], col=rgb(0,0,1,.7), cex=0.5) 
abline(0,0,lty=3, col="grey")
abline(0,1,lty=2, col="grey")
abline(0,-1,lty=2, col="grey")
abline(0,2,lty=2)
abline(0,-2,lty=2)