Comments
Description
Transcript
進化アルゴリズムを用いた分子間相互作用ネットワークの効率的な推定
平成 27 年度 博士論文 進化アルゴリズムを用いた分子間相互作用 ネットワークの効率的な推定手法の開発 古森 朝子 目次 目次 ........................................................................................................................ Ⅰ 第 1 章 諸論 .......................................................................................................... 1 第 1 節 研究の背景・目的 ................................................................................ 1 第2節 第2章 論文の構成 ........................................................................................... 4 ネットワーク推定手法 第1節 緒言 第2節 S-system モデル 第3節 数値最適化 ............................................................................ 6 ....................................................................................................... 6 .................................................................................... 6 ........................................................................................... 8 第 1 項 実数値遺伝的アルゴリズム 第 2 項 目的関数 第4節 S-system モデルのシステムパラメータ最適化手法 ........................ 12 第 1 項 UNDX+MGG 法 第 2 項 SPX+MGG 法 第3項 AREX+JGG 法 第4項 ハイブリッド法 第5節 検証 ..................................................................................................... 22 第 1 項 AREX+JGG 法の有効性の検証 第 2 項 ハイブリッド法の有効性の検証 第3項 新規法の推定可能パラメータ数の検証 第4項 タイムコースデータ数の検証 第6節 第3章 第1節 考察 ..................................................................................................... 36 Random selection method 緒言 ........................................................................ 38 ..................................................................................................... 38 I 第2節 Random selection method のアルゴリズム 第3節 評価方法 第4節 Random selection method の有効性の検証 ......................................... 39 ............................................................................................. 42 ......................................... 43 第 1 項 ネットワークモデル 第 2 項 30 要素ネットワークモデルの推定 第5節 第4章 考察 ..................................................................................................... 52 ノイズを含むデータへの有効性 .......................................................... 54 第1節 緒言 第2節 ノイズを含むタイムコースデータの取得 第3節 ネットワークモデル 第4節 検証 第1項 ..................................................................................................... 54 ....................................... 54 .......................................................................... 55 ..................................................................................................... 58 AREX+JGG 法及びハイブリッド法のノイズを含むタイムコースデータに対す る有効性 第 2 項 許容される最大ノイズの検証 第 3 項 ノイズを含むデータに対する Random selection method の有効性 第5節 考察 ..................................................................................................... 72 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 ................................................................................................................................... 75 第1節 緒言 第2節 ランダムネットワークとスケールフリーネットワーク ..................................................................................................... 75 ................ 76 第 1 項 ランダムネットワーク 第 2 項 スケールフリーネットワーク 第3節 100 要素ネットワークモデルの作成 ................................................ 77 第 1 項 100 要素ランダムネットワークモデル 第 2 項 100 要素スケールフリーネットワークモデル 第4節 100 要素ネットワークモデルの推定 ................................................ 89 第 1 項 100 要素ランダムネットワークモデル推定結果 第 2 項 100 要素スケールフリーネットワークモデル推定結果 第5節 考察 ..................................................................................................... 94 II 第6章 総括 ........................................................................................................ 99 第1節 推定時間 第2節 TP interaction と FP interaction 第3節 今後の展望 参考文献 謝辞 ........................................................................................... 103 ......................................................... 104 ....................................................................................... 107 .............................................................................................................. 108 ........................................................................................................................ Ⅳ III 第1章 緒論 第1章 緒論 第1節 研究の背景・目的 システム生物学は、生命をシステムレベルで理解することを目標とする研究領域であ り、2000 年頃を境に研究が活発化している(Ideker, Galitski, and Hood, 2003; Kitano, 2002a; Kitano, 2002b)。システム生物学では、個々の遺伝子やタンパク質等の生体分子 に着目するのではなく、生命現象を構成要素である生体分子の相互作用システムとして 捉える。近年のシステム生物学研究の牽引力となっているのは、質量分析計をはじめと するオミクスの測定技術の発達だろう。具体的には、ゲノムプロジェクト以降、トラン スクリプトームやプロテオーム、メタボロームなどの各種オミックスデータの測定が可 能となり、網羅的なデータが蓄積されている。これらの膨大なデータによって、分子レ ベルの情報を利用してシステムレベルで解析できる段階へと至った。 システム生物学の研究は大きく以下の 4 つの分野に分けられる(Kitano, 2002a) 。 1.システムの構造:生体分子のネットワーク構造を理解する 2.システムの動的挙動:時間経過により変動するシステムの定量的解析方法やモデ ル・理論の確立 3.システムの制御方法:システムを特定の状態へと誘導する制御理論の確立 4.システムの設計方法:特定の挙動を再現するシステムの設計 システム生物学分野の研究は、生命現象を理解することだけでなく、創薬ターゲットの 探索やオーダーメイド医療などへの応用が期待されている。 本研究は、生体分子の相互作用ネットワーク構造を理解することを主眼としている。 生命現象をシステムとして理解するために、ネットワーク構造を構築することは非常に 1 第1章 緒論 重要である。相互作用ネットワーク構造を理解するために、実験的に得られたタイムコ ースデータから相互作用ネットワークを再構築する研究が盛んに行われている(Bansal et al., 2007; Zhu et al., 2012) 。分子間相互作用ネットワーク推定手法の一つに、実験から 得られた生体分子の濃度のタイムコースデータを、微分方程式モデルである S-system (Savageau, 1976)を用いて再現し、モデルに含まれるパラメータから分子間相互作用 ネットワーク構造を推定する手法がある。この手法は、フィードバックなど、ループ構 造を含む複雑な構造も推定でき、推定されたネットワークの各システムパラメータの感 度解析が行えるなどの利点を持つ。S-system は、一般質量作用則(Generalized Mass Action law; GMA)の生成および分解過程をそれぞれ 1 本のパスで近似したものであり、GMA よりも推定しなければならないパラメータが少ないという利点を持つ。しかしながら、 そのパラメータ数の合計は要素数を n とすると 2n(n+1)個となり、特に要素数が多くな る場合、膨大な数のシステムパラメータを推定する必要がある。また、観測されたタイ ムコースデータを再現できるように微分方程式を推定する問題は逆問題( inverse problem)であり、相互作用ネットワークを一意に決定することができない。富永らは、 これらの問題に対処するために S-system のシステムパラメータの最適化に伝統的な Simple Genetic Algorithm(Simple GA; SGA) (Goldberg, 1989)を適用して遺伝子相互作 用ネットワークを推定する手法を開発した(富永と岡本, 1999; Okamoto et al., 2000) 。し かしながら、SGA は初期収束や進化的停滞といった問題が指摘されている(Baker, 1987; Whitley, 1989)。そこで、上田らは、SGA に代わる数値最適化手法として実数値遺伝的 アルゴリズム(実数値 GA)(Eshelman and Schaffer, 1993)を適用した手法を開発した (Ueda, 2001; Ueda et al., 2002a; Ueda et al., 2002b; 岡本と小野, 2003) 。上田らは、代表的 な実数値 GA の一つである、Unimodal Normal Distribution crossover(UNDX)(Ono, Kobayashi, and Yoshida 1998; 小野ら 1999) と Minimal Generation Gap(MGG) (Satoh, 1996; 佐藤ら, 1997)を用いた UNDX+MGG 法を開発し、5 遺伝子で構成される遺伝子相互作 用ネットワークの推定に成功した。しかしながら、UNDX+MGG 法では、推定可能なパ ラメータ総数は 20 程度に留まっており、実用化には不十分である。そこで、本研究で 2 第1章 緒論 は、より大規模な相互作用ネットワークを効率的に推定するための手法を開発すること を目的とした。 相互作用ネットワーク推定手法を生命システムの解析に用いる場合には、推定可能な 要素数が多ければ多いほど良い。しかしながら、UNDX+MGG 法で推定可能な要素数は 3 つ程度である。そこで、本研究では、まず、30 要素相互作用ネットワークの推定を可 能とすることを目標とした。生命システムにおいて 30 要素程度で形成されるシステム としては、ガンの研究で重要とされている細胞周期チェックポイントの一つである G1/S チェックポイントに関与する遺伝子が挙げられる(Wang et al., 2015) 。従って、30 要素相互作用ネットワークを推定できる手法であれば、このシステムの解析に用いるこ とができる。従来開発されてきた手法を、より大規模な相互作用ネットワークの推定へ と適用するために、本研究では以下の二つのアプローチを採った。 1.新たな実数値 GA を採用した新規法を開発し推定可能なパラメータ数を増大する 2.新規法の適用の仕方を変えることで適用可能な相互作用ネットワークの規模を増 大する まず、一つ目のアプローチのために、新たな実数値 GA として Adaptive Real-coded Ensemble crossover(AREX) (Akimoto et al., 2009)と Just Generation Gap(JGG) (Akimoto et al., 2007; Kobayashi, 2009)を採用し、AREX+JGG 法を開発した(Komori et al., 2012) 。 更に、AREX+JGG 法の収束速度をより早めるために AREX+JGG 法に修正 Powell 法 (Powell, 1964)を導入したハイブリッド法を開発した。これらの新規法は、従来法と 比較して推定可能パラメータ数と推定時間の両方で優れていたが、現在までに推定でき た相互作用ネットワーク構成要素数は 17 要素に留まっている。そこで、ハイブリッド 法を相互作用ネットワーク推定に適用する際のアルゴリズムを開発し、より大規模な相 互作用ネットワークの推定を目指すこととした。 推定可能なパラメータ数が限られる手法を大規模な相互作用ネットワークに適用す 3 第1章 緒論 る方法の一つとして、推定対象としている相互作用ネットワークをクラスタリングする 方法があげられる。この方法は、ネットワーク構成要素をいくつかのクラスタに分け、 各クラスタを推定し大規模な相互作用ネットワークを導出するというアルゴリズムで ある。しかし、適切なクラスタリングを行うためには、各構成要素の濃度のタイムコー スの形状や、各要素に関するデータベース情報などの情報が求められることが多い。そ こで、本研究では元の相互作用ネットワークを構成する要素から、ランダムに要素を選 択し、選ばれた要素で構成されたサブネットワークを推定することで大規模な相互作用 ネットワークを導出する Random selection method(RA 法)を開発した(Komori et al., 2013)。具体的には、ランダムに推定する要素を選択し、選択した要素間の推定に新規 法を適用して、その後、推定した要素間の共通構造を重ねあわせるという手法である。 RA 法では、サブネットワークの構成要素はランダムに決定されるため、人為的な操作 を必要とせず恣意性の低い結果を得ることができる。本研究では、開発した RA 法が新 規法(AREX+JGG 法やハイブリッド法)単独では推定できない相互作用ネットワーク に対して有効であるか精査した。 第2節 論文の構成 第 2 章では、相互作用ネットワーク推定手法で用いられている S-system モデルや実 数値 GA について述べる。先行研究で開発された UNDX+MGG 法や SPX(Simplex cross over)+MGG 法、新規法である AREX+JGG 法やハイブリッド法の原理やアルゴリズム を紹介する。そして、小規模なネットワークモデルを用いて従来法と新規法の性能を比 較する。 第 3 章では、RA 法のアルゴリズムおよび有効性の検証を行う。新規法単独では推定 することのできない 30 要素ネットワークモデルを作成し、RA 法を適用することで、相 互作用を推定できるか検証する。 第 4 章では、各相互作用ネットワーク推定手法のノイズ耐性について検証する。まず、 4 第1章 緒論 小規模なネットワークモデルを用いて UNDX+MGG 法、 SPX+MGG 法、 AREX+JGG 法、 ハイブリッド法のノイズ耐性を調べ、有効性を保てる最大誤差の大きさの検討を行う。 そして、30 要素ネットワークモデルを用いて RA のノイズ耐性を調べた。 第 5 章では、RA 法が様々な特徴を持つような相互作用ネットワークに対しても有効 であるか検証する。検証には、複雑ネットワークの代表的なモデルであるランダムネッ トワークおよびスケールフリーネットワークの特徴を有する 2 つの 100 要素ネットワー クモデルを用いる。そして、我々の身の回りに存在するようなネットワーク構造に対し ても RA 法か有効であるか検証した。 第 6 章では、総合的な考察を行う。 5 第2章 第2章 ネットワーク推定手法 第1節 諸言 ネットワーク推定手法 実験的に得られた生体分子の時系列発現データ(タイムコースデータ)から、分子間 相互作用ネットワークを推定することは、システム応答(アウトプット)からシステム の内部構造(ネットワーク構造)を推定することであるため、数学的には逆問題(inverse problem)であり、不良設定問題(ill-posed problem)に属する。したがって、データを 再現できるようなネットワーク構造は多数あると考えられ、ネットワーク構造を一意に 決定することはできない。このような問題に対処するために、これまで、実験から得ら れたタイムコースデータを微分方程式モデルで再現し、モデルに含まれているパラメー タから相互作用ネットワークを推定する手法が開発されてきた。本研究で開発した手法 では、数理モデルには S-system モデル(Savageau, 1976)を採用し、S-system に含まれ る未知パラメータの最適化には遺伝的アルゴリズム(Genetic Algorithm; GA)を採用し た(北野, 1993; 伊庭, 1994)。本章では、開発してきた手法で用いられている S-system や GA について述べる。 第2節 S-system モデル 観測された生体分子のタイムコースデータを微分方程式で再現するための数理モデ ルとして、一般質量作用則(generalized mass action; GMA)が広く用いられている。一 般質量作用則は次の式(1.1)のように表される。 6 第2章 dxi = dt p q n αik k=1 gijk xj − j=1 ネットワーク推定手法 n hijk βik k=1 xj (i = 1, 2, … , n) (1.1) j=1 ここで、p は状態変数𝑥𝑖 の生成経路の数、q は状態変数𝑥𝑖 の分解経路の数である。𝑛は状 態変数の総数、𝛼𝑖k および𝛽𝑖kはそれぞれ𝑥𝑖 の生成経路 k、分解経路 k の速度定数を表す。 GMA では各生体分子の反応の詳細を記述することができるが、個々の反応経路やメカ ニズムが不明である場合に記述が困難となる。したがって、本研究で扱うような逆問題 を扱う場合には、GMA を近似した S-system(Savageau, 1976)が有効である。S-system は GMA の生成過程、分解過程をそれぞれ 1 本のパスで表現したものであり、式(1.2)の ように表記される。 dxi = αi dt n gij n hij xj − βi j=1 xj (i = 1, 2, … , n) (1.2) j=1 右辺第一項は生成過程、第二項は分解過程を表している。𝑥𝑖 は状態変数、𝑛は状態変数 の総数であり、𝛼𝑖 および𝛽𝑖 は、それぞれ𝑥𝑖 の生成過程、分解過程の見かけの速度定数で ある。gij は状態変数𝑥𝑖 の生成過程に関与する状態変数𝑥𝑗 の相互作用係数、hij は状態変数𝑥𝑖 の分解過程に関与する状態変数𝑥𝑗 の相互作用係数であり、それぞれの符号および値によ って、状態変数間の相互作用の種類を表現することができる。たとえば、gij の符号が正 ならば𝑥𝑗 が𝑥𝑖 の生成過程を促進する相互作用を、負ならば抑制する相互作用をそれぞれ 表している。値が 0 の場合は、𝑥𝑗 は𝑥𝑖 の生成過程に影響を与えない。同様に、hij は𝑥𝑗 が𝑥𝑖 の 分解過程に及ぼす影響であり、hij の符号が正ならば𝑥𝑗 が𝑥𝑖 の分解を促進する相互作用を、 負であるならば抑制する相互作用をそれぞれ表している。値が 0 の場合は、𝑥𝑗 は𝑥𝑖 の分 解過程に影響を与えないことになる。 S-system は 、 生 成 過 程 と 分 解 過 程 を そ れ ぞ れ 一 つ の 項 で 表 現 す る 抽 象 モ デ ル (conceptual modelling)である。したがって、生成項や分解項が複数の経路で構成され 7 第2章 ネットワーク推定手法 ている場合、GMA を近似した表現となる。しかしながら、未知の分子間相互作用ネッ トワークを推定する際、それぞれの生体分子の生成過程や分解過程の詳細な機構が明ら かになっていない場合が多く、GMA による記述は困難である。S-system は、ある反応 種が別の反応種の生成もしくは分解過程にどのような影響を与えているか概念的に記 述しているため、構成要素間の反応形式や化学量論式が明確でないような問題において も記述が可能である。 本研究では、まず、状態変数(𝑥𝑖 )の動的挙動(タイムコース)が与えられ、そして その挙動を説明できるような状態変数間の実数値パラメータの推定を行う。つまり、観 測されたタイムコースデータを再現できるように式(1.2)に含まれる𝛼𝑖 、𝛽𝑖 、gij 、hij を探 索することで分子間相互作用ネットワークの推定を行う。これらの未知パラメータの総 数は、ネットワークを構成する要素数をn個とするとき、2n(n+1)個となり、nが増加す るに従いO(n2 )で増加する。そのため、未知パラメータの数がO(n3 )で増加する GMA と 比べ、大規模な相互作用ネットワークの推定への適用も容易である。 第3節 第1項 数値最適化 実数値遺伝的アルゴリズム GA は、1975 年にジョン・ヘンリー・ホランドによって提案された、生物の進化の過 程を模倣した進化アルゴリズムの一つであり(北野, 1993; 伊庭, 1994)、選択淘汰や突 然変異などの操作を繰り返すことで最適化対象のパラメータの値の近似解を探索する 手法である。GA は多峰性の関数最適化問題に対して優れた実用解を比較的効率よく求 めることができる有力な近似解法であることが知られている(Akimoto et al., 2010)。 GA では、最適化対象となるパラメータ値は個体の遺伝子として表現され、一般的に遺 8 第2章 ネットワーク推定手法 伝子は 2 進数で記述される。しかしながら、解決すべき問題によって 2 進数表現が適切 な形式とは限らないため、最適な形式で遺伝子を表現することが必要である。 遺伝子表現として実数値ベクトルを扱うことを可能とした GA として、実数値 GA (Real-Coded GA)が挙げられる。GA の主な手順を Fig. 2.1 に示す。ある世代(generation) は、集団サイズ(population size)に基づいた個体(individual)の集団(population)か ら構成される。集団サイズは多様性(diversity)を制御する重要なパラメータであり、 通常、世代を通じて一定である。各個体は複数個の遺伝子(gene)の集まりで構成され る染色体(chromosome)によって特徴付けられている。目的関数を適応度(fitness)と し、個体に対して遺伝的操作を行うことで集団を適応度の高いものへと進化させていく。 進化の枠組みとしては、まず、第 t 世代の個体集団の内、交叉に用いる親個体を決める ために複製選択(reproduction)が行われる。選択された親個体に対して、交叉(crossover) を繰り返し適用することで親個体の近傍に子個体を生成する。親個体群と子個体群の中 から次世代に残す個体を生存選択(survival selection)によって選び、個体集団へ戻すこ とで第 t+1 世代集団とする。これを繰り返すことで集団を有望領域へと収束させる。複 製選択と生存選択を合わせて世代交代モデル(generation alternation model)と呼ぶ。実 数値 GA は、実数値ベクトルでの個体表現を可能とした GA であり、現在、実数値 GA の様々な交叉オペレータや世代交代モデルが提案されている(Davis, 1990; Wright, 1991; 坂和と田中, 1995; Kobayashi, 2009) 。 9 第2章 第 t 世代集団 複製選択 親個体群 交叉 子個体群 生存選択 第 t+1 世代集団 Fig. 2.1 GA の手順 10 ネットワーク推定手法 第2章 第2項 ネットワーク推定手法 目的関数 実数値 GA では生成された各個体を目的関数から算出した適応度に準じて評価する。 実数値 GA を適用する問題により様々な目的関数が設定され、適切な目的関数を定める ことが求められる。本研究では、この目的関数として、観測されたタイムコースデータ と実数値 GA を適用した結果求められた数理モデルから計算したタイムコースデータ との相対二乗誤差の和である累積相対二乗誤差を採用した。累積相対二乗誤差 E は、次 の式(1.3)で表される。 𝐷 𝑁 𝑇 𝑒𝑥𝑝 𝐸 = 𝑑=1 𝑖=1 𝑡=1 2 𝑐𝑎𝑙 𝑥𝑑,𝑖,𝑡 − 𝑥𝑑,𝑖,𝑡 ( ) 𝑒𝑥𝑝 𝑥𝑑,𝑖,𝑡 (1.3) ここで、D は観測されたタイムコースデータの総数で、N は観測される状態変数の数、 𝑐𝑎𝑙 T は観測されたタイムコースデータのサンプリングポイント数である。𝑥𝑑,𝑖,𝑡 はタイムコ 𝑐𝑎𝑙 ースデータ d、時刻 t における状態変数𝑥𝑖 の実験値(観測値)を表し、𝑥𝑑,𝑖,𝑡 はタイムコ ースデータ d、時刻 t における状態変数𝑥𝑖 の計算値を表している。実用的には、累積二 乗誤差のタイムコースデータ 1 点当たりの平均値 Eave を算出し、評価値としている。Eave は次の式で表される。 Eave E D N T (1.4) つまり、与えられたタイムコースデータと推定の結果から得られるタイムコースデータ を比較し、各データの 1 点当たりの誤差を指標とし、定めた閾値以下となった時収束し たと判断した。 11 第2章 第4節 第1項 ネットワーク推定手法 S-system モデルのシステムパラメータ最適化手法 従来法と提案法 S-system に含まれる未知パラメータの数は遺伝子数を n 個とするとき 2n(n+1)個であ る。GMA と比較するとモデルに含まれるシステムパラメータは少ないとはいえ、n の 増加に伴い O(n2)で増加する。従って、大規模なネットワークの推定には効率的に未知 パラメータを決定する探索アルゴリズムの開発が求められる。富永らは、simple GA (SGA)(Goldberg, 1989)を基本とする多次元非線形数値最適化法を開発し、遺伝子 相互作用ネットワーク推定問題に適用し効率的なネットワークの推定に成功した(富永 と岡本, 1999; Okamoto et al., 1998)。上田らは、この手法を更に発展させ、実数値 GA を用いた相互作用ネットワーク推定手法を開発した(Ueda, 2001; Ueda et al., 2002a; Ueda et al., 2002b)。 ここでは、先行研究で用いられてきた UNDX+MGG 法および SPX+MGG 法に加え、 新規法である AREX+JGG 法および AREX+JGG 法と修正 Powell 法を組み合わせたハイ ブリッド法の各アルゴリズムについて述べる。 第1項 UNDX+MGG 法 世代交代モデルである MGG(Minimal Generation Gap) (Satoh, 1996; 佐藤, 小野と小 林, 1997)と交叉オペレータである UNDX(Unimodal Normal Distribution crossover) (Ono, Kobayashi, and Yoshida 1998; 小野ら 1999)は相性の良い組み合わせであり、最も有名 な実数値 GA の一つである。上田らは、UNDX+MGG を S-system のパラメータ最適化に 適用した遺伝子相互作用ネットワーク推定手法を開発し、良好な結果を得ている(Ueda, 2001; Ueda et al., 2002a; Ueda et al., 2002b; 岡本と小野, 2003; Nakatsui et al., 2008)。 12 第2章 ネットワーク推定手法 MGG は実数値 GA の代表的な世代交代モデルである。MGG は選択序盤における選択 圧をできるだけ最小に留めることで初期収束を、探索後期においても多種多様な個体を 生存させることで進化的停滞(evolutionary stagnation)を効率的に回避することができ るとされている(Kobayashi, 2009) 。MGG は 2 個体を親として用いる交叉を想定してお り、以下のようなアルゴリズムで行われる。 1. 個体集団から 2 個体をランダムに非復元抽出する(複製選択)。 2. 選択した 2 個体を親として交叉を繰り返し適用して、子個体群を生成する。 3. 全ての子個体について、評価値の計算を行う。 4. 子個体群と親である 2 個体の中から、1 個体は評価値が最良のものを、もう 1 個 体はランキングに基づくルーレット選択により選択する(生存選択) 。 5. 1 で選択した親と、4 で選択した個体とを入れ替え、次世代集団とする。 MGG は子個体に親個体も加えた集団の内、最良の個体を次世代に残すエリート保存戦 略を採用することで、進化的停滞を防ぐ効果を期待している。また、ルーレット選択に より次世代に残す個体を選択することで集団の多様性を維持させている。 UNDX では、Fig. 2.2 のように、両親を結ぶ主軸を主探索成分、主軸に直行する超平 面を副探索成分として正規分布によって子個体を生成する。具体的には、次の式(1.5) のように表現される。 𝑙 𝒄1 = 𝒎 + 𝑧1 𝒆1 + 𝑙 𝑧𝑘 𝒆𝑘 , 𝒄2 = 𝒎 − 𝑧1 𝒆1 − 𝑘=2 𝑧𝑘 𝒆𝑘 𝑘=2 𝒎 = 𝑷1 + 𝑷2 / 2 𝑧1 ~ 𝑁 0, 𝜎12 , 𝑧𝑘 ~ 𝑁 0, 𝜎22 𝑘 = 2, 3, ⋯ , 𝑙 𝜎1 = 𝛼𝑑1 , 𝜎2 = 𝛽𝑑2 / 𝑙 𝒆1 = 𝑷1 + 𝑷2 / 𝑷1 − 𝑷2 𝒆𝑖 ⊥ 𝒆𝑗 𝑖 ≠ 𝑗 𝑖, 𝑗 = 1, 2, ⋯ , 𝑙 13 (1.5) 第2章 P3 σ2 ネットワーク推定手法 σ1 P2 d2 P1 d1 Fig. 2.2 UNDX の概略 ここで、𝑙は次元数、𝒄1 および𝒄2は子個体ベクトル、𝑷1 および𝑷2は親個体ベクトル、𝒎は 𝑷1 と𝑷2 の中点、𝑑1 は𝑷1 と𝑷2 の距離、𝑑2 は第 3 の親𝑷3 と両親を結ぶ主軸との距離、 𝑧1 ~ 𝑁 0, 𝜎12 と𝑧𝑘 ~ 𝑁 0, 𝜎22 は正規乱数を表している。𝜎1 と𝜎2 はそれぞれ主探索方向と 副探索方向の分布を決めるパラメータで、𝛼 = 0.5 と𝛽 = 0.35 が推奨されている。𝒆1 は 主軸方向の単位ベクトルであり、𝒆𝑗 は𝒆1 に垂直かつ線形独立な単位ベクトルである。 UNDX は、目的関数の改善を効率的に行うために複数の変数を同時に変化する必要があ る性質である変数間依存性に対処することを目的として設計された(Kobayashi, 2009)。 第2項 SPX+MGG 法 SPX(Simplex crossover)(Tsutsui, Yamamura, and Higuchi, 1999; Higuchi, Tsutsui, and Yamamura, 2001)は、次元数を l とした場合、l+1 個体の親個体を用いる多親交叉法で ある。SPX では、l+1 個の親個体を頂点とした l 次元単体(simplex)の重心を中心とし て、単体を ε 倍に相似的に拡張した領域の内部に子個体を生成する(Fig. 2.3) 。UNDX の子個体の分布は正規分布に従うが、SPX では子個体は領域内に一様に生成される。 SPX のアルゴリズムは以下のとおりである。 14 第2章 ネットワーク推定手法 1 l Pi l 1 i 0 xk G ( Pk G ) G (k 0) 0 Ck (k 1, , l ) rk 1 ( xk 1 xk C k 1 ) C xl Cl (1.6) 1 rk (u (0,1)) k 1 ここで、G は親個体群 P の重心、C は子個体群で C は子ベクトル、rk は区間[0,1]内の一 様乱数 u(0,1)から返還された乱数である。ε は拡張率で正のパラメータとなる。ε の推奨 値は l 2 である。 SPX+MGG 法では、世代交代モデルに前項で述べた MGG を採用している。MGG で は 2 親交叉を想定して設計された世代交代モデルであるため、SPX のような多親交叉に MGG を用いる場合は、親個体数を np とすると、np-2 個の親個体(つまり l-1 個)は そのまま次世代に残されることとなっている。 P0 G P1 P2 0 ε 1 Fig. 2.3 親個体が 3 個体のときの SPX の子個体生成範囲 15 第2章 第3項 ネットワーク推定手法 AREX+JGG 法 従来法による遺伝子相互作用ネットワークの推定は、推定可能なパラメータ数が限ら れており、計算に膨大な時間を要していた。これまで筆者らは、これらの問題を改善す るために新規実数値 GA として世代交代モデルである JGG(Just Generation Gap) (Kobayashi, 2009; Akimoto et al., 2007)と交叉オペレータである AREX(Adaptive Real-coded Ensemble crossover) (Akimoto, Sakuma, et al., 2009)を適用した AREX+JGG 法を提案した(Komori et al., 2012)。 MGG は二親交叉を対象として設計されているが、二親交叉は変数間でスケールの異 なるような悪スケール性を持ち、変数間依存のある問題において有効に働かないことが 知られている(Akimoto et al., 2010) 。悪スケール性や変数間依存性は実問題において多 く観察され、このような場合には一般的に多親交叉を用いることが望まれる。秋本らは、 これらの問題点を克服し多親交叉を想定した世代交代モデルとして MGG を拡張した JGG を提案した(Akimoto et al., 2007) 。JGG は以下のようなアルゴリズムで行われる。 1. 個体集団から𝑛𝑝 個の個体をランダムに抽出する(複製選択)。 2. 選択した𝑛𝑝 個体を親として交叉を繰り返し適用して、子個体群を生成する。 3. 全ての子個体について、評価値の計算を行う。 4. 子個体群の中から、親個体と同数の μ 個体をランキングに基づき選択する(生存 選択) 。 5. 1 で選択した親と、4 で選択した個体とを入れ替え、次世代集団とする。 JGG では、l を次元数としたとき、親個体数 np=l+1 とすることが推奨されている。JGG では、集団サイズを npop とすると、集団全体の世代交代率が np/npop となり、収束速度が 向上することが期待されている。また、JGG では評価値のランキングが良い個体のみを 生存選択において選択している。 16 第2章 ネットワーク推定手法 多親交叉の実数値 GA としては、前項で述べた SPX や UNDX-m(喜多ら, 2000)が知 られている。UNDX-m は、UNDX を拡張した手法であり、基本的なアルゴリズムは UNDX と変わらないが m 個の親個体の遺伝子を基に子個体が生成される。m の推奨値は 4~5 である。SPX は前項で述べたとおり、l+1 個の親個体が基とされた領域内に一様に子個 体を生成するが、子個体群の分布は親個体群の中心に対して対称ではない。従って、変 数間の依存性が強い問題において性能が低下する可能性がある(Kobayashi, 2009)。ま た、UNDX-m は親の選び方によって子個体生成の分布にバイアスが生じることがある。 小林らは、このような問題へ対処するために、多親交叉である REXstar を REX に基づき 開発した(Kobayashi, 2009)。まず、REX は交叉に参加した親個体数をμ、選択された親 個体を (1 ≤ j ≤ μ) 、その重心を〈 〉 𝜇 ∑ 𝑗 としたとき、以下のように子個体を 生成する。 ci P i , j ( Pj P ) (1.7) j 1 ここで、ε は平均 0、分散σ2 の任意の対称な確率分布であり、σ2 =1/μを満たす。αは拡張 率を表し、子個体生成領域の広がりを決定するパラメータである。REX は統計量の遺 伝を満足する交叉となっている。REX は親個体の重心を交叉中心とし、子個体を親個 体の共分散行列を保存するように生成する。一方で、REXstar では個体集団や親個体群が 大域的最適解を覆っていない状況に対処するために、子個体生成分布の中心を親個体群 の重心から有望な方向へとシフトさせるために大域的降下方向という概念が取り入れ られた。大域的降下方向とは、親個体群の評価値のランキング上位半分の重心と、評価 値のランキング下位半分の重心を用いて、交叉中心をより有望な方へとシフトするよう な概念である。これにより、子個体を有望領域により多く生成しやすくなり、REXstar は稜構造と多峰性を有する関数において優れた性能を示した。秋本らは、REX や REXstar をもとに、悪スケール性や変数間依存性のある関数、更に稜構造関数や多峰性関数の最 17 第2章 ネットワーク推定手法 適化においても効果を発揮する多親交叉法の AREX を開発した(Akimoto, Nagata, et al., 2009)。秋本らは、初期収束の原因は集団が有望領域に移動する際に多様性を失いやす いことにあるとした。秋本らは、個体集団が有望領域にシフトする際に多様性を失うこ とが初期収束の原因であるとして、これらを改善するために拡張率適応(Adaptation of Expansion Rate; AER)(Akimoto, Sakuma, et al., 2009)と交叉中心降下(Crossover Mean Descent; CMD)を REX に組み込んだ。 初期収束を起こしやすい稜構造関数では比較的長く重心をシフトする必要があり、多 峰性関数ではいくつかの谷にまたがり個体が選択されるため重心のシフトが遅い。秋本 らは、これらの点から、実数値 GA は重心を移動している際に初期収束を起こしやすい ことを指摘した。拡張率を大きくして子個体を親個体よりも広い領域で生成した場合、 個体集団の分布を維持しながら拡大することができるが、選択圧がかからなくなるため 個体集団が発散し探索効率が悪くなる。AREX では、探索効率を保ったままで初期収束 を回避するために、拡張率を適応的に制御する手法が組み込まれている。この手法は CMA-ES(evolution strategy with covariance matrix adaptation) (Hansen, 2006)の step-size adaptation (Ostermeier, Gawelczyk, and Hansen, 1994)に着想を得ている。拡張率は以下 のように表現される。 α(g+1) : = max(α(g) √(1 - cα + cα Lcdp /Lavg , 1) (1.8) gは現在の世代番号、cα は学習率を表し、学習率の推奨値は l を次元数としたとき cα= (1/5)l である。Lcdp およびLavg はそれぞれ次のように表される。 𝜇 Lcdp = 𝛼 2 𝜇 − 1 𝑗=1 1 〈𝜀𝑗 〉2𝜇𝛼 − 𝜇 Lavg = 𝛼 2 𝜎 2 𝜇 − 1 18 2 2 𝜇 〈𝜀𝑗 〉𝜇𝛼 𝑗=1 𝜇𝛼 (1.9) 第2章 ネットワーク推定手法 ここで、添字𝑖: 𝜆を𝜆個の子個体の内評価値上位𝑖番目の個体の添字を表すとすると、 𝛼 〈𝜀𝑗 〉𝜇𝛼 = ∑𝜇𝑖=1 𝜀𝑖:𝜆,𝑗 𝜇𝛼 となる。式(1.8)より、適応的な拡張率の制御を実現している。こ れに加え、AREX には子個体生成領域の中心を有望な領域にシフトする仕組みである交 叉中心降下が組み込まれている。この仕組みにより、重心移動を加速し、初期収束を防 ぐとともに探索効率を向上させている。交叉中心降下は REXstar に着想を得ており、式 (1.7)を次のように変更している。 ci = m + α μ 𝜀 j = 1 i,j Pj - P (1.10) μ wj Pi:j .m = j=1 こ こ で 、 添 字 𝑖: 𝜇 は 親 𝜇 個 体 の 内 、 評 価 値 上 位 𝑖 番 目 の 個 体 の 添 字 を 表 し 、 ま た μ wj = 2(μ + 1 - j) / μ(μ + 1)は線形の重みで、∑j = 1 wj = 1である。AREX では、式(1.9)から式 (1.8)を求め、式(1.10)のαに代入する。そして式(1.10)に基づいて子個体を生成する。 これまで、モデルネットワークを用いて AREX+JGG 法の相互作用ネットワーク推定 における有効性を検証しており、従来の手法である UNDX+MGG 法や SPX+MGG 法と 比べて、推定可能なパラメータ数と推定速度の両方で AREX+JGG 法が優れていること を示している(Komori et al., 2012) 。しかしながら、推定の終期に収束速度の停滞が観 察され、推定速度の低下に繋がっていた。 第4項 ハイブリッド法 AREX+JGG 法の推定終期の収束速度の停滞を改善するために、AREX+JGG と修正 Powell 法(Powell, 1964)を組み合わせたハイブリッド法を開発した(Komori et al., 2013)。 19 第2章 ネットワーク推定手法 修正 Powell 法は、目的関数の導関数を必要とせず、m 変数の 2 次式の目的関数の極小 を m 回の繰り返し関数で求められる手法である。修正 Powell 法は局所探索において優 れた性能を示すが、局所最適解が多数存在するような問題においては早々に局所最適解 に陥ってしまう。更に、そこから抜け出す手段を持っていないため、大域的最適解へと 収束する保証はない。そこで、大域的探索において効果的な GA と修正 Powell 法を組 み合わせることで、両者の問題点を補い優れた探索能力を持つ最適化手法となることが 期待される。このことから、岡本らは、SGA と修正 Powell 法を組み合わせた手法を開 発し、SGA 単独での最適化よりも良好な結果を示したことを報告した(Okamoto et al., 1998)。AREX+JGG 法においても、推定後期の大域的最適解の近傍にいると思われる個 体に対し修正 Powell 法を適用することで、大域的最適解へ迅速に収束することが期待 される。 修正 Powell 法のアルゴリズムは次のとおりである。ここで、m 個の未知パラメータ を𝒖 = [𝒖1 , 𝒖2 , … , 𝒖𝑚 ]とし、目的関数𝐹 = 𝐹 𝒖 を極小とする𝒖∗ を求めるとする。𝝃1 , … , 𝝃𝑚 は一次独立な m 次元ベクトルであるとし、座標軸の方向に用いる。 1.出発点を𝒖0 とし、𝑘 = 1, … , 𝑚において、𝐹 𝒖𝑘−1 + 𝛼 𝑘 𝝃𝑘 を最小とする𝛼 𝑘 を計算す る。これにより点𝒖𝑘 = 𝐹 𝒖𝑘−1 + 𝛼 𝑘 𝝃𝑘 を求める。 2.{𝐹(𝒖𝑙−1 ) − 𝐹 𝒖𝑙 }が最大となる整数𝑙 1 ≤ 𝑙 ≤ 𝑚 を求め、Δ = 𝐹(𝒖𝑙−1 ) − 𝐹 𝒖𝑙 とす る。 3.𝐹3 = 𝐹 2𝒖𝑚 − 𝒖0 の値を求め、𝐹1 = 𝐹 𝒖0 、𝐹2 = 𝐹 𝒖𝑚 とする。 4.𝐹1 ≤ 𝐹3 および 1 2 Δ 𝐹1 − 𝐹3 2 ≤ 𝐹1 − 2𝐹2 + 𝐹3 𝐹1 − 𝐹2 + Δ 2 であれば、次の繰 り返し解では古い方向𝝃1 , … , 𝝃𝑚 を用い、次の出発点𝒖0 として𝒖𝑚 を用いる。この式 のいずれも成立しないときには、5 の操作を行う。 5 . 𝝃 = 𝒖𝑚 − 𝒖0 と お き 、 𝐹 𝒖𝑚 + 𝛼𝝃𝑘 を 最 小 と す る 𝛼 を 求 め る 。 そ し て 、 𝝃1 , … , 𝝃𝑙−1 , 𝝃𝑙+1 , … , 𝝃𝑚 , 𝝃を探索の方向として、𝒖𝑚 + 𝛼𝝃𝑘 を次の繰り返し回の出発点 とする。 20 第2章 ネットワーク推定手法 修正 Powell 法において、𝛼を求める際には Rosenbrock 法(Rosenbrock, 1960)を用い た。Rosenbrock 法は、探索方向のベクトルに沿って、極小点を順次絞り込んでいく方法 である。Rosenbrock 法のアルゴリズムは以下のとおりである。 1.探索開始点を𝒖、探索方向を𝝃とし、𝛼と刻み幅ℎの初期値を決定する。 2.目的関数𝐹が𝐹 𝒖 + 𝛼 + ℎ 𝝃 < 𝐹 𝒖 + 𝛼𝝃 を満たせば𝛼 + ℎを新しい𝛼とし、刻み幅 ℎを𝑥倍にする(𝑥 > 1.0)。𝐹 𝒖 + 𝛼 + ℎ 𝝃 < 𝐹 𝒖 + 𝛼𝝃 を満たさなかった場合は、 刻み幅ℎを−𝑦倍にする(0 < 𝑦 < 1) 。 𝛼の最大値と最小値の差が設定した終了条件より小さくなるまで 2 を繰り返す。 これまで、実数値 GA の初期個体すべてに対し修正 Powell 法を適用するようなアル ゴリズムや各世代のエリート個体に対して修正 Powell 法を適用するハイブリッド法が 開発され、その有効性が示されている(Okamoto et al., 1998; Yoshimura et al., 2003)。本 研究で開発したハイブリッド法は、AREX+JGG 法の探索終期の収束停滞の改善を目的 としているため、推定終期に評価値が設定値に達したエリート個体にのみ修正 Powell 法を適用することとした。ハイブリッド法のアルゴリズムは以下の 1~3 のとおりであ る。 1.大域的な探索を AREX+JGG 法を用いて行う。 2.個体集団で最も評価値の良い個体(エリート個体)が定めた評価値となるまで 1 を繰り返す。 3.定めた評価値となったエリート個体に対し修正 Powell 法を適用する。終了条件を 満たした場合、推定を終了する。 21 第2章 第5節 第1項 ネットワーク推定手法 検証 AREX+JGG 法の有効性の検証 AREX+JGG 法の有効性を検証するために、6 要素ネットワークモデルを人工的に作成 した(Fig. 2.4) 。このネットワークモデルにおいて、相互作用は各要素の生成過程に影 響し、各要素は一次崩壊すると仮定した。具体的には、S-system モデル(式(1.2))中の gij のパラメータ値のみ区間[-1.0, 1.0]の値を設定し、hij はi = jのとき 1.0、i ≠ jのとき 0.0 と設定した。また、𝛼i および𝛽𝑖 は 1.0 に設定した。 1 0.7 0.4 2 3 0.3 -0.2 5 4 0.5 0.1 6 Fig. 2.4 6 要素ネットワークモデル →は正の相互作用を表し、―|は負の相互作用を表す 数字は S-system(式(1.2))の gij のパラメータ値を表す 22 第2章 Table 2.1 case1-01 ネットワーク推定手法 各ケーススタディーにおける最適化対象パラメータ 最適化対象 最適化対象とした パラメータ数 S-system の未知パラメータ 20 gij (i = 2, 4, 5, 6, j = 1, 2, … , 6; when i ≠ j) gij (i=1, 2, … , 6, j = 1, 2, …, 6), case1-02 54 hii (i=1, 2, … , 6), αi (i=1, 2, … , 6), βi (i=1, 2, … , 6) gij (i=1, 2, … , 6, j = 1, 2, …, 6), case1-03 84 hij (i=1, 2, … , 6, j = 1, 2, …, 6), αi (i=1, 2, … , 6), βi (i=1, 2, … , 6) 作成した各ネットワークモデルに含まれる S-system のシステムパラメータの内、最 適化対象とするパラメータを変えることで、未知パラメータ数の異なる 3 つのケースス タディーを用意し、従来法である UNDX+MGG 法および新規法である AREX+JGG 法を 用いて推定した(case1-01~case1-03) 。それぞれのケーススタディーに含まれる最適化 対象パラメータは Table 2.1 に示した。 case1-01 の推定における最適化条件を Table 2.2 に示した。これらのパラメータ値は先 行研究における推奨値を用いた。なお、終了条件は平均相対二乗誤差が 1.0%を下回っ たときとした。また、世代数が最大世代数を超えた場合は推定に失敗したものとみなし た。推定に用いたタイムコースデータは、全要素が正常に発現している場合と、各要素 を 1 要素ずつ破壊した場合を想定し用意した。 UNDX+MGG 法および AREX+JGG 法についてそれぞれ 50 試行ずつ推定を行った結 果、両手法ともにすべての試行において最大世代数に達することなく推定を終了するこ とができた。50 試行の世代毎の平均相対二乗誤差の平均値を Fig. 2.5 に示す。Fig. 2.5 から明らかなように、UNDX+MGG 法よりも AREX+JGG 法の方が遥かに少ない世代数 で推定することができている。推定に要した計算時間の平均は、CPU time(CPU: Xeon E5540 2.53GHz, Memory: 8.0GB)で、UNDX+MGG 法は 765.0 秒であったのに対し、 AREX+JGG 法は約 20 倍の速さの 38.6 秒であった。 23 第2章 Table 2.2 ネットワーク推定手法 case1-01 における最適化条件 (a) UNDX+MGG 法 項目 設定 要素数 6 サンプリング数 71 データセット数 7 UNDX パラメータ 主軸方向 0.5 副軸方向 0.3 集団サイズ 300 子個体数 100 親個体数 2 cut off 0.05 最大世代数 100,000 (b) AREX+JGG 法 項目 設定 要素数 6 サンプリング数 71 データセット数 7 集団サイズ 200 子個体数 80 親個体数 21 cut off 0.05 最大世代数 100,000 24 第2章 ネットワーク推定手法 Fig. 2.5 case1-01 における平均相対二乗誤差と世代数の関係 エラーバーは標準偏差を示す 25 第2章 ネットワーク推定手法 case1-02 および case1-03 について、推定を行った結果を Table 2.3 に示す。Table 2.1 に 示したように、case1-02 および case1-03 の最適化対象パラメータ数は、それぞれ 54 個 と 84 個である。また、いずれのケーススタディーにおいても、目標平均相対二乗誤差 が 5.0%を下回った時に推定に成功したこととし、世代数が 100,000 世代に達した場合推 定に失敗したとみなした。なお、case1-02 および case1-03 では、実数値 GA における集 団サイズの条件を変えて検証した。実数値 GA では、関数の景観により必要な集団サイ ズ npop は大きく異なる。UNDX+MGG 法では、これまで用いていたnpop = 300に加え、 次元数(最適化対象パラメータ数)を 𝑙 としたとき、npop = 20l, 30l, 40lとなるように設 定した。AREX+JGG 法では、npop = 20l, 30l, 40l, 50l となるように設定した。それぞれ の集団サイズの設定値で 20 回ずつ試行し、その平均値を Table 2.1 に示している。 UNDX+MGG 法では、case1-02、case1-03 において、いずれの集団サイズの設定値でも 最大世代数 100,000 世代に達するまでに目標平均相対二乗誤差を下回ることができなか った。一方 AREX+JGG 法では、npop = 20lのときを除き、100,000 世代以内に推定を終 了することができた。また、npop = 30lのときに推定に費やす世代数および計算時間が 最小となり効率的に推定できたことが分かった。 26 27 30l 40l 50l 2,520 3,360 4,200 (case1-03) 50l 2,700 84 40l 2,160 (case1-02) 20l 30l 1,620 54 1,380 20l 1,080 30l 40l 2,520 (case1-03) 3,360 1,380 84 20l 2,160 300 40l 1,620 (case1-02) *CPU : Xeon E5540 2.53GHz, Memory : 8.0GB AREX+JGG UNDX+MGG 30l 1,080 336 336 336 336 216 216 216 216 336 336 336 100 216 216 216 100 300 20l 子個体数 集団サイズ 85 85 85 85 55 55 55 55 2 2 2 2 2 2 2 2 親個体数 5.0 5.0 5.0 7.0 5.0 5.0 5.0 6.7 57.3 45.8 33.7 24.6 40.0 13.4 8.9 8.3 (%) 平均相対二乗誤差 case1-02 および case1-03 の結果 実数値 GA のパラメータ 54 パラメータ数(l) 推定対象 Table 2.3 4,682 3,931 3,335 >100,000 3,871 3,759 2,374 >100,000 >100,000 >100,000 >100,000 >100,000 >100,000 >100,000 >100,000 >100,000 平均世代数 推定終了時の 4214.7 3511.8 2801.6 139688.0 2228.2 2177.3 1367.2 59190.8 114232.2 105570.0 113824.0 27033.0 53628.4 74510.4 57680.5 27885.2 [sec]*/試行 平均 CPU time 第2章 ネットワーク推定手法 第2章 第2項 ネットワーク推定手法 ハイブリッド法の有効性の検証 ハイブリッド法の有効性を検証するために、Fig. 2.4 の 6 要素ネットワークモデルの 推定を行った。最適化対象パラメータは case1-01 と同様のパラメータとし、未知パラメ ータ数は 20 個である。この検証における最適化条件を Table 2.4 に示す。同じ条件で 50 試行推定し、各世代におけるエリート個体の平均相対二乗誤差の平均値を Fig. 2.6 に示 した。また、case1-01 の AREX+JGG 法の結果についても同グラフに示す。 Table 2.4 ハイブリッド法の最適化条件 項目 設定 要素数 6 タイムコースデータ数 7 タイムコースデータのサンプリング数 71 ハイブリッド法 Rosenbrock 法 𝑥* 1.618034 𝑦* 0.618034 刻み幅 h* 0.1 終了条件 0.01 修正 Powell 法開始平均相対二乗誤差(%) < 15.0 集団サイズ 200 子個体数 80 親個体数 21 cut off 0.05 最大世代数 100,000 目標平均相対二乗誤差(%) < 1.0 *第 2 章第 4 節第 4 項参照 28 第2章 ネットワーク推定手法 (a) 600 AREX+JGG法 平均相対二乗誤差(%) 500 ハイブリッド法 400 300 200 100 0 0 50 100 150 世代数 200 250 (b) 35 AREX+JGG法 平均相対二乗誤差(%) 30 ハイブリッド法 25 20 15 10 5 0 100 150 200 世代数 Fig. 2.6 AREX+JGG 法とハイブリッド法の推定結果 (a) 6 要素ネットワークモデル推 定における世代数と平均相対二乗誤差 (b) 推定終期の世代数と平均相対二乗誤差 エラーバーは標準偏差を表す 29 第2章 ネットワーク推定手法 推定の結果、AREX+JGG 法よりもハイブリッド法の方が、早い世代で推定すること ができた。推定に要した CPU time では、AREX+JGG 法が 37.0 秒であったのに対し、ハ イブリッド法は 31.2 秒であった。Fig. 2.6(b)からわかるように、ハイブリッド法では AREX+JGG 法で観察されるような推定終期における収束の停滞がほとんどない。この 検証で、AREX+JGG 法に修正 Powell 法を導入したハイブリッド法により、AREX+JGG 法の推定終期の収束停滞を改善できたことを示すことができた。 第3項 新規法の推定可能パラメータ数の検証 開発した AREX+JGG 法とハイブリッド法の推定可能なパラメータ数を検証するため に、Fig. 2.7 のような 17 要素ネットワークモデルを作成し、AREX+JGG 法とハイブリ ッド法でこのネットワークモデルの推定を試みた(case1-04)。この推定における最適化 条件を Table 2.5 に示した。 このとき最適化対象となるパラメータ数は 272 個であった。 推定の結果、どちらの手法でも最大世代数に達することなく推定に成功した。推定に要 した時間は、CPU time で AREX+JGG 法が 101.6 時間、ハイブリッド法が 93.9 時間であ り、ハイブリッド法の方が早い時間で推定することに成功した。 また、17 要素ネットワークモデル(Fig. 2.7)において最適化対象となるパラメータ を αi、βi、gij (i ≠ j)としたとき、AREX+JGG 法とハイブリッド法どちらの手法を用いた 場合でも最大世代数内に目標平均相対誤差に達することができなかった。このときの最 適化対象となるパラメータ数は 306 個であった。 30 第2章 1 1.0 1.0 2 1.0 3 -0.1 4 1.0 -0.1 5 1.0 6 0.7 -0.2 8 7 9 0.3 0.4 -0.1 11 0.3 0.6 12 0.5 14 -0.2 15 0.3 10 -0.2 13 0.4 0.1 16 0.6 17 Fig. 2.7 17 要素ネットワークモデル 31 ネットワーク推定手法 第2章 ネットワーク推定手法 case1-04 におけるハイブリッド法の最適化条件 Table 2.5 項目 設定 最適化対象パラメータ gij (i ≠ j) 要素数 17 タイムコースデータ数 18 タイムコースデータのサンプリングポイント数 71 Rosenbrock 法 ハイブリッド法 𝑥* 1.618034 𝑦* 0.618034 刻み幅 h* 0.1 終了条件 0.01 修正 Powell 法 開始平均相対二乗誤差(%) < 25.0 集団サイズ 8160 子個体数 1088 親個体数 273 cut off 0.05 最大世代数 100,000 目標平均相対二乗誤差(%) < 10.0 *第 2 章第 4 節第 4 項参照 32 第2章 第4項 ネットワーク推定手法 タイムコースデータ数の検証 case1-01~case1-04 における検証では、ネットワーク構成要素数を n とすると(n+1) 個のタイムコースデータが推定に用いられた。実験データを用いる場合、case1-01~ case1-04 で想定したネットワーク中の全要素をそれぞれ破壊した実験のような、十分な 数のデータが得られるとは限らない。そこで、タイムコースデータ数が少ない場合にお いても AREX+JGG 法とハイブリッド法が有効であるか検証した(case1-05)。 Fig. 2.4 の 6 要素ネットワークモデルを用い、最適化対象とした S-system に含まれる パラメータは case1-01 と同様の条件で 20 パラメータとした(Table 2.1)。この検証にお ける最適化の条件を Table 2.6 に示した。本検証では、要素 1~6 の 1 要素破壊実験を想 定したタイムコースデータの内、要素番号が大きなデータから削減し、タイムコースデ ータ数を 7~1 個と変化させた。それぞれの場合について 30 試行ずつ推定し、目標相対 二乗誤差となった解候補すべてに共通する相互作用(共通構造)を抽出した。 検証の結果を Table 2.7 に示した。AREX+JGG 法およびハイブリッド法どちらの手法 においても、タイムコースデータ数が 7、6 個のとき、6 要素ネットワークモデルに含 まれるすべての相互作用を推定することに成功した。しかしながら、タイムコースデー タ数が減る毎に推定に成功した相互作用数も減少し、データ数が 1 個のとき相互作用を 1 つも推定することができなかった。タイムコースデータの数に対する推定に成功した 相互作用数は、AREX+JGG 法およびハイブリッド法いずれの手法でも差はなかった。 また、いずれのタイムコースデータ数においても間違った相互作用は推定されなかった。 本検証から、AREX+JGG 法およびハイブリッド法は、タイムコースデータ数が少な い場合においても相互作用を推定することができるが、すべての相互作用が推定できな くなり効率が悪くなることがわかった。また、あまりにデータ数が少ないと、相互作用 を推定できなくなる可能性が示唆された。 33 第2章 Table 2.6 ネットワーク推定手法 case1-05 における最適化条件 項目 設定 要素数 6 タイムコースデータのサンプリングポイント数 71 Rosenbrock 法 ハイブリッド法 𝑥* 1.618034 𝑦* 0.618034 刻み幅 h* 0.1 終了条件 0.01 修正 Powell 法 開始平均相対二乗誤差(%) < 15.0 集団サイズ 200 子個体数 80 親個体数 21 cut off 0.05 最大世代数 10,000 目標平均相対二乗誤差(%) < 5.0 *第 2 章第 4 節第 4 項参照 34 35 3 2 0 2 1 5 5 3 6 6 4 6 7 4 0 2 2 1 3 3 5 5 4 6 6 4 6 7 *CPU : Xeon E5540 2.53GHz, Memory : 8.0GB ハイブリッド法 AREX+JGG 法 推定に成功した相互作用数(/6) タイムコースデータの数 Table 2.7 case1-05 の推定結果 30 28 30 11 23 30 30 30 29 29 10 25 25 30 成功試行回数(/30) 4.02 59.71 31.02 752.05 366.60 45.20 36.22 3.25 44.48 64.02 952.95 369.64 412.51 50.19 / 試行 平均 CPU time [sec]* 第2章 ネットワーク推定手法 第2章 第6節 ネットワーク推定手法 考察 6 要素ネットワークモデルを用いた case1-01~case1-03 の結果から、AREX+JGG 法が UNDX+MGG 法よりも推定時間および推定可能パラメータ数の両方で優れた手法であ ることが示された(Fig. 2.5、Table 2.3)。AREX+JGG 法は UNDX+MGG 法に比べ、推定 の早い段階で評価値である平均相対二乗誤差が低下していた(Fig. 2.5)。これは、AREX に組み込まれている効率的に有望領域へと交叉中心を移動させる枠組みや、JGG の世代 交代率の高さ(本章第 4 節第 3 項)に起因するものと思われる。また、case1-02 および case1-03 において、最適化対象パラメータ数を l とすると、AREX+JGG 法の集団サイ ズnpop = 30lのときに特に優れた性能を示した。更に、AREX+JGG 法と修正 Powell 法を 組み合わせたハイブリッド法を開発することで、AREX+JGG 法の推定終期の収束停滞 を改善することができた(Fig. 2.6) 。 case1-04 の検証で、これまで、従来法で 3 要素程度の相互作用ネットワークを推定す ることができなかった S-system と実数値 GA の組み合わせ手法を、AREX+JGG 法およ びハイブリッド法により、17 要素程度から成る相互作用ネットワークを推定できる手 法とすることに成功した。しかしながら、生命システムにおける相互作用ネットワーク に適用するためには、より大きな相互作用ネットワークを推定できる必要があり、新規 法に対して更なる工夫を加えることが求められる。 また、case1-05 では、推定に用いるタイムコースデータ数に対する推定に成功した相 互作用数の検証を行った。その結果、タイムコースデータ数が 7、6 個のときはすべて の相互作用を推定することができたが、データ数が減少するごとに推定できた相互作用 の数も減少した。タイムコースデータ数が少ないときは拘束条件が緩やかとなるため、 データを再現できるようなパラメータの組合せが多数存在し、最適解を絞ることができ ず推定できた相互作用数が減少したと考えられる。AREX+JGG 法やハイブリッド法を 用い効率的な推定を行うためには、一定数のタイムコースデータが求められることが示 36 第2章 ネットワーク推定手法 唆された。本検証では、1 要素破壊実験を想定したタイムコースデータを用いて検証が 行われたが、今後、ある要素を過剰発現した場合やパルス刺激を与えた場合、欠損値の ある場合などのタイムコースデータに対する AREX+JGG 法およびハイブリッド法の有 効性の検証が望まれる。 37 第3章 第3章 Random selection method 第1節 諸言 Random selection method AREX+JGG 法およびハイブリット法は従来の S-system と実数値 GA に基づく相互作 用ネットワーク推定手法に比べ、推定規模と推定速度の両方において優れた手法である。 しかしながら、提案法単独では推定可能な要素数は未だ 17 要素程度に留まっている。 そこで、本研究では AREX+JGG 法およびハイブリッド法を用い推定可能規模を増大さ せるためのアルゴリズムを考案した。 最適化手法は変えず大規模な相互作用ネットワークを推定するための比較的容易な 方法として、k-means 法(MacQueen, 1967)などによるクラスタリングの導入が挙げら れる。クラスタリングにより推定対象とする生体分子をいくつかのクラスタに分け、小 規模なネットワークをいくつも推定することで最終的な相互作用ネットワークを決定 する戦略が考えられ得る。しかしながら、相互作用ネットワークを推定する際には、そ の詳細が不明な場合が多く、適切にクラスタリングするのは困難である。そこで、本研 究では対象となるすべての要素からランダムに要素を選択することでサブネットワー クを作成し、そのサブネットワークを複数推定することで全体的な相互作用ネットワー クを推定することを試みた(Random selection method; RA 法) (Komori et al., 2013) 。 38 第3章 第2節 Random selection method Random selection method のアルゴリズム RA 法のアルゴリズムは次の 1~8 の通りである。 1.m 要素のネットワークから k 要素をランダムに選択し、k 要素から成るネットワー クをサブネットワークとする(k < 𝑚)。このとき、選択される要素の頻度を等し く保つために selection rate(SR)を設定する。 2.(k+1)種類のタイムコースデータを用意する。 3.サブネットワークの相互作用をハイブリッド法を用い l 回推定する。このとき設 定した平均相対二乗誤差以下となった個体が示すネットワークをネットワーク候 補群とする。 4.推定されたすべてのネットワーク候補群に共通する相互作用(共通構造)を抽出 する。 5.共通構造を保存する。 6.1~5 の操作を p 回繰り返す。 7.すべての共通構造を組み合わせる。 8.最終的な相互作用ネットワークを決定する。 アルゴリズムの概念図を Fig. 3.1 に示した。 Random selection method では、サブネットワークの要素を選択する際に各要素が均等 に選ばれるように SR が設けられている。SR は次のような式で表される。 SR = p × k m (2.1) 例えば、k=2、m=30、p=30 としたとき、SR=30×2/30 となり、各要素は必ず 2 回 39 第3章 Random selection method ずつ選ばれることとなる。なお、本研究においては、Mersenne Twister(Matsumoto and Nishimura, 1998)を用い乱数を生成することで、ランダムな要素選択を行った。また、 (k+1)種類のタイムコースデータは、全要素が正常に発現している場合と各要素の 1 要 素破壊実験を想定した場合のタイムコースデータを用意した。第 2 章第 5 節第 4 項で示 したように、タイムコースデータ数が少ないとき、推定される相互作用数も少なくなり 推定効率が低くなると考えられる。 40 第3章 29 4 12 ? 16 Random selection method (1) m 要素のネットワークから k 要素をランダムに選択する (2) (k+1)種類のタイムコースデータを与える 11 タイムコースデータ 2 1 3 4 16 5 6 7 8 11 9 (3) ハイブリッド法を用い サブネットワークを l 回推定する 24 14 29 4 12 12 16 17 19 20 21 24 25 26 27 28 29 16 23 23 12 11 13 22 18 29 4 23 12 11 15 18 29 4 11 12 27 29 4 10 11 3 2 7 11 30 (4) 共通構造を抽出する (7) すべての共通 構造を組み合わせる (9) 最終的な相互作用 ネットワークを決定する 16 (5) 共通構造を 保存する 11 (6) 1~5 の操作を p 回繰り返す Fig. 3.1 Random selection method の概略図 41 23 第3章 第3節 Random selection method 評価方法 Random selection method を用い推定された相互作用の種類を true-positive interaction (TP interaction)、false-positive interaction(FP interaction) 、false-negative interaction(FN interaction)の 3 種で定義付けることとした。まず、TP interaction は、random selection method の共通構造から推定した相互作用の内、真のネットワークに存在する相互作用 である。次に、FP interaction は、共通構造から推定した相互作用ではあるものの、真の ネットワークに存在しない相互作用である。そして、FN interaction は、真のネットワー クには存在するが、共通構造から推定できなかった相互作用である。TP interaction と FP interaction は、実問題においては真のネットワーク構造が不明であるため区別できな い点は留意されるべきである。本研究では、ネットワークモデルを用いているため、真 のネットワークに存在するか否か区別ができるものとして定義した。そこで、推定率を 量る指標として、TP 率を用いた。TP 率は次の式で表される。 TP 率 = TPall TPall +FNall n n TPall = ∑ TP , FNall = ∑ FN i=1 (2.2) i=1 ここで、TPi は総数 n の推定された TP interaction の内、i 番目の TP interaction を示し、 FNi は総数 n の推定された FN interaction の内、i 番目の FN interaction を示す。実質上、 TPall+FNall は真のネットワークに存在する相互作用の総数を示すため、TP 率は真のネ ットワークの相互作用の内どの程度推定できているかを表す。真のネットワークに含ま れるすべての相互作用を推定できた場合、TP 率は 1.0 となる。 42 第3章 第4節 第1項 Random selection method Random selection method の有効性の検証 ネットワークモデル RA 法の有効性を検証するために、人工的に作成した 30 要素から成るネットワークモ デルを用いた(Fig. 3.2) 。これまで、AREX+JGG 法およびハイブリッド法単独ではこの ネットワークモデルの推定に成功していない。このネットワークモデルにおいて、相互 作用は各要素の生成過程に影響し、各要素は一次崩壊すると仮定した。具体的には、 S-system モデル(式(1.2))中の gij のパラメータ値のみ区間[-1.0, 1.0]の値を設定し、hij はi = jのとき 1.0、i ≠ jのとき 0.0 と設定した。また、𝛼i および𝛽𝑖 は 1.0 に設定した。こ のとき設定したパラメータ値を Table 3.1 に示した。このネットワークに含まれる相互 作用(gij)の数は 68 である。なお、i ≠ jのときに限ると 38 となる。 Table 3.1 のシステムパラメータを基に、各要素のタイムコースデータを用意した。 S-system の各状態変数の挙動が、各要素のタイムコースとなる。なお、各要素の 1 要素 破壊実験を想定し、破壊対象となる要素を i とすると、状態変数 xi の生成過程の速度定 数 αi を 0.0 とすることでタイムコースを得た。ここで、ネットワークに含まれる要素数 を n とすると、用意したタイムコースデータは(n+1)個となる。 43 第3章 2 3 0.5 0.4 1 1.0 -0.1 1.0 5 4 6 0.3 -0.2 10 11 8 1.0 0.7 -0.2 -0.1 24 0.4 0.3 0.6 0.5 22 21 0.3 0.5 0.4 0.6 20 19 0.3 -0.2 12 16 15 1.0 18 13 0.5 0.2 14 0.6 -0.2 -0.1 9 0.2 0.4 7 1.0 Random selection method 17 0.1 0.2 23 -0.2 26 25 0.4 0.5 0.1 28 27 29 0.6 -0.2 30 Fig. 3.2 30 要素ネットワークモデル →は正の相互作用を表し、―|は負の相互作用を表す 44 第3章 Table 3.1 Random selection method 30 要素ネットワークモデルのシステムパラメータ αi 1.0 βi 1.0 gij g5,1 = 1.0, g6,1 = 1.0, g7,2 = 0.5, g7,3 = 0.4, g8,4 = 0.2, g11,4 = 0.4, g9,5 = 1.0, g9,6 = –0.1, g10,7 = 0.3, g11,7 = –0.2, g13,8 = 0.6, g14,9 = 1.0, g15,10 = 0.2, g16,11 = 0.5, g16,12 = –0.2, g17,13 = 0.5, g1,14 = –0.1, g19,14 = 0.1, g20,15 = 0.7, g24,15 = –0.2, g21,16 = 0.6, g22,16 = 0.5, g8,17 = –0.2, g23,17 = 0.2, g24,18 = –0.1,g24,19 = 0.3, g25,20 = 0.4, g26,21=–0.2, g12,22 = 0.4, g13,23 = 0.1, g27,24 = 0.6, g27,25 = 0.3, g28,25 = 0.5, g20,26 = 0.3, g29,26 = 0.4, g30,27 = 0.6, g26,28 = 0.1,g27,30 = –0.2, other gij = 0.0 hij 1.0(i j), 0.0( i j) 45 第3章 第2項 Random selection method 30 要素ネットワークモデルの推定 まず、30 要素ネットワークモデルを推定するために、RA 法のパラメータ k の適切な 値を検討した。このとき設定した他のパラメータについては Table 3.2 に示した。30 要 素ネットワークモデルを用いて検証した結果、k を 5 としたときに、推定できた相互作 用の数がもっとも多いことがわかった(Table 3.3) 。 そこで、k = 5とし 30 要素ネットワークモデルの推定を試みた。このときの RA 法のパ ラメータを Table 3.4 に示した。SR については、SR = 2とするランダムな選択を 8 回繰 り返すこととした。つまり、各要素が必ず 2 回ずつ選ばれるというランダムな選択を 8 回繰り返すことで、ランダムな選択の総計である p が 96 となるようにした。この推定 における RA 法の概略図を Fig. 3.3 に示した。これらの条件で、30 要素からランダムに 要素を選択した際に、各試行(p)で選ばれた要素番号を Table 3.5 に示した。また、RA 法中のハイブリッド法で用いたパラメータは過去の検証において良好な性能を表した 値を用いた(Table 3.6) (Akimoto, Nagata, et al., 2009; Akimoto, Sakuma, et al., 2009; Komori et al., 2012) 。 RA 法での推定結果を Table 3.7 に示す。推定の結果、最終的に 30 要素ネットワーク モデル中の約 93%の相互作用を推定することができた。 46 第3章 Table 3.2 Random selection method k の最適化のために用いたパラメータ 設定値 l 30 m 30 p 10 Table 3.3 k の最適化の結果 要素数(k) 推定した TP interaction の数(/68) CPU time [h]* 4 22 54.1 5 35 132.4 6 31 327.3 *CPU : Xeon E5540 2.53GHz, Memory : 8.0GB 47 第3章 Random selection method 30 要素ネットワーク推定のために用いたパラメータ Table 3.4 設定値 17 29 9 18 30 7 l 30 m 30 p 96 SR 2×8 30 試行(l)ハイブリッド法で サブネットワークを推定 Random selection 2 6 20 29 4 4 30 要素(m) Fig. 3.3 5 8 10 19 12 16 3 k 60 回 (p) ? 11 21 12 5 要素(k) Random selection method による 30 要素ネットワーク推定の概略図 48 第3章 Table 3.5 Random selection method 各試行で選ばれた要素番号 p 選択された要素番号 p 選択された要素番号 p 選択された要素番号 1 1 8 9 15 17 33 8 14 19 22 25 65 3 9 11 22 25 2 2 5 6 16 22 34 3 5 6 14 30 66 1 8 12 14 26 3 5 7 20 23 28 35 5 13 19 22 23 67 6 17 28 29 30 4 6 17 24 26 27 36 1 13 18 20 28 68 5 10 15 20 27 5 9 14 18 22 27 37 1 8 9 21 23 69 15 18 19 21 25 6 18 20 24 25 29 38 3 9 14 28 30 70 4 8 19 20 30 7 12 14 19 26 29 39 1 5 17 20 29 71 7 21 23 24 27 8 1 8 21 23 25 40 7 20 23 28 29 72 4 7 10 18 24 9 4 7 11 13 15 41 4 10 12 13 15 73 7 14 23 24 25 10 10 12 19 28 30 42 8 12 18 25 26 74 4 8 12 15 16 11 2 3 10 21 30 43 13 19 25 26 27 75 12 16 17 19 22 12 3 4 11 13 16 44 5 6 17 18 22 76 5 6 9 14 26 13 2 13 24 27 28 45 14 15 16 19 24 77 5 11 13 23 28 14 7 11 19 29 30 46 2 3 6 10 27 78 1 2 3 13 30 15 8 13 14 18 25 47 2 11 21 22 24 79 6 15 24 25 27 16 4 6 8 21 28 48 4 7 11 16 30 80 2 8 9 18 22 17 3 21 22 25 26 49 5 9 15 18 25 81 19 20 21 29 30 18 7 9 17 27 29 50 1 11 12 17 25 82 1 3 7 17 29 19 12 17 18 22 24 51 3 5 8 12 21 83 10 18 21 26 28 20 4 10 14 20 26 52 7 9 11 27 30 84 4 10 11 20 27 21 1 9 11 15 16 53 2 18 19 28 30 85 6 7 18 19 22 22 3 5 6 20 23 54 4 6 17 20 27 86 6 12 16 21 22 23 5 10 19 23 30 55 1 15 22 23 28 87 2 7 18 21 27 24 1 2 12 15 16 56 7 8 10 13 23 88 3 12 17 20 23 25 9 12 15 20 29 57 4 6 14 22 29 89 2 3 19 20 26 26 4 10 12 17 27 58 13 16 19 26 29 90 1 4 10 16 24 27 4 7 11 24 27 59 3 20 21 24 26 91 1 5 23 24 30 28 2 11 16 28 29 60 2 10 14 16 24 92 4 10 25 27 28 29 2 7 15 21 24 61 1 2 11 17 29 93 9 11 13 15 26 30 1 3 16 21 25 62 3 5 6 9 13 94 5 8 11 15 28 31 8 9 17 18 26 63 13 14 16 22 23 95 8 9 14 29 30 32 6 10 23 26 30 64 2 12 16 26 28 96 13 14 17 25 29 49 第3章 Table 3.6 Random selection method RA 法による 30 要素ネットワークモデルの最適化条件 項目 設定 αi、βi、gij、hij (i 最適化対象パラメータ 要素数 5 タイムコースデータ数 6 タイムコースデータのサンプリングポイント数 71 Rosenbrock 法 ハイブリッド法 𝑥* 1.618034 𝑦* 0.618034 刻み幅 h* 0.1 終了条件 0.01 修正 Powell 法 開始平均相対二乗誤差(%) < 20.0 集団サイズ 200 子個体数 80 親個体数 21 cut off 0.05 最大世代数 100,000 目標平均相対二乗誤差(%) < 10.0 *第 2 章第 4 節第 4 項参照 50 j) 第3章 Table 3.7 Random selection method 30 要素ネットワークモデル推定結果 推定した相互作用の数(/68) SR p TP FP FN interaction interaction interaction TP 率(%) Total CPU time [h]* 2 12 39 2 29 57.35 358.22 4 24 45 4 23 66.18 741.86 6 36 49 5 19 72.06 1119.14 8 48 53 8 15 77.94 1485.77 10 60 57 11 11 83.82 1883.61 12 72 58 12 10 85.29 2283.77 14 84 60 13 8 88.24 2597.87 16 96 63 14 5 92.65 2967.20 *CPU : Xeon E5540 2.53GHz, Memory : 8.0GB 51 第3章 第5節 Random selection method 考察 RA 法では、もとのネットワークを構成するすべての要素からランダムに要素を選び サブネットワークを形成し、サブネットワーク内の相互作用を推定することで全体的な ネットワークを再現しようとする手法である。本章では、30 要素ネットワークモデル を用い RA 法の相互作用推定における有効性を検証した。 まず、30 要素ネットワークモデルの推定で用いる RA 法のパラメータ k を設定するた めに、k を何点かに振り分け推定を行った。その結果、k=5 のときもっとも多くの TP interaction を推定することができた。また、k=5 と設定した場合よりも k=6 と設定し たときの方が、導出された TP interaction は少なく推定に要した時間も多かった。この ことから、k を大きくすることで必ずしもより多くの TP interaction を推定できるとは限 らないことがわかった。推定対象とするネットワークに応じて最適な k を見つけること が、効率的な推定に繋がることが示唆された。 k=5 と設定し、30 要素ネットワークモデルを推定したところ、SR=16 の段階で約 93% の相互作用を推定することができた。このときの SR と TP 率の関係を Fig. 3.4 に示した。 TP 率が 80%を超えるまで(SR=10)は、SR の増加に対し TP 率は線形に増加している。 しかしながら、TP 率の増加は徐々に飽和状態になっていることから、ある程度の TP 率まで達すると、TP interaction を推定するために増やさなければならない SR 値が徐々 に大きくなり、要する推定時間も増大することが推測される。 本章の結果から、タイムコースデータから相互作用ネットワークを推定する手法とし て、RA 法は有効であることが示された。RA 法を用いることで、AREX+JGG 法やハイ ブリッド法単独では推定することのできない 30 要素から成る相互作用ネットワークに 含まれる相互作用のほとんどを 推定することができた。簡便な処理で推定を始めるこ とができるため、RA 法を用いて大まかなネットワーク構造を決定した後、詳細を検証 していくような場合にも用いることができると考えられる。 52 TP率 第3章 Random selection method 100 90 80 70 60 50 40 30 20 10 0 0 5 10 15 SR Fig. 3.4 30 要素ネットワークモデル推定における TP 率の変化 53 第4章 ノイズを含むデータへの有効性 第4章 ノイズを含むデータへの有効性 第1節 緒言 実問題への適用を目指す上で、提案法への入力値となるタイムコースデータに誤差が 含まれることを考慮しなければならない。 一般的に、遺伝子やタンパク質、代謝物質などのタイムコースデータを得るための実 験では、計測機器や実験者の手技の影響から誤差が含まれる場合がほとんどである。こ れまで、30 要素ネットワークモデル等を用い提案法の有効性を検証してきたが、用い たデータは誤差を考慮していない。そこで、本章では誤差を含むデータを用いる場合で も、AREX+JGG 法およびハイブリッド法、RA 法が有効であるか検証した。 第2節 ノイズを含むタイムコースデータの取得 実験データに含まれる誤差には、計測機器や生物学的ノイズ、実験者の手技等様々な 要因が考えられる。本研究では、誤差が含まれるようなデータを再現するために Langevin 方程式(Gillespie, 2000)を用いた。Langevin 方程式は、Paul Langevin により 定義されたブラウン運動を記述する確率微分方程式である。 本研究では Langevin 方程式と S-system モデル(式(1.2))を組み合わせることにより、 ノイズを含むタイムコースデータを作成した。 n n n n g h g h dxi i x j ij i x j ij i x j ij i x j ij dt j 1 j 1 j 1 j 1 54 (3.1) 第4章 ノイズを含むデータへの有効性 𝜂はホワイトノイズであり、右辺の第 3 項及び第 4 項がノイズ項である。本研究では、 ホワイトノイズは Mox-Muller 法(四辻, 2010)を用い、一様分布に従う確率変数から正 規分布に従う擬似乱数を生成させた。実験データに含まれる誤差には様々な要因が考え られるため、データ点当たりにノイズを含ませることのできるこの方法は、本検証のよ うに入力値としてのデータにおけるノイズを考慮する場合には適していると考えられ る。 第3節 ネットワークモデル 本章での検証では、まず、小規模なネットワークモデルを用い AREX+JGG 法および ハイブリッド法の有効性を検証し、その後 30 要素ネットワークモデルを用い RA 法の 有効性を検証した。 小規模なネットワークモデルとして、4 要素ネットワークモデルと 7 要素ネットワー クモデルを用意した(Fig. 4.1) 。これらのネットワークモデルを基に S-system 式を立式 し、Langevin 方程式を導入することによりノイズを含むタイムコースデータを生成した。 タイムコースデータのサンプリングポイント数は 71 点で、ネットワークに含まれる要 素数 n をとすると、(n + 1)個のタイムコースデータを用意した。(n + 1)個の内 1 つは正 常に全要素が発現している場合を想定し、残り n 個は 1 要素破壊実験を想定している。 1 要素破壊実験のタイムコースデータは、破壊対象となる要素を i とすると、状態変数 xi の生成過程の速度定数 αi を 0.0 とすることで得た。これらのネットワークにおいて、 推定対象とした S-system のパラメータを Table 4.1 に示した。 30 要素ネットワークモデルについては、第 3 章第 2 節に示したネットワークモデル (Fig. 3.1)を基に、Langevin 方程式と組み合わせることでノイズを発生させた。ネット ワークモデルの S-system のシステムパラメータ値を式(3.1)へ代入し、それぞれの状態変 数につき同じ初期値で 3 回計算することでタイムコースを 3 セットずつ作成した。確率 55 第4章 ノイズを含むデータへの有効性 微分方程式であるので、 状態変数 1 つにつき 3 種のタイムコースが得られることになる。 3 セットの平均値を基に標準偏差と相対誤差を算出し、平均値を相互作用ネットワーク 推定のための入力データとして用い、標準偏差と相対誤差をタイムコースデータのノイ ズの大きさの指標とした。 56 第4章 (a) ノイズを含むデータへの有効性 (b) 0.4 1 0.7 0.3 0.5 -0.2 2 1 1.0 3 4 3 2 -0.2 0.5 0.5 0.4 6 5 4 0.7 0.5 7 Fig. 4.1 ネットワークモデルの構造 (a) 4 要素ネットワークモデル (b) 7 要素ネットワークモデル 図中の数字は S-system モデル(式(1.2))における gij の値を表している Table 4.1 推定対象 推定対象 推定対象としない パラメータ数 S-system パラメータ S-system パラメータ αi, βi (i = 1, 2, 3, 4) 4 要素 ネットワーク 28 モデル モデル gij (i = 1, 2, 3, 4, j = 1, 2, 3, 4) hii = 0 (i≠j) hii (i = 1, 2, 3, 4) αi, βi (i = 1, 2, 3, 4, 5, 6, 7) 7 要素 ネットワーク 推定対象とした S-system パラメータ 70 gij (i = 1, 2, 3, 4, 5, 6, 7, j = 1, 2, 3, 4, 5, 6, 7) hii (i = 1, 2, 3, 4, 5, 6, 7) 57 hii = 0 (i≠j) 第4章 第4節 ノイズを含むデータへの有効性 検証 第 1 項 AREX+JGG 法およびハイブリッド法のノイズを含むタイムコースデー タに対する有効性 AREX+JGG 法およびハイブリッド法のノイズ耐性を検証するために、4 要素ネットワ ークモデルと 7 要素ネットワークモデルを用いて実数値 GA に基づくネットワーク推定 法である SPX+MGG 法および UNDX+MGG 法との推定能力の差を比較した。それぞれ の手法に用いられた実数値 GA のためのパラメータを Table 4.2 に示した。 まず、4 要素ネットワークモデルを用い検証した。このとき、入力データとして用い られたタイムコースデータの 1 点当たりの平均標準偏差は約 0.14、平均相対誤差は約 20%であった。各手法における推定の終了条件は、実数値 GA のエリート個体の平均相 対誤差が 25%より小さくなったとき推定成功とし、世代数が 10,000 を超えたとき推定 に失敗したこととした。ハイブリッド法においては、エリート個体の平均相対誤差が 35%より小さくなったときに修正 Powell 法を適用した(Table 4.3(a)) 。4 要素ネットワ ークモデルを各手法 30 試行ずつ推定した結果、UNDX+MGG 法を除く 3 つの手法はす べての試行で評価値が目標平均相対二乗誤差以下となり推定に成功した(Table 4.4) 。 また、4 要素ネットワークモデルに含まれるすべての TP interaction をどの手法でも推定 することができたが、AREX+JGG 法と UNDX+MGG 法で FP interaction が 1 つ推定され た。推定に要した CPU time では、ハイブリッド法がもっとも優れていた。 次に、7 要素ネットワークモデルを用い検証した。入力データとして用いたタイムコ ースデータの 1 点当たりの平均標準偏差は約 0.08、平均相対誤差は約 12%であった。各 手法における推定の終了条件は、実数値 GA のエリート個体の平均相対二乗誤差が 20% より小さくなったとき推定成功とし、世代数が 10,000 を超えたとき推定に失敗したこ ととした。ハイブリッド法においては、エリート個体の平均相対二乗誤差が 30%より小 58 第4章 ノイズを含むデータへの有効性 さくなったときに修正 Powell 法を適用した(Table 4.3(b))。各手法で 30 試行ずつ推定し た結果、UNDX+MGG 法と SPX+MGG 法では設定した最大世代数の内に目標平均相対 二乗誤差に達することができず、すべての試行で推定が失敗した。AREX+JGG 法とハ イブリッド法では、30 試行の内、25 試行で推定に成功した。推定に要した CPU time は、 AREX+JGG 法およびハイブリッド法の間でほとんど差はなかった。AREX+JGG 法では、 8 つの相互作用の内 7 つが TP interaction として推定され、1 つを推定することができな かった。しかし、ハイブリッド法ではすべての TP interaction を推定することに成功し た。 以上の結果から、誤差に対する耐性はハイブリッド法がもっとも優れていると考えら れる。 59 60 28 70 4 7 7 70 7 (AREX+JGG+修正 Powell 法) 28 4 70 28 70 7 4 28 パラメータ数 4 ハイブリッド法 SPX+MGG UNDX+MGG AREX+JGG 要素数 推定対象 2100 840 300 300 300 300 2100 840 集団サイズ 280 112 300 100 300 100 280 112 子個体数 71 29 71 29 2 2 71 29 親個体数 実数値 GA の設定パラメータ Table 4.2 各手法の実数値 GA 設定パラメータ 第4章 ノイズを含むデータへの有効性 第4章 ノイズを含むデータへの有効性 推定のための最適化条件 Table 4.3 (a) 4 要素ネットワークモデル(Fig. 4.1(a)) 項目 設定 最適化対象パラメータ αi、βi、gij、hij (i ≠ j) タイムコースデータのサンプリングポイント数 71 タイムコースデータの平均標準偏差 0.14 タイムコースデータの平均相対誤差 20% Rosenbrock 法 ハイブリッド法 𝑥* 1.618034 𝑦* 0.618034 刻み幅 h* 0.1 終了条件 0.01 修正 Powell 法 開始平均相対二乗誤差(%) < 35.0 cut off 0.05 最大世代数 100,000 目標平均相対二乗誤差(%) < 25.0 *第 2 章第 4 節第 4 項参照 (b) 7 要素ネットワークモデル(Fig. 4.1(b)) 項目 設定 最適化対象パラメータ αi、βi、gij、hij (i ≠ j) タイムコースデータのサンプリングポイント数 71 タイムコースデータの平均標準偏差 0.08 タイムコースデータの平均相対誤差 12% Rosenbrock 法 ハイブリッド法 𝑥* 1.618034 𝑦* 0.618034 刻み幅 h* 0.1 終了条件 0.01 修正 Powell 法 開始平均相対二乗誤差(%) < 30.0 cut off 0.05 最大世代数 100,000 目標平均相対二乗誤差(%) < 20.0 *第 2 章第 4 節第 4 項参照 61 62 8 4 (gij (when i≠j)) 相互作用数 25 0 UNDX+MGG ハイブリッド法 25 AREX+JGG 0 30 ハイブリッド法 SPX+MGG 30 13 UNDX+MGG SPX+MGG 30 AREX+JGG (/30) 成功試行回数 7 0 0 7 7 4 4 0 0 0 0 0 0 0 1 interaction interaction 4 FP TP FN 1 8 8 1 1 0 0 0 interaction 抽出された相互作用数 10727.7 33384.6 33468.3 10694.6 167.10 1258.92 2814.32 249.93 [sec]*/試行 平均 CPU time 4 要素ネットワークモデルおよび 7 要素ネットワークモデルの推定結果 *CPU : Xeon E5540 2.53GHz, Memory : 8.0GB 7 4 要素数 ネットワーク中の Table 4.4 321829.4 1001522 1002909 320838.3 5013.08 37767.52 84429.67 7498.06 [sec]*/30 試行 CPU time 第4章 ノイズを含むデータへの有効性 第4章 第2項 ノイズを含むデータへの有効性 許容される最大ノイズの検証 前項での検証の結果、ノイズに対する耐性はハイブリッド法が最も高く、次にノイズ 耐性が高い手法は AREX+JGG 法であった。本項では、2 つの手法がどの程度の大きさ のノイズまで推定能力を維持できるかを検証した。検証には 4 要素ネットワークモデル (Fig. 4.1(a))を用いた。 入力データとして、1 点当たりの平均標準偏差が約 0.17、平均相対誤差が約 29%のタ イムコースデータと、1 点当たりの平均標準偏差が約 0.19、平均相対誤差が約 32%のタ イムコースデータの 2 つを用意した。ここで、前者のタイムコースデータを用いた検証 を case2-01、後者を用いた検証を case2-02 とする。本検証では、実数値 GA で用いる終 了条件となっている評価値(目標平均相対二乗誤差)の値を 3 種類ずつ設定した。ハイ ブリッド法における修正 Powell 法適用条件は、実数値 GA のエリートの平均相対二乗 誤差が終了条件の平均相対二乗誤差+10%となったときとした。 case2-01 の結果、AREX+JGG 法およびハイブリッド法の両者で、終了条件を【目標 平均相対二乗誤差<40%】としたときすべての TP interaction を推定することに成功した (Table 4.5) 。終了条件を【目標平均相対二乗誤差<35%】とした場合、いずれの試行に おいても設定された世代数の内に推定することができなかった。また、終了条件を【目 標平均相対二乗誤差<45%】としたとき、3 つの TP interaction を推定できたが、1 つ推 定することができなかった。これは、終了条件が【目標平均相対二乗誤差<40%】のと きと比べ、終了条件が緩和されることで、局所解に収束するようになったためであると 考えられる。 次に、case2-02 の結果、終了条件を【目標平均相対二乗誤差<35】、 【目標平均相対二 乗誤差<40%】としたとき、AREX+JGG 法およびハイブリッド法どちらの手法におい ても、すべての試行で推定に失敗し、TP interaction を推定することができなかった。終 了条件を【目標平均相対二乗誤差<45%】とした場合、AREX+JGG 法およびハイブリ ッド法の両者で、2 つの TP interaction を推定することに成功した。このとき、30 試行 63 第4章 ノイズを含むデータへの有効性 すべてで推定に成功したのにもかかわらず推定できなかった相互作用(FN interaction) が 2 つあった。したがって、case2-02 の場合、AREX+JGG 法およびハイブリッド法は、 相互作用のすべてを推定する能力はないとみられる。 以上の結果から、AREX+JGG 法およびハイブリッド法が有効性を示す上で許容され る最大のノイズは、入力データのタイムコースの 1 点当たりの平均相対誤差にすると 29~32%の間であると考えられる。 64 65 0.19 0.17 32% 29% 平均相対誤差 ハイブリッド法 AREX+JGG ハイブリッド法 AREX+JGG 0 30 <40 % <45 % 0 30 <45 % <35 % 0 0 <35 % <40 % 30 29 <45 % <40 % 0 30 <45 % <35 % 28 0 (/30) 成功試行回数 <40 % <35 % 終了条件* ***CPU : Xeon E5540 2.53GHz, Memory : 8.0GB **TP:TP interaction, FP:FP interaction, FN:FN interaction *実数値 GA のエリート個体のタイムコースデータ 1 点当たりの平均相対二乗誤差 case2-02 case2-01 平均標準偏差 Table 4.5 許容誤差の検証結果 2 0 0 2 0 0 3 4 0 3 4 0 TP 0 0 0 0 0 0 0 0 0 0 0 0 FP 相互作用数** 抽出された 2 4 4 2 4 4 1 0 4 1 0 4 FN 249.05 3354.04 3361.56 271.88 3373.96 3362.05 88.2 172.33 3407.84 103.14 210.61 3353.07 [sec]***/試行 平均 CPU time 第4章 ノイズを含むデータへの有効性 第4章 第3項 ノイズを含むデータへの有効性 ノイズを含むデータに対する RA 法の有効性 RA 法のノイズ耐性を検証するため、30 要素ネットワークモデル(Fig. 3.2、Table 3.1) を用いてノイズを含むタイムコースデータを生成した。ここで、ノイズを含むタイムコ ースデータを用いた場合の検証を case3-01、対照実験のためのノイズを含まないタイム コースデータを用いた場合の検証を case3-02 とする。case3-01 および case3-02 で設定し た、RA 法のパラメータを Table 4.6 に示した。パラメータ k については、2~5 の 4 種設 定し、パラメータ p は SR が設定された値を満たすよう定めた。RA 法で適用されてい るハイブリッド法に用いられたパラメータを Table 4.7 および Table 4.8 に示した。 66 第4章 Table 4.6 ノイズを含むデータへの有効性 検証のために用いた RA 法のパラメータ 設定値 Table 4.7 2-5 l 30 m 30 SR 2×4 ハイブリッド法の実数値 GA(AREX+JGG)設定パラメータ 要素数 ハイブリッド法 k 実数値 GA(AREX+JGG)の 推定対象 設定パラメータ パラメータ数 集団サイズ 子個体数 親個体数 2 10 300 40 11 3 18 540 72 19 4 28 840 112 29 5 40 1200 160 41 67 第4章 Table 4.8 ノイズを含むデータへの有効性 推定のための最適化条件 項目 設定 最適化対象パラメータ αi、βi、gij、hij (i ≠ j) タイムコースデータのサンプリングポイント数 71 Rosenbrock 法 ハイブリッド法 𝑥* 1.618034 𝑦* 0.618034 刻み幅 h* 0.1 終了条件 0.01 修正 Powell 法 開始平均相対二乗誤差(%) (目標平均相対二乗誤差+10) cut off 0.05 最大世代数 100,000 *第 2 章第 4 節第 4 項参照 68 第4章 ノイズを含むデータへの有効性 case3-01 で用いた、1 要素破壊実験を想定したデータを含めた 31 種のタイムコースデ ータ中の平均標準偏差は約 0.09、平均相対誤差は約 10%であった。推定の終了条件であ る目標平均相対二乗誤差を 30%とし、エリート個体の平均相対二乗誤差が 30%より小 さくなったとき推定に成功したこととした。推定の結果、k=5、SR=8 のときすべての 相互作用の内およそ 80%の相互作用を推定することに成功した(Table 4.9) 。 case3-02 では、ノイズを含まないタイムコースデータを入力値として用いているため、 推定の終了条件である目標平均相対二乗誤差を 10%とし、エリート個体の平均相対二乗 誤差が 30%より小さくなったとき推定に成功したとした。推定の結果、k=5、SR=8 の ときすべての相互作用の内およそ 78%の相互作用を推定することに成功した(Table 4.10) 。 以上の結果、case3-01 および case3-02 どちらの場合においても、RA 法の性能に違い はほとんどなく、RA 法がノイズを含むデータに対しても有効であると考えられる。 69 70 4 6 8 2 4 6 8 30 45 60 12 24 36 48 4 40 2 2 20 15 8 120 8 6 90 80 4 60 6 2 30 60 SR p 48 36 24 12 60 45 30 15 61 48 33 16 40 31 21 11 (/p) 成功推定回数 *CPU : Xeon E5540 2.53GHz, Memory : 8.0GB 5 4 3 2 (k) 要素数 case3-01 の結果 54 52 46 40 45 42 41 35 44 40 39 34 27 25 24 16 TP interaction 9 8 4 2 9 8 7 2 8 5 4 2 1 0 0 0 FP interaction 14 16 22 28 23 26 27 33 24 28 29 34 41 43 44 52 FN interaction 抽出された相互作用数(/68) Table 4.9 79.41 76.47 67.65 58.82 66.18 61.76 60.29 51.47 64.71 58.82 57.35 50.00 39.71 36.76 35.30 23.50 TP 率(%) 12.61 12.53 12.50 13.12 3.50 3.29 3.28 3.36 2.59 2.56 2.36 2.68 1.65 1.62 1.61 1.59 [h]*/試行 平均 CPU time 605.23 450.99 299.92 157.49 209.82 148.08 98.29 50.42 207.07 153.38 94.60 53.69 197.97 146.00 96.31 47.65 [h]*/ p 試行 CPU time 第4章 ノイズを含むデータへの有効性 71 4 6 8 2 4 6 8 30 45 60 12 24 36 48 4 40 2 2 20 15 8 120 8 6 90 80 4 60 6 2 30 60 SR p 48 36 24 12 60 45 30 15 80 60 40 20 119 99 59 30 (/p) 成功推定回数 *CPU : Xeon E5540 2.53GHz, Memory : 8.0GB 5 4 3 2 (k) 要素数 case3-02 の結果 53 49 45 39 46 44 39 34 46 42 40 38 44 38 36 34 TP interaction 8 5 4 2 8 5 3 0 8 8 6 6 8 7 4 2 FP interaction 15 19 23 29 22 24 29 34 22 26 28 30 22 30 32 34 FN interaction 抽出された相互作用数(/68) Table 4.10 77.94 72.06 66.18 57.35 67.65 64.71 57.35 50.00 67.65 61.76 58.82 55.88 64.71 55.88 52.94 50.00 TP 率(%) 30.95 31.09 30.91 29.85 9.02 9.09 9.21 9.17 1.58 1.70 1.58 1.63 0.10 0.11 0.12 0.09 [h]*/試行 平均 CPU time 1485.77 1119.14 741.86 358.22 541.64 409.21 276.44 137.48 128.99 102.00 63.24 32.59 12.45 9.99 7.39 2.77 [h]*/ p 試行 CPU time 第4章 ノイズを含むデータへの有効性 第4章 第5節 ノイズを含むデータへの有効性 考察 本章では、AREX+JGG 法およびハイブリッド法、そして RA 法のノイズ耐性を検証 した。実問題に適用する場合、入力値として用いる実験データは高い確率で何らかの誤 差を含んでいる。したがって、実用的なネットワーク推定手法かどうか示すために、提 案法が誤差を含んでいるようなデータに対しても有効性を示すか検証することが求め られる。 本章では、まず、4 要素ネットワークモデル(Fig. 4.1(a))と 7 要素ネットワークモデ ル(Fig. 4.1(b))を用い、AREX+JGG 法、UNDX+MGG 法、SPX+MGG 法、ハイブリッ ド法それぞれのノイズ耐性について検証した。4 要素ネットワークモデルを推定した結 果、タイムコースデータに約 20%の平均相対誤差を含んでいる場合でもすべての手法で ネットワークモデルに含まれている相互作用を推定することに成功した。さらに、 AREX+JGG 法および SPX+MGG 法、ハイブリッド法は 30 試行すべての試行で推定に成 功した。しかし、UNDX+MGG 法では 30 試行中 13 試行しか推定を成功させることがで きず、推定に要した時間も長くかかった。これらの手法の交叉法である AREX、UNDX、 SPX のアルゴリズムに着目すると、AREX および SPX は多親交叉であるのに対し、 UNDX は 2 親交叉で子個体を生成する。実数値 GA は親個体を基に子個体を生成し、ほ とんどの場合親個体の周囲に子個体が生成されるため、UNDX は他の手法に比べ探索空 間にばらつきが少なく、最適解に対する収束性も遅くなると考えられる。それが原因と なり、定められた最大世代交代数内に目標平均相対二乗誤差に達することができず、成 功試行回数が半分以下となったと考えられる。また、ハイブリッド法は、AREX+JGG 法をベースとしたアルゴリズムにもかかわらず、AREX+JGG 法の推定結果とは異なる 結果となった。ハイブリッド法ではエリート個体が特定の評価値に達した際にエリート 個体に対して修正 Powell 法が適用される。これは、探索後期に局所収束を図り探索時 間を短縮させる狙いであるが、本検証におけるケースでは最後に局所収束させることで 72 第4章 ノイズを含むデータへの有効性 すべての TP interaction を推定することができたと思われる。したがって本検証で扱っ た問題では大域的最適解の近くに局所的最適解が存在していると考えられ、AREX+JGG 法はその局所的最適解に収束してしまい、4 要素ネットワークモデルでは 1 つの FN interaction が推定され、7 要素ネットワークモデルでは 1 つの TP interaction を推定でき なかったと予測される。 次に、AREX+JGG 法とハイブリッド法の許容誤差を、4 要素ネットワークモデルを用 いて検証したところ、両手法とも平均相対誤差にすると 29%~32%の間に許容誤差があ ると考えられた(Table 4.5) 。また、case2-01 および case2-02 のどちらの検証においても、 入力データに含まれている平均相対誤差に近い値を推定手法の目標平均相対二乗誤差 に設定した場合には、正解のネットワーク構造を推定することができなかった。したが って、誤差を含むようなデータを用いる場合には、推定の終了条件である目標平均相対 二乗誤差を緩やかな条件にしなければならないことがわかった。しかしながら、case01-1 の【目標平均相対二乗誤差<40%】の結果と【目標平均相対二乗誤差<45%】の結果か ら、目標平均相対二乗誤差に余裕を持たせ過ぎると、正解の相互作用を推定できなくな ってしまうことがあった。case2-01 および case2-02 から、目標平均相対二乗誤差を入力 データの平均相対誤差におよそ 10%加えた条件に設定すると最適であることが示唆さ れた。 30 要素ネットワークモデルを用い RA 法のノイズ耐性を検証した。その結果、ノイズ を含まないデータとノイズを含むデータのどちらを用いた場合でも RA 法の推定能力に ほとんど差はなかった。case3-01 および case3-02 どちらの推定でも、k=5 としたときが 最も TP 率が高く、真のネットワークに含まれる相互作用の内約 80%を推定することに 成功した。k=4、k=5 としたとき、case3-01 に比べて case3-02 が推定により多くの CPU time が掛かってしまった原因としては、case3-02 の方が推定条件の目標平均相対二乗誤 差がより厳しく設定されていたことが挙げられる。一方で、k=2、k=3 としたときは、 case3-02 に比べて case3-01 が推定により多くの CPU time がかかっており、k=4、k=5 の場合と逆転している。これは、case3-01 では、要素数が少ないことに加え、誤差を含 73 第4章 ノイズを含むデータへの有効性 むデータを用いたため、各試行で入力データを再現できるような S-system モデルのパ ラメータの組み合わせを見つけることが難しく、推定に成功した回数が少なかったため であると思われる。したがって、RA 法を適用する際には、適切なパラメータ k の値を 設定することが、効率の良い推定に繋がる。 74 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 第5章 異なるネットワークモデルに対する Random selection method の有効性 第1節 緒言 遺伝子制御ネットワークやタンパク質間相互作用ネットワークなどの生体分子によ るネットワークに限らず、World Wide Web(WWW)のハイパーリンク構造や人間関係 のネットワークなど、現実世界には様々なネットワークが存在している(Weng, Bhalla, and Iyengar, 1999; Albert, Jeong, and Barabasi, 1999) 。それらの複雑ネットワークは多様な 構造を有しているが、いくつかの性質を持つことが報告されている。 ポール・エルデシュとアルフレッド・レーニィが 1959 年に各ノードをランダムにリ ンクさせるランダムグラフ理論(Erdős-Rényi モデル; ER モデル)を提唱してから、複 雑ネットワークはランダムネットワークであるとされてきた(Erdös and Rényi, 1959) 。 しかし、1999 年にスモールワールドネットワーク(Watts and Strogatz, 1998)や、スケ ールフリーネットワーク(Barabási and Albert, 1999)が提唱され、現実世界の複雑ネッ トワークは必ずしもランダムネットワークでは表現できないことが明らかとなってき た。 本章では、従来複雑ネットワークの構造として考えられてきたランダムネットワーク と、近年研究が盛んに行われているスケールフリーネットワークの 2 つのネットワーク モデルを用い、両者に対する RA 法の有効性を検証した。 75 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 第2節 第1項 ランダムネットワークとスケールフリーネットワーク ランダムネットワーク ER モデルでは、ある 2 つのノード間を確率 p でエッジが接続され、確率(1-p)でエッ ジが接続しないようにしてネットワークを生成する。各ノード間にエッジが置かれるか どうかは、他のペアからの影響は受けず独立して決定される。確率 p が小さいときネッ トワーク中のエッジの総数が過疎な状態になり、確率 p が大きいときネットワーク中の エッジ数は密な状態となる。p=1 のとき、すべてのノード間にエッジが接続されるた め、完全グラフとなる。ランダムネットワークでは、ノードはランダムに接続されるの で、大部分のノードは同等数のエッジを持っていることとなる。1 つのノードが持つエ ッジの数を次数と呼び、ランダムネットワークにおいては次数が k であるノードの分布 はポアソン分布に従うとされている。 第2項 スケールフリーネットワーク スケールフリーネットワークは、ノードの次数分布がべき乗則に従うスケールフリー 性を有するようなネットワークである(Barabási and Albert, 1999) 。べき乗則の代表的な 例として、ヴィルフレード・パレートが発見した所得の分布が挙げられる。パレートは、 所得の多い人口は少なく、大多数はわずかな所得しかなく、所得の分布がべき乗則に従 うことを示しており、これはパレートの法則と呼ばれる(Hu et al., 2006; Fox Keller, 2005)。この考え方をグラフ理論へと導入し、各ノードが持つエッジ数の度数分布がべ き乗則に従うようなネットワークをスケールフリーネットワークとした。スケールフリ ーネットワークでは、全ノードの中で次数が k であるノードの分布を f(k)とすると、次 数の分布は次の式(5.1)で表される(Komatsu and Namatame, 2012) 。 76 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 f (k ) k (5.1) λ は定数である。 スケールフリーネットワークは、1999 年にアルバート=ラズロ・バラバシとレカ・ アルバートにより提唱され、これまで、現実世界の多くの複雑ネットワークはスケール フリー性を持つことが報告されている(Barabási and Albert, 1999)。たとえば、インター ネット(Faloutsos, Faloutsos, and Faloutsos, 1999)やタンパク質ドメイン間相互作用ネッ トワーク(Wuchty, 2001)においてもスケールフリー性が観察されている。 第3節 第1項 100 要素ネットワークモデルの作成 100 要素ランダムネットワークモデル ランダムネットワークにおける RA 法の有効性を検証するために、100 要素ネットワ ークモデルを作成することとした。ネットワークモデルは、R のネットワーク分析を支 援するパッケージである igraph(Csárdi and Nepusz, 2006)を利用し、ER モデルに基づ き作成した。 ER モデルでは、まず n 個のノードを用意してから、任意の 2 つのノードをランダム に選択しエッジを追加していく。このとき、各ノード間でエッジが張られる確率は p で ある。本検証では、確率 p=0.02 とし 100 要素ランダムネットワークモデルを作成した (Fig. 5.1) 。図中の赤い矢印は正の相互作用を、青い矢印は負の相互作用を表す。作成 したネットワークモデルのそれぞれの相互作用に、ランダムに[-1.0, 1.0]の値を振り分け、 その値を式(1.2)中の gij(i≠j)として S-system モデルを立式した。各相互作用係数の値 を Table 5.1 に示した。作成したネットワークモデル中の相互作用(gij)の総数は 295 77 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 であり、i≠j のときに限ると 195 となる。それぞれの要素が持つ相互作用の次数分布は Fig. 5.2 に示したとおりであり、2~6 つの相互作用を持つ要素がほとんどを占めている。 78 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 Fig. 5.1 100 要素ランダムネットワーク 赤い矢印は正の制御を、青い矢印は負の制御を表す 79 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 Table . 5.1 100 要素ランダムネットワークモデルの相互作用係数 i j 相互作用係数 i j 相互作用係数 i j 相互作用係数 1 19 -0.5 20 51 -1.0 36 28 -1.0 2 40 -0.8 20 62 -0.1 37 63 -0.5 2 56 0.7 20 77 -0.8 37 66 0.5 2 86 -0.8 20 88 -0.2 37 77 -0.4 2 100 0.7 21 38 -0.3 38 39 1.0 3 92 0.7 22 13 -0.3 38 98 0.3 4 19 0.4 22 61 -0.9 39 89 1.0 4 43 0.4 23 37 1.0 40 12 0.4 4 55 0.6 23 39 0.1 40 42 1.0 4 69 -0.3 23 71 -0.7 43 13 0.2 5 62 -0.9 24 35 0.2 43 25 0.8 6 50 0.8 24 74 -1.0 43 42 -0.9 6 77 -0.7 25 79 -0.9 43 44 -0.3 6 91 -0.3 26 4 -1.0 44 7 0.7 8 19 -0.5 26 96 0.8 44 8 0.2 9 57 0.6 27 29 -0.8 44 33 0.1 9 92 -0.3 27 52 -0.6 44 40 -0.8 10 27 1.0 27 100 1.0 45 65 -0.3 11 9 -0.8 28 12 0.6 45 99 0.6 11 49 -0.1 29 1 0.5 47 31 0.1 11 82 -0.9 29 52 0.9 47 64 -0.7 14 24 1.0 29 90 0.4 48 31 -0.7 14 100 -1.0 30 79 -0.2 48 71 0.7 15 41 0.4 31 49 0.9 48 73 0.6 15 98 -0.4 32 3 -0.7 49 44 -0.4 16 51 -0.1 32 82 -0.5 49 82 0.8 17 3 -0.5 33 19 -0.9 49 96 -0.3 100 17 -0.4 33 50 -0.6 50 33 -0.8 17 18 1.0 33 67 0.1 50 52 -0.4 17 37 -0.5 33 82 0.1 50 83 0.6 17 95 0.7 34 21 0.3 51 10 -0.2 19 15 -0.3 100 34 0.7 52 83 -0.9 20 35 -0.4 35 55 -0.6 52 84 0.9 Table 5.1 continued on next page 80 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 Table 5.1 continued i j 相互作用係数 i j 相互作用係数 i j 相互作用係数 54 78 0.3 67 44 -0.4 86 67 0.5 54 92 -0.1 67 98 -0.5 87 39 0.8 54 94 0.3 68 3 -0.1 87 62 -0.6 55 48 0.6 69 18 -1.0 87 72 -0.2 55 76 -0.9 70 7 -1.0 87 88 -0.8 56 13 -0.1 71 14 0.3 88 20 -0.7 56 41 -1.0 71 79 -0.8 88 77 -0.3 57 25 0.5 71 90 0.6 88 89 0.7 58 37 0.4 73 19 0.2 89 82 0.7 58 52 0.9 73 54 -0.4 90 54 -0.9 58 67 0.6 74 50 -1.0 90 73 -0.2 59 24 0.5 75 74 -0.3 90 76 0.3 59 31 -0.1 76 90 0.5 91 54 -0.7 59 46 0.3 78 11 -0.4 92 14 0.1 59 65 -0.1 78 43 -0.6 92 31 0.1 59 87 -0.2 79 98 -0.3 92 66 -0.4 60 8 0.1 80 53 -0.6 93 89 -0.9 60 12 -0.6 80 87 -0.3 94 2 0.7 60 29 -0.1 81 85 0.4 94 10 0.4 60 43 0.5 81 89 0.7 94 56 -0.1 60 49 0.8 82 7 -0.9 94 83 1.0 60 54 0.7 82 51 -0.7 94 86 -0.2 62 27 -0.3 82 65 1.0 94 91 1.0 62 69 -0.9 83 71 -0.4 95 62 0.5 63 55 -0.2 83 93 0.2 95 66 -0.2 63 61 0.6 84 8 0.7 95 98 0.4 64 16 -0.6 84 32 0.3 96 73 0.4 64 22 -0.7 85 9 0.1 97 32 -0.2 64 91 0.3 85 54 0.7 98 55 -0.9 64 95 -0.3 85 63 -0.8 99 42 0.5 66 46 -0.9 86 11 -0.3 66 75 0.8 86 31 -0.4 67 9 0.9 86 50 -0.6 81 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 Fig. 5.2 100 要素ランダムネットワークモデルの次数分布 82 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 第2項 100 要素スケールフリーネットワークモデル スケールフリーネットワークにはいくつかのモデルが提唱されている。本検証では、 スケールフリーネットワークの代表的なモデルである Barabási-Albert モデル(BA モデ ル) (Barabási and Albert, 1999)を用いて 100 要素スケールフリーネットワークモデルを 作成した。 BA モデルはネットワークの成長(growth)と優先的選択(preference attachment)と いう 2 つの大きな特徴を持つ。成長という特徴を取り入れるため、BA モデルでは初期 のネットワークに時間経過ごとに新しくノードが加わることとされている。また、新た に加わったノードが他のノードと接続される際に、次数の高いノードとより接続しやす くなるようにすることを優先的選択と呼ぶ。これらの特徴を取り入れた、BA モデルに よるネットワーク作成方法は次 1~3 の通りである。 1.m0 個ノードから成る完全グラフを作成する 2.作成した完全グラフに、新たなノードを 1 個追加し、既に存在しているノード に対して接続する。このとき、あるノードが新たなノードに接続される確率は次 の式に従う。 (ki ) ki j kj (5.2) 3.2 の操作を指定したエッジ数になるまで繰り返す。 kj はノード i の次数であり、次数の高いノードほど確率 Π(kj)が高くなる。 BA モデルにより 100 要素スケールフリーネットワークモデルを作成するために、ト ポロジジェネレータである BRITE(Medina et al., 2001)を用いた。作成した 100 要素ス ケールフリーネットワークモデルを Fig. 5.3 に示した。図中の赤い矢印は正の相互作用 83 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 を、青い矢印は負の相互作用を表す。それぞれの相互作用には、ランダムに[-1.0, 1.0] の値を振り分けた。振り分けた値を式(1.2)中の gij(i≠j)とし、その値を Table 5.2 に示 した。作成したネットワークモデル中の相互作用(gij)の総数は 287 であり、i≠j のと きに限ると 187 となる。それぞれの要素の次数分布は Fig. 5.4 に示したとおりである。 84 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 Fig. 5.3 100 要素スケールフリーネットワークモデル 赤い矢印は正の制御を、青い矢印は負の制御を表す 85 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 Table 5.2 100 要素スケールフリーネットワークモデルの相互作用係数 i j 相互作用係数 i j 相互作用係数 i j 相互作用係数 1 2 -0.8 17 2 0.2 33 8 -0.6 1 3 0.9 20 1 0.1 34 2 0.1 2 3 0.9 20 5 -0.1 34 1 -0.3 4 2 -0.2 19 13 1.0 39 11 -0.3 4 1 -0.8 19 6 -0.1 39 2 -0.5 5 3 -0.4 22 2 0.2 40 20 -0.1 5 1 -0.5 22 1 -0.6 37 18 -0.2 6 2 -0.3 21 2 -0.4 38 1 -0.7 6 5 0.4 21 5 0.1 38 7 0.1 7 2 0.1 24 17 -0.8 43 39 0.7 7 1 -0.4 24 14 -0.4 43 13 0.6 8 1 -0.8 23 13 -0.3 44 4 0.7 8 5 -0.4 23 16 0.2 44 33 0.6 9 7 -0.3 26 13 -0.1 41 10 0.3 9 2 0.6 26 7 -0.8 41 6 0.8 10 2 -0.4 25 4 0.1 42 40 0.2 10 7 0.4 25 1 -0.9 42 8 -0.6 11 1 0.7 28 1 0.5 47 1 0.2 11 7 -0.5 28 22 -0.9 47 13 0.5 12 2 0.1 27 9 0.9 48 16 -0.3 12 8 -0.1 27 11 -0.6 45 2 -0.8 13 11 0.5 30 24 0.2 45 4 0.2 13 1 0.3 29 5 -0.3 46 44 -0.3 14 7 -0.4 32 24 0.5 46 14 -0.9 14 8 0.4 32 3 0.9 52 46 -0.5 15 7 -0.7 31 5 0.5 52 5 0.9 15 2 0.1 31 1 0.1 51 7 0.8 16 11 0.2 35 24 0.3 51 5 -0.7 16 8 -0.4 35 5 -0.1 50 1 -0.2 18 5 -0.9 36 2 0.8 50 46 0.1 18 16 -0.8 36 10 -0.5 49 25 -0.1 17 4 0.8 33 5 -0.7 49 10 -0.8 Table 5.2 continued on next page 86 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 Table 5.2 continued i j 相互作用係数 i j 相互作用係数 i j 相互作用係数 56 16 -0.9 66 47 0.7 84 72 0.3 56 51 -1.0 67 3 0.7 84 1 0.6 55 7 0.8 67 11 1.0 83 16 -0.5 55 20 -0.2 68 49 0.5 83 56 0.1 54 1 0.3 68 39 -0.1 94 1 -0.6 54 2 -0.8 77 51 -1.0 94 5 0.3 53 3 0.1 77 39 -0.3 93 27 -0.2 53 21 -0.7 78 53 -0.1 93 5 0.7 60 8 -0.9 78 24 -0.7 96 82 -0.9 60 39 0.1 79 24 0.6 96 20 0.9 59 5 -0.1 80 5 -0.5 95 14 0.7 59 22 1.0 80 12 0.3 95 11 0.3 58 10 -0.4 73 58 -0.3 90 2 0.3 58 7 0.7 73 24 0.7 90 87 1.0 57 47 -0.2 74 52 0.2 89 39 0.6 57 22 -0.8 74 50 0.4 89 25 0.2 64 1 -0.3 75 54 0.9 92 38 0.9 64 50 -0.3 75 40 0.1 92 25 -0.8 63 24 0.2 76 60 -0.5 91 56 -0.8 63 12 0.7 76 17 0.4 91 16 -0.1 62 2 0.1 86 10 0.5 99 90 0.7 62 31 -0.5 86 1 -0.6 99 72 0.7 61 7 0.2 85 66 -0.1 100 4 0.7 69 5 0.5 85 51 -1.0 97 16 0.2 70 49 -0.6 88 68 0.7 97 67 0.8 70 69 0.6 88 80 0.6 98 82 0.3 71 16 0.1 87 2 -0.1 98 51 0.5 72 1 -0.7 87 86 0.8 72 30 0.8 82 11 -0.2 65 72 -0.4 82 5 0.3 65 39 -0.3 81 46 -1.0 66 24 0.8 81 47 -0.6 87 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 Fig. 5.4 100 要素スケールフリーネットワークモデルの次数分布 88 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 第4節 第1項 100 要素ネットワークモデルの推定 100 要素ランダムネットワークモデル推定結果 Fig. 5.1 の 100 要素ランダムネットワークモデル推定する際に用いたパラメータを Table 5.3 に示した。パラメータ p については、SR と k のパラメータ値に応じて、式(2.2) を満たすような値とした。RA 法で適用されているハイブリッド法に用いられたパラメ ータを Table 5.4 および Table 5.5 に示した。推定の結果、k=6 としたときが最も TP 率 が高く、真のネットワークの相互作用の内、約 25%の相互作用を推定することができた (Table 5.6)。k=3 としたときが最も TP 率が低く、k=4-5 としたときの TP 率にはほと んど差がなかった。30 要素ネットワークモデルの推定(Table 3.7)に比べ、この推定で は多くの FP interaction が推定された。TP interaction の数は k=6 のときが最も多かった が、FP interaction の数は k=5 が最も多く、次いで k=4 のときであった。 第2節 100 要素スケールフリーネットワークモデル推定結果 推定で用いられた RA 法やハイブリッド法の各設定パラメータは、 Table 5.3、 Table 5.4、 Table 5.5 の通りである。これらの条件で推定した結果、k=6 としたときが最も TP 率 が高く、SR=15 の段階で約 35%の相互作用を推定することが出来た(Table 5.7) 。k=3 としたときが最も TP 率が低く、k の値を大きくするほど TP 率が上昇した。また、100 要素ランダムネットワークモデルの推定と同様に、多くの FP interaction が推定された。 89 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 Table 5.3 検証のために用いた RA 法のパラメータ 設定値 Table 5.4 ハイブリッド法 k 3-6 l 30 m 30 SR 3×6 ハイブリッド法の実数値 GA(AREX+JGG)設定パラメータ 実数値 GA(AREX+JGG)の 要素数 推定対象 (k) パラメータ数 設定パラメータ 集団サイズ 子個体数 親個体数 3 18 540 72 19 4 28 840 112 29 5 40 1200 160 41 6 54 1640 216 55 90 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 推定のための最適化条件 Table 5.5 項目 設定 最適化対象パラメータ αi、βi、gij、hij (i ≠ j) タイムコースデータのサンプリングポイント数 71 Rosenbrock 法 ハイブリッド法 𝑥* 1.618034 𝑦* 0.618034 刻み幅 h* 0.1 終了条件 0.01 修正 Powell 法 開始平均相対二乗誤差(%) < 20.0 cut off 0.05 最大世代数 100,000 目標平均相対二乗誤差(%) < 10.0 *第 2 章第 4 節第 4 項参照 91 92 3 6 75 150 12 15 250 6 100 200 3 50 9 15 300 150 12 6 120 240 3 60 9 15 375 180 12 300 9 15 500 225 12 400 6 200 9 3 100 300 SR p 48 43 30 22 11 38 30 25 18 13 38 34 26 17 10 20 16 12 10 3 TP interaction 116 99 75 56 25 133 109 83 56 28 128 103 84 53 27 113 95 70 51 29 FP interaction 147 152 165 173 184 157 165 170 177 182 157 161 169 178 185 175 179 183 185 192 FN interaction 24.62 22.05 15.38 11.28 5.64 19.49 15.38 12.82 9.23 6.67 19.49 17.44 13.33 8.72 5.13 10.26 8.21 6.15 5.13 1.54 TP 率(%) 100 要素ランダムネットワークモデル推定結果 抽出された相互作用数(/195) *CPU : Xeon E5540 2.53GHz, Memory : 8.0GB 6 5 4 3 (k) 要素数 Table 5.6 56.28 56.30 55.45 54.84 53.78 24.58 24.57 24.64 24.14 23.45 8.50 8.54 8.64 8.66 8.41 2.32 2.34 2.21 2.04 2.05 [h]*/試行 平均 CPU time 14069.45 11259.27 8317.82 5484.32 2689.13 7373.12 5897.80 4435.26 2896.37 1407.25 3186.09 2563.16 1944.55 1298.52 630.62 1160.73 935.90 663.35 408.04 205.33 [h]*/ p 試行 CPU time 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 93 3 6 75 150 12 15 250 6 100 200 3 50 9 15 300 150 12 6 120 240 3 60 9 15 375 180 12 300 9 15 500 225 12 400 6 200 9 3 100 300 SR p 65 52 39 23 10 58 49 38 25 13 45 37 27 23 8 37 34 27 19 10 TP interaction 96 79 61 49 27 83 66 49 35 13 82 60 48 38 17 64 55 36 26 14 FP interaction 122 135 148 164 177 129 138 149 162 174 142 150 160 164 179 150 153 160 168 177 FN interaction 抽出された相互作用数(/187) 34.76 27.81 20.86 12.30 5.35 31.02 26.20 20.32 13.37 6.95 24.06 19.79 14.44 12.30 4.28 19.79 18.18 14.44 10.16 5.35 TP 率(%) 42.28 42.33 42.22 41.80 41.14 17.91 17.97 17.90 17.78 17.62 5.99 6.10 6.05 6.03 6.30 1.67 1.65 1.67 1.70 1.76 [h]*/試行 平均 CPU time 100 要素スケールフリーネットワークモデル推定結果 *CPU : Xeon E5540 2.53GHz, Memory : 8.0GB 6 5 4 3 (k) 要素数 Table 5.7 10568.97 8465.57 6332.39 4179.79 2056.88 5373.87 4313.43 3222.44 2133.19 1056.96 2247.51 1829.00 1360.40 905.02 472.51 834.50 660.76 499.50 339.91 176.39 [h]*/ p 試行 CPU time 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 第5節 考察 本章では、ランダムネットワークとスケールフリーネットワークに焦点を当て、特徴 的な構造を持つネットワークにおいても RA 法が有効であるか検証した。検証の結果、 ランダムネットワークおよびスケールフリーネットワークどちらのネットワークにお いても TP interaction を推定することができ、RA 法が両方のネットワークで推定能力を 有することがわかった。本検証で推定できた TP interaction は真のネットワークの半分 にも満たなかった。しかしながら、SR の値を増加させより多くの組み合わせのサブネ ットワークを推定することで TP 率は増加していくと考えられる。例えば、100 要素ス ケールフリーネットワークモデルの k=6 のときの TP 率の変化に着目すると、TP 率は SR の増加に対して単調増加している(Fig. 5.5) 。100 要素から異なる 6 要素を選ぶ組合 せは 1,192,052,400 通りであり、SR=15 の段階で 250 通りしか検証していない。つまり、 僅か 0.00002%の組合せしか検証していないにもかかわらず、25~35%の相互作用を推 定することができた。SR を増やしていくことで TP 率は今後も増加を続けるだろう。 また、100 要素ランダムネットワークモデルと 100 要素スケールフリーネットワーク モデルの推定では、多くの FP interaction が推定された。FP interaction は実際のネットワ ークには含まれないが、推定により導出された相互作用である。各検証で推定された FP interaction は、迂回経路を有していることが考えられる。ここで、迂回経路を有する FP interaction とは、直接的には接続していないが他の要素を介して間接的に接続されて いるような相互作用を示す。例えば、100 要素ランダムネットワークの推定において k =6、SR=15 としたとき、推定された 116 の FP interaction(Table 5.6)の内、1~3 要素 を経由するような迂回経路の数は 111 であった。本検証で、このように多くの FP interaction が推定された理由の一つは用いたネットワークモデルの頂点間距離にあると 考えられる。頂点間距離とは、異なる 2 ノード間を行き来するための最小のステップ数 である(町田と今野, 2008) 。つまり、あるノードからあるノードへと移動するための最 94 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 短経路を頂点間距離と呼ぶ。たとえば、Fig. 5.5 のような 5 つのノードから成るネット ワークがあるとする。このネットワークに含まれるノードを水色の丸で、エッジを灰色 の実線で示した。このネットワークで A から E に移動する経路として赤、青、黄の 3 つの経路を考える。赤い経路では、{A→D→E}と移動しており、このときステップ数 は 2 となる。同様にして、青い経路では{A→B→C→E}と移動しているので 3 ステッ プ、黄の経路では{A→B→D→C→E}と移動しているので 4 ステップである。このよ うに、A から E まで移動するすべての経路を調べ、最小のステップ数で移動できる経路 のステップ数が頂点間距離となる。Fig. 5.6 の例では、最短経路は赤の経路を辿る時の 経路であり、頂点間距離は 2 となる。100 要素ランダムネットワークモデルと 100 要素 スケールフリーネットワークモデルの頂点間距離の分布を Fig. 5.7 に示した。各要素は、 その他の 99 の要素に対してそれぞれ最短経路を持つため、最短経路数の総計は 9,900 となる。従って、Fig. 5.7 は 9,900 の最短経路の頂点間距離の度数分布である。両ネット ワークモデル共、頂点間距離は 3 と 4 に集中しており、平均頂点間距離はランダムネッ トワークモデルが 3.6 であり、スケールフリーネットワークが 2.2 である。ネットワー クモデルの要素数が 100 であるのに対し、各要素間の距離が 2~4 であるため、要素間 の距離が全体的に短いといえる。一方、第 3 章や第 4 章で用いてきた 30 要素ネットワ ークモデル(Fig. 3.2)の頂点間距離の分布を Fig. 5.8 に示した。平均頂点間距離は 4.8 であり、100 要素ネットワークモデルよりもネットワーク構成要素数が少ないにもかか わらず各要素間の距離は長い。頂点間距離が短いと、直接的な相互作用がない要素から も間接的な影響を受けやすいと考えられる。なぜならば、平均頂点間距離が 2.2 である ならば、わずか 2 もしくは 3 ステップですべての要素へ移動できるため、相互作用も同 様に伝搬しやいからである。したがって、本検証で用いた 100 要素ネットワークモデル は迂回経路が多数あると考えられ、それが原因となり多くの FP interaction が推定された 可能性が高い。 FP interaction は真のネットワークには含まれてはいない相互作用ではあるものの、完 全なる不正解であるとはいえない。実問題においては真のネットワーク構造は未知であ 95 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 るので、TP interaction と FP interaction の区別を付けることはできず、両者を相互作用候 補として扱う必要がある。FP interaction が迂回経路を有する場合、間接的な相互作用で あるとはいえ、2 項間のみに着目すれば制御関係にあり、有用な情報であるといえる。 40 35 TP率 30 25 20 15 10 5 0 0 3 6 9 12 15 18 SR Fig. 5.5 100 要素スケールフリーネットワークモデル推定における TP 率の変化(k=6) B A C D Fig. 5.6 E 頂点間距離の例 96 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 3500 3000 最短経路数 2500 2000 1500 1000 500 0 1 2 3 4 5 6 7 8 9 頂点間距離 (a) 100 要素ランダムネットワークモデルの頂点間距離の分布 4500 4000 3500 最短経路数 3000 2500 2000 1500 1000 500 0 1 2 3 4 5 頂点間距離 (b) 100 要素スケールフリーネットワークモデルの頂点間距離の分布 Fig. 5.7 100 要素ネットワークモデルの頂点間距離の分布 97 第 5 章 異なるネットワークモデルに対する Random selection method の有効性 160 140 120 最短経路数 100 80 60 40 20 0 1 2 3 4 5 6 7 8 9 10 11 12 頂点間距離 Fig. 5.8 30 要素ネットワークモデルの頂点間距離の分布 98 第6章 総括 第 6 章 総括 本研究では、タイムコースデータから要素間相互作用ネットワークを推定する手法を 開発し、その有効性を検証した。 先行研究において、生体分子の挙動を、一般質量微分方程式(GMA)を簡略化した 抽象モデルである S-system で表記し、S-system に含まれるシステムパラメータを実数値 GA を用いて推定する手法が開発されてきた(Ueda, 2001; Ueda et al., 2002a; Ueda et al., 2002b) 。実数値 GA を用いて開発された手法では、代表的な実数値 GA のアルゴリズム である UNDX+MGG や SPX+MGG が用いられていたが、この手法では推定に多大な時 間を要し、推定可能なパラメータ数は 20 程度に限られていた。そこで、新たな実数値 GA として AREX+JGG を取り入れた AREX+JGG 法とハイブリッド法(AREX+JGG+修 正パウエル法)を開発し、第 2 章で、開発した手法が UNDX+MGG 法と比較して計算 時間や推定規模の面で優れていることを示した(Komori et al., 2012; Komori et al., 2013) 。 従来法では 3 要素程度の相互作用ネットワークしか推定できなかった S-system と実数 値 GA に基づく手法を、開発した AREX+JGG 法およびハイブリッド法により 17 要素程 度の相互作用ネットワークを推定可能な手法とすることに成功した(第 2 章第 3 節第 3 項)。しかしながら、実用化のためにはより大規模な相互作用ネットワークに適用でき る必要がある。 そこで、本研究では、推定対象となる相互作用ネットワークを構成する要素から、ラ ンダムに要素を選択し、選ばれた要素で構成されたサブネットワークを推定することで 推定対象の相互作用ネットワークを再構築する Random selection method(RA 法)を開 発し、その有効性を検証した。具体的にはランダムに推定する要素を選択し、選択した 要素間の推定にハイブリッド法を適用して、その後、推定した要素間の共通構造を重ね あわせるという手法である(Komori et al., 2013) 。一般的に、大規模なネットワークを いくつかのクラスタに分けるためには、推定に用いるタイムコースデータの形状や文献 やデータベースの情報などを用いた操作が必要となる場合が多い。しかしながら、RA 99 第6章 総括 法では、必要とされる前処理はランダムに要素を選択するだけであるため、簡便に推定 を開始することができ、人為的な操作も必要とないため恣意性の低い結果を得ることが できる。さらに、実験前にサブネットワークの要素の組み合わせを決めておけば、その 組み合わせの実験のみ行えば良いので、実験にかかる費用や時間、手間を削減すること ができるというメリットもある。 本論文では、まず、第 3 章で RA 法が相互作用ネットワーク推定において有効な手法 であるか 30 要素ネットワークモデル(Fig. 3.2)を用いて検証した。推定された相互作 用は、①推定された相互作用の内、真のネットワークにも存在する相互作用である true-positive interaction(TP interaction)、②推定された相互作用の内、真のネットワーク に存在しない相互作用である false-positive interaction(FP interaction)、③真のネットワ ークに存在するものの推定できなかった相互作用である false-negative interaction(FN interaction)の 3 種で定義付けることとした。まず、30 要素ネットワークモデルを推定 するための適切な RA 法においてサブネットワークの要素数を表すパラメータ k を決定 するために 3 つの推定条件で推定を行ったところ、k=5 のときに、もっとも多くの TP interaction を推定することに成功した(Table 3.3) 。したがって、k=5 として 30 要素ネ ットワークモデルの推定を試みた。RA 法において各要素を均等に選択するためのパラ メータである SR(selection rate)を増加させ、SR=16 のとき、真のネットワークに含ま れている相互作用の内約 93%を推定することに成功した(Table 3.7) 。このとき、要素 のランダムな選択は 96 回(p)行われた。30 要素から異なる 5 要素を選ぶ組合せの総 数は、30C5=142,506 である。すなわち、142,506 ある組合せの内、たかだか 96 の組み合 わせで、ネットワークモデルの約 93%を推定することに成功したのである。今後 SR や p の値を増加させることで更に多くの相互作用を推定することができると考えられる。 この検証から、タイムコースデータから相互作用ネットワークを推定する手法として RA 法は有効であることが示された。AREX+JGG 法やハイブリッド法をより大規模なシ ステムへと適用するために、RA 法は効果的なアルゴリズムであるといえる。 次に、第 4 章では、これまで開発されてきた実数値 GA を用いた 4 つの手法 100 第6章 総括 (UNDX+MGG 法、SPX+MGG 法、AREX+JGG 法、ハイブリッド法)および RA 法の ノイズ耐性を検証した。相互作用ネットワーク推定手法を実問題へと適用する際には、 入力値として用いる実験データには、誤差が含まれている可能性が非常に高い。誤差が 生じる原因としては、測定機器による誤差、生物学的ノイズ等様々ある。本研究では、 各手法のノイズ耐性を検証するために、ブラウン運動を記述する確率微分方程式である Langevin 方程式を S-system モデルに取り入れノイズを発生させることで、タイムコー スデータを作成し、このタイムコースデータに誤差を含む入力データとした。実験の際 に、サンプリングポイントごとにあらゆる要因で発生する誤差を、Langevin 方程式によ り計算時刻ごとにノイズを発生させることで表現した。AREX+JGG 法およびハイブリ ッド法のノイズ耐性の検証には、4 要素ネットワークモデル(Fig. 4.1(a))と 7 要素ネッ ト ワーク モデル (Fig. 4.1(b)) を用い た。こ のと き、先 行研究 で開発 されて い た UNDX+MGG 法および SPX+MGG 法においても同様のネットワークモデルを用いてノ イズ耐性を検証し、 各手法の有効性を比較した。4 要素ネットワークモデルの推定では、 すべての手法において全相互作用を推定することに成功した(Table 4.4) 。しかしなが ら、UNDX+MGG 法および AREX+JGG 法では、FP interaction が一つ推定された。更に、 UNDX+MGG 法では、推定成功率が 50%を下回っていた。推定に要した CPU time(CPU : Xeon E5540 2.53GHz, Memory : 8.0GB)は、ハイブリッド法がもっとも少なかった。7 要 素ネットワークモデルを用いた実験では、UNDX+MGG 法および SPX+MGG 法は、定 められた最大世代数の内では、推定に一度も成功しなかった。AREX+JGG 法とハイブ リッド法の成功試行回数は同等であるが、ハイブリッド法がすべての相互作用を推定し た一方で、AREX+JGG 法では FN interaction が一つ存在した。これらの検証において、 ハイブリッド法はどちらのネットワークモデルにおいてもすべての相互作用を推定す ることができ、さらに計算に要する CPU time も短かった。従って、4 つの手法の中で は、ハイブリッド法が最もノイズ耐性の高い手法であることが示唆された。 UNDX+MGG 法では、推定に成功した試行回数や、推定された相互作用、CPU time の いずれにおいても、他の手法の方が優れており、ノイズ耐性は低いと思われた。また、 101 第6章 総括 第 4 章では 4 要素ネットワークモデルを用い、AREX+JGG 法およびハイブリッド法の 許容できるノイズの大きさについても検証した。タイムコースデータ 1 点あたりの相対 誤差(平均相対誤差)にバリエーションを持たせ推定したところ、平均相対誤差 30% 前後に許容できる最大のノイズがあることがわかった(Table. 4.5) 。RA 法のノイズ耐性 の検証には、30 要素ネットワークモデル(Fig. 3.2)を用いた。ノイズを含むタイムコ ースデータを用いた場合(case3-01)とノイズを含まないタイムコースデータを用いた 場合(case3-02)で、真のネットワークに存在する相互作用の内、推定された相互作用 の割合である TP 率(式(2.2))を比較した。どちらのケーススタディーにおいても、RA 法のパラメータ k=5 としたときがもっとも TP 率が高く、どちらも約 80%の相互作用 を推定することに成功した。このことから、ノイズを含むタイムコースデータを用いた 場合でも、推定能力が低下することはなく、RA 法はノイズを含むタイムコースデータ に対しても有効であることが示された。case3-01 では、k=2~5 の間で、SR の増加に対 する TP 率の増加に差が生じた。さらに、k=2~4、SR=8 のときに着目すると、p 試行 の推定に要した CPU time はいずれの場合も 200 時間前後であるにもかかわらず TP 率に は大きな差がある。この推定では、ノイズを含むデータを用いたためタイムコースデー タの形状による拘束条件が弱いと考えられる。サブネットワークの要素数が少ないとき に推定時間を要したのは、拘束条件が弱い中で入力データを再現できるような S-system パラメータの組合せの探索が困難になり推定に成功した回数が少なかったためである と考えられる。したがって、RA 法を適用する際には、パラメータ k の値を最適な値に 設定することで効率の良い相互作用ネットワークの推定に繋がる。今後、推定対象とな るネットワーク構成要素数やネットワークの構造等による適切なパラメータ k の探索 をすることが必要であると思われる。 複雑ネットワークの構造には、スモールワールド性やスケールフリー性などの特徴を もつものがあり、これまで様々なモデルが提唱されている。第 5 章では、RA 法が特徴 的な構造を持つ相互作用ネットワークにおいても有効であるか検証するため、伝統的な 複雑ネットワークモデルであるランダムネットワークと、近年多くの複雑ネットワーク 102 第6章 総括 がその特徴を持つと報告されているスケールフリーネットワークの二つのネットワー クモデルに対して RA 法を適用した。どちらも構成要素数を 100 とし、ネットワークジ ェネレータを用いてネットワークモデルを作成した(Fig. 5.1、Fig. 5.3) 。どちらの 100 要素ネットワークモデルにおいても相互作用を推定することができ、RA 法は様々なネ ットワーク構造においても有効な手法であることが示唆された。現在のところ、TP 率 は 50%にも達していなかったが、SR の増加に対して TP 率は単調増加しているため、SR を増やすことで TP 率は増加すると思われる。 本研究で考案した RA 法により、AREX+JGG 法やハイブリッド法をより強力な相互 作用ネットワーク推定手法とすることに成功した。RA 法は、AREX+JGG 法やハイブリ ッド法単独では推定することのできなかった規模の相互作用ネットワークの推定に効 果的であることが示された。さらに、ノイズ耐性を有し、特徴的な構造を持つような相 互作用ネットワークに対しても有効であるため、実用的な手法であると考えられる。今 後、トランスクリプトームやプロテオームなどの定量的なタイムコースデータを用い、 実問題への有効性を検証することが望まれる。 第1節 推定時間 RA 法は、大規模な相互作用ネットワークを推定する手法として有効な手法であるこ とが示されたが、推定に膨大な計算時間を要している。100 要素ネットワークモデルの 推定では、サブネットワークの大きさであるパラメータ k の値が大きいほど、SR に対 する TP 率の増加率も大きかった(Table 5.6、Table 5.7)。しかしながら、サブネットワ ークの増加と共に、推定に要する時間も増加している。例えば、100 要素スケールフリ ーネットワークの推定では、k=6、SR=15 としたときに、推定に約 10,569 時間掛かっ ている(Table 5.7)。30 要素ネットワークモデルの推定において、SR の増加に対して TP 率は約 80%になるまで線形に増加していた(Fig. 3.4)。従って、100 要素スケールフ リーネットワークの推定においても、SR の増加に対して TP 率が 80%となるまで線形に 103 第6章 総括 増加すると仮定すると、SR=35 のとき TP 率は 80%を超える。100 要素スケールフリー ネットワークの相互作用の 80%の推定を目指し SR=35 となるまで推定するには、約 24,660 時間(約 1,028 日)という膨大な時間が必要となってしまい、実用的でないよう に思える。しかしながら、本研究で示した CPU time は、すべて 1 コアを用いた場合の 計算時間である。一般的に、大規模な相互作用ネットワークの推定では、大量の CPU コアを搭載したスーパーコンピュータが用いられる。例えば、九州大学のスーパーコン ピュータ(FUJITSU PRIMEHPC FX10)には 12,288 コアが搭載されている(平成 28 年 1 月現在) 。簡単のために、このコアをすべて使用できるとすると、100 要素スケールフ リーネットワークにおいて SR=35 となるまで推定するためには、約 2 時間あれば良い。 1,000 コア程度のコア数であっても、1 日で推定を終えることができる。以上のことか ら、複数のコアを用いることで RA 法でも実用的な時間で大規模な相互作用ネットワー クの推定が可能であることが示唆された。 第2節 TP interaction と FP interaction 第 5 章第 5 節で述べたように、100 要素ネットワークモデルの推定では多くの FP interaction が推定された。実問題においては、真のネットワーク構造が不明であるため、 本研究における検証のように RA 法により推定された相互作用を TP interaction と FP interaction に区別することはできない。TP interaction のみを推定できる手法であること が望ましいが、FP interaction にも有益な情報が含まれているかもしれない。 ここで、推定された FP interaction の中に、迂回経路を有するような無意味ではない相 互作用が含まれていると可能性を考えた。例えば、Fig. 6.1(a)のような構造を持つネッ トワークを推定する場合に、Fig. 6.1(b)の青の点線で示したような要素 1 から要素 3 へ の相互作用が推定されたとする。このとき、要素 1 から要素 3 への相互作用は FP interaction となる。要素 3 は要素 1 とは直接的には相互作用していないが、要素 2 を介 して間接的な影響を受けていると考えられる。つまり、本来は迂回して制御している相 104 第6章 (a) 1 2 3 (b) 1 2 3 Fig. 6.1 迂回経路の例 (a) 真のネットワーク構造 (b) 迂回経路を有する FP interaction 青の破線は FP interaction を表す 105 総括 第6章 総括 互作用が FP interaction として推定されていることとなる。100 要素ネットワークモデル の推定においても、FP interaction の中にこのような迂回経路を有する相互作用が含まれ ているか検証した。100 要素ランダムネットワークモデルの k=6、SR=15 のときに着 目すると、最終的に 116 の FP interaction が推定されている。この 116 の FP interaction の内、Fig. 6.2(a)のように、1 要素もしくは 2 要素を経由した迂回経路を有している相互 作用は 98 であった。さらに、Fig. 6.2(b)のように、3 要素を経由した迂回経路を有して いる相互作用を含めると、推定された FP interaction の 95.7%に当たる 111 の FP interaction が迂回経路を有していた。このような間接的な相互作用は、2 項関係に着目すれば、2 つの要素は制御関係にあるといえる。つまり、100 要素ランダムネットワークの推定で は TP 率は約 25%であったものの、116 の FP interaction の内、111 の FP interaction が制 御関係を表しており、48 の TP interaction と合わせて 159 の制御関係を示唆していたこ ととなる。このように、FP interaction は迂回経路の存在を示唆していることが高いこと (a) 1 2 1 2 3 3 4 (b) 1 2 1 2 1 3 2 3 3 4 4 5 Fig. 6.2 迂回経路 (a) 1 要素・2 要素を介した迂回経路 (b) 1 要素・2 要素・3 要素を介した迂回経路 青の破線は FP interaction を表す 106 第6章 総括 が示され、RA 法で推定された相互作用のほとんどが有用な情報であるといえる。しか しながら、実用化のためには、できる限り FP interaction の数を減らしたり、TP interaction と FP interaction を区別したりするための評価関数を導入するなどの工夫を RA 法へ取り 入れることが今後望まれる。 第3節 今後の展望 RA 法は、これまでの S-system モデルと実数値 GA を基盤とした相互作用ネットワー ク推定手法では推定が困難であった規模の相互作用ネットワークに対して効果的な手 法である。しかしながら、相互作用ネットワークの規模が大きくなればなるほど、SR を増やす必要があり、推定に要する計算時間の増大や、推定結果の解析に手間を要する ことが課題となる。そこで、RA 法による相互作用ネットワークの推定や解析を一つの システムとして実行できるようなソフトウェアを開発することで、より実用的になる。 また、RA 法では複数の CPU を用い並列的に相互作用ネットワークを推定することで、 計算時間が大幅に短縮されると考えらえる。タイムコースデータを入力するだけで、相 互作用ネットワークを並列に推定し、推定された相互作用が出力されるようなシステム を開発することで、より広範に実用化できるような手法となるであろう。 また、S-system モデルは概念モデルであるため、適用範囲が限られておらず、遺伝子 やタンパク質、代謝物質など様々な生体分子の動的挙動を記述できる。したがって、ゲ ノムやトランスクリプトーム、プロテオーム、メタボロームなど多階層オミクス(トラ ンスオミクス)への適用も期待される。従来、オミクスデータを用いた研究は、一つの 階層にのみ着目したものや、特定の分子のみの階層間相互作用に着目したものがほとん どであった。近年、複数の階層の網羅的な多階層相互作用ネットワークを推定しようと いう試みが成されている(Yugi et al., 2014)。本論文で開発した RA 法を多階層オミクス へと応用することで疾患原因のより高度なシステム的理解や薬剤の標的分子探索につ ながることが期待される。 107 参考文献 参考文献 Akimoto, Y., Nagata, Y., Sakuma, J., Ono, I., & Kobayashi, S. (2009). Proposal and Evaluation of Adaptive Real-coded Crossover AREX. Transactions of the Japanese Society for Artificial Intelligence, 24, 446–458. Akimoto, Y., Nagata, Y., Sakuma, J., Ono, I., & Kobayashi, S. (2010). Analysis of The Behavior of MGG and JGG As A Selection Model for Real-coded Genetic Algorithms. Transactions of the Japanese Society for Artificial Intelligence, 25(2), 281–289. Akimoto, Y., Rie, H., Sakuma, J., Ono, I., & Kobayashi, S. (2007). Generation Alternation Model for Real-coded GA Using Multi-parent: Proposal and Evaluation of Just Generation Gap (JGG). SICE Symposium on Decentralized Autonomous Systems 19, 341–346. Akimoto, Y., Sakuma, J., Ono, I., & Kobayashi, S. (2009). Adaptation of expansion rate for real-coded crossovers. Proceedings of the 11th Annual Conference on Genetic and Evolutionary Computation - GECCO ’09, 739–746. Albert, R., Jeong, H., & Barabasi, A. L. (1999). Internet: Diameter of the World-Wide Web. Nature, 401(6749), 130–131. Baker, J. E. (1987). Reducing bias and inefficiency in the selection algorithm. In Proceedings of the 2nd International Conference on Genetic Algorithms, 14–21. 108 参考文献 Bansal, M., Belcastro, V., Ambesi-Impiombato, A., & di Bernardo, D. (2007). How to infer gene networks from expression profiles. Molecular Systems Biology, 3, 78. Barabási, A.-L., & Albert, R. (1999). Emergence of scaling in random networks. Science, 286, 509–512. Csárdi, G., & Nepusz, T. (2006). The igraph software package for complex network research. InterJournal Complex Systems, 1695. Davis, L. D. (1990). Handbook Of Genetic Algorithms. New York: Van Nostrand Reinhold. Erdös, P., & Rényi, A. (1959). On random graphs. Publicationes Mathematicae, 6, 290–297. Eshelman, L., & Schaffer, J. (1993). Real-coded genetic algorithms and interval-schemata. In Foundation of Genetic Algorithms 2, 187 – 202. Faloutsos, M., Faloutsos, P., & Faloutsos, C. (1999). On power-law relationships of the Internet topology. ACM SIGCOMM Computer Communication Review, 29(4), 251–262. Fox Keller, E. (2005). Revisiting “scale-free” networks. BioEssays, 27(10), 1060–1068. Gillespie, D. T. (2000). The chemical Langevin equation. The Journal of Chemical Physics, 113(1), 297. Goldberg, D. E. (1989). Genetic Algorithms in Search, Optimization, and Machine Learning. Addison-Wesley Professional. 109 参考文献 Hansen, N. (2006). The CMA Evolution Strategy : A Comparing Review. Towards a New Evolutionary Computation, 102, 75–102. Higuchi, T., Tsutsui, S., & Yamamura, M. (2001). Simplex Crossover for Real-coded Genetic Algolithms. Transactions of the Japanese Society for Artificial Intelligence, 16(1), 147– 155. Hu, M.-B., Wang, W.-X., Jiang, R., Wu, Q.-S., Wang, B.-H., & Wu, Y.-H. (2006). A Unified Framework for the Pareto Law and Matthew Effect using Scale-Free Networks. European Physical Journal B, 53(2), 5. Ideker, T., Galitski, T., & Hood, L. (2003). A NEW APPROACH TO DECODING LIFE: Systems Biology. Kitano, H. (2002a). Computational systems biology. Nature, 420(6912), 206–210. Kitano, H. (2002b). Systems Biology: A Brief Overview. Science, 295(5560), 1662–1664. Kobayashi, S. (2009). The Frontiers of Real-coded Genetic Algorithms GA. Transactions of the Japanese Society for Artificial Intelligence, 24(1), 147–162. Komatsu, T., & Namatame, A. (2012). Relation between network structure and cascade phenomena, 4(2), 33–40. 110 参考文献 Komori, A., Maki, Y., Nakatsui, M., Ono, I., & Okamoto, M. (2012). Efficient Numerical Optimization Algorithm Based on New Real-Coded Genetic Algorithm, AREX + JGG, and Application to the Inverse Problem in Systems Biology. Applied Mathematics, 03(30), 1463–1470. Komori, A., Maki, Y., Ono, I., & Okamoto, M. (2013). How to Infer the Interactive Large Scale Regulatory Network in “Omic” Studies. Procedia Computer Science, 23, 44–52. Komori, A., Maki, Y., Ono, I., & Okamoto, M. (2015). Investigating noise tolerance in an efficient engine for inferring biological regulatory networks. Journal of Bioinformatics and Computational Biology, 13(3), 1541006. MacQueen, J. (1967). Some methods for classification and analysis of multivariate observations. Proceedings of the Fifth Berkeley Symposium on Mathematical Statistics and Probability, 1(14). Matsumoto, M., & Nishimura, T. (1998). ACM Transactions on Modeling and Computer Simulation (TOMACS) 8.1, 8(1), 3–30. Medina, A., Lakhina, A., Matta, I., & Byers, J. (2001). BRITE: an approach to universal topology generation. In MASCOTS 2001, Proceedings Ninth International Symposium on Modeling, Analysis and Simulation of Computer and Telecommunication Systems, 346– 353. 111 参考文献 Nakatsui, M., Ueda, T., Maki, Y., Ono, I., & Okamoto, M. (2008). Method for inferring and extracting reliable genetic interactions from time-series profile of gene expression. Mathematical Biosciences, 215(1), 105–114. Okamoto, M., Nonaka, T., Ochiai, S., & Tominaga, D. (1998). Nonlinear numerical optimization with use of a hybrid genetic algorithm incorporating the modified Powell method. Applied Mathematics and Computation, 91(1), 63–72. Okamoto, M., Tominaga, D., Koga, N., Maki, Y., Watanabe, S., & Eguchi, Y. (2000). Towards a System for the Inference of Large Scale Genetic Networks: Efficient Numerical Optimization Algorithm Based on Genetic Algorithm for Inverse Problem. IPSJ Symposium Series, 2000(5), 29–37. Ono, I., Kobayashi, S., & Yoshida, K. (1998). Global and multi-objective optimization for lens design by real-coded genetic algorithms. Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, 3482, 110–121. Ostermeier, A., Gawelczyk, A., & Hansen, N. (1994). Step-Size Adaption Based on Non-Local Use of Selection Information. PPSN III Proceedings of the International Conference on Evolutionary Computation. The Third Conference on Parallel Problem Solving from Nature: Parallel Problem Solving from Nature, 866, 189–198. Powell, M. J. D. (1964). An efficient method for finding the minimum of a function of several variables without calculating derivatives. The Computer Journal, 7(2), 155–162. 112 参考文献 Rosenbrock, H. H. (1960). An Automatic Method for Finding the Greatest or Least Value of a Function. The Computer Journal, 3(3), 175–184. Satoh, H. (1996). Minimal Generation Gap Model for GAs Considering Both Exploration and Exploitation. Proceedings of the 4th International Conference on Soft Computing, Neural Nets and Soft Computing, 494–497. Savageau, M. A. (1976). Biochemical Systems Analysis: A Study of Function and Design in Molecular Biology. MA: Addison-Wesley Educational Publishers Inc. Tsutsui, S., Yamamura, M., & Higuchi, T. (1999). Multi-parent Recombination with Simplex Crossover in Real Coded Genetic Algorithms. Proc. of the Genetic and Evolutionary Computation Conference, 1, 657–664. Ueda, T. (2001). Efficient Numerical Optimization Technique Based on Real-Coded Genetic Algorithm fo Inverse Problems. Proceedings of 7th. Int. Symp. on Artificial Life and Robotics (AROB 7th ’02), 453, 451–453. Ueda, T., Koga, N., Ono, I., & Okamoto, M. (2002a). Application of Numerical Optimization Technique Based on Real-coded Genetic Algorithm to Inverse Problem in Biochemical Systems. In GECCO ’02 Proceedings of the Genetic and Evolutionary Computation Conference (p. 701). Ueda, T., Koga, N., Ono, I., & Okamoto, M. (2002b). Development of Efficient Numerical Optimization Method Based on Real-coded Genetic Algorithm: Application to the Estimation of Large Number of Real-valued Parameters. Retrieved February 16, 2013. 113 参考文献 Wang, P., Fiaschi-Taesch, N. M., Vasavada, R. C., Scott, D. K., García-Ocaña, A. and Stewart, A. F. (2015). Diabetes mellitus--advances and challenges in human β-cell proliferation. Nature Reviews Endocrinology, 11(4), 201–212. Watts, D. J., & Strogatz, S. H. (1998). Collective dynamics of “small-world” networks. Nature, 393(6684), 440–442. Weng, G., Bhalla, U. S., & Iyengar, R. (1999). Complexity in biological signaling systems. Science, 284(5411), 92–96. Whitley, D. (1989). The GENITOR algorithm and selection pressure: why rank-based allocation of reproductive trials is best. In Proceedings of the third international conference on Genetic algorithms , 116–121. Wright, A. H. (1991). Genetic Algorithms for Real Parameter Optimization. Foundations of Genetic Algorithms, 205–218. Wuchty, S. (2001). Scale-free behavior in protein domain networks. Molecular Biology and Evolution, 18(9), 1694–1702. Yugi, K., Kubota, H., Toyoshima, Y., Noguchi, R., Kawata, K., Komori, Y., Kuroda, S. (2014). Reconstruction of insulin signal flow from phosphoproteome and metabolome data. Cell Reports, 8(4), 1171–1183. 114 参考文献 Zhu, H., Rao, R. S. P., Zeng, T., & Chen, L. (2012). Reconstructing dynamic gene regulatory networks from sample-based transcriptional data. Nucleic Acids Research, 40(21), 10657– 10667. 伊庭斉志. (1994). 遺伝的アルゴリズムの基礎. オーム社. 岡本正宏, & 小野功. (2003). 実数値GAのバイオ分野への応用 : 大規模遺伝子ネットワ ークの相互作用推定. 人工知能学会誌, 18(5), 502–509. 喜多一, 小野功, & 小林重信. (2000). 実数値GAのための正規分布交叉の多数の親を用 いた拡張法の提案. 計測自動制御学会論文集, 36(10), 875–883. 佐藤浩, 小野功, & 小林重信. (1997). 遺伝的アルゴリズムにおける世代交代モデルの提 案と評価. 人工知能学会誌, 12(5), 734–744. 坂和正敏, & 田中雅博. (1995). 遺伝的アルゴリズム. 朝倉書店. 四辻哲章,. (2010). 計算機シミュレーションのための確率分布乱数生成法. プレアデス 出版. 小野功, 佐藤浩, & 小林重信. (1999). 単峰性正規分布交叉UNDX を用いた実数値GA による関数最適化. 人工知能学会誌, 14(6), 1146–1155. 町田拓也, & 今野紀雄. (2008). 図解入門よくわかる複雑ネットワーク: シミュレーショ ンで見るモデルの性質科学セミナー. 秀和システム. 115 参考文献 富永大介, & 岡本正宏. (1999). 逆問題(inverse problem)解決のための遺伝的アルゴリズ ムを用いた多次元非線形数値最適化手法の開発. 化学工学論文集, 25(2), 220–225. 北野宏明. (1993). 遺伝的アルゴリズム. 産業図書. 116 謝辞 本研究を遂行するにあたり、多くの方々から貴重なご意見、ご指導を賜りましたこと を、心より感謝申し上げます。 九州大学大学院農学研究院生命機能科学部門岡本正宏主幹教授には、本研究実施の機 会を与えていただき深く御礼申し上げます。専門知識の欠けていた私に対し、熱心かつ 温かいご指導、ご鞭撻を賜りました。岡本正宏主幹教授の御高恩を賜り、本論文をまと めるに至ることができました。心より感謝申し上げます。 本論文を審査していただきました九州大学大学院農学研究院生命機能化学部門久原 哲教授ならびに九州大学生体防御医学研究所附属トランスオミクス医学研究センター 久保田浩行教授には、多くの貴重なご助言を賜りました。深く御礼申し上げます。 九州大学大学院農学研究院生命機能科学部門の花井泰三准教授には様々な視点から ご助言、ご意見を賜りました。深く感謝申し上げます。 九州大学大学院農学研究院生命機能科学部門の濱田浩幸助教には、日頃より幅広い知 識と多様な観点から多くのご助言、ご意見を賜りました。深く感謝申し上げます。 牧幸浩先生には、研究全般に渡り多くのご助言、ご意見を賜りました。深く感謝申し 上げます。 東京工業大学大学院総合理工学研究科の小野功准教授には、実数値遺伝的アルゴリズ ムに関する知識や技術について熱心なご指導を賜りました。深く感謝申し上げます。 最後に、常に温かい励ましとご助言を下さった九州大学合成生物学研究室の皆さんに 深く御礼を申し上げるとともに、皆様、そして合成生物学研究室の益々のご発展を祈念 し、謝辞とさせていただきます。 IV