New version!

DESCRIPTION

@@ -1,10 +1,10 @@
 Package: DAtest
 Title: Comparing Differential Abundance methods
-Version: 0.0.0.9000
+Version: 1.0.0
 Authors@R: person("Jakob", "Russel", email = "russel2620@gmail.com", role = c("aut", "cre"))
 Description: What the title says.
 Depends: R (>= 3.2.5)
-Imports: foreach, doSNOW, snow, MASS, glmnet, pROC, ggplot2, cowplot, lme4
+Imports: foreach, doSNOW, snow, MASS, pROC, ggplot2, cowplot, lme4, nlme
 License: GPL (>= 3)
 Encoding: UTF-8
 LazyData: true

NAMESPACE

@@ -8,13 +8,14 @@ export(DA.anc)
 export(DA.aov)
 export(DA.bay)
 export(DA.ds2)
-export(DA.enn)
 export(DA.ere)
 export(DA.erq)
 export(DA.kru)
 export(DA.lao)
 export(DA.lao2)
 export(DA.lim)
+export(DA.lli)
+export(DA.lli2)
 export(DA.llm)
 export(DA.llm2)
 export(DA.lrm)
@@ -23,6 +24,8 @@ export(DA.ltt2)
 export(DA.msf)
 export(DA.neb)
 export(DA.per)
+export(DA.rai)
+export(DA.spe)
 export(DA.ttt)
 export(DA.wil)
 export(DA.zig)
@@ -34,7 +37,6 @@ import(cowplot)
 import(doSNOW)
 import(foreach)
 import(ggplot2)
-import(glmnet)
 import(lme4)
 import(nlme)
 import(pROC)

R/DA.adx.R

@@ -1,12 +1,15 @@
 #' Aldex t.test and wilcox
-
+#' 
+#' @param count_table Matrix or data.frame. Table with taxa/genes/proteins as rows and samples as columns
+#' @param outcome Factor. The outcome of interest. E.g. case and control
+#' @param ... Additional arguments for the aldex function
 #' @export
 
-DA.adx <- function(count_table, outcome, mc.samples = 128, p.adj){
+DA.adx <- function(count_table, outcome, ...){
 
   library(ALDEx2, quietly = TRUE)
 
-  x <- aldex(data.frame(count_table), outcome, mc.samples = mc.samples, verbose = FALSE)
+  x <- aldex(data.frame(count_table), outcome, verbose = FALSE, ...)
   x$Feature <- rownames(x)
   return(x)
 

R/DA.anc.R

@@ -1,15 +1,18 @@
 #' ANCOM
 #'
+#' @param count_table Matrix or data.frame. Table with taxa/genes/proteins as rows and samples as columns
+#' @param outcome Factor. The outcome of interest. E.g. case and control
+#' @param ... Additional arguments for the ANCOM function
 #' @export
 
-DA.anc <- function(count_table, outcome, sig = 0.05, multcorr = 3, tau = 0.02, theta = 0.1, repeated = FALSE){
+DA.anc <- function(count_table, outcome, ...){
 
   library(ancom.R, quietly = TRUE)
 
   input <- as.data.frame(t(count_table))
   input$Group <- outcome
 
-  res <- ANCOM(input, sig, multcorr, tau, theta, repeated) 
+  res <- ANCOM(input, ...) 
 
   df <- data.frame(Feature = rownames(count_table),
                    W = res[[1]],

R/DA.aov.R

@@ -1,11 +1,16 @@
 #' ANOVA
-
+#' 
+#' @param count_table Matrix or data.frame. Table with taxa/genes/proteins as rows and samples as columns
+#' @param outcome Factor. The outcome of interest. E.g. case and control
+#' @param relative Logical. Should count_table be normalized to relative abundances. Default TRUE
+#' @param p.adj Character. P-value adjustment. Default "fdr". See p.adjust for details
+#' @param ... Additional arguments for the aov function
 #' @export
 
-DA.aov <- function(count_table, outcome, p.adj, relative){
+DA.aov <- function(count_table, outcome, relative = TRUE, p.adj = "fdr", ...){
 
   ao <- function(x){
-    tryCatch(as.numeric(summary(aov(x ~ outcome))[[1]][1,5]), error = function(e){NA}) 
+    tryCatch(as.numeric(summary(aov(x ~ outcome, ...))[[1]][1,5]), error = function(e){NA}) 
   }
 
   if(relative){

R/DA.bay.R

@@ -1,11 +1,20 @@
 #' baySeq
-#'
-#' From
+#' 
+#' Implemented as in:
 #' https://microbiomejournal.biomedcentral.com/articles/10.1186/s40168-016-0208-8
-
+#' @param count_table Matrix or data.frame. Table with taxa/genes/proteins as rows and samples as columns
+#' @param outcome Factor. The outcome of interest. E.g. case and control
+#' @param p.adj Character. P-value adjustment. Default "fdr". See p.adjust for details
+#' @param samplesize How large a sample should be taken in estimating the priors? Default 1e5
+#' @param samplingSubset If given, the priors will be sampled only from the subset specified. Default NULL
+#' @param equalDispersions Should we assume equal dispersions of data across all groups in the 'cD' object? Defaults to TRUE
+#' @param estimation Defaults to "QL", indicating quasi-likelihood estimation of priors. Currently, the only other possibilities are "ML", a maximum-likelihood method, and "edgeR", the moderated dispersion estimates produced by the 'edgeR' package
+#' @param zeroML Should parameters from zero data (rows that within a group are all zeros) be estimated using maximum likelihood methods (which will result in zeros in the parameters?
+#' @param consensus If TRUE, creates a consensus distribution rather than a separate distribution for each member of the groups structure in the 'cD' object
+#' @param ... Additional arguments to the getLikelihoods function
 #' @export
 
-DA.bay <- function(count_table, outcome, p.adj){
+DA.bay <- function(count_table, outcome, p.adj = "fdr", samplesize = 1e5, samplingSubset = NULL, equalDispersions = TRUE, estimation = "QL", zeroML = FALSE, consensus = FALSE, ...){
 
   library(baySeq, quietly = TRUE)
 
@@ -14,8 +23,8 @@ DA.bay <- function(count_table, outcome, p.adj){
   libsizes(CD) <- getLibsizes(CD)
   CD@annotation <- data.frame(name=rownames(count_table))
 
-  CD <- getPriors.NB(CD, cl= NULL)
-  CD <- getLikelihoods(CD, cl = NULL)
+  CD <- getPriors.NB(CD, cl= NULL, samplesize=samplesize, samplingSubset=samplingSubset, equalDispersions=equalDispersions, estimation=estimation, zeroML=zeroML, consensus=consensus)
+  CD <- getLikelihoods(CD, cl = NULL, ...)
 
   tc <- topCounts(CD, group = "DE", number=nrow(count_table))
   tc <- tc[,c(1,rev(ncol(tc)-0:4))]

R/DA.ds2.R

@@ -1,12 +1,16 @@
 #' DESeq2
 #'
-#' From:
-#' https://microbiomejournal.biomedcentral.com/articles/10.1186/s40168-016-0208-8
+#' Implemented as in:
+#' https://microbiomejournal.biomedcentral.com/articles/10.1186/s40168-016-0208-8.
 #' Manual geometric means calculated to avoid errors, see https://github.com/joey711/phyloseq/issues/387
-
+#' @param count_table Matrix or data.frame. Table with taxa/genes/proteins as rows and samples as columns
+#' @param outcome Factor. The outcome of interest. E.g. case and control
+#' @param paired Factor. Subject ID for running paired analysis
+#' @param p.adj Character. P-value adjustment. Default "fdr". See p.adjust for details
+#' @param ... Additional arguments for the DESeq function
 #' @export
 
-DA.ds2 <- function(count_table, outcome, paired = NULL, p.adj){
+DA.ds2 <- function(count_table, outcome, paired = NULL, p.adj = "fdr", ...){
 
   library(DESeq2, quietly = TRUE)
   if(is.null(paired)){
@@ -25,7 +29,7 @@ DA.ds2 <- function(count_table, outcome, paired = NULL, p.adj){
   }
   geoMeans = apply(counts(x), 1, gm_mean)
   x = estimateSizeFactors(x, geoMeans = geoMeans)
-  x <- DESeq(x)
+  x <- DESeq(x, ...)
   res <- as.data.frame(results(x)@listData)
   colnames(res)[5] <- "pval"
   res$pval.adj <- p.adjust(res$pval, method = p.adj)

R/DA.ere.R

@@ -1,18 +1,21 @@
-#' EdgeR
+#' EdgeR exact test
 #'
-#' From:
+#' Implemented as in:
 #' https://microbiomejournal.biomedcentral.com/articles/10.1186/s40168-016-0208-8
-
+#' @param count_table Matrix or data.frame. Table with taxa/genes/proteins as rows and samples as columns
+#' @param outcome Factor. The outcome of interest. E.g. case and control
+#' @param p.adj Character. P-value adjustment. Default "fdr". See p.adjust for details
+#' @param ... Additional arguments for the exactTest function
 #' @export
 
-DA.ere <- function(count_table, outcome, p.adj){
+DA.ere <- function(count_table, outcome, p.adj = "fdr", ...){
 
   library(edgeR, quietly = TRUE)
   otu_table <- as.data.frame(count_table)
   x <- DGEList(counts = count_table, group = outcome, genes = data.frame(Feature = row.names(count_table)))
   x <- edgeR::calcNormFactors(x)
   x <- estimateTagwiseDisp(estimateCommonDisp(x))
-  ta <- exactTest(x)[[1]]
+  ta <- exactTest(x, ...)[[1]]
   colnames(ta)[3] <- "pval"
   ta$pval.adj <- p.adjust(ta$pval, method = p.adj)
   ta$Feature <- rownames(ta)