4. Quality Control

Aim:

Assess the quality of raw datasets
Define quality trimming parameters prior running the complete RNAseq gene profiling pipeline

Run RNA-seq QC check

The pipeline requires preparing at least 2 files:

Metadata file (samplesheet.csv) that specifies the name of the samples, location of FASTQ files ('Read 1' and ‘Read 2’), and strandedness (forward, reverse, or auto. Note: auto is used when the strandedness of the data is unknown)
PBS Pro script (launch_nf-core_RNAseq_QC.pbs) with instructions to run the pipeline

Create the metadata file (samplesheet.csv):

Change to the data folder directory:

cd $HOME/workshop/2024-2/session4_RNAseq/data/

Copy the bash script to the working folder

cp /work/training/2024/rnaseq/scripts/create_samplesheet_nf-core_RNAseq.sh $HOME/workshop/2024/rnaseq/data/

Note: you could replace ‘$HOME/workshop/data’ with “.” A dot indicates ‘current directory’ and will copy the file to the directory where you are currently located

View the content of the script:

cat create_samplesheet_nf-core_RNAseq.sh

Example for Single-End data (only ‘Read 1’ is available):

#!/bin/bash -l

#User defined variables

##########################################################

DIR='$HOME/workshop/data/'

INDEX='samplesheet.csv'

##########################################################

#load python module

module load python/3.10.8-gcccore-12.2.0

#fetch the script to create the sample metadata table

wget -L https://raw.githubusercontent.com/nf-core/rnaseq/master/bin/fastq_dir_to_samplesheet.py

chmod +x fastq_dir_to_samplesheet.py

#generate initial sample metadata file

./fastq_dir_to_samplesheet.py $DIR $INDEX \

--strandedness auto \

--read1_extension .fastq.gz

NOTE: when paired end data is available (R1 and R2 FASTQ files for each sample). Then use:

--read1_extension R1.fastq.gz \

--read2_extension R2.fastq.gz

NOTE: check the suffix of the FASTQ files to specify the extension (e.g., R1_001.fastq.gz or R1.fq.gz)

Let’s generate the metadata file by running the following command:

sh create_RNAseq_samplesheet.sh

Check the newly created samplesheet.csv file:

ls -l
cat samplesheet.cvs

sample,fastq_1,fastq_2,strandedness

SRR20622172,/work/training/rnaseq/data/SRR20622172.fastq.gz,,auto

SRR20622173,/work/training/rnaseq/data/SRR20622173.fastq.gz,,auto

SRR20622174,/work/training/rnaseq/data/SRR20622174.fastq.gz,,auto

SRR20622175,/work/training/rnaseq/data/SRR20622175.fastq.gz,,auto

SRR20622176,/work/training/rnaseq/data/SRR20622176.fastq.gz,,auto

SRR20622177,/work/training/rnaseq/data/SRR20622177.fastq.gz,,auto

SRR20622178,/work/training/rnaseq/data/SRR20622178.fastq.gz,,auto

SRR20622179,/work/training/rnaseq/data/SRR20622179.fastq.gz,,auto

SRR20622180,/work/training/rnaseq/data/SRR20622180.fastq.gz,,auto

Copy the PBS Pro script for QC (launch_nf-core_RNAseq_QC.pbs)

Copy and paste the code below to the terminal:

cp $HOME/workshop/2024/rnaseq/data/samplesheet.csv $HOME/workshop/2024/rnaseq/runs/run2_QC
cp $HOME/workshop/scripts/launch_nf-core_RNAseq_QC.pbs $HOME/workshop/2024/rnaseq/runs/run2_QC
cd $HOME/workshop/2024/rnaseq/runs/run2_QC

Line 1: Copy the samplesheet.csv file to the working directory
Line 2: move to the working directory
Line 3: copy the launch_nf-core_RNAseq_QC.pbs submission script to the working directory

View the content of the launch_nf-core_RNAseq_QC.pbs script:

cat launch_nf-core_RNAseq_QC.pbs

#!/bin/bash -l

#PBS -N nfrnaseq_QC

#PBS -l select=1:ncpus=2:mem=4gb

#PBS -l walltime=24:00:00

#work on current directory

cd $PBS_O_WORKDIR

#load java and set up memory settings to run nextflow

module load java

export NXF_OPTS='-Xms1g -Xmx4g'

nextflow run nf-core/rnaseq \

-profile singularity \

-r 3.14.0 \

--input samplesheet.csv \

--outdir results \

--genome GRCm38-local \

--skip_trimming \

--skip_alignment \

--skip_pseudo_alignment

We recommend running the nextflow nf-core/rnaseq pipeline once and then assessing the fastqc results folder to assess if sequence biases are present in the 5'-end and 3'-end ends of the sequences.
Version 3.12.0 allows running the pipeline to do quality assessment only, without any alignment, read counting, or trimming.

Submitting the job

Once you have created the folder for the run, the samplesheet.csv file, and launch.pbs, you are ready to submit the job to the HPC scheduler:

qsub launch_nf-core_RNAseq_QC.pbs

Monitoring the Run

qjobs

to check on the jobs, you are running. Nextflow will launch additional jobs during the run.

You can also check the .nextflow.log file for details on what is going on.

Once the pipeline has finished running - Assess the QC report:

NOTE: To proceed, you need to be on QUT’s WiFi network or signed via VPN.

To browse the working folder in the HPC type in the file finder:

Windows PC

\\hpc-fs\work\training\rnaseq

Mac

smb://hpc-fs/work/training/rnaseq

Evaluate the nucleotide distributions in the 5'-end and 3'-end of the sequenced reads (Read1 and Read2). Look into the “MultiQC” folder and open the provided HTML report.

Items to check:

The overall quality of the experiment and reads. Look at the “Sequence Quality Histogram” plot. For example, if Phred assigns a quality score of 30 to a base, the chances that this base is called incorrectly are 1 in 1000. Phred quality scores are logarithmically linked to error probabilities.

Phred Quality Score	Probability of incorrect base call	Base call accuracy
10	1 in 10	90%
20	1 in 100	99%
30	1 in 1000	99.9%
40	1 in 10,000	99.99%
50	1 in 100,000	99.999%
60	1 in 1,000,000	99.9999%

Assess QC reports (FastQC and MultiQC) to define how many nucleotides should be trimmed from the 5'-end and/or 3-end regions of the FASTQ reads (see Case 3 below).