lefse_analysis_cleaned_prevalence
Loading packages and the second phyloseq object
library(tidyverse)
library(phyloseq)
library(microbiomeMarker)
phy<-readRDS("D:/R/Masha/phy_cleaned_prevalence.rds")
phy## phyloseq-class experiment-level object
## otu_table() OTU Table: [ 1571 taxa and 64 samples ]
## sample_data() Sample Data: [ 64 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 1571 taxa by 6 taxonomic ranks ]
## phy_tree() Phylogenetic Tree: [ 1571 tips and 1569 internal nodes ]This phyloseq object was cleaned by:
Removing NA phylum sequences.
Removing sequences that occur in not more than 3 samples (5% of our samples).
Removing the phyla Cyanobacteria/Chloroplast, Deferribacteres, and Tenericutes which are comprised of low-prevalence features.
This phyloseq has 1571 sequences and 64 samples.
Subset samples with high fiber diet and vitamin A
phy_high<-subset_samples(phy, group %in% c("AH","BH"))
phy_high## phyloseq-class experiment-level object
## otu_table() OTU Table: [ 1571 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 1571 taxa by 6 taxonomic ranks ]
## phy_tree() Phylogenetic Tree: [ 1571 tips and 1569 internal nodes ]table(sample_data(phy_high)$group)##
## AH BH
## 8 8This phyloseq has 1571 sequences and 16 samples.
8 samples with Retinyl Acetate High Fiber Diet.
8 samples with Beta Carotene High Fiber Diet.
Linear discriminant analysis effect size (LEFSe)
For all ranks, “Kingdom” “Phylum” “Class” “Order” “Family” “Genus”
lef_out<-run_lefse(phy_high, group = "group", norm = "CPM",
kw_cutoff = 0.05, lda_cutoff = 2)
lef_out## microbiomeMarker-class inherited from phyloseq-class
## normalization method: [ CPM ]
## microbiome marker identity method: [ lefse ]
## marker_table() Marker Table: [ 34 microbiome markers with 5 variables ]
## otu_table() OTU Table: [ 129 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 129 taxa by 1 taxonomic ranks ]There are 34 microbiome markers.
Microbiome marker table
library(knitr)
dat<-marker_table(lef_out) %>% data.frame() %>% select(1:4)
head(dat)## feature
## marker1 k__Bacteria|p__Firmicutes|c__Erysipelotrichia
## marker2 k__Bacteria|p__Firmicutes|c__Erysipelotrichia|o__Erysipelotrichales
## marker3 k__Bacteria|p__Firmicutes|c__Erysipelotrichia|o__Erysipelotrichales|f__Erysipelotrichaceae
## marker4 k__Bacteria|p__Firmicutes|c__Erysipelotrichia|o__Erysipelotrichales|f__Erysipelotrichaceae|g__Dubosiella
## marker5 k__Bacteria|p__Bacteroidetes|c__Bacteroidia|o__Bacteroidales|f__Muribaculaceae|g__Muribaculaceae_g__
## marker6 k__Bacteria|p__Proteobacteria
## enrich_group ef_lda pvalue
## marker1 AH 5.242020 0.0023220945
## marker2 AH 5.242020 0.0023220945
## marker3 AH 5.242020 0.0023220945
## marker4 AH 5.042356 0.0011313264
## marker5 AH 4.253415 0.0007100757
## marker6 AH 3.490142 0.0032758975dat %>% filter(enrich_group=="BH") %>% head()## feature
## marker16 k__Bacteria|p__Firmicutes|c__Bacilli
## marker17 k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales
## marker18 k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales|f__Lactobacillaceae
## marker19 k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales|f__Lactobacillaceae|g__Lactobacillus
## marker20 k__Bacteria|p__Firmicutes|c__Clostridia
## marker21 k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales
## enrich_group ef_lda pvalue
## marker16 BH 5.027163 0.0007775304
## marker17 BH 5.026887 0.0007775304
## marker18 BH 5.025903 0.0007775304
## marker19 BH 5.004548 0.0007775304
## marker20 BH 4.972364 0.0007775304
## marker21 BH 4.972364 0.0007775304dat %>% kable(align = "c")| feature | enrich_group | ef_lda | pvalue | |
|---|---|---|---|---|
| marker1 | k__Bacteria|p__Firmicutes|c__Erysipelotrichia | AH | 5.242020 | 0.0023221 |
| marker2 | k__Bacteria|p__Firmicutes|c__Erysipelotrichia|o__Erysipelotrichales | AH | 5.242020 | 0.0023221 |
| marker3 | k__Bacteria|p__Firmicutes|c__Erysipelotrichia|o__Erysipelotrichales|f__Erysipelotrichaceae | AH | 5.242020 | 0.0023221 |
| marker4 | k__Bacteria|p__Firmicutes|c__Erysipelotrichia|o__Erysipelotrichales|f__Erysipelotrichaceae|g__Dubosiella | AH | 5.042356 | 0.0011313 |
| marker5 | k__Bacteria|p__Bacteroidetes|c__Bacteroidia|o__Bacteroidales|f__Muribaculaceae|gMuribaculaceae_g | AH | 4.253415 | 0.0007101 |
| marker6 | k__Bacteria|p__Proteobacteria | AH | 3.490142 | 0.0032759 |
| marker7 | k__Bacteria|p__Proteobacteria|c__Betaproteobacteria | AH | 3.489006 | 0.0032759 |
| marker8 | k__Bacteria|p__Proteobacteria|c__Betaproteobacteria|o__Burkholderiales | AH | 3.489006 | 0.0032759 |
| marker9 | k__Bacteria|p__Proteobacteria|c__Betaproteobacteria|o__Burkholderiales|f__Sutterellaceae | AH | 3.489006 | 0.0032759 |
| marker10 | k__Bacteria|p__Proteobacteria|c__Betaproteobacteria|o__Burkholderiales|f__Sutterellaceae|g__Parasutterella | AH | 3.489006 | 0.0032759 |
| marker11 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Ruminococcaceae|g__Sporobacter | AH | 2.366038 | 0.0355562 |
| marker12 | k__Bacteria|p__Actinobacteria|c__Coriobacteriia|o__Eggerthellales | AH | 2.162231 | 0.0086026 |
| marker13 | k__Bacteria|p__Actinobacteria|c__Coriobacteriia|o__Eggerthellales|f__Eggerthellaceae | AH | 2.162231 | 0.0086026 |
| marker14 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Ruminococcaceae|g__Phocea | AH | 2.098340 | 0.0105152 |
| marker15 | k__Bacteria|p__Actinobacteria|c__Coriobacteriia|o__Eggerthellales|f__Eggerthellaceae|g__Adlercreutzia | AH | 2.018143 | 0.0458395 |
| marker16 | k__Bacteria|p__Firmicutes|c__Bacilli | BH | 5.027163 | 0.0007775 |
| marker17 | k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales | BH | 5.026887 | 0.0007775 |
| marker18 | k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales|f__Lactobacillaceae | BH | 5.025904 | 0.0007775 |
| marker19 | k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales|f__Lactobacillaceae|g__Lactobacillus | BH | 5.004548 | 0.0007775 |
| marker20 | k__Bacteria|p__Firmicutes|c__Clostridia | BH | 4.972364 | 0.0007775 |
| marker21 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales | BH | 4.972364 | 0.0007775 |
| marker22 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Lachnospiraceae | BH | 4.810717 | 0.0011313 |
| marker23 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Clostridiales_f__ | BH | 4.394708 | 0.0063229 |
| marker24 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Clostridiales_f__|g__Clostridiales_f___g__ | BH | 4.394708 | 0.0063229 |
| marker25 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Lachnospiraceae|g__Mediterraneibacter | BH | 4.264031 | 0.0045744 |
| marker26 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Lachnospiraceae|g__Hungatella | BH | 3.902328 | 0.0274232 |
| marker27 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Lachnospiraceae|g__Marvinbryantia | BH | 3.729293 | 0.0023221 |
| marker28 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Ruminococcaceae|gRuminococcaceae_g | BH | 3.675561 | 0.0063229 |
| marker29 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Lachnospiraceae|g__Murimonas | BH | 3.253575 | 0.0274232 |
| marker30 | k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales|f__Lactobacillaceae|gLactobacillaceae_g | BH | 3.180085 | 0.0007775 |
| marker31 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Lachnospiraceae|g__Dorea | BH | 2.621998 | 0.0012165 |
| marker32 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Clostridiales_Incertae_Sedis_XIII | BH | 2.423565 | 0.0157143 |
| marker33 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Clostridiales_Incertae_Sedis_XIII|g__Ihubacter | BH | 2.420399 | 0.0208626 |
| marker34 | k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales|f__Lactobacillaceae|g__Bombilactobacillus | BH | 2.141503 | 0.0030757 |
The marker table has 34 rows (for 34 markers)
There are 15 markers for AH group and 19 markers for BH group.
The most features that explain AH group are:
class Erysipelotrichia
order Erysipelotrichales
family Erysipelotrichaceae
genus Dubosiella from the above family, order, and class
The most features that explain BH group are:
class Bacilli
order Lactobacillales
family Lactobacillaceae
genus Lactobacillus from the above family, order, and class
Visualizing marker table as histogram
plot_ef_bar(lef_out)
We have 34 bars as we have 34 markers, they are colored by group.
Visualizing the marker table as cladogram
clade_label_level argument controls the max level of taxa used to label the clade, other level of taxa will be shown on the side.
The default is 4 which is showing the kingdom, Phylum, Class on the cladogram.
# kingdom, Phylum, Class
plot_cladogram(lef_out, color = c("red","blue"), clade_label_level = 4)
# all except genus
plot_cladogram(lef_out, color = c("red","blue"), clade_label_level = 2)
# phylum and kingdom
plot_cladogram(lef_out, color = c("red","blue"), clade_label_level = 5)
Specifying the lefse analysis to the desired taxonomy rank
For class
lef_out<-run_lefse(phy_high, group = "group", norm = "CPM",
taxa_rank = "Class",
kw_cutoff = 0.05, lda_cutoff = 2)
lef_out## microbiomeMarker-class inherited from phyloseq-class
## normalization method: [ CPM ]
## microbiome marker identity method: [ lefse ]
## marker_table() Marker Table: [ 4 microbiome markers with 5 variables ]
## otu_table() OTU Table: [ 10 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 10 taxa by 1 taxonomic ranks ]dat<-marker_table(lef_out) %>% data.frame() %>% select(1:4)
head(dat)## feature enrich_group ef_lda pvalue
## marker1 Erysipelotrichia AH 5.262635 0.0023220945
## marker2 Betaproteobacteria AH 3.867739 0.0032758975
## marker3 Bacilli BH 5.054844 0.0007775304
## marker4 Clostridia BH 4.957244 0.0007775304dat %>% filter(enrich_group=="BH") %>% head()## feature enrich_group ef_lda pvalue
## marker3 Bacilli BH 5.054844 0.0007775304
## marker4 Clostridia BH 4.957244 0.0007775304dat %>% kable(align = "c")| feature | enrich_group | ef_lda | pvalue | |
|---|---|---|---|---|
| marker1 | Erysipelotrichia | AH | 5.262635 | 0.0023221 |
| marker2 | Betaproteobacteria | AH | 3.867739 | 0.0032759 |
| marker3 | Bacilli | BH | 5.054844 | 0.0007775 |
| marker4 | Clostridia | BH | 4.957244 | 0.0007775 |
plot_ef_bar(lef_out)
Only 4 markers (classes) are present, 2 for the AH group and 2 for the BH group.
For order
lef_out<-run_lefse(phy_high, group = "group", norm = "CPM",
taxa_rank = "Order",
kw_cutoff = 0.05, lda_cutoff = 2)
lef_out## microbiomeMarker-class inherited from phyloseq-class
## normalization method: [ CPM ]
## microbiome marker identity method: [ lefse ]
## marker_table() Marker Table: [ 5 microbiome markers with 5 variables ]
## otu_table() OTU Table: [ 14 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 14 taxa by 1 taxonomic ranks ]dat<-marker_table(lef_out) %>% data.frame() %>% select(1:4)
head(dat)## feature enrich_group ef_lda pvalue
## marker1 Erysipelotrichales AH 5.278958 0.0023220945
## marker2 Burkholderiales AH 3.417507 0.0032758975
## marker3 Eggerthellales AH 2.264019 0.0086025506
## marker4 Lactobacillales BH 5.052625 0.0007775304
## marker5 Clostridiales BH 4.951872 0.0007775304dat %>% filter(enrich_group=="BH") %>% head()## feature enrich_group ef_lda pvalue
## marker4 Lactobacillales BH 5.052625 0.0007775304
## marker5 Clostridiales BH 4.951872 0.0007775304dat %>% kable(align = "c")| feature | enrich_group | ef_lda | pvalue | |
|---|---|---|---|---|
| marker1 | Erysipelotrichales | AH | 5.278957 | 0.0023221 |
| marker2 | Burkholderiales | AH | 3.417507 | 0.0032759 |
| marker3 | Eggerthellales | AH | 2.264019 | 0.0086026 |
| marker4 | Lactobacillales | BH | 5.052625 | 0.0007775 |
| marker5 | Clostridiales | BH | 4.951872 | 0.0007775 |
plot_ef_bar(lef_out)
Only 5 markers (orders) are present, 3 for the AH group and 2 for the BH group.
For family
lef_out<-run_lefse(phy_high, group = "group", norm = "CPM",
taxa_rank = "Family",
kw_cutoff = 0.05, lda_cutoff = 2)
lef_out## microbiomeMarker-class inherited from phyloseq-class
## normalization method: [ CPM ]
## microbiome marker identity method: [ lefse ]
## marker_table() Marker Table: [ 8 microbiome markers with 5 variables ]
## otu_table() OTU Table: [ 24 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 24 taxa by 1 taxonomic ranks ]dat<-marker_table(lef_out) %>% data.frame() %>% select(1:4)
head(dat)## feature enrich_group ef_lda pvalue
## marker1 Erysipelotrichaceae AH 5.278959 0.0023220945
## marker2 Sutterellaceae AH 3.420520 0.0032758975
## marker3 Eggerthellaceae AH 2.442457 0.0086025506
## marker4 Lactobacillaceae BH 5.051527 0.0007775304
## marker5 Lachnospiraceae BH 4.782825 0.0011313264
## marker6 Clostridiales_f__ BH 4.382839 0.0063229477dat %>% filter(enrich_group=="BH") %>% head()## feature enrich_group ef_lda pvalue
## marker4 Lactobacillaceae BH 5.051527 0.0007775304
## marker5 Lachnospiraceae BH 4.782825 0.0011313264
## marker6 Clostridiales_f__ BH 4.382839 0.0063229477
## marker7 Clostridiales_Incertae_Sedis_XIII BH 2.538418 0.0157143498
## marker8 Coriobacteriia_o___f__ BH 2.070047 0.0346984836dat %>% kable(align = "c")| feature | enrich_group | ef_lda | pvalue | |
|---|---|---|---|---|
| marker1 | Erysipelotrichaceae | AH | 5.278959 | 0.0023221 |
| marker2 | Sutterellaceae | AH | 3.420520 | 0.0032759 |
| marker3 | Eggerthellaceae | AH | 2.442457 | 0.0086026 |
| marker4 | Lactobacillaceae | BH | 5.051527 | 0.0007775 |
| marker5 | Lachnospiraceae | BH | 4.782825 | 0.0011313 |
| marker6 | Clostridiales_f__ | BH | 4.382839 | 0.0063229 |
| marker7 | Clostridiales_Incertae_Sedis_XIII | BH | 2.538418 | 0.0157143 |
| marker8 | Coriobacteriia_o_f | BH | 2.070047 | 0.0346985 |
plot_ef_bar(lef_out)
Only 8 markers (families) are present, 3 for the AH group and 5 for the BH group.
For genus
lef_out<-run_lefse(phy_high, group = "group", norm = "CPM",
taxa_rank = "Genus",
kw_cutoff = 0.05, lda_cutoff = 2)
lef_out## microbiomeMarker-class inherited from phyloseq-class
## normalization method: [ CPM ]
## microbiome marker identity method: [ lefse ]
## marker_table() Marker Table: [ 17 microbiome markers with 5 variables ]
## otu_table() OTU Table: [ 75 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 75 taxa by 1 taxonomic ranks ]dat<-marker_table(lef_out) %>% data.frame() %>% select(1:4)
head(dat)## feature enrich_group ef_lda pvalue
## marker1 Dubosiella AH 5.065150 0.0011313264
## marker2 Muribaculaceae_g__ AH 4.229787 0.0007100757
## marker3 Parasutterella AH 3.425972 0.0032758975
## marker4 Sporobacter AH 2.372997 0.0355562316
## marker5 Phocea AH 2.113502 0.0105152459
## marker6 Adlercreutzia AH 2.009097 0.0458395302dat %>% filter(enrich_group=="BH") %>% head()## feature enrich_group ef_lda pvalue
## marker7 Lactobacillus BH 5.017786 0.0007775304
## marker8 Clostridiales_f___g__ BH 4.437438 0.0063229477
## marker9 Mediterraneibacter BH 4.295739 0.0045744395
## marker10 Hungatella BH 3.997191 0.0274231544
## marker11 Ruminococcaceae_g__ BH 3.733366 0.0063229477
## marker12 Marvinbryantia BH 3.693001 0.0023220945dat %>% kable(align = "c")| feature | enrich_group | ef_lda | pvalue | |
|---|---|---|---|---|
| marker1 | Dubosiella | AH | 5.065150 | 0.0011313 |
| marker2 | Muribaculaceae_g__ | AH | 4.229787 | 0.0007101 |
| marker3 | Parasutterella | AH | 3.425972 | 0.0032759 |
| marker4 | Sporobacter | AH | 2.372997 | 0.0355562 |
| marker5 | Phocea | AH | 2.113502 | 0.0105152 |
| marker6 | Adlercreutzia | AH | 2.009097 | 0.0458395 |
| marker7 | Lactobacillus | BH | 5.017786 | 0.0007775 |
| marker8 | Clostridiales_f_g | BH | 4.437438 | 0.0063229 |
| marker9 | Mediterraneibacter | BH | 4.295739 | 0.0045744 |
| marker10 | Hungatella | BH | 3.997191 | 0.0274232 |
| marker11 | Ruminococcaceae_g__ | BH | 3.733366 | 0.0063229 |
| marker12 | Marvinbryantia | BH | 3.693001 | 0.0023221 |
| marker13 | Murimonas | BH | 3.268205 | 0.0274232 |
| marker14 | Lactobacillaceae_g__ | BH | 3.151868 | 0.0007775 |
| marker15 | Dorea | BH | 2.621053 | 0.0012165 |
| marker16 | Ihubacter | BH | 2.471276 | 0.0208626 |
| marker17 | Bombilactobacillus | BH | 2.062586 | 0.0030757 |
plot_ef_bar(lef_out)
Only 17 markers (genus) are present, 6 for the AH group and 11 for the BH group.
Focusing on the most prominent features
We decrease the error rate to 0.01 and increase the lda score to 4.
This approach is suggested by the author of the microbiomeMarker package.
For all ranks, “Kingdom” “Phylum” “Class” “Order” “Family” “Genus”
lef_out<-run_lefse(phy_high, group = "group", norm = "CPM",
kw_cutoff = 0.01, lda_cutoff = 4)
lef_out## microbiomeMarker-class inherited from phyloseq-class
## normalization method: [ CPM ]
## microbiome marker identity method: [ lefse ]
## marker_table() Marker Table: [ 15 microbiome markers with 5 variables ]
## otu_table() OTU Table: [ 129 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 129 taxa by 1 taxonomic ranks ]There are 15 microbiome markers.
Microbiome marker table
library(knitr)
dat<-marker_table(lef_out) %>% data.frame() %>% select(1:4)
head(dat)## feature
## marker1 k__Bacteria|p__Firmicutes|c__Erysipelotrichia
## marker2 k__Bacteria|p__Firmicutes|c__Erysipelotrichia|o__Erysipelotrichales
## marker3 k__Bacteria|p__Firmicutes|c__Erysipelotrichia|o__Erysipelotrichales|f__Erysipelotrichaceae
## marker4 k__Bacteria|p__Firmicutes|c__Erysipelotrichia|o__Erysipelotrichales|f__Erysipelotrichaceae|g__Dubosiella
## marker5 k__Bacteria|p__Bacteroidetes|c__Bacteroidia|o__Bacteroidales|f__Muribaculaceae|g__Muribaculaceae_g__
## marker6 k__Bacteria|p__Firmicutes|c__Bacilli
## enrich_group ef_lda pvalue
## marker1 AH 5.320490 0.0023220945
## marker2 AH 5.320490 0.0023220945
## marker3 AH 5.320490 0.0023220945
## marker4 AH 5.075103 0.0011313264
## marker5 AH 4.198036 0.0007100757
## marker6 BH 5.080436 0.0007775304dat %>% filter(enrich_group=="BH") %>% head()## feature
## marker6 k__Bacteria|p__Firmicutes|c__Bacilli
## marker7 k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales
## marker8 k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales|f__Lactobacillaceae
## marker9 k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales|f__Lactobacillaceae|g__Lactobacillus
## marker10 k__Bacteria|p__Firmicutes|c__Clostridia
## marker11 k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales
## enrich_group ef_lda pvalue
## marker6 BH 5.080436 0.0007775304
## marker7 BH 5.080052 0.0007775304
## marker8 BH 5.078891 0.0007775304
## marker9 BH 5.053127 0.0007775304
## marker10 BH 4.966525 0.0007775304
## marker11 BH 4.966525 0.0007775304dat %>% kable(align = "c")| feature | enrich_group | ef_lda | pvalue | |
|---|---|---|---|---|
| marker1 | k__Bacteria|p__Firmicutes|c__Erysipelotrichia | AH | 5.320490 | 0.0023221 |
| marker2 | k__Bacteria|p__Firmicutes|c__Erysipelotrichia|o__Erysipelotrichales | AH | 5.320490 | 0.0023221 |
| marker3 | k__Bacteria|p__Firmicutes|c__Erysipelotrichia|o__Erysipelotrichales|f__Erysipelotrichaceae | AH | 5.320490 | 0.0023221 |
| marker4 | k__Bacteria|p__Firmicutes|c__Erysipelotrichia|o__Erysipelotrichales|f__Erysipelotrichaceae|g__Dubosiella | AH | 5.075103 | 0.0011313 |
| marker5 | k__Bacteria|p__Bacteroidetes|c__Bacteroidia|o__Bacteroidales|f__Muribaculaceae|gMuribaculaceae_g | AH | 4.198035 | 0.0007101 |
| marker6 | k__Bacteria|p__Firmicutes|c__Bacilli | BH | 5.080436 | 0.0007775 |
| marker7 | k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales | BH | 5.080052 | 0.0007775 |
| marker8 | k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales|f__Lactobacillaceae | BH | 5.078891 | 0.0007775 |
| marker9 | k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales|f__Lactobacillaceae|g__Lactobacillus | BH | 5.053127 | 0.0007775 |
| marker10 | k__Bacteria|p__Firmicutes|c__Clostridia | BH | 4.966525 | 0.0007775 |
| marker11 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales | BH | 4.966525 | 0.0007775 |
| marker12 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Lachnospiraceae | BH | 4.822182 | 0.0011313 |
| marker13 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Lachnospiraceae|g__Mediterraneibacter | BH | 4.407640 | 0.0045744 |
| marker14 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Clostridiales_f__ | BH | 4.363845 | 0.0063229 |
| marker15 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Clostridiales_f__|g__Clostridiales_f___g__ | BH | 4.363845 | 0.0063229 |
The marker table has 15 rows (for 15 markers)
There are 5 markers for AH group and 10 markers for BH group.
The most features that explain AH group are:
class Erysipelotrichia
order Erysipelotrichales
family Erysipelotrichaceae
genus Dubosiella from the above family, order, and class
The most features that explain BH group are:
class Bacilli
order Lactobacillales
family Lactobacillaceae
genus Lactobacillus from the above family, order, and class
Visualizing marker table as histogram
plot_ef_bar(lef_out)
We have 15 bars as we have 15 markers.
Visualizing the marker table as cladogram
clade_label_level argument controls the max level of taxa used to label the clade, other level of taxa will be shown on the side.
The default is 4 which is showing the kingdom, Phylum, Class on the cladogram.
# kingdom, Phylum, Class
plot_cladogram(lef_out, color = c("red","blue"), clade_label_level = 4)
# all except genus
plot_cladogram(lef_out, color = c("red","blue"), clade_label_level = 2)
# phylum and kingdom
plot_cladogram(lef_out, color = c("red","blue"), clade_label_level = 5)
The white circles are non significant taxa.
Specifying the lefse analysis to the desired taxonomy rank
For class
lef_out<-run_lefse(phy_high, group = "group", norm = "CPM",
taxa_rank = "Class",
kw_cutoff = 0.01, lda_cutoff = 4)
lef_out## microbiomeMarker-class inherited from phyloseq-class
## normalization method: [ CPM ]
## microbiome marker identity method: [ lefse ]
## marker_table() Marker Table: [ 3 microbiome markers with 5 variables ]
## otu_table() OTU Table: [ 10 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 10 taxa by 1 taxonomic ranks ]dat<-marker_table(lef_out) %>% data.frame() %>% select(1:4)
head(dat)## feature enrich_group ef_lda pvalue
## marker1 Erysipelotrichia AH 5.262635 0.0023220945
## marker2 Bacilli BH 5.054844 0.0007775304
## marker3 Clostridia BH 4.957244 0.0007775304dat %>% filter(enrich_group=="BH") %>% head()## feature enrich_group ef_lda pvalue
## marker2 Bacilli BH 5.054844 0.0007775304
## marker3 Clostridia BH 4.957244 0.0007775304dat %>% kable(align = "c")| feature | enrich_group | ef_lda | pvalue | |
|---|---|---|---|---|
| marker1 | Erysipelotrichia | AH | 5.262635 | 0.0023221 |
| marker2 | Bacilli | BH | 5.054844 | 0.0007775 |
| marker3 | Clostridia | BH | 4.957244 | 0.0007775 |
plot_ef_bar(lef_out)
Only 3 markers (classes) are present, 1 for the AH group and 2 for the BH group.
For order
lef_out<-run_lefse(phy_high, group = "group", norm = "CPM",
taxa_rank = "Order",
kw_cutoff = 0.01, lda_cutoff = 4)
lef_out## microbiomeMarker-class inherited from phyloseq-class
## normalization method: [ CPM ]
## microbiome marker identity method: [ lefse ]
## marker_table() Marker Table: [ 3 microbiome markers with 5 variables ]
## otu_table() OTU Table: [ 14 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 14 taxa by 1 taxonomic ranks ]dat<-marker_table(lef_out) %>% data.frame() %>% select(1:4)
head(dat)## feature enrich_group ef_lda pvalue
## marker1 Erysipelotrichales AH 5.262023 0.0023220945
## marker2 Lactobacillales BH 5.050205 0.0007775304
## marker3 Clostridiales BH 4.951128 0.0007775304dat %>% filter(enrich_group=="BH") %>% head()## feature enrich_group ef_lda pvalue
## marker2 Lactobacillales BH 5.050205 0.0007775304
## marker3 Clostridiales BH 4.951128 0.0007775304dat %>% kable(align = "c")| feature | enrich_group | ef_lda | pvalue | |
|---|---|---|---|---|
| marker1 | Erysipelotrichales | AH | 5.262023 | 0.0023221 |
| marker2 | Lactobacillales | BH | 5.050205 | 0.0007775 |
| marker3 | Clostridiales | BH | 4.951128 | 0.0007775 |
plot_ef_bar(lef_out)
Only 3 markers (orders) are present, 1 for the AH group and 2 for the BH group.
For family
lef_out<-run_lefse(phy_high, group = "group", norm = "CPM",
taxa_rank = "Family",
kw_cutoff = 0.01, lda_cutoff = 4)
lef_out## microbiomeMarker-class inherited from phyloseq-class
## normalization method: [ CPM ]
## microbiome marker identity method: [ lefse ]
## marker_table() Marker Table: [ 4 microbiome markers with 5 variables ]
## otu_table() OTU Table: [ 24 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 24 taxa by 1 taxonomic ranks ]dat<-marker_table(lef_out) %>% data.frame() %>% select(1:4)
head(dat)## feature enrich_group ef_lda pvalue
## marker1 Erysipelotrichaceae AH 5.262038 0.0023220945
## marker2 Lactobacillaceae BH 5.049207 0.0007775304
## marker3 Lachnospiraceae BH 4.780782 0.0011313264
## marker4 Clostridiales_f__ BH 4.388033 0.0063229477dat %>% filter(enrich_group=="BH") %>% head()## feature enrich_group ef_lda pvalue
## marker2 Lactobacillaceae BH 5.049207 0.0007775304
## marker3 Lachnospiraceae BH 4.780782 0.0011313264
## marker4 Clostridiales_f__ BH 4.388033 0.0063229477dat %>% kable(align = "c")| feature | enrich_group | ef_lda | pvalue | |
|---|---|---|---|---|
| marker1 | Erysipelotrichaceae | AH | 5.262038 | 0.0023221 |
| marker2 | Lactobacillaceae | BH | 5.049207 | 0.0007775 |
| marker3 | Lachnospiraceae | BH | 4.780782 | 0.0011313 |
| marker4 | Clostridiales_f__ | BH | 4.388033 | 0.0063229 |
plot_ef_bar(lef_out)
Only 4 markers (families) are present, 1 for the AH group and 3 for the BH group.
For genus
lef_out<-run_lefse(phy_high, group = "group", norm = "CPM",
taxa_rank = "Genus",
kw_cutoff = 0.01, lda_cutoff = 4)
lef_out## microbiomeMarker-class inherited from phyloseq-class
## normalization method: [ CPM ]
## microbiome marker identity method: [ lefse ]
## marker_table() Marker Table: [ 5 microbiome markers with 5 variables ]
## otu_table() OTU Table: [ 75 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 75 taxa by 1 taxonomic ranks ]dat<-marker_table(lef_out) %>% data.frame() %>% select(1:4)
head(dat)## feature enrich_group ef_lda pvalue
## marker1 Dubosiella AH 5.058365 0.0011313264
## marker2 Muribaculaceae_g__ AH 4.274511 0.0007100757
## marker3 Lactobacillus BH 5.036635 0.0007775304
## marker4 Clostridiales_f___g__ BH 4.379020 0.0063229477
## marker5 Mediterraneibacter BH 4.313753 0.0045744395dat %>% filter(enrich_group=="BH") %>% head()## feature enrich_group ef_lda pvalue
## marker3 Lactobacillus BH 5.036635 0.0007775304
## marker4 Clostridiales_f___g__ BH 4.379020 0.0063229477
## marker5 Mediterraneibacter BH 4.313753 0.0045744395dat %>% kable(align = "c")| feature | enrich_group | ef_lda | pvalue | |
|---|---|---|---|---|
| marker1 | Dubosiella | AH | 5.058365 | 0.0011313 |
| marker2 | Muribaculaceae_g__ | AH | 4.274511 | 0.0007101 |
| marker3 | Lactobacillus | BH | 5.036635 | 0.0007775 |
| marker4 | Clostridiales_f_g | BH | 4.379020 | 0.0063229 |
| marker5 | Mediterraneibacter | BH | 4.313753 | 0.0045744 |
plot_ef_bar(lef_out)
Only 5 markers (genus) are present, 2 for the AH group and 3 for the BH group.
Subset samples with low fiber diet and vitamin A
phy_low<-subset_samples(phy, group %in% c("AL","BL"))
phy_low## phyloseq-class experiment-level object
## otu_table() OTU Table: [ 1571 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 1571 taxa by 6 taxonomic ranks ]
## phy_tree() Phylogenetic Tree: [ 1571 tips and 1569 internal nodes ]table(sample_data(phy_low)$group)##
## AL BL
## 8 8This phyloseq has 1571 sequences and 16 samples.
8 samples with Retinyl Acetate Low Fiber Diet.
8 samples with Beta Carotene Low Fiber Diet.
Linear discriminant analysis effect size (LEFSe)
For all ranks, “Kingdom” “Phylum” “Class” “Order” “Family” “Genus”
lef_out<-run_lefse(phy_low, group = "group", norm = "CPM",
kw_cutoff = 0.05, lda_cutoff = 2)
lef_out## microbiomeMarker-class inherited from phyloseq-class
## normalization method: [ CPM ]
## microbiome marker identity method: [ lefse ]
## marker_table() Marker Table: [ 35 microbiome markers with 5 variables ]
## otu_table() OTU Table: [ 123 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 123 taxa by 1 taxonomic ranks ]There are 35 microbiome markers.
Microbiome marker table
library(knitr)
dat<-marker_table(lef_out) %>% data.frame() %>% select(1:4)
head(dat)## feature
## marker1 k__Bacteria|p__Bacteroidetes|c__Bacteroidia|o__Bacteroidales|f__Muribaculaceae|g__Muribaculaceae_g__
## marker2 k__Bacteria|p__Bacteroidetes|c__Bacteroidia|o__Bacteroidales|f__Bacteroidaceae|g__Bacteroides
## marker3 k__Bacteria|p__Bacteroidetes|c__Bacteroidia|o__Bacteroidales|f__Bacteroidaceae
## marker4 k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Ruminococcaceae|g__Intestinimonas
## marker5 k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Ruminococcaceae|g__Phocea
## marker6 k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Ruminococcaceae|g__Anaerotruncus
## enrich_group ef_lda pvalue
## marker1 AL 4.044053 0.020862582
## marker2 AL 3.738189 0.027423154
## marker3 AL 3.729778 0.035691900
## marker4 AL 3.173991 0.035691900
## marker5 AL 2.507475 0.001765096
## marker6 AL 2.384923 0.024860789dat %>% filter(enrich_group=="BL") %>% head()## feature
## marker9 k__Bacteria|p__Firmicutes|c__Bacilli
## marker10 k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales
## marker11 k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales|f__Streptococcaceae
## marker12 k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales|f__Streptococcaceae|g__Lactococcus
## marker13 k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Ruminococcaceae|g__Ruminococcaceae_g__
## marker14 k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Lachnospiraceae|g__Schaedlerella
## enrich_group ef_lda pvalue
## marker9 BL 4.884616 0.027423154
## marker10 BL 4.874638 0.027423154
## marker11 BL 4.863216 0.020862582
## marker12 BL 4.863216 0.020862582
## marker13 BL 4.536076 0.045999367
## marker14 BL 4.256198 0.003275897dat %>% kable(align = "c")| feature | enrich_group | ef_lda | pvalue | |
|---|---|---|---|---|
| marker1 | k__Bacteria|p__Bacteroidetes|c__Bacteroidia|o__Bacteroidales|f__Muribaculaceae|gMuribaculaceae_g | AL | 4.044053 | 0.0208626 |
| marker2 | k__Bacteria|p__Bacteroidetes|c__Bacteroidia|o__Bacteroidales|f__Bacteroidaceae|g__Bacteroides | AL | 3.738189 | 0.0274232 |
| marker3 | k__Bacteria|p__Bacteroidetes|c__Bacteroidia|o__Bacteroidales|f__Bacteroidaceae | AL | 3.729778 | 0.0356919 |
| marker4 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Ruminococcaceae|g__Intestinimonas | AL | 3.173991 | 0.0356919 |
| marker5 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Ruminococcaceae|g__Phocea | AL | 2.507475 | 0.0017651 |
| marker6 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Ruminococcaceae|g__Anaerotruncus | AL | 2.384924 | 0.0248608 |
| marker7 | k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales|f__Lactobacillaceae|gLactobacillaceae_g | AL | 2.314625 | 0.0055780 |
| marker8 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Lachnospiraceae|g__Eisenbergiella | AL | 2.148999 | 0.0272785 |
| marker9 | k__Bacteria|p__Firmicutes|c__Bacilli | BL | 4.884616 | 0.0274232 |
| marker10 | k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales | BL | 4.874638 | 0.0274232 |
| marker11 | k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales|f__Streptococcaceae | BL | 4.863216 | 0.0208626 |
| marker12 | k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales|f__Streptococcaceae|g__Lactococcus | BL | 4.863216 | 0.0208626 |
| marker13 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Ruminococcaceae|gRuminococcaceae_g | BL | 4.536076 | 0.0459994 |
| marker14 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Lachnospiraceae|g__Schaedlerella | BL | 4.256198 | 0.0032759 |
| marker15 | k__Bacteria|p__Actinobacteria | BL | 3.631425 | 0.0086515 |
| marker16 | k__Bacteria|p__Actinobacteria|c__Coriobacteriia | BL | 3.631425 | 0.0086515 |
| marker17 | k__Bacteria|p__Actinobacteria|c__Coriobacteriia|o__Eggerthellales | BL | 3.597609 | 0.0208626 |
| marker18 | k__Bacteria|p__Actinobacteria|c__Coriobacteriia|o__Eggerthellales|f__Eggerthellaceae | BL | 3.597609 | 0.0208626 |
| marker19 | k__Bacteria|p__Actinobacteria|c__Coriobacteriia|o__Eggerthellales|f__Eggerthellaceae|g__Adlercreutzia | BL | 3.597609 | 0.0208626 |
| marker20 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Clostridiaceae_1 | BL | 3.493641 | 0.0272785 |
| marker21 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Clostridiaceae_1|g__Clostridium_sensu_stricto | BL | 3.493641 | 0.0272785 |
| marker22 | k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales|f__Enterococcaceae | BL | 3.428158 | 0.0117187 |
| marker23 | k__Bacteria|p__Firmicutes|c__Bacilli|o__Lactobacillales|f__Enterococcaceae|g__Enterococcus | BL | 3.428158 | 0.0117187 |
| marker24 | k__Bacteria|p__Firmicutes|c__Bacilli|o__Bacillales | BL | 3.243597 | 0.0045744 |
| marker25 | k__Bacteria|p__Firmicutes|c__Bacilli|o__Bacillales|f__Staphylococcaceae | BL | 3.243597 | 0.0045744 |
| marker26 | k__Bacteria|p__Firmicutes|c__Bacilli|o__Bacillales|f__Staphylococcaceae|g__Staphylococcus | BL | 3.243597 | 0.0045744 |
| marker27 | k__Bacteria|p__Firmicutes|c__c__ | BL | 3.234979 | 0.0459994 |
| marker28 | k__Bacteria|p__Firmicutes|c__c__|o__c___o__ | BL | 3.234979 | 0.0459994 |
| marker29 | k__Bacteria|p__Firmicutes|c__c__|o__c___o__|f__c___o___f__ | BL | 3.234979 | 0.0459994 |
| marker30 | k__Bacteria|p__Firmicutes|c__c__|o__c___o__|f__c___o___f__|g__c___o___f___g__ | BL | 3.234979 | 0.0459994 |
| marker31 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Ruminococcaceae|g__Monoglobus | BL | 2.826257 | 0.0458395 |
| marker32 | k__Bacteria|p__Actinobacteria|c__Coriobacteriia|oCoriobacteriia_o | BL | 2.747338 | 0.0117187 |
| marker33 | k__Bacteria|p__Actinobacteria|c__Coriobacteriia|oCoriobacteriia_o|fCoriobacteriia_o_f__ | BL | 2.747338 | 0.0117187 |
| marker34 | k__Bacteria|p__Actinobacteria|c__Coriobacteriia|oCoriobacteriia_o|fCoriobacteriia_o_f__|g__Coriobacteriia_o___f___g__ | BL | 2.747338 | 0.0117187 |
| marker35 | k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Lachnospiraceae|g__Blautia | BL | 2.336873 | 0.0459994 |
The marker table has 35 rows (for 35 markers)
There are 8 markers for AL group and 27 markers for BL group.
The most features that explain AL group are:
Undefined genus from the family Muribaculaceae
genus Bacteroides from the family Bacteroidaceae
family Bacteroidaceae
genus Intestinimonas from the family Ruminococcaceae
genus Phocea from the family Ruminococcaceae
genus Anaerotruncus from the family Ruminococcaceae
The most features that explain BL group are:
class Bacilli
order Lactobacillales
family Streptococcaceae
genus Lactococcus from the above family, order, and class
Visualizing marker table as histogram
plot_ef_bar(lef_out)
We have 35 bars as we have 35 markers, they are colored by group.
Visualizing the marker table as cladogram
# kingdom, Phylum, Class
plot_cladogram(lef_out, color = c("red","blue"), clade_label_level = 4)
# all except genus
plot_cladogram(lef_out, color = c("red","blue"), clade_label_level = 2)
# phylum and kingdom
plot_cladogram(lef_out, color = c("red","blue"), clade_label_level = 5)
Specifying the lefse analysis to the desired taxonomy rank
For class
lef_out<-run_lefse(phy_low, group = "group", norm = "CPM",
taxa_rank = "Class",
kw_cutoff = 0.05, lda_cutoff = 2)
lef_out## microbiomeMarker-class inherited from phyloseq-class
## normalization method: [ CPM ]
## microbiome marker identity method: [ lefse ]
## marker_table() Marker Table: [ 3 microbiome markers with 5 variables ]
## otu_table() OTU Table: [ 10 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 10 taxa by 1 taxonomic ranks ]dat<-marker_table(lef_out) %>% data.frame() %>% select(1:4)
head(dat)## feature enrich_group ef_lda pvalue
## marker1 Bacilli BL 4.873595 0.027423154
## marker2 Coriobacteriia BL 3.660484 0.008651542
## marker3 c__ BL 3.411810 0.045999367dat %>% filter(enrich_group=="BL") %>% head()## feature enrich_group ef_lda pvalue
## marker1 Bacilli BL 4.873595 0.027423154
## marker2 Coriobacteriia BL 3.660484 0.008651542
## marker3 c__ BL 3.411810 0.045999367dat %>% kable(align = "c")| feature | enrich_group | ef_lda | pvalue | |
|---|---|---|---|---|
| marker1 | Bacilli | BL | 4.873595 | 0.0274232 |
| marker2 | Coriobacteriia | BL | 3.660484 | 0.0086515 |
| marker3 | c__ | BL | 3.411810 | 0.0459994 |
plot_ef_bar(lef_out)
Only 3 markers (classes) are present, all for the BL group.
For order
lef_out<-run_lefse(phy_low, group = "group", norm = "CPM",
taxa_rank = "Order",
kw_cutoff = 0.05, lda_cutoff = 2)
lef_out## microbiomeMarker-class inherited from phyloseq-class
## normalization method: [ CPM ]
## microbiome marker identity method: [ lefse ]
## marker_table() Marker Table: [ 5 microbiome markers with 5 variables ]
## otu_table() OTU Table: [ 13 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 13 taxa by 1 taxonomic ranks ]dat<-marker_table(lef_out) %>% data.frame() %>% select(1:4)
head(dat)## feature enrich_group ef_lda pvalue
## marker1 Lactobacillales BL 4.861971 0.027423154
## marker2 Eggerthellales BL 3.552740 0.020862582
## marker3 Bacillales BL 3.278302 0.004574439
## marker4 c___o__ BL 3.216972 0.045999367
## marker5 Coriobacteriia_o__ BL 2.947596 0.011718686dat %>% filter(enrich_group=="BL") %>% head()## feature enrich_group ef_lda pvalue
## marker1 Lactobacillales BL 4.861971 0.027423154
## marker2 Eggerthellales BL 3.552740 0.020862582
## marker3 Bacillales BL 3.278302 0.004574439
## marker4 c___o__ BL 3.216972 0.045999367
## marker5 Coriobacteriia_o__ BL 2.947596 0.011718686dat %>% kable(align = "c")| feature | enrich_group | ef_lda | pvalue | |
|---|---|---|---|---|
| marker1 | Lactobacillales | BL | 4.861971 | 0.0274232 |
| marker2 | Eggerthellales | BL | 3.552740 | 0.0208626 |
| marker3 | Bacillales | BL | 3.278302 | 0.0045744 |
| marker4 | c_o | BL | 3.216972 | 0.0459994 |
| marker5 | Coriobacteriia_o__ | BL | 2.947596 | 0.0117187 |
plot_ef_bar(lef_out)
Only 5 markers (orders) are present, all for the BL group.
For family
lef_out<-run_lefse(phy_low, group = "group", norm = "CPM",
taxa_rank = "Family",
kw_cutoff = 0.05, lda_cutoff = 2)
lef_out## microbiomeMarker-class inherited from phyloseq-class
## normalization method: [ CPM ]
## microbiome marker identity method: [ lefse ]
## marker_table() Marker Table: [ 8 microbiome markers with 5 variables ]
## otu_table() OTU Table: [ 24 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 24 taxa by 1 taxonomic ranks ]dat<-marker_table(lef_out) %>% data.frame() %>% select(1:4)
head(dat)## feature enrich_group ef_lda pvalue
## marker1 Bacteroidaceae AL 3.762836 0.035691900
## marker2 Streptococcaceae BL 4.856424 0.020862582
## marker3 Clostridiaceae_1 BL 3.581284 0.027278512
## marker4 Eggerthellaceae BL 3.546469 0.020862582
## marker5 Enterococcaceae BL 3.420817 0.011718686
## marker6 Staphylococcaceae BL 3.262017 0.004574439dat %>% filter(enrich_group=="BL") %>% head()## feature enrich_group ef_lda pvalue
## marker2 Streptococcaceae BL 4.856424 0.020862582
## marker3 Clostridiaceae_1 BL 3.581284 0.027278512
## marker4 Eggerthellaceae BL 3.546469 0.020862582
## marker5 Enterococcaceae BL 3.420817 0.011718686
## marker6 Staphylococcaceae BL 3.262017 0.004574439
## marker7 c___o___f__ BL 3.219197 0.045999367dat %>% kable(align = "c")| feature | enrich_group | ef_lda | pvalue | |
|---|---|---|---|---|
| marker1 | Bacteroidaceae | AL | 3.762836 | 0.0356919 |
| marker2 | Streptococcaceae | BL | 4.856424 | 0.0208626 |
| marker3 | Clostridiaceae_1 | BL | 3.581284 | 0.0272785 |
| marker4 | Eggerthellaceae | BL | 3.546469 | 0.0208626 |
| marker5 | Enterococcaceae | BL | 3.420818 | 0.0117187 |
| marker6 | Staphylococcaceae | BL | 3.262017 | 0.0045744 |
| marker7 | c_o_f__ | BL | 3.219197 | 0.0459994 |
| marker8 | Coriobacteriia_o_f | BL | 2.953295 | 0.0117187 |
plot_ef_bar(lef_out)
Only 8 markers (families) are present, 1 for the AL group and 7 for the BL group.
For genus
lef_out<-run_lefse(phy_low, group = "group", norm = "CPM",
taxa_rank = "Genus",
kw_cutoff = 0.05, lda_cutoff = 2)
lef_out## microbiomeMarker-class inherited from phyloseq-class
## normalization method: [ CPM ]
## microbiome marker identity method: [ lefse ]
## marker_table() Marker Table: [ 17 microbiome markers with 5 variables ]
## otu_table() OTU Table: [ 70 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 70 taxa by 1 taxonomic ranks ]dat<-marker_table(lef_out) %>% data.frame() %>% select(1:4)
head(dat)## feature enrich_group ef_lda pvalue
## marker1 Muribaculaceae_g__ AL 4.116433 0.020862582
## marker2 Bacteroides AL 3.855283 0.027423154
## marker3 Intestinimonas AL 3.241611 0.035691900
## marker4 Phocea AL 2.536827 0.001765096
## marker5 Anaerotruncus AL 2.376809 0.024860789
## marker6 Lactobacillaceae_g__ AL 2.325802 0.005577994dat %>% filter(enrich_group=="BL") %>% head()## feature enrich_group ef_lda pvalue
## marker7 Lactococcus BL 4.839021 0.020862582
## marker8 Ruminococcaceae_g__ BL 4.477288 0.045999367
## marker9 Schaedlerella BL 4.263494 0.003275897
## marker10 Clostridium_sensu_stricto BL 3.556176 0.027278512
## marker11 Adlercreutzia BL 3.528704 0.020862582
## marker12 Enterococcus BL 3.464011 0.011718686dat %>% kable(align = "c")| feature | enrich_group | ef_lda | pvalue | |
|---|---|---|---|---|
| marker1 | Muribaculaceae_g__ | AL | 4.116433 | 0.0208626 |
| marker2 | Bacteroides | AL | 3.855283 | 0.0274232 |
| marker3 | Intestinimonas | AL | 3.241611 | 0.0356919 |
| marker4 | Phocea | AL | 2.536827 | 0.0017651 |
| marker5 | Anaerotruncus | AL | 2.376809 | 0.0248608 |
| marker6 | Lactobacillaceae_g__ | AL | 2.325802 | 0.0055780 |
| marker7 | Lactococcus | BL | 4.839021 | 0.0208626 |
| marker8 | Ruminococcaceae_g__ | BL | 4.477288 | 0.0459994 |
| marker9 | Schaedlerella | BL | 4.263494 | 0.0032759 |
| marker10 | Clostridium_sensu_stricto | BL | 3.556176 | 0.0272785 |
| marker11 | Adlercreutzia | BL | 3.528704 | 0.0208626 |
| marker12 | Enterococcus | BL | 3.464011 | 0.0117187 |
| marker13 | Staphylococcus | BL | 3.292156 | 0.0045744 |
| marker14 | c_o_f_g | BL | 3.250737 | 0.0459994 |
| marker15 | Monoglobus | BL | 2.866083 | 0.0458395 |
| marker16 | Coriobacteriia_o_f_g__ | BL | 2.827921 | 0.0117187 |
| marker17 | Blautia | BL | 2.282599 | 0.0459994 |
plot_ef_bar(lef_out)
Only 17 markers (genus) are present, 6 for the AL group and 11 for the BL group.
Focusing on the most prominent features
We decrease the error rate to 0.01 and increase the lda score to 4.
This approach is suggested by the author of the microbiomeMarker package.
For all ranks, “Kingdom” “Phylum” “Class” “Order” “Family” “Genus”
lef_out<-run_lefse(phy_low, group = "group", norm = "CPM",
kw_cutoff = 0.01, lda_cutoff = 4)
lef_out## microbiomeMarker-class inherited from phyloseq-class
## normalization method: [ CPM ]
## microbiome marker identity method: [ lefse ]
## marker_table() Marker Table: [ 1 microbiome markers with 5 variables ]
## otu_table() OTU Table: [ 123 taxa and 16 samples ]
## sample_data() Sample Data: [ 16 samples by 4 sample variables ]
## tax_table() Taxonomy Table: [ 123 taxa by 1 taxonomic ranks ]There is only 1 microbiome marker.
Microbiome marker table
library(knitr)
dat<-marker_table(lef_out) %>% data.frame() %>% select(1:4)
head(dat)## feature
## marker1 k__Bacteria|p__Firmicutes|c__Clostridia|o__Clostridiales|f__Lachnospiraceae|g__Schaedlerella
## enrich_group ef_lda pvalue
## marker1 BL 4.209894 0.003275897The marker table has 1 row (for 1 marker)
This marker corresponds to the BL group.
It is:
genus Schaedlerella from family Lachnospiraceae, order Clostridiales, and class Clostridia.
Visualizing marker table as histogram
plot_ef_bar(lef_out)
Visualizing the marker table as cladogram
# kingdom, Phylum, Class
plot_cladogram(lef_out, color = c("red"), clade_label_level = 4)
The white circles are non significant taxa.
Specifying the lefse analysis to the desired taxonomy rank
There is no need to this as we already know this marker.