-------------------------------------------
系統推定に関わる統計的検定について(デモ1)
-------------------------------------------

●ブーツストラップ(Felsenstein 1985)
データからの重複を許す再抽出をするブーツストラップでは,各試行ごとの系統樹探索(search)に時間をかけるのではなく,回数(nreps)を多くとるようにした方がいいでしょう.

Felsenstein, J. 1985. Confidence limits on phylogenies: An approach using the bootstrap. Evolution 39:783-791.

[*** bootstrap nreps=1000 search=heuristic / swap=nni multrees=no ***]
Bootstrap method with heuristic search:
Number of bootstrap replicates = 1000
Starting seed = 499750430
Optimality criterion = parsimony

1000 bootstrap replicates completed
Time used = 2.20 sec

Bootstrap 50% majority-rule consensus tree

/-------------------------------------------------------- Lemur catta(1)
|
| /-------- Homo sapiens(2)
| /--52---+
| | \-------- Pan(3)
| /--100--+
| | \---------------- Gorilla(4)
| /--96---+
| | \------------------------ Pongo(5)
| /--96---+
| | \-------------------------------- Hylobates(6)
| /--97---+
| | | /-------- Macaca fuscata(7)
| | \--------------100--------------+
+--100--+ \-------- M. mulatta(8)
| |
| \------------------------------------------------ Saimiri sciureus(9)
|
\-------------------------------------------------------- Tarsius syrichta(10)

Bipartitions found in one or more trees and frequency of occurrence (bootstrap
support values):

1
1234567890 Freq %
---------------------------
......**.. 1000 100.0%
.***...... 1000 100.0%
.********. 998 99.8%
.*******.. 970 97.0%
.****..... 960 96.0%
.*****.... 958 95.8%
.**....... 523 52.3%
..**...... 460 46.0%

前回の講義でも強調しましたが,上記の合意樹(consensus tree)は系統学的には無意味です.表にされたクレード出現率を元データから推定した系統樹の枝ごとに対応させる必要があります.

●ジャックナイフ(通常)
以下では,全形質の50%を無作為削除する delete-half jackknife を実施します.

[*** jackknife nreps=1000 resample=normal search=heuristic / swap=nni multrees=no ***]
Jackknife method with heuristic search:
Number of jackknife replicates = 1000
Nominal percentage of characters deleted in each replicate = 50
Starting seed = 489163083
Optimality criterion = parsimony

Note: 449 characters are deleted in each replicate; actual deletion
percentage = 50.000

1000 jackknife replicates completed
Time used = 1.16 sec

Jackknife 50% majority-rule consensus tree

/-------------------------------------------------------- Lemur catta(1)
|
| /--------- Homo sapiens(2)
| |
| /---100---+--------- Pan(3)
| | |
| /---97---+ \--------- Gorilla(4)
| | |
| /---97---+ \------------------- Pongo(5)
| | |
| | \---------------------------- Hylobates(6)
| /---95----+
| | | /--------- Macaca fuscata(7)
| | \------------100------------+
+--100---+ \--------- M. mulatta(8)
| |
| \----------------------------------------------- Saimiri sciureus(9)
|
\-------------------------------------------------------- Tarsius syrichta(10)

Bipartitions found in one or more trees and frequency of occurrence (jackknife
support values):

1
1234567890 Freq %
---------------------------
......**.. 1000 100.0%
.***...... 1000 100.0%
.********. 999 99.9%
.****..... 971 97.1%
.*****.... 966 96.6%
.*******.. 951 95.1%
.**....... 497 49.7%
..**...... 491 49.1%


●ジャックナイフ(“Jac”)
各形質が50%の確率で削除される“Jac”オプションでのジャックナイフは,parsimony jackknifing による系統推定(Farris et al. 1996)の方法です.

Farris, J.S., V. Albert, M. Kallersjo, D. Lipscomb, and A.G. Kluge. 1996. Parsimony jackknifing outperforms neighbor-joining. Cladistics 12:99-124.

この方法では,データからの再抽出ごとに抽出される形質数がばらつきます(後述).

[*** jackknife nreps=1000 resample=jac search=heuristic / swap=nni multrees=no ***]
Jackknife method with heuristic search:
Number of jackknife replicates = 1000
Nominal percentage of characters deleted in each replicate = 50
"Jac" resampling method used (actual percentage of characters deleted
varies from replicate to replicate)
Starting seed = 1038032791
Optimality criterion = parsimony

1000 jackknife replicates completed
Time used = 1.15 sec

Jackknife 50% majority-rule consensus tree

/-------------------------------------------------------- Lemur catta(1)
|
| /-------- Homo sapiens(2)
| /--52---+
| | \-------- Pan(3)
| /--100--+
| | \---------------- Gorilla(4)
| /--97---+
| | \------------------------ Pongo(5)
| /--95---+
| | \-------------------------------- Hylobates(6)
| /--96---+
| | | /-------- Macaca fuscata(7)
| | \--------------100--------------+
+--100--+ \-------- M. mulatta(8)
| |
| \------------------------------------------------ Saimiri sciureus(9)
|
\-------------------------------------------------------- Tarsius syrichta(10)

Bipartitions found in one or more trees and frequency of occurrence (jackknife
support values):

1
1234567890 Freq %
---------------------------
......**.. 1000 100.0%
.***...... 999 99.9%
.********. 998 99.8%
.****..... 966 96.6%
.*******.. 958 95.8%
.*****.... 954 95.4%
.**....... 519 51.9%
..**...... 468 46.8%


●PTP検定(Faith & Cranston 1991)
元データの無作為化(randomization)による系統学的情報の検定です.系統樹の全長を検定統計量とし,形質を無作為化することによって全長の帰無分布を構築し,検定を実施します.

Faith, D. P., and P. S. Cranston. 1991. Could a cladogram this short have arisen by chance alone?: On permutation tests for cladistic structure. Cladistics 7:1-28.

[*** permute test=PTP nreps=10000 search=heuristic / swap=nni multrees=no ***]
PTP test with heuristic search:
Starting seed = 803286596
Number of replicates = 10000
All taxa randomized
Optimality criterion = parsimony

10000 permutation test replicates completed
Time used = 26.69 sec

Results of PTP test:

Number of
Tree length replicates
-------------------------
1031* 1
1245 1
1250 1
1251 2
1252 3
1253 2
1254 6
1255 12
1256 8
1257 19
1258 30
1259 27
1260 46
1261 56
1262 93
1263 96
1264 127
1265 154
1266 245
1267 296
1268 294
1269 332
1270 452
1271 490
1272 525
1273 568
1274 595
1275 613
1276 638
1277 619
1278 602
1279 568
1280 479
1281 438
1282 392
1283 305
1284 228
1285 191
1286 142
1287 104
1288 78
1289 50
1290 35
1291 17
1292 8
1293 8
1294 3
1296 1
* = length for original (unpermuted) data

P = 0.000100

●T-PTP検定(Faith 1991)
クレード制約の有無による全長の差を検定統計量とし,上記のPTP検定と同様の形質の無作為化による帰無分布の構築をします.

Faith, D. P. 1991. Cladistic permutation tests for monophyly and nonmonophyly. Systematic Zoology 40:366-375.

[*** permute test=TPTP nreps=10000 search=heuristic / swap=nni multrees=no ***]
T-PTP test for monophyly with heuristic search:
Testing constraint-tree "chg"
Starting seed = 843128832
Number of replicates = 10000
All taxa randomized
Optimality criterion = parsimony

10000 permutation test replicates completed
Time used = 40.92 sec

Results of T-PTP test:

Tree length Number of
difference** replicates
-------------------------
20 1
19* 1
17 3
16 1
15 3
13 4
12 10
11 2
10 13
9 22
8 30
7 38
6 56
5 52
4 86
3 111
2 157
1 181
0 218
-1 284
-2 302
-3 348
-4 411
-5 427
-6 483
-7 451
-8 531
-9 486
-10 458
-11 493
-12 446
-13 461
-14 431
-15 403
-16 345
-17 340
-18 291
-19 268
-20 201
-21 190
-22 162
-23 145
-24 140
-25 104
-26 70
-27 74
-28 55
-29 58
-30 33
-31 24
-32 17
-33 19
-34 21
-35 7
-36 12
-37 2
-38 8
-39 1
-40 1
-41 2
-42 1
-43 1
-44 2
-45 1
-47 1
* = length difference for original (unpermuted) data
** = length of shortest tree not compatible with constraints minus
length of shortest tree compatible with constraints

P = 0.000200

T-PTP検定に対しては論争があることを指摘しておきます:

Swofford, D. L., J. L. Thorne, J. Felsenstein, and B. M. Wiegmann. 1996. The topology-dependent permutation test for monophyly does not test for monophyly. Systematic Biology 45:575-579.

Faith, D. P., and J. W. H. Trueman. 1996. When the topology-dependent permutation test (T-PTP) for monophyly returns significant support for monophyly, should that be equated with (a) rejecting a null hypothesis of nonmonophyly, (b) rejecting a null hypothesis of ''no structure,'' (c) failing to falsify a hypothesis of monophyly, or (d) none of the above? Systematic Biology 45:580-586.

●クレード間差の無作為化検定
T-PTP検定の発展形で,二つのクレード制約に起因する樹長の差異を,無作為化検定します.

[*** permute test=COMPARE2 nreps=10000 search=heuristic / swap=nni multrees=no ***]
Permutation test of two constrained topologies with heuristic search:
Comparing constraint-trees cg and ch
Starting seed = 2136250681
Number of replicates = 10000
All taxa randomized
Optimality criterion = parsimony

10000 permutation test replicates completed
Time used = 39.93 sec

Results of "compare-2" T-PTP test:

Tree length Number of
difference** replicates
-------------------------
35 1
34 1
31 1
30 1
29 3
28 4
27 8
26 11
25 15
24 10
23 20
22 16
21 28
20 35
19 46
18 56
17 83
16 88
15 92
14 133
13 140
12 165
11 221
10 242
9 275
8 308
7 359
6 362
5 397
4 410
3 423
2 448
1 425
0* 460
-1 434
-2 440
-3 445
-4 379
-5 410
-6 349
-7 302
-8 300
-9 246
-10 232
-11 188
-12 196
-13 159
-14 118
-15 107
-16 89
-17 68
-18 52
-19 42
-20 40
-21 29
-22 23
-23 14
-24 13
-25 7
-26 9
-27 11
-28 1
-29 3
-30 1
-31 2
-33 2
-34 2
* = length difference for original (unpermuted) data
** = length of shortest tree compatible with constraint tree "ch"
minus length of shortest tree compatible with constraint tree
"cg"

P = 0.528700

●ILD検定(Incongruence Length Difference Test)
ある形質データの分割された部分が均質であるか否かを無作為分割の反復によって検定する方法である.検定統計量はL[1+2]−(L[1]+L[2])だが,形質全体の樹長L[1+2]は定数であるため,実際にはL[1]+L[2]を検定統計量とし,無作為分割によって構築した帰無分布の上側に棄却域を設定する.

Farris, J.S., M. Kallersjo, A.G. Kluge and C. Bult 1995. Constructing a significance test for incongruence. Systematic Biology, 44(4): 570-572.

[*** HomPart partition=ILDpartition nreps=10000 ***]
Partition-homogeneity test with heuristic search:
Character partition = ILDpartition
Starting seed = 1545647976
Number of replicates = 10000
Optimality criterion = parsimony

10000 partition-homogeneity test replicates completed
Time used = 21.31 sec

Results of partition-homogeneity test:

Sum of Number of
tree lengths replicates
-----------------------------
1018 1
1019 1
1020 7
1021 7
1022 20
1023 50
1024 106
1025 294
1026 579
1027 1026
1028 1680
1029 2489
1030 2321
1031* 1145
1032 140
1033 59
1034 39
1035 14
1036 13
1037 4
1038 3
1039 2
* = sum of lengths for original partition

P value = 1 - (274/10000) = 0.972600

---
なお,上記では説明のため,すべて最節約法のもとでの検定を取り上げた.しかし,目的関数を変更すれば(たとえば尤度関数),上の検定法の多くはそのまま利用できるだろう.


--------------------------------------
●最節約ジャックナイフ法についての補足
--------------------------------------

Farris, J.S. et al. 1996. Parsimony jackknifing outperforms neghbor-joining. Cladistics, 12(2): 99-124.

上記論文は、ブーツストラップではなく、ジャックナイフを用いた方が妥当な系統樹が構築でき、しかも高速で知られている近隣結合法(NJ)よりも巨大データに対してはさらに高速な計算が実現すると主張した論文です。

あるクレードCを支持する形質が元データ行列(全n個)の中にr個あったと仮定します。ブーツストラップは重複を許す形質の無作為再抽出であるため、その形質が抽出される確率はr/n、抽出されない確率は1−r/nです。したがって、あるブーツストラップ反復(全n個)の行列の中にクレードCを支持する形質がまったく含まれない確率は(r/n)^n、Cを支持する形質が少なくともひとつ含まれるという余事象の確率は1−(r/n)^nとなります。n→∞のとき、1−(r/n)^n→1−e^(-r)という有限確定の極限値が得られます。この値は必ず1よりも小さくなります。

ブーツストラップの問題点は、再抽出の確率が形質の総数nに依存する(r/n)という点にあります。したがって、たとえクレードCを支持する形質が同数であったとしても、系統学的情報を持たないノイズ形質(不変形質とか固有派生形質)が多くなれば、支持形質を抽出する確率は小さくなり、結果としてブーツストラップP値は小さくなってしまいます。

Farris et al. は上記論文の中で「Jac」という形質データのジャックナイフ再抽出法を提案しています。この方法では、元データの各形質は全形質数nとは無関係につねにある定数p(<1)の確率で無作為に除去されます。p値が小さいほどジャックナイフ再抽出によるデータ行列の形質数は小さくなると期待され、逆にp値が大きければデータ行列のサイズの期待値は大きくなります。

あるクレードCを支持する形質がn個中r個あったと仮定すると、「Jac」により無作為再抽出されたデータ行列がそのr個をまったく含まない確率はp^rとなります。したがってCの支持形質がジャックナイフ再抽出されたデータ行列の中に少なくともひとつ含まれる確率は1−p^rとなります。この確率が全形質数nとはまったく無関係であることに注意して下さい。ノイズ形質による系統学的情報が攪乱される心配はなくなります。

通常のジャックナイフ法では抽出する形質の個数を固定します。しかし、「Jac」によるジャックナイフ再抽出の特徴は、試行ごとに再抽出データ行列の形質数がばらつくという大きな特徴があります。

「Jac」によるジャックナイフ再抽出によるデータサイズ(形質数)は、二項分布に
従う確率分布をします。各形質の再抽出確率が1−pですので、元データの全形質数
をnとすると、一回のジャックナイフ試行によって抽出された形質数Xは二項分布B
(n,1−p)に従う確率変数となります:

Prob[X=x]=nCx・(1−p)^x・p^(1-x)
[nCxはn個からr個をとる組合せの場合の数]

この式で確率変数Xの範囲はx=0,1,...,nとなります。

このことから、「Jac」を用いたジャックナイフ再抽出における形質数Xの期待値と
分散は、それぞれ:

期待値 E[X]=n(1−p)
分散 V[X]=np(1−p)

となります。たとえば、p=0.5の場合、再抽出形質数は通常のdelete-halfのジャッ
クナイフと同じn/2を期待値としてもちますが、その両側に分散n/4をもってばら
つくことになります。

再抽出形質数と同様に、あるクレードCを支持するr個の形質もまた別の二項分布B
(r,1−p)に従う再抽出数Yを示すことになります:

Prob[Y=y]=rCy・(1−p)^y・p^(1-y)

Yの期待値と分散は次の通り:

期待値 E[Y]=r(1−p)
分散 V[Y]=rp(1−p)

ここで、Yの確率分布は元データの総数nとは無関係であることに注意して下さい。