KEGG数据库包含丰富的功能信息及通路注释，可由此了解生物系统的功能和效用。

而官方已经不再免费提供FTP下载氨基酸以及核酸序列库，并且KAAS和Kofam等官方工具，要么在线注释提交序列过少，要么注释结果会关联到非相关物种，多多少少存在一些问题。

但通过KEGG API检索可以获得一些信息，如：
① KEGG GeneID和KEGG ORTHOLOG的对应关系

KEGG GeneID 与 KEGG ORTHOLOG 对应关系

② NCBI RefSeq与KEGG GeneID的对应关系

RefSeq ID 与 KEGG GeneID 对应关系

所有物种信息

因此可以曲线救国，通过KEGG API检索获得已有基因组注释序列的相关信息，将NCBI proteinid映射到KEGG GeneID，关联至KEGG ORTHOLOG注释信息，从而构建本地数据库进行注释。以下是对构建本地植物KEGG数据库及注释的一次记录。

下载已有基因组蛋白序列进行建库及注释

1.1 下载所有物种列表

curl https://rest.kegg.jp/list/organism > KEGG.Organism.txt

1.2 获取植物物种列表

grep "Plants" KEGG.Organism.txt > KEGG.Plants.txt

植物物种列表：KEGG.Plants.txt

awk -F'\t' '{print "https://www.genome.jp/kegg-bin/show_organism?org=" $1}' KEGG.Plants.txt > KEGG.Plants.link.txt

物种详细信息页面链接：KEGG.Plants.link.txt

1.4 向网页发出GET请求，从返回内容中获取NCBI基因组信息，生成得到对应的蛋白序列下载链接

# 有133个物种，但只有130个物种有蛋白序列信息
缺少的为：Lotus japonicus(lja)  / Oryza sativa japonica (Japanese rice) (dosa) / Bathycoccus prasinos(bpg)
cat KEGG.Plants.link.txt | xargs -n 1 curl -L | grep -wio "https://ftp.ncbi.nlm.nih.gov/genomes[0-9_a-zA-Z\/.\-]*" > KEGG2NCBI.txt

NCBI ftp下载页面链接

cat KEGG2NCBI.txt | awk -F "/" '{print $NF}' | sed "s/$/&_protein.faa.gz/g" > KEGG2NCBI.suffix.txt
paste -d '/' KEGG2NCBI.txt KEGG2NCBI.suffix.txt > KEGG2NCBI.prodown.txt
# 注意NCBI ftp下载链接后面重复了一项，如
拟南芥的基因组页面链接为：https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/001/735/GCF_000001735.4_TAIR10.1
蛋白信息下载链接为：https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/001/735/GCF_000001735.4_TAIR10.1/GCF_000001735.4_TAIR10.1_protein.faa.gz

生成得到下载链接：KEGG2NCBI.prodown.txt

1.5 下载和解压所有蛋白序列（因为文件比较小，可以直接wget下载）

wget -i KEGG2NCBI.prodown.txt -P ./kegg_refseq
pigz -k -d *.faa.gz -p 16

下载所得蛋白序列文件：*.faa.gz & *.faa

1.6 合并蛋白序列，查看数量

cat *.faa > kegg_plants_refseq.faa
grep -c '>' kegg_plants_refseq.faa

kegg_plants_refseq.faa

1.7 diamond建库

diamond makedb --in kegg_plants_refseq.faa --db kp_refseq.dmnd -p 16

1.8 diamond比对注释

diamond blastx -q test.fa -d kp_refseq.dmnd \
--max-hsps 1 --max-target-seqs 1 --sensitive --outfmt 6 --evalue 1e-5 -p 16 \
-o kegg_test_match.out

输出的是熟悉的outfmt6格式

kegg_test_match.out

获取NCBIID与KO号对应关系

2.1 以拟南芥为例，查看ncbiid与KO号对应关系

ncbiid与KO号对应关系

要获取所有植物物种的对应关系，可以直接替换链接中的物种缩写，写个脚本

bash get_ncbi_protein_id.sh

#bash脚本如下：
# 首先，读取input.txt文件中的第二列数据，并将其存储在一个名为"ids"的变量中
ids=$(awk '{print $2}' KEGG.Plants.txt)
# 然后，遍历ids变量中的每一个值，并将其替换为预定义的格式
for id in $ids; do
 # 将每个值替换为http://rest.kegg.jp/conv/<id>/ncbi-proteinid的格式
 # 其中<id>表示原始值
 result=$(echo "http://rest.kegg.jp/conv/${id}/ncbi-proteinid")
 # 将替换结果追加到output.txt文件中
 echo $result >> KEGG.Plants.ncbiproid.temp.txt
done

cat KEGG.Plants.ncbiproid.temp.txt | xargs -n 1 curl -L> KEGG.Plants.ncbiproid.txt

KEGG.Plants.ncbiproid.txt

2.2 以拟南芥为例，查看 KEGG 的基因ID与 KO 注释对应关系

KEGG 的基因ID与 KO 注释对应关系

现在获取KEGG 的基因ID与 KO 注释列表

awk -F'\t' '{print "http://rest.kegg.jp/link/ko/" $2}' KEGG.Plants.txt > KEGG.Plants.Knum.temp.txt
cat KEGG.Plants.Knum.temp.txt | xargs -n 1 curl -L> KEGG.Plants.Knum.txt

2.3 Rstudio下整理KO关联信息，第一列为NCBIID，第二列为Knum

setwd("~/temp/kegg_get")
ncbiid <- read.table("KEGG.Plants.ncbiproid.txt",header = F,sep = "\t")
protid<- ncbiid %>% separate(V1, c("ncbi","protid"), "[:]")
protid <- protid[,c(2,3)]
colnames(protid) <- c("protid","geneid")
head(protid)
knum <- read.table("KEGG.Plants.Knum.txt",header = F,sep = "\t")
koid <- knum %>% separate(V2, c("ko","kid"), "[:]")
koid <- koid[,c(1,3)]
colnames(koid) <- c("geneid","kid")
kplant <- protid %>% left_join(koid)
kplant <- na.omit(kplant)
kplant <- kplant[,c(1,3)]
write.table(kplant,"Plants.NCBI2KEGG.txt",sep = "\t",row.names = F,quote=F)

Plants.NCBI2KEGG.txt

到这里可以把一些中间文件删去，仅保留以下文件，方便以后注释

保留文件

整理关联信息，获得KEGG注释

3.1 提取diamond比对结果，第一列为比对id，第二列为蛋白id，需要去除小数点后面的字符

cat kegg_test_match.out | awk -F "\t" '{split($2,a,"."); print $1 "\t" a[1]}' > kegg_test_match_exct.out

3.2 整理对应关系，获得KEGG注释

python get_kegg_annotation.py

#python脚本：get_kegg_annotation.py，需要提前安装pandas
import pandas as pd
# 读取文件并将它们存储在DataFrame中
df1 = pd.read_csv("kegg_test_match_exct.out", sep="\t", header=None, names=["col1", "col2"]) # "kegg_test_match_exct.out" 需要替换为当次注释的输出结果名称
df2 = pd.read_csv("Plants.NCBI2KEGG.txt", sep="\t", header=None, names=["col3", "col4"])
# 将两个DataFrame合并起来
df_merged = pd.merge(df1, df2, left_on="col2", right_on="col3")
# 将合并后的DataFrame的指定列输出到文件中
df_merged[["col1", "col2", "col4"]].to_csv("kplant_test_out.txt", sep="\t", index=False, header=False)

参考：
https://mp.weixin.qq.com/s/sqg5DGddKN9DPm8uTljW_g
https://cloud.tencent.com/developer/article/1991066
https://www.bioinfo-scrounger.com/archives/230/