Swiss-Prot注释脊椎动物基因组

wget -c ftp://ftp.uniprot.org/pub/databases/uniprot/current_release/knowledgebase/taxonomic_divisions/uniprot_sprot_vertebrates.dat.gz
wget -c ftp://ftp.uniprot.org/pub/databases/uniprot/current_release/knowledgebase/taxonomic_divisions/uniprot_trembl_vertebrates.dat.gz
zcat uniprot_sprot_vertebrates.dat.gz uniprot_trembl_vertebrates.dat.gz > uniprot_vertebrates.dat
awk '{if (/^ /) {gsub(/ /, ""); print} else if (/^AC/) print ">" $2}' uniprot_vertebrates.dat > uniprot_vertebrates.fasta
diamond makedb --in uniprot_vertebrates.fasta -d uniprot_vertebrates
diamond blastp -d uniprot_vertebrates.dmnd -q grass_carp.pep.fasta --evalue 1e-5 > blastp.outfmt6
python -m jcvi.formats.blast best -n 1 blastp.outfmt6
python add_annotation_from_dat.py blastp.outfmt6.best /data/database/UniProt-Plant/uniprot_plants.dat

github获取add_annotation_from_dat.py
之后会输出swiss_annotation.tsv, 输出信息包括如下几列

gene id
uniprot accession
identity
homology species
EnsemblPlants
GO ID
GO component, CC/MF/BP
evidence

twilight处理pangenome结果

github下载脚本(https://github.com/ghoresh11/twilight)
需要的R包data.table,ggplot2,optparse

Rscript classify_genes.R -p gene_presence_absence.Rtab -g groups.tab

必须参数:
-p, –presence_absence:Roary 或 Panaroo输出的gene_presence_absence.Rtab文件
-g, –grouping:一个制表符分隔的分组文件
可选参数:
-o, –out:输出目录名称(默认=“out”)。
-m, –min_size:忽略少于min_size基因组的组(默认 = 10)。
-c, –core_threshold:用于定义每个组内的核心基因的阈值(默认值 = 0.95)。
-r, –rare_threshold:用于定义每组内稀有基因的阈值(默认值 = 0.15)。
-h, –help:打印帮助信息并退出。
输出:
1、classification.tab 基因簇在每个分类中存在的形式(核心、辅助、和稀有)
2、frequencies.csv 基因簇在每个分离中的频率
3、plots 文件夹,包含柱状图和PCA图