原始数据来源于这篇文章https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE50177
1.下载原始数据
https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE50177
https://www.ncbi.nlm.nih.gov/sra?term=SRX339951
点击697.2Mb,然后进入下个页面:
点击SRR957677,这里插一句,下图layout指的是单端或者双端测序。single是单端,paired是双端。
在下一个页面里点击“Data access”,可以看到有一个NCBI的下载地址,copy地址,用wget下载。
wget https://sra-download.ncbi.nlm.nih.gov/traces/sra11/SRR/000935/SRR957678
现在下载的就是sra文件。
2 用sratoolkit从sra文件里提取fastq文件
(在sra文件所在的文件夹里提取,没有加目的文件夹)。--gzip是为了生成压缩的gz格式fastq文件。split-files是把提取的文件分成两份,unpaired的直接去掉。
fastq-dump --gzip --split-files SRR957677 #这是从单个sra文件里提取fastq文件
如果想从好几个sra文件里同时提取fastq文件,可以写一个小脚本运行:
vim fastqdump.sh
#!/bin/bash
for i in SRR* #for循环会遍历你指定的list里的每一个文件,并执行下面的命令
do
echo $i #显示文件名
fastq-dump --gzip --split-files $i #用fastqdump提取fastq文件
done
yanfang@YF-Lenovo:~/Documents$ ./fastqdump.sh#这里我尝试同时提取2个fastq文件,实际也不是同时,是一个一个来的
SRR957679
2019-08-28T17:27:31 fastq-dump.2.8.2 sys: timeout exhausted while reading file within network system module - mbedtls_ssl_read returned -76 ( NET - Reading information from the socket failed )
Read 19909740 spots for SRR957679
Written 19909740 spots for SRR957679
SRR957680
2019-08-28T17:34:48 fastq-dump.2.8.2 sys: timeout exhausted while reading file within network system module - mbedtls_ssl_read returned -76 ( NET - Reading information from the socket failed )
Read 24231941 spots for SRR957680
Written 24231941 spots for SRR957680
3.用安装好的fastqc查看(安装过程略)
直接在终端里调用fastqc
fastqc
点击file,open,选择要查看的fastq文件,然后fastqc会自动分析文件,最后生成一个报告。下面是报告的各项结果:
Basic Statistics:报告整体浏览
Measure | Value |
---|---|
Filename | SRR957678_1.fastq.gz |
File type | Conventional base calls |
Encoding | Sanger / Illumina 1.9 |
Total Sequences | 8828013 |
Sequences flagged as poor quality | 0 |
Sequence length | 50 |
%GC | 46 |
Filename:文件名
File type: 文件类型
Encoding:测序平台的版本和相应的编码版本号,用于计算Phred反推error P时用。
Total Sequences: 输入文本的reads的数量
Sequence length: 测序长度
%GC: GC含量,表示整体序列的GC含量,由于二代测序GC偏好性高,且深度越高,GC含量会越高。
Per base sequence quality:某一位置上所有读段的测序质量评分
(最主要看得数据信息)quality就是Fred值,一条reads某个位置上出错概率为0.01时,quality值就是20,即常说的Q20。就是一个箱线图boxplot,黄色箱子(25%和75%的分数线),红色线(中位数),蓝线是平均数,下面和上面的触须分别表示 10%和 90%的点横坐标reads的碱基位置,最大值即为读长,纵坐标代表质量的好坏(判断的准确性)。如果任何一个位置的下四分位数小于10或者中位数小于25,会显示“警告”;如果任何一个位置的下四分位数小于5或者中位数小于20,会显示“不合格”。这个结果相对来说还是比较好的。
Per tile sequence quality:
图中横轴代表碱基位置,纵轴代表 tile 编号。图中的颜色是从冷色调到暖色调的渐变,冷色调表示这个 tile 在这个位置上的质量值高于所有 tile 在这个位置上的平均质量值,暖色调表示这个 tile 的在这个位置上的质量值比其它 tiles 要差;一个很好的结果,整张图都应该是蓝色,简单来说,就是看图内有无除蓝色外的亮点,有亮点代表低于平均值。当某些tail出现暖色,在后续的分析种把该tail测序结果全部去除。
Per sequence quality scores每条序列平均碱基质量分数
图中横轴为测序质量值,纵轴为 reads 数。红线上的每一个点表示quality值所对应的reads的数量,其面积就是总的reads数。如果最高峰所对应的横坐标质量值小于 27 (错误率 0.2 %) 则会显示“警告”,如果最高峰的质量值小于 20 (错误率 1 %) 则会显示“不合格”。如图所示红线单峰,分值在38左右,所以reads很可靠。
Per base sequence content每个位置的4种碱基组成比例
一个完全随机的文库内每个位置上 4 种碱基的比例应该大致相同,因此图中的四条线应该相互平行且接近25的位置左右。在 reads 开头出现碱基组成偏离往往是建库操作造成的,在reads上加接头的碱基组成不是均一的。会造成明显的碱基组成偏离。如果任何一个位置上的A和T之间或者G和C之间的比例相差10%以上则报“警告”,任何一个位置上的A和T之间或者G和C之间的比例相差 20%以上则报“不合格”。此结果总体上处于25%左右,且A和T比例相等,G和C比例相等,说明质量可以,但在前15个bp位置上严重分离,说明有碱基偏向性。可能有接头的污染。也有可能由于测序平台及测序长度不同,以及测序仪开始状态不稳定经常出现前后波动情况。
Per sequence GC contentGC含量:
横轴表示GC含量,纵轴表示不同GC含量对应的read数,蓝色为理论值,红色是真实值。在一个正常的随机文库中,GC 含量的分布应接近正态分布,且中心的峰值和所测基因组的 GC 含量一致。如果出现不正常的尖峰分布,则说明文库可能有污染, (如果是接头的污染,那么在 overrepresented sequences 那部分结果还会得到提示),或者存在其它形式的偏选;总体上就是看红色的线与蓝色线正态分布趋势是否接近。此图可知道红色线与蓝色线较为接近,质量较好。
Per base N content每个碱基上N的比例
当出现测序仪不能分辨的碱基时会产生N,横轴为碱基分布,纵轴为N比率。如果任何一个位置 N 的比例大于5%则报“警告”,大于20%则报“失败”。此图可知基本无N,皆已测得为ATGC的碱基。测序质量较好。
Sequence Length Distribution Reads的长度分布
测序仪出来的原始 reads 通常是均一长度的,但经过质控软件等处理过的数据则不然;经过质控软件处理过的reads长度则不一样。当 reads 长度不一致时报“警告”,当有长度为 0 的 reads 时则报“不合格”。此图可知为测序仪产出的reads,长度皆为50bp。
Sequence Duplication Levels序列重复的水平
图中横轴代表 reads 的重复次数 ( 1 表示 unique 的序列,2 表示有 2 条完全相同的 reads ...),大于 10 次重复后则按不同的重复次数合并显示。 纵坐标表示各重复次数下的 reads 数占总 reads 的百分比。蓝线展示所有 reads 的重复情况,红线表示在去掉重复以后,原重复水平下的 reads 占去重后 reads 总数的百分比。如果非 unique 的 reads 占总 reads 数的 20 % 以上则报 ”警告“,占总 read 数的 50 % 以上则报 ”不合格“。不合格报错对于此项是正常现象,不需要太过关注。一般测序深度越高,越容易产生一定程度的重复序列。
Overrepresented sequences大量重复出现的序列
这个样品的此项结果为No overrepresented sequences。
我在网上搜到了一个例子:
显示同一条 read 出现次数超过总测序 reads 数的0.1%的统计情况。正常文库内序列的多样性水平很高,不会有同一条 read 大量出现的情况,这部分结果会把大量出现的 reads 列出来,并给出可能来源。如果有任何 read 出现的比例超过总 reads 数的0.1%则报警告,超过总 reads 数的1%则报不合格。如果检测出一条多重复序列,重复次数较多,推测可能是TrueSeq接头序列。
Adapter Content接头含量
显示 reads 中的接头含量,并显示可能的来源。图中横轴为碱基位置,纵轴为含有接头序列的比例。正常的情况下接头的含量应该接近0,如果 reads 中的接头含量过高,说明文库内小片段比例偏高 (这个可以从文库质检报告中看出来),这可能是由于片段选择时选取的长度偏短或者使用切胶的方式回收片段时上样过多致使小片段不能很好的分离等原因造成的;如果接头的含量随着碱基的位置增大而逐渐升高,则表示 reads 中含有接头 (如图所示),这部分接头会影响后续的分析,我们需要截掉 reads 中的接头序列或者将含有接头的 reads 完全删除。如果任何重复 read 超过总 reads 数的5%则报 '警告', 超过总 reads 数的10% 则报 '不合格,由图可知测序是没有接头污染的。如果有接头污染,在序列尾端会出现一个上扬的曲线。
**以上fastqc质控的图是SRR957678的结果。图解摘自两篇文章:
//www.greatytc.com/p/bacb86c78b43
//www.greatytc.com/p/f510dce0ab8c
还有一个英文版的fastqc质量报告解读:https://dnacore.missouri.edu/PDF/FastQC_Manual.pdf