基本概念
单体型:又称单倍体型或单元型,是个体组织中,完全遗传自父母双方中某一亲本的一组等位基因。拥有特定SNP的个体常常在附近某一特定变异位点拥有特定等位基因,这种关系叫做连锁不平衡(linkage disequilibrium , LD),同一染色体上的这一情况即为单体型。
简而言之,单体型就是指一段连锁的核酸序列。在单体型网络图中,不同的单体型就是用序列中的变异位点来区分。因此一般用来构建单体型的材料为:线粒体基因组、X/Y染色体,叶绿体基因组、或者是基因组上一段强连锁的区段。
A图中就是一个单体型网络,每个彩色圆代表一种单体型,圆的大小表示该单体型包含的样本数,两个圆之间的连线表示两个单体型彼此相关,连线上的短线表示从一个单体型变为另一个单体型需要经历的碱基替换个数。另外,每个灰色的圆点表示软件推断可能存在的单体型,但没有样品。
B图为单体型地图,同A图,每个彩色圆代表一种单体型,圆的大小表示该单体型包含的样本数。
**该示意图展示的单体型只在一个群体中出现,但大部分情况一个单体型会出现于多个群体中,因此,多数文章的单体型网络图中,各个单体型圆会用不同颜色填充,也就是常见的饼图,一种颜色表示一个群体(如按地理分布和亚种分布等)。
用PopART构建单体型网络和单体型地图
软件下载
直接官网下载就好,PopART。
数据准备
PopART的输入文件是.nex格式文件,该文件主要包含两个部分内容:DATA和TRAITS,也就是序列信息和分组信息。
DATA准备
重中之重,当然是准备我们要用来建立单体型网络图的序列文件,比如下图文章中用的隐孢子虫gp60基因序列,画图前要先把所有样本的该序列都做alignment(ClustalW、Muscle等等软件都可以)。
拿到序列文件后,就开始构建.nex文件中的DATA部分,把这个复杂的步骤交给DnaSP
打开DnaSP->载入比对好的序列文件.fas->Generate->Haplotype file->储存为.nex文件
TRAITS准备
这个部分如字面意思,告诉软件数据的分组信息,直接人工添加到上一步生成的.nex文件的最下方(用记事本打开.nex文件编辑)。
格式如下图所示:
红框所示部分为样本总分组数;黄框所示部分为样本分组名称;蓝框所示部分为单体型名称和每个分组在该单体型中的样本数量。(注意:分隔符号是Tab键)
画单体型网络图
打开PopART->导入上述.nex文件->选择Network->TCSnetwork
大功告成,可以通过菜单栏各个选项调整群体颜色、字体、比例等。
画单体型地图
选择View->Switch to map view
补充
上面提到,如果我们用来分析的样品不是某基因,而是基因组SNP序列,那么,我们需要找到强连锁区。
#用plink过滤低频位点
plink --vcf xxxx.vcf --maf 0.05 --geno 0.2 --recode vcf-iid -out xxxx-maf0.05 --allow-extra-chr
#用plink筛选连锁区
plink --vcf xxxx-maf0.05.vcf --indep-pairwise 100 50 0.2 -out xxxx-maf0.05-LD --allow-extra-chr --make-bed
#用plink找到强连锁区
plink --noweb --bfile xxxx-maf0.05-LD --blocks no-pheno-req --allow-extra-chr
强连锁区的结果储存在.blocks.det文件中,该文件第一列为染色体/contig名称;第二列为连锁区的开始位置;第三列为连锁区结束位置。找到感兴趣的区域,用beagle软件提出该区域构建DATA文件。
java -Xss5m -Xmx4g -jar YourWayToBeagle/beagle.jar gt=xxxx-maf0.05.vcf out=phased chrom=Chr:start-end
#VCF转换phy
vcf2phylip.py -i phased.vcf.gz
#有了phy文件拿到Win操作系统下,用Bioedit就可以变成fasta格式
上述VCF转换fasta脚本来源:vcf2phylip.py
参考:
单体型网络的构建和展示