(57)【要約】【課題】複数の酵素から構成された代謝ネットワークを動的モジュールと静的モジュールとに適切に分割する。【解決手段】解析対象の代謝ネットワークを構成する代謝物質について得られた濃度時系列データ又は酵素活性時系列データを入力するステップと、上記ステップで入力された代謝物質濃度時系列データから算出した酵素活性時系列データ、或いは上記ステップで入力された酵素活性時系列データを用いて適合関数を最適化することで、上記代謝ネットワークを構成する各酵素を静的モジュール又は動的モジュールに区分けするステップとを含む。【選択図】図２
2008040812.tif 000002
【特許請求の範囲】
【請求項１】
解析対象の代謝ネットワークを構成する代謝物質について得られた濃度時系列データ又は酵素活性時系列データを入力するステップと、
上記ステップで入力された代謝物質濃度時系列データから算出した酵素活性時系列データ、或いは上記ステップで入力された酵素活性時系列データを用いて適合関数を最適化することで、上記代謝ネットワークを構成する各酵素を静的モジュール又は動的モジュールに区分けするステップと
を含む代謝ネットワーク推定方法。
【請求項２】
解析対象の代謝ネットワークに含まれるとともに、当該代謝ネットワークの外側との間に位置する酵素をシステム境界酵素とし、システム境界酵素を除いた当該代謝ネットワークに含まれる酵素をシステム内部酵素とし、
当該システム境界酵素の酵素活性時系列データを代謝物質濃度時系列データから算出し、当該システム内部酵素の酵素活性時系列データを代謝物質濃度の勾配及びシステム境界酵素について算出された活性値を用いて下記式によって算出することを特徴とする請求項１記載の代謝ネットワーク推定方法。
【数１】
2008040812.tif 000012
（式中、Sinternalはシステム内部酵素の化学量論係数行列を示し、Smarginalはシステム境界酵素の化学量論係数行列を示し、diag(-1)は対角行列（行の数＝列の数＝代謝物質の数）を示し、vmarginal(t)はシステム境界酵素活性のベクトルを示し、vinternal(t)はシステム内部酵素活性のベクトルを示し、C'(t)は代謝物質濃度の速度ベクトルを示す。また、上記代謝物質濃度の勾配はC'(t)であり、酵素濃度時系列データの微分によって算出される値である。さらに、上記システム境界酵素について算出された活性値はvmarginal(t)である。）
【請求項３】
上記適合関数は下記式で表され、代謝ネットワークを動的モジュール及び静的モジュールに区分けするパターンについて下記式ｆ値を算出し、当該ｆ値を最小にするパターンを解として出力することを特徴とする請求項１記載の代謝ネットワーク推定方法。
【数２】
2008040812.tif 000013
（式中、Cdata,i,jはj番目のサンプリング点におけるi番目の代謝物質の濃度データを示し、Cestimated,i,jはj番目のサンプリング点におけるi番目の代謝物質の推定濃度データを示し、nmetaboliteは代謝物質の数を示し、nsampling pointはサンプリング点の数を示し、nenzymeはシステム内部酵素の数を示し、nstatic enzymeは静的モジュールに含まれる酵素の数を示し、wは重み係数を示す。）
【請求項４】
解析対象の代謝ネットワークに含まれるとともに、当該代謝ネットワークの外側との間に位置する酵素をシステム境界酵素とし、当該システム境界酵素は必ず動的モジュールに含まれるように上記パターンを定義することを特徴とする請求項３記載の代謝ネットワーク推定方法。
【請求項５】
解析対象の代謝ネットワークを構成する代謝物質について得られた濃度時系列データ又は酵素活性時系列データを入力するデータ入力処理手段と、
上記データ入力処理手段で入力された代謝物質濃度時系列データから算出した酵素活性時系列データ、或いは上記データ入力処理手段で入力された酵素活性時系列データを用いて適合関数を最適化することで、上記代謝ネットワークを構成する各酵素を静的モジュール又は動的モジュールに区分けする動的/静的モジュール区分け演算手段と
を含む代謝ネットワーク推定装置。
【請求項６】
上記動的/静的モジュール区分け演算手段では、解析対象の代謝ネットワークに含まれるとともに、当該代謝ネットワークの外側との間に位置する酵素をシステム境界酵素とし、システム境界酵素を除いた当該代謝ネットワークに含まれる酵素をシステム内部酵素とし、
当該システム境界酵素の酵素活性時系列データを代謝物質濃度時系列データから算出し、当該システム内部酵素の酵素活性時系列データを代謝物質濃度の勾配及びシステム境界酵素について算出された活性値を用いて下記式によって算出することを特徴とする請求項５記載の代謝ネットワーク推定装置。
【数３】
2008040812.tif 000014
（式中、Sinternalはシステム内部酵素の化学量論係数行列を示し、Smarginalはシステム境界酵素の化学量論係数行列を示し、diag(-1)は対角行列（行の数＝列の数＝代謝物質の数）を示し、vmarginal(t)はシステム境界酵素の酵素活性のベクトルを示し、vinternal(t)はシステム内部酵素活性のベクトルを示し、C'(t)は代謝物質濃度の速度ベクトルを示す。また、上記代謝物質濃度の勾配はC'(t)であり、酵素濃度時系列データの微分によって算出される値である。さらに、上記システム境界酵素について算出された活性値はvmarginal(t)である。）
【請求項７】
上記動的/静的モジュール区分け演算手段において上記適合関数は下記式で表され、代謝ネットワークを動的モジュール及び静的モジュールに区分けする全てのパターンについて下記式ｆ値を算出し、当該ｆ値を最小にするパターンを解として出力することを特徴とする請求項５記載の代謝ネットワーク推定装置。
【数４】
2008040812.tif 000015
（式中、Cdata,i,jはj番目のサンプリング点におけるi番目の代謝物質の濃度データを示し、Cestimated,i,jはj番目のサンプリング点におけるi番目の代謝物質の推定濃度データを示し、nmetaboliteは代謝物質の数を示し、nsampling pointはサンプリング点の数を示し、nenzymeはシステム内部酵素の数を示し、nstatic enzymeは静的モジュールに含まれる酵素の数を示し、wは重み係数を示す。）
【請求項８】
解析対象の代謝ネットワークに含まれるとともに、当該代謝ネットワークの外側との間に位置する酵素をシステム境界酵素とし、当該システム境界酵素は必ず動的モジュールに含まれるように上記パターンを定義することを特徴とする請求項７記載の代謝ネットワーク推定装置。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、対象の代謝ネットワークに含まれる酵素群を静的モジュール及び動的モジュールに区分けすることができる代謝ネットワーク推定方法及び代謝ネットワーク推定装置に関する。
【背景技術】
【０００２】
代謝ネットワークは非線形、非平衡及び非定常システムであるため、動的シミュレーションは、そのような複雑なシステムについて詳細且つ現実に即した挙動を分析又は予測する際に特に効果的である。しかしながら、代謝経路の動的モデルを構築する場合、酵素の反応メカニズム、反応速度、酵素及び代謝物質の初期濃度を含む膨大な量の実験的情報が必須となる。包括的な生化学データを取得できる様々な高スループット技術が開発されているが、酵素反応速度を求める実験的手法に関しては、未だ低スループットのままである。近年、酵素速度に関する幾つかのデータベースがインターネット上に公開されている（非特許文献１〜３参照）。しかしながら、多くの場合、これらのデータベース中のパラメータは、適切な代謝モデルを構築するには不十分である。細胞内情報は、公開文献から収集することができるが、実験条件及びターゲットとする生物種が一般的に統一されていない。従って、代謝系の適切な動的モデルを構築するには、膨大な量の実験的労力が必要となる。このため、その実験的労力を削減することができるモデル構築手法は大変有用な方法となり得る。
【０００３】
多くのモデル簡略化手法は動的シミュレーションの単純化を目的として提案されている。なかでもモデル簡略化技術として、動的/静的ハイブリッド（HDS）法が知られている（非特許文献４）。このHDS方法は、代謝系を動的モジュールと静的モジュールとに分割する。動的モジュールに含まれる酵素反応は微分方程式によって表現される。静的モジュールに含まれる酵素の活性は代謝流束解析（MFA）の処理によって算出される。MFAには速度論的な情報が必要でないため、酵素に関する速度論的情報を取得するための実験的労力を削減することができる。
【０００４】
HDS方法の優れた点は、モデル簡略化のスキームが単純であること及び代謝系における代謝物質数の許容性にある。酵素反応間の関係のみがHDS法に採用されるため、モデルを構築する場合、複数の時系列反応が一つの代謝化合物の移行に関連するといった問題を考慮する必要がない。未知の変数（すなわち静的モジュールにおける酵素活性）について、化学量論的係数行列のMoor-Penrose疑似逆行列をMFAに適用するため、未知の変数に対する化学量論的係数行列は正方正則である必要がない。HDS法にはこのような有利な点があるが、初期モデル構築の完成前に動的モジュールと静的モジュールとに区分けする方法は未だ提案されていない。
【０００５】
【非特許文献１】Schomburg, I., Chang, A. & Schomburg, D. (2002) BRENDA, enzyme data and metabolic information, Nucleic Acids Res. 30, 47-9.
【非特許文献２】Sundararaj, S., Guo, A., Habibi-Nazhad, B., Rouani, M., Stothard, P., Ellison, M. & Wishart, D. S. (2004) The CyberCell Database (CCDB): a comprehensive, self-updating, relational database to coordinate and facilitate in silico modeling of Escherichia coli, Nucleic Acids Res. 32, D293-5.
【非特許文献３】Selkov, E., Basmanova, S., Gaasterland, T., Goryanin, I., Gretchkin, Y., Maltsev, N., Nenashev, V., Overbeek, R., Panyushkina, E., Pronevitch, L., Selkov, E., Jr. & Yunus, I. (1996) The metabolic pathway collection from EMP: the enzymes and metabolic pathways database, Nucleic Acids Res. 24, 26-8.
【非特許文献４】Yugi, K., Nakayama, Y., Kinoshita, A. & Tomita, M. (2005) Hybrid dynamic/static method for large-scale simulation of metabolism, Theor Biol Med Model. 2, 42.
【発明の開示】
【発明が解決しようとする課題】
【０００６】
したがって、HDS法を利用して代謝ネットワークを単純化するとしても、代謝ネットワークを動的モジュールと静的モジュールとに適切に区分けする手法の開発が望まれている。そこで、本発明は、このような実情に鑑み、複数の酵素から構成された代謝ネットワークを動的モジュールと静的モジュールとに適切に区分けすることができる代謝ネットワーク推定方法及び推定装置を提供すること目的とする。
【課題を解決するための手段】
【０００７】
上述した目的を達成した本発明に係る代謝ネットワーク推定方法は、解析対象の代謝ネットワークを構成する代謝物質について得られた濃度時系列データ又は酵素活性時系列データを入力するステップと、上記ステップで入力された代謝物質濃度時系列データから算出した酵素活性時系列データ、或いは上記ステップで入力された酵素活性時系列データを用いて適合関数を最適化することで、上記代謝ネットワークを構成する各酵素を静的モジュール又は動的モジュールに区分けするステップとを含む。
【０００８】
また、本発明に係る代謝ネットワーク推定装置は、解析対象の代謝ネットワークを構成する代謝物質について得られた濃度時系列データ又は酵素活性時系列データを入力するデータ入力処理手段と、上記データ入力処理手段で入力された代謝物質濃度時系列データから算出した酵素活性時系列データ、或いは上記データ入力処理手段で入力された酵素活性時系列データを用いて適合関数を最適化することで、上記代謝ネットワークを構成する各酵素を静的モジュール又は動的モジュールに区分けする動的/静的モジュール区分け演算手段とを含む。
【発明の効果】
【０００９】
本発明に係る代謝ネットワーク推定方法及び代謝ネットワーク推定装置によれば、酵素活性（或いは代謝物質濃度）の時系列データを入力値とすることによって、複数の酵素から構成される代謝ネットワークを動的モジュール及び静的モジュールに区分けすることができる。
【発明を実施するための最良の形態】
【００１０】
以下、本発明に係る代謝ネットワーク推定方法及び代謝ネットワーク推定装置を、図面を参照して詳細に説明する。
【００１１】
先ず、本発明を適用した代謝ネットワーク推定装置のハードウェア構成について説明する。代謝ネットワーク推定装置のハードウェア構成は、図１に示すように、中央処理装置（CPU）１０１と、読取り専用メモリー（ROM）１０２と、書き換え可能メモリー（RAM）１０３と、ハードディスクドライブ（HDD）１０４と、ハードディスク（HD）１０５と、ディスプレイ１０６と、インターフェース（I/F）１０７と、キーボード等の入力装置１０８とを備えている。また、各構成部はバス１０９によってそれぞれ接続されている。また、ハードウェア構成としては、図示しないが、例えば、各種インターフェースで接続された外部記憶装置、スキャナ及びプリンタを備えていても良い。また、入力装置１０８としては、キーボードに加えてマウスやタッチパネル式の入力パッド等の装置であっても良い。
【００１２】
また、代謝ネットワーク推定装置は、I/F１０７を介して通信回線網とデータ通信可能であって、通信回線網を経由して各種のデータを入力するものであってもよい。ここで通信回線網は有線及び無線に限定されず、またローカルエリアネットワークでもよいし、インターネットであっても良い。
【００１３】
ここで、CPU１０１は、代謝ネットワーク推定装置の全体の制御を司る。ROM１０２は、後述する処理プログラムを記憶している。RAM１０３は、CPU１０１のワークエリアとして使用される。HDD１０４は、CPU１０１の制御下にHD１０５に対する各種データの書き込み／読み出しを制御する。
【００１４】
次に、代謝ネットワーク推定装置の機能的構成について説明する。代謝ネットワーク推定装の機能的構成は、図２に示すように、解析対象の代謝ネットワークを構成する酵素について得られた代謝物質濃度時系列データ又は酵素活性時系列データを入力するデータ入力処理部２０１と、データ入力処理部２０１で入力された代謝物質濃度時系列データから算出した酵素活性時系列データ、或いはデータ入力処理部２０１で入力された酵素活性時系列データを用いて適合関数を最適化することで、上記代謝ネットワークを構成する各酵素を静的モジュール又は動的モジュールに区分けする動的/静的モジュール区分け演算部２０２と、動的/静的モジュール区分け演算部２０２で区分けした代謝ネットワークを出力する出力処理部２０３とから構成されている。また、代謝ネットワーク推定装置は、既知の酵素活性時系列データ、代謝物質濃度時系列データ及び様々な代謝ネットワークに関するモデル構造データを格納したデータベース２０４を備えていても良い。
【００１５】
次に、代謝ネットワーク推定装置における入力データについて説明する。代謝ネットワーク推定装置における入力データは、解析対象の代謝ネットワークを構成する代謝物質について得られた時系列濃度データ（すなわち、代謝物質濃度時系列データ）、又は解析対象の代謝ネットワークを構成する酵素について得られた時系列活性データ（すなわち、酵素活性時系列データ）である。代謝物質濃度時系列データは、例えば、キャピラリー電気泳動-質量分析装置を使用して得ることができる。キャピラリー電気泳動-質量分析装置によれば、解析対象の代謝ネットワークに含まれる多数の代謝物質に関する時系列濃度データを網羅的に取得することができる。このとき、キャピラリー電気泳動-質量分析装置は、所定の間隔でサンプリングした試料を順次解析にかけることで、代謝物質に関する濃度データを時系列的に取得することができる。代謝ネットワーク推定装置は、キャピラリー電気泳動-質量分析装置から出力された代謝物質濃度時系列データを直接的に入力してもよいし、一旦、データベース２０５に格納してもよい。
【００１６】
本発明に係る代謝ネットワーク推定装置は、従来公知の代謝ネットワークであれば、如何なる代謝ネットワークにも適用することができる。本例では、以下のテストモデルを用いた。２つのテストモデルは、Chassagnoleらによって構築された大腸菌モデル（Chassagnole, C., Noisommit-Rizzi, N., Schmid, J. W., Mauch, K. & Reuss, M. (2002) Dynamic Modeling of the Central Carbon Metabolism of Escherichia coli, Biotechnol Bioeng. 79, 53-73参照）と、Hynneらによって構築された酵母モデル（Hynne, F., Dano, S. & Sorensen, P. G. (2001) Full-scale model of glycolysis in Saccharomyces cerevisiae, Biophys Chem. 94, 121-63参照）である。
【００１７】
大腸菌モデルについては、細胞外グルコース濃度が0.0556mMである定常状態から開始し、グルコースパルスを添加した。このパルス導入の瞬間におけるグルコース濃度は、1.667mMであった。Chassagnoleのオリジナルモデルにおいて、補酵素（ATP/ADP/AMP、NAD(H)、NADP(H)）の時系列は、時間依存性関数によって表現されていた（上記Chassagnole, Cらの論文参照）。しかしながら、本研究では、これら補酵素の濃度は初期値に固定した。酵母モデルについては、供給溶液中のグルコース濃度が2.5mMである定常状態から開始し、グルコース濃度を5.0mMに移行させた。定常状態における両モデル中の代謝物質濃度、すなわち動的シミュレーションのための初期濃度を表１に示した。摂動後の実行時間を、大腸菌については20秒に設定し、酵母については60秒に設定した。計算された時系列代謝物質濃度は、0.1秒間隔で記録され、擬似的実験データとして使用した。
【００１８】
【表１】
2008040812.tif 000003
【００１９】
次に、これら２つのテストモデルを用いた代謝ネットワーク推定装置における処理ステップについて説明する。
【００２０】
酵素活性の推定
先ず、動的/静的ハイブリッド法（HDS法）については、Yugi, K., Nakayama, Y., Kinoshita, A. & Tomita, M. (2005) Hybrid dynamic/static method for large-scale simulation of metabolism, Theor Biol Med Model. 2, 42.に開示されている。HDS法において、代謝ネットワークは静的モジュールと動的モジュールとに区分される。静的モジュールに含まれる酵素活性は以下の式によって計算される。
【００２１】
【数１】
2008040812.tif 000004
【００２２】
式中、vstaticは静的モジュールの酵素活性のベクトルを示し、vdynamicはモジュール境界酵素の酵素活性のベクトルを示し、Sstatic*は静的モジュールに含まれる酵素に関する化学量論係数行列のMoor-Penrose疑似逆行列を示し、Sdynamicは静的モジュールに含まれる酵素に関する化学量論係数行列を示す。また、式中、ｔは時間である。
【００２３】
図３に、動的モジュール、静的モジュール、モジュール境界酵素、システム境界酵素及びシステム内部酵素の関係を示した。当該代謝ネットワークの外側との間に位置する酵素をシステム境界酵素とし、システム境界酵素を除いた当該代謝ネットワークに含まれる酵素をシステム内部酵素とする。また、動的モジュールと静的モジュールとの境界に位置する酵素をモジュール境界酵素とする。モジュール境界酵素及びシステム境界酵素は、動的モジュールに含まれるものとする。１つの酵素がシステム境界酵素であり且つモジュール境界酵素である場合もある。
【００２４】
静的モジュールに含まれる酵素について酵素活性を推定するため、動的モジュールに含まれる酵素の酵素活性を使用する。具体的には、システム境界酵素及びシステム内部酵素の定義とともに、動的遷移状態における代謝系の物質収支方程式は以下の式（２）で与えられる。
【００２５】
【数２】
2008040812.tif 000005
【００２６】
式中、Sは化学量論係数行列を示し、diag(-1)は対角行列（行の数＝列の数＝代謝物質の数）を示し、vmarginal(t)はシステム境界酵素活性のベクトルを示し、vinternal(t)はシステム内部酵素活性のベクトルを示し、C'(t)は代謝物質濃度の速度ベクトルを示す。
【００２７】
ここで、C'(t)及びvmarginalが既知の場合には、システム内部酵素活性は式（２）を変換した式（３）によって推定することができる。
【００２８】
【数３】
2008040812.tif 000006
【００２９】
ある摂動下におけるテストモデルの動的シミュレーションを幾つか実行した。シミュレーションの際には、算出された代謝物質濃度時系列及び酵素活性時系列を記録した。代謝物質濃度時系列は疑似実験データとして見なし、C'(t)及びvmarginalを取得するのに使用した。代謝物質濃度の勾配（C'(t)）は、代謝物質濃度の時系列の時間微分によって取得した。式（３）に代謝物質濃度の勾配とシステム境界酵素の活性とを代入することで、システム内部酵素に関して酵素活性時系列データを推定することができる。
【００３０】
表２には、システム内部酵素の活性について推定した結果をまとめている。表２から、殆ど全ての酵素において誤差を無視できることが判る。大腸菌におけるTKｂ及びTA、並びに酵母におけるAKは、幾つかのサンプリングポイントで０に近い値を示す真の活性を含んでいた（データを示さない）。真の活性は誤差算出式（後述する式（４））における分母として使用するので、これら３つの酵素に関する平均誤差は明らかに大きな値を示すことになる。これらの点を除いた後、これら３つの酵素に関する予測された活性は、十分な正確性を示した。この結果から、酵素活性時系列データは、代謝物質の濃度に関する実験的に測定可能な時系列データから推定できることが示された。この予測は、系における物質収支のみに依拠するためシステム内部酵素に関するいかなる速度情報も必要ではないことが有利な点である。
【００３１】
【表２】
2008040812.tif 000007
【００３２】
なお、表２の誤差は、後述する式（４）によって算出された値である。
【００３３】
酵素活性時系列又は代謝物質濃度時系列における真の値と推定値との間の誤差は、以下の式によって計算された。
【００３４】
【数４】
2008040812.tif 000008
【００３５】
式中、Xtrue(t)は酵素活性又は代謝物質濃度の真の値を示し、XEstimated(t)は酵素活性又は代謝物質濃度の推定値を示し、nはサンプリング点の数を示す。
【００３６】
動的モジュールと静的モジュールとの区分け処理
本発明に係る代謝ネットワーク推定装置では、次に、上述したように算出された酵素活性時系列推定値を用いて、解析対象の代謝ネットワークを動的モジュールと静的モジュールとに区分けする。動的モジュール/静的モジュールとの区分けについての最適パターンは、所定の適合関数を最適化することを目的とした、仮定された区分けの反復修飾によって見いだすことができる。
【００３７】
代謝系における動的/静的酵素の最適組合せを検索するため、遺伝的アルゴリズム（GA）（Willett, P. (1995) Genetic algorithms in molecular recognition and design, Trends Biotechnol. 13, 516-21.）を採用した。本発明に係る代謝ネットワーク推定装置では、GAに対する個体（数値コードセット）は、代謝系に含まれる動的酵素/静的酵素を表現するように定義された。例えば、個体 ”DDSSDD”は、６つの酵素からなる代謝系を示しており、１番目、２番目、５番目及び６番目の酵素が動的モジュールに含まれ、３番目及び４番目の酵素が静的モジュールに含まれる。
【００３８】
GA計算において、静的モジュールの酵素活性時系列は、HDS法の手順と同様に、代謝物質の濃度データから導かれる動的酵素の活性から計算される。その結果、それぞれの代謝物質濃度時系列が、それぞれの酵素活性時系列の数値積分とともに計算される。式5で定義される適合関数は、それぞれの個体について計算され、その後、伝播、交差及び変異を行う。この手順は、式5を最小とする最適解が見つかるまで繰り返される。動的/静的酵素の区分けを行う処理のフローチャートを図４に示す。問題の収束は、生成数に対する適合関数の勾配によって確認される。
【００３９】
【数５】
2008040812.tif 000009
【００４０】
式中、Cdata,i,jはj番目のサンプリング点におけるi番目の代謝物質の濃度データを示し、Cestimated,i,jはj番目のサンプリング点におけるi番目の代謝物質の推定濃度データを示し、nmetaboliteは代謝物質の数を示し、nsampling pointはサンプリング点の数を示し、nenzymeはシステム内部酵素の数を示し、nstatic enzymeは静的モジュールに含まれる酵素の数を示し、wは重み係数を示す。
【００４１】
適合関数（式5）における第二項は、代謝系に含まれる静的酵素の比を評価しており、静的モジュールにおける酵素の数を調整するために加えられている。調整用のパラメータである重み係数（w）が、モデル誤差と静的酵素比との間の重要性均衡を制御する目的で第二項に乗算されている。
【００４２】
大腸菌及び酵母におけるテスト代謝系についての動的/静的酵素の最適区分けを、式5のもとに推定した。３つの様々な重み係数の値をそれぞれ採用した（w＝1.00、0.10及び0.01）。表３に動的/静的酵素の最適区分けを示す。予測されたように、重み係数が低い値へと推移するのに沿って、誤差は低減し、静的酵素の割合は減少した。
【００４３】
【表３】
2008040812.tif 000010
【００４４】
表３中、Dは動的モジュールへの区分けを意味し、Sは静的モジュールへの区分けを意味する。また、システム境界酵素は、全て動的モジュールに区分けされるため、表３に挙げていない。
【００４５】
次に、大腸菌及び酵母における全動的モデルを、予測された動的/静的の最適区分け（表３）を用いてハイブリッドモデルに変換した。表４には、もとの全動的モデルと変換後のハイブリッドモデルとの間における、代謝物質濃度の誤差を示している。
【００４６】
動的モジュールと静的モジュールの区分け工程、すなわち、酵素活性時系列を示す所定の曲線を数値積分する工程において、１つのサンプリング点における計算された静的酵素の活性は、次のサンプリング点に影響しない。一方、HDS法、すなわち連立微分方程式の初期値問題においては、１つのサンプリング点における計算された静的酵素の活性は、次のステップにおける計算に影響する。したがって、動的モジュールと静的モジュールの区分け工程において算出された誤差は、取得されたハイブリッドモデルにおける誤差と必ずしも等価ではない。したがって、これら２種の計算における誤差間の比較をすることが必要である。
【００４７】
【表４】
2008040812.tif 000011
【００４８】
表４に示す平均誤差は表３に示したものと類似している。この結果は、上述した動的モジュールと静的モジュールの区分け工程がハイブリッドモデルの構築に使用可能であることを示している。
【００４９】
重み係数1.00では、約半数の酵素が静的モジュール（表３）に区分けされた。これにより、速度データを取得するための実験的仕事量を大幅に削減できる。この重み誤差1.00における平均誤差は、大腸菌について14.5％であり、酵母について17.3％であった（表４）。これらの誤差サイズは、代謝物質濃度の実験的測定における精度を考慮すると実用的に受け入れられるものである。より精緻なモデルが要求される場合には、より低い値の重み係数を使用することができる。重み係数を0.10にセットしたとしても、両テスト系について静的酵素の比率は30％程度に維持されていた（表４）。
【００５０】
また、式（５）で示す適合関数における重み係数（W）はチューニングパラメータである。重み計数の最適値は推測的に与えられるものではないため、その値をどのようにして規定するか考慮する。図５には、重み係数値と誤差との関係、並びに重み係数値と静的酵素の比率との関係を示している。２つのテスト系においては、非常に類似した結果が得られた。これは、両モデルが共に主要炭素代謝を扱っているためであると考えられる。誤差サイズと静的酵素の比率との間のトレードオフの関係からすると、重み係数は約0.1の近傍であることが望ましい。
【００５１】
以上のように、本発明に係る代謝ネットワーク推定装置によれば、多くのフィードバック・ループを有し、且つ生細胞の恒常性を維持するように厳格に制御された主要炭素代謝ネットワークを動的モジュールと静的モジュールとに区分けすることができた。したがって、本発明に係る代謝ネットワーク推定装置は、あらゆる種類の代謝系に対しても広く適用可能であると考えられる。
【００５２】
なお、本明細書では、上述した２つのテストモデルにおける代謝物質及び酵素並びに反応を略語で示したが、略語は以下の物質及び酵素並びに反応を意味している。
【００５３】
代謝物質
2PG：      2-phosphoglycerate
3PG：      3-phosphoglycerate
6PG：      6-phosphogluconate
ACA：      acetaldehyde, intracellular
ACAx：      acetaldehyde, extracellular
ADP：      adenosine diphosphate
AMP：      adenosine monophosphate
ATP：      adenosine triphosphate
CNo：      cyanide, mixed flow
CNx：      cyanide, extracellular
DHAP：      dihydroxyacetone phosphate
E4P：      erythrose 4-phosphate
EtOH：      ethanol, intracellular
EtOHx：      ethanol, extracellular
F6P：      fructose 6-phosphate
FDP：      fructose 1,6-bisphosphate
G1P：      glucose 1-phosphate
G6P：      glucose 6-phosphate
GAP：      glyceraldehyde 3-phosphate
Glco：      glucose, mixed flow
Glcx：      glucose, extracellular
Glyc：      glycerol, intracellular
Glycx：      glycerol, extracellular
NAD：      nicotinamide adenine dinucleotide, oxidized
NADH：      nicotinamide adenine dinucleotide, reduced
NADP：      nicotinamide adenine dinucleotide phosphate, oxidized
NADPH：      nicotinamide adenine dinucleotide phosphate, reduced
PEP：      phosphoenolpyruvate
PGP：      1,3-bisphospho glycerate
Pyr：      pyruvate
Rib5P：      ribose 5-phosphate
Ribu5P：      ribulose 5-phosphate
Sed7P：      sedoheptulose 7-phosphate
Xyl5P：      xylulose 5-phosphate
【００５４】
酵素/反応
ADH：      acetaldehyde dehydrogenase
AK：      adenylate kinase
ALDO：      aldolase
consum：      ATP consumption
difACA：      diffusion of acetaldehyde
difEtOH：      diffusion of EtOH
difGlyc：      diffusion of glycerol
ENO：      enolase
G6PDH：      glucose 6-phosphate dehydrogenase
GAPDH：      glyceraldehyde 3-phosphate dehydrogenase
GlcTrans：      glucose transporter
HK：      hexokinase
lpGlyc：      lumped glycerol formation reaction
lpPEP：      lumped PEP formation reaction
PDC：      pyruvate decarboxylase
PFK：      phosphofructokinase
PGDH：      6-phosphogluconate dehydrogenase
PGI：      glucose-6-phosphate isomerase
PGK：      phosphoglycerate kinase
PGluMu：      phosphoglycerate mutase
PGM：      phosphoglucomutase
PK：      pyruvate kinase
R5PI：      ribose-phosphate isomerase
Ru5P：      ribulose-phosphate epimerase
TA：      transaldolase
TIS：      triosephosphate isomerase
TKa：      transketolase, reaction a
TKb：      transketolase, reaction b
【図面の簡単な説明】
【００５５】
【図１】本発明を適用した代謝ネットワーク推定装置のハードウェア構成の一例を示すブロック図である。
【図２】本発明を適用した代謝ネットワーク推定装置の機能的構成の一例を示すブロック図である。
【図３】動的モジュール、静的モジュール、モジュール境界酵素、システム境界酵素及びシステム内部酵素の関係を示す特性図である。
【図４】本発明に係る代謝ネットワーク推定装置の処理ステップを示すフローチャートである。
【図５】適合関数における重み係数（W）を変更した際の誤差への影響及び静的酵素の比率への影響を検討した結果を示す特性図である。
【符号の説明】
【００５６】
101…CPU、102…ROM、103…RAM、104…HDD、105…HD、106…ディスプレイ、107…I/F、108…入力装置、109…バス、201…データ入力処理部、202…動的/静的モジュール区分け演算部、203…出力処理部、204…データベース
【図１】
2008040812.tif 000016
【図２】
2008040812.tif 000017
【図３】
2008040812.tif 000018
【図４】
2008040812.tif 000019
【図５】
2008040812.tif 000020