<center><img src="https://i.imgur.com/rPIZUIq.png" alt="drawing" width="700"/></center>
# ACEIDHA: Quality control short-read reads and assembly
In this hands-on exercise, you will work on sequence data using the Galaxy platform. Galaxy EU is an open, web-based platform for accessible, reproducible, and transparent computational biomedical research, and usegalaxy.eu is a European Galaxy server for life science data.
> For this course, use **https://usegalaxy.eu/join-training/aceidha_bioinf22** link. This gives a dedicated computer powers for our workshop.
> After the course, your data and the server is accessible at **https://usegalaxy.eu. **
We will first analyse some of the genomes of the bacterium *Campylobacter jejuni* presented in this paper: https://pubmed.ncbi.nlm.nih.gov/33489938/. The genomes were generated using paired end (PE) *short-read technology* on the Illumina HiSeq or MiSeq platform (see [video](https://youtu.be/fCd6B5HRaZ8)). Our goal is to check how closely related the isolates are; we dont know which isolates we picked....
For each isolate, two `fastq` files were produced; one forward and one reverse (R1, R2). These files should already be available to you in your `History` as uploading these files were part of the pre-work needed to be done for this excerise.
First, we will inspect and correct the quality of the reads using the software **Fastp** (https://github.com/OpenGene/fastp) and compile a report for all reads using **MultiQC**. The reads will thereafter be assembled into bacterial genomes using a assembly pipeline called **Shovill** (https://github.com/tseemann/shovill).
All needed software (Fastp; MultiQC; Shovill) is available on the Galaxy platform, so no need for software download!
You can access the **usegalaxy.eu** at [https://usegalaxy.eu/join-training/aceidha_bioinf22]
</div>
</html>
<br/>
:::info
:information_source: **Note:** Names of files, folders, commands, and other functionality elements will be highlighted through the exercises in the following way: `file_name_example.txt`. User specific (and thus variable) character strings as for example user_name will be written between `<` and `>` as `<user_name>` to indicate the meaning. There are also `TASKS` written throughout the document; think about these and provide answers as you work through the HackMD.
:::
Suggested reading/useful pages for Sequence technology and Quality Control:
[**Next-generation sequencing** technologies: An overview *by Hu et al 2020* ](https://www.sciencedirect.com/science/article/pii/S0198885921000628?casa_token=EyzT_nHHEK4AAAAA:4bQYEpNGgYoCk49lZnlZ9fz91WXQZ4jmsK9_SAQZaRUTKf_x3GcQRYl6P0rxESMJmaXc7c9E7A#s0010)
[**Glossary** of commonly used terms and catch-phrases in the field of bacterial whole genome sequencing](https://www.clinicalmicrobiologyandinfection.com/action/showFullTableImage?isHtml=true&tableId=tbl1&pii=S1198743X18302064)
[**FastP**: an ultra-fast all-in-one FASTQ preprocessor. Chen et al. 2018](https://academic.oup.com/bioinformatics/article/34/17/i884/5093234?login=true)
[**MultiQC**: summarize analysis results for multiple tools and samples in a single report. by Ewels et al (2016)](https://doi.org/10.1093/bioinformatics/btw354)
[**Shovill** by Torsten Seeman](https://github.com/tseemann/shovill)
Suggested videos:
[Illumina MiSeq Workflow](https://www.youtube.com/watch?v=JA6mofeuntk)(short video)
[Illumina sequencing ](https://www.youtube.com/watch?v=6jf_6STEnI4)(rather long)
[MultiQC interpretation ](https://www.youtube.com/watch?v=qPbIlO_KWN0)
[Theory of assembly from MIT](https://www.youtube.com/watch?v=ZYW2AeDE6wU)
## Import data
<html>
<div style="background-color: #f2f3f4">
</div>
</html>
This can be done by many different ways.
**One simple way is to `Import History` See preassignment
If you have done that, your good to go - jump to point: **III**
If you have data on your own computer, you can follow these instructions:
**I.** In the **Tools** menu, press `Upload Data` [1]. This should redirect you to **your own computer**
<center> Figure 1 <img src="https://i.imgur.com/7tG7aEf.png" alt="drawing" width="800"/></center>
- [x] **II.** Browse your computer by clicking (1)(Fig 2)for Campylobacter reads downloaded from Google Drive. You should have received a link per e-mail. Upload all files with the `.fastq.gz` (Fig 3) extension; these are the fastq files. Set `Type` (set all) to `fastqsanger.gz`. This will then start uploading the reads to Galaxy server (2). It should take around 60 minutes, so hopefully you have done this already.
<center> Figure 2<img src="https://i.imgur.com/veLPrDC.png" alt="drawing" width="800"/></center>
<center> Figure 3<img src="https://i.imgur.com/fEz3TUN.png" alt="drawing" width="800"/></center>
**III.** When done, the imported data should appear in your current `History`. Name it something sensible, like `<Campylobacter_fastQC>`
<br/>
**IV.** Inspect the datatype (format). Preview the imported datasets by clicking on the file names [1] and by pressing Edit [2].
*TASK:
If you need to change the format, go to `Edit dataset` (the pen) and change the dataformat to `fastqssanger.gz` by pressing `Detect datatype` under `Datatypes`.*
<br/>
<center><img src="https://i.imgur.com/Tf5S6lb.png" alt="drawing" width="800"/></center>Inspect the fastq-file by pressing `'View'`[3]. Do so and inspect one selected fastq-file.
*TASK:
Locate a single entry in the fastq-file. Locate the actual genome sequence. What does the three other lines represent?*
## Quality control of short-read Illumina data
**V.** Now we will check the quality and trim the reads. There are a number of software available for this, but here we will utilize **Fastp** (https://academic.oup.com/bioinformatics/article/34/17/i884/5093234?login=true))
<html>
<div style="background-color: #f2f3f4">
**ABOUT**: Fastp is a quality control tool for high throughput sequence data, and provides comprehensive quality profiling for both before and after filtering data (quality curves, base contents, KMER, Q20/Q30, GC Ratio, duplication, adapter contents...). It has a range of features, such as filter out bad reads (too low quality, too short, or too many N...) and/or
cut low quality bases for per read in its 5' and 3' by evaluating the mean quality from a sliding window (like Trimmomatic but faster).
***Features - it can do all this, but you choose!***
* Filter out bad reads (too low quality, too short, or too many N...)
* Cut low quality bases for per read in its 5' and 3' by evaluating the mean quality from a sliding window (like Trimmomatic but faster)
* Trim all reads in front and tail
* Cut adapters. Adapter sequences can be automatically detected, which means you don't have to input the adapter sequences to trim them.
* Correct mismatched base pairs in overlapped regions of paired end reads, if one base is with high quality while the other is with ultra-low quality
* Trim polyG in 3' ends, which is commonly seen in NovaSeq/NextSeq data. Trim polyX in 3' ends to remove unwanted polyX tailing (i.e. polyA tailing for mRNA-Seq data)
* Report JSON format result for further interpreting
* Visualize quality control and filtering results on a single HTML page (like FASTQC but faster and more informative)
* Split the output to multiple files (0001.R1.gz, 0002.R1.gz...) to support parallel processing. Two modes can be used, limiting the total split file number, or limitting the lines of each split file (Not enabled in this Galaxy tool)
* Support long reads (data from PacBio / Nanopore devices)
However, with so many different features, more understanding is needed from its user.
</div>
</html>
- [x] **VI.** Go to the **Tools** section, and scroll down to **GENOMIC FILE MANIPULATION** -- **FASTA/FASTQ** [1], and click on **FastP** [2]. Choose `paired reads` [3]. Pick your input data and choose `Multiple datasets` [4] to include all fastq files and mark these blue. Enable `Overrepresented Sequence Analysis` and under `Output Options` choose YES to `Output JSON report`. Leave the remaining settings on default (i.e. what is already there). Press `Execute`.
<center> Figure 4 <img src="https://i.imgur.com/lQaTgBr.png" alt="drawing" width="800"/></center>
- [ ] **VII.** **Rename your files!** Fastp will produce four files for each isolate: two output reads (1 and 2, quality controlled and trimmed reads), one html report (3) and one json report (not shown). Now you need to rename them according to isolate name - it will get very messy after a while if not.
To rename your files, press Edit (arrow), rename. Be careful to name them correctly.
<center> Figure 4 <img src="https://i.imgur.com/WyhuPSB.png" alt="drawing" width="800"/></center>
- [ ] **VIII.** After renaming, we can inspect the `Fastp` results from one random read. Inspect the html-webpage by pressing the eye symbbol or downloading it and opening up in a web-browser.
**TASK: Answer the following questions**
1. What files was used as input?
3. How many basepairs are there in the fastq files analyzed?
4. Of what length is most of the reads (i.e. sequence length distribution)?
5. Considering that this is *Campylobacter jejuni*, which has a genome size of 1.6Mbp, what is the average coverage achieved by a single end sequence run?
**IX.** The Fastp results are for each fastq files, and in theory you could just inspect them one by one. However, when you have many files, it becomes a bit trickier. Its then useful to compile the quality assessment files using a tool called **MultiQC** [Ewels, P., Magnusson, M., Lundin, S., & Käller, M. (2016). MultiQC: summarize analysis results for multiple tools and samples in a single report. Bioinformatics, 32(19), 3047–3048.](https://doi.org/10.1093/bioinformatics/btw354)
<html>
<div style="background-color: #f2f3f4">
ABOUT: MultiQC aggregates results from bioinformatics analyses across many samples into a single report also viewable on a web-page. It's a general use tool, perfect for summarizing the output from numerous bioinformatics tools. Here we will summarize the data from the four FastQC analyzis.
</div>
</html>
Try to find the **MultiQC** software (see picture for hints if needed). Choose correct tool from which you wish **MultiQC** to aggregate the results from (hint: it starts with an "F" and ends with "astp"), and then choose all **JSON-files** from the **Fastp** analysis you just did (choose for all isolates). Execute.
Inspect the report from the **MultiQC** by downloading the webpage (**Webpage**) and opening up in a web-browser. It will have to be copied and inserted into a web-browswer to work properly. The text reports you can view in usegalaxy.com by clicking **Stats** in history.
***Are your reads ok?***
***What are duplicates and why do we have them?***
***What is the most likely reason for the poor per base sequence content?***
***Do you have any adapters left?***
Need help looking at MultiQC? [Video](https://www.youtube.com/watch?v=qPbIlO_KWN0 )
<br/>
<center>Figure 4<img src="https://i.imgur.com/gxolTZS.png" alt="drawing" width="800"/></center>
<br/>
#### Extra if time allows
**IXb.** Another quality control software is `FastQC`. Use the input reads and run `FastQC`, then `MultiQC`.
TASK
* What is the differences between Fastp and FastQC results and why do you get it?
## Assembly
**X**: Now its time to make assemblies! This will be done using an assembly pipeline **[Shovill](https://github.com/tseemann/shovill)**.
<html>
<div style="background-color: #f2f3f4">
ABOUT:
Main steps in **Shovill**:
1. Estimate genome size and read length from reads
2. Reduce FASTQ files to a sensible depth
3. Trim adapters from reads
4. Conservatively correct sequencing errors in reads
5. Pre-overlap paired-end reads
6. Assemble with SPAdes/SKESA/Megahit (you choose, SPAdes default) with modified kmer range and PE + long SE reads
7. Correct minor assembly errors by mapping reads back to contigs
8. Remove contigs that are too short, too low coverage, or pure homopolymers
9. Produce final FASTA with nicer names and parseable annotations
</div>
</html>
Find `Shovill` from `Tools`, and choose `Multiple Dataset` mode. Choose correct forward and reverse `*corrected and renamed reads*`, and leave `Trim reads` at default which is OFF. Under `Advance` options, write 1600000 under `Estimated genome size` and unclick `Disable post-assembly correction` (you want that). Run.
Hint: It might be that you have to autodetect filetype for all your `fastq` files output from `Fastp` anyway. Sometimes Galaxy is weird.
## Workflow
**XI**: Create your first workflow using the workflow editor. Follow and complete Steps 1.1 -- 1.3 of the [Galaxy Workflow Hands-On.](https://training.galaxyproject.org/training-material/topics/galaxy-interface/tutorials/workflow-editor/tutorial.html#introduction)
**XII**: Back to your own History. You have done each step individually, but you could make a workflow. In Galaxy, you have the option of turning your `History` into a workflow. In your `History` bar, press settings, and choose `Extract Workflow`. Name your Workflow something sensible, and view it under `Workflow`. We will use it in the next exercise (short, long and hybrid assembly-work). Optional: Follow the [hands-on on extracting a workflow from your history](https://training.galaxyproject.org/training-material/topics/galaxy-interface/tutorials/history-to-workflow/tutorial.html) from the Galaxy training materials.