QC reports for the pollen seq data

# , include=TRUE, warning=FALSE, echo=TRUE, error=FALSE
knitr::opts_knit$set(root.dir=normalizePath('../../'))
knitr::opts_chunk$set(warning=FALSE, message=FALSE, error=FALSE, echo=TRUE)
library(fastqcr)
library(tidyverse)

Load QC report Data

qcdir1 <- "largedata/fqqc/BM-MB-1st/"
qcdir2 <- "largedata/fqqc/BM-MB-2nd/"

Aggregating Reports

QC reports for the first round of pollen sequencing

They will report 10 modules and 7 columns for each fq file. The column names are: - sample: sample names
- module: fastqc modules
- status: fastqc module status for each sample
- tot.seq: total sequences (i.e.: the number of reads)
- seq.length: sequence length
- pct.gc: percentage of GC content
- pct.dup: percentage of duplicate reads

QC reports for the 2nd round of pollen sequencing

Binding data together and get WARN and FAIL modules:

Column names:

• module: fastqc modules
  
• nb_samples: the number of samples tested
  
• nb_pass, nb_fail, nb_warn: the number of samples that passed, failed and warned, respectively.
  
• failed, warned: the name of samples that failed and warned, respectively.

qc <- rbind(qc1, qc2)
qc %>%
    select(sample, module, status) %>%    
    filter(status %in% c("WARN", "FAIL")) %>%
    arrange(sample)

## Summary of qc
summary(qc)

General Statistics

Column names:

• pct.dup: the percentage of duplicate reads,
  
• pct.gc: the percentage of GC content,
  
• tot.seq: total sequences or the number of reads and
  
• seq.length: sequence length or the length of reads.

mean(as.numeric(st$pct.dup))

[1] 45.9713

Note, we should have 24 files, but BW350-1_L5_I20006.R2.clean.fastq.gz failed to unzip.

Inspecting Problems

Per module problems:

qc_fails(qc, "module")

#qc_warns(qc, "module")
#qc_problems(qc, "module")

Per sample problems

#qc_warns(qc, "sample")
#qc_problems(qc, "sample")
qc_fails(qc, "sample")

Building an HTML Report

qc_report(qc.path=qcdir1, result.file = "largedata/fqqc",
          experiment = "Pollen BM")


qc <- qc_read(list.files(qcdir1, full.names = T)[3])
qc_plot(qc, "Summary")
qc_plot(qc, "Basic Statistics")
qc_plot(qc, "Kmer Content")

qc_plot(qc, "per_tile_sequence_quality")

qc_plot(qc, "Per sequence quality scores")

qc_plot(qc, "Per base sequence content")

qc_plot(qc, "Per base N content")


# Demo file
qc.file <- system.file("fastqc_results", "S1_fastqc.zip",  package = "fastqcr")

# Read all modules
qc <- qc_read(qc.file)
names(qc)

# Plot per sequence GC content
qc_plot(qc, "Per sequence GC content")

# Per base sequence quality
qc_plot(qc, "per_base sequence quality")

# Per sequence quality scores
qc_plot(qc, "Per sequence quality scores")

# Per base sequence content
qc_plot(qc, "Per base sequence content")

# Sequence duplication levels
qc_plot(qc, "Sequence duplication levels")