---
title: "Example 16S Annotation Workflow"
author: "Nathan D. Olson"
date: "`r Sys.Date()`"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Vignette Title}
  %\VignetteEngine{knitr::rmarkdown}
  \usepackage[utf8]{inputenc}
---

## Overview
The metagenomeFeatures package and associated annotation packages (e.g. greengene13.5MgDb) can be used to and taxonomic annotations to MRexperiment objects.  The metagenomeSeq package includes the MRexperiment-class, and has a number of methods for analyzing marker-gene metagenome datasets.  This vignette demonstrate how to add taxonomic annotation to a MRexperiment-class object using metagenomeFeatures.   

To generate a MRexperiment object with taxonomic information you need;  

* assayData - matrix with count data for each OTU and sample  
* phenoData - an AnnotatedDataFrame with sample metadata  
* featureData - a metagenomeAnnotation object with taxonomic information for each OTU  

For MRexperiment objects, the assayData and phenoData need to have the sample sample order similarly assayData and featureData need to have the sample OTU order.

__Required Packages__
```{r message=FALSE}
library(metagenomeFeatures)
library(metagenomeSeq)
```

## Loading assayData
A file with count data from the msd16s study is included in the metagenomeFeatures package.
```{r}
dataDirectory <- system.file("extdata", package = "metagenomeFeatures")
assay_data <- load_meta(file.path(dataDirectory, "msd16s_counts.csv"),sep = ",")
```

## Loading phenoData
```{r}
pheno_data <- AnnotatedDataFrame(read.csv(file.path(dataDirectory, 
                                                    "msd16s_sample_data.csv"), 
                                          row.names = 1, stringsAsFactors = F))
```

## featureData
### Generating a metagenomeAnnotation-class object
The `msd16s_query_df` is a dataframe with greengenes taxonomic ids for the msd16s data.
```{r}
data("msd16s_query_df")
head(msd16s_query_df)
```

A subset of the greengenes database, `msd16s_MgDb` is used to generate a MgDb-object for annotating the msd16s data. This database is only used here for demonstration purposes, the `gg13.5MgDb` database in the greengenes13.5MgDb package can also be used to generate the metagenomeAnnotation object.
```{r echo=FALSE, message=FALSE}
library(Biostrings)
metadata <- list(ACCESSION_DATE = "3/31/2015",
                 URL = "https://greengenes.microbio.me",
                 DB_TYPE_NAME = "GreenGenes-MgDb-msd16s",
                 DB_TYPE_VALUE = "MgDb",
                 DB_SCHEMA_VERSION = "1.0")
msd16s_MgDb <- new("MgDb",
                seq = readDNAStringSet(file.path(dataDirectory,"msd16s_MgDb_seq.fasta.gz")),
                taxa = file.path(dataDirectory,"msd16s_MgDb_taxa.sqlite"),
                metadata = metadata)
```

```{r}
msd16s_feature_data <- annotate(msd16s_MgDb, query_df = msd16s_query_df)

## Using the greengenes13.5MgDb annotation package
# library(greengenes13.5MgDb)
# msd16s_feature_data <- annotate(gg13.5MgDb, query_df = msd16s_query_df)
```

### Adding featureNames and ordering features
Prepare the metagneomeAnnotation object for use in a MRexperiment by defining featureName and ordering the OTUs to that they agree with the assayData. Note that the OTU column in `msd16s_feature_data` are factors and are first converted to numeric before ordering.
```{r}
OTU_numbers <- as.numeric(as.character(msd16s_feature_data$OTU))
featureNames(msd16s_feature_data) <- OTU_numbers
ord <-  sort(OTU_numbers)
msd16s_feature_data <- msd16s_feature_data[order(OTU_numbers),]
```


## Creating MRexperiment object
```{r}
obj = newMRexperiment(assay_data$counts,
                      phenoData=pheno_data,
                      featureData = msd16s_feature_data)
obj
```




```{r}
sessionInfo()
```
