研究の成果

ソフトウェアの成果

研究の初期段階においては、与えられた系統樹の尤度を計算するための方程式 の構造を、そのままKL1でコーディングしたが、この結果得られた結果は、き わめて悲観的なものだった。再帰的に尤度を計算するアルゴリズムに検討を加 え、再帰的な計算の途中で生じる中間結果をベクターとして保存するようにし たところ、劇的な高速化を達成できた。

プログラムの尤度曲面の変化の様子を視覚化するためのプログラム contour/1 をPIM上に実装し、その後Sunワークステーション上に移植した。 このプログラムを用いると、系統樹の任意の2本の枝を少しずつ変化させたと き、尤度曲面がどのような振る舞いをするかを直感的に観察することができ、 その性質を研究するのに有用である。だだし、系統樹作成法における尤度とい う尺度の解析的な研究は今後必要である。contour/1 の負荷分散の実験は、 PIM/m上のみで行った。多量の浮動少数点演算をともなう計算において、我々 が期待するような負荷分散の結果が得られた。

実用的なサイズのアミノ酸配列を扱えるよう改良した枝長最適化プログラム、 traverse/3 をKL1で書いたが、現在のところPIM上では動作していない。原因 は、プログラムの稼働中にシステムがメモリ不足に陥り、OSそのものが停止し てしまうからである。ゴール・プライオリティの変更やプログラムの修正を行 い、要求メモリ量の圧縮を試みたが、現在のところごく少量のデータを用いた 場合を除いて、すべて失敗している。

他方、新たに構築したPVMによる分散環境では、traverse/3 は動作に成功し ている。理由は、後者の方が1CPU当たりのメモリ資源の大きさが圧倒的に大き いからであると考えられる。このプログラムを用いて、距離行列法て推定した 系統樹の枝長を再推定するなどの研究を行った。

traverse/3 に与えるためのトポロジーを生成する topology/2 をKL1で書いた。これは、いわば完全探索のための仮説空間生成器である。 traverse/3 と topology/2 を結合して、go/0 というプログラムを書き、分散環境に実装した。

我々はデータ解析の仕事において必要なプログラムを、基本的にKL1を用いて 書いている。その結果、実用的なデータ解析プログラムの資産が蓄積しつつあ る。例えば、データベースの検索、データの抽出、ドットマトリックスの作成、 そして文献データベースの生成プログラム等である。

開発したプログラム(述語)の一覧は以下のようになる。括弧内はモジュール名 とKL1プログラムの行数を表す。

ancestor/3

Infers an ancestral site for a pair of descendants (ancestor:345).

beautify/2

Beautifies a printed format of a nested list structure, especially for expressing structure of a database (beautify:122).

b/4

Generates a binomial distribution (probability:87).

c/3

Takes a binomial coefficient (probability:87).

c2b/2

Converts a cluster-like topology format to a branch-like one (tree:442).

clike/8

Computes a conditional likelihood value (traverse:1924).

cofactor3/3

Calculates a cofactor for a given 3x3 matrix (matrix:103).

consensus/2

Generates the consensus sequence for a given set of sequences (consensus:96).

connect/3

Finds the descendant nodes for a given node (tree:442).

connect/4

Finds the descendant nodes for a given node in a unrooted tree (tree:442).

convert_tree/2

Converts a distance (Number of substitution) to an internal expression P (traverse:1924).

curvature/5

Draws a likelihood surface for given two nucleotide sequences in a tree (curvature:597).

dX/3

Computes a unbiased dx for given sequences (dX:215).

det2/2

Calculates a determinant of a given 2x2 matrix (matrix:103).

det3/2

Calculates a determinant of a given 3x3 matrix (matrix:103).

divide/3

Divides a list into a given number of segments equally (divide:104).

dot_matrix/4

Generates a dot matrix (dot_matrix:125).

extract_author/2

Extracts authors from PIR database (extract_author:70).

extract_fields/3

Extracts selected fields from PIR database (extract_fields:235).

factorial/2

Takes an approximate factorial to a given floating point number (float:59).

floatadd/3

Adds floating point data having arbitrary accuracy (add:185).

format/2

Makes formatted lists from blocked lines (format:149).

fx/8

Computes Fx (fx:306).

get_seq

Reads aligned sequences, and makes a character state format (get_seq:554).

int:factorial/2

Takes an exact factorial to a given integer (int:18).

interleave/2

Converts an interleaved format to a blocked format (interleave:118).

invert/2

Inverts an order of a given list (invert:23).

lines2string/2

Converts given lines to a string (lines2string:29).

join/3

Joins two lists (list:71).

join_lists/2

Joins lists (list:71).

get_flat_list/2

gets a flat list from a nested list (list:71).

log_like/5

Computes log likelihood for a given tree (log_like:324).

make_bib/2

Makes a bibliography from PIR database (make_bib:952).

make_destbl/2

Makes a description table for annotating a given tree (make_destbl:276).

make_gap_matrix/2

Makes a minimum gap matrix (mgm) (make_gap_matrix:408).

make_records/2

Makes records from PIR database (make_records:254).

make_seq/2

Extract sequences from PIR database, and converts them to FASTA format (make_seq:493).

minor3/2

Calculates a minor of a given 3x3 matrix (matrix:103).

minorm/4

Generates a minor matrix of a given matrix (matrix:103).

newton/7

Solves a given equation by the newton method (newton:616).

p/3

Calculates number of permutations (probability:87).

plot_like/3

Plots likelihood values (plot_like:575).

plot_xyz/1

Plots 3-D coordinates from a given list(termio:171).

print_sequences/2

Prints sequence from a list of integer-presented sequence data (print_sequences:93).

putt/1

Puts a Prolog term to standard output (termio:171).

putts/1

Puts a Prolog terms to standard output (termio:171).

read_terms/2

Reads Prolog terms from a given file (termio:171).

remove_gap_site/2

Removes gap sites from given sequences, and generates input file for traverse/3 (remove_gap_site:188).

search_str/3

Searches a string (search_str:158).

starling/2

Calculates a factorial by Starling's formula (probability:87).

str2lst/2

Converts a string to a list (string:52).

synonymous/3

Computes numbers of synonymous and nonsynonymous sites of a given sequence (synonymous:372).

topology/2

Generates a topology space for a given number of sequences (tree:442).

traverse/4

Optimizes branch lengths of a given trees according to likelihood values (traverse:1924).

このうち、述語traverse/4の出力例を図--2に示す。

  
図 --2: 述語traverse/4の出力例。グラフィックデータの表示にはTreeview (Page R., 1996)を用いている。

研究上の成果

contour/1 による、人工的なデータを用いて描いた尤度曲面から、場合によっ てはこの曲面がきわめて扁平で非対称な形をしていることがわかった。実際に は枝長の最適化アルゴリズムとは言っても、曲面の微分値が0になる点を、EM アルゴリズムで探索しているに過ぎないので、プログラムがどのような探索パ スをたどるかによって、結果が変わってくる可能性がある。少なくともこの尤 度曲面を描いたデータセットからは、枝長の長い方からではなく、短い方から 探索を開始すべきであるということが示唆される。なぜなら、EMアルゴリズム に与える初期値が最大尤度に対応する値より大きい場合、探索パスが最大尤度 点に達することができないとすると、推定された枝長はより大きな誤りを含む からである。もちろん、探索パスが最大尤度点近傍を通り越したとき後戻りで きるようにアルゴリズムを書くことはできる。しかし、この場合も本質的な解 決になるわけではない。今度は枝長の短い方から探索を開始した場合、同じ問 題にぶつかってしまうからである。このように考えると、どこから尤度曲面上の 探索を開始するかということは、厳密に尤度曲面の頂上を決めようとするとき 重要な問題となってくることがわかった。

traverse/3 は、1回目の枝長の最適化で、ほぼ最大尤度を得ることができる。 むしろ、2回目以降は微々たる改良でしかない。現在のアルゴリズムでは、改 良された尤度の差を見て探索の打ち切りを決めているのだが、この1回目の改 良された尤度の値を見て、これ以降の処理を続けるかどうかを判断できるかも しれない。つまり、多くの系統樹を吟味し、1回目の改良で得られた尤度と真 の尤度との間の関係を統計的に調べるのである。。また、トポロジーに無理の ある系統樹の尤度はなかなか改良されないと予想できる。つまりこのような系 統樹の尤度曲面は、きわめて「平坦」となる。このことが統計的に有意である ことを示せれば、改良のための繰り返し回数は、処理の打ち切りのためのよい 指標となろう。

なお、本研究は以下の学会等で発表された。

• OOta S., Saitou N., and Kunifuji S.: Inference of the Phylogenetic Relationship Among Motor Proteins Using a Newly Developed Program for the Maximum Likelihood Method. International Symposium on Network and Evolution of Molecular Information (NEMI), Tokyo, April.

• OOta S., Saitou N., and Kunifuji S.: Inference of the Phylogenetic Relationship Among Motor Proteins Using a Newly Developed Program for the Maximum Likelihood Method. The Fourth Annual Meeting of the Society for Molecular Biology and Evolution (SMBE4), Tucson, June.

• 太田 聡史、斎藤 成也：最尤推定プログラムによるモータータンパク質間の進 化的関係の推定。第19回日本分子生物学・第69日本生化学会大会合同年会(名 古屋)，８月。

• 太田 聡史、斎藤 成也：モータータンパク質の進化的関係の推定。日本遺伝学 会第６８回大会(名古屋)，１０月。