exportInteractive

NAMESPACE

@@ -16,6 +16,7 @@ export(doAdonis)
 export(doBoruta)
 export(drawEGGNOG)
 export(drawPATRIC)
+export(exportInteractive)
 export(genomeHeatmap)
 export(loadInStrain)
 export(loadMIDAS)
@@ -27,6 +28,7 @@ export(plotGenes)
 export(plotHeatmap)
 export(plotKEGGPathway)
 export(plotPCA)
+export(plotTree)
 export(setAnnotation)
 export(strainClusterHeatmap)
 export(summariseAbundance)

R/checkEGGNOG.R

@@ -80,9 +80,9 @@ drawEGGNOG <- function(annot_file, geneIDs, candPlot) {
 #' @importFrom data.table fread
 #' @export
 #' 
-checkEGGNOG <- function(annot_file, ret="all", checkIDs=NULL) {
+checkEGGNOG <- function(annot_file, ret="all", checkIDs=NULL, fill=TRUE) {
   ann <- data.table::fread(annot_file,
-                           skip=4,sep="\t")
+                           skip=4,sep="\t",fill=fill)
   if (ret!="all") {
     parsed <- ann[,c("#query",ret), with=FALSE] |>
       tidyr::pivot_longer(-1) |>
@@ -211,7 +211,8 @@ checkEGGNOG <- function(annot_file, ret="all", checkIDs=NULL) {
 summariseAbundance <- function(stana, sp, anno, how=mean) {
   geneDf <- stana@genes[[sp]]
   merged <- list()
-  for (i in anno$value |> unique()) {
+  annoval <- anno$value |> unique()
+  for (i in annoval) {
     candID <- (anno |> dplyr::filter(anno$value==i))$ID
     ints <- intersect(row.names(geneDf), candID)
     if (length(ints)>0) {

R/exportInteractive.R

@@ -0,0 +1,69 @@
+#' exportInteractive
+#' 
+#' @param stana stana object
+#' @param out output directory
+#' export the stana results for the inspection in interactive
+#' application
+#' @export
+#' @return output the results to specified directory
+exportInteractive <- function(stana, out=".", db="uhgg", calcko=TRUE,
+	calctree=TRUE,
+	species=NULL) {
+	if (stana@type!="MIDAS2") {stop("This feature is for MIDAS2 only")}
+	if (is.null(species)) {
+		species <- stana@ids	
+	}
+	if (calcko) {
+		for (sp in species) {
+	        if (is.null(stana@kos[[sp]])) {
+	        	qqcat("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 {
+	            qqcat("Using pre-computed KO table\n")
+	        }
+	    }
+	}
+
+	for (tre in names(stana@treeList)) {
+		if (is.null(stana@treeList[[tre]])) {
+			qqcat("No tree for @{tre}\n")
+			if (calctree) {
+				qqcat("Calculating ... \n")
+				stana <- consensusSeqMIDAS2(stana, tre, tree=TRUE,
+					max_sites=50) # temp	
+			}
+		}
+	}
+	treLen <- sum(sapply(stana@treeList, function(x) !is.null(x)))
+	koLen <- length(stana@kos)
+	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"))
+	qqcat("Tree number: @{treLen}, KO number: @{koLen}\n")
+	qqcat("Exporting ... \n")
+	if (!dir.exists(out)) {
+		dir.create(out)
+	}
+	dir.create(paste0(out,"/tree"))
+	for (i in names(stana@treeList)) {
+		tre <- stana@treeList[[i]]
+		if (!is.null(tre)) {
+			ape::write.tree(tre, paste0(out,"/tree/",i,".cons.tree"))
+		}
+	}
+	dir.create(paste0(out,"/KOs"))
+	for (i in names(stana@kos)) {
+		ko <- stana@kos[[i]]
+		write.table(ko, paste0(out, "/KOs/", i, ".txt"), sep="\t", quote=FALSE)
+	}
+	write.table(meta, paste0(out, "metadata.tsv"), sep="\t", quote=FALSE)
+	return(stana)
+}