convert mirlyn object back to phyloseq

[PedroDaPos]

Hi! First, thank you for creating this package, it is really helpful. I was wondering if I could ask for your help with the following. I just rarefied my data using the function mirl, and now I would like to transform the resulting mirl object back to a phyloseq object. There are a few steps in my downstream analysis that I can not do using the mirl object (i.e. permanova test, beta diversity with unifrac distances, etc). Would you be able to explain to me how to transform the mirl object back to a phyloseq object? I am wondering if i could get an average of the iterations for each observed ASV, if possible.

[mailasb]

Hello again. Were you able to resolve this issue? I completed the entire Mirlyn workflow with my data, but I still have this question. The article suggests multiple rarefaction as a normalization tool. Can I obtain a phyloseq object from the mirl object rarefied 1000 times at the same chosen depth, or would the mirl object be like multiple phyloseq objects? In other words, can I deviate from the Mirlyn workflow to perform other analyses using the mirl object? Thanks in advance. Maila

[PedroDaPos]

Hi Maila, unfortunately I was not able to resolve this. I believe Pat Schloss has some videos on youtube where he shows how to do repeated rarefactions using R. Please let me know if you find a solution.

[escamero]

We will be looking into adding the functionality for this later this month hopefully - apologies for the delay in this we have both just been incredibly busy with other projects and research.

In the meantime, a solution that might work (I have not tried it at a larger scale yet) although is not the prettiest is below. Basically just creating average counts of all the repetitions for each sample, and then feeding that into the otu table position of a phyloseq object:

mirl <- mirl(example, rep = 10)
mirl_otu <- vector("list", length(mirl))

for (i in 1:length(mirl)){
  colnames(mirl[[i]]@otu_table) <- paste0(colnames(mirl[[i]]@otu_table), "_",i)
  (mirl_otu[[i]] <- mirl[[i]]@otu_table)
}

mirl_rep_df <- do.call(cbind, mirl_otu)

example_id <- example@sam_data$Id

average_counts <- vector("list", length(example_id))

for (i in 1:length(example_id)){
  sample_df <- mirl_rep_df[,grepl(example_id[[i]], colnames(mirl_rep_df))]
  sample_average <- data.frame(rowMeans(sample_df))
  colnames(sample_average) <- example_id[[i]]
  average_counts[[i]] <- sample_average
}
average_count_df <- do.call(cbind, average_counts)

mirl_phyloseq <- phyloseqobject
mirl_phyloseq@otu_table@.Data <- average_count_df

[mailasb]

Hello everyone. It seems that everything works fine up to the second-to-last line. When you write: mirl_phyloseq <- phyloseqobject, does it mean that I have to rename the phyloseq object (in the example it would be the "example" object) as mirl_phyloseq?"

[escamero]

Hi - yes sorry, the second last line is just creating another phyloseq object and should be as follows for the example data:

mirl_phyloseq <- example

You could also just do the last line like:

example@otu_table@.Data <- average_count_df

but then this will get rid of your raw count data which is why I had the extra bit of creating a new phyloseq object just for the mirl count data.

[mailasb]

Excellent! "The only issue encountered was that I had to transform the data.frame "average_count_df" into a matrix in order to run: mirl_phyloseq@otu_table@.Data <- average_count_df since the class of "mirl_phyloseq@otu_table@.Data" is of type 'matrix' 'array', but it works great!" Thank you so much!