转载自链接：https://www.jianshu.com/p/3cfeb49c821a
古老WGD检测有两种方法，一种是共线性分析，另一种则是根据Ks分布图。其中Ks定义为平均每个同义位点上的同义置换数，与其对应的还有一个Ka，指的是平均每个非同义位点上的非同义置换数。
如果没有WGD或是大片段重复，那么基因组中的旁系同源基因的同义置换符合指数分布(exponential distribution), 反之，Ks分布图中就会出现一个由于WGD导致的正态分布峰(normal distributed peak). 而古老WGD的年龄则可通过分析这些峰中的同源置换数目来预测(Tiley et al., 2018)。
wgd是一个分析流程软件，通过调用各种工具软件来实现的。
WGD安装使用 参考地址 wgd githubs
wgd进行全基因组复制分析
conda create -n wgd python=3.7 diamond blast mcl muscle mafft prank paml fasttree cmake libpng mpi=1.0=mpich
conda activate wgd
i-adhore下载地址
安装i-adhore-3.0.01，
先用which conda查看你的conda的路径，修改prefix，或者是正常安装，之后添加环境变量也行。
tar -zxvf i-adhore-3.0.01.tar.gz
cd i-adhore-3.0.01
mkdir -p build && cd build
cmake .. -DCMAKE_INSTALL_PREFIX=~/miniconda3/envs/wgd/
make -j 4
make install
目前版本wgd要求python>=3.6
从github下载wgd的zip文件
安装wgd
unzip wgd.zip
cd wgd
cd wgd
pip install .
dmd:功能和mcl一样，新推出的功能，比对基于diamond，速度比blast快。
kde: 对Ks分布进行KDE拟合
ksd : Ks分布构建
mcl:All-vs-ALl的BLASP比对 + MCL分类分析.
mix: Ks分布的混合建模.
pre: 对CDS文件进行预处理
syn: 调用I-ADHoRe 3.0利用GFF文件进行共线性分析
viz: 绘制柱状图和密度图
wf1: 全基因组旁系同源基因组(paranome)的Ks标准分析流程，调用mcl, ksd和syn
wf2: one-vs-one 同源基因(ortholog)的Ks标准分析流程，调用wcl和ksd
WGD使用示例
输入文件准备
物种基因组的snail.cds.fa文件、snail.gff3文件
一步到位分析法
直接mcl+ksd+syn(单物种自身分析)
wgd wf1 -gff snail.gff3 -n 48 snail.cds.fa snail_wgd_out
参数：-n 指定进程数量48， snail_wgd_out是输出目录 gff3文件和cds文件
2个物种比较分析（wcl+ksd）
wgd wf2 -n 48 snail.fasta whale.fasta snail_vs_whale_ks_out
参数：没有做syn,所以不需要gff3.可以后续自己做syn
可视化
wgd kde wgd_ksd/sample.cds.fasta.ks.tsv #Ks直方图和拟合曲线(自带的作图算法，很丑)
getKaKs.R脚本，只需要修改工作目录和输入的tsv文件名。
Rscript getKaKs.R内容如下
setwd("./wgd_ksd")
library("ggplot2")
#读取wgd-ksd产出的文件
data_kaks <- read.table("snail.cds.fa.ks.tsv",header = TRUE,fill=TRUE,na.strings = "")
#删除包含NA的行
#data_kaks <- na.omit(data_kaks)
##Ks
p0 <- ggplot(data=data_kaks, aes(x=Ks)) +
geom_density(alpha=0.4,color="red")+ xlab('Ks value')+ylab('Density')+
xlim(0,15) #设置x的最大值
p0
ggsave('wgd_Ks.pdf',dpi=300)
##Ka
p1 <- ggplot(data=data_kaks, aes(x=Ka)) +
geom_density(alpha=0.4,color="green")+ xlab('Ka value')+ylab('Density')+xlim(0,2)
#labs(title = "Distribution of Ka distances", size = 1.5)+guides()
p1
ggsave('wgd_Ka.pdf',dpi=300)
##Ka/Ks
data_kaks <- na.omit(data_kaks)
data_kaks$Ka.Ks <- data_kaks$Ka/data_kaks$Ks
p2 <- ggplot(data=data_kaks, aes(x=Ka.Ks)) +
geom_density(alpha=0.4,color="orange")+ xlab('Ka/Ks value')+ylab('Density')+xlim(0,2)
p2
ggsave('wgd_KaKs.pdf',dpi=300)
##计算复制时间divergence date根据公式T(divergence date)=Ks/(2*r) r=7*10^-9
###r=7*10^-9 参考https://science.sciencemag.org/content/327/5961/92
r <- 7*10^-9
data_kaks$Time <- (data_kaks$Ks/(2*r))/10^6
p3 <- ggplot(data=data_kaks, aes(x=Time)) +
geom_density(alpha=0.4,color = "blue")+xlab('Million years ago (MYA)')+ylab('Density')+xlim(0,500)
p3
ggsave('Age distribution of wgd.pdf',dpi=300)
##拼图
p4 <- cowplot::plot_grid(p0,p1,p2,p3,nrow=2,labels = c("a","b","c","d")) #nrow 指定2行
p4
ggsave("wgd_kaks_mya.pdf",dpi=300)
后续还有共线性的可视化（待续）
wgd这个程序开发者，真是为用户操碎了心。脚本都不用写，一句话直接出结果。难点就是依赖软件的安装，安装依赖软件可能会报一堆错误。
分步分析法
查找基因组内部同源基因
wgd mcl -n 20 --cds --mcl -s sample.cds.fasta -o sample_cds.out
# -n: 线程
# --cds: 输入为cds
# --mcl: mcl聚类
输出结果在sample_cds.out，有两个结果输出
blast.tsv: BLASTP的outfmt6输出结果
blast.tsv.mcl: MCL聚类结果，每一行可以认为是一个基因家族(gene family)
使用wgd ksd构建Ks分布
wgd ksd -n 80 sample_cds.out/sample.cds.fasta.blast.tsv.mcl sample.cds.fasta
输出结果在默认在wgd_ksd目录下
ks.tsv: 每个基因家族中基因对的各项统计，其中包括Ka和Ks
ks.svg: Ks分布
使用wgd syn进行共线性分析
wgd syn -feature gene --gene_attribute ID \
-ks wgd_ksd/sample.cds.fasta.ks.tsv \
sample.gff3 sample_cds.out/sample.cds.fasta.blast.tsv.mcl
#--feature: 从第三列选择特征
#--gene_attribute: 从第九列提取编号
Ks的下游分析通常还包括对Ks分布的统计建模，这可以使用wgd kde进行核密度拟合(Kernel density estimate, KDE)或用wgd mix建立高斯混合模型(Gaussian mixture models)
wgd安装中报错解决方案
blast中makeblastdb 命令提示缺少libssl.so.1.0.0
error while loading shared libraries: libssl.so.1.0.0: cannot open shared object file: No such file or directory
解决办法：conda install -c bioconda rgi openssl=1.0
