获取数据的来源

获取序列信息的数据库有很多，首先介绍从NCBI获取数据：

想从NCBI上获取数据，学会使用Entrez是必不可少的！什么是Entrez？

Entrez 是美国国家生物技术信息中心所提供的在线资源检索器。该资源将GenBank序列与其原始文献出处链接在一起。 Entrez 是由NCBI主持的一个数据库检索系统。

Entrez中核酸数据库为：GenBank,EMBL,DDBJ

蛋白质数据库为：Swiss-Prot,PIR,PFR,PDB

PubMed

基因组和染色体图谱资料

NCBI上是如何整合数据的呢？

LOCUS   AF086833     18959 bp      cRNA linear         VRL      13-FEB-2012

DEFINITION    Ebola      virus - Mayinga, Zaire, 1976, complete genome.

ACCESSION   AF086833

VERSION         AF086833.2   GI:10141003

LOCUS：这一行里包括基因座的名字,核酸序列长度,分子的类别,拓扑类型,原核生物的基因拓扑类型都是线性的,最后是更新日期。

DEFINITION：对序列的简短定义。

ACCESSION：就是在搜索条中输入的那个数据库编号,也叫做检索号,每条记录的检索号在数据库中是唯一且不变的。即使数据提交者改变了数据内容,Accession 也不会变。

Version：版本号和 Locus,Accession 长得差不多。版本号的格式是“检索号点上一个数字”。版本号于 1999 年 2 月由三大数据库采纳使用。主要用于识别数据库中一条单一的特定核苷酸序列。在数据库中,如果某条序列发生了改变,即使是单碱基的改变,它的版本号都将增加,而它的 Accession 也就是检索号保持不变。比如,版本号由 U12345.1 变为 U12345.2,而检索号依然是 U12345。版本号后面还有个 GI 号。

GI号：与前面的版本号系统是平行运行的。当一条序列改变后,它将被赋予一个新的 GI 号,同时它的版本号将增加。

Locus 是一个同学的真实姓名,而 Accession 是这个同学的学号。同一个人在不同的学校里会有不同的学号,而名字只有一个。基因也是一样,同一个基因在不同的数据库中会有不同的检索号,而基因的名字只有一个。

如何安装和使用Entrez？ 参见：Entrez安装和使用说明

使用示例：

#查看其序列信息（前10行）                                                                                   

efetch -db=nuccore -format=gb -id=AF086833 | head

#导入Genbank format                                                                                              

efetch -db=nuccore -format=gb -id=AF086833 > AF086833.gb

#导入Fasta format                                                                                                    

efetch -db=nuccore -format=fasta -id=AF086833 > AF086833.fa

#查看序列的前三个碱基(区域可选)                                                                             

efetch -db=nuccore -format=gb -id=AF086833 -seq_start=1 -seq_stop=3

#查看序列的反向互补序列                                                                                      

efetch -db=nuccore -format=fasta -id=AF086833 -seq_start=1 -seq_stop=5 -strand=1                                                                                                                   

efetch -db=nuccore -format=fasta -id=AF086833 -seq_start=1 -seq_stop=5 -strand=2

---

如何用Entrez搜索？

#在核酸与蛋白质数据库中搜索‘PRJNA257197’                                                

esearch -db nucleotide -query PRJNA257197                                                 

esearch -db protein -query PRJNA257197

#结果如下：

#将搜索到的核酸与蛋白质序列全部抓取                                                          

esearch -db nucleotide -query PRJNA257197 | efetch -format=fasta >genomes.fa                                                                                                            

esearch -db protein -query PRJNA257197 | efetch -format=fasta > proteins.fa

如何从Short Read Archive（SRA）获取数据？

利用SRA Toolkit，使用与安装见以下链接：

安装与使用说明

#将sra数据转换成fastq格式                                                                                      

fastq-dump

#下载sam格式数据                                                                                                         

sam-dump

#产生一个SRR1553607.fastq                                                                                     

fastq-dump SRR1553607

fastq -dump --split-files SRR1553607                                                                             # --split-files：将每一个read读到单独的文件中。文件将接收到读取相应的后缀数。针对pairend 测序结果，将原文件转换成两个fastq文件，里面的reads是成对的，文件名其中一个是*_1.fastq,一个是*_2.fastq                                                           得到：              

SRR1553607_1.fastq                                                                               

SRR1553607_2.fastq

fastq-dump --split-files -X 10000 SRR1553607                                                           

#提取需要的子集（前10000个读长）

#简要查看数据信息                                                                                                      

sra-stat --xml --quick SRR1553607

#查看一个数据的读长统计                                                                                             

sra-stat --xml --statistics SRR4237168

如何批量获取SRA数据？

以 bioproject number PRJNA257197为例：

#首先搜索数据库并下载runinfo并导入到一个csv文本中：                                 

esearch -db sra -query PRJNA257197 | efetch -format runinfo > info.csv

#查看第一列信息：                                                                                                        

cat info.csv|cut -d ',' -f 1|head

#将第一列SRR开头的id写入一个文本：                                                                    

cat info.csv | cut -f 1 -d ',' | grep 'SRR' | head > ids.txt

#获取ids.txt中id的前10000个读长                                                                                 

cat ids.txt | xargs -n 1 fastq-dump -X 10000 --split-files $1

#还可以通过筛选日期来得到想要的序列                                                                      

cat info.csv | grep '2014-08-19' | cut -f 1 -d ',' | grep SRR | head -10 > ids.txt

下面说说如何查看测序数据信息

#根据发布日期获取测序数据信息：                                                                        

esearch -db sra -query '"2015/05"[Publication Date]'|efetch -format runinfo>2015-05.cs

#收集SRR开头为编号的数据：                                                                                  

cat sra/*.csv | grep 'SRR' > allruns.csv

下面利用csvcut查看csv文件：

#统计一共测序了多少个碱基：                                                                                   

cat allruns.csv | csvcut -c 5 | datamash sum                                       

1 4.966817753558e+15    #如果报错可删除对应的行

#查看测序物种并统计频次：                                                                                       

cat allruns.csv|csvcut -c 29|sort|uniq -c|sort -nr|head -25

#对使用的测序平台进行统计：                                                                                    

cat allruns.csv|csvcut -c 19|sort|uniq -c|sort -rn|head -25

#查看测序仪型号：                                                                                                         

cat allruns.csv|csvcut -c 20|sort|uniq -c|sort -rn|head -25

#查看哪个数据容量最大：                                                                                            

cat allruns.csv|csvcut -c 1,8|datamash -f -t ',' mean 2 sum 2 max 2