In fbreitwieser/pavian: Visualize and analyze metagenomics classification results

Most functions that are available in the Pavian shiny interface can also be used from the command line.

library(pavian)
sample_data <- pavian::read_sample_data(system.file("shinyapp","example-data","brain-biopsies",package="pavian"))
reports <- pavian::read_reports(sample_data$ReportFilePath, sample_data$Name)
merged_reports <- pavian::merge_reports2(reports, col_names = sample_data$Name)

pavian::summarize_reports(reports)

tax_data <- merged_reports[[1]]
clade_reads <- merged_reports[[2]]
taxon_reads <- merged_reports[[3]]

colSums(clade_reads,na.rm = T)
colSums(taxon_reads,na.rm = T)
sel_rows <- pavian::filter_taxa(tax_data,
                                rm_clades = c("Chordata", "artificial sequences", "unclassified"),
                                taxRank = "S")
summary(sel_rows)
filtered_clade_reads <- pavian::filter_cladeReads(clade_reads, tax_data, c("Chordata", "artificial sequences", "unclassified"))

tax_data1 <- tax_data[sel_rows,]
filtered_clade_reads1 <- filtered_clade_reads[sel_rows, ]
taxon_reads1 <- taxon_reads[sel_rows, ]

head(cbind(tax_data1[,1:3],clade_reads[sel_rows, ])[order(-apply(filtered_clade_reads1,1,max, na.rm=T)),])
normalized_clade_reads <- normalize(filtered_clade_reads1)
normalized_taxon_reads <- normalize(taxon_reads[sel_rows,], sums = colSums(filtered_clade_reads1,na.rm = T))
head(cbind(tax_data1[,1:3],max=apply(cbind(normalized_clade_reads),1,max, na.rm=T), normalized_clade_reads)[order(-apply(cbind(normalized_clade_reads),1,max, na.rm=T)),])

reads_zscore <- robust_zscore(100*cbind(normalized_clade_reads,normalized_taxon_reads), 0.001)
clade_reads_zscore <- reads_zscore[,1:20]
reads_zscore_df <- cbind(tax_data1[,1:3],max=apply(clade_reads_zscore,1,max, na.rm=T), clade_reads_zscore)[order(-apply(clade_reads_zscore,1,max, na.rm=T)),]
## Calculate z-score from the clade reads