systemPipeR 1.24.6
Users want to provide here background information about the design of their VAR-Seq project.
This report describes the analysis of a VAR-Seq project studying the genetic differences among several strains … from organism ….
Typically, users want to specify here all information relevant for the analysis of their NGS study. This includes detailed descriptions of FASTQ files, experimental design, reference genome, gene annotations, etc.
systemPipeRdata package is a helper package to generate a fully populated systemPipeR workflow environment in the current working directory with a single command. All the instruction for generating the workflow template are provide in the systemPipeRdata vignette here.
After building and loading the workflow environment generated by genWorkenvir
from systemPipeRdata all data inputs are stored in
a data/
directory and all analysis results will be written to a separate
results/
directory, while the systemPipeVARseq.Rmd
script and the targets
file are expected to be located in the parent directory. The R session is expected to run from this parent
directory. Additional parameter files are stored under param/
.
To work with real data, users want to organize their own data similarly
and substitute all test data for their own data. To rerun an established
workflow on new data, the initial targets
file along with the corresponding
FASTQ files are usually the only inputs the user needs to provide.
For more details, please consult the documentation
here. More information about the targets
files from systemPipeR can be found here.
Now open the R markdown script systemPipeVARseq.Rmd
in your R IDE (e.g. vim-r or RStudio) and
run the workflow as outlined below.
Here pair-end workflow example is provided. Please refer to the main vignette
systemPipeR.Rmd
for running the workflow with single-end data.
The systemPipeR
package needs to be loaded to perform the analysis steps shown in
this report (H Backman and Girke 2016).
library(systemPipeR)
If you are running on a single machine, use following code as an example to check if some tools used in this workflow are in your environment PATH. No warning message should be shown if all tools are installed.
targets
fileThe targets
file defines all FASTQ files and sample comparisons of the analysis workflow.
targetspath <- system.file("extdata", "targetsPE.txt", package = "systemPipeR")
targets <- read.delim(targetspath, comment.char = "#")
targets[1:4, 1:4]
## FileName1 FileName2
## 1 ./data/SRR446027_1.fastq.gz ./data/SRR446027_2.fastq.gz
## 2 ./data/SRR446028_1.fastq.gz ./data/SRR446028_2.fastq.gz
## 3 ./data/SRR446029_1.fastq.gz ./data/SRR446029_2.fastq.gz
## 4 ./data/SRR446030_1.fastq.gz ./data/SRR446030_2.fastq.gz
## SampleName Factor
## 1 M1A M1
## 2 M1B M1
## 3 A1A A1
## 4 A1B A1
The following removes reads with low quality base calls (here a certain pattern) from all FASTQ files.
targetsPE <- system.file("extdata", "targetsPE.txt", package = "systemPipeR")
dir_path <- system.file("extdata/cwl/preprocessReads/trim-pe",
package = "systemPipeR")
trim <- loadWorkflow(targets = targetsPE, wf_file = "trim-pe.cwl",
input_file = "trim-pe.yml", dir_path = dir_path)
trim <- renderWF(trim, inputvars = c(FileName1 = "_FASTQ_PATH1_",
FileName2 = "_FASTQ_PATH2_", SampleName = "_SampleName_"))
trim
output(trim)[1:2]
preprocessReads(args = trim, Fct = "trimLRPatterns(Rpattern='GCCCGGGTAA', subject=fq)",
batchsize = 1e+05, overwrite = TRUE, compress = TRUE)
writeTargetsout(x = trim, file = "targets_trimPE.txt", step = 1,
new_col = c("FileName1", "FileName2"), new_col_output_index = c(1,
2), overwrite = TRUE)
The following seeFastq
and seeFastqPlot
functions generate and plot a series of
useful quality statistics for a set of FASTQ files including per cycle quality box
plots, base proportions, base-level quality trends, relative k-mer
diversity, length and occurrence distribution of reads, number of reads
above quality cutoffs and mean quality distribution. The results are
written to a PDF file named fastqReport.pdf
. Use the output from previous step
(fastq trimming) as the demonstration here to generate fastq report.
fqlist <- seeFastq(fastq = infile1(trim), batchsize = 1e+05,
klength = 8)
pdf("./results/fastqReport.pdf", height = 18, width = 4 * length(fqlist))
seeFastqPlot(fqlist)
dev.off()