Skip to content

How to start

Jump to bottom
Daniel Fischer edited this page Aug 4, 2022 · 11 revisions

Download the pipeline

The easiest way to download the pipeline is to clone the GitHub repository like this:

ADD CLONE COMMAND

Requirements

In the current configuration, the pipeline is only able to run on systems that support singularity. Other systems (inlcuding locally installed tools) are under development.

Preparations

Once the pipeline is cloned to the local system, there are a few configuration files that need to be prepared. The repository contains default/standard files that need to be adjusted for the concrete user case. Further, a sample sheet file needs to be created.

Sample sheet

Creating the sample sheet could be laborious, considering you have lots of files. If file names are systematically given (e.g. like this: sampleName_S12_L001_R1_001.fastq.gz), you could create a samplesheet like this:

ls | cut -d'_' -f1-3 | uniq > rawsample
cat rawsample | cut -d'_' -f1-2 > sample_name
cat rawsample | cut -d'_' -f3 | > lane
sed -e 's/$/_R1_001.fastq.gz/' rawsample > read1
sed -e 's/$/_R2_001.fastq.gz/' rawsample > read2
paste -d'\t' rawsample sample_name lane read1 read2 > samplesheet.tsv
sed -i '1s/^/rawsample\tsample_name\tlane\tread1\tread2\n/' samplesheet.tsv
rm rawsample sample_name lane read1 read2

Further, the pipeline currently only works for paired-end data. The support for single-end data will be added 'upon personal demand'...

Pipeline targets

Like every Snakamake pipeline, the default target rule is such that all modules of the pipeline will be executed. However, if needed, also the individual modules can be executed, as described in the following.

Preparations

To run the preparations module, just add the corresponding keyword to the pipeline command:

 bash run_WGS-Pipeline.sh preparations -np

In that case, a couple of rules are execute:

  1. Download of existing variants In case an url to a vcf with known variants is provided in the configuration file, this file is downloaded and registered to the pipeline

  2. Download the reference genome In case an url to a reference genome is provided in the configuration file, this file is downloaded and registered to the pipeline

  3. Create the bwa index for the reference genome The bwa inded for the reference genome is created

  4. Create the fasta index for the reference genome The fai is created for the reference genome

Trimming

The trimming module for the raw data can be initialised with that command

bash run_WGS-Pipeline.sh trimming

The trimming parameters can be adjusted in the configuration file under params and trimmomatic:

Here, the default values are:

threads: 20
illuminaclip: "/Trimmomatic/adapters/TruSeq3-PE-2.fa:2:30:3:1:true"
leading: 20 
trailing: 20  
slidingwindow: "3:15"  
avgqual: 20
minlen: 35

Here, the module takes the raw reads as they are provided in the sample sheet file and creates then the following output:

R1="%s/FASTQ/TRIMMED/{rawsamples}_R1.fastq.gz" % (config["project-folder"]),
R1S="%s/FASTQ/TRIMMED/{rawsamples}_R1-singleton.fastq.gz" % (config["project-folder"]),
R2="%s/FASTQ/TRIMMED/{rawsamples}_R2.fastq.gz" % (config["project-folder"]),
R2S="%s/FASTQ/TRIMMED/{rawsamples}_R2-singleton.fastq.gz" % (config["project-folder"]),
S="%s/FASTQ/TRIMMED/{rawsamples}.summary" % (config["project-folder"])

Here, R1 is the trimmed read 1 and R2 the corresponding paired-end read. The files with ending S (R1S and R2S) are the corresponding singletons, meaning the reads that do not have a proper pair anymore after trimming. And finally, the file labeled S contains the summary file from trimmomatic.

Quality checks

The quality check module can be executed via

bash run_WGS-Pipeline.sh qc

Here, it will run FastQC for all read files (including the trimmed ones, so the qc-module will also trigger the trimming module) as well as MultiQC.

Alignment

Here, all the lane-wise files will be aligned and then merged. THe command is

bash run_WGS-Pipeline.sh qc

In detail, the follwoing steps will be processed

  1. Align the data (bwa). Here, the different lanes are aligned to the reference genome

  2. Sort and index the alignments (samtools). The lane-wise alignments are sorted and indexed

  3. Mark the duplicates in the BAM files Duplicates are marked in the BAM files

  4. Merge the lane-wise BAM-files into sample-wise BAMs These are then the final BAM files that are used for variant calling.

Variant calling

THe module to call variants consists also out of several rules. First, the module is started like this:

bash run_WGS-Pipeline.sh variant_calling

and it triggers a few steps:

  1. Perform the base recalibration (PICARD)
  2. Analyzeis of the covariates (PICARD)
  3. Running the haplotype caller (GATK)
  4. Calculate Depth of Coverage (GATK)
  5. Combine gvcfs (GATK)
  6. Apply GenotypeGVCFs (GATK)
  7. Train Variant Recalibrator
  8. ApplyVQSR, filtering step

Prepare the output

To prepare the output, you can run the module

bash run_WGS-Pipeline.sh prepare_output
```

Clone this wiki locally