自家製BLAST用DBから必要な配列エントリ取得
BLASTでヒットがあったエントリの配列データを切り出してきたい時はままあるかと。公共データベースのエントリならtogowsのRESTで取得するなり、ググるなりNCBIで検索するなりで配列データにありつけるかと思います。ところが、自分で作成したDBに対してBLASTした場合にはそうもいかないことがあるわけで。その場合どうするか、なのですが、とくに新規のコマンドをインストールしなくてもできるいい方法があります。
それにはまず、BLASTのDBを作成するコマンドmakeblastdbを実行するときに-parse_seqidsというオプションも足して実行します。例えば、自分の手持ちのDBがoreno.aaというFASTA形式のファイルだとすると
#!/bin/sh
makeblastdb -in oreno.aa -dbtype prot -hash_index -parse_seqids
な感じで。 そして、BLASTヒットがorenogene1だったとすると
#!/bin/sh
blastdbcmd -db oreno.aa -entry orenogene1
とすればFASTA形式でそのエントリの配列が取得できます。ただ、FASTAのヘッダ行のラベルが|(パイプ)などで長くなっていると、行頭からスペースが現れるまでの文字列を指定する必要があるようです。例えば、
>gi|5524211|gb|AAD44166.1| cytochrome b
のようになっている場合です。その場合は以下のように指定することになります。
#!/bin/sh
blastdbcmd -db oreno.aa -entry 'gi|5524211|gb|AAD44166.1|'