Closed kenbyrrrd closed 1 year ago
Update: I have now changed taxa CSV files to look identical to american gut project data. Still returns a fatal error in r. At this point, I am at a loss as what to do.
Filtered data for OTUs under 1% and associated taxa names, then was able to run network
Hi, I am trying to use the SpecEasi package on a phyloseq object. I successfully get the package to run on the amgut2.filt.phy properly, but when I apply the package to my data, it runs for a while then R force quits, and restarts the session stating it encountered a fatal error.
I am unsure what I am doing incorrectly, as the phyloseq object is created with no error. I have tried creating a phyloseq object multiple ways, including using CSV (shown below) and QZA files directly from QIIME2. I have also tried transposing my OTU table data, although I believe this is not an issue for phyloseq.
I have attached all my files here (except QZA) and posted both code snippets below. Thanks for the help in advance!
CODE for working with CSV:
QIIME2 artifacts import to R for phyloseq
setwd("~/Grad_Project/SchillerPark/Network_Analysis/Network_code")
install packages
library(tidyverse)
qiime2r
if (!requireNamespace("devtools", quietly = TRUE)){install.packages("devtools")} devtools::install_github("jbisanz/qiime2R") library(qiime2R) library(devtools) library(plotly) install_github("jbisanz/MicrobeR") library(magrittr) library(phyloseq) library(ggplot2) library(vegan) install_github("zdk123/SpiecEasi") library(SpiecEasi)
SP_Meta1 <- read.csv("SchillerProject_Meta_R.csv", row.names = 1) SP_OTU1 <- read.csv("OTU-Table-Schiller.csv", row.names = 1) SP_taxa1 <- read.csv("Taxonomy-Table-Schiller.csv", row.names = 1)
Save unaltered data by changing the name for analysis
SP_Meta <- SP_Meta1 SP_OTU <- SP_OTU1 SP_taxa <- SP_taxa1
Set up phyloseq
otu_mat <- as.matrix(SP_OTU) tax_mat <- as.matrix(SP_taxa)
OTU = otu_table(otu_mat, taxa_are_rows = TRUE) TAX = tax_table(tax_mat) META = sample_data(SP_Meta) phylo <- phyloseq(OTU, TAX, META)
network analysis
SP_mb <- spiec.easi(phylo, method='mb', lambda.min.ratio=1e-2, nlambda=20)
CODE for QZA:
QIIME2 artifacts import to R for phyloseq
setwd("~/Grad_Project/SchillerPark/Network_Analysis/Network_code")
install packages
library(tidyverse)
qiime2r
if (!requireNamespace("devtools", quietly = TRUE)){install.packages("devtools")} devtools::install_github("jbisanz/qiime2R") library(qiime2R) library(devtools) library(plotly) install_github("jbisanz/MicrobeR") library(magrittr) library(phyloseq) library(ggplot2) library(vegan) install_github("zdk123/SpiecEasi") library(SpiecEasi) install_github("hallucigenia-sparsa/seqtime") library(seqtime)
test datasets
data(amgut1.filt) data(amgut2.filt.phy)
Load phyloseq
initial_phyloseq_obj <- qza_to_phyloseq(features = "Filtered_Schiller_table.qza", tree = "Filtered_Schiller_rootedtree.qza", taxonomy = "Filtered_Schiller_taxonomy.qza", metadata = "SchillerProject_Meta_R.tsv") otus <- otu_table(initial_phyloseq_obj) otus
taxa <- tax_table(initial_phyloseq_obj)
check test data structure for comparison to personal data
sample_names(amgut2.filt.phy) rank_names(amgut2.filt.phy) sample_variables(amgut2.filt.phy)
Check personal data structure
sample_names(initial_phyloseq_obj) rank_names(initial_phyloseq_obj) sample_variables(initial_phyloseq_obj)
transpose OTU data
updated_otus <- t(otu_table(initial_phyloseq_obj)) phyloseq_obj <- phyloseq(updated_otus, taxa)
Normalization is handled inside of speceasi function so no need to run below code
normalize to the same threshold as QIIME
set.seed(111995)
Schiller_phylo = rarefy_even_depth(Schiller_phylo, sample.size = 19496)
network analysis
se.mb.schiller <- spiec.easi(phyloseq_obj, method='mb', lambda.min.ratio=1e-2, nlambda=20)
Taxonomy-Table-Schiller.csv OTU-Table-Schiller.csv SchillerProject_Meta_R.csv