NGS分析入门:测序原始数据操作

我们在完成了测序之后的第一步工具就是把测序的结果拷贝到指定的地方去,有很多种方法,比如:

  • 使用scp命令进行拷贝。其格式为:
    scp [options] username1@source_host:directory1/filename1 username2@destination_host:directory2/filename2
  • 使用ftp,http, smb工具进行拷贝。比如mac os中的Fetch工具。这些工具多半有图型化界面,比如容易掌握。也可以使用命令行比如wget进行拷贝。
  • 从GEO/UCSC等在线资源下载数据。

我们需要知道的几种数据源:

一般来说,原始的测序结果都会保存在fastq文件中。当我们下载到一个sra文件时,需要把它转换成fastq,需要用到工具:fastq-dump。这个是SRA toolkit当中的一个工具。

一个标准的fastq文件其格式如下:

@SEQ_ID
GATTTGGGGTTCAAAGCAGTATCGATCAAATAGTAAATCCATTTGTTCAACTCACAGTTT
+
!''*((((***+))%%%++)(%%%%).1***-+*''))**55CCF>>>>>>CCCCCCC65

其文件第一行为序列名,第二行为序列,第三行是序列终止标记,第四行为序列相对应每个位置的测序质量。
对于Illumina测序的序列名,比如:

@HWUSI-EAS100R:6:73:941:1973#0/1

它用”:”分格,分别代表下列信息:

HWUSI-EAS100R the unique instrument name
6 flowcell lane
73 tile number within the flowcell lane
941 ‘x’-coordinate of the cluster within the tile
1973 ‘y’-coordinate of the cluster within the tile
#0 index number for a multiplexed sample (0 for no indexing)
/1 the member of a pair, /1 or /2 (paired-end or mate-pair reads only)

 

而1.8版本的会变成

@EAS139:136:FC706VJ:2:2104:15343:197393 1:Y:18:ATCACG

,其含义为:

EAS139 the unique instrument name
136 the run id
FC706VJ the flowcell id
2 flowcell lane
2104 tile number within the flowcell lane
15343 ‘x’-coordinate of the cluster within the tile
197393 ‘y’-coordinate of the cluster within the tile
1 the member of a pair, 1 or 2 (paired-end or mate-pair reads only)
Y Y if the read fails filter (read is bad), N otherwise
18 0 when none of the control bits are on, otherwise it is an even number
ATCACG index sequence

在这里需要特别指出的是字母Y/N在1.8中的作用,它表明了read的结果是否可以接受。可以通过这个标记对fastq文件进行筛选。

对于记录质量的那一行,其记录是的phred质量:phred=-10*log(p/(1-p), base=10)。值越大表明其质量越好。当p-value等于0.05时,phred值应该是13左右。但是我们看到的那一行,并不是数字。它其实是数字,只是在广本中把它当成一个字符显示了出来。但是不同的fastq文件,对质量的记录是有区别的,

  SSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSS.....................................................
  ..........................XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX......................
  ...............................IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII......................
  .................................JJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJ......................
  LLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLL....................................................
  !"#$%&'()*+,-./0123456789:;<=>?@ABCDEFGHIJKLMNOPQRSTUVWXYZ[\]^_`abcdefghijklmnopqrstuvwxyz{|}~
  |                         |    |        |                              |                     |
 33                        59   64       73                            104                   126
  0........................26...31.......40                                
                           -5....0........9.............................40 
                                 0........9.............................40 
                                    3.....9.............................40 
  0........................26...31........41                               
 
 S - Sanger        Phred+33,  raw reads typically (0, 40)
 X - Solexa        Solexa+64, raw reads typically (-5, 40)
 I - Illumina 1.3+ Phred+64,  raw reads typically (0, 40)
 J - Illumina 1.5+ Phred+64,  raw reads typically (3, 40)
     with 0=unused, 1=unused, 2=Read Segment Quality Control Indicator (bold) 
     (Note: See discussion above).
 L - Illumina 1.8+ Phred+33,  raw reads typically (0, 41)

有时候我们拿到的数据是经过压缩的,我们需要了解两个工具,解压/压缩工具,以及md5工具。通常的工具有tar, gzip/gunzip, zip/unzip等等。md5的工具就是md5。它的作用在于检查文件的传输的过程中是否产生了变化,如果md5的结果不一致,表明数据可能有所损失。对于这些工具的使用,请google。

10 thoughts on “NGS分析入门:测序原始数据操作

  1. 你好,我做HIV-1 的illumina的深度测序,初次接触,对数据处理一点不懂 怎么办呢,能否告知您的联系方式 谢谢

  2. 博主你好,首先感谢下楼主分享经验,我刚接触这个领域,看了你的文章受益匪浅。
    我前几天做demultiplexing,output fastq中有好多的undetermined(Miseq数据,用bcl2fastq跑的),其中大部分是forward barcode 显示的是adapter的一部分,这个是什么原因呢?

  3. 楼主: 我想问一下NIH roadmap project 上下载的bed wig数据可以直接解释出来么,我搜了好几天,但是没有结构。然后我想用matlab做可行么?

    1. 它们都是文本文件,除非你下载的是bigwig或者bigbed文件。 它们的描述文件在:http://genome.ucsc.edu/FAQ/FAQformat.html. 你完全可以使用matlab做。但是我对matlab不熟悉,所以不知道具体用什么包。在R/Bioconductor中,使用rtracklayer或者trackViewer就可以很好的import这些数据。

发表评论

电子邮件地址不会被公开。 必填项已用*标注