一、质控

前面我们从GEO下好了SRA数据并转换为fastq文件，现在需要对fastq文件进行质控，这里用的软件为fastqc。首先建好文件夹用来存放数据

mkdir 01raw 02clean 03alin 04count

fastqc

1.png

R1=/mnt/d/bioinfo/project/rna/01raw/SRR15859344_1.fastq.gz
R2=/mnt/d/bioinfo/project/rna/01raw/SRR15859344_2.fastq.gz
fastqc $R1 $R2

这步会产生一个html文件，这就是fastqc的结果，我们可以双击在网页中打开进行查看。（详见下面介绍的fastqc报告解读）

2.png

fastp

用于数据过滤

3.png

fastp -i $R1 -I $R2 -o ../02clean/SRR15859344_1.clean.fastq.gz -O ../02clean/SRR15859344_2.clean.fastq.gz &

经过质控后的fastq文件就可以进行后面的比对分析了。

二、fastqc报告解读

本文参考：

https://blog.csdn.net/qq_44520665/article/details/113779792

1. Basic Statistics

Basic statistics是该fastq一些基本信息:

• Filename:文件名

• File type: 文件类型

• Encoding：测序平台的版本和相应的编码版本号，用于计算Phred反推error P时用

• Total Sequences: 输入文本的reads的数量

• Sequence length: 测序长度

• %GC: GC含量，表示整体序列的GC含量，由于二代测序GC偏好性高，且深度越高，GC含量会越高。

  
  
  4.png

2. Per base sequence quality

• 横轴为read长度，纵轴为质量得分，Q-score = -10✖lg（error P）
• 柱状表示该位置所有序列的测序质量的统计，柱状（黄色）是25%~75%区间质量分布，error bar（触须）是10%~90%区间质量分布，蓝线表示平均数，红色为中位数。
• 一般要求所有位置的10%小于20，即最多允许该位置10%的序列低于Q20，即90%的序列的碱基质量都大于Q20，即90%的序列碱基错误率不超过99%。当任何碱基质量低于10，或者任何中位数低于25时报WARN,需注意；当任何碱基质量低于5或者任何中位数低于20报FAIL。

5.png

碱基质量值Q

• 碱基质量值，Q，即每个碱基的正确识别率，是衡量测序质量的重要标准

• Q值通过测序Phred值计算而得，公式为：Q-score = -10 ✖ lg P。

• Phred值：不正确的碱基识别率，在碱基识别过程中通过一种概率模型计算得到，该模型可准确预测碱基识别的错误率。

• 碱基质量值越高表明碱基识别越可靠，准确度越高。

• Q20与Q30的含义：
  【1】Q20为每100个碱基中会有一个识别错，即正确识别率为2个9，99%，当Phred = 20 时，碱基识别出错率为1/100，碱基识别正确率为99%，Q-score = -10 ✖ lg 10-2=20
  【2】Q30为每1000个碱基中会有一个识别错，正确识别率为3个9，99.9%，当Phred = 30 时，碱基识别出错率为1/1000，碱基识别正确率为99.9%，Q-score = -10 ✖ lg 10-3=30

• Q30 > 90%,即碱基质量值 ≥ Q30的碱基所占百分比 ≥ 90%。

3. Per Sequence Quality Scores

• 每条reads的quality的均值的分布
• 横轴表示Q值，纵轴表示每个值对应的read数目，当测序结果主要集中在高分中，证明测序质量良好

• 当峰值小于27（错误率0.2%）时报"WARN"，当峰值小于20（错误率1%）时报"FAIL"

  
  
  6.png

4. Per Base Sequence Content

• 对所有reads的每一个位置，统计ATCG四种碱基（正常情况）的分布，用于检查是否有AT,GC分离现象
• 横轴为碱基长度分布，纵轴表示百分比，图中4条线分别代表A，C，T，G在每个位置上的平均含量
• 由于测序平台及测序长度不同，以及测序仪开始状态不稳定经常出现前后波动情况
• 好的样本中四条线应该平行且接近。当部分位置碱基的比例出现bias时，即四条线在某些位置纷乱交织，往往提示我们有overrepresented sequence的污染。当所有位置的碱基比例一致的表现出bias时，即四条线平行但分开，往往代表文库有bias (建库过程或本身特点)，或者是测序中的系统误差。
• 在碱基含量分布图，前几个碱基可能会出现较大波动，这是由于随机引物扩增偏差原因造成的

• 当任一位置的A/T比例与G/C比例相差超过10%，报"WARN"；当任一位置的A/T比例与G/C比例相差超过20%，报"FAIL"

  
  
  7.png

5. Per Sequence GC Content

• 统计reads的平均GC含量的分布
• 红线是实际情况，蓝线是理论分布（正态分布，均值不一定在50%，而是由平均GC含量推断的）
• 曲线形状的偏差往往是由于文库的污染或是部分reads构成的子集有偏差（overrepresented reads）。形状接近正态但偏离理论分布的情况提示我们可能有系统偏差
• 偏离理论分布的reads超过15%时，报"WARN"；偏离理论分布的reads超过30%时，报"FAIL"。

8.png

6. Per base N content

• 当出现测序仪不能分辨的碱基时会产生N，横轴为碱基分布，纵轴为N比率

• 当任一位置N的比率超过5%报WARN，超过20%报FAIL

  
  
  9.png

7. Sequence Length Distribution

• reads长度的分布

• 理论上每次测序仪测出的read长度时一致的，但是由于建库等因素通常会导致一些小片段，reads长度不一致时报"WARN"；当有长度为0的read时报“FAIL

  
  
  10.png

8. Sequence Duplication Levels

• 统计序列完全一致的reads的频率，横坐标是duplication的次数，纵坐标是duplicated reads的数目，以unique reads的总数作为100%。
• 上图的情况中，相当于unique reads数目～20%的reads是观察到两个重复的，～3%是观察到三次重复的，依此类推
• 一般测序深度越高，越容易产生一定程度的重复序列
• 如果原始数据很大（事实往往如此），做这样的统计将非常慢，所以fastqc中用fq数据的前200,000条reads统计其在全部数据中的重复情况。重复数目大于等于10的reads被合并统计，大于75bp的reads只取50bp（不知道怎么选的）进行比较。但由于reads越长越不容易完全相同（由测序错误导致），所以其重复程度仍有可能被低估。

• 当非unique的reads占总数的比例大于20%时，报"WARN"；当非unique的reads占总数的比例大于50%时，报"FAIL“

  
  
  11.png

9. Overrepresented sequences

• 如果有某个序列大量出现，就叫做over-represented。
• fastqc的标准是占全部reads的0.1%以上。和上面的duplicate analysis一样，为了计算方便，只取了fq数据的前200,000条reads进行统计，所以有可能over-represented reads不在里面。而且大于75bp的reads也是只取50bp。如果命令行中加入了-c contaminant file，出现的over-represented sequence会从contaminant_file里面找匹配的hit（至少20bp且最多一个mismatch），可以给我们一些线索。

• 当发现超过总reads数0.1%的reads时报”WARN“，当发现超过总reads数1%的reads时报”FAIL“。

  
  
  12.png

10. Adapter Content

• 横轴表示碱基位置，纵轴表示百分比
• 当fastqc分析时没有选择参数-a adapter list时，默认使用图例中的4种通用adapter序列进行统计。
• 若有adapter残留，后续必须去接头

13.png

RNA-seq入门

RNA-seq入门（一）数据下载和格式转换

RNA-seq入门（二）质控及fastqc报告解读

RNA-seq入门（三）比对、输出表达矩阵