Update Long_Read_Assembly

Long_Read_Assembly/Dockerfile

+FROM rocker/r-ver:3.5.3
+
+ENV PATH /opt/biotools/bin:$PATH
+ENV ROOTSYS /opt/biotools/root
+ENV LD_LIBRARY_PATH '$LD_LIBRARY_PATH:$ROOTSYS/lib'
+
+RUN apt-get update
+RUN apt-get install -yq tzdata
+RUN apt-get install -y locales
+RUN locale-gen "en_US.UTF-8"
+RUN export LC_ALL=en_US.UTF-8
+RUN export LANG=en_US.UTF-8
+
+RUN apt-get install -y curl wget apt-utils
+RUN apt-get install -y gcc fort77 aptitude
+RUN aptitude install -y g++ xorg-dev libreadline-dev  gfortran
+RUN apt-get install -y libssl-dev libxml2-dev libpcre3-dev liblzma-dev libbz2-dev libcurl4-openssl-dev liblapack3 git nano graphviz python3 python3-pip
+
+RUN apt-get install -y  autotools-dev automake cmake grep sed dpkg fuse zip build-essential pkg-config bzip2 ca-certificates libglib2.0-0 libxext6 libsm6 libxrender1 mercurial subversion zlib1g-dev libncurses5-dev libncursesw5-dev
+
+RUN apt-get clean
+
+RUN if [ ! -d "/opt/biotools" ];then mkdir /opt/biotools; fi
+RUN if [ ! -d "/opt/biotools/bin" ];then mkdir /opt/biotools/bin; fi
+RUN chmod 777 -R /opt/biotools/
+ENV PATH /opt/biotools/bin:$PATH
+
+RUN pip3 install snakemake==5.4.0 oyaml
+RUN pip3 install multiqc==1.7
+RUN fic=$(find /usr/ -name stacks.py) && sed -i 's/out_dict\[s_name\] = cdict/out_dict\[content\[0\]\] = cdict/' $fic
+
+RUN Rscript -e 'install.packages("yaml",Ncpus=8,repos="https://cloud.r-project.org/")'
+RUN Rscript -e 'install.packages("DT",Ncpus=8,repos="https://cloud.r-project.org/")'
+RUN Rscript -e 'install.packages("shiny",Ncpus=8,repos="https://cloud.r-project.org/")'
+RUN Rscript -e 'install.packages("shinydashboard",Ncpus=8,repos="https://cloud.r-project.org/")'
+RUN Rscript -e 'install.packages("shinyjs",Ncpus=8,repos="https://cloud.r-project.org/")'
+RUN Rscript -e 'install.packages("shinyFiles",Ncpus=8,repos="https://cloud.r-project.org/")'
+
+
+COPY files /workflow
+COPY sagApp /sagApp
+
+ENV PATH /opt/biotools/minimap2-2.17_x64-linux:$PATH 
+RUN cd /opt/biotools \
+ && wget https://github.com/lh3/minimap2/releases/download/v2.17/minimap2-2.17_x64-linux.tar.bz2 \
+ && tar -xvjf minimap2-2.17_x64-linux.tar.bz2 \
+ && rm minimap2-2.17_x64-linux.tar.bz2
+
+ENV PATH /opt/biotools/miniasm-0.3:$PATH 
+RUN cd /opt/biotools \
+ && wget -O miniasm.tar.gz https://github.com/lh3/miniasm/archive/v0.3.tar.gz \
+ && tar -xvzf miniasm.tar.gz \
+ && cd miniasm-0.3 \
+ && make -j 10 \
+ && cd .. rm miniasm.tar.gz
+
+ENV PATH /opt/biotools/racon-v1.4.10/bin:$PATH 
+RUN cd /opt/biotools \
+ && wget -O racon.tar.gz https://github.com/lbcb-sci/racon/releases/download/1.4.10/racon-v1.4.10.tar.gz \
+ && tar -xvzf racon.tar.gz \
+ && cd racon-v1.4.10 \
+ && cmake -DCMAKE_BUILD_TYPE=Release \
+ && make -j 10 \
+ && cd .. \
+ && rm racon.tar.gz
+
+RUN wget -O htslib.tar.bz2 https://github.com/samtools/htslib/releases/download/1.9/htslib-1.9.tar.bz2 \
+ && tar -xvjf htslib.tar.bz2 \
+ && cd htslib-1.9 \
+ && ./configure \
+ && make -j 10 \
+ && make install
+
+RUN cd /opt/biotools \
+ && wget https://github.com/samtools/samtools/releases/download/1.9/samtools-1.9.tar.bz2 \
+ && tar -xvjf samtools-1.9.tar.bz2 \
+ && cd samtools-1.9 \
+ && ./configure && make -j 10 \
+ && cd .. \
+ && mv samtools-1.9/samtools bin/samtools \
+ && rm -r samtools-1.9 samtools-1.9.tar.bz2
+
+RUN apt-get -y install tabix
+
+RUN pip3 install medaka==0.10.1
+
+RUN apt-get install -y zlib1g-dev pkg-config libfreetype6-dev libpng-dev python-matplotlib python-setuptools
+
+RUN cd /opt/biotools \
+ && wget https://github.com/ablab/quast/releases/download/quast_5.0.2/quast-5.0.2.tar.gz \
+ && tar -zxvf quast-5.0.2.tar.gz \
+ && cd quast-5.0.2/ \
+ && python3 ./setup.py install_full
+
+EXPOSE 3838
+CMD ["Rscript", "-e", "setwd('/sagApp/'); shiny::runApp('/sagApp/app.R',port=3838 , host='0.0.0.0')"]

Long_Read_Assembly/files/Snakefile

@@ -57,6 +57,11 @@ def medaka_inputs():
 	inputs["assembly_fasta"] = rules.racon.output.assembly_corrected
 	return inputs
 
+def quast_inputs():
+	inputs = dict()
+	inputs["assembly"] = rules.medaka.output.consensus_assembly
+	return inputs
+
 def prepare_report_inputs():
     inputs = list()
     for step in STEPS:
@@ -106,6 +111,9 @@ def step_outputs(step):
 
 	if (step == "polishing"):
 		outputs = rules.medaka.output
+
+	if (step == "assembly_quality"):
+		outputs = rules.quast.output
 		
 	if (step == "all"):
 		outputs = list(rules.multiqc.output)
@@ -197,19 +205,20 @@ rule racon:
 		config["racon_threads"]
 	shell: 
 		"racon "
+		"-u "
 		"-t {threads} "
 		"{input.reads} "
 		"{input.overlaps} "
 		"{input.assembly} "
 		"> {output.assembly_corrected} "
-		"|& tee {log}"
+		"2> >(tee {log} >&2)"
 
 
 rule medaka:
 	input:
 		**medaka_inputs(),
 	output: 
-		consensus_assembly = config["results_dir"] + "/" + config["medaka_output_dir"] + "/consensus_assembly.fasta",
+		consensus_assembly = config["results_dir"] + "/" + config["medaka_output_dir"] + "/consensus.fasta",
 	params: 
 		output_dir = config["results_dir"] + "/" + config["medaka_output_dir"]+ "/",
 		model = config["medaka_model"],
@@ -227,6 +236,24 @@ rule medaka:
 		"|& tee {log}"
 
 
+rule quast:
+	input:
+		**quast_inputs(),
+	output: 
+		report = config["results_dir"] + "/" + config["quast_output_dir"] + "/report.tsv",
+	params: 
+		output_dir = config["results_dir"] + "/" + config["quast_output_dir"]+ "/",
+	log: 
+		config["results_dir"] + "/logs/quast/quast_log.txt"
+	threads: 
+		config["quast_threads"]
+	shell: 
+		"quast.py "
+		"{input.assembly} "
+		"-o {params.output_dir} "
+		"--threads {threads}"
+
+
 
 rule prepare_report:
     input:

Long_Read_Assembly/files/params.total.yml

@@ -19,6 +19,9 @@ params:
   medaka_output_dir: medaka
   medaka_threads: 4
   medaka_model: r941_min_high
+  assembly_quality: quast
+  quast_output_dir: quast
+  quast_threads: 4
 samples: []
 groups: []
 steps:
@@ -47,6 +50,11 @@ steps:
   tools:
   - medaka
   default: medaka
+- title: Assembly Quality
+  name: assembly_quality
+  tools:
+  - quast
+  default: quast
 params_info:
   results_dir:
     type: output_dir
@@ -84,6 +92,11 @@ params_info:
     rule: medaka
     type: select
     label: Medaka model
+  quast_threads:
+    tool: quast
+    rule: quast
+    type: numeric
+    label: Threads to use
 prepare_report_scripts: []
 prepare_report_outputs: {}
 outputs:
@@ -113,11 +126,17 @@ outputs:
   medaka:
     medaka:
     - name: consensus_assembly
-      file: consensus_assembly.fasta
+      file: consensus.fasta
       description: Consensus assembly file
+  quast:
+    quast:
+    - name: report
+      file: report.tsv
+      description: Report (tabulated file)
 multiqc:
   minimap2_overlap_self: custom
   miniasm: custom
   minimap2_reference: custom
   racon: custom
   medaka: custom
+  quast: quast

Long_Read_Assembly/files/singularity.def

@@ -101,6 +101,8 @@ singularity run --app appName this_container.sif
 	cd racon-v1.4.10
 	cmake -DCMAKE_BUILD_TYPE=Release
 	make -j 10
+	cd ..
+	rm racon.tar.gz
 	echo 'export PATH="/opt/biotools/racon-v1.4.10/bin:$PATH"' >>$SINGULARITY_ENVIRONMENT
 
 	wget -O htslib.tar.bz2 https://github.com/samtools/htslib/releases/download/1.9/htslib-1.9.tar.bz2
@@ -123,3 +125,11 @@ singularity run --app appName this_container.sif
 
 	pip3 install medaka==0.10.1
 
+	apt-get install -y zlib1g-dev pkg-config libfreetype6-dev libpng-dev python-matplotlib python-setuptools
+
+	cd /opt/biotools
+	wget https://github.com/ablab/quast/releases/download/quast_5.0.2/quast-5.0.2.tar.gz
+	tar -zxvf quast-5.0.2.tar.gz
+	cd quast-5.0.2/
+	python3 ./setup.py install_full
+

Long_Read_Assembly/sagApp/R/menugauche.R

@@ -12,12 +12,14 @@ MenuGauche = sidebarMenu(id="sidebarmenu",
 
 	menuItem("Polishing", tabName="polishing", icon=icon("pencil", lib="font-awesome"), newtab=FALSE),
 
+	menuItem("Assembly Quality", tabName="assembly_quality", icon=icon("pencil", lib="font-awesome"), newtab=FALSE),
+
 	menuItem("Rule Graph", tabName="RULEGRAPH", icon=icon("gear", lib="font-awesome"), newtab=FALSE),
 
 	tags$br(),
 
 	numericInput("cores", label = "Threads available", min = 1, max = 24, step = 1, width =  "auto", value = 16),
-selectInput("force_from", label = "Start again from a step : ", selected = "none", choices = list('none'='none','Find overlaps'='find_overlaps','Assembly'='assembly','Mapping'='mapping','Correction'='correction','Polishing'='polishing',"All"="all")),	tags$br(),
+selectInput("force_from", label = "Start again from a step : ", selected = "none", choices = list('none'='none','Find overlaps'='find_overlaps','Assembly'='assembly','Mapping'='mapping','Correction'='correction','Polishing'='polishing','Assembly Quality'='assembly_quality',"All"="all")),	tags$br(),
 	actionButton("RunPipeline", "Run pipeline",  icon("play"), class="btn btn-info"),
 
 	actionButton("StopPipeline", "Stop pipeline",  icon("stop"), class="btn btn-secondary"),

Long_Read_Assembly/sagApp/app.R

@@ -20,6 +20,7 @@ source("./pages/pages_def_assembly.R", local=T)
 source("./pages/pages_def_mapping.R", local=T)
 source("./pages/pages_def_correction.R", local=T)
 source("./pages/pages_def_polishing.R", local=T)
+source("./pages/pages_def_assembly_quality.R", local=T)
 tabRULEGRAPH = fluidPage(box(title = 'Rule Graph :', width = 12, status = 'primary', collapsible = TRUE, solidHeader = TRUE, uiOutput('RULEGRAPH_svg'),actionButton('refresh_rg', 'Refresh',  icon('sync'), class='btn btn-info')))
 tabReport = fluidPage(box(title = 'Report :', width = 12, status = 'primary', collapsible = TRUE, solidHeader = TRUE, uiOutput('report_html')))
 tabRUN = fluidPage(box(title = 'Run :', width = 12 , status = 'primary', collapsible = TRUE, solidHeader = TRUE, uiOutput('run_out',style = 'overflow-y: scroll; height: 600px')),actionButton("unlock", "Unlock the directory in case of previous failure"))
@@ -57,7 +58,9 @@ tags$head(style),
 
 	tabItem(tabName = "correction",         tabcorrection),
 
-	tabItem(tabName = "polishing",         tabpolishing)
+	tabItem(tabName = "polishing",         tabpolishing),
+
+	tabItem(tabName = "assembly_quality",         tabassembly_quality)
 
 	,tabItem(tabName = "RULEGRAPH",         tabRULEGRAPH)
 	,tabItem(tabName = "Report",         tabReport)

Long_Read_Assembly/sagApp/pages/pages_def_assembly_quality.R

+tabassembly_quality = fluidPage(
+
+box(title = "Parameters :", width = 12, status = "primary", collapsible = TRUE, solidHeader = TRUE,
+
+	hidden(textInput("selectassembly_quality", label = "", value="quast")),box(title = "QUAST", width = 12, status = "success", collapsible = TRUE, solidHeader = TRUE,
+		numericInput("quast_threads", label = "Threads to use", min = 1, max = 32, step = 1, width =  "auto", value = 4),
+
+		p("QUAST: Quality Assessment Tool for Genome Assemblies"),
+
+		p("Website : ",a(href="http://bioinf.spbau.ru/quast","http://bioinf.spbau.ru/quast",target="_blank")),
+
+		p("Documentation : ",a(href="http://quast.bioinf.spbau.ru/manual.html","http://quast.bioinf.spbau.ru/manual.html",target="_blank")),
+
+		p("Paper : ",a(href="https://doi.org/10.1093/bioinformatics/btt086","https://doi.org/10.1093/bioinformatics/btt086",target="_blank"))
+
+	)))
+
+

Long_Read_Assembly/sagApp/server/opt_global.R

@@ -71,6 +71,14 @@ save_params <- function(path_param){
 			res = paste0(res, paste("medaka_model:", paste0('"', input$medaka_model, '"'), "\n", sep = " "))
 		}	
 
+	# 	Page : assembly_quality
+		res = paste0(res , paste("assembly_quality:", paste0('"', input$selectassembly_quality, '"'), "\n", sep = " "))
+		if(!is.na(as.numeric(input$quast_threads))) {
+			res = paste0(res, paste("quast_threads:", input$quast_threads, "\n", sep = " "))
+		} else {
+			res = paste0(res, paste("quast_threads:", paste0('"', input$quast_threads, '"'), "\n", sep = " "))
+		}	
+
 	a = yaml.load_file("/workflow/params.total.yml")
 	p = a[["params"]]
 	a["params"] = NULL