(57)【要約】【課題】  互いにスプライシング・バリアントの関係にある配列を含む多数の生体高分子配列データから，複数の同一配列のエクソンが同一の入力配列に現れる場合を考慮しつつ，個々のエクソンに相当する部分配列を，効率よく抽出する。【解決手段】  エクソン配列が満たすべき条件を用いて，複数の同一配列のエクソンが同一の入力配列に現れる場合を考慮したエクソン配列のモデルを構築し，そのモデルによって定義されるエクソン配列を，入力として与えられた生体高分子配列から抽出する。本発明の方法は，入力された生体高分子配列に基づきsuffix treeを構築し，配列中の文字数や配列数に依存しない回数の，suffix tree上の深さ優先探査および配列中の各文字の位置に関する処理によって，入力として与えられた生体高分子配列の長さの総和に対し，線形時間で抽出処理を完了できる。【選択図】  図１
2005316605.tif 000002
【特許請求の範囲】
【請求項１】
文字列入力手段から複数の文字列が入力されたときに，演算部において，
入力された文字列のうち少なくとも1つの文字列に，1箇所または2箇所以上出現する部分文字列であって，長さが与えられた整数u以上であり，該部分文字列の右または左に一文字追加して得られる新しい文字列が出現する入力された文字列中の位置の集合が，文字を追加する前の該部分文字列が出現する位置の集合と異なる集合となる部分文字列を，MMと呼ぶとき，
(1)入力された文字列のうち少なくとも1つの文字列に，1箇所または2箇所以上出現する部分文字列であって，
(2)長さが前記整数u以上であり，
(3)該部分文字列が，該部分文字列を完全に含むMM以外のどのMMとも，どの入力された文字列上でも文字を共有することがなく，
上記3つの条件(1)〜(3)のみを満足するほかの文字列の一部となっていない，EGと呼ばれる文字列を抽出する処理を実行することを特徴とする，文字列解析方法。
【請求項２】
請求項1記載の文字列解析方法において，前記演算部は，
入力された文字列のうち少なくとも1つの文字列に，1箇所または2箇所以上出現する部分文字列であって，長さが与えられた整数u以上であり，該部分文字列の右に一文字追加して得られる新しい文字列が出現する入力された文字列中の位置の集合が，文字を追加する前の該部分文字列が出現する位置の集合と異なる集合となる部分文字列である，右MMと呼ぶ部分文字列を抽出する工程，
入力された文字列のうち少なくとも1つの文字列に，1箇所または2箇所以上出現する部分文字列であって，長さが前記整数u以上であり，該部分文字列の左に一文字追加して得られる新しい文字列が出現する入力された文字列中の位置の集合が，文字を追加する前の該部分文字列が出現する位置の集合と異なる集合となる文字列である，左MMと呼ぶ文字列を抽出する工程または，前記MMと呼ばれる部分文字列を抽出する工程，
MMであるか，またはある入力された文字列上においてMMの次の文字から始まる右MMであって，右MMである文字列を真のprefixとして含まない文字列である，右EG-holderと呼ばれる文字列を抽出する工程，
を実行することを特徴とする文字列解析方法。
【請求項３】
文字列入力手段から複数の文字列が入力されたときに， 
入力された文字列のうち少なくとも1つの文字列に，1箇所または2箇所以上出現する部分文字列であって，長さが与えられた整数u以上であり，該部分文字列の右または左に一文字追加して得られる新しい文字列が出現する入力された文字列中の位置の集合が，文字を追加する前の該部分文字列が出現する位置の集合と異なる集合となるMMと呼ばれる部分文字列を抽出する方法であって，
演算部が，
入力された文字列のうち少なくとも1つの文字列に，1箇所または2箇所以上出現する部分文字列であって，長さが前記整数u以上であり，該部分文字列の右に一文字追加して得られる新しい文字列が出現する入力された文字列中の位置の集合が，文字を追加する前の該部分文字列が出現する位置の集合と異なる集合となる，右MMと呼ばれる部分文字列を抽出する処理を実行することを特徴とする，文字列解析方法。
【請求項４】
文字列入力手段から複数の文字列が入力されたときに，演算部が，
入力された文字列のうち少なくとも1つの文字列に，1箇所または2箇所以上出現する部分文字列であって，長さが与えられた整数u以上であり，該部分文字列の右に一文字追加して得られる新しい文字列が出現する入力された文字列中の位置の集合が，文字を追加する前の該部分文字列が出現する位置の集合と異なる集合となる部分文字列である，右MMと呼ぶ部分文字列を抽出する工程，
入力された文字列のうち少なくとも1つの文字列に，1箇所または2箇所以上出現する部分文字列であって，長さが前記整数u以上であり，該部分文字列の左に一文字追加して得られる新しい文字列が出現する位置の集合が，文字を追加する前の該部分文字列が出現する位置の集合と異なる集合となる部分文字列である左MMと呼ばれる文字列を抽出する工程を実行し，
入力された文字列のうち少なくとも1つの文字列に，1箇所または2箇所以上出現する部分文字列であって，長さが前記整数u以上であり，該部分文字列の右または左に一文字追加して得られる新しい文字列が出現する入力された文字列中の位置の集合が文字を追加する前の該部分文字列が出現する位置の集合と異なる集合となる部分文字列をMMと呼ぶとき，
MMであるか，またはある入力された文字列上においてMMの次の文字から始まる右MMであって，右MMである文字列を真のprefixとして含まない文字列である，右EG-holderと呼ばれる部分文字列を抽出する処理を実行することを特徴とする文字列解析方法。
【請求項５】
文字列入力手段から複数の文字列が入力されたときに，
入力された文字列のうち少なくとも1つの文字列に，1箇所または2箇所以上出現する部分文字列であって，長さが与えられた整数u以上であり，該部分文字列の右に一文字追加して得られる新しい文字列が出現する入力された文字列中の位置の集合が，文字を追加する前の該部分文字列が出現する位置の集合と異なる集合となるとき該部分文字列を右MMと呼び，
入力された文字列のうち少なくとも1つの文字列に，1箇所または2箇所以上出現する部分文字列であって，長さが前記整数u以上であり，該部分文字列の右または左に一文字追加して得られる新しい文字列が出現する入力された文字列中の位置の集合が，文字を追加する前の該部分文字列が出現する位置の集合と異なる集合となるとき該部分文字列をMMと呼ぶとき，
演算部において，
MMであるか，またはある入力された文字列上においてMMの次の文字から始まる右MMであって，右MMである文字列を真のprefixとして含まない文字列である，右EG-holderと呼ばれる部分文字列を抽出する工程を実行する文字列解析方法であって，
入力された文字列におけるMMの位置を計算するために，他の右MMを真のsuffixとしない右MMである，極小右MMと呼ばれる文字列の位置を計算する処理を実行することを特徴とする文字列解析方法。
【請求項６】
文字列入力手段から複数の文字列が入力されたときに，
入力された文字列のうち少なくとも1つの文字列に，1箇所または2箇所以上出現する部分文字列であって，長さが与えられた整数u以上であり，該部分文字列の右に一文字追加して得られる新しい文字列が出現する入力された文字列中の位置の集合が，文字を追加する前の該部分文字列が出現する位置の集合と異なる集合となるとき該部分文字列を右MMと呼び，
入力された文字列のうち少なくとも1つの文字列に，1箇所または2箇所以上出現する部分文字列であって，長さが前記整数u以上であり，該部分文字列の右または左に一文字追加して得られる新しい文字列が出現する入力された文字列中の位置の集合が，文字を追加する前の該部分文字列が出現する位置の集合と異なる集合となるとき該部分文字列をMMと呼ぶとき，
演算部において，
MMであるか，またはある入力された文字列上においてMMの次の文字から始まる右MMであって，右MMである文字列を真のprefixとして含まない文字列である，右EG-holderと呼ばれる部分文字列を抽出する工程を実行する文字列解析方法であって，
入力された文字列におけるMMの位置を計算するために，他の右MMを真のprefixとしない右MMである，prefix極小右MMと呼ばれる文字列の位置を計算する処理を実行することを特徴とする文字列解析方法。
【請求項７】
請求項１記載の文字列解析方法において，
入力された文字列のうち少なくとも1つの文字列に，1箇所または2箇所以上出現する部分文字列であって，長さが与えられた整数u以上であり，該部分文字列の右に一文字追加して得られる新しい文字列が出現する入力された文字列中の位置の集合が，文字を追加する前の該部分文字列が出現する位置の集合と異なる集合となる部分文字列を右MMと呼ぶとき，
前記演算部は，
前記MMであるか，またはある入力された文字列上においてMMの次の文字から始まる右MMであって，右MMである文字列を真のprefixとして含まない文字列である，右EG-holderのprefixだけを，抽出すべき部分文字列EGの候補とすることを特徴とする文字列解析方法。
【請求項８】
請求項５に記載の文字列解析方法において，ある右MMに，該右MMの真のsuffixである別の右MMが存在するか否かを調べる際，suffix tree上で該右MMをパスラベルのprefixとするノードの，親ノードのsuffix linkを参照することを特徴とする文字列解析方法。
【請求項９】
請求項１に記載の文字列解析方法において，入力文字列上でEGの境界とすべき位置が与えられた場合に，EGが該位置を含まないようにEGを分割することを特徴とする文字列解析方法。
【請求項１０】
請求項１〜９のいずれか１項記載の文字列解析方法をコンピュータに実行させるためのプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は，複数の文字列，特に生体高分子配列に共通に存在する部分文字列を抽出するための文字列解析方法に関し，生体高分子配列のスプライシングパターン解析方法に適用して好適な方法に関する。
【背景技術】
【０００２】
国際共同プロジェクト及び米国ベンチャー企業により，2000年6月にヒトゲノムの文字配列決定の完了が宣言された。ゲノム配列の解析が進む一方で，発現している遺伝子について調べるために，mRNAの解析が行なわれている。mRNAは，遺伝子が発現する際，ゲノムDNAから生成されるRNA分子で，遺伝子の機能発現の過程で不可欠な物質である。mRNA分子は分解しやすいため，逆転写によりmRNAよりも安定な物質であるcDNAに転換し解析されることが多い。cDNAを配列決定して得られるcDNA配列には様々な利用価値があるが，そのひとつが，遺伝子がどのようなエクソンにより構成されているかを明らかにすることである。
【０００３】
同一遺伝子由来のmRNAで，スプライシングを受ける領域が異なるものを，互いにスプライシング・バリアントであると言う。ヒトの場合，生体内のタンパク質は10万種類以上存在するといわれる一方，遺伝子数は3万〜4万のみといわれており，スプライシング・バリアントがこの差を埋めていると考えられている。こうしたスプライシング・バリアントの解析は，遺伝子の機能発現を解析するにあたり不可避であり，生命現象の解明や，ゲノム創薬のために不可欠である。同一遺伝子に由来するあらゆるスプライシング・バリアントの配列を網羅する十分な量のcDNA配列を収集し解析することは，mRNAのエクソン構造を明らかにする有力な手段となる。cDNA配列など転写産物由来の配列のみに基づきスプライシング・バリアントの解析を行うためには，複数の配列から個々のエクソンに相当する部分配列をいかにして抽出するかが課題となる。ところが，与えられた配列のあらゆる部分配列を網羅的に列挙し，エクソン配列か否かを調べるアプローチをとる場合は，各入力配列に対してその配列長の2乗のオーダーの部分配列が存在し，複数の入力配列間でそれらの比較を行うと，少なくとも，与えられた配列の数をさらに乗じたオーダーの計算時間が必要となる。従って，処理すべきデータ量が増加すると，処理に必要な時間が急激に増加し，実用的なシステムを作ることは困難である。その一方で，米国公共機関のデータベースに蓄積されている，cDNAの一部を配列決定して得られたESTと呼ばれる配列の数は，配列決定技術の進歩に伴い急速に増大しつつあり，ヒトだけで450万配列を超えている。こうした膨大な配列データを解析するために，高速な計算方法が必要である。
配列の比較解析を行う従来技術には，以下に述べるものがある。
【０００４】
非特許文献１のホモロジー検索の方法はcDNA配列の解析にも使われている。もし，cDNA配列以外にゲノム配列も利用可能ならば，ゲノム配列の部分文字列で， cDNA配列にも存在する部分文字列がエクソンであると，高い精度で予想することができる。ただし，対応するゲノム配列が無いcDNA配列を解析する場合には，転写産物由来配列間の相互比較を行う必要がある。ホモロジー検索の方法では，総当りで配列間の相互比較を行うと最悪の場合配列数の2乗のオーダーの計算時間が必要となり，配列数が増加するとスプライシングパターンの解析を行うのは，計算量的に困難である。
【０００５】
複数の配列に共通する部分配列の抽出を効率よく行う手法として，非特許文献２のDelcherらの手法が知られている。Delcherらの手法は，主に2つの近縁種のゲノムを比較することを目的としている。この手法では，まず次の(D1)-(D4)の条件を満足するMUM(Maximal Unique Match)と呼ばれる部分配列を探索する。
(D1) MUMは，長さがユーザに指定されたパラメータ以上の配列長をもつ部分配列である。
(D2) MUMは，2つの配列の，共通の部分配列である。
(D3) MUMは，2つのそれぞれの配列内に，１つしか含まれない。
(D4) MUMは，別のMUMの部分配列にならない。すなわち，前記(D1)-(D3)を満足する配列を，可能な限り延長して得られた配列のみがMUMである。
【０００６】
Delcherらの方法をcDNA配列に適用し，MUMとしてアラインされる部分や大きな挿入や欠失が入る位置を同定すれば，MUMはエクソン配列に対応し，大きな挿入や欠失は，選択的なエクソンに対応すると考えることができ，エクソン構造を同定するためのツールとして使うことができる。しかし，Delcherらの方法は２配列のみを対象としており，配列数が3つ以上になったときには適用することができない。また，同一配列のエクソンが同一配列中に現れないことは，一般には保障されない。Delcherらは，非特許文献３において，suffix tree実装のコンパクト化，片方の配列のみをsuffix treeに格納するだけで済む方法の開発を行っているが，３つ以上の配列の同時比較ができない問題は，未解決である。
【０００７】
３つ以上の文字列から共通する部分配列を抽出する方法としては，非特許文献４のchapter 9に記載の，longest common substring problemに対する線形時間解法が知られている。この方法を用いると，任意の数の文字列について，指定された数以上の文字列に共通の部分文字列の中で最も長い部分文字列を，文字総数の線形時間で抽出することができる。ただし，抽出できるのが共通の部分文字列の中で最も長いものだけであり，しかも，互いに重なっている文字列を抽出してしまう場合もあり得るため，この方法で抽出した配列がエクソンに相当する配列と考えることはできない。
【０００８】
Hohlらは，主にゲノム配列のアラインメントに使用することを想定し，３つ以上の配列に共通する部分配列を抽出し，抽出された部分配列を整列し，その結果に基づき複数配列のアラインメントを高速に生成する方法を開発した(非特許文献５)。しかし，Hohlらの方法は，全配列に共通の配列のみを対象としており，一部の配列にのみ共通する配列を抽出することができないため，スプライシング・バリアント配列に適用した場合，選択的スプライシングが行なわれるエクソンを抽出できないという問題がある。
【０００９】
【非特許文献１】Altschul, S.F. et al., Nucleic Acid Research, 25:3389-3402, 1997
【非特許文献２】Delcher et al., Nucleic Acids Research, 1999, 27(11):2369-2376
【非特許文献３】Delcher et al., Nucleic Acids Research, 2002, 30(11):2478-2483
【非特許文献４】Gusfield, D., Algorithms on Strings, Trees, and Sequences: Computer Science and Computational Biology, Cambridge University Press, New York
【非特許文献５】Hohl, M., Kurtz, S., and Ohlebusch, E., Efficient multiple genome alignment, Bioinformatics 18 Suppl.1, S312-S320, 2002
【発明の開示】
【発明が解決しようとする課題】
【００１０】
本発明は，転写産物由来の複数の配列に基づき，これらの配列に含まれるエクソン配列を抽出することを目的としている。
【００１１】
転写産物由来の配列だけに基づき，個々のエクソンの配列を完全に同定できれば理想だが，現実にはそれが常に可能なわけではない。ある遺伝子が4つのエクソンA,B,C,Dを持ち，それらのエクソンはこの順序でゲノム上に並んでいるとする。仮に，図14のように，BとCがスプライシングの過程を経て常に同時にmRNAに現れるならば，cDNA配列だけを見ている限り，B,Cという2つの別のエクソンがあることや，cDNA配列上のどこにB,Cの境界があるかを知ることは不可能である。このほか，図15のように，エクソンAが直後にBかCを必ず伴うとし，BとCには必ずAが先行する場合を考える。BとCの先頭の配列が一致するならば，A,B,Cという３つのエクソンの存在を知ることができても，どこがA,B,Cの境界なのかを知ることはできない。そこで本発明では，常に隣接して共存するエクソンを連結した配列の部分配列であって，他のエクソンに由来する塩基を含まないことが確実な最長の配列を同定することを目指す。こうした部分配列を，以下ではExon Group (EG)と呼ぶ。ところで，例えば数文字程度の，長さが極端に短い部分配列は入力配列上のいたるところに現れ，複数の入力配列に共通する部分配列となり得るが，それらはエクソン配列とみなすべきでない。よって，EGはある程度以上長い部分配列であるべきである。
【００１２】
EGを抽出する方法を述べるために，上記で直感的に説明したEGを，厳密に定義する必要がある。その準備のために，３つの概念を定義する。なお，以下の説明は，塩基配列だけでなくアミノ酸配列など他の種類の文字列を扱う場合でもそのまま当てはまる。そこで，以下ではより広い問題も扱えることを明確化するために，配列の代わりに文字列，塩基の代わりに文字という言葉を用いる。
【００１３】
一つ目の定義は，「文字列が重なる」ことの，本明細書における厳密な定義である。本明細書では，図13に示したように，次の4つの条件を満足する文字列t,t',t''が存在するときに，かつそのときに限り，「文字列s,s'が重なる」という。
(1) t,t',t''が，いずれも，長さが1以上の文字列である。
(2) s=tt'，すなわち，sはtとt'を連結した文字列である。
(3) s'=t't''，すなわち，s'はt'とt''を連結した文字列である。
(4) ある入力文字列が，tt't''を部分文字列として含む。
なお，(1)-(4)が満足されるとき，「sはs'と右で重なる」，「s'はsと左で重なる」という。
【００１４】
二つ目の定義は，文字列sが出現する入力された文字列中での位置の集合Q(s)である。Q(s)は，次の式により定義される。
Q(s)={(i,j)|sは，与えられたi番目の文字列のj番目の塩基から始まる文字列}
また，Q(s)と整数kに対し，以下の演算を定義する。
Q(s)+k={(i,j+k)|(i,j)∈Q(s)}，Q(s)−k=Q(s)+(−k)
なお，Q(s)の要素(i,j)を，本明細書では，sのアピアランスと呼ぶ。
【００１５】
三つ目の定義は，複数のエクソン配列の連結になっている可能性のある文字列maximal match (MM)の定義である。 (M1)-(M4)の条件を満たす文字列を，本明細書では，以下ではMMと呼ぶ。
(M1) 1つ以上の文字列の部分文字列である。
(M2) 長さu以上。ただし，uとは本発明の方法で抽出したい最も短いエクソンの長さで，本発明の方法の利用者がパラメータとして与える整数である。
(M3) 入力された文字列中に出現する任意の文字aに対して，Q(m)≠Q(ma)である。ここに，maは文字列mの末尾に文字aを付加した文字列である。
(M4) 入力された文字列中に出現する任意の文字aに対して，Q(m)≠Q(am)+1である。ここに，amは文字列mの先頭に文字aを付加した文字列である。
【００１６】
図３に，MMの例を示す。
EG 305は，MM 303を部分文字列として完全に含んでいたり，重なっていたりしない文字列であるべきである。そこで，EGを次の条件を満たす文字列と定義する。
(E1) 少なくとも1つの入力文字列の部分文字列。
(E2) 長さu以上。
(E3) MMと重ならない。
(E4) MMを真の部分文字列としない。なお，ある文字列の真の部分文字列とは，その文字列自身以外の部分文字列のことである。
(E5) (E1)-(E4)を満たす文字列の，真の部分文字列にならない。
【００１７】
巨大なデータを効率よく処理するためには，処理時間を可能な限り少なく抑えることが望ましい。入力文字列を読むために，文字総数の線形時間が必要だから，処理時間を線形時間に抑えれば，理論上最も効率のよい計算時間オーダーを実現したことになる。以上を踏まえ，本発明の課題は，複数の生体高分子配列に相当する文字列が与えられたときに，(E1)-(E5)を満足する文字列すなわちEGを，与えられた入力文字列の総文字数の線形時間で抽出する方法を提供することである。
【課題を解決するための手段】
【００１８】
本発明では，入力として複数の文字列301が与えられたとき，後述の右MM 302，後述の左MMまたはMM 303，後述の右EG-holder 304を順次抽出し，それらを利用してEG 305を抽出する。なお，左MMの抽出(ステップS503)と，MMの抽出(ステップS504)は，どちらか片方だけを行えばよい。
【００１９】
以下，本発明の特徴を述べる。
本発明の方法は，複数の入力文字列が与えられたときに，EGを抽出する方法であって，入力文字列の文字数の総和に対し線形時間で処理を完了できる。
【００２０】
本発明の方法は，EGを抽出するために，右MMを抽出する工程と，MMを抽出する工程または前述の左MMを抽出する工程と，右EG-holderを抽出する工程とを有し，それぞれの工程で，入力文字列を格納したsuffix treeを用いることを特徴とする。
【００２１】
本発明の方法はまた，MMを抽出する工程において，右MMを抽出しその中からMMの条件を満足する部分文字列を選択し，抽出することを特徴とする。
【００２２】
本発明の方法はまた，右EG-holderを抽出する工程を有することを特徴とする。
【００２３】
本発明の方法はまた，右EG-holderを抽出する工程において，入力文字列上のMMの位置を解析するために，他の右MMを真のsuffixとしない右MMの位置のみを計算するか，または，他の右MMを真のprefixとしない右MMの位置のみを計算することを特徴とする。
【００２４】
本発明の方法はまた，他の右MMを真のsuffixとしない右MMを求める際，該右MMに相当するsuffix tree上のノードがリーフノードである場合に，該ノードの親ノードのsuffix linkを使用することを特徴とする。
【００２５】
本発明の方法はまた，右EG-holderのprefixのみを，EGの候補とみなすことを特徴とする。
【発明の効果】
【００２６】
本発明によれば，cDNA配列のような生体高分子配列から，EGを，与えられた配列長の総和の線形時間で抽出することができる。
【発明を実施するための最良の形態】
【００２７】
以下，本発明の実施の形態について説明する。最初に，本明細書で使用する記号と用語及び概念を定義する。なお，本発明の方法は，塩基配列だけでなく，アミノ酸配列など他の種類の文字列に対しても適用可能である。そこで，以下では，より一般的な問題を扱うことが可能であることを明確にするために，配列の変わりに文字列，塩基の代わりに文字という言葉を用いる。
【００２８】
はじめに，すでに公知である概念について，本明細書で用いる記号および厳密な意味を定義する。
●空文字列
長さが0の文字列。本明細書では，以下，空文字列をεと表記する。
●文字列の連結
文字列s,tを連結した文字列をstと表記する。
【００２９】
●|s|
文字列sの長さ。
●s[j] (0≦j＜|s|)
文字列sの，j番目の文字。なお，本明細書では，文字列の先頭の文字を０番目と数える。
【００３０】
●s[i..j]
文字列sの，i番目の文字からj番目の文字までの，i番目の文字とj番目の文字をいずれも含む部分文字列。i＞jのときは，s[i..j]=εと定義する。
●Prefix
文字列sについて，s[0..j] (j≦|s|)を，sのprefixという。
【００３１】
●Suffix
文字列sについて，s[j..|s|−1] (j≧0)を，sのsuffixという。
●|A|
集合Aの要素数。
【００３２】
●⊆,⊇,⊂,⊃
本明細書では，2つの集合A,Bについて，「A⊆B」とは，Aの全ての要素がBの要素であることを表し，「A⊇B」とはB⊆Aを表す。「A⊂B」とは，A⊆B かつA≠Bのことであり，「A⊃B」とはB⊂A を表す。
【００３３】
●O(f(n))
関数g(n)に対しg(n)=O(f(n))であるとは，ある定数Cが存在して，十分大きなnに対しg(n)≦Cf(n)が成立することである。また，ある量が「O(f(n))である」とは，その量がnの関数g(n)以下であり，g(n)=O(f(n))であることを意味する。
【００３４】
●Ω(f(n))
g(n)= Ω(f(n)) であるとは，ある定数Cが存在して，十分大きなnに対しg(n)≧Cf(n)が成立することである。また，ある量が「Ω(f(n))である」とは，その量がnの関数g(n)以上であり，g(n)= Ω(f(n))であることを意味する。
【００３５】
●Suffix tree
文字列P1,P2,…,Pnのsuffix tree とは，P1,…,Pn に存在しない文字$について，$を各文字列の末尾に付加して得られるP1$,…,Pn$に現れる全てのsuffix を格納したtree状のデータ構造で，以下の性質をもつ。なお，2つの文字列「ATATG」，「TTAGTA」を格納したsuffix tree 401を図4に図示した。
(S1) ルートノード402を持つ有向木である。
(S2) P1$, …, Pn$に存在する互いに異なるsuffixの数に等しい数のリーフを持つ。本発明では，それらのリーフには，文字列番号iと文字列中の位置jのペア(i,j)の集合{(i,j)}が，対応付けられ，その集合が格納されているデータ構造に，そのデータ構造を指すポインタ等により，O(1)でアクセスが可能であることを前提とする。
(S3) それぞれのエッジ404には，P1$,..., Pn$のうち，どれか１つの文字列の部分文字列がラベル405として付加されている。
(S4) 任意のノード407を起点とするエッジに，同じ文字で始まるラベル405を持つエッジのペアは存在しない。
(S5) ルートノード402から文字列番号iと文字列中の位置jが割り当てられたリーフ403へのパス上にあるエッジのラベルを，このパス上で出会う順に連結すると，入力文字列iのj番目の文字で始まるsuffixの末尾に$を付加した文字列となる。
(S6) 任意のノードvを起点として伸びるエッジの数は1でない。
【００３６】
Suffix treeに格納されている文字列の文字数の合計をnとするとき，リーフ403の数がn以下なので，リーフ以外のノードの数はn-1個以下である。従って，suffix tree全体のノード数は2n-1個以下で，O(n)である。なお，ルートノードから，ラベルが’$’であるエッジを一回通って到達できるリーフには，空集合φが対応付けられる。
【００３７】
●パスラベル
Suffix tree中のノード407のパスラベルとは，ルートノード402からこのノード407にいたるパス408上のエッジのラベルを，連結して得られる文字列である。
●Suffix link
vを，suffix treeの，リーフでないノードとする。ノードvのパスラベルが，ある文字aとある文字列sによりasと表現できるとき，vからパスラベルをsとするノードへのポインタをsuffix linkと呼ぶ。本明細書では，vのsuffix linkを， suffixlink(v)と表記する。vのパスラベルがasならば，suffixlink(v)が指し示すべきノードv(s)は必ず存在することが知られている(非特許文献４, chapter 6)。図4の409は，suffix linkの例である。なお，誤解を招く恐れのないときは，「suffixlink(v)」という表記を，suffixlink(v)が指すノードと同一視する場合がある。
【００３８】
次に，本明細書で特に用いる記号・用語および，本明細書固有の概念について，以下で定義する。
●入力文字列（図３の301参照）
本発明の方法に入力として与えられる，塩基配列やアミノ酸配列などの，生体高分子配列。
●N
本発明の方法を適用すべき入力文字列の文字数の総和。
●k
本発明の方法を適用すべき入力文字列の数。
●u
本発明の方法に与えられるパラメータで，エクソン配列と考えられる最も短い文字列の長さを指定する整数である。
【００３９】
●p(v)
ノードvのパスラベル。
●p’(v)
ノードvのパスラベルp(v)の末尾に$があれば，その$を除去して得られる文字列。それ以外の場合は，p(v)と同一の文字列。
【００４０】
●v(s)
パスラベルが文字列sのノード。
●v’(s)
sが出現する任意の入力文字列Si上で，sがSiのsuffixとしてのみ出現する場合は，v(s$)。それ以外の場合，v(s)。
【００４１】
●Σ
入力文字列中に存在するすべての文字種の集合。
●Si
i番目の入力文字列。ただし，0≦i＜kである。
【００４２】
●MM（図３の303参照）
上述の条件(M1)-(M4)を満たす文字列。
●EG（図３の305参照）
上述の条件(E1)-(E5)を満たす文字列。
●右MM（図３の302参照）
【００４３】
次の条件(R1)-(R3)を満たす文字列を，右MMという。
(R1) 1つ以上の入力文字列の部分文字列である。
(R2) 長さはu以上である。
(R3) 右MM rは，任意の文字aに対してQ(r)≠Q(ra)を満たす。すなわち，右MMは，右側へ延長することはできない。
【００４４】
なお，任意の右MM rに対し，Q(sr)+|s|=Q(r)を満たす空文字列かも知れない最長の文字列sをとることができ，srは(M1)-(M4)を満たすので，MMである。従って，次の性質(R4)が成り立つ。
(R4) 任意の右MM rについて，rをsuffixとするMM mで，Q(m)+|m|-|r|=Q(r)を満たすものが存在する。
【００４５】
また，右MM rがあるMM mのprefixであると仮定する。このとき，任意の文字aについてQ(am)+1⊂Q(m)より，mのあるアピアランス(i,j)について，(i,j-1)はamのアピアランスではない。従って，Si[j-1]≠aだから，(i,j-1)はarのアピアランスではない。ところが，(i,j)はmのアピアランスだから，そのprefixであるrのアピアランスでもある。従ってQ(ar)+1≠Q(r)。よって，rは(M4)を満たすので，MMである。つまり，右MMについて次の性質(R5)が成り立つ。
(R5) MMのprefixである右MMは，MMである。
【００４６】
●左MM
次の条件(L1)-(L3)を満たす文字列を，左MMという。
(L1) 1つ以上の入力文字列の部分文字列である。
(L2) 長さはu以上である。
(L3) 左MM lは，任意の文字aに対してQ(l)≠Q(al)+1 を満たす。すなわち，左MMは，左側へ延長することはできない。
【００４７】
●右EG-holder（図３の304参照）
次の条件(H1)-(H3)を満たす文字列を，右EG-holderという。
(H1) 右MMである。
(H2) 次の(H2a)，(H2b)のうち，どちらかの条件を満たす。
(H2a) MMである。
(H2b) あるMM mが存在し，mhが少なくとも１つの入力文字列の部分文字列となるような，文字列hである。
(H3) 右MMを，真のprefixとして持たない。
【００４８】
上述したように，本発明では右MM，左MM又はMM，右EG-holderを順次抽出し，それらの抽出結果を利用してEGを抽出する。図５は，本発明による処理の全体を示すフローチャートである。以下では，図５を参照して，右MM，左MM，MM，右EG-holder，EGを抽出する方法を説明する。
【００４９】
■Suffix tree Tの構築(ステップS501)
本発明の方法では，まず，入力文字列に基づきsuffix tree Tを構築する。Suffix tree T は，Ukkonen のアルゴリズム(非特許文献４, chapter 6)を用いて，O(N)の処理時間で，全suffix linkを設定する工程を含めて構築可能である。本発明の方法は，以降，suffix tree T上の深さ優先探査を繰り返し実行する。
【００５０】
■右MMの抽出(ステップS502)
Suffix tree Tを深さ優先探査し，ノードvに出会ったとき，文字列p’(v)の長さがu以上で，かつvの親ノードからvに至るエッジのラベルが「$」でないとき，vに「右MMである」と印をつける。なお，パスラベル長|p(v)|は，深さ優先探査の際にルートノードからvまでのパスにあるエッジ長の和を保持することで，1ノード当たり定数時間で計算でき，最後に通過したエッジのラベルの末尾が$のときは|p’(v)|=|p(v)|−1，それ以外のときは|p’(v)|=|p(v)|である。
【００５１】
この処理が完了した時点で，「右MMである」と印がついたノードvについてp’(v)は右MMであり，かつそれらだけが右MMである。その理由を述べる。
(R1)について:
suffix tree T の任意のノードのパスラベルp(v)は，少なくとも１つの入力文字列Siの末尾に$を付加した文字列Si$のsuffixのprefixである。従って，p(v)の末尾に$が存在すればその$を除去した文字列p’(v)は，Siのsuffixのprefixであり，少なくとも１つの文字列の部分文字列である。よって，右MMであると印がついたノードvについて，p’(v)は，(R1)を満たす。
(R2)について:
本発明の方法で「右MMである」と印がつくノードvは，|p’(v)|≧uなるノードだけである。従って，こうしたノードvについて，p’(v)は(R2)を満たす。
(R3)について:
「右MM」であると印がついたノードvについて，どのように文字a∈Σを選んでも，p’(v)のあるアピアランス(i,j)が存在し，(i,j)はp’(v)aのアピアランスではない。なぜなら，そのようなアピアランス(i,j)が存在しないなら，ノードvから出るエッジのラベルの最初の文字は全てaであり，suffix treeの条件(S4)に反するからである。従って，p’(v)は(R3)を満足する。
【００５２】
さらに，本発明の方法で右MMの印がついたノードvのp’(v)でなければ，(R1)-(R3)のどれかの条件を満足しないことを示す。
【００５３】
まず，suffix tree Tのノードのパスラベルのprefixとして現れない文字列は，どの入力文字列の部分文字列でもないから，(R1)に反する。
【００５４】
次に，ノードのパスラベルであって$を含まないprefixとなっている文字列rについて，ノードv’(r)が存在しないと仮定する。rが入力文字列のsuffixとしてのみ出現するならば，v’(r)=v(r$)が存在し仮定に反するため，rはある入力文字列Siのsuffix以外の位置に出現する。すなわち，あるアピアランス(i,j)∈Q(r)についてj+|r|＜|Si|である。このとき，a=Si[j+|r|]とおく。aはSi中の文字だからa∈Σである。rがある入力文字列Siのsuffixである場合，raをprefixとするパスラベルと，r$をprefixとするパスラベルが両方存在するため，v’(r)=v(r)が存在するから仮定に矛盾する。よってrはどの入力文字列のsuffixでもない。aでない文字bについて，rbを部分文字列とする入力文字列の存在を仮定すると，raをprefixとするパスラベルと，rbをprefixとするパスラベルが両方存在するため，v’(r)=v(r)が存在するから仮定に矛盾し，rbを部分文字列とする入力文字列の存在の存在は否定される。従って，rは入力文字列中で右側に常に文字a∈Σを伴って現れるからQ(r)=Q(ra)である。よって，rは(R3)を満足しない。
【００５５】
さらに，本発明の方法は，|p’(v)|＜uなるノードvには印をつけないが，このようなノードvについて文字列p’(v)は (R2)を満足しない。
【００５６】
以上で検討した以外に，「右MMである」と印がついたノードvのp’(v)以外の文字列は存在しない。つまり，本発明の方法で，すべての右MM rに対して，v’(r)に，過不足なく「右MMである」と印をつけることができる。
【００５７】
■左MMの抽出(ステップS503)
ステップS501でsuffix tree Tを構築し，ステップS502で右MMを抽出したのと同様の方法で，左MMを抽出する。なお，ステップS504によりMMを抽出する場合，左MMを抽出するステップS503は行わなくてもよい。
【００５８】
ステップS503では，左MMを抽出するために，すべての入力文字列を逆順に見た文字列について，suffix tree T'を構築する。例えば，入力文字列がATATGとTTAGTAの場合は，これらを逆順にして得られるGTATAとATGATTについて，suffix tree T'を構築する。そして，「suffix tree T」を「suffix tree T'」に，「右MM」を「左MM」に読み替えて，本発明の右MM抽出方法を実行する。この方法で，左MMを過不足無く抽出できるのは，ステップS502で過不足なく右MMが抽出されることから明らかである。
【００５９】
■MMの抽出(ステップS504)
MMを抽出する方法について述べる。なお，ステップS503により左MMの抽出を行った場合，MMを抽出するステップS504は行わなくてもよい。
【００６０】
全てのMMは，(M1)-(M3)より(R1)-(R3)を満たし，右MMである。任意の右MM rについてv’(r)に既に上記の方法で「右MMである」と印がついているので，これらのrが (M4)を満足するときv’(r)に「MMである」と印をつければよい。(M4)を満足する右MMに相当するノードを同定するには，非特許文献４のP144〜P145に記載の“left diverse node”同定方法を，次のように修正して適用すればよい。
(1) 非特許文献４における“left diverse”であるノードvの定義を，「p’(v)が，(M4)を満足するノードv」に変更する。
(2) リーフノードvについては，非特許文献４の方法のように，「リーフは常に”left diverse”でない」とは判断しない。代わりに，ある文字a∈Σが存在し，p’(v)の任意のアピアランス(i,j)についてj≧1かつSi[j-1]=aが成立すれば，vの”left character”をaとし，このようなaが存在しなければ，vは”left diverse”であると判断する。
(3) リーフでないノードvについては，vの子ノードに”left diverse”であるノードがなく，かつ全ての子ノードの”left character”があるa∈Σに一意に決まるならば，v自身も”left diverse”ではなくvの”left character”をaとし，それ以外の場合は，vは”left diverse”であるとする。
(4) “left diverse”であると判断されたノードに「右MMである」と印があれば，「MMである」と印をつける。
【００６１】
この方法の計算時間は，リーフに対応付けられているアピアランスの総数と，suffix tree Tのノード数のオーダーである。前者は，入力文字列長の総和Nに等しく，後者は2N-1以下だから，全体の処理時間はO(N)である。
【００６２】
■MM末端位置の抽出(ステップS505)
ステップS502で抽出した右MMの情報と，ステップS503で抽出した左MMの情報またはステップS504で抽出したMMの情報を用いて，ステップS505で2次元配列Head，Tailを以下のように設定する。なお，Head[i,j], Tail[i,j]は，0≦i＜k，0≦j＜|Si|なる(i,j)に対して値を持つ。
Head[i,j] :
(i,j)をアピアランスとするMMが存在すれば1。それ以外の場合は0。（図７参照）
Tail[i,j] :
(i,j-|m|+1)をアピアランスとするMM mが存在すれば1。それ以外の場合は0。（図７参照）
【００６３】
まず，Tailの値を設定する方法を２つ述べる。一つ目の方法は，以下の(t1)-(t3)を実行することである。
(t1) 2次元配列Tail 702の全ての要素を0で初期化する。
(t2) リストminimal_right_MMを空リストで初期化する。
(t3) suffix tree Tを深さ優先探査し，出会ったノードvに対し(t3a)-(t3c)の処理を行なう。
(t3a) p’(v)=asなる文字a∈Σと文字列sについて，v’(s)を同定する。vがリーフでない場合，sは$を含まないのでv’(s)=v(s)=suffixlink(v)である。vがリーフの場合は，次のような処理を行う。まず，深さ優先探索時に常に，その時点で処理対象となっているノードの親ノードを記憶しておく。リーフvの親ノードをwとするとき，wがリーフでないことは明らかである。w’=suffixlink(w)とおく。w-v間のエッジのラベルである文字列に沿って，w’からsuffix tree上を移動すると，ノードw’’=v(s)に到達する。w’’に直接至るエッジのラベルが一文字の文字列”$”であればv(s)はw’’の親ノード，それ以外の場合はv’(s)=w’’である。
(t3b) v’(s)に「右MMである」と印がついていない場合，リストminimal_right_MMの末尾に，|p’(v)|を追加する。なお，リストminimal_right_MMに追加した値は，ノードvをルートノードとする部分木の処理が終了した時点で削除する。この削除の操作は，単にminimal_right_MMの末尾の要素を削除するだけで実現できる。
(t3c) vがリーフのとき，vに対応付けられた任意のアピアランス(i,j)およびリストminimal_right_MMの任意の要素kについて，Tail[i,j+k-1] 702を1に設定する。
【００６４】
この方法では，Tailの値を設定するにあたり，MMの位置を直接計算する代わりに，他の右MMをsuffixとして持たない右MM(以下，極小右MMと呼ぶ)である文字列rの終了位置を計算している。すなわち，任意の極小右MM rに対し，rの任意のアピアランス(i,j)について，Tail[i,j+|r|-1] 702に，1を書き込む。全ての極小右MMの終了位置の集合と，全てのMMの終了位置の集合は等しい。その理由は，全ての右MM rには(R4)よりrをsuffix としQ(r)+|r|=Q(m)+|m|であるMM mが存在する一方，全てのMMは極小右MMをsuffixとして持つからである。したがって，上記(t1)-(t3)の方法で，全てのMMの終了位置に対応するTail[i,j] 702を1に設定することができる。
【００６５】
この方法(t1)-(t3)を処理するために要する計算時間について検討する。まず，(t1)の処理時間は入力文字列長の総和Nに対しO(N)である。(t2)の処理時間はO(1)。(t3)で， (t3a)，(t3c)以外の処理は，Tのノード数がO(N)で各ノードについての処理時間がO(1)だから，全体でO(N)である。一方，同一の終了位置を持つ異なる２つの極小右MMが存在すると仮定すれば，長い方が短い方をsuffixとするため，２つとも極小右MMとの仮定に矛盾。よって，同一の終了位置を持つ極小右MMはひとつであり，(t3c)の処理で，同一の(i,j)に対してTail[i,j]が2度以上更新されることはない。Tailの要素数はNだから，(t3c)の処理時間もO(N)である。
【００６６】
(t3a)の処理に要する処理時間の総和が，O(N)であることを示す。vがリーフでない場合，v’(s)=suffixlink(v)はO(1)で同定できるから，リーフでないノード全体ではO(N)である。vがリーフの場合，w’’からv’(s)を同定する処理は各リーフについてO(1)で，全体でO(N)である。w’’を同定する処理時間について検討する。リーフvからリーフw’’を同定する処理は，各リーフvについて1度ずつ行われるが，その処理時間の合計は，すべての0≦i＜kについて，文字列Si$の最初のsuffixから|Si|番目のsuffixまで順番にリーフを同定していく時間以下である。そこで，後者の処理時間の上限を考察する。変数xを，suffix tree T上のノードを表す変数とする。そして，Siのj番目のsuffixに$を付与した文字列をパスラベルとするリーフvから，j+1番目のsuffixに$を付与した文字列をパスラベルとするリーフw’’へと辿る過程を考える。xの値は，注目しているノードに順次更新されるとする。なお，w’からw’’までxが変化する過程を，以下ではダウンウォークと呼ぶ。xをvからwに更新するために要する処理時間は，深さ優先探索の過程で常に親ノードを記憶しているからO(1)であり，wからw’への更新もsuffixlinkを用いればO(1)である。w’からw’’へのパスをたどる時間は，ノードを1つ通過するごとにO(1)だから，このパス上のノード数nに対しO(n)である。ダウンウォークに要する処理時間の，文字列Siについて処理全体の和をとると，O(|Si|)であることを述べる。ルートノードからxへ至るパス上にあるノードの数（ルートノードおよびxを含む）を，xのnode depthと定義し，その増減について考察する。xのnode depthは，xがw’からw’’へ移る過程で，ノードをひとつ通過するごとに１増加する。一方，xがリーフvからその親ノードwへ更新されるときには1減少する。さらに，xがw’へ更新されるときに，node depthの減少量の最大値は1である。なぜなら，一般にノードv1,v2がv2=suffixlink(v1)の関係にあるとき|p(v1)|=|p(v2)|+1であり，ルートノードからwへのパス上にパスラベル長が等しいノードは存在しないので，ルートノードからwに至るパス上のルートノード以外の任意のノードは，ルートノードからw’に至るパス上の互いに異なるノードへsuffix linkを持つから，w’のnode depthは(wのnode depth)-1以上となるためである。従って，Siのsuffixに相当するリーフを順に辿る過程で，xのnode depthは最大で2×|Si|減少する。一方，Siの末尾に$を付加した文字列のsuffixの長さの最大値は|Si|+1だから，Si$のsuffixに相当するリーフのnode depthは常に|Si|+2以下である。従って，Siに関する全ダウンウォークにおいてxが更新される回数の和をnとすれば，n−2×|Si|≦|Si|+2よりn≦3|Si|+2=O(|Si|)である。従って，|Si|のsuffixに相当するリーフを辿るために必要な処理時間が，O(|Si|)であることがわかった。このことから，全リーフを辿るために要する処理時間がO(Σ|Si|)=O(N)であることが直ちにわかる。
【００６７】
二つ目の方法は，Tailの値を設定するにあたり，MMの位置を直接計算する代わりに，他の右MMをprefixとして持たない右MM rの位置を計算し，そのような任意のrについて，(i,j)がrのアピアランスなら，Tail[i,j+|r|-1] 702に，1を書き込むことである。このような，他の右MMをprefixとしない右MM r (以下，prefix極小右MMと呼ぶ)の終端の集合と，全てのMMの終端の集合は等しい。その理由を述べる。まず(R4)より，全ての右MM rには(R5)よりrをsuffix としQ(r)=Q(m)+|m|-|r|であるMM mが存在する。従って，次式(1)が成り立ち，prefix極小右MMの終端全体の集合は，MMの終端の集合に含まれる。
【００６８】
【数１】
2005316605.tif 000003
【００６９】
ただし，式(1)で，pmrMMとは，「prefix極小右MM」の集合を表す。一方，mの長さu以上の任意のsuffix sは右MMである。なぜなら，sはmの部分文字列だから(R1)を満たし，長さu以上だから(R2)を満たし，次の議論により(R3)も満たす。まず，mがある入力文字列のsuffixであるとき，sもその文字列のsuffixであり，sが(R3)を満たすのは明らか。mがどの文字列のsuffixでもないとき，異なる文字a,bについて，それぞれあるアピアランス(i,j), (i’,j’)が存在し，(i,j)∈Q(m), Si[j+|m|]=a, (i’,j’)∈Q(m), Si’[j’+|m|]=bである。すなわち，mは右側にaを伴って現れる場合と，bを伴って現れる場合がある。sはmの|m|-|s|番目の文字から始まるsuffixだから(i,j+|m|-|s|)∈Q(s), Si[(j+|m|-|s|)+|s]]=Si[j+|m|]=a, (i’,j’+|m|-|s|)∈Q(s), Si’[(j’+|m|-|s|)+|s]]=Si’[j’+|m|]=bが成り立つ。すなわち，sも右側にaを伴って現れる場合と，bを伴って現れる場合があるため，任意の文字a∈Σに対しQ(s)≠Q(sa)だから，sは(R3)を満たす。以上に加え，MM mはそれ自身が長さu以上の右MMであるからmには長さu以上のsuffix sは必ず存在し，そのうち長さがちょうどuのsuffix rは上記の議論から右MMであるが，そのprefixの長さはu-1以下だから，右MMを真のprefixに持たない。すなわち，rはprefix極小右MMであり，すべてのMMはprefix極小MMをsuffixに持つことがわかる。従って，式(2)が成り立つ。すなわち，prefix極小右MMの終端全体の集合は，MMの終端の集合を含む。
【００７０】
【数２】
2005316605.tif 000004
【００７１】
式(1)，式(2)から式(3)が成立し，prefix極小右MMの終端全体の集合は，MMの終端の集合に等しいことがわかる。
【００７２】
【数３】
2005316605.tif 000005
【００７３】
以上を踏まえ，次の(t’1)-(t’3)の方法で，全ての右MMの末端に対応するTail[i,j] 702を1に設定することができる。
(t’1) 2次元配列Tail 701の全ての要素を0で初期化。
(t’2) 変数kを-1で初期化する。
(t’3) suffix tree Tを深さ優先探査し，出会ったノードvに対し(t’3a)，(t’3b)の処理を行う。
(t’3a) vに「右MMである」と印がついており，かつk<0の場合，kに|p’(v)|を代入する。kは，vのvを含む部分木の深さ優先探査が終わった時点で，−1に戻す。
(t’3b) vがリーフのとき， k>0ならば，vに対応付けられている任意のアピアランス(i,j)についてTail[i,j+k-1] 701を1に設定する。
【００７４】
この方法により，入力文字列iのj番目の文字で終了するprefix極小右MM rが存在するときに，かつそのときに限り，Tail[i,j+|r|-1]に1を書き込む。逆に，このようなi,j以外では，Tail[i,j]は1にならない。深さ優先探査で，同一のリーフを2回以上処理することはないから，この方法で，同一の(i,j)に対し2回以上Tail[i,j]が更新されることはない。ゆえに，(t’1)-(t’3)はO(N)で完了できる。
【００７５】
なお，極小右MMやprefix極小右MMの位置のみを計算するのではなく，全てのMM mについて，それらが出現するすべてアピアランス(i,j)についてTail[i,j+|m|-1]=1とする方法では，O(N)の処理時間では完了できない場合がある。Nが文字列数kで割り切れ，かつN≧3k(u+k)と仮定し，文字列集合{Si=T^(u+i)AT^(N/k-(u+i+1))} (0≦i＜k)について考察する（図６）。ただし，「a^n」は文字aのn回の繰り返しを表すものとする。この文字列集合において，計k-1個の，Tがu〜u+k-2個連続する文字列T^u, ..., T^(u+k-2)のそれぞれが，Si[0..u+i-1]=T^(u+i), Si[u+i]=A, S(i+1)[0..u+i-1]=T^(u+i), S(i+1)[u+i]=Tより，(R1)-(R3)を満たすので，これらは右MMである。ただし，S(i+1)はi+1番目の入力文字列を表す。ここで，Sj(0≦j＜k)ではSj[u+j]=Aの後にN/k-(u+j+1)個のTが並ぶことから，T^(u+i)はSjにN/k-(u+j+1)-(u+i)+1=N/k-2u-i-j 回以上出現する(N≧3k(u+k),i≦k,j≦kより，N/k-2u-i-j≧k(u+k)≧0)。従って，考察している文字列集合全体でのMMのアピアランスの数は，式(4)より(k-1)N/3以上だから，Ω(Nk)である。以上により，すべてのMMのアピアランスを列挙して終端を探索する方法では，O(N)の処理時間でTailの値を設定できるとは限らないことがわかる。
【００７６】
【数４】
2005316605.tif 000006
【００７７】
次に，Head[i,j] 701を設定する方法を示す。まず，ステップS503で，左MMの抽出が完了している場合に， Head[i,j]を設定するためには，suffix tree Tをsuffix tree T'，Tail[i,j]をHead[i,|Si|-1-j]と読み替えて，Tailを設定する上記の方法を適用すればよい。この方法で2次元配列Headを正しく設定できることは，上記の方法がTailを正しく設定することから自明であり，O(N)の時間で処理を完了できる。
【００７８】
一方，ステップS503ではなくステップS504を実行した場合には，suffix tree Tを深さ優先探査し，「MMである」と印がついたノードvで，ルートノードからvに至るパス上に別の「MMである」と印のついたノードが存在しないノードvを発見した場合，vの部分木内では，出会った任意のリーフに対応付けられている任意のアピアランス(i,j)について，Head[i,j] 701を1に設定すればよい。この方法で，同一のアピアランス(i,j)に対し，2回以上Head[i,j] 701が更新されることはないから，O(N)の処理時間で完了できる。ところで，この方法により他のMM mを真のprefixとして持つMM m'の位置は計算されないが，mはm'のprefixでありQ(m)⊇Q(m')だから，m'が入力文字列iのj番目の文字で始まる部分文字列ならば，mもそうであり，Head[i,j]は正しく1に設定される。この方法で，入力文字列iのj番目から始まる任意の文字列がMMでない場合に，Head[i,j]が1にならないことは明らかである。
【００７９】
■右EG-holderの抽出(ステップS506)
以下に述べる方法(h1)-(h3)によって，p’(v)が右EG-holderであるノードvに「右EG-holderである」と印をつけることができる。 (h1)-(h2)で2次元配列Tvの値を設定し，(h3)で「右EG-holderである」と印をつける。
(h1) 2次元配列Tvを用意し，Tvの全要素をヌルポインタで初期化する。なお，TvもHead，Tailと同様に，入力文字列iと該文字列の長さよりも小さな非負整数jに対して値が定義される。さらに，変数shortest_right_MMを用意し，ヌルポインタを代入する。
(h2) suffix tree Tを深さ優先探査し，出会ったノードvについて，(h2a)と(h2b)の処理を行う。
(h2a) ノードvに「右MMである」と印があり，shortest_right_MMがヌルポインタであれば，shortest_right_MMにvへのポインタを代入する。なお，vをルートノードとする部分木の処理が終了した段階で，shortest_right_MMの値はヌルポインタに戻す。
(h2b) ノードvがリーフであった場合には，このリーフのパスラベルに対応付けられている任意のアピアランス(i,j)について，変数shortest_right_MMの値をTv[i,j]に代入する。
(h3) (i,j) (0≦i＜k，0≦j＜|Si|)が，次の(h3a)，(h3b)の条件がともに満たすとき，Tv[i,j]が指し示すノードに「右EG-holder」と印を付ける。
(h3a) Head[i,j]が1であるか，j≧1かつTail[i,j-1]が1。
(h3b) Tv[i,j]がヌルポインタでない。
【００８０】
この方法(h1)-(h3)で，p’(v)が右EG-holderであるノードvに，かつそれらのノードvに限り，「右EG-holderである」と印がつくことを示す。まず，(h2)によって，Tv[i,j]の値は，以下の値に設定される。
(1) (i,j)をアピアランスとする右MMが１つ以上存在するなら，それらの中でprefixに他の右MMをもたない右MMをrとするとき，v’(r)のポインタ。
(2) (i,j)をアピアランスとする右MMが存在しないなら，ヌルポインタ。
【００８１】
これを踏まえ，(h1)-(h3)の処理により，「右EG-holderである」と印がついたノードvについて，p’(v)は右EG-holderであることを示す。
(H1)について:
「右MMである」と印があるノードにしか「右EG-holderである」との印をつけないので，p’(v)は(H1)を満たす。
(H3)について:
ノードvに「右EG-holderである」と印がつくならば，p’(v)のあるアピアランス(i,j)について，Tv[i,j]=v’(p(v))。ここで，p’(v)の真のprefixで，右MMである文字列rの存在を仮定すると，Tv[i,j]の定義に矛盾。よって，p’(v)が(H3)を満足することがわかる。
(H2)について:
本発明の方法は，Head[i,j]が1であるか，j≧1かつTail[i,j-1]が1であるようなTv[i,j]の指すノードにしか「右EG-holderである」との印をつけない。従って，「右EG-holderである」との印がついたノードvに対して，p’(v)は(i,j)をアピアランスとするMM mのprefixである右MMか，入力文字列SiにおいてあるMMの次の文字から始まる右MMである。前者の場合，rは(R5)よりMMである。従って，(H2)も満たされる。
【００８２】
次に，ある右EG-holder hについて，v’(h)に右EG-holderの印がつかないと仮定する。まず，Tv[i,j]にv’(h)へのポインタが入らないと仮定する。このとき，hは右MMでないか，hの真のprefixである右MMが存在する。それぞれ，(H1)，(H3)よりhが右EG-holderであることに矛盾しTv[i,j]にはv’(h)へのポインタが入らなければならないことがわかる。それでもhに「右EG-holderである」との印がつかないということは，hの任意のアピアランス(i,j)について，Head[i,j]=0 かつTail[i,j-1]=0(j≧1のとき)ということである。このとき，任意のMM mについて，(i,j)と(i,j-|m|)はいずれもmのアピアランスではないが，そのようなhは(H2a)，(H2b)をいずれも満たさないため，やはりhが右EG-holderであることに矛盾する。従って，v’(h)に「右EG-holderである」と印がつかないと仮定したことが誤りであり，本発明の方法は右EG-holderに相当するノードに過不足なく「右EG-holderである」と印をつけられることがわかる。
【００８３】
■EGの抽出(ステップS507)
EGは，次の(e1)-(e3)の方法によって，過不足なく抽出することができる。その理由は後述する。本発明の方法は，EGを抽出するために，これまでのステップでで計算したHead, Tail, Tvの値を利用する。
(e1) suffix tree T内の全てのノードvについて，変数EGlength(v)を設ける。そして，深さ優先探査を行ない，EGlength(v)=|p’(v)|と初期化する。
(e2) 全ての入力文字列Siについて，次の(e2a)-(e2d)までの処理を行う。
(e2a) 変数c,jをc=1,j=|Si|-1で初期化する。
(e2b) Tv[i,j]がヌルポインタでなければ，Tv[i,j]が指し示すノードをvについて変数EGlength(v)の値を，cの値703とそれまでのEGlength(v)の値のうち，小さい方の値で置き換える。
(e2c) Head[i,j]が1であるか，j≧1かつTail[i,j-1]が1であれば，cに1を代入し，それ以外のときはcに1を加える。
(e2d) jの値から1を減じる。減じた後でjが0以上ならば，(e2b)-(e2d)までの処理を再び行う。jが0より小さければ，入力文字列Siに関する処理は終了である。
(e3) Suffix tree T上を深さ優先探査し，ノードvが，次の条件(e3a)，(e3b)をともに満たすとき，vに「EGである」と印をつける。
(e3a) vに「右EG-holderである」と印がついている。
(e3b) EGlength(v)≧u
【００８４】
以上の方法で，「EGである」と印がついたノードvの長さEGlength(v)のprefixはEGであり，かつそれらだけがEGであることを説明する。その説明に先立ち，(E1)-(E5)の条件が次の(E’1)-(E’5)の条件と等価であることを示す。
(E’1) ある右EG-holderのprefixである。
(E’2) 長さがu以上である。
(E’3) 右でも左でも，MMと重ならない。
(E’4) MMを真の部分文字列としない。
(E’5) (E’1)-(E’4)を満足する文字列の真のprefixにならない。
【００８５】
(E1)-(E5)の条件が(E’1)-(E’5)と等価であることを示す準備として，(E’1)と(E’2)を満足する文字列eが，次の(E’6)を満足することを示す。
(E’6) eをprefixとし，Q(e)=Q(h)なるEG-holder hが，ただひとつ存在する。
【００８６】
(E’6)が成り立つことを示すために，互いに異なる右EG-holder h,h’が存在すると仮定する。esがh,h’の共通のprefixである最長の文字列となる文字列sについて検討する。hとh’は互いに異なるから，少なくとも一方は|es|より長い。一般性を失わず，|h|＞|es|と仮定できる。h中で，prefix esの次の文字をaとする。esがhとh’に共通する最長のprefixだから，esaはh’のprefixでない。よってQ(esa)∩Q(h’)=φだが，Q(es)⊇Q(h’)≠φよりQ(esa)≠Q(es)。一方，a以外の任意の文字bについて，Q(esb)∩Q(h)＝φだが，Q(es)⊇Q(h)≠φよりQ(esb)≠Q(es)。従って，esは(R3)を満たす。Q(es)⊇Q(h)≠φよりesは(R1)を満たし，|es|≧|e|≧uより，(R2)も満たす。従って，esは右MMである。|h|＞|es|より|es|は|h|の真のprefixだが，これはhが右EG-holderであり(H3)を満たすことに矛盾する。ゆえに，互いに異なる右EG-holder h,h’の存在を仮定したのが誤りで，eをprefixとする右EG-holderは一意である。
【００８７】
一方，eは(E’1)を満たすからある右EG-holder hのprefixであり，Q(e)⊇Q(h)である。ここで，Q(e)⊃Q(h)を仮定し，ある文字列sをQ(e)=Q(es)⊃Q(h)なる最長の文字列とする。このとき，esはQ(es)=Q(e)≠φより(R1)を，|es|≧|e|≧uより(R2)を満たし，sの定義より(R3)も満たすから右MMである。Q(es)⊃Q(h)より|es|＞|h|だからesはhの真のprefixだが，hが右EG-holderであり(H3)を満足することに矛盾する。従って，Q(e)⊃Q(h)の仮定が誤りでQ(e)⊆Q(h)。eはhのprefixだからQ(e)⊇Q(h)より，Q(e)=Q(h)。以上で，(E’6)が正しいことがわかる。
【００８８】
まず，(E’1)-(E’5)を満たす文字列eが，(E1)-(E5)を満たすことを示す。
文字列eは(E’1)-(E’4)を満たすので， (E1)-(E4)も満たすことは明らか。
【００８９】
(E5)も満たすことを示す。e’=ses’≠eが，(E1)-(E4)を満たすと仮定し，矛盾を導く。eには(E’6)より，eをprefixとする右EG-holder hが一意に存在し，hは(H1)より右MMだから，(R4)より，Q(h)+|h|=Q(m)+|m|なるMM mでhをsuffixとするものが存在。
【００９０】
ここで，|es’|＞|h|を仮定すると，Q(es’)⊆Q(e)=Q(h)=Q(m)+|m|-|h|より，ses’が存在する入力文字列上では，図９に示すようにes’のprefixとしてhが存在し，そのhをsuffixとするmが存在する。|es’|＞|h|だから，e’=ses’はmと左で重なり(E3)に反するか，mを真の部分文字列として含み(E4)に反し，矛盾が生じる。従って，|es’|≦|h|であり，Q(es’)⊆Q(e)=Q(h)よりes’はhのprefixである。さらに，Q(e)⊇Q(es’)⊇Q(h)=Q(e)よりQ(es’)=Q(h)である。
【００９１】
次に，|s|≧1を仮定する。このとき，hがMMなら，Q(e’)+|s|⊆Q(es’)=Q(h)よりある入力文字列Si上に部分文字列shが存在するから，e’=ses’は図10に示すように|es’|＜|h|ならばMM hと右で重なり，|es’|=|h|ならばhを真の部分文字列として持つため，e’が(E3)と(E4)を満たすことに矛盾するから，hはMMでない。しかし，hは (H2)を満たすから，(H2b)を満たす。従って，あるMM m’について，ある入力文字列上に文字列m’hが存在する。ここで，tをQ(tes’)+|t|=Q(es’)=Q(h)を満足する最長の文字列とする(図11)。|t|≧|s|を仮定すると，Q(tes’)+|t|⊆Q(ses’)+|s|⊆Q(es’)=Q(h)=Q(tes’)+|t|であり，図１１に示すように，ある (i,j)∈Q(ses’)+|s|について(i,j)∈Q(m’h)+|m’|が成立する。このとき，e’=ses’中のes’がm’h中のhのprefixだから，e’はm’と左で重なり(E3)に矛盾するか，m’を真の部分文字列として含み(E4)に矛盾する。よって|t|＜|s|。ところが，tはQ(tes’)+|t|=Q(es’)=Q(h)を満たす最長の文字列なので，任意の(i,j)∈Q(tes’)+|t|について(i,j-|t|)∈Q(t)および(i,j)∈Q(es’)=Q(h)でありSi[j-|t| .. j+|h|-1]=thより(i,j)∈Q(th)+|t|だからQ(tes’)+|t|⊆Q(th)+|t|である一方，任意のa∈ΣについてQ(ates’)+|t|+1⊂Q(tes’)+|t|⊆Q(th) +|t|だから，tはQ(th)+|t|=Q(h)を満たす最長の文字列でもある。従って，thは(M4)を満たしMMであるが，|t|＜|s|より，e’はthと右で重なり(E3)に矛盾するか，thを真の部分文字列として含み(E4)に矛盾する(|es’|=|h|のとき)。結局，|s|≧1の場合は矛盾が生じるので，|s|=0であり，e’=es’。|es’|＞|h|およびQ(es’)=Q(h)より，es’は右EG-holder hのprefixであり，(E’1)を満足する。さらに，e’=es’は，(E1)-(E4)を満たすと仮定したから(E’2)-(E’4)も満たす。ところが，e≠e’よりs’≠εであり，eはe’の真のprefixである。これは，eが(E’5)を満たすことに矛盾。結局，(E1)-(E4)を満足するe’=ses’≠eの存在を仮定したのが誤りであり，eが(E5)を満たすことがわかる。
【００９２】
次に，文字列eが(E1)-(E5)を満足すれば，(E’1)-(E’5)を満足することを示す。まず，(E’1)について検討する。sをQ(es)=Q(e)なる最長の文字列とする。このesがEG-holderであることを示す。eは(E1)を満足するからQ(es)=Q(e)≠φ。よってesは(R1)を満たす。|es|≧|e|≧uより，esは(R2)を満たす。定義よりsはQ(e)=Q(es)なる最長の文字列だから，(es)は(R3)を満たす。よって，esは右MMであり，(H1)を満たす。
【００９３】
esの真のprefixである右MM rの存在を仮定する。rは(R3)を満たすから，esのprefixであるrに続く文字をaとすれば，Q(r)⊃Q(ra)⊇Q(es)=Q(e)より，|r|＜|e|。(R4)より，rをsuffixとするあるMM mが存在しQ(m)+|m|=Q(r)+|r|。ゆえに，eはmと左で重なり(E3)に矛盾するか，mを真の部分文字列として含み(E4)に矛盾する(図12)。よってこうしたrの存在を仮定したのが誤りで，esは右MMを真のprefixとして持たないことがわかる。すなわち，esは(H3)を満たす。
【００９４】
さらに，esが(H2)を満たすことを示す。esが(H2b)を満足する場合は明らかに(H2)は満足される。esが(H2b)を満足しない場合を考察する。esがMMでないと仮定する。esをprefixとする入力文字列があればesはMMであり仮定に反するから，esはどの入力文字列のprefixでもない。ところで，esは右MMであったから，esがMMでないということは(M4)が満たされないことを意味し，ある文字aについてQ(aes)+1=Q(es)≠φが成り立つ。このとき，aeは，以下に示すように，(E1)-(E4)を満足する。
(E1)について:
eは(E’1)を満たすからQ(aes)+1=Q(es)=Q(e)≠φ。ゆえに，Q(ae)⊇Q(aes)≠φ。
(E2)について:
|ae|=|e|+1＞u。
(E3)について:
aeが，あるMMと重なると仮定する。eは(E3)を満たすからMMと重ならないので，aeはaをsuffixとするMMと左で重なるか，eをprefixとするMMと右で重なるかのいずれかである。esは(H2b)を満足しないので，任意のMM mに対し，mesを部分文字列とする入力文字列は存在しない。Q(es)=Q(e)より，meを部分文字列とする入力文字列も存在しない。ゆえに，aeは左ではMMと重ならない。一方，eをprefixとするMM mの存在を仮定すると，Q(es)=Q(e)より，|m|＜|es|ならばesが(H3)を満たすことに矛盾し，|m|≧|es|ならば(R5)よりesはMMだが，esがMMでないことに矛盾する。よって，いずれの場合も矛盾が生じるので，aeはMMと重ならないことがわかる。
(E4)について:
aeの真の部分文字列でMM mの存在を仮定する。eが(E4)を満たすから，eの真の部分文字列にMMはない。さらに，eがMMとすれば，s=εならばesがMMでないことに矛盾し，s≠εならばesの真のprefixに右MMがないことに矛盾。よって，mはaeの真のprefixである。m=aと仮定するとQ(aes)+1=Q(es)より，esは常に左側にMM aを伴って入力文字列中に出現するから，(H2b)に矛盾。従って，m≠a。ところが，aeの真のprefixであるmが2文字以上の文字列ならば，eはmと左で重なり，eが(E3)を満たすことに矛盾。以上より，aeの真の部分文字列であるmの存在を仮定すれば必ず矛盾が生じるので，aeが(E4)を満たすことがわかる。
【００９５】
以上より，aeは(E1)-(E4)を満たすが，aeはeを真の部分文字列として含むので，eが(E5)を満たすことに矛盾。従って，esがMMでないとした仮定が誤りである。つまり，esは(H2b)または(H2a)を満たすから(H2)を満たし，右EG-holderである。よって，eは(E’1)を満たす。
【００９６】
さらに，esは|es|≧|e|≧uより(E’2)を満たし，(E3)と(E4)を満足するから(E’3)と(E’4)を満たす。eを真のprefixとし，(E’1)-(E’4)を満たす文字列e’の存在を仮定すると，e’は(E’1)-(E’4)より(E1)-(E4)をそれぞれ満たすが，これはeが(E5)を満たすことに矛盾。よってこのようなe’は存在せず，eが(E’5)を満たすことがわかる。
【００９７】
以上で，条件(E1)-(E5)が条件(E’1)-(E’5)と等価であることを示した。従って，本発明の方法で，条件(E1)-(E5)を満たす文字列すなわちEGを過不足無く抽出できることを示すには，(E’1)-(E’5)の条件を満たす文字列の集合と，「EGである」と印がつくノードvのパスラベルp(v)の長さEGlength(v)のprefix eの集合が，一致することを示せばよい。
【００９８】
まず，eが(E’1)-(E’5)を満たすことを示す。
(E’1)について:
EGlength(v)はp(v)から$を除いた文字列の長さに初期化されて以降，値が増加することはないため，p(v)の長さEGlength(v)のprefixは$を含まない。よって，eは(E’1)を満たす。
(E’2)について:
条件(e3b)より，「EGである」と印がついたノードvについて，EGlength(v)≧uだから，eは(E’2)を満たす。
(E’3)，(E’4)について:
ステップ(e2b)において，変数cの値703（図７参照）は，(i,j)をアピアランスとする部分文字列で，入力文字列iにおいてj番目の文字以外には，任意のMMの先頭の文字と終端の次の文字を含まない，最も長い文字列の長さとなっている。一方，任意の右EG-holder hの任意のアピアランス(i,j)について，v=Tv[i,j]とすると，hが右EG-holderだから，v=v’(h)であり，さらに，ステップ(e2b)でEGlength(v)はその時点のEGlength(v)と変数cの値703のうち，小さいほうの値に更新される。また，EGlength(v)の初期値は，|p’(v)|=|h|である。従って，(e2)が完了した時点で，EGlength(v)は，次の条件(e4a)および(e4b)を満足する最大の整数nとなる。
(e4a) n≦|h|
(e4b) hの，長さnのprefixで，MMと重なる文字列や，MMを部分文字列とする文字列は存在しない。
ゆえに，文字列eは(E’3)と(E’4)を満たす。
(E’5)について:
EGlength(v)は，eをprefixとする右EG-holderのprefixで，(E’3)と(E’4)を満たす最長のものの長さである。eの長さはEGlength(v)だから，(E’5)も満たされる。
【００９９】
逆に，本発明の方法で「EGである」と印がついたノードvのEGlength(v)のprefix以外の文字列は，EGでないことを示す。
【０１００】
まず，右EG-holderのprefix以外の文字列は排除されるが，そうした文字列は(E’1)を満たさない。条件(e3b)により，長さがu未満の文字列は排除されるが，これらは(E’2)を満たさない。さらに，EGlength(v)は右EG-holderであるp’(v)のprefixで(E’3)と(E’4)を満たす最も長い文字列の長さだから，vのパスラベルのprefixでEGlength(v)より長い文字列は(E’3)または(E’4)に矛盾する。さらに，EGlength(v)より短い文字列は(E’5)に矛盾する。
【０１０１】
以上により，本発明の方法でEGであると判断されない文字列は，EGではないことがわかった。
【０１０２】
■本発明の方法の処理時間
本発明の方法の処理時間の，入力文字列の文字の総数Nに対する振る舞いについて検討する。本発明の処理は，入力文字列数や入力文字列の文字数に依存しない回数のsuffix tree Tに対する深さ優先探査と，それぞれの深さ優先探査における各ノードに対する定数時間の処理と，全要素数がNである3つの2次元配列Head, Tail, Tvの各要素に対する入力文字列数や入力文字列の文字数に依存しない回数の定数時間の処理と，(t3)において各リーフを探査するO(N)時間で完了できる処理の，固定回数の組み合わせである。従って，本発明の方法全体の計算時間は，O(N)である。
【０１０３】
■本発明の方法を実現する装置
本発明は，前記方法を実行するための装置も提供する。図８に，装置の構成の一例を示す。該装置は，主記憶806に上記方法を実行するプログラム805を格納し，さらにsuffix tree Tや入力文字列を格納する。プログラム805は，中央演算装置801により実行される。計算結果は，ディスプレイ802を通じて表示されるか，補助記憶装置807に格納されるか，またはそれら両方の処理が行なわれる。ユーザからの入力はキーボード803およびポインティングデバイス804を用いて行なわれる。本発明の装置は，インターネットやイントラネット等のネットワークによって他の装置と通信可能に接続されていてもよい。入力文字列は，例えばファイルの形で与えられ，CD-R等の記録媒体に記録されたファイル，あるいはネットワークを介して受信したファイルを読み込むことによって主記憶806に取り込まれる。もちろん，キーボードから入力文字列を入力することも可能である。本明細書では，入力文字列を装置の主記憶806に取り込むための手段を総称して文字列入力手段という。
【０１０４】
本発明の方法で得られたEGをディスプレイ802に表示する場合，見易さや解析の容易さの観点から，図1の例のように，EGを入力文字列そのものまたは該入力文字列を象徴的に示す線分や矩形101上の，EGである部分文字列に相当する位置104に，色を変えたり文字や記号で示すなど視覚的に見やすい方法で表示することが好ましい。また，入力文字列をそのものまたは該入力文字列を象徴的に示す線分や矩形101とともに表示するか否かに関わらず，EGの入力文字列上の位置102を表示してもよい。また，同じ入力文字列に由来するEGを表示する際には，各EGの開始位置の昇順とすることが好ましい。昇順の表示は，各入力文字列iについて変数jを，0から，入力文字列iの長さから1を減じた値まで増加させ，Tv[i,j]が「EGである」と印がついたノードvへのポインタであったとき，p(v)のprefixで長さがEGlength(v)であるEGを表示することで，実現できる。EGを補助記憶装置807に格納する際にも，同じ入力文字列から得られたEGをともに格納する場合には，上記の表示方法で開始位置の昇順の表示を得たのと同様の方法で，EGの開始位置の昇順を得てその順に格納することが好ましい。
【０１０５】
また， EGを表示する場合には，色や記号，文字列103などを用いて複数の入力文字列に存在する同一のEGを識別できるように表示することが好ましい。また，補助記憶装置807に記録する際も，複数のEGが同一であると識別可能なように番号あるいはEGごとに一意である数値または文字列を同時に記録することが好ましい。
【０１０６】
■本発明の方法の利用者による，エクソン境界の指定
本発明の方法では，cDNA配列の相互比較によりエクソンの境界を探索しているが，この方法で発見できないエクソン境界が，何らかの手段により知られている場合があり得る。そこで，本発明の方法の利用者が既知のエクソン境界の位置を入力することにより，その境界でEGを分割できることが望ましい。そのためには，「既知のエクソン境界が，入力文字列Siのj番目の文字とj-1番目の文字の間にある」と入力された場合に，ステップS505においてHeadとTailを設定した直後に，以下の操作を行えばよい。
(b1) Head[i,j]を１に設定する。
(b2) j≧1ならば，Tail[i,j-1]を1に設定する。
このようにしてからステップS506以降を実行すると，入力文字列Siのj番目の文字とj-1番目の文字の間にMMの境界があるときと同様に処理が行われ，EGが分割される。なお，既知のエクソン境界が複数箇所存在する場合には，それぞれの境界に対し(b1)-(b2)の処理を繰り返し実行すればよい。
【０１０７】
■本発明の方法の，cDNA配列以外への適用
本発明は，mRNA配列のスプライシングパターン解析を目的としているが，比較ゲノム解析にも応用が可能である。入力文字列をmRNA配列の代わりにゲノム配列とし，得られたEGがゲノム内の保存領域と解釈すれば，本発明の方法により，ゲノム配列の保存領域の位置と，複数のゲノム配列中の保存領域の対応関係が明らかにでき，比較ゲノム解析の有力な手段となる。
【０１０８】
この他，文字列に使用される文字の種類が有限であれば，本発明の方法は塩基配列に依存せず一般の文字列を処理することができるため，アミノ酸配列のような塩基配列以外の配列や文字列に対しても，そのまま適用可能である。
【実施例１】
【０１０９】
本発明の方法を実装したシステムを作成し，本発明の方法で塩基配列のスプライシングパターンを解析できることを実証した。以下では，uの値を15とした場合の結果について述べる。
【０１１０】
はじめに，乱数を用いて，長さが30文字の4つの塩基配列s1,s2,s3,s4を用意した。そして，これらを組み合わせ，長さが90の4つの入力文字列s1s2s3，s1s2s4，s1s3s1, s2s3s4を作成し，本発明の方法を実装したソフトウェアを用いて，これらがs1,s2,s3,s4の組み合わせであることを認識できるか試みた。3つ目の文字列に，s1が２度現れることに注目されたい。パラメータuの値は，15とした。なお，ここで用いたs1,s2,s3,s4は以下の配列である。
s1= TTCAACAAAGACGGAAGTGTCCTAAATAGG
s2= GTGCTGACAGTGCTGTTAGAACTACAGGCT
s3= GAAGAAAGGTAACGCATATAGTGCGACGAA
s4= GGTGGCATGCCATGGACGCATACTCCGTAA
【０１１１】
処理の結果，図１に示す4つのEGが得られた。この結果において，s1とs4は，EG1, EG4として抽出されたが，EG2，EG3は，それぞれ先頭の1塩基が欠けた配列として抽出された。この原因は，s2,s3がいずれもGで始まるため，s1およびs2の末尾にGを加えた配列がMMとなり，s2,s3の先頭が1塩基削られた配列が右EG-holderとなったためである。このように，末尾での配列の偶然の一致により，EGが実際のエクソンより短くなる場合がありうる。しかし，完全な一致でなかったとはいえ，EG2，EG3はそれぞれs2,s3にほぼ一致した。つまり，本発明の方法により，入力文字列が長さ約30塩基の4配列の組み合わせであると認識できたといえる。さらに，本発明の方法では，３つ目の配列に，s1がEG1として２箇所現れていることも認識することができた。
【実施例２】
【０１１２】
本発明の方法を実装したシステム用いて，米国公共機関のデータベースRefSeqに登録されている，実際の遺伝子のスプライシング・バリアント配列のうち，ヒトのウィルムス腫瘍遺伝子WT1の配列を取得して解析し，得られたスプライシングパターンをデータベースの記述と比較した。使用した配列は，アクセッション番号がNM_004906, NM_152857, NM_152858の3つの配列である。データベースには，NM_152857の576番塩基以降とNM_152858の675番塩基以降の配列が選択的3'端配列であることと，NM_152857に選択的5'端配列が存在することが記載されている。この3配列に対し，本発明の方法を実装したシステムを適用した結果は，図２の通りである。データベースに記載されていたNM_152857, NM_152858の選択的3'配列は，開始位置を含めて正確にEG5として認識された。さらに，NM_152857の選択的5'端配列に関しては，データベースに位置の記述はなかったが，先頭から115番目の塩基までの配列が選択的5'端配列であり，NM_004906, NM_152858の先頭から213番目の塩基までの5'端配列に置き換わっていることがわかった。しかし，NM_004906の選択的3’配列が，1566塩基目までと1592塩基以降の２つに分断されEG4とEG6となり，それらに挟まれた25塩基が空白区間となった。この空白区間の25塩基は全てTであった。(M1)-(M4)より，長さ15塩基以上のTの連続は全てMMとなり，EGがこれらのMMを真の部分文字列としないように選ばれたのが，選択的3’配列がEG4とEG6に分断された原因であった。このことから，本発明の方法を適用する場合には，リピート配列に関する注意が必要だとわかる。
【０１１３】
なお，3配列ともに末尾にEGに含まれない空白区間が見られたが，これらはpolyA配列であり，エクソン配列ではない。この区間では17塩基以上Aが連続しているが，15個以上のAの連続はMMとなるためpoly-A配列の箇所にはMMの末端が密集し，EGは存在しなかった。
【０１１４】
なお，この実験で用いた３つの入力配列の総塩基数は5,529bpで，計算時間はCPUのクロック周波数が1.7GHzであるパーソナルコンピュータを用いて0.1秒程度であった。
【図面の簡単な説明】
【０１１５】
【図１】４つのエクソンの組み合わせになっている４つの入力配列を，本発明の方法により解析した結果を表示した例を示す図。
【図２】本発明の方法により，ウィルムス腫瘍遺伝子のスプライシング・バリアントからエクソン配列を抽出した結果を表示した例を示す図。
【図３】入力文字列，右MM，MM，右EG-holder，EGの説明図。
【図４】Suffix treeの説明図。
【図５】本発明の方法における，処理全体のフローチャート。
【図６】MMのアピアランスの数の総和がΩ(Nk)となる例を示す図。
【図７】図5中，ステップS505の右EG-holder抽出，および，ステップS506のEG抽出方法の説明図。
【図８】本発明の方法を実現する装置の一例の説明図。
【図９】条件(E’1)-(E’5)より条件(E5)が導かれることの説明において，|es’|≦|h|を導くために検討した文字列の関係を説明する図。
【図１０】条件(E’1)-(E’5)より条件(E5)が導かれることの説明において，あるMM m’について，ある入力文字列上に文字列m’hが存在することを導くために検討した文字列の関係を説明する図。
【図１１】条件(E’1)-(E’5)より条件(E5)が導かれることの説明において，あるMM m’について，|t|＜|s|を導くために検討した文字列の関係を説明する図。
【図１２】条件(E1)-(E5)より条件(E’1)を導くために検討した文字列の関係を説明する図。
【図１３】文字列s,s'が重なることの定義を説明する図。
【図１４】複数のエクソンが存在することを，転写産物由来の配列のみに基づき認識するのが困難な例を説明する図。
【図１５】複数のエクソンの境界を，転写産物由来の配列のみに基づき同定するのが困難な例を説明する図。
【符号の説明】
【０１１６】
101: 入力文字列を象徴的に示す矩形
102: EGの入力文字列上の位置を表す数字
103: 複数の入力文字列上にあるEGを識別するための文字列の例
104: 101上で，EGである部分文字列に相当する位置
301: 入力文字列の例
302: 301の入力文字列集合に対する右MMを列挙したもの
303: 301の入力文字列集合に対するMMを列挙したもの
304: 301の入力文字列集合に対する右EG-holderを列挙したもの
305: 301の入力文字列集合に対するEGを列挙したもの
401: 文字列ATATGとTTAGTAから構築されたsuffix tree
402: Suffix tree 401のルートノード
403: Suffix tree 401のリーフの１つ（「i, j」と書かれているリーフは，文字
列iのj番目から始まるsuffixに対応）
404: Suffix tree 401の，エッジの１つ
405: Suffix tree 401の，エッジのラベルの１つ
407: Suffix tree 401の，ノードの１つ
408: Suffix tree 401の，ルートノード402からあるノードにいたるパス
409: Suffix tree 401のsuffix linkの１つ
701: 2次元配列Headの，入力文字列2についての要素を取り出したもの
702: 2次元配列Tailの，入力文字列2についての要素を取り出したもの
703: EG抽出方法において用いられる変数cのステップ(e2b)での値
801: 本発明の装置の中央演算装置
802: 本発明の装置のディスプレイ
803: 本発明の装置のキーボード
804: 本発明の装置のポインティングデバイス
805: 本発明の方法を実行するためのプログラム
806: 本発明の装置の主記憶装置
807: 本発明の装置の補助記憶装置
【図１】
2005316605.tif 000007
【図２】
2005316605.tif 000008
【図３】
2005316605.tif 000009
【図４】
2005316605.tif 000010
【図５】
2005316605.tif 000011
【図６】
2005316605.tif 000012
【図７】
2005316605.tif 000013
【図８】
2005316605.tif 000014
【図９】
2005316605.tif 000015
【図１０】
2005316605.tif 000016
【図１１】
2005316605.tif 000017
【図１２】
2005316605.tif 000018
【図１３】
2005316605.tif 000019
【図１４】
2005316605.tif 000020
【図１５】
2005316605.tif 000021
【配列表】
2005316605000001.app