JP2005346340A

JP2005346340A - フラグメントによる配列クラスタリング・アラインメント方法

Info

Publication number: JP2005346340A
Application number: JP2004164387A
Authority: JP
Inventors: Tomohiro Yasuda; 知弘安田
Original assignee: Hitachi Ltd
Current assignee: Hitachi Ltd
Priority date: 2004-06-02
Filing date: 2004-06-02
Publication date: 2005-12-15

Abstract

【課題】複数の配列を比較した結果である、複数の配列に共通する部分配列や、特定の塩基配列にのみ存在する部分配列の情報（フラグメント情報）が与えられたときに、配列を互いに共通部分が多い配列のグループごとに分類するクラスタリング処理と、グループ内で、配列の各部位の対応関係を明らかにするアラインメント処理を行なう。
【解決手段】与えられたフラグメント情報に基づき、それらのフラグメントをノードとしてもつ有向グラフを構築する。この有向グラフから閉路を削除した後、最小重み経路を探索する公知の手法を用いて最小重み経路を同定する。この経路に対応するフラグメント列をもとに、互いに共通部分の多い配列を同定し、さらに、前記最小重み経路を利用して、各配列をフラグメント単位でマルチプルアラインメントする。
【選択図】図１

Description

複数の塩基配列に共通に存在する部分文字列を抽出した結果に基づき、配列のクラスタリングおよびアラインメントを行なうための方法に関する。

国際共同プロジェクト及び米国ベンチャー企業により、２０００年６月にヒトゲノムの文字配列決定の完了が宣言されたほか、大腸菌や酵母菌、線虫、ショウジョウバエ、マウスなど、多くの生物のゲノム配列が続々と決定された。ゲノム配列において、機能上重要な部位は、種間で保存性が高いことが知られている。これは、機能上重要な部位では、中立的な突然変異が起こりにくいためである。したがって、近縁種で相同性のある部位の配列を比較することにより、ゲノム配列中の機能上重要な部位を同定することが可能である。

ゲノム配列と並び、ｃＤＮＡ配列の解析も重要と考えられている。ｃＤＮＡとは、生体内のｍＲＮＡを逆転写して得られるＤＮＡであり、その配列を解析することは、生体内で発現している遺伝子の配列を解析することに他ならない。ｃＤＮＡ配列の一部を配列決定したＥＳＴと呼ばれる配列は、米国の公共機関のデータベースに２千万件以上蓄積されている折、ヒトの全長ｃＤＮＡ配列を、全遺伝子の２／３以上を網羅すると言われる２万配列以上をアノテーション付きで公開するデータベースも構築された（今西等の文献“Imanishi, T., et al., Integrative Annotation of 21,037 Human Genes Validated by Full-Length cDNA Clones, PLoS Biol. Vol 2, Issue 6, 2004.”）。これらの膨大なｃＤＮＡ配列データを用いて遺伝子の機能を調べたり、同一遺伝子に由来するｃＤＮＡ配列を同定したりする際にも、ｃＤＮＡ配列の比較解析が有効である。

多数の配列の同時相互比較を、最適化問題と見て最適解を探索する処理は、計算量的に極めて困難であり、短時間で効率よく最適解を得る方法は存在しないと言われている。そのため、最適解を得ることは保証されないものの、高速処理が可能で実用的に十分な比較結果を得るための手法が開発されてきた。それらの手法のうち、配列の類似性が極めて高い領域をまず同定し、その結果に基づいて配列全体の比較を行なう手法が、ゲノム配列やｃＤＮＡ配列の比較においては有効である。このような配列の比較解析を行なう従来技術には、以下に述べるものがある。

非特許文献１のホモロジー検索の方法は、ひとつの問い合わせ配列と複数の配列を格納したデータベースが与えられたとき、問い合わせ配列に類似する領域をもつほかの配列を探索する技術であり、配列の比較解析のツールとして利用できる。ただし、この方法は、２配列間の類似する領域、または、ある１配列に類似する多数の配列上の領域を同定することができるのみである。多数の配列に共通する類似配列の探索や、多数の配列間で対応する部分の同定、多数のｃＤＮＡ配列中で同一遺伝子に由来する配列の同定といった問題を、直接的に処理することはできない。

２配列の比較を、効率よく行なう手法として、非特許文献２のＤｅｌｃｈｅｒらの手法が知られている。この手法では、まず与えられた２配列に共通するＭＵＭ（Maximal Unique Match）と呼ばれる部分配列を抽出する。そして、得られたＭＵＭを、数列中の最長の単調増加部分列を抽出する方法を用いて整列させ、２配列の共通部分が配列全体でどう分布しているかを明らかにする。共通部分間のギャップは、同様の手法を再帰的に適用することで埋める。Delcherらは、非特許文献３においてこの手法の改良を行なっているが、比較できる配列の数は２にとどまっており、３つ以上の配列を同時に比較することはできない。

Ｈｏｈｌらは、複数のゲノム配列を高速にアラインメントすることを目的とし、３つ以上の配列に共通する部分配列を抽出し、抽出された部分配列を整列し、その結果に基づき複数配列のアラインメントを高速に生成する方法を開発した（非特許文献４）。しかし、Ｈｏｈｌらの方法は、全配列に共通の配列のみを対象としており、一部の配列にのみ共通する配列を抽出することができない。そのため、複数のゲノム配列の全てではなく一部のみの配列に共通する部分配列や、ｃＤＮＡ配列中の選択的スプライシングが行なわれるエクソン配列を見落とすという問題がある。

一方で、３つ以上の配列から、共通する部分配列を抽出する試みも行なわれている（非特許文献５）。しかし、非特許文献５は、ＵＥＢ（共通する部分配列）を抽出する手法を提案するだけであって、それだけでは、配列を共通部分が多い配列ごとのグループに分類するクラスタリングや、配列間で対応する部分を明らかにするアラインメントを行なうための方法が提供されるわけではない。

Altschul, S.F., Madden, T.L., Schaffer, A.A., Zhang, J., Zhang, Z., Miller, W. and Lipman, D.J., Gapped BLAST and PSI-BLAST: a new generation of protein database search programs, Nucleic Acid Research, 25:3389-3402, 1997

Delcher, A.L., Kasif, S., Fleischmann, R.D., Peterson, J., White, O. and Salzberg, S.L., Alignment of whole genomes, Nucleic Acids Research, 27(11):2369-2376, 1999. Delcher, A.L., Phillippy, A., Carlton, J., and Salzberg, S.L., Fast algorithms for large-scale genome alignment and comparison, Nucleic Acids Research, 30(11):2478-2483, 2002. Hohl, M., Kurtz, S., and Ohlebusch, E., Efficient multiple genome alignment, Bioinformatics 18 Suppl.1, S312-S320, 2002. Yasuda, T., Kimura, K., Nishikawa, T., Towards Splicing Pattern Detection Based on cDNA Sequences, Genome Informatics 14:422-423, 2003.

本発明は、塩基配列のフラグメント情報が与えられた場合に、そのフラグメント情報に基づき配列のクラスタリングおよびアラインメントを行なう方法を提供することを目的としている。

図１は、本発明のフラグメント列解析方法を模式的に示す図であり、処理対象となる三つの塩基配列の塩基、フラグメント、フラグメント列、および、処理結果の例を示す。

本発明では、入力として、複数のフラグメントの列が与えられることを想定する。本発明の処理対象となる塩基配列は、ゲノム配列やｃＤＮＡ配列といった塩基配列であり、分析の対象となる塩基配列が、当初からフラグメント情報を持ったものとして与えられるわけでは無い。しかし、背景技術で説明したように、ゲノム配列やｃＤＮＡ配列といった塩基配列をフラグメント情報に変換する手法は種々の提案があるので、これらのいずれかにより処理されたフラグメント情報に基づき配列のクラスタリングおよびアラインメントを行なう方法を提供するものである。

本発明では、複数の塩基配列に共通するか、または、特定の塩基配列のみに存在する部分配列をフラグメントと呼び、ある塩基配列に現れるフラグメントを、該塩基配列に現れる順に並べたものをフラグメント列と呼ぶ。

塩基配列１はフラグメントｆ_１，ｆ_２およびｆ_３からなるフラグメントの列から成り、塩基配列２はフラグメントｆ_１およびｆ_３からなるフラグメントの列から成り、塩基配列３はフラグメントｆ_４，ｆ_５およびもう１つのｆ４からなるフラグメントの列から成るものとする。このような複数の塩基配列１−３のフラグメント情報が与えられたとき、それぞれの塩基配列からフラグメントを抽出して、各塩基配列におけるフラグメントの並び、すなわち、フラグメント列を導出する。塩基配列１のフラグメント列Ｆ１は＜ｆ_１，ｆ_２，ｆ_３＞であり、塩基配列２のフラグメント列Ｆ２は＜ｆ_１，ｆ_３＞であり、塩基配列３のフラグメント列Ｆ３は＜ｆ_４，ｆ_５，ｆ_４＞である。

これらのフラグメント列Ｆ１−Ｆ３を相互比較し、互いに共通部分が多いグループ１０３_１、１０３_２に分類するクラスタリングと、クラスタリングされたフラグメント列間のアラインメントを行なって複数のフラグメント列間で共通する部分の対応関係１０４_１、１０４_２を明らかにする。

上記の処理を効率良く行なう方法を提供することが、本発明の課題である。

上記の課題を解決するために、次の（１）−（４）の工程からなる方法を実行する。
（１）複数のフラグメント列のフラグメント情報が与えられたとき、常に隣り合って現れるフラグメントを統合する。
（２）フラグメント列の集合をＤＡＧ（directed acyclic graph）と呼ばれるデータ構造に変換し、ＤＡＧ上の２点間の最短路を同定する方法を用いて、互いに共通部分の多いフラグメント列を同定する。
（３）ＤＡＧ上の最短路を探索する際に、前もって与えられたスコア関数により各フラグメントのスコアを計算する。あるフラグメントｆのスコアを計算するスコア関数としては、例えば、フラグメント長をＬ（ｆ）、該フラグメントの全フラグメント列中における総出現回数をＮ（ｆ）とするとき、ｓｃｏｒｅ（ｆ）＝Ｌ（ｆ）Ｎ（ｆ）、ｓｃｏｒｅ（ｆ）＝Ｌ（ｆ）、または、ｓｃｏｒｅ（ｆ）＝Ｌ（ｆ）ｌｏｇ（Ｎ（ｆ））を用いる。
（４）ＤＡＧ上の最短路に対応するフラグメント列が、全フラグメント列中の一部のフラグメント列からなる、互いに共通部分の多いグループを特徴付けているとみなし、該フラグメント列を用いて、該グループに属するフラグメント列を同定する。
（５）前記グループに属するフラグメント列間のマルチプルアラインメントを行なう。

本発明によれば、塩基配列に基づき作られたフラグメント列に基づき、同一遺伝子に由来するｃＤＮＡ配列や近縁種のゲノム配列のようなグループごとに、フラグメント列を分類し、マルチプルアラインメントにより配列中の対応する箇所を同定することが可能になる。

以下、本発明の実施の形態について説明する。最初に、本明細書で使用する記号と用語及び概念を定義する。
（１）フラグメント：解析の対象である塩基配列の一部分。ある塩基配列に固有の場合や、複数の塩基配列に共通する場合がある。複数の配列に共通する場合には、主に完全一致の場合を想定しているが、少数の塩基の置換や短いギャップを許すなどバリエーションも考えられる。非特許文献５で説明されているＵＥＢがフラグメントの例であるが、他の種類のフラグメントに対しても、本発明の方法は有効である。
（２）フラグメント列：フラグメントの順序つきの列。フラグメントｆ_１，ｆ_２，．．．，ｆ_ｎがこの順に並んだ列を＜ｆ_１，ｆ_２，．．．，ｆ_ｎ＞と表記する。図１に例を挙げた。なお、本明細書ではｉ番目のフラグメント列をＦｉで表し、フラグメント列の数をＫとすれば、１≦ｉ≦Ｋとする。
（３）Ｏ（ｆ（ｎ））：関数ｇ（ｎ）に対しｇ（ｎ）＝Ｏ（ｆ（ｎ））であるとは、ある定数ｃが存在して、十分大きなｎに対しｇ（ｎ）≦ｃｆ（ｎ）が成立することである。また、ある量が「Ｏ（ｆ（ｎ））である」とは、その量がｎの関数ｇ（ｎ）以下であり、ｇ（ｎ）＝Ｏ（ｆ（ｎ））であることを意味する。

以下、図面を参照しつつ、本発明の方法について詳述する。

図２は、本発明による処理の全体を示すフローチャートである。

ステップＳ２０１：配列のバリエーションにより分断されたフラグメントの統合
本発明では、複数のフラグメント列が与えられたとき、まず、各フラグメント列について、それぞれのフラグメントが単一のフラグメントと見なすのが好ましいが、ＳＮＰやマイクロサテライトにより２つ以上に分断されているフラグメントに対しては、統合して一つのフラグメントに置き換える処理を行なう。そのため、あるフラグメントｆまたはｆ'を含むすべてのフラグメント列において、フラグメントｆの直後のフラグメントがｆ’である場合、フラグメントｆとフラグメントｆ’を統合し新しいフラグメントｆ’’で置き換える。

図３は分断されているフラグメントを統合する例の説明図である。上段に示すように、塩基配列１ではフラグメントｆの直後のフラグメントｆ’が塩基Ａで分断され、塩基配列２ではフラグメントｆの直後のフラグメントｆ’が塩基無しで分断され、塩基配列３はフラグメントｆの直後のフラグメントｆ’が塩基Ｇで分断されている例である。ここで、３０１は分断されたフラグメント間に存在するギャップであり、ＳＮＰと考えられるものの例である。統合後のフラグメントｆ’’のフラグメント長Ｌ（ｆ’’）は、Ｌ（ｆ’’）＝Ｌ（ｆ）＋Ｌ（ｆ’）で定義する。統合後のフラグメントｆ’’を下段に示す。

フラグメントｆと直後のフラグメントｆ’の、元の塩基配列中での距離が与えられれば、それらの中間のギャップに相当する領域３０１の性質を解析することができる。図３の例のように、フラグメントｆと直後のフラグメントｆ’の距離が常に１以下でありギャップ３０１に存在する塩基が一意でない場合には、このギャップ３０１はＳＮＰまたはシーケンシングエラーである。また、フラグメントｆと直後のフラグメントｆ’の距離が、常にある定数Ａ，Ｂと整数ｎを用いてｎ＊Ａ＋Ｂと表現できれば、フラグメントｆと直後のフラグメントｆ’の間にタンデムリピートの存在が示唆される。ギャップ３０１の配列をｓとすれば、ｓを塩基数でＡだけずらした配列とｓ自身を比較することにより、実際にマイクロサテライト等のタンデムリピートが存在するか否か判定できる。

ステップＳ２０２：全フラグメント列に基づく有向グラフＧの構築
全てのフラグメント列に基づき、有向グラフＧを構築する。図４はフラグメント列に基づき構築される有向グラフＧを、図１の例から構築して示す図である。有向グラフＧは、ノード４０１の集合Ｖと、ノード間を結ぶ有向辺（エッジ４０２）の集合Ｅから成る。このとき、Ｅ⊆Ｖ×Ｖである。

まず、ノードの集合Ｖは次のように構築する。

Ｖ＝｛ｆ｜ｆをフラグメント列の要素にもつフラグメント列が存在｝∪｛ｓ，ｔ｝
ただし、ｓおよびｔは、それぞれ、配列の先頭と終端を象徴的に表すシンボルである。

一方、エッジ（ノード間を結ぶ有向辺）の集合Ｅは次のように構築する。

Ｅ＝｛（ｆ，ｆ’）｜あるフラグメント列で、ｆはｆ’に先行（ｆの直後にｆ’が存在しなくても良い）｝
∪｛（ｓ，ｆ）｜ｆはフラグメント列中の任意のフラグメント｝
∪｛（ｔ，ｆ）｜ｆはフラグメント列中の任意のフラグメント｝
−｛（ｆ，ｆ’）｜あるフラグメント列で、フラグメントｆ’はフラグメントｆに先行（ｆ’の直後にｆが存在しなくても良い）｝
エッジ構築後、出入りするエッジのないノードは、削除してもよい。集合Ｅと集合Ｖを構築する順番を入れ替え、エッジの無いノードは初めから作らないことにしてもよい。

各エッジ（ｆ，ｆ’）に、重みｗ４０３を設け、ｗ（ｆ，ｆ’）＝−ｓｃｏｒｅ（ｆ’）で定義する。ただし、ｗ（ｆ，ｔ）＝０と定める。さらに、任意の（ｆ，ｆ’）∈Ｅに対し変数Ｄ（ｆ，ｆ’）４０４を準備し、次の式を満足するように初期化する。

Ｄ（ｆ，ｆ’）＝（全フラグメント中で、ｆがｆ’に先行する回数）
ステップＳ２０３：有向グラフＧにおける閉路の除去
有向グラフＧ上を、深さ優先探査し、コルメン等の文献”Cormen, T.H., Leiserson, C.E., Rivest, R.L., Clifford, S.,Introduction to algorithms (second edition), chapter 22 and 24, MIT Press, Cambridge, MA, 2001.”で述べられている後退辺(文献中ではback edge)を探索する。なお、後退辺とは、有向グラフＧ上の深さ優先探索中に、あるノードｕから出てその祖先ノードｖへ入るエッジのことである。後退辺であるエッジが見つかったら、そのエッジを削除する。後退辺であるエッジを削除しても、ｖからｕに至るパスが存在するので、有向グラフＧの連結成分の数は不変である。深さ優先探査が終了した時点で、グラフＧは、閉路の無い有向グラフになっている。このようなグラフはＤＡＧ（directed acyclic graph）と呼ばれている（前出：コルメン等の文献）。

ステップＳ２０４：有向グラフＧにおける重み最小経路の探索
有向グラフＧを、前記コルメン等の文献に記載の方法ＴＯＰＯＬＯＧＩＣＡＬ−ＳＯＲＴ（Ｇ）を用いてトポロジカルソートする。そして、同じく、コルメン等の文献に記載の方法ＤＡＧ−ＳＨＯＲＴＥＳＴ−ＰＡＴＨＳ（Ｇ，ｗ，ｓ）を用いて、配列の先頭ｓから配列の終端ｔに至る最小重みの経路Ｐｏを探索する。ただし、ここで言う経路とは、有向グラフＧのノードを順序をつけて並べたものであって、隣合うノード間には先行するノードから後のノードへのエッジが存在するものである。ある経路Ｐの重みを、Ｐ中でノードｖ１の次のノードをｖ２とするとき、式（１）で定義する。ただし、ｓｃｏｒｅ（ｔ）＝０と定める。Ｐｏとは、式（１）を最大化する経路である。

ステップＳ２０５：Ｐｏを用いた同一グループに属する配列の判定
式（１）を最大化する経路Ｐｏから、配列の先頭ｓ、配列の終端ｔを除いて得られるフラグメント列をＦｏとする。ラグメント列Ｆｏに現れる任意のフラグメントｆについて、フラグメントｆがあるフラグメント列Ｆｉにｎ（ｆ，Ｆｉ）回現れるならば、フラグメント列Ｆｉ毎に設けられる変数ｏｖｅｒｌａｐ（ｉ）に、Ｌ（ｆ）ｎ（ｆ，Ｆｉ）を加える。実装を簡単にするために、単にＬ（ｆ）を加える方法もある。なお、ｏｖｅｒｌａｐ（ｉ）の初期値は０とする。

この処理が終わった段階で、ｏｖｅｒｌａｐ（ｉ）の値は、式（２）となる。

一方、フラグメント列Ｆｉが抽出された元の塩基配列Ｓｉの長さを｜Ｓｉ｜とする。ここで、ｏｖｅｒｌａｐ（ｓ）／｜Ｓｉ｜が与えられたパラメータＴを上回る場合には、フラグメント列Ｆｉはフラグメント列Ｆｏと十分に良く重なると判断する。フラグメント列Ｆｏと十分によく重なると判断されたフラグメント列Ｆｉを、互いに共通部分が多いフラグメント列のグループに属するものとする。このグループを、以下グループＣと呼ぶ。

ステップＳ２０６：グループＣのフラグメント列のアラインメント
グループＣに属すると判定されたフラグメント列を、マルチプルアラインメントする。そのためには、トンプソン等の文献"Thompson, J.D., Plewniak, F., and Poch, O., A comprehensive comparison of multiple sequence alignment programs, Nucleic Acids Research, 27: 2682-2690, 1999."に挙げられている各種のシステムで活用されている方法を使用すればよい。

また、アラインメントの質を若干落としつつも、高速に動作し、かつ実用上十分な質の結果が得られる方法として、本発明では次の方法を提供する。

フラグメント列Ｆｏは、ステップＳ２０５にて得られたグループＣに属するフラグメント列と共通部分が多いと判定されているフラグメント列である。そこでフラグメント列Ｆｏに、グループＣに属するフラグメント列を、ダルビン等の文献”Durbin, R., Eddy, S., Krogh, A., Mitchison, G., Biological sequence analysis: Probabilistic models of proteins and nucleic acids, chapter 2, Cambridge University Press, Cambridge, UK, 1998.”に記載の大域アラインメントアルゴリズムを用いてアラインメントする。その結果得られるコンセンサス配列を新たにフラグメント列Ｆｏと見なして、グループＣの未処理で残っているフラグメント列に、同様の処理を行なう。これを、グループＣから未処理のフラグメント列が無くなるまで、繰り返し実行する。

ステップＳ２０７：有向グラフＧを、未処理フラグメント列に対し再構築
有向グラフＧから、グループＣのフラグメント列を削除する。その方法を２つ述べる。１つ目は、グループＣに属するフラグメント列を除いたフラグメント列に対して、有向グラフＧを作り直すことである。この場合、ステップＳ２０３も改めて実行し、有向グラフＧ中の閉路を除去する必要がある。２つ目は、以下の（１）−（２）の処理を行なう方法である。
（１）フラグメントｆ，ｆ’を、式（１）を、最大化する経路Ｐｏ中で隣接する２つのノードとする。
（２）エッジ（ｆ，ｆ’）が存在する場合には、変数Ｄ（ｆ，ｆ’）から１を減じる。変数Ｄ（ｆ，ｆ’）が０になったら、エッジ（ｆ，ｆ’）を削除する。その際、ノードｆに出入りするエッジが無くなれば、ノードｆをＧから削除してもよい。同様に、ノードｆ’に出入りするエッジが無くなれば、ノードｆを有向グラフＧから削除してもよい。

この方法で、エッジが削除されることはあっても新たに導入されることはないので、有向グラフＧに新たな閉路が生じることはなく、有向グラフＧはＤＡＧのままに保たれる。

繰り返し
新たに再構築された有向グラフＧを用いて、ステップＳ２０４−Ｓ２０７を繰り返す。残っているフラグメント列が存在する間、この繰り返しを行なう。

本発明の方法の処理時間
Ｌｆをフラグメント列の長さ（フラグメント列のフラグメント数）の最大値、
Ｌｃをグループのコンセンサス配列の長さの最大値、
Ｍを最大のグループサイズ、
Ｎをクラスタ数、
Ｆをフラグメントの数とする。

このとき、ステップＳ２０１に要する処理時間は、Ｏ（Ｆ）、
ステップＳ２０２に要する処理時間は、Ｏ（｜Ｖ｜＋｜Ｅ｜）＝Ｏ（｜Ｖ｜Ｌｆ^２）＝Ｏ（ＦＬｆ^２）、
ステップＳ２０３に要する処理時間は、Ｏ（｜Ｅ｜）＝Ｏ（｜Ｖ｜Ｌｆ^２）＝Ｏ（ＦＬｆ^２）、
ステップＳ２０４を一回実行するために要する時間は、Ｏ（｜Ｖ｜＋｜Ｅ｜）＝Ｏ（ＦＬｆ^２）、
ステップＳ２０５を一回実行するために要する処理時間は、Ｏ（ＦＬｃ＋Ｆ）＝Ｏ（ＦＬｃ）、
ステップＳ２０６を一回実行するために要する処理時間は、Ｏ（ＭＬｃＬｆ）、
ステップＳ２０７を一回実行するために要する処理時間は、Ｏ（ＭＬｆ）である。

ステップＳ２０４−Ｓ２０７は、Ｎ回繰り返す。

ただし、本発明の方法全体で、処理するフラグメント列数はＦであるため、
ステップＳ２０６を実行するために要する処理時間は、Ｏ（ＦＬｃＬｆ）、
ステップＳ２０７を実行するために要する処理時間は、Ｏ（ＦＬｆ）である。

よって、本発明の方法全体の処理時間は、漸近的時間計算量で表すと、式（３）の通りとなる。

本発明の方法は、フラグメント数が極端に多くなると処理に時間がかかるが、フラグメント数は塩基数よりも少ないため、実用的には十分短い時間で処理が可能である。

本発明の方法を実現する装置
本発明は，前記方法を実行するための装置も提供する。図５に、本発明の方法を実現する装置の一例の説明図を示す。該装置は，主記憶装置５０６に上記方法を実行するプログラム５０５を格納し，さらに有向グラフＧやフラグメント列を格納する。プログラム５０５は，中央演算装置５０１により実行される。計算結果は、ディスプレイ５０２を通じて表示されるか、補助記憶装置５０７に格納されるか、またはそれら両方の処理が行なわれる。

ユーザが直接フラグメント列やパラメータを与える場合には、入力はキーボード５０３およびポインティングデバイス５０４を用いて行なわれる。本発明の装置は、インターネットやイントラネット等のネットワークによって他の装置と通信可能に接続されていても良い。フラグメント列は、例えばファイルの形で与えられ、ＣＤ−Ｒ等の記録媒体に記録されたファイル、あるいはネットワークを介して受信したファイルを読み込むことによって主記憶装置５０６に取り込まれる。本明細書では、フラグメント列を装置の主記憶装置５０６に取り込むための手段を総称してフラグメント列入力手段という。

本発明の方法で得られた、フラグメント列のクラスタリングおよびマルチプルアラインメントの結果を表示する場合、見易さや解析の容易さの観点から、図１の例の下部のように、クラスタリングの結果得られたグループを象徴的に表すボックス１０３により区別し、マルチプルアラインメントにより対応するとわかったフラグメント１０４を、図１のように揃えて表示することが好ましい。

（実施例１）
本発明の方法を実装したシステムを作成し、本発明の方法で塩基配列のスプライシングパターンを解析できることを実証した。以下では、パラメータＴの値を０．２５、ｓｃｏｒｅ（ｆ）＝Ｌ（ｆ）Ｎ（ｆ）とした場合の結果について述べる。また、有向グラフＧの再構築には、１つ目の方法を用い、マルチプルアラインメントには本明細書で述べた高速手法を用いた。

米国公共機関のデータベースＲｅｆＳｅｑに、２００４年５月６日の時点で登録されているヒト遺伝子のｃＤＮＡ配列のうち、”transcript variant”との表記があるもの全て（６３５７配列）を収集し、その中から無作為に選択した１００〜２０００配列について、非特許文献５の方法で得られたＵＥＢをフラグメントとして、フラグメント列を構成し、本発明の方法によりクラスタリングおよびアラインメントを行った。プルット等の文献”Pruitt, K.D., and Maglott, D.R., RefSeq and LocusLink: NCBI gene-centered resources, Nucleic Acids Research 29:137-140, 2001.”におけるＬｏｃｕｓＬｉｎｋにおいて、同一遺伝子座の配列を同一クラスタと見なしたデータを正解データとし、本発明の方法によるクラスタリング結果の評価を行なった。正解データに含まれるクラスタの集合をＣ０、本発明の方法により生成されたクラスタの集合をＣ１とするとき、式（４）で表されるジャッカードインデックス−Jaccard index−Ｊ（Ｃ１，Ｃ２）を、評価指標として用いた。

図６はＬｏｃｕｓＬｉｎｋを正解データとした場合の本発明のクラスタリング結果の精度を示すグラフである。

ゲノム配列を用いずｃＤＮＡ配列だけでクラスタリングを行なう場合、ＳＮＰやシーケンシングエラーによる配列のバリエーションや、ファミリー遺伝子や共通ドメインにより別遺伝子でも配列類似性が存在するといった問題がある。完全なクラスタリングは困難であり、ジャッカードインデックスで０．９を超える結果は、十分な精度であると言える。

図７は本発明の実施に要する処理時間を示すグラフである。２０００配列の処理に要した処理時間は７４．８秒で、本発明が実用的な処理速度を持つことがわかる。ただし、フラグメント列の数が増加するに従い、処理時間が急激に増大するため、多数の配列を扱う場合には注意が必要である。なお、この実施例においては、クロック周波数が３ＧＨｚのＣＰＵを持ち、２ＧＢの主記憶を搭載した計算機を使用した。

図８は、図1で考え方を説明した分析結果のディスプレイ５０２上の表示例を示す図である。フラグメントｆ１，ｆ２，---，ｆｎを横軸に表示し、フラグメント列Ｆ１，Ｆ２，---を縦軸に表示する。各フラグメント列について存在が確認できたフラグメントの位置を見やすいＸＹマトリクス表示で示す。さらに、クラスタリングの結果得られたフラグメント列のグループごとに分離して表示する。この結果、クラスタリングとマルチプルアラインメントの結果が、明確に表示できることになる。さらに、ＸＹマトリクス表示の下段に、各フラグメントの塩基配列を表示すれば、解析者にとって参考となる。

本発明を、ｃＤＮＡ配列やゲノム配列の解析に用いることで、ｉｎｓｉｌｉｃｏアノテーションが進み、創薬、診断等の技術開発に欠かせない有用遺伝子の発見が加速されるものと期待される。

本発明の方法が処理対象としている、塩基配列、フラグメント列、および、処理結果の例を示す図。本発明の方法における、処理全体のフローチャート。分断されたフラグメントを統合する例の説明図。フラグメント列に基づき構築される有向グラフＧの例を示す図。本発明の方法を実現する装置の一例の説明図。ＬｏｃｕｓＬｉｎｋを正解データとした場合の本発明のクラスタリング結果の精度を示すグラフ。本発明の実施に要する処理時間を示すグラフ。図1で考え方を説明した分析結果のディスプレイ５０２上の表示例を示す図。

符号の説明

１０３…互いに共通部分の多いフラグメント列のグループを表す矩形、１０４…互いに共通部分の多いフラグメント列のグループ中で複数のフラグメント列間で対応するフラグメントを表す矩形、３０１…分断されたフラグメント間に存在するギャップ、４０１…フラグメントから生成された、有向グラフのノード、４０２…フラグメント列から構築される有向グラフのエッジ、４０３…フラグメント列から構築される有向グラフの重みｗ、４０４…フラグメント列から構築される有向グラフの変数Ｄの値、５０１…中央演算装置、５０２…ディスプレイ、５０３…キーボード、５０４…ポインティングデバイス、５０５…本発明の方法を実行するためのプログラム、５０６…主記憶装置、５０７…補助記憶装置。

Claims

与えられたフラグメント情報に基づき、それらのフラグメントをノードとしてもつ有向グラフを構築すること、
該有向グラフから閉路を削除した後、最小重み経路を同定すること、
該経路に対応するフラグメント列をもとに、互いに共通部分の多い配列を同定すること、よりなることを特徴とするフラグメント列解析方法。
前記フラグメント情報において、あるフラグメントｆ１またはｆ２を含むすべてのフラグメント列において、フラグメントｆ１の直後にフラグメントｆ２が続く場合に、フラグメントｆ１とフラグメントｆ２を統合し、新たなフラグメントとする請求項１記載のフラグメント列解析方法。
前記フラグメント情報において、あるフラグメントｆ１またはｆ２を含むすべてのフラグメント列において、フラグメントｆ１の直後にフラグメントｆ２が続く場合に、フラグメントｆ１とフラグメントｆ２の元の塩基配列中での距離が１以下であり、フラグメントｆ１とフラグメントｆ２の間にある塩基が一意でない場合、フラグメントｆ１とフラグメントｆ２の間にＳＮＰが存在すると判定する請求項１記載のフラグメント列解析方法。
前記フラグメント情報において、あるフラグメントｆ１またはｆ２を含むすべてのフラグメント列において、フラグメントｆ１の直後にフラグメントｆ２が続く場合に、フラグメントｆ１とフラグメントｆ２の元の塩基配列中の距離がある定数Ａ、Ｂと整数ｎによりｎ＊Ａ＋Ｂと常に表現できる場合には、フラグメントｆ１とフラグメントｆ２の間にある部分配列ｓを塩基数でＡずらした配列と部分配列ｓを比較することによって、タンデムリピートの存在を判定する請求項１記載のフラグメント列解析方法。
前記フラグメント情報において、与えられたフラグメント列中のフラグメントをノードとし、あるフラグメント列で、あるフラグメントｆ１があるフラグメントｆ２に先行するとき、ノードｆ１とノードｆ２の間にエッジを設ける請求項１記載のフラグメント列解析方法。
前記フラグメント情報において、互いに共通部分を多く持つとされたフラグメント列のマルチプルアラインメントを行なう請求項１記載のフラグメント列解析方法。
前記フラグメント情報において、前記最小重み経路を、フラグメント列のマルチプルアラインメントにおける初期値とする請求項１記載のフラグメント列解析方法。
与えられたフラグメント情報に基づき、それらのフラグメントをノードとしてもつ有向グラフを構築すること、
該有向グラフから閉路を削除した後、最小重み経路を同定すること、
該経路に対応するフラグメント列をもとに、互いに共通部分の多い配列を同定すること、よりなることを特徴とするフラグメント列解析プログラム。
与えられたフラグメント情報に基づき、それらのフラグメントをノードとしてもつ有向グラフを構築すること、
該有向グラフから閉路を削除した後、最小重み経路を同定すること、
該経路に対応するフラグメント列をもとに、互いに共通部分の多い配列を同定すること、よりなる処理を経て解析されたフラグメント列解析結果を
フラグメントｆ１，ｆ２，---，ｆｎと、フラグメント列Ｆ１，Ｆ２，---とをＸＹマトリクスで表示するとともに、各フラグメント列について存在が確認できたフラグメントについてＸＹマトリクス交点の位置にフラグメントの存在を示す表示を行なうことを特徴とするフラグメント列解析結果の表示方法。