WO2017122790A1

WO2017122790A1 - 情報処理方法及び装置

Info

Publication number: WO2017122790A1
Application number: PCT/JP2017/001056
Authority: WO
Inventors: 正朗長崎; 嵩矩長谷川; 要小島; 洋介河合
Original assignee: Tohoku University
Current assignee: Tohoku University
Priority date: 2016-01-15
Filing date: 2017-01-13
Publication date: 2017-07-20
Anticipated expiration: 2018-07-15
Also published as: JP2019053340A

Abstract

本方法は、ゲノム中の指定部分に存在する一塩基多型の集合についてジェノタイプとフェノタイプとの関係を表す統計量を観測データに基づいて算出する工程と、観測データにおいてフェノタイプの値を入れ替えた場合における統計量ｓを算出する処理を、算出した統計量が統計量ｓより高くなる回数の割合が第１の閾値を超えるか若しくは第２の閾値を下回るという事象を検出するまで、又は上記事象が検出されなければ予め指定された回数実行するまで、実行する処理実行工程と、上記割合を一塩基多型の集合に対するｐ値として出力する工程とを含む。ここで第１の閾値及び第２の閾値は、統計量ｓより高い統計量が算出されるケースの出現分布及び予め指定された信頼度に基づき、予め指定された有意水準を表す値から動的に決定されるものである。

Description

情報処理方法及び装置

　本発明は、Gene-based Testに関する。

　ＧＷＡＳ（Genome-Wide Association Studies）解析は、フェノタイプと変異の関係を調べる方法で、ＳＮＰ（一塩基多型。Single-Nucleotide Polymorphism）を中心に調べるSNP-based testと、特定遺伝子上に関連するＳＮＰが集合しているという仮定の下でフェノタイプ（形質）とＳＮＰ集合の関連性を検定するGene-based testとがある。

　SNP-based testの場合、特定のフェノタイプとＳＮＰとの関連性を検出するためには、そのＳＮＰが集団中で共通であることが要件となる。しかしながら、集団中のレアなＳＮＰに関しては、実際に関連性があったとしても、統計的に有意な水準で関連性を検出することが出来ない。そのため、特定のフェノタイプとレアなＳＮＰとの関連性を明確に検出するためには、遺伝子上に存在するレアなＳＮＰを集合として観測することでフェノタイプとの関係を見るGene-based testで行うことがよいとされている。

　しかしながら、既存のGene-based testに用いられる統計手法を用いて正確なｐ値（帰無仮説が正しいという条件の下で、実際のデータから計算された統計量よりも大きな統計量が観測される確率）を得るためには、permutation test（並べ替え検定とも呼ぶ）を行うことになる。このpermutation testには、多大な計算コストを要することが分かっている。

　より具体的には、全ゲノム中の遺伝子部分（エクソン）に存在するレアバリアント（レアなＳＮＰ）とフェノタイプの関係を見ようとすると、通常の統計学的な立場からは、そのレアバリアントとフェノタイプとの関係を統計的に処理した場合に、所定の閾値（例えば０．０１や０．０５）よりも小さいｐ値を有するかどうかを判断することになる。

　このｐ値を求める統計的処理を全ゲノムの中に存在する遺伝子部分毎に行うことになり、この遺伝子は全ゲノム中に２００００個ほど存在している。そのため、多重検定という概念を用いると、遺伝子毎の統計学的処理によりその関係性の確からしさを確認するためには、所定の閾値／遺伝子数を基準とすることになる。そうすると、Gene-based testを行う場合のようにpermutation testを実行する場合、有意であるかどうかを見極めるために要する試行回数は最低でも１／（所定の閾値／遺伝子数)となる。

　例えば、所定の閾値＝０．０１、遺伝子数＝２００００であれば、permutation testをそれぞれの遺伝子に対して少なくとも２００万回行う必要があり、遺伝子全体ではその２００００倍である４００億回行う必要がある。さらに、複数のフェノタイプについて関連性を計算する場合には、その分permutation testの試行回数は増加することになる。

　このため、permutation testの試行回数を削減するための技術が幾つか提案されている。

　例えば、以下に示す非特許文献１では、adaptive permutation testとして、最終的に計算されるｐ値が有意ではないと分かったところでテストを中止することで試行回数を削減する技術が開示されている。

　さらに、以下に示す非特許文献２では、Ｍ回実行するpermutation testの途中におけるｐ値が、閾値αを超える場合にはテストを途中で打ち切ること、また、ｐ値が、閾値α＋（ｐ値の分布の標準偏差σ_Mの６倍）を超えるか又は閾値α－（ｐ値の分布の標準偏差σ_Mの６倍）を下回る場合にテストを途中で打ち切ることが開示されている。

　これらの技術でもある程度の高速化がなされるが、特に非特許文献２では、閾値αを中心として「標準偏差σ_Mの６倍」の２倍という固定的な範囲をテスト継続の条件としているが、このような範囲が適切な範囲であるかどうかは不明である。よって、無駄にテストを継続してしまう若しくは中断してしまう場合もある。

PLINK/SEQ, PSEQ association analysis [online] [2015年10月1日検索] インターネット＜https://atgu.mgh.harvard.edu/plinkseq/assoc.shtml＞ SNP&Variation Suite v8.4.1. Manual Collapsing Methods for DNA-Sequence Analysis [Online] [2015年10月16日検索]　＜http://doc.goldenhelix.com/SVS/latest/svsmanual/collapsing_methods.html＞

　従って、本発明の目的は、一側面によれば、ユーザ指定の信頼度を反映した精度を実現しつつpermutation testを高速化させる技術を提供することである。

　本発明に係る方法は、（Ａ）ゲノム中の指定部分に存在する一塩基多型の集合についてジェノタイプとフェノタイプとの関係を表す統計量を観測データに基づいて算出するステップと、（Ｂ）観測データにおいてフェノタイプの値を入れ替えた場合における統計量を算出する処理を、算出した統計量が観測データに基づき算出した統計量より高くなる回数の割合が第１の閾値を超えるか若しくは第２の閾値を下回るという事象を検出するまで、又は上記事象が検出されなければ予め指定された回数実行するまで、実行する処理実行ステップと、（Ｃ）上記割合を一塩基多型の集合に対するｐ値として出力するステップと、を含む。そして、第１の閾値及び第２の閾値は、観測データに基づき算出した統計量より高い統計量が算出されるケースの出現分布及び予め指定された信頼度に基づき、予め指定された有意水準を表す値から動的に決定されるものである。

図１は、実施の形態に係る情報処理装置の構成例を示す図である。図２は、本実施の形態に係る処理フローを示す図である。図３は、観測データの一例を示す図である。図４は、本実施の形態に係る処理フローを示す図である。図５は、置換処理の処理フローを示す図である。図６は、コンピュータの機能ブロック図である。

　本発明の一実施の形態に係る情報処理装置の一例を図１に示す。

　本実施の形態に係る情報処理装置１００は、入力部１１０と、観測データ格納部１２０と、第１データ格納部１３０と、ｐ値算出部１４０と、第２データ格納部１５０と、抽出部１６０と、出力部１７０とを有する。

　入力部１１０は、ユーザから、各種パラメータ、処理すべき観測データ及びＳＮＰの集合（ＳＮＰｓｅｔと呼ぶ）の指定などを受け付け、第１データ格納部１３０などに格納する。入力部１１０は、情報処理装置１００がネットワークを介して接続されている他のコンピュータから、それらのデータを受信する場合もある。観測データ格納部１２０は、複数の被観測者の各々について、フェノタイプのデータと、ゲノムのデータとを対応付けて格納している。これらのデータについては予め格納されている場合もあれば、入力部１１０から入力される場合もある。

　ｐ値算出部１４０は、各ＳＮＰｓｅｔについてｐ値を算出し、算出されたｐ値をＳＮＰｓｅｔに対応付けて第２データ格納部１５０に格納する。このため、ｐ値算出部１４０は、統計量算出部１４１と、置換処理部１４２とを有する。統計量算出部１４１は、予め定められた統計量を算出する。置換処理部１４２は、ブートストラップ法などの方法に従って観測データにおける各被観測者のフェノタイプの値を入れ替える処理等を実行する。

　抽出部１６０は、第２データ格納部１５０から、所定の条件を満たすｐ値が算出されたＳＮＰｓｅｔを抽出し、出力部１７０に出力する。出力部１７０は、ユーザ又は他のコンピュータに対して、抽出部１６０の抽出結果を出力する。なお、出力部１７０は、例えば第１データ格納部１３０に格納されているデータに基づき、第２データ格納部１５０に格納されているデータ等の処理途中のデータを出力するようにしても良い。

　次に、図２乃至図５を用いて、このような情報処理装置１００が行う処理について詳細に説明する。

　入力部１１０は、例えばユーザから、各種パラメータの入力を受け付け、第１データ格納部１３０に格納する（図２：ステップＳ１）。本実施の形態では、１つのＳＮＰｓｅｔについて統計量の最大計算回数Ｂと、第１の信頼度α1と、第２の信頼度α2と、以下で説明する置換処理において算出する統計量の計算回数Ｔと、第１の信頼度α1に対応する有意水準ｐ_sigとを入力する。なお、ｐ_sigについては、多重検定を考慮して設定される。すなわち、ｐ_sig＝（１－α1）／ＳＮＰｓｅｔ数となる。また、Ｂは、少なくとも１／ｐ_sigとなる。さらに、第２の信頼度α2に対応する有意水準ｐ_a2を入力するようにしても良い。また、第１の信頼度α1と第２の信頼度α2については同一である場合もある。

　また、入力部１１０は、例えばユーザから、処理すべき観測データ及びＳＮＰｓｅｔの指定を受け付け、第１データ格納部１３０又は観測データ格納部１２０に格納する（ステップＳ３）。例えば、観測データ格納部１２０において、どこからどこまでの観測データを用いるか、例えば遺伝子の位置及びその遺伝子におけるＳＮＰの位置等のデータを指定する。

　観測データ格納部１２０に格納される観測データの一例を図３に示す。図３では、被観測者毎に、フェノタイプ（ここではｘ及びｙ。例えば身長体重など。）と、ゲノムのデータとが対応付けられて格納されている。例えば、ゲノムには、複数の遺伝子が含まれているが、ここでは遺伝子Ｈに着目しており、その中で変異がある複数の遺伝子座を処理すべきＳＮＰｓｅｔとして指定する。図３の例では変異がある２つのＳＮＰが指定される。なお、フェノタイプについても、多数の種類のフェノタイプが観測データに含まれる場合には、いずれのフェノタイプについて処理するのかについて指定する場合もある。

　さらに、遺伝子毎にＳＮＰｓｅｔを指定するのではなく、ある程度の遺伝子のまとまりについてＳＮＰｓｅｔを指定するようにしても良い。また、ここでは、複数のＳＮＰｓｅｔが指定されるものとする。

　そうすると、ｐ値算出部１４０は、処理すべきＳＮＰｓｅｔのうち、未処理のＳＮＰｓｅｔを１つ特定する（ステップＳ５）。例えば、処理すべきＳＮＰｓｅｔのリストを用意して上から順番に特定する。

　また、ｐ値算出部１４０の統計量算出部１４１は、特定されたＳＮＰｓｅｔについて、観測データそれ自身に関する所定の統計量ｓを算出する（ステップＳ７）。本実施の形態に係る統計量ｓは、ゲノム中の指定遺伝子等に存在する一塩基多型の集合に関してジェノタイプとフェノタイプとの関係を表す統計量である。より具体的には、例えばＳＫＡＴ（Sequential Kernel Association Test）の統計量ｓを算出することとする。

　ＳＫＡＴの統計量ｓは、以下のように表される。
ｓ＝（ｙ－μ）^TＧＷＧ^T（ｙ－μ）　　（１）
　ここでｎを被観測者数、ｍをＳＮＰｓｅｔに含まれるＳＮＰの数として、ｙは、観測されたフェノタイプのｎ次元ベクトルを表し、μは、帰無仮説（すなわち着目しているフェノタイプはジェノタイプとの関係を有さないという仮説）の下におけるフェノタイプの予測平均のｎ次元ベクトルである。また、行列Ｇは、以下のように表される。
Ｇ＝（ｇ₁，・・・，ｇ_i，・・・ｇ_m）^T　　（２）
　ｇ_iは、ｉ番目のＳＮＰについてｎ人の被観測者のジェノタイプを含むｎ次元のベクトルである。なお、ｇ_iは、離散値（例えば０，１，２）でも連続値（例えば０以上２以下）でも良い。
　行列Ｗは、以下のように表される。
Ｗ＝diag（ｗ₁，・・・，ｗ_j，・・・ｗ_m）　（３）
　Ｗは、ｍ×ｍの対角行列であり、ｗ_jは、ｊ番目のＳＮＰについての重み値を表している。

　本実施の形態で採用可能な統計量は、ＳＫＡＴの統計量だけに限定されず、例えば、Burden Testの統計量やc-alpha Test統計量などを用いることができる。これらの統計量については、例えば、Seunggeng Lee et al, "Rare-Variant Association Analysis: Study Designs and Statistical Tests", The American journal of Human Genetics 95, 5-23, July 3, 2014、Benjamin et al. “Testing for an Unusual Distribution of Rare Variants” PLoS Genet 7(3), 2011などを参照のこと。これらの文献を本明細書に取り込む。

　次に、ｐ値算出部１４０は、集合Ｓ_BcにＮＵＬＬを代入し、統計量の計算回数を表す変数Ｂ_cに０を代入する（ステップＳ９）。

　そして、ｐ値算出部１４０は、Ｂ_cがユーザ指定のＢより小さいか否かを判断する（ステップＳ１１）。この条件が満たされる場合には、処理は端子Ａを介して図４の処理に移行する。

　一方、ステップＳ１１の条件が満たされない場合、すなわち計算回数Ｂ_cがユーザ指定の最大計算回数Ｂ以上となった場合には、ｐ値算出部１４０は、ステップＳ５で特定されたＳＮＰｓｅｔについての推定ｐ値として、図４の処理で計算されるｐ_Bcを、当該ＳＮＰｓｅｔのデータ（例えばＳＮＰｓｅｔに含まれるＳＮＰの位置データ）に対応付けて第２データ格納部１５０に格納する（ステップＳ１３）。

　その後、ｐ値算出部１４０は、未処理のＳＮＰｓｅｔが存在するか否かを判断する（ステップＳ１５）。存在する場合には、処理はステップＳ５に戻る。一方、未処理のＳＮＰｓｅｔが存在しない場合には、抽出部１６０は、第２データ格納部１５０において、推定ｐ値がユーザが指定したｐ_sig未満であるＳＮＰｓｅｔを有意水準を超えるＳＮＰｓｅｔとして抽出して、抽出したＳＮＰｓｅｔを出力部１７０に出力する。出力部１７０は、抽出部１６０の抽出結果を、出力装置（例えば表示装置、印刷装置、他のコンピュータなど）に出力する（ステップＳ１７）。

　図２に示した処理は、従来とおおよそ変わらない処理である。処理の高速化については図４及び図５を用いて詳しく述べる。

　端子Ａを介して図２のステップＳ１１から移行してくると、置換処理部１４２は、置換処理を実行する（図４：ステップＳ２１）。この置換処理については、図５を用いて詳しく述べる。

　置換処理部１４２は、カウンタｔを０に初期化し、集合Ｓ_TにＮＵＬＬを代入する（図５：ステップＳ４１）。そして、置換処理部１４２は、観測データ格納部１２０に格納されている観測データにおいてフェノタイプの入れ替えを実行する（ステップＳ４３）。各被観測者について、復元抽出又は非復元抽出によってフェノタイプの値を抽出して、抽出したフェノタイプの値に入れ替えて、フェノタイプ入れ替え後の観測データを生成する。すなわち、処理対象の被観測者のフェノタイプの値を母集団として、母集団から標本を選ぶ場合に１つの個体を選んだ後にその個体を母集団の中に戻す方法が復元抽出であり、一度選んだ個体を母集団の中に戻さずに標本を選ぶ方法が非復元抽出である。図３の例では、被観測者１００１のフェノタイプの値を、被観測者１００３のフェノタイプの値として設定するといった処理を行う。

　そして、置換処理部１４２は、統計量算出部１４１に、フェノタイプ入れ替え後の観測データに関する所定の統計量ｓ_tを算出させる（ステップＳ４５）。統計量算出部１４１は、処理対象の観測データは異なるが、算出方法はステップＳ７で述べた方法と同じである。

　その後、置換処理部１４２は、算出された統計量ｓ_tを、集合Ｓ_Tに入れる（ステップＳ４７）。また、置換処理部１４２は、カウンタｔを１インクリメントする（ステップＳ４９）。

　そして、置換処理部１４２は、カウンタｔの値が、ユーザ指定の定数Ｔを超えたか否かを判断する（ステップＳ５１）。カウンタｔの値が、ユーザ指定の定数Ｔ以下であれば、処理はステップＳ４３に戻る。一方、カウンタｔの値が、ユーザ指定の定数Ｔを超えた場合には、処理は呼び出し元の処理に戻る。なお、Ｔは１であってもよいが、ある程度の数値である方が、処理の効率化が図られる。

　このような処理を実行すれば、Ｔ個の統計量ｓ_tが集合Ｓ_Tに入れられることになる。

　図４の処理の説明に戻って、ｐ値算出部１４０は、集合Ｓ_Bcに、集合Ｓ_Bc及び集合Ｓ_Tを代入する（ステップＳ２３）。また、ｐ値算出部１４０は、変数Ｂ_cに、集合Ｓ_Bcに含まれる統計量の数を代入する（ステップＳ２５）。さらに、ｐ値算出部１４０は、変数ｂに、集合Ｓ_Bc内において統計量ｓよりも高い統計量ｓ_tの数を計数して代入する（ステップＳ２７）。

　そして、ｐ値算出部１４０は、ｂ／Ｂ_cを算出して、推定ｐ値を表す変数ｐ_Bcに代入する（ステップＳ２９）。

　さらに、ｐ値算出部１４０は、ガウス分布Ｎ（０，ｐ_Bc（１－ｐ_Bc）／Ｂ_c）から、信頼度α1の信頼区間の幅ｄ_aを算出する（ステップＳ３１）。

　本実施の形態では、推定ｐ値がユーザが指定したｐ_sig未満（又は以下）であるＳＮＰｓｅｔを抽出することになるので、明らかにｐ_sigを下回るＳＮＰｓｅｔ及び明らかにｐ_sigを上回るようなＳＮＰｓｅｔについては、途中で処理を打ち切っても問題が生じない。しかしながら、どのような条件を満たした場合に、途中で処理を打ち切るかという問題が別に生ずる。

　具体的には、ｐ_sigを中心として一定の幅に入るｐ_Bcが計算された場合には処理を継続し、ｐ_sigを中心として一定の幅に入らないｐ_Bcが計算された場合には処理を途中で打ち切ることとする。そこで、本実施の形態では、ユーザが指定した信頼度α1に応じた精度が出るように、当該信頼度α1に応じた幅を設定する。

　より具体的には、統計量ｓ未満の統計量が算出されると１となり、それ以外であれば０となるバイナリのランダムな変数Ｘを考え、計算された統計量についてのＸの期待値Ｅｘｐ［Ｘ］及び分散Ｖａｒ［Ｘ］は、以下のように表される。
Ｅｘｐ［Ｘ］＝ｐ_Bc＝ｂ／Ｂ_c
Ｖａｒ［Ｘ］＝Ｅｘｐ［Ｘ²］－Ｅｘｐ［Ｘ］²＝ｐ_Bc（１－ｐ_Bc）

　大数の法則から期待値の分散は以下のようになる。
Ｖａｒ［ｐ_Bc］＝ｐ_Bc（１－ｐ_Bc）／Ｂ_c　　（４）

　さらに中心極限定理から、ｐ_Bcはガウス分布Ｎ（０，ｐ_Bc（１－ｐ_Bc）／Ｂ_c）（平均０で、分散は（４）式のガウス分布）に対応するものと考える。そうすると、このガウス分布における信頼度α1に対応する信頼区間の幅ｄ_aを、信頼度α1に応じた幅として用いることが適当であることが分かる。

　ガウス分布Ｎ（０，ｐ_Bc（１－ｐ_Bc）／Ｂ_c）は、統計量ｓ_tの計算回数が増加すれば動的に変化する。そして、この分布Ｎ（０，ｐ_Bc（１－ｐ_Bc）／Ｂ_c）が分かれば、ユーザ指定の信頼度α1に対応する信頼区間を算出する関数は既知であるから、その幅ｄ_aの詳細な計算方法については説明を省略する。但し、一般的に、信頼区間の下限及び上限は、以下のように表される。
下限：標本平均－ｔ×標本標準偏差／（標本の数）^0.5
上限：標本平均＋ｔ×標本標準偏差／（標本の数）^0.5
よって、信頼区間の幅は、ｔ×標本標準偏差／（標本の数）^0.5となる。但し、ｔは、既知のｔ分布表から得られる標本の数－１と信頼度とに対応する値である。

　一方、一時的にｂ＝０又はｂ＝Ｂ_cとなる状態が生じ得る。このような場合には、分散が０になるので、上で述べたような基準にて処理を途中で終了させるのは妥当ではない。そこで、本実施の形態では、負の二項分布を採用する。ここで、Ｙを試行回数を表す正の整数のランダム変数として、ｐ_{sig_single}＝ｐ_sig＊ＳＮＰｓｅｔ数とする。そうすると、真のｐ値は、最も高くてｐ_sigであり、最も低くてｐ_{sig_single}となる。

　よって、ｂ＝０の場合には、Ｂ_cの確率は負の二項分布ＮＢ（Ｙ＝Ｂ_c；Ｂ_c，１－ｐ_sig）と表され、この確率が、最も高くてｐ_Bc＝ｐ_sigで信頼度α2となる有意水準ｐ_a2（＝１－α2）より小さければ、処理を途中で停止して、ｐ_Bc＝１／Ｂ_cと設定する。

　すなわち、ｐ_Bcが０である場合には、以下の事象を検出すると、処理を停止する。
　予め指定された信頼度α2に対応する有意水準ｐ_a2＞（算出した統計量ｓ_tが観測データに基づき算出した統計量ｓより低くなる確率を（１－ｐ_sig）として、Ｂ_c回だけ実行すると、Ｂ_c回だけ、算出した統計量ｓ_tが観測データに基づき算出した統計量ｓより低くなる確率）

　また、ｂ＝Ｂ_cの場合には、Ｂ_cの確率は負の二項分布ＮＢ（Ｙ＝Ｂ_c；Ｂ_c，ｐ_{sig_single}）と表され、この確率が、最も低くてｐ_Bc＝ｐ_{sig_single}で有意水準ｐ_a2より小さければ、処理を途中で停止して、ｐ_Bc＝１と設定する。

　すなわち、ｐ_Bc＝１である場合は、以下の事象を検出すると、処理を停止する。
　予め指定された信頼度α2に対応する有意水準ｐ_a2＞（算出した統計量ｓ_tが観測データに基づき算出した統計量ｓより高くなる確率をｐ_{sig_single}として、Ｂ_c回だけ実行すると、Ｂ_c回だけ、算出した統計量ｓ_tが観測データに基づき算出した統計量ｓより高くなる確率）

　このような考え方から、ｐ値算出部１４０は、ｐ_Bc＝０且つＮＢ（Ｙ＝Ｂ_c；Ｂ_c，１－ｐ_sig）＜ｐ_a2という条件を満たすか否かを判断する（ステップＳ３３）。このような条件を満たす場合には、ｐ_Bc＝１／Ｂ_cとして処理を途中で停止しても良いので、ｐ値算出部１４０は、ｐ_Bcに対して１／Ｂ_cと設定して（ステップＳ３４）、処理は端子Ｃを介して図２のステップＳ１３に戻る。

　一方、ステップＳ３３の条件を満たさない場合には、ｐ値算出部１４０は、ｐ_Bc＝１且つＮＢ（Ｙ＝Ｂ_c；Ｂ_c，ｐ_{sig_single}＝ｐ_sig×ＳＮＰｓｅｔ数）＜ｐ_a2という条件を満たすか否かを判断する（ステップＳ３５）。このような条件を満たす場合には、ｐ_Bc＝１として処理を途中で停止しても良い。よって、処理は端子Ｃを介して図２のステップＳ１３に戻る。

　一方、ステップＳ３５の条件を満たさない場合には、ｐ値算出部１４０は、ｐ_Bc＞ｐ_sig＋ｄ_a／２又はｐ_Bc＜ｐ_sig－ｄ_a／２という条件を満たしているか否かを判断する（ステップＳ３９）。すなわち、ｐ_Bcが、ｐ_sigを中心に、信頼区間の幅ｄ_aに入るか否かを判断している。この条件を満たす場合には、処理を途中で停止して良いので、端子Ｃを介して図２のステップＳ１３に戻る。なお、この条件に含まれる不等式は、その関係を維持しつつ変形することが可能である。一方、この条件を満たさない場合には、処理は端子Ｂを介して図２のステップＳ１１に戻る。

　このように、ｐ_Bcが、ｐ_sigを中心に、ユーザ指定の信頼度が反映されて決定される信頼区間の幅ｄ_aに入る場合には、入念に置換処理を実行して正確なｐ値を算出するようにする。ｄ_aは、統計量ｓ_tの算出状況に応じて変化する分布Ｎ（０，ｐ_Bc（１－ｐ_Bc）／Ｂ_c）に基づき動的に変化する。

　ステップＳ３３及びＳ３５の例外を除けば、ｐ_Bcが、ｐ_sigを中心に、信頼区間の幅ｄ_aに入らなければ、その時点で次のＳＮＰｓｅｔについての処理に移行するので、それだけ置換処理の実行回数を削減できる。従って、全体の処理速度を高速化できるようになる。

　以上述べた実施の形態は一例であって、本発明はこれに限定されるものではない。例えば、図１に示した機能ブロック図は一例であって、プログラムモジュール構成とは一致しない場合もある。処理フローについては、処理結果が変わらない限り、処理順番を入れ替えたり、処理を並列に実行するようにしても良い。本実施の形態の場合、特にＳＮＰｓｅｔ毎に、異なるコンピュータ又はプロセッサコアで並列化しても良い。

　また、重要なＳＮＰｓｅｔについては、ｐ_Bc＜ｐ_sigが成り立つ限り置換処理を繰り返し実行することで正確なｐ値を算出するようにしても良い。重要か否かについては、フェノタイプの入れ替え無しで算出された統計量ｓの値に基づき決定するようにしても良い。特に、フェノタイプの入れ替え無しで算出された統計量ｓの値が、所定の閾値以上となる場合には、本処理の対象外とするようにしても良い。

　さらに、図１の情報処理装置１００のように１台のコンピュータで実装するのではなく、クライアント端末から１又は複数の情報処理装置１００にネットワークを介してアクセスして処理を行わせるようにしても良い。複数の情報処理装置１００を用いる場合には、上で述べた機能を分担すればよい。

　なお、上で述べた情報処理装置１００は、コンピュータ装置であって、図６に示すように、メモリ２５０１とＣＰＵ（Central Processing Unit）群２５０３とハードディスク・ドライブ（ＨＤＤ：Hard Disk Drive）２５０５と表示装置２５０９に接続される表示制御部２５０７とリムーバブル・ディスク２５１１用のドライブ装置２５１３と入力装置２５１５とネットワークに接続するための通信制御部２５１７とがバス２５１９で接続されている。オペレーティング・システム（ＯＳ：Operating System）及び本実施例における処理を実施するためのアプリケーション・プログラムは、ＨＤＤ２５０５に格納されており、ＣＰＵ群２５０３により実行される際にはＨＤＤ２５０５からメモリ２５０１に読み出される。ＣＰＵ群２５０３は、アプリケーション・プログラムの処理内容に応じて表示制御部２５０７、通信制御部２５１７、ドライブ装置２５１３を制御して、所定の動作を行わせる。また、処理途中のデータについては、主としてメモリ２５０１に格納されるが、ＨＤＤ２５０５に格納されるようにしてもよい。本技術の実施例では、上で述べた処理を実施するためのアプリケーション・プログラムはコンピュータ読み取り可能なリムーバブル・ディスク２５１１に格納されて頒布され、ドライブ装置２５１３からＨＤＤ２５０５にインストールされる。インターネットなどのネットワーク及び通信制御部２５１７を経由して、ＨＤＤ２５０５にインストールされる場合もある。このようなコンピュータ装置は、上で述べたＣＰＵ群２５０３、メモリ２５０１などのハードウエアとＯＳ及びアプリケーション・プログラムなどのプログラムとが有機的に協働することにより、上で述べたような各種機能を実現する。

　以上述べた本実施の形態をまとめると、以下のようになる。

　本実施の形態に係る方法（検定方法、統計処理方法、テスト方法等とも呼ぶ）は、（Ａ）ゲノム中の指定部分に存在する一塩基多型の集合についてジェノタイプとフェノタイプとの関係を表す統計量を観測データに基づいて算出するステップと、（Ｂ）観測データにおいてフェノタイプの値を入れ替えた場合における統計量を算出する処理を、算出した統計量が観測データに基づき算出した統計量より高くなる回数の割合が第１の閾値を超えるか若しくは第２の閾値を下回るという事象を検出するまで、又は当該事象が検出されなければ予め指定された回数実行するまで、実行する処理実行ステップと、（Ｃ）上記割合を一塩基多型の集合に対するｐ値として出力するステップとを含む。そして、上で述べた第１の閾値及び第２の閾値は、観測データに基づき算出した統計量より高い統計量が算出されるケースの出現分布及び予め指定された信頼度に基づき、予め指定された有意水準を表す値から動的に決定されるものである。

　このように、予め指定された信頼度を反映させた精度が実現されるように第１及び第２の閾値が設定され、第１の閾値を超えるか又は第２の閾値を下回る場合には処理が途中で打ち切られるので、permutation testの高速化も実現されている。さらに、第１及び第２の閾値は、観測データに基づき算出した統計量より高い統計量が算出されるケースの出現分布にも基づくものであるため、第１及び第２の閾値は統計量の算出状況に応じて動的に変化する。すなわち、本実施の形態は、統計量の算出方法に依存せずに適用可能である。

　なお、上で述べた第１の閾値が、予め指定された有意水準を表す値に対して、上記出現分布において予め指定された信頼度に対応する信頼区間の幅の１／２の値を加算した値であり、上で述べた第２の閾値が、予め指定された有意水準を表す値から、上記出現分布において予め指定された信頼度に対応する信頼区間の幅の１／２の値を引いた値である場合もある。このように、予め指定された有意水準を表す値を中心に、「予め指定された信頼度に対応する信頼区間の幅」を設定することで、精度からして処理継続が適切な一塩基多型の集合を適切に特定できるようになる。

　上で述べた処理実行ステップにおいて、（ｂ１）上記割合が０である場合には、（予め指定された第２の信頼度に対応する第２の有意水準）＞（算出した統計量が観測データに基づき算出した統計量より低くなる確率を（１－予め指定された有意水準）として、上記処理を現在の実行回数だけ実行すると、上記処理の現在の実行回数だけ、算出した統計量が観測データに基づき算出した統計量より低くなる確率）となる第２の事象、又は（ｂ２）上記割合が１である場合には、（上記第２の有意水準）＞（算出した統計量が観測データに基づき算出した統計量より高くなる確率を予め指定された有意水準を表す値に対して処理すべき一塩基多型の集合の数を乗じた値として、上記処理を現在の実行回数だけ実行すると、上記処理の現在の実行回数だけ、算出した統計量が観測データに基づき算出した統計量より高くなる確率）となる第３の事象を検出するまで上記処理を実行するようにしても良い。この場合、上記第２の事象を検出すると、１／処理の現在の実行回数をｐ値として出力し、上記第３の事象を検出すると、１をｐ値として出力するようにしても良い。上記割合が０である場合や１である場合等、特別なケースにも対処できる。

　また、本実施の形態に係る方法は、（Ｄ）上記ｐ値に基づき、一塩基多型の集合の採否を決定するステップをさらに含むようにしても良い。これによって有意な一塩基多型の集合を得ることができるようになる。

　また、上で述べた処理実行ステップが、上記処理を所定回数実行する毎に、上記事象の発生の有無を判断するステップを含むようにしても良い。このようにすれば、より処理の効率化が図られる。

　なお、上で述べたような処理をプロセッサ又はコンピュータに実行させるためのプログラムを作成することができ、当該プログラムは、例えばフレキシブル・ディスク、ＣＤ－ＲＯＭなどの光ディスク、光磁気ディスク、半導体メモリ（例えばＲＯＭ）、ハードディスク等のコンピュータ読み取り可能な記憶媒体又は記憶装置に格納される。なお、処理途中のデータについては、ＲＡＭ等の記憶装置に一時保管される。

Claims

　ゲノム中の指定部分に存在する一塩基多型の集合についてジェノタイプとフェノタイプとの関係を表す統計量を観測データに基づいて算出するステップと、
　前記観測データにおいて前記フェノタイプの値を入れ替えた場合における統計量を算出する処理を、算出した前記統計量が前記観測データに基づき算出した前記統計量より高くなる回数の割合が第１の閾値を超えるか若しくは第２の閾値を下回るという事象を検出するまで、又は前記事象が検出されなければ予め指定された回数実行するまで、実行する処理実行ステップと、
　前記割合を前記一塩基多型の集合に対するｐ値として出力するステップと、
　を、コンピュータに実行させ、
　前記第１の閾値及び前記第２の閾値は、前記観測データに基づき算出した前記統計量より高い統計量が算出されるケースの出現分布及び予め指定された信頼度に基づき、予め指定された有意水準を表す値から動的に決定される
　ことを特徴とするプログラム。
　前記第１の閾値が、前記予め指定された有意水準を表す値に対して、前記出現分布において前記予め指定された信頼度に対応する信頼区間の幅の１／２の値を加算した値であり、
　前記第２の閾値が、前記予め指定された有意水準を表す値から、前記出現分布において前記予め指定された信頼度に対応する信頼区間の幅の１／２の値を引いた値である
　請求項１記載のプログラム。
　前記処理実行ステップにおいて、
　前記割合が０である場合には、（予め指定された第２の信頼度に対応する第２の有意水準）＞（算出した前記統計量が前記観測データに基づき算出した前記統計量より低くなる確率を（１－前記予め指定された有意水準）として、前記処理を現在の実行回数だけ実行すると、前記処理の現在の実行回数だけ、算出した前記統計量が前記観測データに基づき算出した前記統計量より低くなる確率）となる第２の事象、又は
　前記割合が１である場合には、（前記第２の有意水準）＞（算出した前記統計量が前記観測データに基づき算出した前記統計量より高くなる確率を前記予め指定された有意水準を表す値に対して処理すべき前記一塩基多型の集合の数を乗じた値として、前記処理を現在の実行回数だけ実行すると、前記処理の現在の実行回数だけ、算出した前記統計量が前記観測データに基づき算出した前記統計量より高くなる確率）となる第３の事象
　を検出するまで前記処理を実行して、
　前記第２の事象を検出すると、１／前記処理の現在の実行回数をｐ値として出力し、
　前記第３の事象を検出すると、１をｐ値として出力する
　請求項１又は２記載のプログラム。
　前記ｐ値に基づき、前記一塩基多型の集合の採否を決定するステップ
　をさらに前記コンピュータに実行させるための請求項１乃至３のいずれか１つ記載のプログラム。
　前記処理実行ステップが、
　前記処理を所定回数実行する毎に、前記事象の発生の有無を判断するステップを
　含む請求項１乃至４のいずれか１つ記載のプログラム。
　ゲノム中の指定部分に存在する一塩基多型の集合についてジェノタイプとフェノタイプとの関係を表す統計量を観測データに基づいて算出するステップと、
　前記観測データにおいて前記フェノタイプの値を入れ替えた場合における統計量を算出する処理を、算出した前記統計量が前記観測データに基づき算出した前記統計量より高くなる回数の割合が第１の閾値を超えるか若しくは第２の閾値を下回るという事象を検出するまで、又は前記事象が検出されなければ予め指定された回数実行するまで、実行する処理実行ステップと、
　前記割合を前記一塩基多型の集合に対するｐ値として出力するステップと、
　を、コンピュータが実行し、
　前記第１の閾値及び前記第２の閾値は、前記観測データに基づき算出した前記統計量より高い統計量が算出されるケースの出現分布及び予め指定された信頼度に基づき、予め指定された有意水準を表す値から動的に決定される
　ことを特徴とする情報処理方法。
　ゲノム中の指定部分に存在する一塩基多型の集合についてジェノタイプとフェノタイプとの関係を表す統計量を観測データに基づいて算出する算出手段と、
　前記観測データにおいて前記フェノタイプの値を入れ替えた場合における統計量を算出する処理を、算出した前記統計量が前記観測データに基づき算出した前記統計量より高くなる回数の割合が第１の閾値を超えるか若しくは第２の閾値を下回るという事象を検出するまで、又は前記事象が検出されなければ予め指定された回数実行するまで、実行する処理実行手段と、
　前記割合を前記一塩基多型の集合に対するｐ値として出力する手段と、
　を有し、
　前記第１の閾値及び前記第２の閾値は、前記観測データに基づき算出した前記統計量より高い統計量が算出されるケースの出現分布及び予め指定された信頼度に基づき、予め指定された有意水準を表す値から動的に決定される
　ことを特徴とする情報処理装置。