学校网课笔记~
今天的网课主要介绍怎么把原始的fastq文件根据fastqc结果进行去接头、以及过滤掉低质量的reads。
主讲人介绍目前主要使用的去接头/过滤的软件有三种：
- TrimGalore (+ CutAdapt) https://github.com/FelixKrueger/TrimGalore/blob/master/Docs/Trim_Galore_User_Guide.md
- Trimmomatic
http://www.usadellab.org/cms/?page=trimmomatic - Fastx-toolkit
http://hannonlab.cshl.edu/fastx_toolkit/
主讲人推荐的是第一个软件。TrimGalore软件里各种参数详解在这里：here
先看一下原始fastq文件的QC结果：
可以看到这里有adapter是需要去掉的，有些read的质量也不太好。还是用脚本来写命令行（服务器运行）：
#!/bin/bash
#SBATCH --job-name=trimming_example # Job name
#SBATCH --mail-type=END,FAIL # Mail events (NONE, BEGIN, END, FAIL, ALL)
#SBATCH --mail-user=***************** # Where to send mail
#SBATCH --ntasks=1 # Run on a single CPU
#SBATCH --mem=4gb # Job memory request
#SBATCH --time=02:00:00 # Time limit hrs:min:sec
#SBATCH --output=/gpfs/home/ID/log_files/trimming_%j.log # Standard output and error log
#SBATCH -p cpu_short
module load trimgalore/0.5.0
module load python/cpu/2.7.15-ES
module load fastqc/0.11.7
# Trimming
trim_galore --paired --q 30 --gzip --fastqc /gpfs/home/ID/download/SRR1523671_1.fastq.gz /gpfs/home/ID/download/SRR1523671_2.fastq.gz -o /gpfs/home/ID/trimmed_files/
#QC again
fastqc -o QC/ /gpfs/home/ID/trimmed_files/*
上面需要注意的参数可能就是--paired和--q 30这两个参数，因为我的文件是双端测序的结果，所以要注明是--paired，并且要跟着两个fastq文件。还有就是你要过滤的read质量，这里有个表：
一般过滤的标准是q 30 或者是q 20，一般不用40,50,60，如果标准设置过高，你有可能就把一些重要的信息给过滤掉了
运行后会生成一个txt文件（每一个fastq生成一个文件），是你过滤fastq后的一个summary：
SUMMARISING RUN PARAMETERS
==========================
Input filename: /gpfs/home/ID/download/SRR1523671_1.fastq.gz
Trimming mode: paired-end
Trim Galore version: 0.5.0
Cutadapt version: 1.8.2
Quality Phred score cutoff: 30
Quality encoding type selected: ASCII+33
Adapter sequence: 'AGATCGGAAGAGC' (Illumina TruSeq, Sanger iPCR; auto-detected) #这里会注明你用的是哪一种adaptor
Maximum trimming error rate: 0.1 (default)
Minimum required adapter overlap (stringency): 1 bp
Minimum required sequence length for both reads before a sequence pair gets removed: 20 bp
Running FastQC on the data once trimming has completed
Output file will be GZIP compressed
This is cutadapt 1.8.2 with Python 2.7.15 #这个软件是要在python环境下运行的
Command line parameters: -f fastq -e 0.1 -q 30 -O 1 -a AGATCGGAAGAGC /gpfs/home/ID/download/SRR1523671_1.fas
tq.gz
Trimming 1 adapter with at most 10.0% errors in single-end mode ...
Finished in 316.18 s (15 us/read; 4.00 M reads/minute).
=== Summary ===
Total reads processed: 21,069,824 #所有的reads
Reads with adapters: 5,833,844 (27.7%) #有adaptor的reads
Reads written (passing filters): 21,069,824 (100.0%)
Total basepairs processed: 2,128,052,224 bp
Quality-trimmed: 387,450,070 bp (18.2%)
Total written (filtered): 1,719,160,219 bp (80.8%)
=== Adapter 1 ===
Sequence: AGATCGGAAGAGC; Type: regular 3'; Length: 13; Trimmed: 5833844 times.
No. of allowed errors:
0-9 bp: 0; 10-13 bp: 1
Bases preceding removed adapters:
A: 22.2%
C: 34.9%
G: 22.0%
T: 20.3%
none/other: 0.6%
Overview of removed sequences #这里是一长串的记录，有哪些reads被去掉了
length count expect max.err error counts
1 3396604 5267456.0 0 3396604
2 1146056 1316864.0 0 1146056
3 310513 329216.0 0 310513
4 109696 82304.0 0 109696
......
RUN STATISTICS FOR INPUT FILE: /gpfs/home/ID/download/SRR1523671_1.fastq.gz
=============================================
21069824 sequences processed in total
过滤后，再做一次QC，做个对比：
