Information

Version: 1.0.0
Title: Robust Phylogenetic Tree-based Microbiome Association Test using Repeatedly Measured Data for Composition Bias (mTMAT)
R package Dependancies: ggplot2 raster stringr ape grid mmmgee edgeR parallel mmmgee statVisual Biobase pvca msir

URL: http://healthstat.snu.ac.kr/software/mtmat/
Maintainer: Kangjin Kim <rekki@channing.harvard.edu>

Usage

Obtain mTMAT result

mTMAT(Data_Read,Meta_Table,Phenotype,ID,Time,Cov,Tree,Taxonomy_Table,Taxonomy_level="Genus",ncore=1,str="independence",CompositionalBias_correction="Default",mTMAT_type="IM",taxon_filter_zero_prop=0.25,Total_Readcount=NULL)

Input

  • Data_Read: data.frame of metagenome data whose rows corresponds to OTUs. Its first column indicates unique OTU ID corresponding to Tree object. The other columns has absolute abundances of OTUs for each samples. Column names should corresponds to first column of Meta_Table object. Example file: http://viva1109.duckdns.org/paper/mTMAT/temporal_spatial_otuTable.csv
  • Meta_Table: data.frame containing phenotype and covariates data whose rows corresponds to the samples. First column has to represent sample ID. First column should be a subset of column names of Data_Read object. Example file: http://viva1109.duckdns.org/paper/mTMAT/temporal_spatial_meta.csv
  • Phenotype: character object that represents phenotype of interest. This has to be cantained in the column names of Meta_Table object.
  • ID: character object that represents the subject on which data is repeatedly measured. This has to be cantained in the column names of Meta_Table object.
  • Time: character object that represents Time variable that represents the order of repeatedly measured data. This has to be cantained in the column names of Meta_Table object.
  • Cov: character vector that represents the names of covaraites to adjust. This has to be cantained in the column names of Meta_Table object.
  • Phenotype: Tree object whose tip labels corresponds to unique OTU ID of Data_Read object.
  • Tree: Tree object whose tip labels contains the unique OTU IDs of Data_Read object. Example file: http://viva1109.duckdns.org/paper/mTMAT/temporal_spatial_tree.rds
  • Taxonomy_Table: data.frame object that contains taxonomy informaion. First column has to contain the unique OTU IDs of Data_Read. From second to seventh columns contain Kingdom, Phylum, Class, Order, Family, Genus and Species taxonomy corresponding to each OTUs. Example file: http://viva1109.duckdns.org/paper/mTMAT/temporal_spatial_taxonomy.csv
  • Taxonomy_level: character object that specify the taxonomy level on which the analysis will be conducted. Can be either "Genus", "Family", "Order" or "Class". Upper than Order level is not recommended. (Dafault: "Genus") Example file: http://viva1109.duckdns.org/paper/mTMAT/temporal_spatial_taxonomy.csv
  • ncore: Number of core to use. (Dafault: 1). Option with >1 value is not available for windows OS system.
  • str: Working correaltion structure for sandwitch estimator. Can be either "independence", "exchangeable", "ar1", or "unstructured". (Default: "independence")
  • CompositionalBias_correction: Transformation applied to Gij for additional compositional bias adjustment. Can be either "Default", "ALR" or "CLR". (Default: "Default")
  • mTMAT_type: Type of TMAT statistic. Can be either "M" for TMATM or "IM" for TMATIM.
  • taxon_filter_zero_prop: Criteria that is applied to OTU filtering. OTUs which have more excessive proportion of zero count in samples than this criteria will be discared from analysis. (Default: 0.25)
  • Total_Readcount: If NULL, the program automatically calcualte it. You might want to specify this if you have OTU-filtered dataset. The values should match to the samples corresponding to the columns of Data_Read.

Output

mTMAT returns an object of class list containing the following compoenents.

  • resultTable: data.frame with the following columns.
  • $Upper Taxonoy : The one-rank higher taxonomy of Taxonomy_level.
  • $Target Taxonomy: The taxonomy specified in input object Taxonomy_level.
  • $mTMAT_beta: mTMAT beta.
  • $mTMAT_pvalue: mTMAT p-values.

  • $mTMAT_FDR: False discovery rate for each taxa.

  • output_raw: list object containing mTMAT results for plotting.

Plot mTMAT result

plot_mTMAT(target_taxon=NULL,taxon_index=NULL,output_mTMAT,phenotype_level=c("0","1"),x_label="Time",Cex_node=7,Cex_nodetxt=1.09, ...)

Input

  • target_taxon: character object specifying name of the taxon to plot.
  • taxon_index: numeric objec specifying the row number of the taxon to plot in resultTable object from mTMAT function.
  • output_mTMAT: Output object from mTMAT function.
  • phenotype_level: character vector that represents the label of phenotype. (Default: c("0","1"))
  • x_label: character object that represents the x-axis of distribution plot. This one corresponds to Time object in the input of mTMAT function.
  • Cex_node: numeric object that represents the size of nodes.
  • Cex_nodetxt: numeric object that represents the size of text inside of nodes.

Output

  • Distribution plot for the chosen taxon

Example

source("http://viva1109.duckdns.org/RFunctions/mTMAT.R",encoding = "UTF-8")

Data load

path_data<-"http://viva1109.duckdns.org/paper/mTMAT"
phenofile<-file.path(path_data,"temporal_spatial_meta.csv")
otufile<-file.path(path_data,"temporal_spatial_otuTable.csv")
taxonomyfile<-file.path(path_data,"temporal_spatial_taxonomy.csv")
treefile<-file.path(path_data,"temporal_spatial_tree.rds")

#Check the file paths
print(c(phenofile,otufile,taxonomyfile,treefile))
## [1] "http://viva1109.duckdns.org/paper/mTMAT/temporal_spatial_meta.csv"    
## [2] "http://viva1109.duckdns.org/paper/mTMAT/temporal_spatial_otuTable.csv"
## [3] "http://viva1109.duckdns.org/paper/mTMAT/temporal_spatial_taxonomy.csv"
## [4] "http://viva1109.duckdns.org/paper/mTMAT/temporal_spatial_tree.rds"
Metagenome file
#Metagenome file.
data_read<-read.csv(otufile,stringsAsFactors=F,check.names = F)
data_input<-data_read[,-1]
rownames(data_input)<-data_read[,1]
data_input[1:3,1:6]
##         1000301298 1000301308 1000301318 1000301328 1000301338 1000301488
## 4330849          0          0          0          0          0          0
## 4400869          0          0          0          0          0          0
## 4457086          0          0          0          0          0          0
Taxonomy information file
#This taxonomy information corresponds to the OTUs of data_input.
taxonomy_Table<-read.csv(taxonomyfile,stringsAsFactors=F)
head(taxonomy_Table,3)
##   X.OTU.ID  Kingdom Phylum   Class   Order  Family   Genus Species
## 1  4330849 Bacteria  P:SR1   P:SR1   P:SR1   P:SR1   P:SR1   P:SR1
## 2  4400869 Bacteria  P:SR1   P:SR1   P:SR1   P:SR1   P:SR1   P:SR1
## 3  4457086 Bacteria P:GN02 C:BD1-5 C:BD1-5 C:BD1-5 C:BD1-5 C:BD1-5
Phenotype and covariates file
#Phenotype file. Note that SampleID column should be a subset of the second to column names of data_input.
metadata<-read.csv(phenofile,stringsAsFactors=F,check.names = F)
head(metadata,3)
##     SampleID SubjectID Preg log_GDColl Time Birthweight_kg Length_at_birth
## 1 1000301298     10003 TRUE   5.288267    1           3.26              48
## 2 1000301308     10003 TRUE   5.323010    2           3.26              48
## 3 1000301318     10003 TRUE   5.356586    3           3.26              48
##   Number_of_babies            Indication         Hypertensive.Disorder
## 1                1 Hypertensive_disorder Mild_gestational_hypertension
## 2                1 Hypertensive_disorder Mild_gestational_hypertension
## 3                1 Hypertensive_disorder Mild_gestational_hypertension
##   SES_household_income
## 1          30000_50000
## 2          30000_50000
## 3          30000_50000
Phylogenetic tree
#Tree file.
download.file(treefile, "./temporal_spatial_tree.rds")
tree_raw<-readRDS("./temporal_spatial_tree.rds")
tree_raw
## 
## Phylogenetic tree with 1146 tips and 1144 internal nodes.
## 
## Tip labels:
##   4330849, 4400869, 4457086, 4457085, 655496, 4440970, ...
## 
## Unrooted; includes branch lengths.

Genus level analysis

Run mTMAT

output_mTMAT<-mTMAT(Data_Read=data_input,Meta_Table=metadata,Phenotype="Preg",ID="SubjectID",Time="log_GDColl",Cov=c("log_GDColl","Birthweight_kg"),Tree=tree_raw,Taxonomy_Table=taxonomy_Table)

# Significant result
ind_sig<-which(output_mTMAT$resultTable$mTMAT_FDR<0.05)
print(output_mTMAT$resultTable[ind_sig,])
##                  Family           Genus mTMAT_beta mTMAT_pvalue  mTMAT_FDR
## 1    F:Lactobacillaceae   Lactobacillus  1.2731392  0.001085959 0.02171917
## 2   F:[Tissierellaceae]    Anaerococcus -0.9723689  0.011065940 0.02213188
## 3      F:Prevotellaceae      Prevotella -0.6481359  0.024802292 0.04509508
## 4   F:[Tissierellaceae]   Peptoniphilus -0.8691580  0.007801490 0.02213188
## 5     F:Veillonellaceae       Dialister -0.5644144  0.007460551 0.02213188
## 9  F:Campylobacteraceae   Campylobacter -0.6524956  0.007299725 0.02213188
## 11  F:[Tissierellaceae]      Finegoldia -0.7217598  0.004227707 0.02213188
## 14      O:Clostridiales O:Clostridiales -0.7279819  0.004936155 0.02213188
## 15 F:Porphyromonadaceae   Porphyromonas -0.5812502  0.006163016 0.02213188
## 17   F:Streptococcaceae   Streptococcus -0.7355147  0.010312580 0.02213188
## 19   F:Actinomycetaceae     Varibaculum -0.4764613  0.009127836 0.02213188

Generate distribution plots

y_level<-c("Unpregnant","Pregnant")
x_axis_label<-"log of days from first visit"
plot_mTMAT(taxon_index=1,output_mTMAT=output_mTMAT,phenotype_level=y_level,x_label=x_axis_label)

plot_mTMAT(target_taxon="Anaerococcus",output_mTMAT=output_mTMAT)

plot_mTMAT(target_taxon="Streptococcus",output_mTMAT=output_mTMAT,phenotype_level=y_level,x_label=x_axis_label)

plot_mTMAT(taxon_index=19,output_mTMAT=output_mTMAT,phenotype_level=y_level,x_label=x_axis_label)

Family level analysis

Run mTMAT

output_mTMAT_family<-mTMAT(Data_Read=data_input,Meta_Table=metadata,Phenotype="Preg",ID="SubjectID",Time="log_GDColl",Cov=c("log_GDColl","Birthweight_kg"),Tree=tree_raw,Taxonomy_Table=taxonomy_Table,Taxonomy_level="Family")
# 
ind_sig<-which(output_mTMAT_family$resultTable$mTMAT_FDR<0.05)
# 
print(output_mTMAT_family$resultTable[ind_sig,])
##                  Order               Family mTMAT_beta mTMAT_pvalue  mTMAT_FDR
## 1      O:Clostridiales  F:[Tissierellaceae] -0.9515440  0.023160039 0.04340401
## 2    O:Lactobacillales   F:Lactobacillaceae  1.2731392  0.001085959 0.01520342
## 3      O:Bacteroidales     F:Prevotellaceae -0.6481359  0.024802292 0.04340401
## 5      O:Clostridiales    F:Veillonellaceae -0.5644144  0.007460551 0.02088954
## 7  O:Campylobacterales F:Campylobacteraceae -0.6524956  0.007299725 0.02088954
## 11     O:Bacteroidales F:Porphyromonadaceae -0.5812502  0.006163016 0.02088954
## 13   O:Lactobacillales   F:Streptococcaceae -0.7355147  0.010312580 0.02406269
## 14     O:Clostridiales      O:Clostridiales -0.7279819  0.004936155 0.02088954

Generate distribution plots

plot_mTMAT(taxon_index=ind_sig[1],output_mTMAT=output_mTMAT_family,phenotype_level=y_level,x_label=x_axis_label)

plot_mTMAT(taxon_index=ind_sig[2],output_mTMAT=output_mTMAT_family,phenotype_level=y_level,x_label=x_axis_label)