Question

Heatmap of mothur output with Phyloseq

0

Entering edit mode

4.5 years ago

dpc ▴ 250

Hi community!!! I want to visualise my data output analysed by mothur pipeline with Phyloseq package. I have gone through the phyloseq tutorial but still confused as I have no idea with R (Although I've successfully installed the necessary packages). I don't know which outputs from Mothur should I use as input for phyloseq and how to process them before running with phyloseq. Can anyone please help me. I am really struggling with this.

Thanks and regards, DC7

phyloseq mothur 16s metagenome • 2.3k views

ADD COMMENT • link 4.5 years ago by dpc ▴ 250

0

Entering edit mode

Hi DC,

here is a tutorial on how to build a phyloseq object from mothur files.

Hope this helps!

ADD REPLY • link 4.5 years ago by andres.firrincieli 3.8k

0

Entering edit mode

Thanks andres!!! Can you please help me with this error? I have imported mothur files with these codes:

 physeq = import_mothur(mothur_shared_file = "stability.paired.trim.contigs.unique.good.filter.unique.precluster.pick.pick.opti_mcc.shared",
                               mothur_constaxonomy_file = "stability.paired.trim.contigs.unique.good.filter.unique.precluster.pick.pick.opti_mcc.0.03.cons.taxonomy")

But, after head(phyloseq) , I am getting the following ERROR:

Error in x[seq_len(n)] : object of type 'S4' is not subsettable

Thanks and Regards, DC7

ADD REPLY • link 4.5 years ago by dpc ▴ 250

0

Entering edit mode

you can't use head() on a S4 object. To inspect the OTU table simply use otu_table(physeq). For the taxonomy tax_table(physeq).

Source

ADD REPLY • link 4.5 years ago by andres.firrincieli 3.8k

score 0 · Answer 1 · 2020-06-22

0

Entering edit mode

4.5 years ago

dpc ▴ 250

Thanks a lot andres.firrincieli. Can you please check out this query of mine? Link to question

ADD COMMENT • link 4.5 years ago by dpc ▴ 250

0

Entering edit mode

Hi andres.firrincieli, can I draw your attention for an isuue I am facing? I am following the tutorial you have provided. I am here at the "a stacked barplot for the class level" of the tutorial. I accordingly used the codes as in tutorial:

    # melt to long format (for ggploting) 
# prune out phyla below 1% in each sample
# selecting the taxa at the level: Phylum

mdata_class <- mouse_data %>%
  tax_glom(taxrank = "Class") %>%                     # agglomerate at phylum level
  transform_sample_counts(function(x) {x/sum(x)} ) %>% # Transform to rel. abundance
  psmelt() %>%                                         # Melt to long format
  filter(Abundance > 0.01) %>%                         # Filter out low abundance taxa
  arrange(Class)                                      # Sort data frame alphabetically by phylum

str(mdata_class)

And the expected outcome supposed be like:

    ## 'data.frame':    77 obs. of  10 variables:
##  $ OTU        : chr  "Otu0059" "Otu0010" "Otu0010" "Otu0010" ...
##  $ Sample     : chr  "F3D142" "F3D144" "F3D3" "F3D6" ...
##  $ Abundance  : num  0.0119 0.0949 0.0765 0.0747 0.0734 ...
##  $ group      : Factor w/ 19 levels "F3D0","F3D1",..: 4 6 14 16 11 7 9 2 4 18 ...
##  $ dpw        : int  142 144 3 6 149 145 147 1 142 8 ...
##  $ time       : chr  "late" "late" "early" "early" ...
##  $ body_weight: num  40.2 41.6 15.5 20.5 40.6 ...
##  $ Kingdom    : Factor w/ 1 level "Bacteria": 1 1 1 1 1 1 1 1 1 1 ...
##  $ Phylum     : Factor w/ 9 levels "Actinobacteria",..: 1 5 5 5 5 5 5 5 5 5 ...
##  $ Class      : Factor w/ 15 levels "Actinobacteria",..: 1 2 2 2 2 2 2 2 2 2 ...

Here clearly number of different kingdom, phylums, classes is shown which is 1, 9, 15 respectively. But in my dataset I am using the same code, but the output doesn't show number of different classes. As a result I am unable to use colours for them:

 str(mdata_phylum)
'data.frame':   37 obs. of  9 variables:
 $ OTU       : chr  "Otu00009" "Otu00009" "Otu00009" "Otu00009" ...
 $ Sample    : chr  "DRR046797" "DRR046803" "DRR046802" "DRR046805" ...
 $ Abundance : num  0.1336 0.1205 0.0908 0.0899 0.0893 ...
 $ sample    : chr  "DRR046797" "DRR046803" "DRR046802" "DRR046805" ...
 $ type      : chr  "control" "obese" "obese" "obese" ...
 $ time      : chr  "early" "late" "late" "early" ...
 $ bodyweight: num  20 45 45 78 40 45 45 45 32 45 ...
 $ Kingdom   : chr  "k__Bacteria" "k__Bacteria" "k__Bacteria" "k__Bacteria" ...
 $ Phylum    : chr  "p__Actinobacteria" "p__Actinobacteria" "p__Actinobacteria" "p__Actinobacteria" ...

Can you please help me out?

Thanks and Regards, DC7

ADD REPLY • link 4.5 years ago by dpc ▴ 250

0

Entering edit mode

What is the output of head(stability.paired.trim.contigs.unique.good.filter.unique.precluster.pick.pick.opti_mcc.0.03.cons.taxonomy)

and tax_table(physeq)

ADD REPLY • link 4.5 years ago by andres.firrincieli 3.8k

score 0 · Answer 2 · 2020-06-23

I have changed my script from earlier. I think the following information may help:

#ASSIGN VARIABLES FOR IMPORTED DATA
sharedfile = "stability.paired.trim.contigs.unique.good.filter.unique.precluster.pick.pick.opti_mcc.shared"
taxfile = "stability.paired.trim.contigs.unique.good.filter.unique.precluster.pick.pick.opti_mcc.0.03.cons.taxonomy"
treefile <- "stability.paired.trim.contigs.unique.good.filter.unique.precluster.pick.pick.opti_mcc.0.03.rep.otu.phylip.tre" 
mapfile = "metadata.xlsx"

#IMPORT MOTHUR_DATA
mothur_data <- import_mothur(mothur_shared_file = sharedfile, 
                             mothur_constaxonomy_file = taxfile,
                             mothur_tree_file = treefile)

#IMPORT SAMPLE METADATA INTO R
map <- read_xlsx(mapfile)

As you asked for tax_table(physeq) which should be tax_table(mothur_data) here (a portion is given, full table is not being accepted in biostars):

>tax_table(mothur_data)
Taxonomy Table:     [23170 taxa by 7 taxonomic ranks]:
         Rank1         Rank2                      Rank3                           
Otu22126 "k__Bacteria" "p__Bacteroidetes"         "c__Bacteroidia"                
Otu13776 "k__Bacteria" "p__Bacteroidetes"         "c__Bacteroidia"                
Otu13459 "k__Bacteria" "p__Bacteroidetes"         "c__Bacteroidia"                
Otu03511 "k__Bacteria" "p__Bacteroidetes"         "c__Bacteroidia"                
Otu18225 "k__Bacteria" "p__Bacteroidetes"         "c__Bacteroidia"                
Otu12185 "k__Bacteria" "p__Bacteroidetes"         "c__Bacteroidia"                
Otu08666 "k__Bacteria" "p__Bacteroidetes"         "c__Bacteroidia"                
Otu00807 "k__Bacteria" "p__Bacteroidetes"         "c__Bacteroidia"                
Otu18214 "k__Bacteria" "p__Bacteroidetes"         "c__Bacteroidia"                
Otu18675 "k__Bacteria" "p__Bacteroidetes"         "c__Bacteroidia"                

         Rank4                                 Rank5                                
Otu22126 "o__Bacteroidales"                    "f__Porphyromonadaceae"              
Otu13776 "o__Bacteroidales"                    "f__Porphyromonadaceae"              
Otu13459 "o__Bacteroidales"                    "f__Porphyromonadaceae"              
Otu03511 "o__Bacteroidales"                    "f__Porphyromonadaceae"              
Otu18225 "o__Bacteroidales"                    "f__Porphyromonadaceae"              
Otu12185 "o__Bacteroidales"                    "f__Porphyromonadaceae"              
Otu08666 "o__Bacteroidales"                    "f__Porphyromonadaceae"              
Otu00807 "o__Bacteroidales"                    "f__Porphyromonadaceae"              
Otu18214 "o__Bacteroidales"                    "f__Porphyromonadaceae"              
Otu18675 "o__Bacteroidales"                    "f__Porphyromonadaceae"              

         Rank6                                   Rank7                                         
Otu22126 "g__Parabacteroides"                    "s__distasonis"                               
Otu13776 "g__Parabacteroides"                    "s__distasonis"                               
Otu13459 "g__Parabacteroides"                    "s__distasonis"                               
Otu03511 "g__Parabacteroides"                    "s__distasonis"                               
Otu18225 "g__Parabacteroides"                    "s__distasonis"                               
Otu12185 "g__Parabacteroides"                    "s__distasonis"                               
Otu08666 "g__Parabacteroides"                    "s__distasonis"                               
Otu00807 "g__Parabacteroides"                    "s__distasonis"                               
Otu18214 "g__Parabacteroides"                    "s__distasonis"                               
Otu18675 "g__Parabacteroides"                    "s__distasonis"                               
     [ reached getOption("max.print") -- omitted 23028 rows ]

And here is another one you asked for:

head("stability.paired.trim.contigs.unique.good.filter.unique.precluster.pick.pick.opti_mcc.0.03.cons.taxonomy")
[1] "stability.paired.trim.contigs.unique.good.filter.unique.precluster.pick.pick.opti_mcc.0.03.cons.taxonomy"

Thanks and Regards, DC7