知財求人 - 知財ポータルサイト「IP Force」
(19)【発行国】日本国特許庁(JP)
(12)【公報種別】公開特許公報(A)
(11)【公開番号】P2024141008
(43)【公開日】2024-10-10
(54)【発明の名称】情報処理装置、方法及びプログラム
(51)【国際特許分類】
G16B 20/00 20190101AFI20241003BHJP
G16B 40/20 20190101ALI20241003BHJP
【ＦＩ】
G16B20/00
G16B40/20
【審査請求】未請求
【請求項の数】10
【出願形態】ＯＬ
(21)【出願番号】P 2023052425
(22)【出願日】2023-03-28
(71)【出願人】
【識別番号】000208891
【氏名又は名称】ＫＤＤＩ株式会社
(74)【代理人】
【識別番号】100092772
【弁理士】
【氏名又は名称】阪本 清孝
(74)【代理人】
【識別番号】100119688
【弁理士】
【氏名又は名称】田邉 壽二
(72)【発明者】
【氏名】中村 みゆき
(72)【発明者】
【氏名】安田 圭志
(72)【発明者】
【氏名】永田 雅俊
(57)【要約】
【課題】効果的なエピゲノム解析に寄与することができる情報処理装置を提供する。
【解決手段】ＤＮＡの複数の修飾サイトのそれぞれについて、修飾レベルデータから目的形質を推定する項目応答理論のパラメータを入力として読み込み、目的形質の所与の範囲（Ｓ４１）について、項目応答理論において所与の情報量が確保（Ｓ４３）されるように、前記複数の修飾サイトの中から選別を行う（Ｓ４４，Ｓ４５，Ｓ４６）第１処理を実行する
【選択図】図１１
【解決手段】ＤＮＡの複数の修飾サイトのそれぞれについて、修飾レベルデータから目的形質を推定する項目応答理論のパラメータを入力として読み込み、目的形質の所与の範囲（Ｓ４１）について、項目応答理論において所与の情報量が確保（Ｓ４３）されるように、前記複数の修飾サイトの中から選別を行う（Ｓ４４，Ｓ４５，Ｓ４６）第１処理を実行する
【選択図】図１１
【特許請求の範囲】
【請求項1】
ＤＮＡの複数の修飾サイトのそれぞれについて、修飾レベルデータから目的形質を推定する項目応答理論のパラメータを入力として読み込み、目的形質の所与の範囲について、項目応答理論において所与の情報量が確保されるように、前記複数の修飾サイトの中から選別を行う第１処理を実行することを特徴とする情報処理装置。
【請求項2】
前記第１処理では、前記複数の修飾サイトのそれぞれについて、前記パラメータにより目的形質の関数としてのテスト情報量を求め、当該テスト情報量の和が前記所与の範囲において前記所与の情報量を確保するように、前記選別を行うことを特徴とする請求項１に記載の情報処理装置。
【請求項3】
前記第１処理では、目的形質を第１軸とし、テスト情報量を第２軸とする空間において、前記所与の範囲及び前記所与の情報量による空間範囲を探索範囲として定め、当該探索範囲との重複がより大きくなるテスト情報量に対応する修飾サイトを優先して、前記選別することを特徴とする請求項２に記載の情報処理装置。
【請求項4】
前記第１処理では、前記重複の大きさが上位側に属すると判定される修飾サイトが複数あった場合に、前記所与の範囲内におけるテスト情報量の最大値が大きくなる修飾サイトを優先して、前記選別することを特徴とする請求項３に記載の情報処理装置。
【請求項5】
前記第１処理では、前記重複の大きさが上位側に属すると判定される修飾サイトが複数あった場合に、前記所与の範囲外におけるテスト情報量が大きくなる修飾サイトを優先して、前記選別することを特徴とする請求項３に記載の情報処理装置。
【請求項6】
前記第１処理を、互いに異なる複数の探索範囲について実行し、当該異なる探索範囲ごとに、前記選別された結果を提示することを特徴とする請求項３に記載の情報処理装置。
【請求項7】
前記第１処理では、前記複数の修飾サイトが選別される個数の入力を受け付け、当該個数に渡って前記選別されたテスト情報量の和が前記所与の範囲を確保することをどの程度だけ実現できているかの情報を提示することを特徴とする請求項３に記載の情報処理装置。
【請求項8】
前記項目応答理論のパラメータは項目困難度パラメータ及び項目識別力パラメータを含み、
前記第１処理では、前記複数の修飾サイトのうち、対応する項目困難度パラメータが前記所与の範囲内にあるものを予め候補として絞り込み、当該候補の中から前記選別することを特徴とする請求項１に記載の情報処理装置。
前記第１処理では、前記複数の修飾サイトのうち、対応する項目困難度パラメータが前記所与の範囲内にあるものを予め候補として絞り込み、当該候補の中から前記選別することを特徴とする請求項１に記載の情報処理装置。
【請求項9】
ＤＮＡの複数の修飾サイトのそれぞれについて、修飾レベルデータから目的形質を推定する項目応答理論のパラメータを入力として読み込み、目的形質の所与の範囲について、項目応答理論において所与の情報量が確保されるように 、前記複数の修飾サイトの中から選別を行う第１処理をコンピュータが実行することを特徴とする情報処理方法。
【請求項10】
コンピュータを請求項１ないし８のいずれかに記載の情報処理装置として機能させることを特徴とするプログラム。
【発明の詳細な説明】
【技術分野】
【0001】
本発明は、エピゲノムデータに関する処理を行う情報処理装置、方法及びプログラムに関する。
【背景技術】
【0002】
アデニン（A）、グアニン（G）、チミン（T）、シトシン（C）の４塩基から構成されるDNA（デオキシリボ核酸）はゲノムを構成し、生物の生涯において不変である。一方、DNA塩基配列の内、CとGが連続するCpGサイト（CG配列）と呼ばれる箇所に存在するCに関しては、環境などの影響を受け、メチル化や脱メチル化が生じることが知られ、エピゲノムとしてゲノムの様々な働きに関係している。近年、DNA解析技術の発展により、DNAの塩基配列のみならず、Cにおける、メチル化の有無までも計測可能になっている。
【0003】
すなわち、ヒトを含む哺乳類のゲノムDNAにおけるCG配列はメチル化修飾を受けうる。ゲノムDNAにおける特定の領域における化学修飾DNAメチル化レベルは、加齢や環境刺激により後天的に変化するため、生物学的年齢などの指標の推定に用いられる。図１及び図２はこれらの模式図である。
【0004】
図１はゲノムとエピゲノムとの模式図であり、例EX1に示すようにゲノムは各個体に固有のものとして生涯においてその配列が不変であり、それぞれの個体において現在もX年後も同じであるが、例EX2に示すように、エピゲノムはゲノムにおけるCpGサイトの状態として可変であり、例EX3にも示すように、環境や経年経過で様々な影響を受けて変化しうる。
【0005】
この影響は、図２に示されるように、メチル化率の形で数値化することができる。図２では、左側に示されるように４個の細胞の４個のDNAにおいて、２つのCpGサイトCpG1,CpG2があって図示されるようにメチル化状態の有無が存在する場合に、各CpGサイトのメチル化率を平均値等による代表値として、右側に示される通り算出されること（CpG1は4個中4個で100%、CpG2は4個中1個で25%）が例として示されている。（ただし、図２は模式例であり、実際は血液などのサンプルから得られる膨大な数のDNAのそれぞれにおける膨大な数のCpGサイトより、メチル化率を算出することができる。）
【0006】
従来技術において、年齢等の影響がメチル化率に現れているものとしてメチル化率から年齢等を推定する際に使われるゲノム領域（CpGサイト）は、過去の文献に基づいた経験上よく知られている特定のもの、あるいはゲノム規模で解析されたものからElastic net（弾性ネット）などの正則化回帰で選択されてきた。
【0007】
こうした従来技術として例えば、特許文献１～４（それぞれ「DNAメチル化に基づく年齢推定マーカと計算モデル」３件及び「DNAメチル化を用いた体内検査値推定マーカと生物学的年齢計算モデル」）は、これらの選択された領域とそれらの領域のDNAメチル化レベルの情報を元に年齢を推定するモデルが提案されている。また、特許文献５（「健康診断に基づく年齢推定装置、年齢学習方法及び年齢学習プログラム」）では、DNAメチル化と血液検査項目を用いて死亡リスクを計算することを提案しており、DNAメチル化からの指標の計算にはElastic netにより得られたLevine2018のモデル（非特許文献１）を使用している。
【先行技術文献】
【特許文献】
【0008】
【特許文献1】米国特許出願公開第２０１６／０２２２４４８号明細書
【特許文献2】米国特許出願公開第２０２０／０１９０５６８号明細書
【特許文献3】国際公開第２０１８／１３９８２６号
【特許文献4】国際公開第２０２０／０７６９８３号
【特許文献5】特許第６９０１１６９号公報
【特許文献6】特願２０２２－１９４００４号
【非特許文献】
【0009】
【非特許文献1】Levine, Morgan E. et al (2018). "An epigenetic biomarker of aging for lifespan and healthspan". Aging 10(4)
【発明の概要】
【発明が解決しようとする課題】
【0010】
従来技術には次のような問題があった。
（１） 分析する際に、DNAメチル化測定結果と分析対象の環境要因に関するラベル（例えば年齢、喫煙履歴、運動履歴など）が必要であり、形質ラベルの無いデータは分析に用いることができない。
（２） 分析結果として、分析対象の形質とDNAメチル化箇所の関連性が回帰係数等の形で得られるが、ここからは「関連性の強さ」のみしか分からない。
（１） 分析する際に、DNAメチル化測定結果と分析対象の環境要因に関するラベル（例えば年齢、喫煙履歴、運動履歴など）が必要であり、形質ラベルの無いデータは分析に用いることができない。
（２） 分析結果として、分析対象の形質とDNAメチル化箇所の関連性が回帰係数等の形で得られるが、ここからは「関連性の強さ」のみしか分からない。
【0011】
このような問題に対して、より効果的なエピゲノム解析を実現すべく、本出願人は特許文献６として項目応答理論を用いたメチル化DNA解析手法を開発した。具体的に、当該解析手法ではゲノム規模で解析されたDNAメチル化レベルから目的とする形質（例えば、実年齢や生物学的年齢、累積喫煙量など）を推定する場合に項目応答理論（IRT）を用いる。すなわち、形質とDNAメチル化レベルとの関連付ける際にIRTを用いることで、モデル構築に必要な検体（サンプル）数を必要最小限に抑えることができる。また、個々のゲノム領域についての項目困難度パラメータ、項目識別力パラメータを得ることが可能であり、これらのパラメータを利用して当該ゲノム領域の特性を知ることができる。
【0012】
このように、特許文献６の解析手法はIRTの適用により効果的なゲノム解析が可能となるものであるが、当該解析手法をさらに効果的に運用するためには、次のような問題に対処する必要がある。
【0013】
すなわち、DNAメチル化から生物学年齢のような何らかの健康指標（目的形質）を推定値として算出する際には、ゲノム中の複数の領域を用いる。候補となる領域は多くの場合、数百から数千万箇所と複数あるが、実際の検査時にすべての領域に関してDNAメチル化の値を測定するにはコストがかかる。検査コストを抑制して可能な限り少ない領域のDNAメチル化レベルの測定結果から精度高く目的の指標の推定をできるようにすることで、既存のゲノムワイドのデータから効率的にゲノム領域を選択することで、低コスト且つ高精度な検査パネル設計などを実現する手法が求められる。
【0014】
上記従来技術の課題に鑑み、本発明は、効果的なエピゲノム解析に寄与することができる情報処理装置、方法及びプログラムを提供することを目的とする。
【課題を解決するための手段】
【0015】
上記目的を達成するため、本発明は情報処理装置であって、ＤＮＡの複数の修飾サイトのそれぞれについて、修飾レベルデータから目的形質を推定する項目応答理論のパラメータを入力として読み込み、目的形質の所与の範囲について、項目応答理論において所与の情報量が確保されるように、前記複数の修飾サイトの中から選別を行う第１処理を実行することを特徴とする。また、当該装置に対応する方法及びプログラムであることを特徴とする。
【発明の効果】
【0016】
本発明によれば、所与の情報量が確保されるように複数の修飾サイトの中から選別することにより、効果的なエピゲノム解析に寄与することができる。
【図面の簡単な説明】
【0017】
【図1】ゲノムとエピゲノムとの模式図である。
【図2】メチル化率の模式図である。
【図3】一実施形態に係る情報処理装置の機能ブロック図である。
【図4】一実施形態に係るIRTパラメータ算出部の機能ブロック図である。
【図5】一実施形態に係るサブセット決定部の機能ブロック図である。
【図6】一実施形態に係る情報処理装置の全体概要的な動作のフローチャートである。
【図8】IRTにおける項目パラメータに応じたロジスティック関数の模式例である。
【図9】識別力パラメータが負の場合のロジスティック関数グラフの例である。
【図10】本実施形態の有効性を確認した実験データ例を２つのグラフとして示す図である。
【図12】本実施形態のIRTで利用する各関数の模式例をグラフとして示す図である。
【図13】範囲設定部がIIAを設定する模式例を示す図である。
【図15】変形例による選別の模式例を示す図である。
【図16】変形例として、探索範囲を狭めた際にマーカがどれだけ減らせるかを提示することの模式例を示す図である。
【図17】探索範囲の縦軸の値を一定値以外にも設定してよいことの模式例を示す図である。
【図18】一般的なコンピュータにおけるハードウェア構成の例を示す図である。
【発明を実施するための形態】
【0018】
図３は、一実施形態に係る情報処理装置100の機能ブロック図である。情報処理装置100は、IRTパラメータ算出部40、サブセット決定部50及び目的形質推定部33を備え、その全体的な動作として図示する通り、（１）入力としてデータベース構築のための情報をIRTパラメータ算出部40において読み込み、多数に渡るIRTパラメータとしての項目パラメータai,biを中間出力として算出し、（２）次いで、この多数のパラメータai,biと管理者等によって入力される事前情報などをサブセット決定部50で読み込み、当該多数（セット全体）の中からサブセットとして選別された項目パラメータai,biを出力し、（３）さらに、この選別された項目パラメータai,biを用いて目的形質推定部33が対象者のDNA修飾レベルデータ（メチル化率）から、当該対象者の潜在特性としての目的形質を推定する。
【0019】
【0020】
図４に示される通り、IRTパラメータ算出部40は、エピゲノムデータセット入力部11、エピゲノムデータ解析部12、ラベル入力部13、対応付け保存部14及び第１データベース15を含む第１構築部10と、エピゲノムデータセット入力解析部21及び第２データベース22を含む第２構築部20と、前処理部31及び項目パラメータ推定部32を含む解析部30と、を備える。図５に示される通り、サブセット決定部50は、範囲設定部51及びマーカ選別部52を備える。
【0021】
図６は、一実施形態に係る情報処理装置100の全体概要的な動作のフローチャートである。図示される通り、ステップS1～S5のうちステップS1～S3がIRTパラメータ算出部40による処理内容であり、ステップS4がサブセット決定部50による処理内容であり、ステップS5が目的形質推定部33による処理である。各ステップの処理内容は以下の通りである。
【0022】
ステップS1では、第１構築部10において第１データベース15を構築する。ステップS2では、第２構築部20において第２データベース22を構築する。ステップS3では、これら構築された第１データベース15及び第２データベース22を参照して解析部30が解析処理を実施することで、IRTパラメータとしての項目パラメータai,biを推定する。ステップS4では、サブセット決定部50が事前情報なども参照して解析処理を実施することにより、ステップS3で推定された多数の項目パラメータai,biの中から選別されるサブセットとして、目的形質を効率的に推定するのに用いることができる少数の項目パラメータai,biを決定する。ステップS5では目的形質推定部33が、対象者について、そのDNA修飾レベルデータの入力を受け付け、この対象者の目的特性を推定して出力して、図６のフローを終了する。
【0023】
当該構築される第１データベース15には、複数のサンプル（被験者）について、そのサンプルIDに、DNA修飾レベルデータ及びラベルを紐づけて記録した情報が記録される。当該構築される第２データベース22には、複数のサンプルについて、そのサンプルIDに、DNA修飾レベルデータを紐づけて記録した情報が記録される。第１データベース15のサンプルと第２データベース22のサンプルとは別途のものであってよい。
【0024】
なお、ステップS1～S3では第１データベース15及び第２データベース22を構築して、これを参照して項目パラメータai,biを推定するが、これらの処理は本出願人による特許文献６と同様の処理を用いることができる。
【0025】
以下、図３の各ステップにおいて実行されるIRTパラメータ算出部40及びサブセット決定部50の各部と目的形質推定部の処理の詳細を説明する。
【0026】
＜ステップS1＞ ステップS1では、複数のサンプル（健康関連の情報（ゲノム関連情報及びラベル）の取得対象となる複数の被験者）についての情報を取得して第１データベース15を構築する。このステップS1は以下の手順S11～S14として構成することができる。なお、手順S11,S12と手順S13とはその先後を問わず実施でき、並行して実施してもよい。これら手順S11,S12と手順S13とを終えた後に手順S14を実施することができる。
【0027】
（手順S11） エピゲノムデータセット入力部11が各サンプルのゲノム関連情報を取得し、これをエピゲノムデータ解析部12へと出力する。
【0028】
ゲノム関連情報は、一般的な次世代シーケンサーの出力フォーマットであるFASTQ形式または参照配列にマッピング済みのBAM形式及びそれに準ずる形式で、ユーザID、年齢（実年齢）、性別、人種を含むユーザ情報を含む。エピゲノムデータセット入力部11はこれらのデータを受け取り、フォーマットのチェックを行ったうえでエピゲノムデータ解析部12へ出力する。
【0029】
ここで、データ入力がFASTQ形式の場合は、冗長な配列情報やシークエンス結果の精度が低いデータを除外するため、データに含まれるクオリティー値に基づいたチェックを行うとともに、シーケンスリードのトリミングを行った上でヒトゲノム参照配列への位置づけによりBAMファイルとなる。なお、メチル化測定手法はデジタルなものからアナログなものまで種々のものが存在しており、最終的に修飾レベル（DNAメチル化率）の情報が得られれば任意の手法を用いるようにしてよく、このようなファイル形式を経る手法を用いることは一例に過ぎない。
【0030】
（手順S12） 手順S11で得た各サンプルのゲノム関連情報をエピゲノムデータ解析部12が解析することで、個々のエピゲノムの箇所における修飾レベルを算出して対応付け保存部14へと出力する。
【0031】
ここで、領域ごとのエピゲノムによる修飾レベルは、統計的検定による有意確率とともに算出することができる。なお、これら手順S11,12の一連の処理のためのツール類は一般に公開されているため、手動による調整を行なってもよい。これにより、各々の対象者のゲノム配列における、DNAメチル化のようなエピゲノムによるDNA修飾レベルのデータとなる。
【0032】
こうして手順S12において、前述の図２で模式的に示した通りの領域（CpGサイト）ごとのDNA修飾レベルの代表値（平均値など）としてのDNAメチル化率が得られる。なお、ヒトでは2,800万程のCpGサイトがある。図１，２等の模試例では各々のCpGサイトを、識別子としての添え字(i=1,2,…,n,…)を付してCpGi(i=1,2,…,n,…)として示しているが、以降でも同様にCpGサイトを添え字番号で区別して説明する。なお、後述する通り、情報処理装置100において推定し、そのサブセットを選別するIRTの項目パラメータai,bi(i=1,2,…,n,…)とは、共通の添え字iを用いていることで表されている通り、当該i番目のCpGサイトであるCpGi(i=1,2,…,n,…)に対応するパラメータである。
【0033】
図１の例EX2にも示されている通り、DNAメチル化等のエピゲノムは遺伝子スイッチの役割を果たすものであることが知られている。高メチル化状態はDNAが折りたたまれており遺伝子を読み取ることができない状態（不活性）であるのに対し、逆に低メチル化状態はDNAにアクセスして遺伝子を読み取ることができる状態（活性）である。なお、ゲノム領域の種類によっては、このように高メチル化が不活性で低メチル化が活性である、ということが該当しないものも相当数が存在する。このため、より一般には、DNAメチル化は遺伝子の活性と何らかの関連があるものと考えればよい。
【0034】
図１の例EX3にも示されている通り、エピゲノム（DNAメチル化）は、各個人の身体環境をモニタリングした情報を含んでいる。そして、図２の模式例に示されるように個別の細胞（及びこの細胞が構成する器官）ごとにDNAメチル化の状態は異なりうるが、多数のDNAにおける各サイトのメチル化率の代表値の情報を手順S12において取得することが可能である。このメチル化率は各個人の身体環境が何らかの形で反映された数値となっている。
【0035】
（手順S13） ラベル入力部13が各サンプル（エピゲノムデータセット入力部11においてゲノム関連情報を受け付けたのと共通の各サンプル）についてのラベルを取得し、これを対応付け保存部14へと出力する。ラベルとしては各サンプルの環境要因に関するラベルとして、例えば年齢、喫煙履歴、運動履歴あるいはその他の心身の健康に影響しうる履歴や実績等などの入力を受け付けることができる。以下では年齢（実年齢）をラベルとして受け付けた場合を例として説明するが、年齢以外のラベルに関しても同様に、情報処理装置100において分析処理を行うことが可能である。なお、喫煙履歴は後述するPack-yearとして喫煙累積量の形で受け付けるようにしてもよく、同様に、運動履歴に関しても、運動累積量（歩いた距離の総和など）や平均運動時間（例えば１日あるいは１週間ごとに、走る等の積極的な運動を行った時間など）の形で受け付けるようにしてもよい。心身の健康に影響しうる履歴や実績として、BMI（既知のように「体重kg」÷「身長mの２乗」で定義され、肥満度を表す体格指数とされるボディマス指数）、累積睡眠時間、平均睡眠時間、大気汚染物質への累積暴露期間、心理尺度などをラベルに用いてもよい。心理尺度としては例えば、PTSD診断尺度などの、複数の質問項目に対する回答の結果に所定の関数を適用して、心理状態が良いか悪いかを数値で表現したものを用いるようにしてよい。
【0036】
（手順S14） 以上の手順S12,13で得たDNA修飾レベルデータとラベルとを、対応付け保存部14が対応するサンプルIDと紐づけたうえで、第１データベース15へと保存することによって第１データベース15を構築する。
【0037】
＜ステップS2＞ ステップS2では、複数のサンプル（健康関連の情報としてゲノム関連情報の取得対象となる複数の被験者）についての情報をエピゲノムデータセット入力解析部21において取得して各サンプルのDNA修飾レベルデータを得ることで、第２データベース22を構築する。
【0038】
エピゲノムデータセット入力解析部21の処理は、ステップS1において手順S11,S12として説明したエピゲノムデータセット入力部11及びエピゲノムデータセット解析部12による処理と共通であり、入力データとして用いるサンプルのみが異なるため、重複する説明を省略する。こうして、第２データベース22には、エピゲノムデータセット入力解析部21で得られた、各サンプルについてのサンプルIDと、DNA修飾レベルデータとが紐づけて保存される。
【0039】
＜ステップS3＞ ステップS3では解析部30が、以上の通り構築された第１データベース15及び第２データベース22を参照し分析した結果のパラメータを用いて、IRTパラメータとしての項目パラメータai,biの推定結果を得る。
【0040】
図７は、一実施形態に係る図６のステップS3の詳細を示すフローチャートであり、一実施形態に係る分析部30の処理を示すものである。図７のフローに従って解析部30では具体的に、前処理部31及び項目パラメータ推定部32が以下の処理を行う。
【0041】
＜前処理部31…ステップS101～S106＞
前処理部31は、第１データベース15を参照することで、各CpGサイトi（i=1,2,…,I）（ここで、CpGサイトの総数をIとする）のうち、ラベルと相関があると判定されるものを項目パラメータ推定部32へと出力する。
前処理部31は、第１データベース15を参照することで、各CpGサイトi（i=1,2,…,I）（ここで、CpGサイトの総数をIとする）のうち、ラベルと相関があると判定されるものを項目パラメータ推定部32へと出力する。
【0042】
図７の各ステップのうちステップS101～S106が前処理部31による処理であり、各CpGサイトi（i=1,2,…,I）のうちステップS104を経由したものが、ラベルと相関があると判定されたものとして項目パラメータ推定部32へと出力される。
【0043】
ステップS101では各CpGサイトを表すインデクスiを初期値i=1に設定してからステップS102へと進む。ステップS102では、第１データベース15の全サンプルs（s=1,2,…,M）（Mは第１データベース15に記録されているサンプル総数とする）でのラベルL(s)（例えば年齢）と、当該CpGサイトiについてのサンプルsのメチル化率β(i,s)との相関係数r(i)を求めてからステップS103へと進む。
【0044】
ステップS103では、ステップS102で求めた相関係数r(i)の絶対値|r(i)|が予め設定しておく閾値を超えるか否かを判定することにより、当該CpGサイトi（そのメチル化率βi）とラベルLとに正または負の相関があるか否かを判定し、相関があると判定されればステップS104へと進み、相関がないと判定された場合はステップS105へ進む。
【0045】
ステップS104では、当該CpGサイトiを次の項目パラメータ推定部32で利用するものとして設定してからステップS105へ進む。（従って、ステップS103で否定判定となったCpGサイトiに関しては、当該ステップS104がスキップされることから、次の項目パラメータ推定部32で利用するものとしては設定されないものとなる。）ステップS105では、i=Iであるか、すなわち、全てのCpGサイトiについて以上の相関の有無の判定が完了しているかを判定し、完了していればステップS201へと進み、未完了であればステップS106へ進み、インデクスiを次の値i+1に更新してからステップS102へ戻ることで、未完了のCpGサイトiについての相関有無の判定処理を継続する。
【0046】
＜項目パラメータ推定部32…ステップS201＞
項目パラメータ推定部32はステップS201において、上記の前処理部31でラベルに対して相関があると判定されたCpGサイトiのそれぞれについて、第２データベース22を参照することによって項目パラメータai,biを推定してサブセット決定部50へと出力する。本実施形態の項目パラメータ推定部32では特に、教育測定の分野等でテスト分析等の手法として用いられている項目応答理論（IRT）を利用することで、項目パラメータai,biを推定する。
項目パラメータ推定部32はステップS201において、上記の前処理部31でラベルに対して相関があると判定されたCpGサイトiのそれぞれについて、第２データベース22を参照することによって項目パラメータai,biを推定してサブセット決定部50へと出力する。本実施形態の項目パラメータ推定部32では特に、教育測定の分野等でテスト分析等の手法として用いられている項目応答理論（IRT）を利用することで、項目パラメータai,biを推定する。
【0047】
教育分野でのIRTにおいては、テスト（試験問題）が複数（i∈I）の問題（項目）から構成されているとし、各項目の特性をロジスティック関数でモデル化する。IRTのモデルとして様々なモデルが提案されているが、ここでは、２パラメタロジスティックモデル（2PL）について説明する。2PLでは次式(1)により、項目をモデル化している。
【0048】
【数1】
【0049】
ここで、xi,sは被験者(s∈S)が、項目iを解答した際の正誤情報を表し、正解の場合は１、不正解の場合は０の値をとる。また、θsは被験者sの潜在特性(能力)を表し、ai,biは項目iにおける項目パラメータで、識別力パラメータと、困難度パラメータをそれぞれ表している。上式では、ロジスティック関数により、能力θsの被験者が、ある項目iを解答した場合に正解する確率を表現している。通常、項目パラメータai,biは、複数の被験者による、実際のテストの解答結果から推定され、また、θsについては、推定された項目パラメータを用いて推定される。
【0050】
なお、式(1)は正解確率を表しているが、これと同一の式において各被験者sの潜在能力θsの取りうる値を変数θとして、潜在能力θにおける正解確率を意味するものとして、項目パラメータai,biにより項目iについての情報項目関数Pi(θ)を以下の式(1')のように定義できる。
【0051】
【数2】
【0052】
図８はIRTにおける項目パラメータに応じたロジスティック関数（情報項目関数）の模式例であり、グラフG1,G2,G3は、識別力パラメータaiを固定したうえで困難度パラメータbiをそれぞれ大、中、小としたものであり、グラフG4,G5,G6は、困難度パラメータbiを固定したうえで識別力パラメータaiをそれぞれ大、中、小としたものである。これらから見て取ることができるように、困難度パラメータbiは試験問題の困難度（難易度）に対応して、ロジスティック関数のグラフの水平位置を定めるものであり、識別力パラメータaiは試験問題の識別力に対応して、ロジスティック関数のグラフの傾斜具合を定めるものである。（傾斜具合が大きい極限としてステップ関数状となり、中間ではS字形状となり、傾斜具合が小さい極限として値0.5の定数関数となる。）
【0053】
なお、図８の例は全て識別力パラメータaiが正であり、能力θsの高い被験者sほど試験問題の正解率が高い場合を示しているが、識別力パラメータaiが負の場合のロジスティック関数グラフG7の例を図９に示す通り、ai <0の場合は能力θsの高い被験者sほど試験問題の正解率が低くなり、試験問題自体が不適切であった等の状況が想定される。
【0054】
本実施形態の項目パラメータ推定部32では、（既に共通のインデクスiで説明している通り、）測定されたCpGサイトiを、IRTにおける項目iとみなし、各CpGサイトiのメチル化の有無を、IRTにおける解答の正誤情報として処理を行う。ここで、第２データベース22に記録されているCpGサイトのメチル化測定においては、前述したように複数の細胞のメチル化の有無の平均（メチル化率）として得られることが多く、一般的である。上記教育分野等での２PLを適用する場合の入力は、０（不正解）か１（正解）である必要があるため、このような場合は、項目パラメータ推定部32においてIRTの項目パラメータ推定ai,biの処理にかける際に、予め設定しておく閾値以上のメチル化率を１（メチル化があると判定される結果として、試験問題の場合の正解に対応）、閾値未満の場合を０（メチル化がないと判定される結果として、試験問題の場合の不正解に対応）とするなどして、二値化の処理を行えばよい。（前述した例とは異なり、第２データベース22の記録がメチル化率有無平均のメチル化率としてではなく、メチル化の有無で２値として用意されている場合は、この閾値による２値化の処理を省略してよい。）
【0055】
ただし、メチル化に関しては、試験問題の場合に正解の方が良い評価となる（能力θsが高い評価となる）こととは異なり、必ずしも、メチル化があることが良く、メチル化がないことが悪いことを意味するとは限らない。
【0056】
すなわち、通常の学力テストとしての試験問題においては、能力が高いほど、テストの点数が良くなり、また、能力が高いほど、難しい問題に対する正答率も高くなるという前提条件がおけるが、DNAメチル化分析においては、それぞれのCpGサイトが、分析対象のラベルとそもそも関連性が全く無い場合も多く、事前のCpGサイトのフィルタリングが重要となる。
【0057】
本実施形態では、前処理部31がこのフィルタリングを実現している。図１０は、本実施形態の有効性を確認した実験データ例を２つのグラフG11,G12として示す図であり、グラフG11は、前処理部31で相関ありと判定されフィルタリングされた各CpGサイトiについて、項目パラメータ推定部32で推定した項目パラメータai,biでのロジスティック関数を描いたものであり、ラベルLとして年齢を推定する（（後述する目的形質推定部33において推定する））ために、抽出された100のCpGサイトに対する項目特性曲線を示している。なお、項目パラメータaiの値は正負のいずれも存在しうるが、グラフG11では便宜上、項目パラメータaiが全て正であるものとして（すなわち、項目パラメータaiが負となる場合は左右反転させて）描いている。（２つめのグラフG12については後述する。）
【0058】
本実施形態では以上のようにメチル化率を２値化したものを教育分野のIRTの場合の正解／不正解に見立て、ラベルとして与える潜在特性について、IRTの項目パラメータai,biを項目パラメータ推定部32において推定する。当該推定するための数値計算処理については、既存のIRTと同様の任意の計算手法を用いてよく、前処理部31で選択された（相関が高い）CpGサイトをまとめた上で、周辺最尤推定により、項目パラメータを計算すればよい。例えば、相関が高いとされた項目数が100の場合、100回の周辺最尤推定を行うのではなく、100項目のデータを一括処理して、1回の周辺最尤推定で、100個のaiとbiを推定すればよい。（なお、このパラメータ推定の際は、第２データベース22を参照するので、教育分野のIRTの場合と同様に、潜在特性としてのラベルの値は第２データベース22に記録されておらず未知の状態であり、いわゆる教師無し学習により、パラメータ推定結果と共に、当該未知であった各サンプルの潜在特性を推定した結果を得ることができるが、本実施形態ではこの結果の情報は利用しなくてもよい。）
【0059】
なお、各々のCpGサイトiに関しては、そのメチル化率の検査結果から目的形質を推定するマーカの役割を果たすものであるため、文脈に応じて適宜、「マーカ」や「検査項目」等とも呼ぶこととする。
【0060】
＜ステップS4…サブセット決定部50＞
ステップS4ではサブセット決定部50が、以上のIRTパラメータ決定部50によってラベルに相関があると判定された複数のCpGサイトiのそれぞれについて算出された項目パラメータai,biの中から、対象者の目的形質を効果的に推定できる項目パラメータai,biをサブセットとして決定する。
ステップS4ではサブセット決定部50が、以上のIRTパラメータ決定部50によってラベルに相関があると判定された複数のCpGサイトiのそれぞれについて算出された項目パラメータai,biの中から、対象者の目的形質を効果的に推定できる項目パラメータai,biをサブセットとして決定する。
【0061】
【0062】
ここで、図１１の各ステップの詳細を説明するための前提事項（サブセット決定部50の処理における考え方）を説明する。本実施形態のサブセット決定部50では、IRTを用いることで、目的形質とDNAメチル化率に関して個々のゲノム領域の情報項目関数Pi(θ)（2PLについて前述の式(1')）から、メチル化率から推定される潜在特性（目的形質）の確度を意味するテスト情報量（テスト情報関数）Ii(θ)が得られること（テスト情報量に関する以下の式(2)）とそれらが加算的に取り扱えること（テスト情報量の加算性に関する以下の式(3)）を利用して、特定検体のおおよその潜在特性（形質の度合いに関するおおよその度合い）を推定するために最低限のマーカのサブセットを選択するのがサブセット決定部50の役割である。この際、テスト情報量からその平方根の逆数として算出できるによる誤差（以下の式(4)）も適宜、考慮する。
【0063】
【数3】
【0064】
図１２は、以上の式(1)～(4)による、本実施形態のIRTで利用する各関数（いずれも目的形質θの関数）の模式例をグラフとして示す図である。例EX1に示すように、項目iの各々について、式(1')の情報項目関数Pi(θ)から式(2)のテスト情報量Ii(θ)を算出できる。（また、当該加算しない単独のテスト情報量Ii(θ)について、式(3),(4)によりその誤差も算出できる。）
【0065】
また、例EX2は式(3)のテスト情報量の加算性の例を示す図であり、上段部分に示すある項目i=1のテスト情報関数（θが小さい箇所に１つのピークを有する）と中段部分に示すある項目i=2のテスト情報関数（θが大きい箇所に１つのピークを有する）とを加算することで、下段部分に示すテスト情報関数（θが小さい箇所と大きい箇所とで２つのピークを有する）が得られることを示している。また、例EX2に示される通り、加算前や加算後のいずれのテスト情報関数についても、その誤差関数を算出することができる。例EX2では２つのテスト情報関数を加算しているが、式(3)に一般的に示される通り、３つ以上も全く同様に加算することができる。（なお、式(3)では加算性の例として項目i=1,2,…,ioptのように1から連番で順にiopt番目まで加算しているが、実際に加算する項目は例えばi=3,5,10,19,…のように、連番に限らない任意のものでよい。）
【0066】
＜ステップS41…範囲設定部51＞
ステップS41では、範囲設定部51が、管理者等によって入力される事前情報をもとに、対象者の目的形質を効果的に推定するための注目対象となる範囲を、テスト情報量のグラフ範囲（図１２に例示したような、目的形質の横軸θと情報量の縦軸Iとで定義されるグラフ範囲）内に設定する。
ステップS41では、範囲設定部51が、管理者等によって入力される事前情報をもとに、対象者の目的形質を効果的に推定するための注目対象となる範囲を、テスト情報量のグラフ範囲（図１２に例示したような、目的形質の横軸θと情報量の縦軸Iとで定義されるグラフ範囲）内に設定する。
【0067】
すなわち、本実施形態によるマーカのサブセットの設計を行う管理者等による入力を受け付けることで、事前情報として目的形質（横軸θ）の探索範囲「θ下限≦θ≦θ上限」を決める。例えば、目的形質が生物学的年齢であり、老若を問わずに任意のユーザをマーカによる推定対象として設定したい用途の場合であれば、「0歳≦θ≦100歳」のように広範囲に目的形質の範囲を設定すればよく、特定年代層として例えば20代のユーザのみに限定してマーカによる推定対象として設定したい用途の場合であれば、「20歳≦θ＜30歳」のように特定範囲に限定して目的形質の範囲を設定すればよい。
【0068】
次に、同様にして管理者等による入力を受け付けることで許容される誤差を決め、そこから式(4)により逆算で検査に必要な最低限のテスト情報量I最低限を算出し、I＞I最低限が成立するためにテスト情報量の関数グラフ上で覆いつくされている必要がある範囲として、縦軸Iの範囲「0≦I≦I最低限」を定める。（すなわち、図１４の例EX43として後述するように、テスト情報量Iの関数グラフがこの縦軸範囲「0≦I≦I最低限」を覆いつくしていることは、テスト情報量Iの関数値が「I＞I最低限」であることと、同値である。）例えば、誤差の値として一定値Eが許容範囲であれば、式(4)により、テスト情報IについてI＞1／E2であることによってこの誤差許容範囲が達成されることから、I最低限の値としてI最低限=1／E2を用いて、縦軸Iの範囲を設定できる。
【0069】
あるいは例えば、IRTパラメータ算出部40において相関ありとして得られた全てのマーカi（CpGサイトi）についてそのパラメータai,biより式(2)でテスト情報量Ii(θ)を求め、横軸範囲「θ下限≦θ≦θ上限」内におけるそれらの最大値Imaxの一定割合r（例えば80%としてr=0.8）として、以下の式(5a),(5b)のようにI最低限を決定してもよい。
【0070】
【数4】
【0071】
以上の通り、範囲設定部51において２次元(θ,I)グラフの空間範囲内に設定する範囲を、マーカによる検査が持つテスト情報量が最低限満たさなければならない注目領域としてIIA（Interest information area）と呼ぶこととする。図１３に、範囲設定部51がIIAを設定する模式例を示す。
【0072】
＜ステップS42からステップS47まで…マーカ選別部52＞
ステップS42からステップS47まではマーカ選別部52によるマーカ選別処理である。以下、これら各ステップを説明するに際しては、動作主体であるマーカ選別部52の言及は原則として省略するものとし、その処理内容のみを説明する。
ステップS42からステップS47まではマーカ選別部52によるマーカ選別処理である。以下、これら各ステップを説明するに際しては、動作主体であるマーカ選別部52の言及は原則として省略するものとし、その処理内容のみを説明する。
【0073】
＜ステップS42＞
ステップS42では、困難度パラメータbiが探索範囲IIA内にある（すなわち、探索範囲IIAの横軸範囲「θ下限≦θ≦θ上限」内にある）マーカのプールを、最終的に選別されるマーカの候補として作成する。項目iのマーカをMiと表記することにより、当該作成されるプールPLは以下のように表記することができる。
PL={Mi|θ下限≦bi≦θ上限}
ステップS42では、困難度パラメータbiが探索範囲IIA内にある（すなわち、探索範囲IIAの横軸範囲「θ下限≦θ≦θ上限」内にある）マーカのプールを、最終的に選別されるマーカの候補として作成する。項目iのマーカをMiと表記することにより、当該作成されるプールPLは以下のように表記することができる。
PL={Mi|θ下限≦bi≦θ上限}
【0074】
図８や図１２の例EX1から理解されるように、当該作成されるプールPLは、そのテスト情報量Ii(θ)のピーク位置bi（2PLの情報項目関数Pi(θ)の値が1/2となる位置）が、探索範囲IIA内にあるものとなる。
【0075】
なお、このステップS42は、マーカ選別部52が最終的な選別結果として得るマーカの候補を予め絞り込む役割を有するものであるため、省略されてもよい。（なお、省略された場合は、サブセット決定部50で得たマーカMiの全部が当該プールを構成しているものとして扱えばよい。）
【0076】
＜ステップS43＞
ステップS43では、最終的に選別されるサブセット（SBとする）の初期値として空集合φを設定（SB=φ）し、且つ、探索範囲IIAを網羅するのに最低限必要となるテスト情報量（IIAミニマムとする）の初期値をこの探索範囲IIAに一致するように設定（IIAミニマム=IIA）してから、ステップS44へと進む。
ステップS43では、最終的に選別されるサブセット（SBとする）の初期値として空集合φを設定（SB=φ）し、且つ、探索範囲IIAを網羅するのに最低限必要となるテスト情報量（IIAミニマムとする）の初期値をこの探索範囲IIAに一致するように設定（IIAミニマム=IIA）してから、ステップS44へと進む。
【0077】
なお、以降のステップS44,S45,S46において当該サブセットSBに属する要素が逐次的に追加されていき、当該最低限必要なテスト情報量（の残り）IIAミニマムがゼロとなる（探索範囲IIAが網羅され尽くす）ことにより、最終的なサブセットSBが決定されることとなる。
【0078】
以下では説明のため、繰り返しのステップS44,S45,S46の繰り返し回数をk(k=1,2,3,…)とする。ステップS46からステップS44へ戻る都度、この繰り返し回数kの値は１ずつ増える。
【0079】
＜ステップS44＞
ステップS44では、当該k回目の時点でのプールPL（=PL(k)とする）の中から、当該k回目の時点での残りの最低限必要なテスト情報量IIAミニマム（=IIAミニマム[k](θ)とする）をカバーするテスト情報量（目的形質θの関数としての各項目jのテスト情報量Ij(θ)と、目的形質θの関数であるテスト情報量IIAミニマム[k](θ)と、のグラフ上での重複面積）が最も大きいマーカ（Mj(k)とする）を決定して、サブセットSBに当該決定されたマーカMj(k)を追加してから、ステップS45へ進む。
ステップS44では、当該k回目の時点でのプールPL（=PL(k)とする）の中から、当該k回目の時点での残りの最低限必要なテスト情報量IIAミニマム（=IIAミニマム[k](θ)とする）をカバーするテスト情報量（目的形質θの関数としての各項目jのテスト情報量Ij(θ)と、目的形質θの関数であるテスト情報量IIAミニマム[k](θ)と、のグラフ上での重複面積）が最も大きいマーカ（Mj(k)とする）を決定して、サブセットSBに当該決定されたマーカMj(k)を追加してから、ステップS45へ進む。
【0080】
このマーカMj(k)はプールPL(k)内から、以下の式(6a),(6b)の通り探索範囲IIA内での２つのテスト情報量関数（IIAミニマム[k](θ)とIj(θ)）重複面積overlap(j)が最大となるものとして決定することができる。なお、式(6)にて{PL(k)}は、プールPL(k)に属するマーカMjのインデクスjの集合とする。
【0081】
【数5】
【0082】
＜ステップS45＞
ステップS45では、当該k回目のステップS44で新たに選別されたマーカMj(k)が追加されたサブセットSB（k回目のものである点を明示するためSB[k]とする）においてテスト情報量の和I(θ)（k回目のものである点を明示するためI(θ)[k]とする）を以下の式(7a)の通り算出し、選別途中における当該k回目のテスト情報量IIAミニマム[k](θ)も同様に以下の式(7b)の通り更新して、さらに、プールPL(k)も後述の通り更新してから、ステップS46へと進む。
ステップS45では、当該k回目のステップS44で新たに選別されたマーカMj(k)が追加されたサブセットSB（k回目のものである点を明示するためSB[k]とする）においてテスト情報量の和I(θ)（k回目のものである点を明示するためI(θ)[k]とする）を以下の式(7a)の通り算出し、選別途中における当該k回目のテスト情報量IIAミニマム[k](θ)も同様に以下の式(7b)の通り更新して、さらに、プールPL(k)も後述の通り更新してから、ステップS46へと進む。
【0083】
【数6】
【0084】
なお、式(7a)にて、{SB[k]}は、サブセットSB[k]に属するマーカMjのインデクスjの集合である。また、式(7b)では、左辺と右辺とのいずれにも関数であるテスト情報量IIAミニマム[k](θ)が現れるが、これはプログラミングの慣行的な表記における変数への値の代入の場合のように、（例えば、「a=a+1」の表記により、変数aの現在の値（右辺）に1を加算したものを新たな変数aの値（左辺）とすることを意味するように、）右辺が更新前のテスト情報量IIAミニマム[k](θ)を意味し、左辺が更新後のテスト情報量IIAミニマム[k](θ)を意味するものとする。（すなわち、この更新後の左辺のテスト情報量IIAミニマム[k](θ)は、次のk+1回目のステップS44において、IIAミニマム[k+1](θ)として利用されるものとなる。）
【0085】
ステップS45では同様に、マーカMj(k)が選別されたことに伴い、当該マーカMj(k)が今後（k+1回目以降）に重複して２回以上選別されることがないように、プールPL(k)から当該マーカMj(k)を除外することによってプールPL(k)を更新する。
【0086】
ステップS46では、当該k回目において前記式(7a)で求まっているテスト情報量の和I(θ)[k]が、探索範囲IIA内の全体（横軸範囲「θ下限≦θ≦θ上限」の全体）において許容範囲（すなわち、θの関数としての探索範囲IIA）を上回っているかを、以下の式(8a)（またはこれと同値な式(8b)）が真となるか否かにより判定し、上回っていればステップS47へ進み、上回っていなければステップS44へと戻り、次のk+1回目について以上のステップS44,S45,S46を同様に繰り返す。
I(θ)[k]≧IIA …(8a)
IIAミニマム[k](θ)=0 …(8b)
I(θ)[k]≧IIA …(8a)
IIAミニマム[k](θ)=0 …(8b)
【0087】
なお、範囲「θ下限≦θ≦θ上限」において式(8b)が成立し、IIAミニマム[k](θ)の値が常にゼロとなることは、当該k回目の時点での残りの最低限必要なテスト情報量がゼロであることから、当初の探索範囲IIAはk回に渡って選別された和I(θ)[k]によって覆いつくされていること、すなわち、式(8a)が成立することを意味するため、これらの式(8a),(8b)は同値である。
【0088】
＜ステップS47＞
ステップS47では、当該k回目のステップS46にて真となる判定が得られることで、当初の探索範囲IIAを覆いつくすことがわかったサブセットSB[k]を、マーカの最終的な選別結果として出力し、図１１のフローを終了する。
ステップS47では、当該k回目のステップS46にて真となる判定が得られることで、当初の探索範囲IIAを覆いつくすことがわかったサブセットSB[k]を、マーカの最終的な選別結果として出力し、図１１のフローを終了する。
【0089】
図１４は、マーカ選別部52による以上の図１１のステップS42～S47の処理内容の模式例をテスト情報量のグラフ（横軸θ、縦軸I(θ)）上で例EX41～EX43として示す図である。例EX40はステップS42,S43の例として、選別開始時の初期値として最小限必要なテスト情報量IIA[ミニマム]を当初のIIAに一致するようにして設定していることと、且つ、選別対象の候補であるプールPLに属するマーカの例としてマーカM1,M2,…,M7,…があるものとし、これらに対応するテスト情報関数I1,I2,…,I7,…があることが示されている。
【0090】
例EX41ではk=1回目のステップS44～S46において、最小限必要なテスト情報量IIAミニマムとマーカMjのテスト情報量Ijとの重複が最大となるマーカとしてj=3（j(1)=3）のマーカM3が選別されてサブセットSB[1]＝{M3}が構成されてそのテスト情報量の和I(θ)[1]=I3(θ)が定まり、重複箇所を減算して更新することにより最小限必要なテスト情報量IIAミニマム[1]が以下のように得られる例が示されている。
IIAミニマム[1]=IIA-I[1](θ)=IIA-I3(θ)
なお、上記の通りk=1回目のステップS45で更新されたIIAミニマム[1]が、次のk=2回目でのステップS44におけるIIAミニマム[2]として、最大重複判定の対象として利用されることとなる。
IIAミニマム[1]=IIA-I[1](θ)=IIA-I3(θ)
なお、上記の通りk=1回目のステップS45で更新されたIIAミニマム[1]が、次のk=2回目でのステップS44におけるIIAミニマム[2]として、最大重複判定の対象として利用されることとなる。
【0091】
同様に、例EX42ではk=2回目のステップS44～S46において、最小限必要なテスト情報量IIAミニマム[2]とマーカMjのテスト情報量Ijとの重複が最大となるマーカとしてj=5（j(2)=5）のマーカM5が選別されてサブセットSB[2]＝{M3,M5}が構成されてそのテスト情報量の和I(θ)[2]=I3(θ)+I5(θ)が定まり、重複箇所を減算して更新することにより最小限必要なテスト情報量IIAミニマム[2]が以下のように得られる例が示されている。
IIAミニマム[2]=IIA-I[2](θ)=IIA-I3(θ)-I5(θ)
なお、上記の通りk=2回目のステップS45で更新されたIIAミニマム[2]が、次のk=3回目でのステップS44におけるIIAミニマム[3]として、最大重複判定の対象として利用されることとなる。
IIAミニマム[2]=IIA-I[2](θ)=IIA-I3(θ)-I5(θ)
なお、上記の通りk=2回目のステップS45で更新されたIIAミニマム[2]が、次のk=3回目でのステップS44におけるIIAミニマム[3]として、最大重複判定の対象として利用されることとなる。
【0092】
こうして、各k回目の処理で最大重複するものをサブセットに追加して、テスト情報量の和や、残っている最小限必要となるテスト情報量IIAミニマム[k]を更新していくことで、例EX43に示されるように、最終的にサブセット内のテスト情報量の和が当初のIIAを覆いつくした時点で、サブセットの選別が終了し、最終的なサブセット選別結果を得ることができる。
【0093】
【0094】
＜ステップS5…目的形質推定部33＞
ステップS5では目的形質推定部33が、対象者（第１データベース15や第２データベース22に記録されているサンプルとは別途の新規の対象者）について、そのDNA修飾レベルデータ（エピゲノムデータセット入力解析部21で出力されるのと同じ形式でのDNA修飾レベルデータ）の入力を受け付け、この対象者の目的特性を推定して出力し、図６のフローを終了する。この際、サブセット決定部50でサブセットとして決定されたマーカ集合を用いて、目的特性を推定する。
ステップS5では目的形質推定部33が、対象者（第１データベース15や第２データベース22に記録されているサンプルとは別途の新規の対象者）について、そのDNA修飾レベルデータ（エピゲノムデータセット入力解析部21で出力されるのと同じ形式でのDNA修飾レベルデータ）の入力を受け付け、この対象者の目的特性を推定して出力し、図６のフローを終了する。この際、サブセット決定部50でサブセットとして決定されたマーカ集合を用いて、目的特性を推定する。
【0095】
すなわち、サブセットに属するマーカMjについて項目パラメータ推定部32で推定された項目パラメータaj,bjを用いて最尤推定などを用いることで、このマーカ集合{Mi}から推定されるものとして、目的特性を推定することができる。例えば、前処理部31で相関があるとされてIRTパラメータ算出部40においてその項目パラメータai,biが算出されたCpGサイトi（マーカMi）が全部で1000個であり、この1000個のマーカ集合全体{Mi}の中からサブセットとして10個のみのマーカ集合{Mj}が選別された場合であれば、既存のIRTの場合と同様に最尤推定などの手法で、当該選別された10個のマーカ集合{Mj}（CpGサイトjの集合）に関して共通の１つの値として、対象者の潜在特性θを推定することができる。
【0096】
すなわち、この対象者sのCpGサイトj（マーカMj）のメチル化率βjを、（項目パラメータ推定部32でのパラメータ推定処理時のように閾値で２値化するのではなく、）0≦βj≦1の範囲内での値として、式(1)（式(1)の逆関数）に入力することによって、対応する目的形質θsを推定することができる。
【0097】
以上、本発明の実施形態によれば、サブセット決定部50において情報量の多いマーカを優先的に選別する形で、段階的に選別されるマーカを増やしていくことにより、必要最低限のマーカセットを自動選別することが可能となる。当該選別された必要最低限のマーカにより対象者の検査（目的形質推定部33の処理に相当）を行うことで、検査コストを抑制でき、（例えば、検査コストのうちマーカ数に比例するコストがあったとすると、マーカが1000個から10個に選別された場合は、当該比例コストを1/100に抑制でき、）且つ、マーカ全体でIIAを覆いつくす情報量の多いマーカが利用されていることから、検査精度も確保することが可能となる。
【0098】
すなわち、マーカの選別はIIAを覆いつくすように行われることから、IIAの横軸範囲において、その縦軸の高さで表現される情報量が確保されており、当該情報量に応じた検査精度を確保することができる。また、マーカの選別はIIAをより効果的に覆うマーカを優先して選別するようにしているので、少数のマーカによってコスト効率の高い検査が可能となる。
【0099】
以下、本発明の実施形態に関する種々の補足例、追加例、代替例などについて説明する。
【0100】
（１） 本発明の実施形態において、目的形質推定部33で推定される対象者の目的形質θは、対象者の健康状態等が反映され数値化されたものである。従って、健康状態の維持・向上のための行動変容、あるいは少なくとも意識付けを促す数値指標として目的形質を推定することができるため、国連が主導する持続可能な開発目標（ＳＤＧｓ）の目標３「あらゆる年齢のすべての人々の健康的な生活を確保し、福祉を推進する」に貢献することが可能となる。
【0101】
（２） 本発明の実施形態に係る情報処理装置100のうち、サブセット決定部50によりサブセットを絞り込む処理以外のIRTパラメータ算出部40と目的形質推定部33との処理は、本出願人による特許文献６の手法と同様の処理を用いることができる。従って、情報処理装置100をサブセット決定部50のみで構成されるものとして用意し、IRTパラメータ算出部40と目的形質推定部33の処理は情報処理装置100とは別途の外部サーバ等で行われる処理として、外部に委ねるものとしてもよい。例えば、IRTパラメータ算出部40で算出されるべき目的形質と相関のある項目パラメータai,biは、予め用意されているものとして、サブセット決定部50で読み込むようにしてよい。
【0102】
サブセット決定部50を少なくとも含んで構成される情報処理装置100によって実現できる用途として例えば、以下を挙げることができる。
●目的形質（潜在特性）の検査パネルの設計サービス
●ゲノム規模で広く取られたDNAメチル化率の測定データと目的形質とのデータを入力として、その中から、目的形質の推定に必要な最小ゲノム領域（CpGサイト）のセットを出力するサービス
●事前情報をもとに検査に必要なゲノム領域を選択するサービス
●生物学的年齢など、DNAメチル化進行度に相関する形質の推定サービス
●目的形質（潜在特性）の検査パネルの設計サービス
●ゲノム規模で広く取られたDNAメチル化率の測定データと目的形質とのデータを入力として、その中から、目的形質の推定に必要な最小ゲノム領域（CpGサイト）のセットを出力するサービス
●事前情報をもとに検査に必要なゲノム領域を選択するサービス
●生物学的年齢など、DNAメチル化進行度に相関する形質の推定サービス
【0103】
（３） 変形例に係る実施形態として、マーカを選ぶ過程でIIAとの加算分の重複面積が同値のマーカが複数存在する場合、次の２つの方針のどちらかに基づき同値マーカ内から１つの選択のマーカを選択するようにしてよい。
【0104】
＜方針１＞ユーザが「探索範囲内での推定精度を重視したい場合」
上記同値マーカの中から、探索範囲内においてIIAとの重複にこだわらず（縦軸方向に）より大きなI(θ)を持つマーカを１つ選択する。
＜方針２＞ユーザが「探索範囲の周辺の情報を得たい場合」
上記同値マーカの中から、探索領域の両端より外側に関して（横軸方向に）より大きなI(θ)を持つマーカを１つ選択する。
上記同値マーカの中から、探索範囲内においてIIAとの重複にこだわらず（縦軸方向に）より大きなI(θ)を持つマーカを１つ選択する。
＜方針２＞ユーザが「探索範囲の周辺の情報を得たい場合」
上記同値マーカの中から、探索領域の両端より外側に関して（横軸方向に）より大きなI(θ)を持つマーカを１つ選択する。
【0105】
すなわち、この変形例は図１１のステップS44の変形例として、式(6a)のように重複面積overlap(j)が最大となるものを選別することに代えて、閾値判定で上位側に属すると判定されるもの（例えば、最大重複面積の80%以上の重複面積を有するマーカを上位側と判定）が複数存在する場合に、当該複数の中から、方針１または方針２の基準により、１つのマーカを選別するものである。ステップS44においてプールPL(k)の中から重複面積が上位側に属すると判定されたものをPL(k)[上位]とし、そのインデクス集合を{PL(k)[上位]}とすると、方針１は以下の式(8-1)及び式(8-1-1)により、方針２は以下の式(8-2)、式(8-2-1)及び式(8-2-2)により、マーカMj(k)を選別することができる。
【0106】
【数7】
【0107】
なお、式(8-2-1),(8-2-2)において、θ下限拡張＜θ下限＜θ上限＜θ上限拡張であり、探索範囲IIAの横軸範囲「θ下限≦θ≦θ上限」に対して所定の低い側の範囲「θ下限拡張≦θ≦θ下限」と所定の高い側の範囲「θ上限≦θ≦θ上限拡張」とを、探索範囲の両側より外側（θが低い側及び高い側）での面積area低い側(j)及びarea高い側(j)をそれぞれ評価するために、事前情報などにより予め設定しておけばよい。探索範囲の外側の面積に関して、低い側と高い側の両方ではなく、片方のみを用いるようにしてもよい。
【0108】
なお、式(8-2)等では方針２における探索範囲の外側でより大きな値を持つことを面積で評価しているが、面積に代えて、以下の式(8-3),(8-3-1)のように、当該外側範囲での情報量の最大値で評価するようにしてもよい。
【0109】
【数8】
【0110】
図１５は、当該変形例による選別の模式例を示す図であり、例えばk=1の初回で、IIAに対しての重複面積が上位側に属すると判定された複数のマーカとして２つのマーカM10及びM20（それぞれ情報量I10及びI20が対応）があった場合に、方針１ではIIA範囲内での情報量の最大値が大きいことからマーカM10が選別されるのに対して、方針２ではIIA範囲外での値（重複面積または最大値で評価される値）が大きいことからマーカM20が選別される。
【0111】
（４） 変形例に係る実施形態として、探索範囲IIA（その横軸範囲）をどの程度狭めた場合にマーカをどれだけ減らせるかを提示するようにしてもよい。
【0112】
すなわち、例えば図１６に模式的に示されるように、探索範囲IIAとして、最も広い探索範囲IIA-1と、中間の探索範囲IIA-2と、最も狭い探索範囲IIA-3とを事前知識などにより予め設定しておき、それぞれの探索範囲についてサブセット決定部50の処理を適用することで、最も広い探索範囲IIA-1の場合の最も多いマーカ数n1（及び選別された各マーカ）と、中間の探索範囲IIA-2の場合の中間のマーカ数n2（及び選別された各マーカ）と、最も狭い探索範囲IIA-3の場合の最も少ないマーカ数n3（及び選別された各マーカ）とを、結果として出力するようにしてもよい。
【0113】
（５） 上記（４）と同様に、変形例に係る実施形態として、探索範囲をその横軸範囲に限らず、許容誤差が反映されるその縦軸範囲も含めて動的に変更した場合に必要となるマーカ個数（及び選別された各マーカ）を、提示するようにしてもよい。この場合も、図１６の模式例と同様に、動的に変更させた探索範囲IIAのそれぞれについて、サブセット決定部50の処理を適用して、結果のマーカ数（及び選別された各マーカ）を出力すればよい。
【0114】
（６） 変形例に係る実施形態として、最終的なマーカの個数を先に定めた場合に、探索範囲でどれだけの精度を保証できるかを結果として提示するようにしてもよい。
【0115】
前掲の図１４はこの変形例の模式例も示しており、例えば最終的なマーカの個数を２個として設定した場合には、例EX42に示されるように、選別された２個のマーカM3及びM5について、その情報量のグラフI3+I5を探索範囲IIAと重ねて表示することで、ユーザに対して視覚的に提示することにより、ユーザはマーカ２個の場合には探索範囲IIAにおいて、概ね半分程度の情報量が確保されることを把握することができる。
【0116】
（７） 探索範囲IIAについて、横軸範囲「θ上限≦θ≦θ上限拡張」に対して縦軸範囲は一定で設定する場合を説明したが、図１７の模式例に示される通り、事前情報などに基づいて縦軸範囲を目的形質θ（横軸）に対して変動しうるものとして設定してもよい。図１７の例EX51はこれまで説明した一定値の例であり、例EX52はθが高くなるほど情報量を減らして設定する例であり、例EX53はθが高くなるほど情報量を増やして設定する例であり、例EX54は中間のθで情報量の最大値を設定する例である。
【0117】
（８） 特許文献６と同様に、目的形質推定部33では、予め第１データベース15に記録される年齢等のラベルに関して、推定した目的形質（潜在特性や能力を意味する）θs（及び／又は推定の際に用いた項目パラメータai,bi）を、以下のような種々の用途や内容に関して出力するようにしてよい。
【0118】
＜潜在特性の予測サービス＞
対象者の検体から、前述の前処理部31により選択されたCpGサイトiのメチル化率βiを計測し、そこから計算された潜在特性を提示するサービスを実現できる。例えば、ラベルとして年齢（実年齢）をもとに構築されたモデルであれば、潜在特性として「DNAメチル化年齢」のようなものが推定できる。すなわち、この潜在特性としてのDNAメチル化年齢が意味するのは、加齢に伴うDNAメチル化の変化がどの程度進んでいるかを表す、実年齢とは別途の概念としての、生物学的な年齢のようなものであるものと考えることができる。（ただし、一般には実年齢と当該生物学的な年齢とは、対象者によって程度の違いはあっても、ある程度は関連があるものと考えられる。）また、ラベルとして喫煙履歴（後述するPack-year等の時間積分値）を対象にすれば、この潜在特性が意味するものは、「喫煙によるDNAメチル化影響度」のようなものであると考えることができる。同様に、ラベルとして前述した各種の運動履歴や心身の健康に影響しうる履歴や実績を用いる場合、対応する潜在特性が意味するものは、「当該履歴や実績によるDNAメチル化影響度」であると考えることができる。
対象者の検体から、前述の前処理部31により選択されたCpGサイトiのメチル化率βiを計測し、そこから計算された潜在特性を提示するサービスを実現できる。例えば、ラベルとして年齢（実年齢）をもとに構築されたモデルであれば、潜在特性として「DNAメチル化年齢」のようなものが推定できる。すなわち、この潜在特性としてのDNAメチル化年齢が意味するのは、加齢に伴うDNAメチル化の変化がどの程度進んでいるかを表す、実年齢とは別途の概念としての、生物学的な年齢のようなものであるものと考えることができる。（ただし、一般には実年齢と当該生物学的な年齢とは、対象者によって程度の違いはあっても、ある程度は関連があるものと考えられる。）また、ラベルとして喫煙履歴（後述するPack-year等の時間積分値）を対象にすれば、この潜在特性が意味するものは、「喫煙によるDNAメチル化影響度」のようなものであると考えることができる。同様に、ラベルとして前述した各種の運動履歴や心身の健康に影響しうる履歴や実績を用いる場合、対応する潜在特性が意味するものは、「当該履歴や実績によるDNAメチル化影響度」であると考えることができる。
【0119】
＜CpGサイト分析サービス＞
複数の利用者（サンプル）の検体を用いて、本実施形態の分析を行うことにより、各CpGサイトの分析を行うサービスを実現できる。例えば、ラベルとして年齢を対象にした場合では、推定された項目パラメータを用いて、以下のような情報A,B等を提示することが可能となる。
複数の利用者（サンプル）の検体を用いて、本実施形態の分析を行うことにより、各CpGサイトの分析を行うサービスを実現できる。例えば、ラベルとして年齢を対象にした場合では、推定された項目パラメータを用いて、以下のような情報A,B等を提示することが可能となる。
【0120】
●情報A…識別力パラメータが正の値となるCpGサイトに関して、以下の各パラメータを以下の通りの意味付けを有するものとして提示することができる。
困難度パラメータ → メチル化が起こりやすい年齢
識別力パラメータが大 → 加齢に伴い、急速にメチル化が進むCpGサイト
識別力パラメータが小 → 加齢に伴い、緩やかにメチル化が進むCpGサイト
困難度パラメータ → メチル化が起こりやすい年齢
識別力パラメータが大 → 加齢に伴い、急速にメチル化が進むCpGサイト
識別力パラメータが小 → 加齢に伴い、緩やかにメチル化が進むCpGサイト
【0121】
●情報B…識別力パラメータが負の値となるCpGサイトに関して、以下の各パラメータを以下の通りの意味付けを有するものとして提示することができる。
困難度パラメータ → 脱メチル化が起こりやすい年齢
識別力パラメータの絶対値が大 → 加齢に伴い、急速に脱メチル化が進むCpGサイト
識別力パラメータの絶対値が小 → 加齢に伴い、緩やかに脱メチル化が進むCpGサイト
困難度パラメータ → 脱メチル化が起こりやすい年齢
識別力パラメータの絶対値が大 → 加齢に伴い、急速に脱メチル化が進むCpGサイト
識別力パラメータの絶対値が小 → 加齢に伴い、緩やかに脱メチル化が進むCpGサイト
【0122】
情報A,Bのいずれについても、分析される対象者について潜在特性として推定された年齢との関係において提示するようにしてよい。すなわち、図８等で示したロジスティック関数のグラフ内での対象者の年齢の位置づけも示すことで、例えば情報Aであれば対象者自身が既にメチル化がどの程度進んでいるのか、あるいは進んでいないのか、という情報を提示するようにしてもよい。
【0123】
例えば前述した図１０のグラフG11は年齢を対象としたメチル化分析に用いた100のCpGサイトの項目特性曲線を表しており、CpGサイトごとに曲線が水平方向にばらついていることから、様々な困難度のCpGサイトが得られており、また、識別力（曲線の傾き）のバリエーションも存在することが確認できる。
【0124】
同様に、分析対象となるラベルとして年齢ではなく、1日の喫煙量×喫煙年数などで定量化されるPack-yearを用いた場合は、以下のような情報C,Dを提示することができる。前述と同様に、ロジスティック関数のグラフ上で対象者の位置づけと関連付ける形でこうした情報C,Dを提示してもよい。
【0125】
●情報C…識別力パラメータが正の値となるCpGサイトに関して、以下の各パラメータを以下の通りの意味付けを有するものとして提示することができる。
困難度パラメータ → メチル化が起こりやすい喫煙量
識別力パラメータが大 → 喫煙に伴い、急速にメチル化が進むCpGサイト
識別力パラメータが小 → 喫煙に伴い、緩やかにメチル化が進むCpGサイト
困難度パラメータ → メチル化が起こりやすい喫煙量
識別力パラメータが大 → 喫煙に伴い、急速にメチル化が進むCpGサイト
識別力パラメータが小 → 喫煙に伴い、緩やかにメチル化が進むCpGサイト
【0126】
●情報D…識別力パラメータが負の値となるCpGサイトに関して、以下の各パラメータを以下の通りの意味付けを有するものとして提示することができる。
困難度パラメータ → 脱メチル化が起こりやすい年齢
識別力パラメータの絶対値が大 → 喫煙に伴い、急速に脱メチル化が進むCpGサイト
識別力パラメータの絶対値が小 → 喫煙に伴い、緩やかに脱メチル化が進むCpGサイト
困難度パラメータ → 脱メチル化が起こりやすい年齢
識別力パラメータの絶対値が大 → 喫煙に伴い、急速に脱メチル化が進むCpGサイト
識別力パラメータの絶対値が小 → 喫煙に伴い、緩やかに脱メチル化が進むCpGサイト
【0127】
以上の通り、本実施形態では新規にIRTモデルを導入したことにより、DNAメチル化分析手法において以下のような利点A1,A2を提供することができる。
【0128】
●利点A1 前処理部31の処理（CpGサイトの絞り込み処理）のみで、ラベル付きデータを必要とするため、後続の項目パラメータ推定部32では、ラベル無しデータも利用可能となる。通常、エピゲノム計測の費用は高額であるため、別の研究目的で計測されたデータ等で、研究対象のラベルが無いデータの再利用性が高まることの経済的効果は大きい。（例えば、加齢の研究で収集されたエピゲノムデータのため喫煙歴の情報が収集されていないような場合であっても、喫煙の影響の分析に使うことが可能となる。）
【0129】
すなわち、ラベルを必要とするのは第１データベース15のみであり、第２データベース22にはラベルを必要としない。（なお、図２に線L16で示される以上の実施形態に対する変形例として、図２に点線L15で示されるように、項目パラメータ推定部22では第２データベース22を参照することに代えて、第１データベース15をラベル情報を無視したうえで参照することも可能である。）
【0130】
図１０のグラフG12は、年齢を対象にしたデータを用いた実験結果であり、フィルタリング時に用いるラベル付きデータの数（横軸）と、後続のIRTにより推定される潜在特性と暦年齢との相関（縦軸）を求めたグラフである。ここでは、項目パラメータ推定のために利用するラベル無しデータの数は、1138名分一定とした。本グラフより、10名のラベル付きデータでは相関係数が0.55程度と低いが、20名程度以上のラベル付きデータがあれば、相関係数が0.85付近に飽和し、精度良くフィルタリングが行えることが示されている。
【0131】
●利点A2 IRTを用いることにより、各CpGサイトに対して、困難度、識別力パラメータを明らかにすることができる。例えば、加齢の研究においては、Elastic Net等の従来法では、各CpGサイトの加齢への関連性の強さのみしか分からないのに対して、本発明手法では、情報A,B,C,D等として前述したように、例えば困難度パラメータにより、「何歳ごろにメチル化が進みやすいCpGサイトであるか」や、識別力パラメータにより「CpGサイトのメチル化の進み具合が、緩やかであるか急峻であるか」などの情報が得られる。このような情報は、CpGサイト近傍の遺伝子と関連性を見る際に非常に有益な情報である。例えば、加齢における身体への影響を調査する目的では、高齢でメチル化が進むCpGサイト近傍の遺伝子を重点的に分析するなどの目的に利用可能である。
【0132】
（９） 前処理部11では、図７のステップS101～S105として、CpGサイトiのうちラベルと相関があると判定されるものを選出したが、相関に限らずより一般に、任意の既存手法によって推定することで、あるいは、事前知識として既知の関係を参照して利用することで、ラベルと関連があると判定されるものを選出するようにしてよい。例えば、ラベルが年齢（実年齢）である場合には、非特許文献１等の既存手法で、当該年齢と関連性があると判定されるCpGサイトを特定することができる。前処理部11では当該既存手法で予め得られている結果を参照して用いることで、この年齢と関連性のあるCpGサイトを用いるようにしてよい。
【0133】
（１０） 図１８は、一般的なコンピュータ装置70におけるハードウェア構成の例を示す図である。情報処理装置100は、このような構成を有する１台以上のコンピュータ装置70として実現可能である。なお、２台以上のコンピュータ装置70で情報処理装置100を実現する場合、ネットワーク経由で処理に必要な情報の送受を行うようにしてよい。コンピュータ装置70は、所定命令を実行するCPU（中央演算装置）71、CPU71の実行命令の一部又は全部をCPU71に代わって又はCPU71と連携して実行する専用プロセッサとしてのGPU（グラフィックス演算装置）72、CPU71（及びGPU72）にワークエリアを提供する主記憶装置としてのRAM73、補助記憶装置としてのROM74、通信インタフェース75、ディスプレイ76、マウス、キーボード、タッチパネル等によりユーザ入力を受け付ける入力インタフェース77と、これらの間でデータを授受するためのバスBSと、を備える。
【0134】
情報処理装置100の各機能部は、各部の機能に対応する所定のプログラムをROM74から読み込んで実行するCPU71及び／又はGPU72によって実現することができる。なお、CPU71及びGPU72は共に、演算装置（プロセッサ）の一種である。ここで、表示関連の処理が行われる場合にはさらに、ディスプレイ76が連動して動作し、データ送受信に関する通信関連の処理が行われる場合にはさらに通信インタフェース75が連動して動作する。第１データベース15及び第２データベース22の記憶機能は、補助記憶装置としてのROM74として実現してよい。
【符号の説明】
【0135】
100…情報処理装置、10…第１構築部、20…第２構築部、30…解析部、40…IRTパラメータ算出部、50…サブセット決定部
11…エピゲノムデータセット入力部、12…エピゲノムデータ解析部、13…ラベル入力部、14…対応付け保存部、15…第１データベース
21…エピゲノムデータセット入力解析部、22…第２データベース
31…前処理部、32…項目パラメータ推定部、33…目的形質推定部
51…範囲設定部、52…マーカ選別部
11…エピゲノムデータセット入力部、12…エピゲノムデータ解析部、13…ラベル入力部、14…対応付け保存部、15…第１データベース
21…エピゲノムデータセット入力解析部、22…第２データベース
31…前処理部、32…項目パラメータ推定部、33…目的形質推定部
51…範囲設定部、52…マーカ選別部
【図1】
【図2】
【図3】
【図4】
【図5】
【図6】
【図7】
【図8】
【図9】
【図10】
【図11】
【図12】
【図13】
【図14】
【図15】
【図16】
【図17】
【図18】