diff phyloseq_net.r @ 0:af6d9ad14a0f draft

"planemo upload for repository https://github.com/QFAB-Bioinformatics/metaDEGalaxy/tree/master/phyloseq_net commit 8bd68662b72404f6291e9628327dcb109b5fa55e"
author qfabrepo
date Mon, 14 Sep 2020 08:13:43 +0000
parents
children
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/phyloseq_net.r	Mon Sep 14 08:13:43 2020 +0000
@@ -0,0 +1,147 @@
+library('getopt')
+library('data.table')
+suppressPackageStartupMessages(library('phyloseq'))
+suppressPackageStartupMessages(library('DESeq2'))
+
+Sys.setenv("DISPLAY"=":1")
+
+options(warn= -1)
+option_specification = matrix(c(
+   'infile','i',2,'character',
+   'metafile','m',2,'character',
+       'biom','b',2,'character',
+     'obsfile','o',2,'character',
+     'norm','n',2,'logical',
+     'xcolumn','x',2,'numeric',
+     'lcolumn','l',2,'numeric',
+     'outdir','d',2,'character',
+   'htmlfile','h',2,'character'
+),byrow=TRUE,ncol=4);
+
+
+options <- getopt(option_specification)
+options(bitmapType="cairo")
+
+matrix.format<-function(x) {
+    m<-as.matrix(x[,-1])
+    rownames(m)<-x[,1]
+    m
+}
+
+
+gm_mean = function(x, na.rm=TRUE){ exp(sum(log(x[x > 0]), na.rm=na.rm) / length(x))}
+
+
+tax_col_norm <- c("Kingdom","Phylum","Class","Order","Family","Genus","Species")
+tax_col_extra <- c("None","Kingdom","Phylum","Class","Order","Family","Genus","Species") 
+
+tax_col_norm_otu <- c("OTUID","Kingdom","Phylum","Class","Order","Family","Genus","Species")
+tax_col_extra_otu <- c("OTUID","None","Kingdom","Phylum","Class","Order","Family","Genus","Species")
+
+if (!is.null(options$outdir)) {
+  # Create the directory
+  dir.create(options$outdir,FALSE)
+}
+
+is.biom<-options$biom
+
+
+pdffile <- gsub("[ ]+", "", paste(options$outdir,"/pdffile.pdf"))
+pngfile_net <- gsub("[ ]+", "", paste(options$outdir,"/net.png"))
+htmlfile <- gsub("[ ]+", "", paste(options$htmlfile))
+
+if(is.biom=="set_biom"){
+
+galaxy_biom <- import_biom(options$infile)
+galaxy_map <- import_qiime_sample_data(options$metafile)
+
+number.of.tax.rank<-length(colnames(tax_table(galaxy_biom)))
+
+	if(number.of.tax.rank == 7){
+	    colnames(tax_table(galaxy_biom)) <- tax_col_norm
+	} else {
+	    colnames(tax_table(galaxy_biom)) <- tax_col_extra
+	}
+
+physeq_galaxy <- merge_phyloseq(galaxy_biom,galaxy_map)
+
+} else {
+
+count.table<-read.table(options$infile,header=T,sep="\t",comment.char="",stringsAsFactors = F)
+meta.table<-read.table(options$metafile,header=T,sep="\t",comment.char="",stringsAsFactors = F)
+tax.table<-read.table(options$obsfile,header=T,sep="\t",comment.char="",stringsAsFactors = F)
+
+
+colnames(count.table)<-gsub("^X","",colnames(count.table))
+colnames(meta.table)<-gsub("^X.","",colnames(meta.table))
+colnames(tax.table)<-gsub("^X.","",colnames(tax.table))
+
+count.table.formatted<-matrix.format(count.table)
+OTU<-otu_table(count.table.formatted,taxa_are_rows = TRUE)
+
+tax.table.new<-as.data.frame(cbind(tax.table[,1],t(as.data.table(strsplit(tax.table[,2],";")))))
+
+
+if(length(colnames(tax.table.new)) != length(tax_col_extra_otu)){
+colnames(tax.table.new)<-tax_col_norm_otu
+}else{
+colnames(tax.table.new)<-tax_col_extra_otu
+}
+
+tax.table.formatted<-matrix.format(tax.table.new)
+
+TAX<-tax_table(tax.table.formatted)
+
+physeq_galaxy <- phyloseq(OTU, TAX)
+
+
+galaxy_map<-meta.table
+
+rownames(galaxy_map)<-meta.table[,1]
+
+sampledata<-sample_data(as.data.frame(galaxy_map,row.names=sample_names(galaxy_map),stringsAsFactos=F))
+
+sample_data(physeq_galaxy)<-sampledata
+
+}
+
+
+x.selectedColumn<-colnames(galaxy_map)[options$xcolumn]
+l.selectedColumn<-colnames(galaxy_map)[options$lcolumn]
+
+### normalisation
+if(is.null(options$norm) || options$norm =="false"){
+suppressMessages(raw.count.deseq2.obj<-phyloseq_to_deseq2(physeq_galaxy,as.formula(paste('~',x.selectedColumn,sep=""))))
+geoMeans = apply(counts(raw.count.deseq2.obj), 1, gm_mean)
+
+deseq.obj = estimateSizeFactors(raw.count.deseq2.obj, geoMeans = geoMeans)
+deseq.obj.norm<-otu_table(as.matrix(counts(deseq.obj,normalized=T)),taxa_are_rows=TRUE)
+
+otu_table(physeq_galaxy)<-deseq.obj.norm
+}
+
+
+# Produce PDF file
+
+pdf(pdffile); 
+plot_net(physeq_galaxy,point_label=x.selectedColumn,color=l.selectedColumn)
+garbage<-dev.off();
+
+#Cairo(pngfile_net, type="png", bg="white",pointsize=12,dpi=100,units="in",width=6,height=6) 
+png(pngfile_net,units="in",width=6,height=6,pointsize=12,res=100,bg="white")
+plot_net(physeq_galaxy,point_label=x.selectedColumn,color=l.selectedColumn)
+garbage<-dev.off() 
+
+# Produce the HTML file
+htmlfile_handle <- file(htmlfile)
+html_output = c('<html><body>',
+	            '<table align="center">',
+		        '<tr>',
+		        '<td valign="middle" style="vertical-align:middle;">',
+                '<a href="pdffile.pdf"><img src="net.png"/></a>',
+		        '</td>',
+		        '</tr>',
+		        '</table>',
+                '</html></body>');
+writeLines(html_output, htmlfile_handle);
+close(htmlfile_handle);