structure-2.3---中文使用手册


2023年12月14日发(作者:天天斗地主)

Structure 2.3中文使用手册

Jonathan K. Pritcharda

Xiaoquan Wena

Daniel Falushb 1 2 3

a芝加哥大学人类遗传学系

b牛津大学统计学系

软件来自

/

2010年2月2日

1我们在Structure项目中的其他的同事有Peter Donnelly、Matthew Stephens和Melissa

Hubisz。

2开发这个程序的第一版时作者(JP、MS、PD)在牛津大学统计系。

3关于Structure的讨论和问题请发给在线的论坛上:***********************************。 在邮递问题之前请查对这个文档并搜索以前的讨论。

1 引言

程序Structure使用由不连锁的标记组成的基因型数据实施基于模型的聚类方法来推断体结构。这种方法由普里查德(Pritchard)、斯蒂芬斯(Stephens)和唐纳利(Donnelly)(2000a)在一篇文章中引入,由Falush、斯蒂芬斯(Stephens)和普里查德(Pritchard)(2003a,2007)在续篇中进行了扩展。我们的方法的应用包括证明体结构的存在,鉴定不同的遗传体,把个体归到体,以及鉴定移居者和掺和的个体。

简言之,我们假定有K个体(这里K可能是未知的)的一个模型,每个体在每个位点上由一组等位基因频率来刻画。样本内的个体被(按照概率)分配到体,或共同分配到两个或更多个体,如果它们的基因型表明它们是混和的。假定在体内,位点处于哈迪-温伯格平衡和连锁平衡。不精确地讲,个体被按达到这一点那样的方法指定到体。

我们的模型不假定一个特别的突变过程,并且它可以应用于大多数通常使用的遗传标记,包括微卫星(microsatellites)、SNP和RFLP。模型假定在亚体内标记不处于连锁不平衡(LD),因此我们不能处理极其靠近的标记。从2.0版开始,我们现在能够处理弱连锁的标记。

虽然这里实现的计算方法是相当强有力的,但是为了保证明智的答案,在运行程序的过程中还是需要谨慎。例如,不可能从理论上确定合适的运行长度(时间),这需要用户自己做一些实验。这份资料描述软件的使用和解释,并补充发表的文章,这些文章提供了对方法的更正式的描述和评价。

1.1 概述

软件包Structure由几个部分组成。程序的计算部分用C语言编写。我们发布源码和用于各种平台(目前有苹果机,Windows,Linux,Sun)的可执行文件。C可执行文件读取用户提供的一个数据文件。还有一个Java前端为用户提供各种有帮助的工具,包括对输出的简单的处理。你也可以从命令行调用Structure而不是使用前端。

这份资料包括关于怎样格式化数据文件、怎样选择合适的模型、以及怎样解释结果的信息。它也有关于使用两种界面(命令行和前端)的细节以及各种用户定义的参数的汇总。

1.2 在2.3版中有哪些更新?

2.3版(2009年4月发布)引入了新的模型用于改进数据集结构的推论,其中(1)数据对于通常的结构模型来说信息不够,不足以提供准确的推论,但是(2)抽样的地点与体归属关系(population membership)相关。在这种情形下,通过明确利用抽样地点信息,我们使结构得到改善,经常允许性能提高很多(Hubisz et al., 2009)。我们希望在下几个月释放更进一步的改进。

表1:实例数据文件。这里MARKERNAMES = 1, LABEL = 1, POPDATA = 1, NUMINDS = 7,

NUMLOCI = 5, MISSING = -9, POPFLAG = 0, LOCDATA = 0, PHENOTYPE = 0,

EXTRACOLS = 0。第2列显示个体的地理取样位置。我们也可以把数据存储为每个个体一行(ONEROWPERIND = 1),在这种情况下第一行为“George 1 -9 -9 145 -9 66 64 0 0 92 94”。

乔治

乔治

保拉

保拉

马修

马修

鲍勃

鲍勃

Anja

Anja

彼得

彼得

卡斯坦

卡斯坦

1

1

1

1

2

2

2

2

1

1

1

1

2

2

Loc_a

-9

-9

106

106

110

110

108

-9

112

114

-9

110

108

110

Loc_b

145

-9

142

148

145

148

142

142

142

142

145

145

145

145

Loc_c

66

64

68

64

-9

66

64

-9

-9

66

66

-9

62

64

Loc_d

0

0

1

0

0

1

1

0

1

1

0

1

0

1

Loc_e

92

94

92

94

92

-9

94

94

-9

94

-9

-9

-9

92

2 数据文件的格式

基因型数据的格式显示在表2中(表1显示一个例子)。基本上,整个数据集被作为一个矩阵安排在单个文件里,其中个体的数据在行里,位点在列里。用户能对格式做出若干选择,大多数这些数据(除基因型外!)是可选择的。

对于一个二倍体生物,每个个体的数据可以是作为连续的2行被储存,其中每个位点在一列,或者在一行中,其中每个位点在连续的两列。除非你打算使用连锁模型(见下面),否则单个个体的等位基因的次序并不重要。预基因型(pre-genotype)数据列(见下面)对每个体记录两次。(更一般地,对于n倍体生物来说,每个个体的数据被储存在n个连续的行中,除非ONEROWPERIND选项被使用。)

2.1 数据文件的组成部分:

输入文件的要素如下所列。如果给出,它们一定按以下顺序,然而大多数是可选的并且可以被完全删除。用户必须指明哪些数据被给出,或者在前端里(front end),或者(当从命令行运行Structure时)在一个单独的文件mainparams里。同时,用户也要指定个体和位点的数目。

2.2 行

1. 标记名称(可选择;字符串) 文件的第一行可以包含数据集里的每个标记的标识符的一个列表。这一行包含整数或字母的L个字符串,其中L是位点的数目。

2. 隐性等位基因(仅用于有显性的标记数据;整数)SNP或者微卫星数据一般将不包括这一行。但是如果选项RECESSIVEALLELES被设置为1,则程序要求有这一行来表明每个标记上哪个等位基因(如果有的话)是隐性的。关于更多的信息请参阅第4.1节。该选项用于象AFLP那样的数据,以及用于多倍体的情形,其中基因型可能是含糊的。

3. 标记之间的距离(可选择;实数)文件里的下一行是一个标记之间距离的集合,供有连锁的位点使用。这些应该是遗传距离(例如,厘摩),或者是这种距离的一些替代,基于(例如)物理距离。如果标记距离(粗略地)与重组率成正比,则距离

的实际单位不是那么重要 。前端从数据估计一个合适的尺度,但是命令行版本的用户必须在文件extraparams里设置LOG10RMIN、LOG10RMAX和LOG10RSTART。标记必须按照连锁中的图谱次序排列。当连续的标记来自不同的连锁(例如,不同的染体)时,这应该用数值-1注明。第一个标记也被赋值为-1。所有其他的距离都是非负的。这一行包含L个实数。

4. 连锁相信息(可选择;仅用于二倍体数据;在范围[0, 1]内的实数)。这只供连锁模型使用。这是L个概率的一行,出现在每个个体的基因型数据之后。如果连锁相是完全知道的,或者没有连锁相信息可用,则这些行是不必要的。当有来自家系数据的部分连锁相信息,或者当来自雄性的单倍体X染体数据和二倍体常染体数据被一起输入时,它们可能是有用的。对于连锁相信息有两种可选择的表示:(1)个体的两行数据被假设为分别与父本的和母本的相对应。连锁相行表明当前标记上的排序正确的概率(设置MARKOVPHASE = 0);(2)连锁相行表明与以前的等位基因有关的一个等位基因的连锁相是正确的概率(设置MARKOVPHASE = 1)。第一项应该填入0.5,以便把这行填写到L项。例如下列数据输入表示来自一个男性的信息,有5个连锁相未知的常染体微卫星位点,后面是3个X染体位点,使用母本/父本相模型:

102 156 165 101 143 105 104 101

100 148 163 101 143 -9 -9 -9

0.5 0.5 0.5 0.5 0.5 1.0 1.0 1.0

其中-9表示“缺失数据”,这里缺失是由第二X染体缺乏造成的,0.5表明常染体位点的连锁相是未知的,1.0表明X染体位点由母本遗传的概率为1.0,因此其连锁相是已知的。相同的信息可以用markovphase模型来描述。这样的话输入文件将读为:

102 156 165 101 143 105 104 101

100 148 163 101 143 -9 -9 -9

0.5 0.5 0.5 0.5 0.5 0.5 1.0 1.0

这里,2 1.0 s 表明那个第1 和第2,其次和第3 个X染体位点彼此完全同相。注意站点以站点产量在这些2 模式下将不同。在第一例子中,Structure将输出母亲和父亲染体的任务可能发生的事件。在第2 个情况下,它将输出在输入文件里列举的每等位基因的可能发生的事件。

5. 个体/ 基因型数据(必需的)取样的每一个个体的数据象在下面描述的那样安排成一行或多行。

2.3 个体/基因型数据

个体数据的每一行包含下列要素。这些形成数据文件里的列。

1. Label(标签)(可选择;字符串) 一串整数或者字母,用来指明样本中的每个个体。

2. PopData(可选择;整数)一个整数,指明一个用户定义的体,从其中获得个体(例如这些整数可以指明个体取样的地理位置)。在默认的模型中,这个信息不被聚类算法使用,但是能用来帮助组织输出(例如,将来自相同的预定义体的个体彼此紧挨着绘图)。

3. PopFlag(可选择;0或者1)一个布尔标签,表明使用学习样本时是否使用PopData(见USEPOPINFO,在下面)。(注:布尔(Boolean)变量(标签)是取值为TRUE或FALSE的变量,在这里分别用整数1(使用PopData)和0(不使用PopData)表示。)

4. LocData(可选择;整数)一个整数,为每个个体指明一个用户定义的取样地点(或者其他特性,例如一个分享的表现型)。当LOCPRIOR模型被打开时,这个信息用来帮助聚类。如果你仅仅希望使用LOCPRIOR模型的PopData,那么你可以省略LocData列,并设置LOCISPOP = 1(这告诉程序使用PopData来设置地点)。

5. Phenotype(可选择;整数) 一个整数,为每个个体指明一个所关心的表现型的值(表中的(i))。(表现型信息实际上没有用于Structure。这里用来允许与关联作图程序STRAT有一个平滑的接口。)

6. Extra Columns(可选择;字符串) 用户把被程序忽略的附加数据包括在输入文件里可能是方便的。这些数据就在这里输入,可以是由整数或字符组成的串。

7. Genotype Data(必需的;整数) 一个给定位点上的每个等位基因应该由一个独特的整数来编码(例如微卫星重复得分)。

2.4 缺失的基因型数据

缺失数据应该用没在数据中的其他地方出现过的一个数字来标明(按照惯例经常使用

-9)。这个数字也可以用于有单倍体和二倍体数据混合的地方(例如男性中的X和常染体位点)。缺失数据值是与描述数据集特性的其它参数一起被设置的。

2.5 格式化的错误。

我们已经进行了相当仔细的错误检查,以保证数据集的格式正确,并且程序将试图提供一些关于存在的任何问题的性质的提示。前端要求在每行的结束回车,不允许在行内回车;Structure的命令行版本以与处理空格或制表符(Tab)同样的方式处理回车。

可能出现的一个问题是,在将数据导入Structure之前用来组装数据的编辑程序可能引入隐藏的格式化字符,经常在行的末尾,或者在文件的末尾。前端能自动除去大多数这些错误,但是当数据文件好像处于正确的格式时,这类问题可能对错误负责。如果你正在把数据导入到一个Unix系统,dos2unix功能可能对彻底清理这些错误有帮助。

3 用户的建模决策

3.1 祖先模型

个体的祖先有4个主要模型:(1) 非混合模型(个体离散地来自一个体或者另一个体);(2)混合模型(每个个体从K个体中的每一个抽取他/她的基因组的一部分);(3)连锁模型(象混合模型一样,但是连锁的位点更可能来自相同的体);(4)有先验信息的模型(允许Structure使用关于取样地点的信息:或者帮助用弱的数据进行的聚类,发现迁移者,或者预定义一些体)。关于模型1、2 、4的详情见Pritchard等(2000a)和Hubisz

等(2009),关于模型3的详情见Falush等(2003a)。

1. 非混合模型。每个体完全来自K个体之一。输出报告个体i来自体k的后验概率。每个体的先验概率是1 / K。这个模型适合于研究完全离散的体,并且经常比混合模型在检测微妙的结构方面更强有力。

2. 混合模型。个体可能具有混合的祖先。这可以表述为个体i从体k中的祖先那里继承了他的/她的基因组的一部分。输出记录这些比例的后验平均估计值。以祖先向量q(i)为条件,每个等位基因的起源是独立的。

我们推荐这个模型作为大多数分析的起始点。这是处理真实体的大多数复杂性的一个相当灵活的模型。混合是真实数据的一个普通特征,如果你使用非混合模型,你或许不会发现它。混合模型也能以一种自然的方式处理混合的区域(hybrid zones)。

表2:数据文件的格式,为两行的格式。大多数这些组成部分是可选的(欲了解详细信息,参见正文)。Ml是标记l的标识符。rl表明哪个等位基因,如果有的话,在每个标记上是隐性的(仅针对显性的基因型数据)。Di,i+1是标记i和i + 1之间的距离。ID(i)是个体i的标签,g(i)是个体i的一个预先定义的体索引(PopData);f(i)是一个被用来合并学习样品的标签(PopFlag);l(i)是个体i的取样地点(LocData);(i)可以储存个体i的表现型;y1(i), ..., yn(i)用于储存额外的数据(这些数据会被程序忽略);(xli,1, xli,2)储存个体i在位点l上的基因型。pi(l)是个体i中的标记l的连锁相的信息。

3. 连锁模型。这实质上是将混合模型推广,来处理“混合连锁不平衡”,即,在最近混和的体中的连锁标记之间出现的相关性。Falush等(2003a)描述了该模型和更详细的计算。

基本的模型是,过去的t个世代,有一次混合事件,将K个体混合了。如果你考虑

单个染体,它由一系列“块(chunk)”组成,这些“块”是从混合时的祖先那里作为离散的单位遗传来的。出现混合LD是因为连锁的等位基因经常在相同的块上,因此来自相同的祖先体。

块的大小被假设为独立的指数随机变量,具有平均长度1/t(以摩尔根为单位)。在实践中我们估计“重组率”r,所用的数据对应于从现在的块切换到新的块的比率。个体i里的每个块以概率qk(i)独立地来自体k,其中qk(i)是那个个体的祖先来自体k的比例。

总起来,新模型保留了混合模型的主要要素,但是在单个块上的全部等位基因必须来自相同的体。新的MCMC算法结合了可能的块大小和断点。它对于每个体报告总的祖先,考虑连锁,并且也能报告染体的每一点儿的起源的可能性,如果用户想要的话。

当使用连锁的位点来研究混合的体时,这个新模型表现得比原先的混合模型更好。它得到对祖先向量的更准确的估计,并且能从数据中抽出更多的信息。这对混合作图应该是有用的。该模型不是用于处理非常紧密连锁的标记之间的背景LD的。

显然,这个模型是大多数混合体的复杂现实的大大的简化。不过,混合的主要的效应是在连锁的标记之间建立长远的相关性,因此我们这里的目的是在一个相当简单的模型中将那个特征包括进来。

计算比混合模型的要慢一点,特别对于大的K和不知道连锁相的数据。不过,它们对于数千个位点和个体以及多个体来说还是切实可行的。如果有关于标记的相对位置的信息(通常是一张遗传图谱),则只能使用该模型。

4. 使用先验的体信息。Structure的默认模式只使用遗传学的信息来了解体结构。不过,经常有可以与聚类相关的附加信息(例如,取样的个体的物理特性或者取样的地理位置)。目前,Structure可以用3种方式使用这种信息:

• LOCPRIOR模型:利用取样位置作为先验信息来辅助聚类——用于结构信号比较弱的数据集。有一些数据集,其中有真实的体结构(例如,取样位置之间的显著的FST),但是信号太弱,标准的Structure模型不能发现。对于标记很少、个体很少或者非常弱的Structure,经常是这样的情况。

在这种情形下,为了提高性能,Hubisz等(2009)发展了新模型,利用地点信息来辅助聚类。对于这样的数据集,其中结构的信号太弱以致使用标准的Structure模型不能被发现,新模型经常能提供体结构和个体祖先的准确的推断。

简言之,LOCPRIOR模型的基本原理如下。通常,Structure假定个体的所有部分都大约是先验等可能的。因为可能的部分的数目非常巨大,对于Structure来说,需要信息非常丰富的数据来断定个体的任何特定的部分被聚类到具有强的统计支持。相反,LOCPRIOR模型认为实际上,来自相同的取样位置的个体经常来自相同的体。因此,建立LOCPRIOR模型以期望取样的位置可能关于祖先是信息丰富的。如果数据表明位置是信息丰富的,那么LOCPRIOR模型允许Structure使用这种信息。

Hubisz等(2009)发展了一对LOCPRIOR模型:一种用于没有混合的情况,一种用于有混合的情况。在两种情况中,内在的模型(以及似然函数)与标准版本相同。关键的差别是允许structure使用地点信息来帮助聚类(即,通过修改先验信息来得到与位置有关的更偏爱的聚类解决方案)。

LOCPRIOR模型具有合乎需要的特性:(i)当不存在结构时,它们不倾向于发现结构;(ii)当个体的祖先与取样位置不相关时,他们能够忽视取样的信息;(iii)当体结构的信号非常强大时,旧模型和新模型基本上给出相同的答案。因此,我们建议在大多数数据数量非常有限的情形下使用新模型,特别是当标准的Structure模型不提供一个Structure的清晰信号时。但是,因为现在已经积累了标准的Structure模型的很多经验,我们建议对于信息非常丰富的数据集将基本模型作为默认(Hubisz

等等,2009)。

为了运行LOCPRIOR模型,用户必须首先为每个个体指定“取样地点”,作为一个整数编码。即,我们假定样品是在一组分离的位置收集的,并且我们不使用关于地点的任何空间信息。(我们认识到,在一些研究中,每个个体可能在一个不同的地点收集,因此将个体塞进一套更小的分离的地点可能不是对数据的理想的代表。)

“地点”也可以代表一个表现型、生态型(ecotype)或者民族团体(ethnic group)。

地点被键入到输入文件中,要么在PopData列(设置LOCISPOP = 1)中,要么作为一个单独的LocData列(参阅第2.3节)。为了使用LOCPRIOR 模型,你必须首先指定或者用混合模型用非混合的模型。如果你使用的是图形用户界面版本,则勾选“use sampling locations as prio”(用取样位置作为先验信息)框。如果你使用的是命令行版本,则设置LOCPRIOR = 1。(注意,LOCPRIOR与连锁模型不兼容。)

我们迄今的经验是当不存在结构时,LOCPRIOR模型不偏向于检测到假的结构。你可以把相同的诊断用于是否有真的结构,当你没使用LOCPRIOR时。另外查看r的值可能有帮助,它确定由位置携带的信息的数量。r的值接近1,或者<1,表明位

置是信息丰富的。r的更大的值表明或者没有体结构,或者结构不依赖位置。

USEPOPINFO模型:使用取样位置来对移居者或者杂交种进行检验——供信息非常丰富的数据数据集使用。在一些数据集里,用户可能发现预确定的组(例如取样位置)几乎正好与结构聚类相对应,除了少数似乎被错误归类的个体以外。Pritchard等(2000a)提出了正式的Bayesian检验,用于评价是否在这个样品内的任何个体是他们认为的体的移民,或者具有新近的移民祖先。

注意这个模型假定被预先规定的体通常是正确的。它采用十分强大的数据来克服先验的错误分类。在使用USEPOPINFO模型之前,你也应该在没有体信息的情况下运行程序,以保证预确定的体与遗传学的信息粗略一致。

为了使用这模型,把USEPOPINFO设置为1,并且选择MIGRPRIOR的一个值(在Pritchard等(2000a)中它是v)。你可以在0.001到0.1的范围内为v选择一个值。

每个个体的预确定的体被设置在输入数据文件中(见PopData)。用这种方式,在输入文件里被分配到体k的个体在Structure算法中将被分配到k。因此,被预先规定的体应该是在1和MAXPOPS (K)之间的整数。如果任何个体的PopData超出这个范围,它们的q将按正常的方式被更新(即没有先验的体信息,根据将被使用的模型,如果USEPOPINFO被关上的话。)

USEPOPINFO模型:预先指定一些个体的起源的体来帮助未知起源的个体的祖先估计。使用USEPOPINFO模型的第二个方法是定义“学习样本”(learning samples),它被预定义为来自特定的。然后用Structure来聚类剩下的个体。注意:在前端里,这个选项使用“Update allele frequencies using only individuals with POPFLAG=1”选项被打开,位于“Advanced Tab”标签下。

学习样品是利用数据文件里的PopFlag列实现的。预先规定的体被用于那些个体,它们的PopFlag = 1(并且它们的PopData在(1...K)中)。对于PopFlag = 0的个体,PopData值被忽略。如果数据文件里没有PopFlag列,那么当USEPOPINFO被开启时,PopFlag被为全部个体设置为1。具有PopFlag = 0的或者PopData不在(1...K)中的个体的祖先,根据混合或者没有混合的模型被更新,象由用户指定的那样。如上所述,如果有很少的个体没有预先规定的体,将设置成一个明智的值来可能是有帮助的。

USEPOPINFO的应用可能在几个方面有帮助。例如,可能有一些个体的来源是已知的,我们的目标是对未知来源的另外的个体进行归类。例如,我们可能从一

已知品种(编号为1 . . .K)的狗中收集数据,然后使用Structure为未知的(也许是杂交种)起源的另外的狗估计祖先。通过预先设置体数目,我们可以保证Structure聚类对应于预先确定的品种,这使输出更可解释,并且能改进推论的准确性。(当然,如果两个预先确定的品种在遗传上是相同的,那么未知起源的狗可能被推断为具有混合的祖先。

USEPOPINFO的另一种用途是用于这样一种情况:用户想要只使用个体的一个子集来更新等位基因频率。通常,Structure分析使用全部可得到的个体来更新等位基因频率估计值。但是有一些情况,在那里你可能想对于一些个体估计祖先,没有那些个体会影响等位基因频率的估计。例如你可以有学习样品的一个标准的收集,然后周期性地你想要为新的一批基因型化的个体估计祖先。使用默认的选项,个体的祖先估计(稍微)取决于它们所在的批次。通过使用PFROMPOPFLAGONLY,你可以保证等位基因频率估计值只依赖于PopFlag = 1的那些样品。在不同的情况下,Murgia等(2006)想要确定一套无性系的狗瘤的起源。那些瘤如此紧密有关以至于使用的缺省设置时瘤形成它们自己的一类。通过使用PFROMPOPFLAGONLY,Murgia等迫使瘤与其他canid聚类分在一组。

意见:我们建议首先运行Structure的基本的版本,以便证实被预先规定的标签确实的确符合实际的遗传学体。其次,当使用学习样品时,通过设置比0大的MIGRPRIOR来允许一些错误的分类可能是明智的。

3.2 等位基因频率模型

对于等位基因频率有两个基本的模型。一个模型假定每个体内的等位基因频率是独立的,从一个分布中抽取,这个分布由参数指定。那是用于Pritchard等(2000a)种的原先的模型。通常我们设置 = 1;这是缺省设置。

Falush等(2003a)实施了一个模型,具有相关的等位基因频率。这个模型标明不同体中的频率很可能是相似的(或许由于迁移或者由于共有的祖先)。更详细的资料如下。

独立的模型对于很多数据集表现不错。粗略地说,这最先说我们期望在不同的体中的等位基因频率彼此不同。相关的频率模型说它们实际上可能十分相似。对于亲缘关系近的体,这经常改进聚类,但是可能增加过高估计的K的危险(如下)。如果一个体与其他体分歧较大,则当那个体被除去时,相关的模型有时可以取得更好的推论。

估计: 固定 = 1对于大多数数据是一个好主意,但是在一些情况下,例如SNP数据,

其中大多数次要的等位基因是稀少的,这时候较小的数值可能工作得更好。对于这个原因,你可以让程序为你的数据估计。你可能想要这样做一次,或许对于K = 1来说,然后将固定在被估计的值上,因为在试图同时国际太多的假设参数(, , F)时对于非识别性(non-identifiability)好像有一些问题。

相关的等位基因频率模型: 如同Falush等(2003a)描述的那样,相关的频率模型使用一个(多维的)矢量,PA,它记录假设的“祖先”体中的等位基因频率。假定在我们的样品中代表的K个体每个都已经经历过与这些祖先频率的独立的漂移,速率分别用参数F1, F2,

F3, ..., FK表示。除归因于有点不同的模型的差别和估计的差别外,被估计的Fk值应该数量上类似于FST值。此外,对于具有许多混合的数据要准确地估计Fk很难。

PA被假设为具有Dirichlet先验,具有与上面的体频率使用的相同的形式:

pAl· ~ D(1, 2, . . . ,

Jl), (1)

对每个l独立。然后,体k中的频率的先验为

pkl.~D(PAl11Fk1Fk,PAl2,FkFk,PAlJl1Fk), (2)

Fk对每个k和l独立。在这个模型里,F与遗传学距离FST有密切的关系。按照FST的标准的参数化方法,每个体中的期望频率由总的平均频率给出,当等位基因的总频率为p时,跨越亚体的频率的方差为p(1 − p)FST。这里的模型几乎一样,除了我们对模型稍微做了推广以外,通过允许每个体以一个不通的速率(Fk)漂离祖先体,如同体具有不同的大小时可能被期望的那样。我们也试图估计“祖先频率”,而不是使用平均的频率。

我们将独立的先验(prior)放于Fk上,与平均数为0.01、标准差为0.05的γ分布成正比(但是有Pr[Fk  1] = 0)。γ先验分布的参数可以由用户修改。一些实验表明,0.01的先验平均值对应于非常低细分的水平,对于独立频率模型的数据经常导致好的表现。在其他的问题中(其中体之间的差别更加明显),好像数据通常压倒了这个Fk的先验。

3.3 程序要运行多长时间?

程序从一个随机的配置启动,从那里采取一系列步骤穿过参数空间,每个步骤(只)依赖于前一个步骤的参数值。这个程序在运行期间引起不同的点上的Markov链的状态之间的相关性。希望是通过运转模拟足够久,相关性将可以被忽视。

有两个问题要担心:(1) burnin长度:在收集数据使启动配置的影响减到最小之前模

拟要运行多久,(2)在burnin以得到准确的参数估计之后模拟要运行多久。

要选择合适的burnin长度,看看由这个程序打印的归纳统计量的值是真的有帮助的(例如(,F,在体之间的分歧距离Di,j,以及似然),以便了解它们是否已经收敛。通常10000——100000的burnin非常足够了。

要选择适当的运行长度,你需要在每个K上做几次运行,也许长度不同,看看你是否得到一致的答案。通常,利用10000——100000步运行你能得到参数(P和Q)的好的估计,但是Pr(X|K)的准确的估计可能需要更长时间的运行。实际上,你的运行时间的长度可能决定于你的计算机速度和耐心。如果你正处理极其大的数据集,并且被运行时间阻止,你可以试着修剪运行的长度和标记/个体的数量,至少为探索的分析。

前端提供了几个主要参数的时间序列曲线。在burnin阶段结束之前你应该看看这些曲线,以便了解这些曲线是否看起来达到了平衡。如果在burnin阶段结束时数值仍然在增加或者减少,你需要增加burnin长度。

如果在整个运行期间(即,不只是在burnin期间)的估计值变化非常大,你可以通过增大ALPHAPROPSD来得到对Pr(X|K)的更准确的估计,这改进了在那种形势下的混合。(见在第5节的一个有关的问题)。

4 缺失数据,无效的等位基因和显性标记

当不断改进Q和P时,程序忽略缺失的基因型数据。当在一个特别的位点有漏缺数据的可能性与个体在那里有什么等位基因无关时,这种方法是正确的。当具有漏缺数据的个体的Q的估计不那么准确时,没有特别的原因阻止这样的个体参加分析,除非他们根本几乎没有数据。

当数据以系统的方式遗漏时,出现一个严重的问题,如同用无效等位基因那样。这些不适合假设的模型,即使没有体结构,也能够导致明显的违背哈迪-温伯格。人们不应该期望假设的模型对这类破坏是稳健的。但是如果无效的等位基因可能是一个重要的问题的话,则显性标记模型(下面)可以被使用。

在样本中有多名家庭成员也会破坏模型假定。这有时会导致K的过高估计,特别对于相关的频率模型(Falush等,2003a),但是当K固定时,这对将个体分配给体的影响很小。

4.1 显性标记、无效等位基因和多倍体基因型

对一些类型的遗传学标记(例如AFLP)来说,区分全部基因型是不可能的。其它类型的标记可能导致模棱两可的基因型,如果由于附近序列的变化导致PCR产物不能扩增,一部分等位基因为“无效”。从2.2版开始,我们实现了一个模型,处理与显性标记相关的基因型的模糊性。

总之,我们假定在任何特定的位点可能有对全部其他等位基因(例如A)为隐性的单个的等位基因,而全部其他的标记是共显性的。因此AB和BB将作为“表现型”B出现在未加工的基因型数据中,AC和CC将被记录为C,而BC将被记录为BC。当有模糊性时,模型在可能的基因型上求和。全部的细节在Falush等(2007)里给出。

为了执行这些计算,必须告诉算法每个位点上的哪个等位基因(如果有的话)是隐性的。这通过设置RECESSIVEALLELES=1来进行,并且在输入文件顶上包括一行单L整数,在标记名称和图谱距离的(可选的)行之间,表明在数据集里的L个位点的每个上面的隐性等位基因。如果一个给定的位点上的全部标记是共显性的,那么那个位点上的隐性值必须被调整成MISSING(缺失的)数据值。相反,如果隐性等位基因从未在纯合状态被观察到,但是你认为它可能存在(例如可能有无效的等位基因),那么就把隐性值设置成在那个位点没被观察到的等位基因(而不是MISSING!) .

编码基因型数据:如果表现型是不含糊的,那么它被在Structure输入文件里按照它本来的样子编码。如果它是含糊的,那么它被作为显性等位基因的纯合体编码。例如,表现型A

被编码为AA,B被编码为BB,BC被编码为BC,等等。如果标记是其他方面为二倍体的一个个体中的单倍体(例如男性中的X染体),那么第2个等位基因被象以前一样编码为MISSING(缺失)。当A是隐性的时,基因型AB、AC等等在输入文件里是不合法的。

当RECESSIVEALLELES被用来处理无效的等位基因时,看起来是无效的纯合体(homozygote null)的基因型应该作为隐性等位基因的纯合体而不是作为缺失数据被输入。在实践中可能不确定是否一个失败的基因型真的归因于纯合的无效等位基因。Structure应该对这些编码为缺失的数据是稳健的,除非无效等位基因在一个位点上的频率很高。

在多倍体(PLOIDY>2)中形势更复杂,因为甚至对共显性标记都可能有基因型的含糊。在杂合体中准确地识别出基因型经常是困难的。例如在三倍体中,表现型AB可能是AAB或者ABB。如果Structure在RECESSIVEALLELES=0的条件下运行,那么就假定没有含糊。

对于多倍体,当RECESSIVEALLELES=1时,Structure允许数据包含具有基因型模糊和

不具有基因型模糊的位点。如果一些位点不含糊那么设置代码NOTAMBIGUOUS为一个整数,这个整数不与数据内的的任何等位基因相匹配,并且不等于MISSING(缺失)。然后在输入文件顶上的隐性等位基因的行里为不含糊的位点放置NOTAMBIGUOUS代码。如果不是那样,而是在一个特定的位点上等位基因全部是共显性的,但是有关于每个的数目(例如为在四倍体里的微卫星)含糊,那么就把隐性等位基因代码设置为MISSING。最后,如果有隐性等位基因,并且还有关于每个等位基因的数目的含糊性,则设置隐性等位基因代码来表明哪个等位基因是隐性的。存在拷贝数含糊性的等位基因的编码与存在显性标记的那些相似。因此,举例来说,在四倍体中,观察到3个共显性位点B、C和D,这应该被编码为B C

D D或者等效地B B C D或者任何包括3个等位基因中的每一个的其他组合。它不应该被编码为B C D (MISSING),因为这表明该特定的个体在所指的位点是三倍体。如果在这个位点上存在一个隐性等位基因A,它也不能被编码为B C D A。

Pr(K)的估计: 当RECESSIVEALLELES被用于二倍体时,Markov链的每个步骤上的似然值是通过在可能的基因型上求和来计算的。为了便于编码,当要么PLOIDY>2要么使用了连锁模型时,我们以当前推算的(imputed)基因型为条件。这减小似然值,并且好像大大地扩大似然值的方差。有限的经验表明在后一种情况下这导致对K的估计效果变差,你应该把K的这种估计看做是不可靠的。

5 K(体数目)的估计

在描述这个程序的我们的文章里,我们指出这个问题应该被小心对待,由于两个原因:(1)要获得对Pr(X|K)的准确估计在计算上是困难的,我们的方法仅仅提供一个专门的(ad

hoc)近似;(2)K的生物学解释可能不是简单的。

在我们的经验里我们发现真正的困难在于第2个问题。我们的用于估计K的程序一般在具有少量离散的体的数据集中计算效果较好。不过,很多现实世界的数据集并不准确地符合Structure模型(例如,由于通过距离或者近交而产生的隔离)。在那些情况里对于什么是K的正确值可能没有一个自然的答案。

或许由于这种原因,在真实的数据中我们的模型选择标准的值随着增加的K而继续增加是不稀有的。那么集中于捕获数据中的大多数结构的K的值通常是讲得通的,这在生物学上似乎是合理的。

5.1估计K的步骤

1. (命令行版本)在文件extraparams里把COMPUTEPROBS和INFERALPHA设置为1。(前端版本)确保允许改变。

2. 对不同的MAXPOPS (K)值运行MCMC方案。最后它将输出一行“Estimated Ln

Prob of Data”。这是ln Pr(X|K)的估计。你应该对每个K独立地运行几次,以便证实不同运行得到的估计值是一致的。如果与不同的K获得的估计值的变异性相比,一个给定的K的不同运行的变异性是显著的,那么你可能需要使用更长的运行或者更长的burnin时期。如果lnPr(X|K)看起来是双峰的(bimodal)或者多峰的(multimodal),则MCMC方案可能到不同的答案。你可以对此进行验证,通过比较在单个K上的不同运行的Q。(参看Pritchard et al. (2000a)的数据集2A(Data

Set 2A),也见下面有关多峰性(Multimodality)的部分,)。

3. 计算K的后验概率。例如,对于论文中的数据集2A(这里K是2),我们得到

K

1

2

3

4

5

ln Pr(X|K)

-4356

-3983

-3982

-3983

-4006

我们一开始可以假定一个关于K = {1, ..., 5}的均匀先验分布。然后根据贝叶斯定理,Pr(K = 2)由下式给出:

e3983 (3)

43563983398239834006eeeee如果我们将该式简化为下面的公式,计算就会更容易

e10.21 (4)

e374e1e0e1e245.2 轻微的违背模型可能导致过高估计K

当存在真正的体结构时,这导致不连锁的位点之间的LD,以及违背哈迪温伯格比例。

粗略地说,这是被Structure算法使用的信号。但是模型的一些违背也能导致哈迪温伯格或连锁不平衡。这些包括近交和基因型鉴定错误(例如偶然的、未被发现的无效的等位基因)。即使在没有体结构的情况下,对于K >1,这些类型的因素也可能导致弱的统计信号。

从2版本开始,我们提出相关的等位基因频率模型(correlated allele frequency model)应该被用作默认,因为它在困难的问题上经常实现更好的执行,但是用户应该意识到,在这样的设置中可能更容易过高估计K,与独立的频率模型下相比(Falush et al. (2003a))。

下一节讨论怎样确定推断的结构是否是真实的。

5.3 关于选择K的非正式提示;结构是真实的吗?

有两个非正式的提示,可能有助于选择K。第一个是,对于比合适的值(有效零)更小的K,Pr(K)常常是非常小的,对于更大的K,则有或多或少的高原,如同上面显示的数据集2A的例子中那样。在这种情形中(其中K的几个值给出log Pr(X|K)的相似的估计下),似乎这些估计中最小的常常是正确的。

对于我们通过“或多或少的高原”所表示的东西,要提供一个坚固的规则有点难。对于小数据集来说,这可能意味着log Pr(X|K)的值在5-10的范围内,但是Daniel Falush写道“在非常大的数据集中,K = 3和K = 4之间的差别可能是50,但是如果K = 3和K = 2之间的差别是5

000,那么我将肯定选择K = 3”。想要使用更正式的标准(这种标准将这一点纳入了考虑)的读者可能对Evanno等(2005)的方法感兴趣。

我们认为考虑这一点的一种明智的方法是就模型选择而言。即,我们可能不总是能知道K的真值,但是我们应该致力于捕获数据里的主要结构的K的最小的值。

第二提示是,如果真的有单独的体,那个,通常有许多有关的值的信息,一旦Markov

链收敛,通常将相对恒定(范围经常为0.2或更少)。不过,如果没有任何真正的结构,在运行过程中通常变化很大。

这一点的一个必然的结果是当没有体结构时,你将通常将看到分配给每个体的样本的比例是大致匀称的(每个体中1/K),大多数个体将被公平地混和。如果一些个体被强烈地分配到一体或者另一个,以及如果分配给每组的比例不对称,那么这是你有真正的体结构的强的迹象。

假定你有两个清楚的体,但是你试图决定是否这些中之一是更进一步再分(例如,Pr(X|K = 3)的值类似于P(X|K = 2),或者也许比P(X|K = 2)还大一点)。那么,你能尝试的一

件事情是只使用你怀疑可能被再分的体内的个体来运行Structure,看看是否有一个如上所述的强信号。

总之,你应该对根据小的Pr(K)的差别推断的体结构持怀疑态度,如果(1)对于分派没有清楚的生物学解释,(2)对全部体的分派大致匀称,没有个体被强烈地分配。

5.4 通过距离数据的隔离

通过距离的隔离指的是这样的想法:个体可能跨越一些地区呈空间分布,带有本地分散的。在这种形势下,等位基因频率跨越地区逐渐变化。基础的Structure模型对来自这种情况的数据不很适合。当这发生时,推断的K 的值,以及在每组中的相应的等位基因频率可能相当任意。取决于取样的计划,大多数个体可能在多个组中具有混合的成员身份。即,算法将尝试使用K的不同组分的加权平均数来对跨越地区的等位基因频率建模。在这样的形势下,结果的解释可能具有挑战性。

6 背景LD和其他miscellania

6.1 序列数据,紧密连锁的SNP和单体型数据

Structure模型假定位点在体内是独立的(即,在体内不处于LD)。序列数据或者来自非重组区域的(比如Y染体或者mtDNA)的数据很可能违反这个假定。

如果你有序列数据或来自多个独立区域的密集的SNP数据,那么尽管数据不完全适合模型,Structure实际上可能表演得想当好。粗略地说,这将发生,倘若跨越不同的区域有足够的独立性,以至于区域内的LD不在数据中占优势。当有足够的独立区域时,区域内的依赖性(dependence)的主要代价将是Structure在特别的个体的分派中低估不确定性。

例如,Falush等(2003b)把Structure用于来自H. pylori的MLST(多位点序列)数据,以了解H. pylori的体结构和迁移历史。在那种情况下,在区域内有足够的重组以至于体结构的信号超过了背景LD。(关于MLST数据的更多情况,也见第10节)。在人类的应用中,Conrad等(2006)发现来自36个连锁的区域的3000个SNP生产明智(但是嘈杂)的答案,在一个全世界的样本中,基本上与基于微卫星的以前的结果一致[见他们的Supplementary

Methods Figure SM2]。

然而,如果数据被一个或者少数非重组的或在低重组的区域主导,那么,Structure可能被严重地误导。例如,如果数据只由Y染体数据组成,那么估计的结构大概将反映出关于Y染体树的某些事情,而非体结构本身。使用这样的数据的影响很可能是:(1)算法低估祖先估计中的不确定性的程度,在最坏的情况下,可能是有偏的或者不准确的;(2)K的估计不可能表演得好。如果你有Y或者mtDNA数据加上许多核标记,一个安全和有效的解决办法是重新编码来自每个连锁区域的单体型,以至于单体型被描述为一个具有n等位基因的单个位点。如果有许多单体型,则可以把相关的单体型归类到一起。

注意连锁模型不一定比(非)混合模型对于处理这些问题更好。连锁模型不是设计来处理体内的背景LD的,并且很可能被类似地干扰。

6.2 多峰性

Structure算法在参数空间中的一个随机的地方开始,然后朝着参数空间的一个峰收敛。(在这里的上下文内,峰可以被看作,松弛地讲,具有高的后验概率的一个聚类解。)当不使用先验的标签时,在K类的编号过程中没有固有的意思,因此有K!个对称的峰,对应于对类标签的排列。理论上,Structure可以在这些峰之间切换,但是这通常对真正的数据集不发生(Pritchard et al., 2000a)。为准备用于出版的图,诺厄·罗森堡(Noah Rosenberg)的实验室有一个有帮助的程序,CLUMPP,在对数据绘图之前跨越不同的运行把类标签排列成行(第10节)。

除了这些对称的峰之外,一些数据集可能还有另外的非对称的峰。Structure的当前的实施通常不在实际长度的运行过程中在这些之间穿过。这表明不同的运行可能产生显著不同的答案,并且更长时间的运行或许将不能修正这一点。

这主要是对于非常复杂的数据集的一个问题,具有大的K值,比如K >5或者K >10(但是见Pritchard等(2000a)中的数据集2A的例子)。你可以检查Q的结果,以了解这是否发生了。Rosenberg等(2001)提出了这类情形的仔细的分析,对一个数据集,其中估计的K大约是19。

6.3 当大多数个体被混合时,估计混合比例

如果亲本的体的代表非常少,估计混合比例可能特别具有挑战性。在Pritchard等(2000b)中对于模拟数据有一个这样的例子。数据假定是来自大多数个体有某种程度的欧洲祖先的美国黑人体的一个样本。对那些数据来说,估计的祖先比例与真实(模拟)值高度相关,但

是祖先的实际比例是有偏的。那个例子也是我们用真正的数据的更新近的经历的代表。

这发生因为在缺乏任何程度非掺和的个体的情况下,可能有一些不可识别性(nonidentifiability),在那里有可能把等位基因频率推得更开,把混合比例挤压到一起(反过来也是这样),获得几乎一样的模型拟合。当有强烈非对称的混合时,使用POPALPHAS

= 1(每个体单独的)能帮助一点,但是不解决基本的问题。因此,在这些情形下的混合性的估计应该被谨慎对待。

7 从命令行运行Structure

有许多由用户设置的程序参数。这些在两个文件中(mainparams和extraparams),每当程序执行的时候,这些文件被读取。mainparams指定数据文件的输入形式和最基本的运行参数。extraparams指定多种程序选项。你需要设置mainparams里的所有值,而以extraparams里的默认值开始运行或许没有问题。注意默认模型假定混合,并且不使用用户定义的PopData。

在这两个文件中每个参数都以大写字母打印,在参数前面有“#define”(它们在整个这份资料中也全部以大写字母打印。) 值是紧跟在参数的名字之后立即被确定的(例如“#define

NUMREPS 1000”把MCMC的重复数量设置为1000)。

继每个参数定义之后,有一个简短的注释(用“//”标明),描述这个参数。这包括期望的值的类型。这些包括: “(str)”,用于字符串(用于输入和输出文件的名字);“(int)”,用于整数;“(d)”,用于双精度数(即,一个实数,比如3.14;“(B)”,用于布尔数(Boolean)(即,参数取值为TRUE或FALSE,通过分别把这个设置为1或0)。

程序对参数的次序不敏感,因此你可以重排它们或者增加说明,等等。一次给定的运行所使用的全部参数的值被打印在输出文件的末尾。

7.1 程序参数

在这一节里我们列举可以被用户确定的所有参数。这些被根据用在Structure的命令行版本的参数文件排序。

7.2 文件mainparams里的参数

在运行程序之前,用户需要确定所有这些参数。这些参数的一些(LABEL, POPDATA,

POPFLAG, PHENOTYPE, EXTRACOLS)表明是否特殊类型的数据存在于输入文件中;这些在第2节中描述。

基本的程序参数

MAXPOPS (整数) 为程序的一次特定的运行假设的体的数目。Pritchard等(2000a)称之为K。有时(取决于数据的性质)有可以被使用的K的自然值,否则K可以通过在K的不同值检查模型的配合被估计(参阅第5节)。

BURNIN (整数) 在开始数据收集之前的burnin时期的长度。(参阅第3.3节)

NUMREPS (整数) 在burnin之后的MCMC重复的数目。(参阅第3.3节)

输入/输出文件

INFILE (字符串) 输入数据文件的名字。最大长度30个字符(或者也许更少,取决于操作系统)。

OUTFILE (字符串) 程序输出文件的名字( 后缀“_1”, “_2”, ...,“_m”(用于中间的结果)以及“_f”(最终结果)被添加到这个名字)。有这些名字的现有的文件将被覆盖。名字的最大的长度为30个字符(或者也许更少,取决于操作系统)。

数据文件格式

NUMINDS (整数) 数据文件里的个体的数目。

NUMLOCI (整数) 数据文件里的位点的数目。

PLOIDY (整数) 生物体的倍性。默认是2(二倍体)。

MISSING (整数) 为遗漏基因型数据给的价值。必须是一个整数,并且不可在数据集里的其他地方出现过。默认是-9。

ONEROWPERIND (布尔值) 每个个体的数据被安排在一单个的行中。例如,对于二倍体数据来说,这将意味着每个位点的两个等位基因在相同的行中按照连续的次序,而不是被安排在相同的纵行里的连续两行中。关于输入格式的细节请参阅第2节。

LABEL (布尔值) 输入文件包含每个个体的标签(名字)。1 = 是的;0 = 不。

POPDATA (布尔值) 输入文件为每个个体包含一个用户定义的体来源。1 = 是的;0 =

不。

POPFLAG (布尔值) 输入文件包含一个指示变量,它表明当USEPOPINFO==1(见下面)时是否使用popinfo。1 = 是的;0 = 不。

LOCDATA (布尔值) 输入文件为每个个体包含一个用户定义的取样位置。1 = 是的;0 =

不。用于LOCPRIOR模型。可以设置LOCISPOP = 1以便在LOCPRIOR模型里改为使用POPDATA。

PHENOTYPE (布尔值) 输入文件包含一个表现型信息的列。1 = 是的;0 = 不。

EXTRACOLS (整数) 在基因型数据开始之前,在表现型之后的附加的数据列的数目。这些被程序忽略。0 = 没有额外的列。

MARKERNAMES (布尔值) 数据文件的顶行包含对应于使用的标记的一个L名字的列表。

RECESSIVEALLELES (布尔值) 数据文件的下一行包含一个L整数的列表,表明每个位点上的等位基因是隐性的。把这设置为1意味着使用的是显性标记模型。

MAPDISTANCES (布尔值) 数据文件的下一行(或者如果MARKERNAMES==0,则为第一行)包含相邻的位点之间的mapdistances的一个列表。

高级数据文件选项

PHASED (布尔值) 供连锁模型使用。表明数据处于正确的连锁相。如果(LINKAGE=1,

PHASED=0),则PHASEINFO可以被使用——这是输入文件里的一个额外的行,给出连锁相的概率。当PHASEINFO = 0时每个值被设置为0.5,暗示没有连锁相的信息。当连锁模型被用于多倍体时,要求PHASED=1。

PHASEINFO (布尔值) 每个个体的基因型数据的行后面有一行关于单倍型连锁相的信息。这只供连锁模型使用。有关详细情形,参阅第2节和3.1节。

MARKOVPHASE (布尔值) 连锁相信息跟随一个Markov模型。欲了解详细信息,参见第2.2节和9.6节。

NOTAMBIGUOUS (整数) 当RECESSIVEALLELES=1时,供多倍体使用。定义代码,表明一个标记上的基因型数据是不含糊的。不可与数据里的MISSING或者任何等位基因值相匹配。

7.3 文件extraparams里的参数

这些选项允许用户以各种各样的方式来精炼模型,并且做复杂的分析。默认值首先或许是好的。对于Boolean选项,键入1用于“Yes”,或者“Use this option”;0用于“No”或者“Don’t

use this option”。

程序选项

NOADMIX (布尔值) 假定模型没有混合(Pritchard et al., 2000a)。(每个个体被假设为完全来自K个体之一)。在输出中,程序不是如同在混合的情况下打印Q的平均值,而是打印每个个体来自每个体的后验概率。1 = 没有混合;0 = 有混合的模型。

LINKAGE (布尔值) 使用连锁模型。参阅第3.1节。RLOG10START设置每单位距离的重组率r的初始值。RLOG10MIN和RLOG10MAX确定用于log10r的最小和最大允许值。RLOG10PROPSD设置每次更新中log10r的建议的改变的大小。前端对这些做出一些猜测,但是用户需要注意,确信那些值对于特定的应用是合理的。

USEPOPINFO (布尔值) 使用先验的体信息来分派个体到(cluster)。也见MIGRPRIOR

和GENSBACK。必须有POPDATA = 1。

LOCPRIOR (布尔值) 使用地点信息来提高结构信息弱的数据的性能。

FREQSCORR (双精度) 使用“F model”,其中等位基因频率跨越体相关 (Falush et al.,

2003a)。更具体地说,不是假定一个先验,在其中每个体中的等位基因频率是来自一个均匀的Dirichlet分布的独立抽样,而是从一个分布开始,这个分布以样本的平均等位基因频率为中心。这个模型对于亲缘关系非常近的体更真实(其中我们预期等位基因

频率跨越体相似),并且可以产生更好的聚类(第3.2节)。Fk的先验是用FPRIORMEAN

和FPRIORSD确定的。当FREQSCORR被开启时,可能有一种过高估计K的趋势。

ONEFST (布尔值) 对全部体假定相同的Fk的值(类似于Wright的传统的FST)。这不为大多数数据推荐,因为在实践中你或许期望每个体中的不同的分歧(divergence)水平。当K = 2时,有时可能难以分别估计FST的二个值(但是见Harter等(2004))。当你试图估计K时,你应该把相同的模型用于全部K(我们建议ONEFST = 0)。

INFERALPHA (布尔值) 根据数据推断模型参数的值;否则被固定在值ALPHA上,它是由用户选择的。在NOADMIX模型下这个选项被忽略。(祖先向量Q的先验是具有参

数(, , ..., )的Dirichlet。小的暗示大多数个体基本上来自一个体或者另一个体,而>1 暗示大多数个体被混和。)

POPALPHAS (布尔值) 为每个体推断一个单独的。多数情况下不推荐,但是对于具有非对称混合的情况可能是有用的。

ALPHA (双精度) 混合程度的Dirichlet参数()(如果INFERALPHA ==1,这个就是初始值)。

INFERLAMBDA (布尔值) 为推断一个适当的值。不为大多数分析推荐。

POPSPECIFICLAMBDA (布尔值) 为每个体推断一个单独的。

LAMBDA (双精度) 参数化等位基因频率的先验,对于大多数数据默认值1好像工作得相当好。如果大多数标记上的频率严重向低/高频率偏斜,更小的值可能导致潜在更好的表现。它好像对估计没有对其它超参数(和F)的估计那么好。

先验

这些值用来参数化假定的概率模型。多数情况下缺省设置应该相当明智,你可能不想要担心这些。

FPRIORMEAN,FPRIORSD (双精度) 见FREQSCORR。Fk的先验被取具有平均数FPRIORMEAN和标准偏差FPRIORSD的γ。我们的缺省设置将许多权重放于F的小值上。我们发现这使得算法对细微的结构敏感,但是增加了过高估计K的风险(Falush et al.,

2003a)。

UNIFPRIORALPHA (布尔值) ,ALPHAMAX (双精度) 对假设一个均匀分布的先验,运行在0和ALPHAMAX之间。这个模型好像工作很好;其它的模型(当UNIFPRIORALPHA

= 0时)是取作为具有γ分布的先验,具有平均数ALPHAPRIORA×ALPHAPRIORB,方差ALPHAPRIORA×ALPHAPRIORB2。

LOG10RMIN,LOG10RMAX,LOG10PROPSD,LOG10RSTART (双精度) 当连锁模型被使用时,切换速率r被取为在对数尺度上具有一个均匀的先验,在LOG10RMIN和LOG10RMAX之间。这些值需要由用户确定,以便就被使用的图距单位的尺度而言有意义。

使用先验的体信息(USEPOPINFO)

GENSBACK (整数) 这对应于G(Pritchard et al., 2000a) 。当对个体使用先验的体信息(USEPOPINFO = 1)时,程序测试是否每个体在最后G代有一个移民祖先,其中G = 0对应于个体本身作为一个移民。为了有适当的功效,G应该被设置得相当小(比如2),除非数据的信息非常丰富。

MIGRPRIOR (双精度) 必须是在[0,1]中。这是Pritchard等(2000a)中的v。明智的值可能在0.001—0.1的范围内。

PFROMPOPFLAGONLY (布尔值) 这个选项是2.0版新有的,使仅仅使用预先指定的个体的子集来更新等位基因频率P成为可能。为了使用这个选项,包括一个POPFLAG列,对应该用来更新P的个体设置POPFLAG=1,对不应该用来更新P的个体的设置POPFLAG = 0。USEPOPINFO打开或者不打开,都可以使用这个选项。

这个选项将是有用的,例如,如果你有来自已知体的一组准参考个体,然后你想要估计一些未知的个体的祖先。使用这个选项,每个未知个体的q估计值只依赖于参考个体组,而不依赖于样本中的其它未知个体。这个特点有时是合乎需要的。

使用地点信息的LOCPRIOR模型

LOCISPOP (布尔值) 当LOCPRIOR模型打开时,这个选项指导程序使用输入文件中的PopData列作为地点数据。当LOCISPOP = 0时,程序需要一个LocData列来使用LOCPRIOR。

LOCPRIORINIT (双精度) LOCPRIOR参数r的初始值,确定体的信息有多丰富(citepHubiszEtAl09)。我们发现LOCPRIORINIT = 1有助于取得好的收敛。

MAXLOCPRIOR (双精度) r的范围是来自(0,MAXLOCPRIOR)。我们建议MAXLOCPRIOR = 20。

输出选项

PRINTNET (布尔值) 打印聚类之间的“净核苷酸距离”。体A和B之间的这个距离DAB被计算为

DA,BJl(HAHB)1L(l)(l)ˆAˆ, (5)

1pp,jB,jLl1j12(l)ˆx其中p,j是体x中位点l上的等位基因j的等位基因频率的后验平均估计,L是位点的数目,

Jl是位点l上的等位基因的数目,其中

Jl1L(l)2ˆxHx1p,j, (6)

Ll1j1净核苷酸距离是分别来自体A和B的一对等位基因不同的平均概率,小于体内平均杂合性。或许更直觉地,这可以被认为是来自不同体的等位基因之间成对差异的平均数量,超出每个体内存在的变异的数量。该距离有合适的特性,即相似的体之间的距离接近0,并且尤其,DAA = 0。注意到距离是对称的,以至于DAB = DBA。这个距离适于画体的树来帮助直观化聚类之间的差异的水平(Falush et al., 2003b)。

PRINTKLD (布尔值) [不赞成]这个选项不再可用。

PRINTLAMBDA (布尔值) 打印的当前到屏幕。

PRINTQSUM (布尔值) 打印当前的Q估计的摘要到屏幕;这打印PopData的每个值的平均数。

SITEBYSITE (布尔值) (连锁模型)为数据里的每种基因型打印分派概率的完整的摘要。这被打印到一个单独的文件,具有后缀“ss”。这个文件可能是大的!

PRINTQHAT (布尔值) 当这个被打开时,Q的点估计不仅被打印到主要的结果文件,而且被打印到一个单独的后缀为“q”的文件。为了运行姊妹程序STRAT,需要这个文件。

UPDATEFREQ (整数) 对屏幕打印更新的频率。如果这= 0则自动设置。

PRINTLIKES (布尔值) 在每个迭代过程中把似然函数的当前值打印到屏幕。

INTERMEDSAVE (整数) 如果你渴望在运行结束之前看见初步的结果,你可以让程序在MCMC运行期间每隔一段时间将结果打印到文件。总共打印INTERMEDSAVE个这样的文件,在继BURNIN的完成之后相等的间隔。通过设置为0关掉这个选项。使用OUTFILE名字创造的这些文件的名字。

ECHODATA (布尔值) 把数据集的概要打印到屏幕和输出文件。(打印输入文件最顶行和最底行的开头和结尾,以便允许用户检查它已经被正确地读取。)

ANCESTDIST (布尔值) 收集关于每个体的Q的分布的信息,也只估计平均值。当这被打开时,输出文件包括每个q(i)的概率区间的左端和右端。(概率区间是置信区间的贝叶斯类似

物。)打印的值显示概率区间的中间100p%,其中p是0.0到1.0的范围内的一个数字,用ANCESTPINT设置。Q的分布通过记录在0和1之间的许多箱子的每一个点击(hit)的数目被估计,以便形成某种直方图。这些箱子的宽度一样大小,是使用NUMBOXES确定的。

杂项

COMPUTEPROB (布尔值) 在每次更新时打印数据的对数似然,并估计给定K和模型的数据的概率(参阅第5节)。这被用于估计K,并且对于burnin(老化?)是否足够长也是一个有用的诊断。关掉这个选项的主要原因将是加速程序(~ 10–15%)。

ADMBURNIN (整数) (当RECOMBINE=1时使用)当使用连锁模型时,在大多数情形下强烈推荐利用混合物模型的短的burnin(比方说500次迭代)。没有这样的burnin,连锁模型经常生产独特的结果。设置ADMBURNIN

ALPHAPROPSD (双精度) 的Metropolis-Hastings更新步骤涉及从一个具有平均数#和标准差ALPHAPROPSD>0的正态分布中挑选一个值。ALPHAPROPSD的值不影响Markov链的渐近的行为,但是可能对收敛的速率具有显著的影响。如果有许多关于的信息,则小的ALPHAPROPSD值对于获得合理的接受速率更可取。如果没有很多关于#的信息,则更大的值产生更好的混合。

STARTATPOPINFO (布尔值) 使用给定的体作为体起源的起始条件。(需要POPDATA==1)。这个选项提供Markov链正在正确地收敛的一个检查,如果你期望推断的结构与输入标签相配,而它没有。这个选项假定输入文件里的PopData在1和k之间,其中k  MAXPOPS。PopData不在这个范围内的个体随机起始。

RANDOMIZE (布尔值) 为每次运行使用一个不同的随机数种子,取自系统时钟。(也见SEED)

SEED(整数) 如果RANDOMIZE==0,则模拟种子被预置为SEED。这允许运行被准确地重复。如果RANDOMIZE==1,则在SEED里指定的任何值被忽略。注意即使当RANDOMIZE==1时,注意那个==1,程序输出仍然表明开始的种子值,这样的话,如果需要,可以重复特定的运行。

METROFREQ (整数) 在混合物模型下使用Metropolis-Hastings步骤来更新Q的频率。当这个选项被使用时,为每个q(i)选择一个新的提议值q(i)。这个提议值是从先验(即q(i) ~ D(,

, ..., ))中抽样的。有这个更新的基本原理是当很小时,它可能改进混合,通过使个体在体之间跳跃更容易。Metropolis-Hastings移动在每次METROFREQ迭代中使用一次。如果METROFREQ被调整到0,则它从未被使用。

REPORTHITRATE (布尔值) 报告q(i)的Metropolis更新的接受速率(见METROFREQ)。

7.4 命令行转换成参数值

为了简化批运行以及使涉及结构的模拟更容易,我们已经增加了命令行标签(flag),更新某些参数的值,覆盖在mainparams里设置的值。这些如下:

-m(mainparams) 读取一个不同的参数输入文件而不是mainparams。

-e(extraparams) 读取一个不同的参数输入文件而不是extraparams。

-s(stratparams) 读取一个不同的参数输入文件而不是stratparams。(供姊妹程序STRAT使用,用于关联作图)

-K(MAXPOPS) 改变体的数目。

-L(NUMLOCI) 改变位点的数目。

-N(NUMINDS) 改变个体的数目。

-i (输入文件) 从一个不同的输入文件中读取数据。

-o(输出文件) 把结果打印到一个不同的输出文件。

-D(SEED) 使用值SEED初始化随机数的产生。注意使用这个选项时那RANDOMIZE

MUST必须被设置为0。)

因此,为了覆盖预设参数值之一,我们调用Structure然后使用相关的标签,继之以的新参数值。标签和新值用空格隔开。标签可以按任意顺序使用。

例如,为了把假定的体数目改为5,并且把输出指向一个名称为output5的文件,我们可以调用Structure如下:

./structure -K 5 -o output5

8 前端

这个部分提供一些一般的指导,以及关于使用前端的一点建议。一般的主题在上面讨论

了,你可以通过查看第7节得到关于各种参数选项的详细资料。

8.1 下载和安装。

首先,从网页中下载合适的程序文件。对于不同的平台有单独的版本(目前有Windows,Sun,Linux和Mac OS X)。

Windows文件是一个可执行的安装文件。双击图标开始安装。将引导你进行安装。通过双击Structure图标运行程序。

在Unix或者苹果机系统上,把文件放进一个临时文件夹。然后,解压缩该文件(“gzip –dc

<文件名> | tar xvf - ”),这里<文件名>是下载的文件的名字。通过键入“./install”来运行安装脚本。安装成功以后,将创建一个Structure开始脚本,这个脚本也可以移动到程序的一个标准文件夹中,例如,/usr/local/bin/。为了启动前端,仅仅执行这个开始脚本就可以了。

除了Windows OS,我们不再与Structure软件包一起发布Java虚拟机(从Structure 2.2

版开始)。安装Structure之前需要由Sun Microsystem开发的一种Java运行期环境(Java Runtime

Environment,JRE,版本>1.5.0)。与各种操作系统的兼容的JRE可以免费从

/download下载。在那个网站上可以到JRE的安装指导。

图1:该例子显示一个项目的组成部分。项目数据(Project Data)是数据文件;项目信息(Project

Information)指定数据文件格式。模拟摘要(Simulation Summary)作为这个项目的一部分运行的全部MCMC模拟的摘要。参数集(Parameter Sets)由使用不同的参数设置的3组MCMC运行组成:ps1,ps2,ps3;这些中的每一个显示设置(Settings),以及一个利用这些参数值运行的完成的MCMC的结果的列表。用户可以点击任何这些来查看细节。

8.2 概述

前端把数据分析组织到“项目”里。每个项目连接到一个单个的数据文件。当创建一个项目时,用户也提供信息,确切说明怎样读取数据文件(位点的数目,个体的数目,等等)。

这些是数据文件的特征,并且在这个项目内总是相同的。

每个项目也包含一或更多个“参数集”。这些允许用户指定MCMC运行的细节,包括重复的次数,老化(burnin)长度,等等,也指定分析的模型(例如,是否允许混合,等位基因频率的模型,等等)。然后用户可以在K的选定值上对一个给定的参数集运行Markov链。图1显示了一个称为“popdata”的项目的组成部分的例子。

程序然后可以使用这些参数值被运行。前端储存结果的各种各样的归纳,包括许多图表,在下面描述。

8.3 建立项目

首先你需要建立一个输入文件。这在第2节介绍了。

现在,点击File(文件)  New Project(新项目)。这打开一个向导来导入数据(图2)。数据被从指定的输入文件复制到为项目选择的工作目录中。

该向导由4个框组成:

1. 指定项目目录,项目名称和输入数据文件(图2)

2. 指定数据文件的基本特性(个体的数目,数据的倍性(对二倍体生物体键入‘2’),位点的数目,用来表明漏缺数据的值。点击,“Show data file format”来获得数据文件里的长度和行数的摘要(图3)。

图2:导入数据(第1步)。用户为项目指定目录(这里为数据),项目目录的名字(样本项目;这是在数据内的一个目录),以及要被程序样本数据读取的数据文件。

3. (Rows(行))指定(如果有的话)存在哪些可选的行数据:标记名称行;标记间距离的行;以及每个体之后的一个连锁相数据行。如果每个个体的数据储存在单个行内,而不是按照标准格式每个个体二行的话,还要点击“single line”。

4. (Columns(列))指定有哪些可选的列数据:个体的身份证明(标签);来源的体(POPDATA);USEPOPINFO标签(flag——说明当使用先验的体信息模型时为某些个体使用POPDATA信息的标签;表现型数据(用于关联作图(Pritchard et

al., 2000b));应该被Structure忽略的基因型数据之前的其他额外的数据列。

当你已经完成这些步骤时,你将得到数据格式的一个摘要;如果这看起来是正确的,则点击“proceed”。程序现在将试图加载数据文件并且创建新项目(图4)。

图3:导入数据(第2步)——指定数据文件的特性。在这里,以及向导的接下来的两个框(未示出)中,用户指定数据文件的特性(位点的数目,个体的数目,数据的类型,等等)。

图4:在成功地导入数据之后,前端加载数据文件,并创建新项目。

8.4 配置参数集

一旦你已经成功加载了一个数据文件,你就可以开始运行Structure了。你将创建一个或多个“参数集”;这些代表了你做出的关于如何分析数据的选择的全部列表。我们已经输入了一系列缺省设置,并且这些是开始的好地方。你或许将想要对每个参数集运行Structure多次,在K的不同的值上,前端被建立使这变得容易。

图5:指定一个新的参数集:设置运行长度。

转到“Parameter Set”(参数集)下的下拉菜单。你可以建立一个新的参数集,修改一个存在的参数集,或者删除一个参数集。点击“New”。你现在看见一个对话框,有4个标签(图5)。点击这些中的每个:

Run Length(运行长度)。关于这个的讨论请参阅第3.3节。注意前端提供了一些主要参数的时间序列曲线,帮助你确定运行长度是否看起来是足够的。

Ancestry Model(祖先模型)。参阅AncestryModels节。对于大多数数据集,混合模型是一个起始的好地方。连锁模型将失活,除非你输入关于标记的连锁信息。第7节提供了一些详细的选项的一些额外的细节,包括“Use Population Information”选项下的GENSBACK和MIGRPRIOR。注意连锁模型是相对计算密集的。

Allele Frequency Model(等位基因频率模型)。参阅第3.2节。我们建议使用相关的频率模型和独立的频率模型。相关的频率模型有更好的功效发现细微的体结构,但是K的后验概率可能稍微向上有偏。相关的频率模型按FST被参数化,每个体具有单独的参数(细节在第7.3节)。在最初的调查期间或许不必要推断。

Advanced(高级)。关闭计算后验概率(用于估计K)的功能使程序显著加速。你也可以让程序输出每个个体的祖先的后验的可信区域(credible region)(见ANCESTDIST,第7.3节)。在第7.3节的STARTATPOPINFO下更详细地描述了“Initialize at POPINFO”。

8.5 运行模拟

既然你有了一个参数集,你可以开始运行程序,通过转到[Parameter Set](参数集) 

[Run](运行)。你将被要求设置体的数目(K)。你也可以在相同的地方停止模拟([Parameter

Set](参数集) [Stop](停止))。

文本数据将被打印到屏幕底部的操纵台(图6)。你也可以查看各种的关键归纳统计量的实时时间序列曲线:FST,,似然值,等等,(图7)。

图6:运行时间输出被显示在底部的操纵台里

图7:FST的时间序列曲线

一旦你有不止一个参数集,你需要确切说明你想要为新的MCMC运行使用哪一个参数集。任何时候,一个参数集被称为“active”(活动的)(见左边的窗口)。你可以切换活动参数集,通过转到[Parameter Set](参数集) [Parameter Set List](参数集列表),并且加亮适当的选择,或者通过双击项目树里相应的参数集提示。

8.6 批运行

你可以安排一系列Structure运行,通过转到[Project]  [Start a Job]。这打开一个调度程序,允许你挑选:(1)参数集(使用控制键+鼠标来选择多个参数集);(2)K的值;(3)每

个的运行的次数。参阅图8。

8.7 从前端导出参数文件

你可以使用前端使写基于文本的参数文件自动化,用于Structure的命令行版本。转到[Project] “generating ”。这个选项可能是有用的,例如,帮助你在一组计算上建立大量的运行。

8.8 从自命令行程序导入结果

来自Structure程序的命令行版本的结果可以被导入到前端中,通过转到[File]  [Load

]。你将被要求提供2个文件:Structure结果文件(必需的;这个文件通常有后缀“_f”),一个文件包含对屏幕的运行时间Structure输出(可选的)。流动通过利用重定向到一个文件(例如,Structure> )的输出在操纵台中运行Structure,可以获得后者。成功加载文件后,你应该可以在前端中读取条形图,三角形图,以及各种的时间序列曲线(如果你提供运行时间文件则只有后者)。

图8:建立一批运行任务

8.9 分析结果

Summaries(摘要)。你可以查看目前完成的全部运行的详细摘要细节一张表,通过转到[View]  [Simulation Summary]。

Text results(文本结果)。你可以查看每个运行的全文结果,通过转到左边的窗口,然后点击适当的运行。全文输出将出现在右边的窗口里。解释文本输出的细节在第9节给出。

ˆ(每个个体的估计的成Plots of Ancestry estimates(祖先估计值的曲线)。我们提供Q员系数(membership coefficient),在每个聚类中)的两种类型的曲线图。当你在左边的窗口里点击相应的运行时,第一个图像自动出现。数据集里的每个个体都表示为一条单个的垂直的线,它被剖分成K个彩的片段,代表个体在K个推断的聚类的每一个中的估计的成员比例(membership fraction)。

个体的祖先的第2个图像把每个个体绘图到一个三角形(图10)。这类图对于直观化K =

3的数据是有用的(Pritchard et al., 2000a)。这是一种探索高维数据的有趣工具,但是条形图通常更容易解释。

Plots of summary statistics(归纳统计量的绘图)。我们也提供若干有趣的归纳统计量的绘图。[[Plotting]  [Data Plotting]包含运行过程中某些归纳统计量的时间序列曲线。一个例子如图7所示;注意到在运行一开始的时候有一股短的阶段,其中值在达到其稳定分布之前大大地增加。你应该检查这些曲线,以确信归纳统计量似乎可以在老化(burnin)结束之前稳定。

还有FST和的直方图(图8.9)。这些是对这些参数的后验分布的估计。

图9:Q的估计值的归纳绘图。每个个体表示为一个垂直线,垂直线被分解成K个彩的片段,片段的长度与K个推断的聚类的每一个成正比。数字(1..4)对应于预先规定的体。

图10:Q的三角形图。每个个体都表示为一个有的点。颜对应于先验的体标签。个体的估计的祖先矢量由总计为1的K个分量组成。当K = 3时,祖先矢量可以被绘图到一个三角形上,象显示的那样。对一个给定的点来说,3个分量中的每一个由到三角形的一边的距离给出。因此在角落之一的个体被完全分配到一个体或者另一个体。对K >3来说,我们是这样来表示数据的:允许用户一次挑选出两个推断的聚类,然后把所有其它的聚类归到一起。

9 解释文本输出

这一节描述在运行期间被打印到操纵台和输出文件的数据。前端也提供额外地数据图,在下面描述。

图11:Fst的直方图

9.1 运行期间的屏幕输出

Rep#: Alpha Corr

4500: 0.009 6.1

4525: 0.010 9.1

4550: 0.008 7.3

D1,2 Ln Like Est Ln P(D)

0.831 -2730 -2776

0.806 -2716 -2776

0.792 -2734 -2775

上面的例子显示了在运行期间打印到屏幕的输出的一部分。这里,“rep”给出到目前为止的MCMC迭代的次数,包括老化(burnin);“Alpha”是的当前值;“Corr”是相关系数f(l)的当前平均数(跨越位点)(见FREQSCORR);“D1,2”是体1和2之间的分歧(divergence)的度量(见PRINTNET);“Ln Like”给定P和Q的当前值时数据的对数似然;“Est Ln P(D)”是ln(P(X|K))的当前估计值(从burnin阶段结束起的全部迭代上平均)。这些列中的一些可能不被打印,取决于选择的程序选项。如果K变大,为了使输出更容易阅读,可以设置PRINTNET=0。

9.2 Q的打印资料

Label (%Miss) Pop:

1 17 (0) 2 :

2 :

2 :

1 :

3 :

Inferred clusters (and 90% probability intervals)

0.977 0.023 (0.829,1.000) (0.000,0.171)

0.997 0.003 (0.988,1.000) (0.000,0.012)

0.833 0.167 (0.003,1.000) (0.000,0.997)

0.005 0.995 (0.000,0.020) (0.980,1.000)

0.006 0.994 (0.000,0.016) (0.984,1.000)

2 1219 (7)

3 1223 (0)

4 1329 (7)

5 15 (0)

当程序被运行而没有使用先验的体信息时,Q的结果被以上面显示的格式给出(这里K被设置为2)。这阅读如下。从第1行读:个体的标签(取自数据文件)= 17;这个个体的漏缺数据的百分比 = 0%;用户分配的体 = 2;在聚类1和2中估计的成员分别 = 0.977和0.023(这些是q(17)的平均值);在q1(17)和q2(17)上的90%的概率区间分别是(0.829,1.000)和(0.000,0.171)。(注意只有ANCESTDIST被打开时才显示概率区间。)

9.3 当使用先验的体信息时Q的打印资料

当程序使用先验的体信息运行时,打印资料有点不同。这里我们已经以K = 3,

USEPOPINFO = 1,以及GENSBACK = 1运行了相同的数据:

Label (%Miss) Pop :

1 17 (0) 2 : 0.998 | Pop 1: 0.000 0.001 | Pop 3:

0.000 0.000 | Pop 3:

0.000 0.000 | Pop 3:

0.000 0.000 | Pop 3:

0.000 0.001 | Pop 2:

0.000 0.000

0.000 0.000

0.004 0.383

0.000 0.000

0.000 0.000

2 1219 (7) 2 : 1.000 | Pop 1:

3 1223 (0) 2 : 0.612 | Pop 1:

4 1329 (7) 1 : 1.000 | Pop 2:

5 15 (0) 3 : 0.999 | Pop 1:

通常,结果的第一列(继“Pop : ”之后)显示有问题的个体被正确地分配到给定体的后验概率。随后的列显示它来自其它体(或者在其它体中有祖先)的概率。对于其它体的每一个有GENSBACK+1个条目,显示个体来自那个体、具有来自那个体的父母、祖父母、曾祖父母,...等等,的概率(按这个顺序)。

例如,从第3行读:个体1223(它有0%的漏缺数据)实际上以0.612的概率来自假定的体(2)。这个个体在体1中有新近的祖先的后验概率(大约)为零,但是它可能在体3中有新近的祖先(该个体来自体3,或者有来自体3的单个的父母的概率分别是0.004和0.383)。

9.4 等位基因频率分歧的打印资料

聚类之间的等位基因频率分歧(divergence)(净核苷酸距离),利用P的点估计值数计算。

1

1 -

2

0.0357

- 2 0.0357

相同的聚类中的个体之间的平均距离(期望的杂合性):

cluster 1 :

cluster 2 :

0.7845

0.7686

这个例子显示K = 2个聚类之间的成对净距离的矩阵(上面),以及每个聚类内部的(期望)杂合性(下面)。关于怎样计算这些距离的细节见第7.3节中的PRINTNET。注意成对矩阵是对称的(即Dij = Dji)。这些距离替换了2.2版以下的Structure中的库尔贝克-莱布勒分歧(Kullback-Leibler divergence)。

你可能发现基于净核苷酸距离画出树来表示聚类之间的距离是有帮助的。例子被显示在Falush等(2003b)中。

9.5 估计的等位基因频率(P)的打印资料

Locus 5

3 alleles

19.0% missing data

2 0.511 0.821 0.656

3 0.444 0.171 0.317

1 0.045 0.008 0.027

这个例子显示了位点5上估计的等位基因频率(P)的打印资料。列2、3和4显示在列1中列出的等位基因估计的频率(分别在聚类1、2和3中)。

9.6 连锁模型的位点接位点输出

当SITEBYSITE选项被选择时,有单独的输出文件,具有后缀“_ss”,包含每个个体在每个位点上的每个等位基因拷贝的来源分派的后验体。对大的数据集来说,这个文件可能有好几兆。

每一行显示一个个体的一个位点的分派概率(assignment probability)。该行的前两列注明个体的数目(从1到NUMINDS)和位点的数目(分布在从1到NUMLOCI的范围内),按照它们在数据文件里出现的次序。

后验分派概率的格式取决于参数组合。如果LINKAGE=0或者PHASED=1,那么输出的前面K行给出这个位点上的第一个等位基因拷贝来自体1…K的概率。对于二倍体或者多倍体数据,随后的等位基因拷贝的类似的概率在另外的列内显示。

如果使用连锁模型(LINKAGE=1)并且数据没被完全确定连锁相(PHASED=0),则每个位点上的等位基因拷贝的后验分派概率可能是强烈相互依赖的(co-dependent)。Structure因此输出二个等位基因拷贝的联合分派概率,意味着每个位点有K2个条目(注意这个选项对PLOIDY 2无效)。

如果MARKOVPHASE = 1,则前面K列给出数据文件中的第一个等位基因拷贝在体1中、第2个等位基因拷贝在体1…K中的概率,随后的列与体2…K中的第一个等位基因拷

贝的概率相关。如果MARKOVPHASE = 0,则不是指到数据文件中的第1和第2个列举的等位基因拷贝,概率指的是母亲和父亲链(strand)的起源的体。如果没有连锁相的信息(PHASEINFO = 0),则后验概率基体理论上应该是对称的,因此母亲的等位基因在体k1中而父亲的等位基因在k2中的概率将等于母亲的等位基因在体k2而父亲的等位基因在体k1中的概率。实际上,因为用MCMC来估计该矩阵,如果NUMREPS是小的,将有显著的与对称性的不符合。

例如,假定下面是一个没有连锁相信息的二倍体个体的两个位点的位点接位点(site-by-site)输出,具有MARKOVPHASE = 0。

1 1 0.001 0.000 0.008 0.000 0.000

0.001 0.007 0.001 0.982

1 2 0.001 0.000 0.008 0.000 0.000

0.001 0.008 0.001 0.982

然后为了计算第一个位点上的母亲和父亲等位基因拷贝的分派概率,数字被累加如下:

locus 1 pop1 pop2 pop3 origin of maternal(X)

chromosome

pop1

pop2

pop3

0.001 0.000 0.008

0.000 0.000 0.001

0.008 0.000 0.982

0.009

0.001

0.990

origin of paternal 0.009 0.000 0.991

chromosome (missing)

在这个例子中,数据来自一个雄性的X染体,因此实际上第2个等位基因拷贝是缺失的。

注意该格式被从2.1版本简化,在那里结果被放在与其余输出相同的文件里。标签和标记名字不再被打印,输出用小数形式而不是科学计数法打印每个数字。这些变化是为了紧凑。

10 供Structure使用的其他资源

10.1 Structure结果的绘图

CLUMPP和distruct是诺厄·罗森堡(Noah Rosenberg)的实验室为制作Q矩阵的好的图编制的两个程序。前端生产相似的图,但是这两个程序提供对图表输出的很多更精细的控制。见

/

10.2 将细菌的MLST数据导入Structure格式

由Xavier Didelot和Daniel Falush开发的软件xfma2struct按照扩展的Fasta格式提取单倍体序列数据,并且把它们转变成Structure格式。见ClonalFrame网站:

/

11 怎样引用这个程序

对基本方法的合适的引用是Pritchard et al. (2000a)。Falush等(2003a)的文章是在2.0版中实现的连锁模型和相关等位基因频率模型的合适的参考文献。含糊的基因型数据,例如显性标记(2.2版新增)的方法由Falush等(2007)描述。小数据集的信息性先验的模型(2.3版新增)由Hubisz等(2009)描述。

12 书目

参考文献

Beaumont, M., Gottelli, D., Barratt, E. M., Kitchener, A. C., Daniels, M. J., Pritchard, J. K., and

Bruford, M. W. (2001). Genetic diversity and introgression in the Scottish wildcat. Molecular

Ecology, 10:319–336.

Conrad, D., Jakobsson, M., Coop, G., Wen, X., Wall, J., Rosenberg, N., and Pritchard, J. (2006).

A worldwide survey of haplotype variation and linkage disequilibrium in the human genome.

Nature Genetics, 38:1251–1260.

Evanno, G., Regnaut, S., and Goudet, J. (2005). Detecting the number of clusters of individuals

using the software STRUCTURE: a simulation study. Mol. Ecol., 14:2611–2620.

Falush, D., Stephens, M., and Pritchard, J. K. (2003a). Inference of population structure:

Extensions to linked loci and correlated allele frequencies. Genetics, 164:1567–1587.

Falush, D., Stephens, M., and Pritchard, J. K. (2007). Inference of population structure using

multilocus genotype data: dominant markers and null alleles. Molecular Ecology Notes,

7:574–578.

Falush, D., Wirth, T., Linz, B., Pritchard, J. K., Stephens, M., and 13 others (2003b). Traces of

human migrations in Helicobacter pylori populations. Science, 299:1582–1585.

Harter, A., Gardner, K., Falush, D., Lentz, D., Bye, R., and Rieseberg, L. (2004). Origin of extant

domesticated sunflowers in eastern North America. Nature, 430:201–205.

Hubisz, M., Falush, D., Stephens, M., and Pritchard, J. (2009). Inferring weak population structure

with the assistance of sample group information. Molecular Ecology Resources, In Press.

Murgia, C., Pritchard, J. K., Kim, S., Fassati, A., andWeiss., R. (2006). Clonal origin and

evolution of a transmissible cancer. Cell, 126:477–487.

Pritchard, J. K., Stephens, M., and Donnelly, P. (2000a). Inference of population structure using

multilocus genotype data. Genetics, 155:945–959.

Pritchard, J. K., Stephens, M., Rosenberg, N. A., and Donnelly, P. (2000b). Association mapping

in structured populations. Am. J. Hum. Genet., 67:170–181.

Rosenberg, N. A., Burke, T., Elo, K., Feldman, M. W., Freidlin, P. J., Groenen, M. A., Hillel, J.,

Maki-Tanila, A., Tixier-Boichard, M., Vignal, A., Wimmers, K., and Weigend, S. (2001).

Empirical evaluation of genetic clustering methods using multilocus genotypes from 20

chicken breeds. Genetics, 159:699–713.


本文发布于:2024-09-22 13:26:44,感谢您对本站的认可!

本文链接:https://www.17tex.com/fanyi/284.html

版权声明:本站内容均来自互联网,仅供演示用,请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系,我们将在24小时内删除。

标签:数据   群体   模型   个体   使用   可能
留言与评论(共有 0 条评论)
   
验证码:
Copyright ©2019-2024 Comsenz Inc.Powered by © 易纺专利技术学习网 豫ICP备2022007602号 豫公网安备41160202000603 站长QQ:729038198 关于我们 投诉建议