2005-12-16

Arlequinを用いた遺伝解析実習

集団遺伝学実習 Arlequin

全７限

このシリーズでの取り扱い範囲

遺伝子(もしくは、ある領域)の遺伝的多様性の評価(記事はこちら)
- Heterozygosity、ハプロタイプ間距離、ハプロタイプ関係樹・ハプロタイプ関係ネットワーク、Mismatch distribution、ハプロタイプ頻度推定
連鎖不平衡の評価(記事はこちら)
- ペアワイズ連鎖不平衡係数(D', $r^2$ 、LD検定(最尤法))
中立仮説検定(記事はこちら)
- 観測遺伝子多様性の中立仮説との適合度検定
集団の構造化解析(記事はこちら)
- 遺伝的異質集団の組み合わせに分解

集団遺伝学。その基礎。資料整備が間に合わないので、当座は、こちら(National Biological Information Infrastructure,USA)で。

2005-12-16

第４限　Haplotypic dataとGenotypic data、そのGametic phase既知と未知

集団遺伝学実習 Arlequin

第３限にて、Haplotypic dataを用いて、Diversity indicesの解析とその結果の解釈を行った。第４限では、Genotypic dataに拡張する。

！注意！　Diversity indicesのうち、Standard diversity、および、Theta(Hom)はhaplotypic dataとPhase既知のGenotypic dataにのみ適用すること。これらは、Phase未知のジェノタイプサンプルの記載の仕方により結果が異なる。２SNPのダブルヘテロのデータを



sample AA

       CC

と書くか

sample AC

sample CA

と書くか

で、結果が異なる。これらは集団中のホモ個体率・ヘテロ個体率を用いて計算しているからである。また、MINIMUM SPANNING TREE/NETWORK、Mismatch distributionは出力されない。Haplotype frequencyはEMアルゴリズムの結果が出力される。Theta(Hom)を除くMolecular diversityのデータが得られるのは、それらは、ハプロタイプフェーズに依存しない指標だからである。

サンプルデータの作成
- 第３限で用いたhaplotypic dataと同じ染色体構成となるようなGametic phase 既知のGenotypic dataを作成して解析を実行せよ。そうすることで、Gametic dataの解析とHaplotypic dataの解析の異動を確認する。Genotypic data-Phase既知の場合には、ハプロタイプベースの解析はHaplotypic dataのそれと同一になるはずであり、Genotypic data-Phase未知の場合には、推定haplotypeがhaplotypic dataのそれと異なるので、ハプロタイプベースの解析が異なることを、確認する。実際には、Standard diversity indices, Molecular diversity indices、Mismatch distributionはジェノタイプ情報、ハプロタイプ情報にのみ依存する。次の２つの内容を比較せよ。上はhaplotyic data解析で用いたSampleAのデータ、下はそれと同じhaplotic構成となるようなGenotypic dataの１例である。SampleSizeが半分(102=204/2)となっていること、また、上のA1,...A5はホモとして人数が半分に、A6,A7はあわせて4本の染色体からA6ホモ1人、A6-A7のヘテロ１人を作っている。Phase既知の場合はGameticPhase=1,未知の場合はGameticPhase=0とする



SampleName="SampleA"

SampleSize=204

SampleData={

A1  98 CTTGGA

A2  52 GTTGGG

A3  34 GTTGGA

A4  10 GTTTGG

A5  6 GATGAA

A6  3 GTCGGA

A7  1 CATGGG

}



SampleName="SampleA"

SampleSize=102

SampleData={

A1  49 CTTGGA

       CTTGGA

A2  26 GTTGGG

       GTTGGG

A3  17 GTTGGA

       GTTGGA

A4  5 GTTTGG

      GTTTGG

A5  3 GATGAA

      GATGAA

A6  1 GTCGGA

      GTCGGA

A7  1 CATGGG

      GTCGGA

}

Genotypic data Phase既知の場合とhaplotypic dataの違い
- Standard diversity indices



===============================

== Standard diversity indices : (SampleA)

===============================

<中略>

Sum of square freqs.     :  0.3270

Gene diversity           :  0.6763   0.0224

No. of Heterozigotes     : 1

No. of Homozygotes       : 101

Obs.Heter.               : 0.0098

Exp.Heter.               : 0.6763

No. of Heterozigotes, No. of Homozygotesの２項が追加されている
- Genotypic サンプルのうちヘテロ個体数・ホモ個体数であり、Genotypic data固有の情報
Obs.Heter., Exp.Heter.の２項が追加されている
- Obs.Heter.は入力データにおけるヘテロ個体の割合であり、Exp.Heter.はHardy-Weinberg平衡を仮定したときのヘテロ個体の割合の期待値である。Gene diversityの値と一致していることを確認せよ。

Molecular diversity indices
- 距離行列・MINIMUM SPANNING TREE/NETWORK・Thetas
  - ハプロタイプ名にh1,...h7が振られている
  - それ以外はすべて同一である

Inter-haplotypic distance matrix (s.d. above diagonal):

h1 h2 h3 h4 h5 h6 h7

h1 1.1547 0.9129 1.2247 1.2247 1.1547 1.1547
h2 2.0000 0.9129 0.9129 1.2247 1.1547 1.1547
h3 1.0000 1.0000 1.1547 1.1547 0.9129 1.2247
h4 3.0000 1.0000 2.0000 1.1547 1.2247 1.2247
h5 3.0000 3.0000 2.0000 4.0000 1.2247 1.2247
h6 2.0000 2.0000 1.0000 3.0000 3.0000 1.1547
h7 2.0000 2.0000 3.0000 3.0000 3.0000 4.0000

List of Haplotypes:

h1 : CTTGGA
h2 : GTTGGG
h3 : GTTGGA
h4 : GTTTGG
h5 : GATGAA
h6 : GTCGGA
h7 : CATGGG

Mismatch distribution
- haplotypic dataのときと異なるのは、以下の抜粋項目中、P(Sim. Ssd ....), P(Sim. Rag ....)のみであることを確認せよ。その値は、下記とも異なっているはずである
- また、同一入力ファイルにて、再度実行し、その値場合とも比較せよ。やはりPのみが異なるはずである
- Estimated Parameters...、および、Sum of Squared deviation, Harpending's R...は、決定的解析結果(Determinisitic) (入力データと、定められた計算式とから算術的に計算される値)であるために、同一である。他方、２つのP値は、Bootstrapサンプリング(乱数を用いて複数回試行した結果を用いて得られた確率分布を基準にして得られる値(確率的解析結果(Stochastic))であるために、解析のために異なる値が得られている。確率的解析結果の再現性を確認するためには、解析にて用いる擬似乱数列の初期引数(シード(seed,種))を特定して再実行すればよい。ただしArlequinの実行オプションには、シードの選択はないようである



==========================

== Mismatch distribution : (SampleA)

==========================
Estimated parameters from the sudden expension model:

2005-12-16

第１限　インストールと起動、設定

集団遺伝学実習 Arlequin

以下にはWindows XPにおけるインストール。Mac,Linuxへのインストールについては、オフィシャルページの解説を適宜参考にすること

ホームページ
ダウンロード
- bin + jre版"Arlequin20jre_zip.exe"をダウンロードし、ダブルクリック。解凍されて現れる２つの.exeファイル"jre117-win32.exe""Arlequin20_zip.exe"を順次ダブルクリック。
- Arlequin 2001 update (patch)"arlpatch2001.zipをダウンロードし、unzipし、できた２ファイル"arlecore2001.exe""arlequin2001.jar"を"ArlequinFolder"にコピーし、同フォルダに存在している"arlecore.exe"を"arlecore2000.exe"と改名、"arlequin.jar"ファイルを"arlequin2000.jar"と改名した上で、コピーしてきた"arlecore2001.exe"を"arlecore.exe"に、"arlequin2001.jar"を"arlequin.jar"に改名する
arlequin.exeをダブルクリックして起動
設定
- サンプル入力ファイル("hoge.arq")の読み込み
  - arlequin.exeをダブルクリックすると、メインウィンドウが立ち上がる
  - メニューバーより、File->Open projectを選択すると、ファイル読み込みウィンドウが立ち上がる。右下のAdd to list...ボタンを押すとファイル選択ウィンドウが立ち上がるので、サンプルファイルを指定する。サンプルファイルは、アプリケーションのダウンロード・解凍の後に.../ArlequinFolder/datafiles/hoge/以下に出来上がった、拡張子が"arp"のファイルである。List of recent Projectsに選択ファイルが追加されるので、それをクリックして左下のOKボタンを押す。実行画面が立ち上がる
- 実行画面から、設定
  - Configurationタブを選び、中段の２つのBrowseボタンを用いて実行結果の出力のための基本アプリケーションを選ぶ。上のBrowseボタンで選ぶのは、htmファイルを開くアプリケーションであるブラウザ(Internet Explorer、FireFoxなど)。下のボタンで選ぶのはテキストエディタである(NotePad、秀丸など)
試験実行
- ウィンドウ上部の"Run"ボタンを押すと、実行されて、結果がブラウザで表示される。細かい設定などは、気にせず、とりあえず動くことをここまでの目標とする

2005-12-16

第２限　入力ファイルと全解析に共通する設定

集団遺伝学実習 Arlequin

Arlequinの入力ファイルフォーマット・出力フォーマットは大量データのコマンドラインを意識したつくりになっていないので、そのつもりで。

SNP解析を前提にする。それ以外の場合は、SNP解析での経験をもとに、マニュアルを参照のこと。データを構成する基幹情報は以下のとおり(例示用データにおける数を挙げてある)

集団数 Nsample = 2 (SampleA,SampleB)
SNP数 Nsnp = 5
各サンプリング集団の検体数 SampleSize = (SampleA 6人、SampleB 4人)とする
検体がDiploidかHaploidか
- Diploidの場合には、Phase既知か未知か
検体数表示か頻度表示か
そのほかに、SNPのアレルの表示の仕方の指定をする

<例１>Phase未知のSNPデータが得られており、SNPのアレルを"0""1"で表現している場合

ジェノタイプは"00""01""11"の３通りで、不明コールを"99"で表すものとする。

今、SampleAの第1番検体の検体IDがA_1,その５SNPのジェノタイプが"00","01","00","11","99"だったとする。その個人のジェノタイプ情報は２行に渡って記載される。Phaseは不明なので、便宜的に１行目と2行目に個々のSNPのつ２アレルが割り振られる

集団遺伝学的アプローチであると、集団中に認められる同一ジェノタイプ保有個体・同一ハプロタイプ本数は集計して提示することがデフォルトとして適当であるが、連鎖不平衡マッピングの立場からは、個々の情報をそのまま扱うことの方が適当なことも多い。そのような場合には、個人のジェノタイプをデータとして示しつつ、そのようなデータの検体を１個あるとした上で、同一パターンのデータを複数回記載するのも手である。以下はそのようにした場合である。左端にデータ名、ついでその観測数、そしてその後にSNPのジェノタイプが２行に渡って記載されている



A_1 1  00019

       01019

A_2 1  00019

       01019

A_3 1  10010

       11010

A_4 1  00011

       01011

A_5 1  00010

       01011

A_6 1  00000

       00010

これを10人分にすると



A_1 1  00019

       01019

A_2 1  00019

       01019

A_3 1  10010

       11010

A_4 1  00011

       01011

A_5 1  00010

       01011

A_6 1  00000

       00010

B_1 1  00011

       01011

B_2 1  09000

       09010

B_3 1  00011

       01011

B_4 1  10001

       10111

ただし、SampleAとSampleBとの区別を入れないといけないので、Sample単位で名称を与え、その人数情報を付加した上で、帰属サンプルの情報を{}でくくる。また、検体の遺伝型データのほかにその書式を指定する欄も付加する必要があるので、それと区別するために、これらの冒頭に"[Data] Samples?"を追加する。また、SampleAとSampleBとをまとめて解析するために、グループ化するために、末尾にStructureと題して始まる記述をしSampleA SampleBを名称"Test"のGroupに属することを宣言する。以下のような記述となる



[Data]

[[Samples]]

SampleName="SampleA"

SampleSize=6

SampleData={

A_1 1  00019

       01019

A_2 1  00019

       01019

A_3 1  10010

       11010

A_4 1  00011

       01011

A_5 1  00010

       01011

A_6 1  00000

       00010

}

SampleName="SampleB"

SampleSize=4

SampleData={

B_1 1  00011

       01011

B_2 1  09000

       09010

B_3 1  00011

       01011

B_4 1  10001

       10111

}

[[Structure]]
StructureName = "Test"

NbGroups = 1

Group =  {   

"SampleA"

"SampleB"

}

これで全データが記載された。あとは、その書式をアプリケーションに提示する部分の記述である。書式についての情報は"[Profile]"で始まる部分に書き、タイトル名(Title)、サンプル数(NbSamples)、Diploidデータかhaploidデータか(GenotypicData)、Diploidの場合には、Phase未知か既知か(GameticPhase)、不明アレルをどう表しているか(MissingData)、同一データパターンの数を観測数で表すか頻度で表すか(Frequency)、ローカス間の区切りは何か(LocusSeparator)、SNPアレルの表示の仕方(DataType)を指定する。

１サンプルがDiploidデータなので、GenotypicData=1(ハプロイドデータの場合には0)、ハプロタイプ推定をしていないので、GameticPhase=0(個人の２本のハプロタイプデータになっていれば1)、SNPとSNPとの間になにも区切り文字が入っていないので、LocusSeparator=NONE(タブ区切りならTAB,半角スペースならWHITESPACE)、アレルは0,1なので、下の例ではDataType=RFLPとしたが、STANDARDとしても動く。STANDARDはもっと一般的なアレル名(HLAのDNAベース表記(DR*0401など)の文字列も受け付けるオプションである。もちろんこの場合には、ローカス間に区切り文字を入れる必要がある。塩基（ATGCにするならば、DNA、マイクロサテライトの繰り返し数にする場合には、MICROSAT）。

上記の記載パターンに即した[Profile]の書き方は次のとおり



[Profile]

Title="Sample1"

NbSamples=2

GenotypicData=1

GameticPhase=0

MissingData="9"

DataType=RFLP

LocusSeparator=NONE


[Data]

[[Samples]]

SampleName="GroupA"

SampleSize=6

SampleData={

A_1 1  00019

       01019

A_2 1  00019

       01019

A_3 1  10010

       11010

A_4 1  00011

       01011

A_5 1  00010

       01011

A_6 1  00000

       00010

}

SampleName="GroupB"

SampleSize=4

SampleData={

B_1 1  00011

       01011

B_2 1  09000

       09010

B_3 1  00011

       01011

B_4 1  10001

       10111

}

[[Structure]]
StructureName = "Test"

NbGroups = 1

Group =  {   

"GroupA"

"GroupB"

}

上記のファイルを"hoge.arq"として保存する。第１　インストールと起動、設定に記載した方法で実行(デフォルト設定)するとブラウザが立ち上がり、以下の出力がなされる。実行設定の出力である

////////////////////////////////////////////////////////////////////
RUN NUMBER 1 (15/12/05 at 16:09:33)
////////////////////////////////////////////////////////////////////

Project information:

2005-12-16

第６限　Neutralityテスト

集団遺伝学実習 Arlequin

分子進化の中立仮説にもとづいたモデルからの予想値と実測データの乖離の程度を評価するのが、Neutralityテストである。テストは、実測データのあるものの分布が、モデルにからの予想の分布との差について検定することにより行われる。アレルに関する分布を用いる場合と、アレルを構成する個々の多型箇所に関する分布を用いる場合との２種類がある。アレルに関する分布を与えるモデルは、Infinite-alleles model(記事はこちら、また、アレル単位での多型・遺伝的多様性の指標であるHomozygosity,Heterozygosityに関する記事はこちら) であり、アレルを構成する個々の多型箇所およびそれらの相互関係についての分布を与えるモデルが、Infinite-sites model(記事はこちら)である。
Arlequinが提供するNeutralityテストは、次の４つ。アレルベースのものと個々の多型箇所ベースのものとに分類される
- Infinite-alleles modelベースのテスト
  - Ewens-Watterson neutrality tests
  - Chakraborty's test of population amalgamation
- Infinite-sites modelベースのテスト
  - Tajima's test
  - Fu's Fs test
Infinite-alleles modelベースのテスト
- Ewens-Watterson neutrality test
  - Infinite-alleles modelにおいては、変異の新規発生とその遺伝的浮動との影響から、複数のアレルが異なるアレル頻度を有する状態で定常になる。そのようなときにサンプリングし、そのアレル数を観測すると、観測アレル数は、Theta( $4N￥mu$ )とサンプル数とで決まる分布を取ることが知られている。観測データのモデル適合度をPermutationベースで評価したものがEwens-Watterson neutrality testである
- Chakraborty's test of population agalgamation
  - 均一な集団に認められるアレル数と不均一集団に認められるアレル数では、後者の方が多くなる。このことから、観測データのHomozygosityから推定されるTheta(Theta(Hom))を用いて、均一集団仮説からの逸脱の程度を検定したものである
Infinite-sites modelベースのテスト
- Tajima's test
  - 第３限で示したとおり、Thetaの推定には、複数の方法がある。それぞれのTheta算出方法には、それぞれモデル・仮定がある。今、観測データに対してThetaを推定すると複数のThetaの値は(普通)一致しないが、その理由には、それぞれのThetaの背景となっているモデルの違いがある。その点を利用して、ThetaS(対象領域の多型箇所数から算出)とThetaPi(染色体ペアの異塩基数から算出)との異同をもとにNeutralityからのずれを評価したもの
- Fu's Fs
  - Tajima's test と同様に観測データから得られるThetaPiを用いて、そのThetaの値の下で、集団に認められるべきアレル数と実際に観測されたアレル数との差異をもとにNeutralityを評価する
出力例
- 第３限に用いたSampleAについて上記４種類のNeutralityを選択して実行した結果を以下に示す。
- Tajima's testでは、Molecular diversityオプションで計算させたTheta(Pi)と同じ値がMean No. of pairwise difference (Pi)として再掲されているともに、Thesa(S)と同じ値がObs. Theta(S)に再掲されている。Theta(S)の値をPermutationしていることが示されている。
  - Simulationから算出したP値は、P(D simul < D obs)に記載されている



==========================================

== Tajima's test of selective neutrality : (SampleA)

==========================================
Reference: Tajima, F. 1989a.

           Tajima, F., 1996.

Sample size                           :      204

No. of sites with substitutions (S)   :        6

Mean No. of pairwise differences (Pi) :  1.17782

Distance method                       : Pairwise difference (no Gamma correction, indels not taken into account)
Tajima's D                            :  0.31759
P(D random < D obs)                   : -0.35930   (Beta distribution aproximation)
No. of simulations                    :     1000

Obs. Theta(S)                         :  1.01818

Mean Theta(S)                         :  0.83134

S.D. Theta(S)                         :  0.23665

Mean D                                :  0.16136

S.D. D                                :  1.09841
P(D simul < D obs)                    :  0.60800

Ewens-Watterson testsではObserved F valueとしてStandard diversity indices のSum of square freqs.が示されている
- 検定結果としてのP値が Slatkin's Exact P.Value で示されている



==================================================

== Ewens-Watterson tests of selective neutrality : (SampleA)

==================================================
Reference: Ewens, W.J. 1972.

           Watterson, G., 1975.

           Stewart, F. M. 1977.

           Slatkin, M. 1994b.

           Slatkin , M., 1996.

Original haplotype definition was used for the tests
No. of genes in sample       :   204.00000

No. of haplotypes in sample  :           7

Observed F value             :     0.32704

Expected F value             :     0.42068

No. of simulated samples     :        1000

Watterson F P.Value          :     0.32600

Slatkin's Exact P.Value      :     0.33400

2005-12-16

第３限　配列の違いの評価(Diversity indices)

集団遺伝学実習 Arlequin

サンプルデータ
- ２群(GroupA,GroupB)について、全６SNPが作るハプロタイプが観測された。GroupAでは102人=204本、GroupBでは69人=138本。GroupAには７種類のハプロタイプが認められ、A1...A7と名前をつけた。それぞれの観測本数は、98,52,...3,1本ずつである。入力ファイル(xxx.arqは以下のとおり。第２限のxxxx.arqファイルと次の点が異なることに留意せよ。GenotypicDataオプションが0(haplotypic)になり、GameticPhaseオプションは削られ、DataTypeがDNAと変更になっている。不明アレルは今回のデータにはないが、実行オプションとしてはMissingData="N"と変更してある



[Profile]

Title="Hap1"

NbSamples=2

GenotypicData=0

MissingData="N"

DataType=DNA

LocusSeparator=NONE


[Data]

[[Samples]]

SampleName="SampleA"

SampleSize=204

SampleData={

A1  98 CTTGGA

A2  52 GTTGGG

A3  34 GTTGGA

A4  10 GTTTGG

A5  6 GATGAA

A6  3 GTCGGA

A7  1 CATGGG

}

SampleName="SampleB"

SampleSize=138

SampleData={

B1  55 CTTGGA

B2  43 GTTGGG

B3  24 GTTGGA

B4  10 GTCGGA

B5  6 GATGAA

}

[[Structure]]
StructureName = "Test"

NbGroups = 1

Group =  {   

"SampleA"

"SampleB"

}

Diversity indicesの実行プログラムの設定
- このカテゴリのすべての解析を実行することとする
  - Molecular diversity
    - Standard diversity indices,Molecular diversity,Compute minimum spanning network among haplotypes,Print distance matrix,Theta(Hom),Theta(S),Theta(K),Theta(Pi)をチェックし、Molecular distanceはデフォルトのPairwise differecen,Gamma a:0.0とする
  - Mismatch distributionにチェックを要れ、Molecular distanceはデフォルトのPairwise difference, Number of bootstrap replicates: もデフォルトの100とする
  - Haplotype frequencies
    - Gene frequency estimation,Estimate allele frequencies at all loci,Search for shared haplotypes between populationsにチェックを入れる
実行
- Runボタンを押すとウィンドウのMessage欄に解析の進行を示す文字列が現れ、終了すると、ブラウザが立ち上がって結果が表示される。結果は別記事(こちら)に示すとおり。上から順番に概説をする
Settings used for Calculations
- 実行オプションを記録した部分
Checking for haplotypes shared among populations:
- 今回の例では２サンプル(SampleAとSampleB)で相互に一致するハプロタイプを提示し、その本数と頻度とを出力
これ以降はSampleAとSampleBについて別個に解析した結果が、SampleAについての全結果,SampleBについての全結果の順番で記載されている
Standard diversity indices
- 対象領域の多様性の強さを0から1の値で示している
- Gene diversityがその解(その他の出力はそのための計算過程で必要とする数値)
  - SampleAでは７ハプロタイプが観測されている。今、１ハプロタイプしかないとするとこのハプロタイプ領域のdiversityは０(まったく多様性がない)。他方、無限大のハプロタイプ数が存在していれば、diversityは最大値の１で与えられる。また、ハプロタイプの種類数が同じく７だとしても７種類のハプロタイプのうち１つだけが極端に頻度が高く、残り６種類の頻度はほとんどゼロだとした場合には、diversityは、ただ１種類のハプロタイプが観測された場合に近く、７種類のハプロタイプ頻度が等しいときに、diversityの値が大きくなるような指標である。実際には、ハプロタイプをサンプリングした母集団においてランダムメイティングを仮定したときに、ヘテロ個体の普遍推定頻度をもって、diversityとしている $Gene diversity = ￥frac{n}{n-1}(1-Sum of square freqs.)$
- Homozygosity Heterozygosityについての記事はこちら
Molecular diversity indices
- 観測配列集合の遠近関係を示している
  - ７ハプロタイプが作るハプロタイプペア 21通り(7x6/2)について、遠近関係を数値化し、それに対して、ある手法を用いて７ハプロタイプ全体の遠近相関関係を提示している
  - Inter-haplotypic distance matrix ハプロタイプペア間の距離の算出(距離行列の作成)
    - 解析オプションにてPairwise difference法を採用した。同法は２ハプロタイプ配列を比較し相互に異なる箇所数を数える方法である。A1=CTTGGAとA2=GTTGGGは第１、第６塩基の２塩基が異なるので、その距離は２。この方法はすべての塩基の違いが点突然変異によるとした場合の距離に相当し、組み換えは考慮されていない。Pairwise difference法以外にも複数のハプロタイプペア間距離計算方法がオプションで選択できるがいずれも、組み換えを想定していない。進化系統樹の解析にあたっては、適当にオプションを変更すること
- MINIMUM SPANNING TREE between 距離行列からハプロタイプ集合の遠近関係の抽出
  - 全ハプロタイプペアについて相互距離の算出がなされたが、それだけでは「親子関係」に類似した関係がわからない。親子関係は「木」で表されるので、「木」を選ぶ必要がある、そのアルゴリズムが、最小木アルゴリズムである(記事は[ http://d.hatena.ne.jp/ryamada22/20051215/1134593910:title=こちら])。簡単に言うと、すべてのハプロタイプを含み、かつ、サイクルを持たない木の形状をしているもののうち、もっとも枝の長さの総和が短いもののことを最小木という。すべてのハプロタイプがつながることを目指し、サイクル(閉路)を許さないアルゴリズムであるので、距離行列では短い距離になっているハプロタイプ同士の距離が結果として長くなる場合もある。たとえば、A2とA7の距離は２だが、下のMINIMUM SPANNING TREEでは、A2-A3-A1-A7とたどる必要があり、その距離は木の上では、４となっている
  - 木は、あるハプロタイプとあるハプロタイプを結ぶエッジ(辺)のありなし、と、ある場合にはその長さがわかればよく、それを表したのが、



   OTU 1     OTU 2    Connection length

   =====     =====    =================

      A1        A3             1.00000

      A3        A2             1.00000

      A2        A4             1.00000

      A3        A6             1.00000

      A3        A5             2.00000

      A1        A7             2.00000

の部分である。また、それを、Newick書式を含むNEXUS書式にしたのが NEXUS notation for MSTの部分である。Newick書式とその描図については、こちらを参照。

MINIMUM SPANNING NETWORK
- MINIMUM SPANNING TREEはサイクルを許さないが、サイクルを許したのが、MINIMUM SPANNING NETWORKである。上の例で、A2-A7間は距離２であったので、A2-A7間に長さ２のエッジを加える。このエッジを加えても、どのハプロタイプペア間のNETWORK上の距離も、距離行列に示された距離より短くなっていないことに留意。新たに生じたパスについて確認するとA3-A7は距離行列上、距離３、TREE上、距離３、NETWORK上(新パス)で距離３、A4-A7も同様に、３、５、３。A6-A7は、４、４、４、である
- これを出力ファイルでは、次のようにして示している



     OTU    List of alternative links

     ===    =========================

      A7       A2 (2.00000)

グラフの用語(ノード・エッジ・木・グラフ・サイクル・パス、などについては、こちらを参照)

Mean number of pairwise differences
- 204染色体が作る204x203/2ペアについて距離の標本平均値と標本標準偏差(同一ハプロタイプ間の距離は0)
- 上記ハプロタイプ距離行列・MINIMUM SPANNING TREE/NETWORKは７種類のハプロタイプについての解析結果であったが、実際には、204本の染色体の情報が７種類のハプロタイプに代表されている。204本のサンプル全体ついての情報は、各々のハプロタイプの数を考慮に入れて解析することが必要。それについての解析結果
Nucleotide diversity
- 同じく、204x203/2ペアについての異なる塩基数を単位距離あたりに直したもの。具体的には、Mean number of pairwise differencesを配列距離(この場合はSNP数)で割った値。これが意味を持つのは、非多型部分も含めて配列を与えた場合
Thetas
- Theta, $￥Theta$ は集団中の平均変異率を反映した数値。集団遺伝学においては、平均変異率uに対し、集団個体数N(diploid)にあって、theta = 4Nu で定義される
Thetaの値の推定方法は複数のものが知られ、Arlequinでは、観測データそのものが与えるホモ個体比率から算出する方法(Theta(Hom))、解析対象範囲に認められた多型箇所数から推定する方法(Theta(S))、解析対象範囲の配列種類数から推定する方法(Theta(k))、Pairwise difference数平均値から推定する方法(Theta(pi))を選べる。Nucleotide diversityと同様、今回の入力データ(SNP位置のみの情報を入力)では不適
それぞれ、標準偏差もしくは信頼区間が与えられる

Mismatch distribution
- 母集団サイズの急増が最近起きたか、起きたとすれば、どのくらいのサイズからどのくらいのサイズへ、どのくらいの世代時間で起きたかを推定し、その母集団サイズ急増モデルへの当て嵌まりのよさを検定したもの。推定の基礎となる原理は、母集団の急増が起きている間には組み換えが起こらないことから、集団中に認められるハプロタイプの違いの分布には、急増の影響が認められるという仮定である。たとえば、時間の経過とともに、変異が蓄積して行く。AAAAAA->AGAAAA->AGGAAA->AGGAAGという変化が起きる一方で、AAAAAA->GAAAAA->GAAGAA->GAAGGAという変化が起きたとする。これに対して、組み換えが起きて出来上がった、７種類のハプロタイプAAAAAA,GAAAAG,AGAAGA,GAAAGG,AGGAAA,GAAGGA,AGGAAGがあるとすると、前者に見られる特徴を持って急増したとみなす。
- 推定値として、Theta0が急増開始前の平衡状態時の集団サイズを反映したTheta値の推定値を、Theta1が急増の結果、大きくなった集団サイズを反映したTheta値の推定値を表す（上記、Thetasのところでも記した式 theta = 4Nu を思い出すこと)。Tauが母集団サイズ急増期間を表す推定値(Tau=2ut; ただしtは時間)。また、急増モデルの当て嵌まりの評価はP(Sim. Ssd>= Obs. Ssd)に表されている。急増によって生じていると考えている観測データの「偏り」がSum of Squared deviationとして計算される。一方、モデルの変数推定値のもとでBootstrapサンプリングを行って計算したSum of Squared deviation値の分布を求め、その分布における観測値の偏りの程度をPとして与える。Sum of Squared deviationのほかに、raggedness indexについて同様にP値化したものも与えられている。出力ファイルの本項の出力は、変数推定とそのBootstrapの過程に関する出力である
Haplotype frequencies estimation
- 今回の入力データはすでにハプロイドデータであり、ここの数値はモーメント法で出力されている(はず)
Thetaとdiversityについては、別記事も参考

2005-12-16

第５限　連鎖不平衡判定と連鎖不平衡係数の計算

集団遺伝学実習 Arlequin

Arlequinの連鎖不平衡解析は遅いので原則、用いない。しかしながら、出力が丁寧なので、原理の学習という意味で連鎖不平衡のみなぞる。

EMアルゴリズムにてハプロタイプ頻度を推定し、それに基づいて連鎖不平衡の検定、および連鎖不平衡係数を計算している
連鎖不平衡に関する記載は、同じくSNPの連鎖不平衡を解析するアプリケーションであるHaploviewの出力との対比も含めて、こちらから
Hardy-Weinberg平衡とその検定
- HHW平衡に関する解説はこちら
- 個々のSNPにおけるHWE検定のほか、Phase既知の場合には、ハプロタイプ単位での検定も行う。仔細に出力を確認していないが、ヘテロ個体の観測期待頻度の値が正しいかどうか、疑義が残っているので、使用を現時点では勧めない。
入力ファイル(２SNP、連鎖不平衡あり)のファイルは末尾の通り



[Profile]

Title="Gen1"

NbSamples=1

GenotypicData=1

GameticPhase=0

MissingData="N"

DataType=DNA

LocusSeparator=NONE


[Data]

[[Samples]]

SampleName="SampleA"

SampleSize=100

SampleData={

A1  80 AA

       AA

A4  5 AA

       CA

A5  5  AA

       CC

A6  4  AC

       CC

A7  4  CA

       CA

A8  2  CA

       CC

}

2005-12-16

第７限　Genetic structure

集団遺伝学実習 Arlequin 構造化

集団は、均一集団の寄せ集めになっている(ことが多い)。そのよせ集まり具合の評価。最近、論文でよく使われている"Structure"はこちら。Structureでは、観測データから個人をグループ分けして、グループ間の遠近関係を出す。こちらの手法は、観測データとともに、個体の所属関係を引数として与えた上で、「所属集団」間の遠近関係を算出する。

実装されている解析は、複数集団が均一か否かを解析するAMOVA、複数集団中の集団ペアについて $F_{ST}$ を算出・検定、また、 $F_{ST}$ にもとづき、集団ペア間の遺伝的距離を推定する方法、それに伴い、集団サイズの推定、集団の分化についての推定・検定など。

記事『Population subdivision と $F_{ST}$ (FST), AMOVA』はこちら。とくにAMOVAについてはこちら。
入力ファイル
- Haplotypic dataの場合
  - Arlequinをダウンロードするとついてくるサンプルデータ\Amova\amovahap.arp　と　そのハプロタイプ定義ファイル \Amova\56hapdef.txt を用いて、その出力を確認する。
  - 入力ファイルには、解析したい構造が記されている。サンプルは本記事の末尾の通り
  - 出力ファイルの該当箇所は次の通り