Initial commit

.Rbuildignore

+^.*\.Rproj$
+^\.Rproj\.user$

.gitignore

+.Rproj.user
+.Rhistory
+.RData
+.Ruserdata

DESCRIPTION

+Package: ROBITools2
+Type: Package
+Title: What the Package Does (Title Case)
+Version: 0.1.0
+Author: Who wrote it
+Maintainer: The package maintainer <yourself@somewhere.net>
+Description: More about what it does (maybe more than one line)
+    Use four spaces when indenting paragraphs within the Description.
+License: What license is it under?
+Encoding: UTF-8
+LazyData: true
+RoxygenNote: 6.0.1
+Imports: tidyverse
+Suggests: vegan,
+          roxygen2,
+          knitr,
+          rmarkdown
+Collate: 
+    'color.R'
+    'utils.R'
+    'robiobject.R'
+    'dplyr.R'
+    'fasta.R'
+    'robimetabar.R'
+    'taxonomy.R'
+    'vegan.R'

NAMESPACE

+# Generated by roxygen2: do not edit by hand
+
+S3method("$",robimetabar)
+S3method("$<-",robimetabar)
+S3method("[",robiobject)
+S3method("[[",robimetabar)
+S3method("[[<-",robimetabar)
+S3method("motus<-",robimetabar)
+S3method("samples<-",robimetabar)
+S3method(as_tibble,robiobject)
+S3method(decostand,"NULL")
+S3method(decostand,default)
+S3method(decostand,robimetabar)
+S3method(dim,robimetabar)
+S3method(dimnames,robidata)
+S3method(dimnames,robimetabar)
+S3method(filter,robiobject)
+S3method(motus,robimetabar)
+S3method(names,robimetabar)
+S3method(samples,robimetabar)
+export("motus<-")
+export("samples<-")
+export(Fasta.reader)
+export(build.robimetabar)
+export(colors.tol)
+export(decostand)
+export(ecofind)
+export(family)
+export(genus)
+export(gradient.tol1)
+export(gradient.tol2)
+export(is.robiobject)
+export(is.subcladeof)
+export(kingdom)
+export(lowest.common.ancestor)
+export(make.robiobject)
+export(motus)
+export(ntaxa)
+export(parent)
+export(rainbow.tol)
+export(rank.list)
+export(read.fasta)
+export(robiclass)
+export(samples)
+export(scientificname)
+export(species)
+export(superkingdom)
+export(taxid.list)
+export(taxid.max)
+export(taxid.min)
+export(taxid.validate)
+export(taxon.at.rank)
+export(taxon.path)
+export(taxonomic.rank)
+export(unmake.robiobject)
+import(dplyr)
+import(stringr)
+import(tibble)

R/color.R

+#' Qualitative color schemes by Paul Tol
+#'
+#' Provide a set of color palette covenient even for blind color people
+#'
+#' @param n the minimal number of color you are interested in
+#' @param alpha the transparency value between $0$ full transparency
+#'        and $1$ for opaque colors.
+#' @return a character vector containing hexadecimal color codes
+#'
+#' @examples
+#'     x = colors.tol(10,0.5)
+#'     x
+#'     plot(1:10,1:10,col=x,cex=1:10,pch=16)
+#'
+#' @see <http://www.sron.nl/~pault/>
+#' @export
+colors.tol = function(n,alpha=1) {
+
+  color.tol = list(
+    "1" = c("#4477AA"),
+
+    "2" = c("#4477AA", "#CC6677"),
+
+    "3" = c("#4477AA", "#DDCC77", "#CC6677"),
+
+    "4" = c("#4477AA", "#117733", "#DDCC77", "#CC6677"),
+
+    "5" = c("#332288", "#88CCEE", "#117733", "#DDCC77",
+            "#CC6677"),
+
+    "6" = c("#332288", "#88CCEE", "#117733", "#DDCC77",
+            "#CC6677","#AA4499"),
+
+    "7" = c("#332288", "#88CCEE", "#44AA99", "#117733",
+            "#DDCC77", "#CC6677","#AA4499"),
+
+    "8" = c("#332288", "#88CCEE", "#44AA99", "#117733",
+            "#999933", "#DDCC77", "#CC6677","#AA4499"),
+
+    "9" = c("#332288", "#88CCEE", "#44AA99", "#117733",
+            "#999933", "#DDCC77", "#CC6677", "#882255",
+            "#AA4499"),
+
+    "10" = c("#332288", "#88CCEE", "#44AA99", "#117733",
+             "#999933", "#DDCC77", "#661100", "#CC6677",
+             "#882255", "#AA4499"),
+
+    "11" = c("#332288", "#6699CC", "#88CCEE", "#44AA99",
+             "#117733", "#999933", "#DDCC77", "#661100",
+             "#CC6677", "#882255", "#AA4499"),
+
+    "12" = c("#332288", "#6699CC", "#88CCEE", "#44AA99",
+             "#117733", "#999933", "#DDCC77", "#661100",
+             "#CC6677", "#AA4466", "#882255", "#AA4499"),
+
+    "14" = c("#882E72", "#B178A6", "#D6C1DE", "#1965B0",
+             "#5289C7", "#7BAFDE", "#4EB265", "#90C987",
+             "#CAE0AB", "#F7EE55", "#F6C141", "#F1932D",
+             "#E8601C", "#DC050C"),
+
+    "15" = c("#114477", "#4477AA", "#77AADD", "#117755",
+             "#44AA88", "#99CCBB", "#777711", "#AAAA44",
+             "#DDDD77", "#771111", "#AA4444", "#DD7777",
+             "#771144", "#AA4477", "#DD77AA"),
+
+    "18" = c("#771155", "#AA4488", "#CC99BB", "#114477",
+             "#4477AA", "#77AADD", "#117777", "#44AAAA",
+             "#77CCCC", "#777711", "#AAAA44", "#DDDD77",
+             "#774411", "#AA7744", "#DDAA77", "#771122",
+             "#AA4455", "#DD7788"),
+
+    "21" = c("#771155", "#AA4488", "#CC99BB", "#114477",
+             "#4477AA", "#77AADD", "#117777", "#44AAAA",
+             "#77CCCC", "#117744", "#44AA77", "#88CCAA",
+             "#777711", "#AAAA44", "#DDDD77", "#774411",
+             "#AA7744", "#DDAA77", "#771122", "#AA4455",
+             "#DD7788")
+  )
+
+  if (n > 21)
+    stop("Only color schemas with at maximum 21 colors are handled")
+
+  limits = as.integer(names(color.tol))
+
+  n.tol = as.character(limits[n <= limits][1])
+  tol   = color.tol[[n.tol]]
+
+  hex.alpha = toupper(as.hexmode(floor(alpha * 255)))
+  tol.alpha = paste(tol,hex.alpha,sep="")
+
+  return(tol.alpha)
+}
+
+
+#' Colors for variations around 0 (blue-red)
+#' @see <http://www.sron.nl/~pault/>
+#' @export
+
+gradient.tol1 = function(x,alpha=1) {
+  rcol = function(x) 0.237 - 2.13*x + 26.92*x**2 - 65.5*x**3 + 63.5*x**4 - 22.36*x**5
+  gcol = function(x) ((0.572 + 1.524*x - 1.811*x**2)/(1 - 0.291*x + 0.1574*x**2))**2
+  bcol = function(x) 1/(1.579 - 4.03*x + 12.92*x**2 - 31.4*x**3 + 48.6*x**4 - 23.36*x**5)
+
+  return(rgb(rcol(x),gcol(x),bcol(x),alpha))
+}
+
+#' Gradient colors,
+#' @see <http://www.sron.nl/~pault/>
+#'
+#' @export
+
+gradient.tol2 = function(x,alpha=1) {
+  rcol = function(x) (1 - 0.392*(1 + erf((x - 0.869)/ 0.255)))
+  gcol = function(x) (1.021 - 0.456*(1 + erf((x - 0.527)/ 0.376)))
+  bcol = function(x) (1 - 0.493*(1 + erf((x - 0.272)/ 0.309)))
+
+  return(rgb(rcol(x),gcol(x),bcol(x),alpha))
+}
+
+#' Rainbow gradient colors,
+#' @see <http://www.sron.nl/~pault/>
+#'
+#' @export
+
+rainbow.tol = function(x,alpha=1) {
+  rcol = function(x) (0.472-0.567*x+4.05*x**2)/(1.+8.72*x-19.17*x**2+14.1*x**3)
+  gcol = function(x) 0.108932-1.22635*x+27.284*x**2-98.577*x**3+163.3*x**4-131.395*x**5+40.634*x**6
+  bcol = function(x) 1./(1.97+3.54*x-68.5*x**2+243*x**3-297*x**4+125*x**5)
+
+  return(rgb(rcol(x),gcol(x),bcol(x),alpha))
+}

R/dplyr.R

+#' @import dplyr
+#' @include robiobject.R
+NULL
+
+#' @author Eric Coissac
+#' @export
+filter.robiobject = function(.data, ...) {
+  attrib=attributes(.data)
+  message("filter.robiobject : ",
+          deparse(attrib))
+  rep = NextMethod("filter")
+
+  return(.set_robiattribs(rep,
+                          attrib = attrib,
+                          set.class=TRUE))
+
+}
\ No newline at end of file

R/fasta.R

+#' @include robiseq.R
+#' @import stringr
+#'
+NULL
+
+#' @author Eric Coissac
+#' @export
+Fasta.reader = function(con,chunk.size=100000) {
+  previous.sequence = NULL
+  previous.header   = NULL
+  headers     = character(0)
+  ids         = character(0)
+  sequences   = character(0)
+  finnished   = FALSE
+  current.seq = 0
+
+  read.chunk = function() {
+
+    # I read chunk.size lines from the connection
+    buffer.lines = readLines(con = con,n = chunk.size)
+    n.lines      = length(buffer.lines)
+
+    # If I actually read less than chunk.size this means that
+    # the file is ended
+    finnished <<- n.lines < chunk.size
+
+    # I match the header lines
+    header.idx = stringr::str_detect(buffer.lines,
+                                     "^>")
+
+    seq.idx     = 0
+    seq.factor  = integer(n.lines)
+
+    for (i in seq_len(n.lines)){
+      if (header.idx[i]) seq.idx = seq.idx + 1
+      seq.factor[i]=seq.idx
+    }
+
+    # I put the header line on one side
+    header.lines = buffer.lines[header.idx]
+
+    # The sequence lines on the other
+    header.idx = ! header.idx
+    sequence.lines =  buffer.lines[header.idx]
+
+    # Sequence lines bellonging the same sequence are pasted
+    seq.factor     =  seq.factor[header.idx]
+    sequences <<- tapply(X = sequence.lines,
+                         INDEX = seq.factor,
+                         FUN = paste,
+                         collapse="",
+                         simplify = TRUE)
+
+    # We have already read the beginning of the sequence
+    if (! is.null(previous.sequence) && seq.factor[1]==0) {
+      sequences[1] <<- paste(previous.sequence,sequences[1])
+      header.lines = append(header.lines,
+                            previous.header,
+                            after = 0)
+    }
+
+    if (! finnished) {
+      # The file is not finnished, I save the last sequence
+      n = length(sequences)
+      previous.header   <<- header.lines[n]
+      previous.sequence <<- sequences[n]
+    }
+    else {
+      # The file is finnished, store NULL
+      previous.header   <<- NULL
+      previous.sequence <<- NULL
+    }
+
+    # Divide the header in ids and titles
+    header.lines = stringr::str_split(string = header.lines,
+                                      pattern = " ",
+                                      n = 2,simplify = TRUE)
+
+    ids     <<- str_sub(header.lines[,1],2,-1)
+    headers <<- header.lines[,2]
+
+    current.seq <<- 0
+  }
+
+  fasta.iterator = function() {
+    if (current.seq==length(headers)) {
+      if (finnished)
+        return(NULL)
+      else
+        read.chunk()
+    }
+
+    current.seq <<- current.seq + 1
+
+    return(list(id=ids[current.seq],
+                header=headers[current.seq],
+                sequence=sequences[current.seq]))
+  }
+
+  if (is.character(con))
+    con = file(con[1])
+
+  return(fasta.iterator)
+}
+
+#' @author Eric Coissac
+#' @export
+read.fasta = function(con,n=-1L,chunk.size=100000) {
+  previous.lines = character(0)
+  continue = TRUE
+
+  while(continue) {
+    str_detect(string, pattern)
+  }
+
+}
+

R/robimetabar.R

+#' @import dplyr
+#' @import tibble
+NULL
+
+
+#' @author Eric Coissac
+#' @export
+build.robimetabar = function(reads,samples,motus,
+                             layers = vector(mode = "list",0),
+                             distances = vector(mode = "list",0),
+                             motu.ids="id",
+                             sample.ids="id",
+                             verbose = TRUE,
+                             check.dims = TRUE
+                 ) {
+
+  # Reads must be a matrix of something cohecible
+  # to a matrix
+
+  if (! is.matrix(reads))
+    reads = as.matrix(reads)
+
+  #
+  # If no motus table provided, a fake one is
+  # created from the col names of reads or
+  # col positions if no col names are available
+  #
+
+  if (missing(motus)) {
+    if (verbose) warning("missing MOTU descriptions")
+    motu.names = colnames(reads)
+    if (is.null(motu.names)){
+      if (verbose) warning("no column names on the reads matrix")
+      motu.names = paste0("motu_",seq_len(ncol(reads)))
+    }
+    motus = tibble(id=motu.names)
+    motu.ids = "id"
+  }
+
+  #
+  # If no samples table provided, a fake one is
+  # created from the row names of reads or
+  # row positions if no row names are available
+  #
+
+  if (missing(samples)) {
+    if (verbose) warning("missing sample descriptions")
+    sample.names = rownames(reads)
+    if (is.null(sample.names)){
+      if (verbose) warning("no row names on the reads matrix")
+      sample.names = paste0("sample_",seq_len(nrow(reads)))
+    }
+    samples = tibble(id=sample.names)
+    sample.ids = "id"
+  }
+
+  #
+  # motus and samples must be cohercible to tibble
+  #
+
+  motus = as.tibble(motus)
+  samples = as.tibble(samples)
+
+  #
+  # If no name of column is provided for motu  ids
+  # send a warning
+  #
+  if (missing(motu.ids)) {
+    if (verbose) warning("missing MOTU ids column name use the default id column")
+  }
+
+  #
+  # Checks if the motu.ids column exists in motus
+  # otherwise solve the problem by using colnames
+  # or column index
+  #
+  motu.names = motus[[motu.ids]]
+  if (is.null(motu.names)) {
+    if (verbose) warning("no MOTU id column in the motu descriptions")
+    motu.names = colnames(reads)
+    if (is.null(motu.names)) {
+      if (verbose) warning("no column names on the reads matrix")
+      motu.names = paste0("motu_",seq_len(ncol(reads)))
+      motus[[motus.ids]]=motu.names
+    }
+  }
+
+  #
+  # If no name of column is provided for sample ids
+  # send a warning
+  #
+  if (missing(sample.ids)) {
+    if (verbose)
+      warning("missing MOTU ids column name use the default id column")
+  }
+
+  #
+  # Checks if the sample.ids column exists in samples
+  # otherwise solve the problem by using rownames
+  # or row index
+  #
+  sample.names = samples[[sample.ids]]
+  if (is.null(motu.names)) {
+    if (verbose) warning("no sample ids column in the sample descriptions")
+    sample.names = rownames(reads)
+    if (is.null(sample.names)) {
+      if (verbose) warning("no row names on the reads matrix")
+      sample.names = paste0("sample_",seq_len(nrow(reads)))
+      samples[[sample.ids]]=sample.names
+    }
+  }
+
+  rownames(reads) = samples[[sample.ids]]
+  colnames(reads) = motus[[motu.ids]]
+
+  metabar = new.env(parent = globalenv(),size = 20)
+  make.robiobject(metabar,robiclass = "robimetabar")
+
+  assign(x = "reads",value = reads,
+         envir = metabar)
+  assign(x = "samples",value = samples,
+         envir = metabar)
+  assign(x = "motus",value = motus,
+         envir = metabar)
+  assign(x = "layers",value = layers,
+         envir = metabar)
+  assign(x = "distances",value = distances,
+         envir = metabar)
+  assign(x = "motu.ids",value = motu.ids,
+         envir = metabar)
+  assign(x = "sample.ids",value = sample.ids,
+         envir = metabar)
+
+  return(metabar)
+
+}
+
+
+update.robimetabar = function(metabar,
+                              reads,
+                              samples,motus,
+                              layers,
+                              distances,
+                              motu.ids,
+                              sample.ids,
+                              verbose = TRUE,
+                              check.dims=TRUE
+                 ) {
+
+  if (! missing(reads))
+    assign(x = "reads",value = reads,
+           envir = metabar)
+  if (! missing(samples))
+    assign(x = "samples",value = samples,
+           envir = metabar)
+  if (! missing(motus))
+    assign(x = "motus",value = motus,
+           envir = metabar)
+  if (! missing(layers))
+    assign(x = "layers",value = layers,
+           envir = metabar)
+  if (! missing(distances))
+    assign(x = "distances",value = distances,
+           envir = metabar)
+  if (! missing(motu.ids))
+    assign(x = "motu.ids",value = motu.ids,
+           envir = metabar)
+  if (! missing(sample.ids))
+    assign(x = "sample.ids",value = sample.ids,
+           envir = metabar)
+
+  return(metabar)
+}
+
+copy.robimetabar = function(metabar,
+                              reads,
+                              samples,motus,
+                              layers,
+                              distances,
+                              motu.ids,
+                              sample.ids,
+                              verbose = TRUE,
+                              check.dims=TRUE
+                 ) {
+  if (missing(reads))
+    reads = get("reads",metabar)
+  if (missing(samples))
+    samples = get("samples",metabar)
+  if (missing(motus))
+    motus = get("motus",metabar)
+  if (missing(layers))
+    layers = get("layers",metabar)
+  if (missing(distances))
+    distances = get("distances",metabar)
+  if (missing(motu.ids))
+    motu.ids = get("motu.ids",metabar)
+  if (missing(sample.ids))
+    sample.ids = get("sample.ids",metabar)
+
+  build.robimetabar(reads = reads,
+                    samples = samples,
+                    motus = motus,
+                    layers = layers,
+                    distances = distances,
+                    motu.ids = motu.ids,
+                    sample.ids = sample.ids,
+                    verbose = verbose,
+                    check.dims = check.dims)
+
+}
+
+
+#' @author Eric Coissac
+#' @export
+dim.robimetabar = function(x)
+{eval(quote(dim(reads)),envir = x)}
+
+#' @author Eric Coissac
+#' @export
+dimnames.robimetabar = function(x) {
+  list(samples = eval(quote(samples[[sample.ids]]),envir = x),
+       motus   = eval(quote(motus[[motu.ids]]),envir = x))
+}
+
+#' @author Eric Coissac
+#' @export
+`[[.robimetabar` = function(x, ...) {
+
+  args = list(...)
+  if (length(args) > 1)
+    warning("Only first argument passed to [[ will be used")
+
+  layer = args[[1]]
+
+  if (length(layer) > 1)
+    error("Only a single layer can be returned by [[")
+
+  if (is.numeric(layer)) {
+    if (layer==1)
+      return(eval(quote(reads),envir = x))
+
+    layer = eval(quote(layers[[layer-1]]),
+                 envir = list(layer=layer),