基于SSR数据的群体结构分析和DeltaK值可视化

1.ssr做群体结构数据准备
1.1数据转换。

1.2原始数据是不符合structure的格式要求的，分享给大家我的python脚本做数据格式转换。

$cat dat2structure.py
#!/usr/bin/env python3

import sys

dat = sys.argv[1]

with open(dat,'r') as f:
	num=1
	for LN in f:
		if num ==1:
			print(LN.strip())
			num-=1
		else:
			spl = LN.strip().split()
			cnt = len(spl)
			#print(' '.join(spl[0:2]))
			print(' '.join(spl[0:3]),' '.join([spl[i] for i in range(3,cnt,2)]))
			print(' '.join(spl[0:3]),' '.join([spl[i] for i in range(4,cnt,2)]))

1.3转化后的数据格式，由于SSR数据是共显性标记，每两行为一个样本的SSR标记值：

2.structure 命令行计算群体结构
mainparams主要参数及说明。

extraparams使用默认即可。

for K in {1..10};do (nohup /biodata/02.software/structure/structure -i structure.new.txt -m mainparams -e /biodata/02.software/structure/extraparams -K ${K} -o kk_run1_k${K} &);done # 运行>=2次做后续DeltaK 折线图

3.利用structureHarvester做K值评估
我选择下载structureHarvester.py本地版命令行进行统计，也可在线分析，不用自己再另行绘图。

python /home/hang/software/structureHarvester-master/structureHarvester.py --dir=test/Results/ --out=./LK --evanno

4.统计绘制DeltaK折线图，找最佳K值。

Rscript /home/mmm/script/GWAS/structure_DeltaK.R evanno.txt notgu

下面是我的绘图代码

library(ggplot2)
args = commandArgs(T)
ipt = args[1] # evanno.txt
out = args[2]

dat = read.table(ipt)


func = 'DeltaK = mean(|L"(K)|) / sd(L(K))'
ggplot(dat,aes(x=V1,y=V7)) + geom_point(color='blue',size=rel(2)) + geom_line(color='blue') + 
  scale_x_continuous(breaks = 1:10,labels = c('',2:9,'')) +
  xlab('K') + ylab('Delta K') + 
  geom_hline(yintercept = 0,color="gray") +
  annotate("text", x=5, y=max(dat$V7,na.rm = T)*1.1, parse=F, label=func) +
  theme_classic()
ggsave(paste0(out,'_DeltaK.png'),width = 8,height = 6,dpi = 300)
ggsave(paste0(out,'_DeltaK.pdf'),width = 8,height = 6)

5.CLUMPP
使用CLUMPP对structure分析的重复运算结果进行重复抽样分析。得到最佳K值的Q-matrix结果。
在当前目录创建好 paramfile 配置文件（如下：红框内是你需要对应自己数据修改的参数）

生成popfile文件参考命令：

ls *_k8_f|while read i;do grep -A 206 "Inferred clusters" $i|grep -v "Inferred clusters"|awk '{print $1":",$6,$7,$8,$9,$10,$11,$12,$13,1}END{print ""}' >> notgu.popfile;done

运行下列命令即可得到最终K值 Q矩阵：
nohup /biodata/02.software/CLUMPP_Linux64.1.1.2/CLUMPP &
6.群体结构图可视化
可选择将CLUMPP的结果传递给distruct，进行structure图形的绘制。
这里我使用R自行绘制。（由于我的数据是真实项目数据，这里做了部分信息涂抹处理）