转录组注释离不开Blast+,需要比对到NCBI nt nr等数据库。

1. 软件下载:ftp://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/LATEST/     如:ncbi-blast-2.7.1+-x64-linux.tar.gz,解压后可用。

2. 使用

 (1) 建立数据库索引

        NCBI的nt和nr有直接使用的数据库,不需要再做处理,当然也可以下载fasta格式数据库再做makeblastdb。

 $ makeblastdb -in input.fasta -dbtype molecule_type -title database_title -parse_seqids -out database_name -logfile File_Name

   -dbtype    nucl为核酸,prot为蛋白

   -title 给数据库起个名

   -parse_seqids 推荐加上

   -out 后接数据库名,以后blast+搜索时要用到的-db的参数

   -logfile 日志文件,如果没有默认输出到屏幕

 (2) 作比对

  核酸序列比对核酸数据库:
          $ blastn -query xx.fasta -db yyy -out zzz -evalue 1e-5 -max_target_seqs 10 -outfmt 6 -num_threads 12
                不同物种之间e-5,同物种间e-10,线程数最多12,输出格式6为列表式,可以自己设置如:
                                   -outfmt '6 qseqid qlen sseqid sgi slen pident length mismatch gapopen qstart qend sstart send evalue bitscore staxid ssciname'
  核酸序列比对蛋白质数据库:
        $ blastx -query -out -db -outfmt 6 -evalue 1e-5 -max_target_seqs 10 -num_threads 48
             不同物种之间e-5,同物种间e-10,线程数最多12,输出格式6为列表式,可以自己设置如:
                                   -outfmt '6 qseqid qlen sseqid sgi slen pident length mismatch gapopen qstart qend sstart send evalue bitscore staxid ssciname'
                                                                                                                                                                                                                                                        BLASTn tabular output format 6
                                                                      qseqid sseqid pident length mismatch gapopen qstart qend sstart send evalue bitscore
 1.  qseqid  query (e.g., gene) sequence id
 2.  sseqid  subject (e.g., reference geome) sequence id
 3.  pident  percentage of identical matches
 4.  length  alignment length
 5.  mismatch  number of mismatches
 6.  gapopen  number of gap openings
 7.  qstart  start of alignment in query
 8.  qend  end of alignment in query
 9.  sstart  start of alignment in subject
 10.  send  end of alignment in subject
 11.  evalue  expect value
 12.  bitscore  bit score
 
3. 数据库下载地址
         https://www.ncbi.nlm.nih.gov/public/使用Aspera下载速度非常快。
   Nt/Nr数据库:ftp://ftp.ncbi.nlm.nih.gov/blast/db/(直接使用)    ftp://ftp.ncbi.nlm.nih.gov/blast/db/FASTA/(需建索引) 
   KOG/COG 数据库:ftp://ftp.ncbi.nlm.nih.gov/pub/COG/KOG/
   Swiss-Prot/TrEMBL/UniRef90数据库:http://www.uniprot.org/downloads

4. 从blast结果中取出每个query搜到的evalue最小的结果(参考:https://blog.csdn.net/weixin_40099163/article/details/83215747
   blast的-outfmt 6格式默认最上面的比对到的evalue最小,因此可以利用awk根据第一列去重,默认会保留最上面的一条记录,即evalue最小值。
   awk '!a[$1]++{print}' 1.blast > uniq.blast