exporting finalization

DESCRIPTION

@@ -7,6 +7,6 @@ Description: Analysis toolkit for intra-species diversity from metagenomics.
 Authors@R: person("Noriaki", "Sato", email = "nori@hgc.jp", role = c("cre", "aut"))
 Depends: ggplot2, ggstar, ggraph, igraph
 Imports: GetoptLong, BiocFileCache, RCurl, vegan, methods, data.table, phangorn, RColorBrewer, ggtree, circlize, ComplexHeatmap, ggkegg, ape, dplyr, exactRankTests, ggblend, ggh4x, scales, tidygraph, ggplotify, ggtreeExtra, ggnewscale, scico, MKmisc, NMF, pillar, BiocStyle, cowplot, patchwork, reshape2, ggrepel, Boruta, tidyr, stringr
-Suggests: simplifyEnrichment, knitr, rmarkdown, NNLM
+Suggests: simplifyEnrichment, knitr, rmarkdown, NNLM, shiny
 RoxygenNote: 7.3.1
 VignetteBuilder: knitr

R/exportInteractive.R

@@ -11,17 +11,138 @@
 #' @param calctree calculate consensus tree
 #' @param species candidate species, default to all the species
 #' @param dataset_name dataset name
+#' @param notRun do not run the shiny app immediately
+#' @param clear when exporting clear the gene and snp slot to reduce file size
 #' @export
 #' @return output the results to specified directory
 exportInteractive <- function(stana, out=".", db="uhgg", calcko=FALSE,
+	calctree=FALSE, dataset_name=NULL, clear=TRUE,
+	species=NULL, notRun=FALSE) {
+	if (is.null(dataset_name)) {
+		dataset_name <- gsub(":", "-", gsub(" ", "-", as.character(Sys.time())))
+	}
+	dir.create(paste0(out,"/data"))
+	if (length(stana@names)==0) {
+		stana@names <- stana@ids %>% setNames(stana@ids)
+	}
+	if (is.null(species)) {
+		species <- stana@ids	
+	}
+	if (calcko) {
+		for (sp in species) {
+	        if (is.null(stana@kos[[sp]])) {
+	        	if (!is.null(stana@eggNOG[[sp]])) {
+		        	cat_subtle("# Summarizing abundances ", sp, "\n")
+		            ko_tbl <- summariseAbundance(stana,sp = sp,
+		                        checkEGGNOG(annot_file=stana@eggNOG[[sp]],
+		                        	"KEGG_ko"),
+		                        how="mean")
+		            stana@kos[[sp]] <- ko_tbl	        		
+	        	}
+	        } else {
+	            cat_subtle("# Using pre-computed KO table\n")
+	        }
+	    }
+	} else {
+		for (sp in species) {
+			if (is.null(stana@treeList[[sp]])) {
+				if (is.null(stana@genes[[sp]])) {
+					cat_subtle("# Need at least the gene copy number table if tree is not available\n")
+				}
+			}
+	        if (is.null(stana@kos[[sp]])) {
+	        	## Insert gene table instead
+	        	stana@kos[[sp]] <- stana@genes[[sp]]
+	        }
+	    }
+	}
+
+	for (tre in species) {
+		if (is.null(stana@treeList[[tre]])) {
+			cat_subtle("# No tree for ", tre, "\n")
+			if (calctree) {
+				cat_subtle("Calculating ... \n")
+				stana <- consensusSeqMIDAS2(stana, tre)
+                stana <- inferAndPlotTree(stana, tre)
+			}
+		}
+	}
+	if (length(stana@treeList)!=0) {
+		treLen <- sum(sapply(stana@treeList, function(x) !is.null(x)))	
+	} else {
+		treLen <- 0
+	}
+	koLen <- length(stana@kos)
+
+    all_samples_in_dataset1 <- unique(lapply(stana@treeList, function(x) {x$tip.label}) %>% unlist())
+    all_samples_in_dataset2 <- unique(lapply(stana@kos, function(x) {colnames(x)}) %>% unlist())
+    all_samples_in_dataset <- union(all_samples_in_dataset1, all_samples_in_dataset2)
+
+	## If no metadata is available in meta slot:
+    if (dim(stana@meta)[1]==0) {
+        if (length(stana@cl)==0) {
+            meta <- data.frame(all_samples_in_dataset) %>% `colnames<-`(c("samples"))
+            meta[["label"]] <- meta[,1]
+            row.names(meta) <- meta[,1]
+        } else {
+            meta <- NULL
+            for (i in names(stana@cl)) {
+                meta <- rbind(meta,
+                    cbind(stana@cl[[i]],
+                        rep(i, length(stana@cl[[i]]))))
+            }
+            meta <- meta |> data.frame() |> `colnames<-`(c("label","group"))
+            row.names(meta) <- meta[["label"]]            
+        }
+    } else {
+        meta <- stana@meta
+    }
+
+	cat_subtle("# Tree number: ", treLen, " KO (or gene) number: ", koLen,  "\n")
+	cat_subtle("# Exporting ... \n")
+	if (clear) {
+		stana@genes <- list()
+		stana@snps <- list()
+		stana@snpsInfo <- list()
+	}
+	save(file=paste0(out, "/data/", dataset_name,".rda"), stana, compress="xz")
+    ## Copy the main script and run the app
+    file.copy(system.file("extdata", "app_stana.R", package = "stana"), out)
+    if (requireNamespace("shiny")) {
+    	if (!notRun) {
+            shiny::runApp(paste0(out,"/app_stana.R"))
+    	}
+    }
+	return(stana)
+}
+
+#' exportInteractiveRawFiles
+#' 
+#' export the current stana object to interactive application 
+#' for the convenient analysis and visualization
+#' if `cl` and `meta` slot is filled, `meta` slot is exported.
+#' This outputs the raw files in TSV, not stana object itself.
+#' 
+#' @param stana stana object of type MIDAS2
+#' @param out output directory
+#' @param db db used to profile 'uhgg' or 'gtdb'
+#' @param calcko calculate KO abundance
+#' @param calctree calculate consensus tree
+#' @param species candidate species, default to all the species
+#' @param dataset_name dataset name
+#' @noRd
+#' @return output the results to specified directory
+exportInteractiveRawFiles <- function(stana, out=".", db="uhgg", calcko=FALSE,
 	calctree=FALSE, dataset_name=NULL,
 	species=NULL) {
 	if (is.null(dataset_name)) {
 		dataset_name <- gsub(":", "-", gsub(" ", "-", as.character(Sys.time())))
 	}
 	dir.create(paste0(out,"/data"))
 	dir.create(paste0(out,"/data/",dataset_name))
-
+	if (length(stana@names)==0) {
+		stana@names <- stana@ids %>% setNames(stana@ids)
+	}
 	# if (stana@type!="MIDAS2") {stop("This feature is for MIDAS2 only")}
 	if (is.null(species)) {
 		species <- stana@ids	

R/loadMIDAS2Refine.R

@@ -217,6 +217,7 @@ loadMIDAS2 <- function(midas_merge_dir,
     stana@snps <- snpsList
     stana@genes <- geneList
     stana@ids <- union(names(geneList),names(snpsList))
+    stana@names <- loadDic()[[db]][stana@ids] %>% setNames(stana@ids)
     stana <- initializeStana(stana,cl)
 
 

R/stana.R

@@ -2,6 +2,8 @@
 #' @import methods
 #' @importFrom utils object.size
 #' @importFrom grDevices colorRampPalette
+#' @slot ids identifiers distinguishing species
+#' @slot names name of species (named vector)
 #' @slot snpsInfo snp info
 #' @slot snpsDepth snp depth
 #' @slot snpsSummary snp summary
@@ -19,6 +21,7 @@ setClass("stana", slots=list(
                             mergeDir="character",
                             db="character",
                             ids="character",
+                            names="character",
                             comparisonTable="data.table",
                             genomeWideCompare="data.table",
                             strainClusters="data.table",