Nr数据库子库构建

# blast
wget https://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/LATEST/ncbi-blast-2.12.0+-x64-linux.tar.gz
tar -zxvf ncbi-blast-2.12.0+-x64-linux.tar.gz
echo 'PATH=your_path/blast/ncbi-blast-2.12.0+/bin:$PATH' >> ~/.bashrc
source ~/.bashrc

# aspera
wget -c https://d3gcli72yxqn2z.cloudfront.net/connect_latest/v4/bin/ibm-aspera-connect_4.1.0.46-linux_x86_64.tar.gz
tar -zxvf ibm-aspera-connect_4.1.0.46-linux_x86_64.tar.gz
# 自动安装
sh ibm-aspera-connect_4.1.0.46-linux_x86_64.tar.gz
echo 'PATH=~/.aspera/connect/bin/:$PATH' >> ~/.bashrc
source ~/.bashrc

ascp -i ~/.aspera/connect/etc/asperaweb_id_dsa.openssh -l 200M -k 1 --overwrite=diff -QT anonftp@ftp.ncbi.nlm.nih.gov:blast/db/FASTA/nr.gz ./

# 使用参数说明：
-i   # Use public key authentication and specify the private key file. 
     # 该密匙保存位置一般为：~/.aspera/connect/etc/
-l   # Set the target transfer rate in Kbps. Default: 10000
     # 即设置最大传输速率
-k   # 可选参数为{0|1|2|3}，即断点续传，一般设置为1
     # Enable resuming partially transferred files. (Default: 0). Please note that this must be specified for your first transfer; otherwise, it will not work for subsequent transfers.
     # “-k 1的具体说明”：-k 1 If the transfer stops because of connection loss or other issues the k option tells the transfer to resume from where it left off rather than restarting the entire transfer over. 
     ### 其他参数说明：
     # 0 Always retransfer the entire file.
     # 1 Check file attributes and resume if the current and original attributes match.
     # 2 Check file attributes and do a sparse file checksum; resume if the current and original attributes/checksums match.
     # 3 Check file attributes and do a full file checksum; resume if the current and original attributes/checksums match.
--overwrite=diff   # 默认设置即为diff，表示如果目前下载的文件与之前下载的文件不一致，则会overwrite
-Q   # Enable fair (-Q) or trickle (-QQ) transfer policy. Use the -l and -m to set the target and minimum rates.
-T   # Disable encryption for maximum throughput.

for i in {00..55}
do
   echo "ascp -i ~/.aspera/connect/etc/asperaweb_id_dsa.openssh -l 200M -k 1 --overwrite=diff -QT anonftp@ftp.ncbi.nlm.nih.gov:/blast/db/nr.${i}.tar.gz ./" >> aspera.preformatted.nr.sh
   echo "ascp -i ~/.aspera/connect/etc/asperaweb_id_dsa.openssh -l 200M -k 1 --overwrite=diff -QT anonftp@ftp.ncbi.nlm.nih.gov:/blast/db/nr.${i}.tar.gz.md5 ./" >> aspera.preformatted.nr.sh
done
ParaFly -c aspera.preformatted.nr.sh -CPU 1 &

# 检查下载的数据库是否完整
for i in {00..55};do md5sum -c nr.${i}.tar.gz.md5;done

perl ~/biosoft/blast/ncbi-blast-2.12.0+/bin/update_blastdb.pl --blastdb_version 5 --decompress nr --source gcp --passive

# 使用参数说明：
--blastdb_version  # Specify which BLAST database version to download (default: 5). 
                   # Supported values: 4, 5
--decompress       # Downloads, decompresses the archives in the current working directory, and deletes the downloaded archive to save disk space, while preserving the archive checksum files (default: false).
--source           # to download BLAST databases from (default: auto-detect closest location). Supported values: ncbi, aws, or gcp.
--passive          # Use passive FTP, useful when behind a firewall or working in the cloud (default: true). 
                   # To disable passive FTP, configure this option as follows: --passive no

# 网址：https://ftp.ncbi.nlm.nih.gov/blast/db/
# 自行编写shell脚本进行下载
for i in {00..55}
do
  echo "wget -c https://ftp.ncbi.nlm.nih.gov/blast/db/v5/nr.${i}.tar.gz ./" >> preformatted.wget.nr.sh
  echo "wget -c https://ftp.ncbi.nlm.nih.gov/blast/db/v5/nr.${i}.tar.gz.md5 ./" >> preformatted.wget.nr.sh
done

ParaFly -c preformatted.wget.nr.sh -CPU 10 2>err &

# taxonkit
conda install -c bioconda taxonkit
# csvtk
conda install csvtk

grep -P "\|\s+[pP]lant\w*\s*\|" ~/.taxonkit/names.dmp

echo Embryophyta | taxonkit name2taxid --show-rank

wget -c https://ftp.ncbi.nlm.nih.gov/pub/taxonomy/taxdump.tar.gz
wget -c https://ftp.ncbi.nlm.nih.gov/pub/taxonomy/taxdump.tar.gz.md5
md5sum -c taxdump.tar.gz.md5 # 检查下载数据库是否完整
mkdir ~/.taxonkit
tar zxf taxdump.tar.gz -C ~/.taxonkit

wget -c https://ftp.ncbi.nlm.nih.gov/pub/taxonomy/accession2taxid/prot.accession2taxid.gz
wget -c https://ftp.ncbi.nlm.nih.gov/pub/taxonomy/accession2taxid/prot.accession2taxid.gz.md5
md5sum -c prot.accession2taxid.gz.md5

# 获取对应类群的taxids
taxonkit list --ids 3193 --indent "" > plants.taxid.txt

# 获取对应类群的accession ID
pigz -dc prot.accession2taxid.gz | csvtk grep -t -f taxid -P plants.taxid.txt | csvtk cut -t -f accession.version,taxid | sed 1d > plant.acc2taxid.txt && cut -f 1 plant.acc2taxid.txt > plant.acc.txt &

get_species_taxids.sh -t 3193 > plants.blast.taxids.txt
# 得到上述taxids list之后，结合csvtk即可得到accession ID

# 解压下载好的preformatted nr database
for i in {00..55}
do 
  tar -zxvf nr.${i}.tar.gz -C ../nr_preformatted_20211122_decompress
done

blastdb_aliastool -seqidlist sequence.seq -db nr -out nr_plant -title nr_plant
# sequence.seq为对应类群的accession list

blastdbcmd -db nr -entry_batch sequence.seq -out - | pigz -c > blastdbcmd.nr.3193.fa.gz

# 提取preformatted nr数据库中的所有序列
blastdbcmd -db nr -dbtype prot -entry all -outfmt "%f" -out - | pigz -c > nr.fa.gz

time cat <(echo) <(pigz -dc nr.fa.gz) \
    | perl -e 'BEGIN{ $/ = "\n>"; <>; } while(<>){s/>$//;  $i = index $_, "\n"; $h = substr $_, 0, $i; $s = substr $_, $i+1; if ($h !~ />/) { print ">$_"; next; }; $h = ">$h"; while($h =~ />([^ ]+ .+?) ?(?=>|$)/g){ $h1 = $1; $h1 =~ s/^\W+//; print ">$h1\n$s";} } ' \
    | seqkit grep -f 3193.acc.txt -o nr.3193.fa.gz

seqkit stats blastdbcmd.nr.3193.fa.gz
file           format  type       num_seqs        sum_len  min_len  avg_len  max_len
blastdbcmd.nr.3193.fa.gz  FASTA   Protein  19,326,820  7,670,152,256        6    396.9   23,040

seqkit stats nr.3193.fa.gz
file           format  type       num_seqs        sum_len  min_len  avg_len  max_len
nr.3193.fa.gz  FASTA   Protein  19,326,820  7,670,152,256        6    396.9   23,040