The Human Genome

ヒトゲノム配列TOPページへ


0.序論  
1.DNA供給源と塩基配列解析法  
2.ゲノムアセンブリ戦略と特徴  
3.遺伝子予測と注釈
4.ゲノムの構造  
5.ゲノムの進化
6.ゲノム全域の配列変異検査
7.ヒトゲノムにおいて予測される蛋白質コード遺伝子の概観  
8.結論


2.ゲノムアセンブリ戦略と特徴

要約

ここでは、ゲノムをアセンブリするために用いた2つの手法について述べる。第一の方法では、全配列の読み取り結果とGenBankの断片化されたデータとを、コンピュータにより組み合わせ、一つの独立した偏りのないゲノムの概観を作成する。第二の手法では、マッピング情報に基づいて、全ての断片を、ある領域または染色体にまとめ、次に、まとめたデータを断片化して、コンピュータによるアセンブリを行なう。これらの方法では、正しい順序と正しい方向性でアセンブリされたDNA塩基配列が、基本的には同じように再構築される。第二の手法では、わずかに大きい配列包括度(sequence coverage)(大きいほどギャップが少なくなる)が得られ、ゲノム解析段階で用いる主要な配列となった。また次に、このアセンブリ過程の完全さと正確さについて論じ、主に独立したBACごとのアプローチ(BAC-by-BAC approach)により再構築された国際ゲノムプロジェクトのゲノム配列との比較を行う。我々が行ったアセンブリは、ヒト染色体の真正染色質領域を効率よくカバーしている。ゲノムの90%以上が、10万bp以上のアセンブリ骨格中にあり、25%が1000万bp以上のアセンブリ骨格の中にあった。

ショットガン法によるアセンブリは遡行問題(inverse problem)をはらむ典型的な例である。すなわち、解析標的とした配列から無作為にサンプリングした1セットの読み取り配列が与えられた時、標的配列の中の順序と位置を再構築せよというわけである。キイロショウジョウバエのために開発したゲノムアセンブリアルゴリズムは、今や約25倍も大きいヒトゲノム解析向けに拡張されている。セレラ社のアセンブリは、1セットのコンティグ(連結断片)からなり、これが配列骨格に順に並べられ方向付けされた後、既知のマーカーによって染色体の位置にマッピングされるのである。コンティグは、読み取り配列断片の集まりであり、配列重複部分があるために、ゲノムの隣接区間のコンセンサス配列を再構築する土台になる。メイト対は、アセンブリ戦略の中心的な構成要素であり、連続したコンティグ間のギャップサイズがかなり精密にわかっているアセンブリ骨格を作成するのに用いる。これは、一方が1つのコンティグの中にあり、もう一方が別のコンティグの中にある一対の読み取り配列が、二つのコンティグ間の方向性と距離を暗示していることを観察することにより達成できる(図3)。最終的には、我々のアセンブリでは、報告した最終セットのアセンブリ骨格に、読み取った配列を全て組み込んだわけではない。組み込まなかった読み取り配列は「チャフ(chaff、おじゃまむし)」と名付けた。チャフの典型例は、高反復領域内から出てきた読み取り配列、多数のゲノムプロジェクトで認められていた他生物由来の配列で様々な経路で持ち込まれたもの、および品質の低い配列やトリミングされていないベクターの配列である。




図3  全ゲノムアセンブリの解剖図。オーバーラップした断片化bactigフラグメント(赤線)と、5人からセレラ社内で得られた読み取り配列(黒線)を組み合わせ、コンティグとコンセンサス配列を作る(緑線)。メイト対の情報を使って、コンティグをつなげてアセンブリ骨格(赤)を作る。アセンブリ骨格は、STS(青星)物理的マップ情報により、ゲノム(灰色線)にマッピングされる。

2.1. アセンブリデータセット

アセンブリには、二つの独立したデータセットを用いた。第一のデータセットは、セレラ社で作成した平均長543bpの、2727万個の読み取り配列からなるランダムショットガン・データセットである。これは主に、5人のドナー由来のDNAサンプルから作成した16個のライブラリーのメイト対読み取り配列で構成される。挿入配列サイズが2、10、50kbpのものを用いた。1つのライブラリーに由来するメイト対が、ゲノム内の既に配列決定された区域にどのように位 置しているかを知ることで、各ライブラリーにおける挿入配列サイズの範囲の特徴がわかり、その平均値と標準偏差を決定することができた。表1に、読み取り配列の数、配列カバー倍数、データセットにより達成したクローン配列カバー倍数(clone coverage)の詳細を示した。クローン配列カバー倍数は、両端からの配列を持つ各クローンの挿入片全体を考慮したものだが、クローニングされたDNA中のゲノムの配列カバー倍数といえる。ひいては、ゲノムの物理的DNAカバー倍数(DNA coverage)を計る尺度となる。ゲノムサイズを2.9Gbpと仮定すると、セレラ社のトリミングした配列はゲノムの5.1倍であり、クローン配列カバー倍数は2、10、50kbpライブラリーのそれぞれ3.42、16.40、18.84倍、総計38.7倍であった。

第二のデータセットは、公的資金によって展開された国際ヒトゲノムプロジェクト(PFP)から得たもので、主にBACクローンを塩基配列解析に用いて決定したものである(30)。アセンブリに入力したこのBACデータは、2000年9月1日にGenBankからダウンロードして得たもので(表2)、全長4443.3 Mbpの配列であった。各BACのデータは、その完成度によって4つのフェーズに分類される。フェーズ0のデータは、BACクローンをごく軽いショットガン法にかけて得たものであり、一般的にはアセンブリされていない読み取り配列セットで、全ゲノムの1倍未満を特徴としている。フェーズ1のデータは、BACコンティグまたは"bactig"と呼んでいるコンティグが順番に並んでいない(unordered)集合体であり、フェーズ2のデータは順番に並んだ(ordered)bactigの集合体である。フェーズ3のデータは、完全なBAC配列である。過去2年間、PFPは、各BACクローンの3−4倍の配列カバー倍数の軽ショットガン法データから得られるフェーズ1データの作成に専念してきた。時間的にはより早くできるが品質と完成度は低いデータの作成になってしまっている。

表2 アセンブリに入力したGenBankのデータ


Center Statistics Completion phase sequence
0 1 and 2 3

Whitehead Institute/ Number of accession records 2,825 6,533 363
MIT Center for Number of contigs 243,786 138,023 363
Genome Research, Total base pairs 194,490,158 1,083,848,245 48,829,358
USA Total vector masked (bp) 1,553,597 875,618 2,202
Total contaminant masked (bp) 13,654,482 4,417,055 98,028
Average contig length (bp) 798 7,853 134,516
Washington University, Number of accession records 19 3,232 1,300
USA Number of contigs 2,127 61,812 1,300
Total base pairs 1,195,732 561,171,788 164,214,395
Total vector masked (bp) 21,604 270,942 8,287
Total contaminant masked (bp) 22,469 1,476,141 469,487
Average contig length (bp) 562 9,079 126,319
Baylor College of Number of accession records 0 1,626 363
Medicine, USA Number of contigs 0 44,861 363
Total base pairs 0 265,547,066 49,017,104
Total vector masked (bp) 0 218,769 4,960
Total contaminant masked (bp) 0 1,784,700 485,137
Average contig length (bp) 0 5,919 135,033
Production Sequencing Number of accession records 135 2,043 754
Facility, DOE Joint Number of contigs 7,052 34,938 754
Genome Institute, Total base pairs 8,680,214 294,249,631 60,975,328
USA Total vector masked (bp) 22,644 162,651 7,274
Total contaminant masked (bp) 665,818 4,642,372 118,387
Average contig length (bp) 1,231 8,422 80,867
The Institute of Physical Number of accession records 0 1,149 300
and Chemical Number of contigs 0 25,772 300
Research (RIKEN), Total base pairs 0 182,812,275 20,093,926
Japan Total vector masked (bp) 0 203,792 2,371
Total contaminant masked (bp) 0 308,426 27,781
Average contig length (bp) 0 7,093 66,978
Sanger Centre, UK Number of accession records 0 4,538 2,599
Number of contigs 0 74,324 2,599
Total base pairs 0 689,059,692 246,118,000
Total vector masked (bp) 0 427,326 25,054
Total contaminant masked (bp) 0 2,066,305 374,561
Average contig length (bp) 0 9,271 94,697
Others* Number of accession records 42 1,894 3,458
Number of contigs 5,978 29,898 3,458
Total base pairs 5,564,879 283,358,877 246,474,157
Total vector masked (bp) 57,448 279,477 32,136
Total contaminant masked (bp) 575,366 1,616,665 1,791,849
Average contig length (bp) 931 9,478 71,277
All centers combineddagger Number of accession records 3,021 21,015 9,137
Number of contigs 258,943 409,628 9,137
Total base pairs 209,930,983 3,360,047,574 835,722,268
Total vector masked (bp) 1,655,293 2,438,575 82,284
Total contaminant masked (bp) 14,918,135 16,311,664 3,365,230
Average contig length (bp) 811 8,203 91,466

* Other centers contributing at least 0.1% of the sequence include: Chinese National Human Genome Center; Genomanalyse Gesellschaft fuer Biotechnologische Forschung mbH; Genome Therapeutics Corporation; GENOSCOPE; Chinese Academy of Sciences; Institute of Molecular Biotechnology; Keio University School of Medicine; Lawrence Livermore National Laboratory; Cold Spring Harbor Laboratory; Los Alamos National Laboratory; Max-Planck Institut fuer Molekulare, Genetik; Japan Science and Technology Corporation; Stanford University; The Institute for Genomic Research; The Institute of Physical and Chemical Research, Gene Bank; The University of Oklahoma; University of Texas Southwestern Medical Center, University of Washington.
dagger The 4,405,700,825 bases contributed by all centers were shredded into faux reads resulting in 2.96× coverage of the genome.


我々は、以下の3つのデータセットに対して、BLASTアルゴリズムを用いて、bactig配列の夾雑物をスクリーニングした。(i)Univec coreのベクター配列(38)について:98%の同一性を持つ配列の配列末端で25bpが一致するか、配列内部で30bpが一致するか、(ii)GenBankのハイスループットゲノム配列部門(High Throughput Genomic(HTG)Sequence division)のヒト以外の生物由来の配列(39)について: 98%同一性を持つ配列で200 bpが一致するか、(iii)霊長類ウイルスやヒトウイルの侵入遺伝子を含まないGenBankの非反復塩基配列について:98%同一性を持つ配列で200 bpが一致するか、でふるい分けた。コンティグの末端50 bp以内に25 bp以上のベクター配列を認め次第、ベクターの部分まで末端を削除した。これらの基準の下で、我々は、夾雑物やベクターと考えられる配列を、フェーズ3のデータから2.6 Mbp分、フェーズ1と2のデータから61.0 Mbp分、フェーズ0のデータから16.1 Mbp分削除した(表2)。これにより、計4363.7 MbpのPFP配列データ(20%が完全状態、75%が粗稿状態(フェーズ1と2)、5%が一回だけの読み取り配列(フェーズ0))が残った。さらに、104,018対のBAC 末端配列メイト対をダウンロードし、二つのアセンブリ処理データセットに加えた(18)

2.2. アセンブリ戦略

2種類のアセンブリ法を用いた。第一の手法は、セレラ社のデータとPFPデータをあわせて合成ショットガンデータとして用いる全ゲノムアセンブリプロセス(whole-genome assembly process)である。第二の手法は、最初にセレラ社とPFPのデータを、大きな染色体分節に局在している組に分けて、各組ごとに最初からショットガンアセンブリを行う区画化アセンブリプロセス(compartmentalized assembly process)である。図4にこれらの工程の全体的な流れを示す。




図4  セレラ社の2分岐アセンブリ戦略。長円形はそれぞれのラベルで示される機能を実行するコンピュータプロセス(訳注:本文中では例えばマッチング機などと記載)を示す。長円形間の矢印上のラベルは、プロセスによって作り出される、または使われるオブジェクトの性質を示す。この図は本文中の議論をまとめたものである。使われている用語と成句は、本文中に定義されている。

全ゲノムアセンブリでは、PFPデータを、最初に、完全にbactigの2倍分カバーする550 bpの読み取り配列の合成ショットガンデータセットに分解または断片化した。これにより、BACデータセットのなかに配列重複があるために、結果的にゲノムを2.96倍カバーするのに十分な、1605万個の「人造(faux)」読み取り配列ができた。しかしPFPアセンブリプロセスに固有の偏りは組み込んでいない。次に、4332万個の読み取り配列を合わせたデータセット(ゲノムの8倍)と全ての関連したメイト対情報に対して、我々の全ゲノムアセンブリアルゴリズムを適用し、ゲノムの再構築を行った。ゲノム内のBACの位置も、bactigのアセンブリも、このプロセスでは用いなかった。bactigうち2.13%のアセンブリが間違っているという有力な証拠を見いだしたため(40)、bactigは読み取り配列断片長にまで断片化した。その上、PFPの出した物理地図上に正しく配置されていないBACがいくつかあったこと、ならびに、おそらくサンプル追跡ミスの結果(後述)、少なくとも2.2%のBACが、そのBACのものではない配列データを含んでいるという有力な証拠を見いだしたこと(41)から、BAC位置情報は無視した。要するに、我々は、外部で作成されたデータからはちょっとばかり良い配列カバー倍数を得ただけであって、メイト対情報やアセンブリされたbactig情報、ゲノムの位置情報は使わず、本当に最初から全ゲノムアセンブリを行ったのである。

区画化ショットガンアセンブリ(compartmentalized shotgun assembly, CSA)では、セレラ社とPFPのデータを、自信を持って決定できる最も大きい染色体分節と思われるもの、すなわち「コンポーネント(区画)」に分割し、この区画化したサブセットに対してショットガンアセンブリを適用した。この場合、そのコンポーネント個別の、最初からのアセンブリであるということを確実にするため、bactigデータを再び人造読み取り配列断片長に断片化した。データをこの方法で分けることによって、全体的なコンピュータ作業が削減され、染色体間の重複の影響が改善された。また、結果としては、全ゲノムアセンブリ法とは相互に独立したゲノムの再構築を行うことになり、一致しているかどうか2つのアセンブリを比較することができた。異なるゲノム領域が混らないようにするには、コンポーネント分割の質(quality)が重要である。われわれは、このコンポーネントを、(i)各BACから得た最長の配列骨格からと、(ii)セレラ社データセット特有のアセンブリ骨格から作成した。BACアセンブリは、bactigとこれらのbactigにマッピングした5倍数のセレラ社データをインプットデータとして用い、結合アセンブラ(combining assembler)を使って得たものである。この作業は、中間段階として行った。所定の配列に対してアセンブリ骨格が正確で完全であればあるほど、配列重複部位とメイト対情報に基づき、これらの骨格を一層正確な連続したコンポーネントとして並べることができるからという理由による。我々はさらに正確さを増すために、コンポーネント中のアセンブリ骨格の並び方を視覚的に探究し整理した。最終的なCSAアセンブリでは、区画化以外の件は全て無視し、区画化した該当セレラ社データと、区画化した該当bactigデータを断片化した人造読み取り配列断片に対して、我々の全ゲノムアセンブリアルゴリズムを適用した結果、各コンポーネントで個別に当初段階からの配列再構成を行うことができた。

2.3. 全ゲノムアセンブリ

ヒトゲノムの全ゲノムアセンブリ(WGA)に用いたアルゴリズムは、キイロショウジョウバエゲノムの配列決定[詳細な報告は(28)]に用いた方法を拡大した方法である。

WGAアセンブラは、5つの主要な段階、すなわち、ふるい機(Screener)、オーバーラップ検出装機(Overlapper)、ユニティグ検出機(Unitigger)、アセンブリ骨格形成機(Scaffolder)、反復配列解析機(Repeat Resolver)で構成されるパイプラインからなる。ふるい機(Screener)は、6-bp未満のエレメントで構成されるマイクロサテライト反復配列を全て発見してマークし、Alu、Line、リボゾームDNAなど散在する既知の反復配列を全てふるい落とす。残りのマークされた領域はオーバーラップ(塩基配列一致)の有無が検索される。ふるい落とされた領域はオーバーラップ検索をされることはないが、オーバーラップの一部でふるいに引っかからないがつき合わせれば一致するセグメントを含む可能性もある。

オーバーラップ検出機(Overlapper)は、1つひとつの読み取り配列を他の全読み取り配列と比較し、最短でも40bpの、端から端まで完全なオーバーラップ(塩基配列一致)がある配列と、つき合わせ配列どうし中の差異が6%以下しかないものを検出する。全てのデータは、綿密にベクターをトリミングしているので、オーバーラップ検出機は、完全にオーバーラップする配列があるぞあるぞと主張できる。全てのオーバーラップセットを計算するのに、4ギガバイトのRAMを搭載した4プロセッサーAlpha SMP1台では、約10,000CPU時間かかる。このような機械を40台並行して作動させたところ、実所要時間は4〜5日となった。

上述のように計算処理した1つひとつのオーバーラップは、統計学的には1017分の1の確率の事象であり、従って偶発的なものではない。アセンブリの組み合わせを困難にしているものは、実際にゲノムがオーバーラップしている領域からサンプリングされたオーバーラップ配列は多いが(すなわち、その配列はアセンブリすべきものである)、実際にはもっと多くのオーバーラップが、前述のスクリーニングで除外されなかった低反復配列数エレメントを含む2つの異なる配列由来のものである(組合わせると誤りになる)という点である。我々は前者を「真のオーバーラップ」、後者を「反復配列によって引き起こされる間違ったオーバーラップ(repeat-induced overlaps)」と呼ぶ。アセンブラは、特に、工程初期に、この間違ったオーバーラップを選別しないよう注意しなければならない。

我々は、この目的をユニティグ検出機で達成した。先ず、他の読み取り配列アセンブリに比べて、これは一種類しかなく文句なく採択できると思われる読み取り配列アセンブリを見つける。これらのサブアセンブリから構成されるコンティグ群をunitigと呼ぶ(ユニークにアセンブリしたコンティグの意(uniquely assembled contigs))。正式には、これらのunitigは、全てのオーバーラップを示すグラフの中で、より確実に正しいと思われる中間段階のサブグラフである(42)。経験的にはこれらのアセンブリの多くは正しいが(すなわち、真のオーバーラップしか含まない)、残念なことにいくつかは、実際は複数の反復塩基配列由来の配列の寄せ集めで、単一サブアセンブリ(体)になるまで分断しすぎたものとなってしまう。しかし、過剰分断したunitigは、平均カバー倍数の程度(coverage depth)が高くなりすぎて、全体の配列カバー倍数レベルとは一致しないため、容易に見分けることができる。我々は、unitigが、特有のDNA塩基配列だけか2つ以上のコピーを含む反復配列で構成されているかを判定し、オッズ比の対数を算出する単純な統計学的識別プログラムを開発した。この識別プログラムは、閾値を十分に厳しく設定すれば、正しいと確信できるunitigのサブセットを同定できる。また、あまり厳しくない閾値では、正しくアセンブリされた可能性の高い、残りのunitigのサブセットを同定できる。これらのうちで、一貫してアセンブリの足場となり得るものを選別していくと(後述参照)、選別されたunitigはほぼ正しいと考えられる。これらの2セットをあわせて、U-unitigsと呼ぶ。ヒト第22番染色体を6倍カバー倍数シュミレーションショットガン法でアセンブリした結果から、経験的に、我々が得たU-unitigsは長さ2 kbp以上の特有のDNA塩基配列の98%をカバーすることがわかった。我々はさらに、一つのU-unitigの末端で反復配列が始まる境界を同定できた。これを梃子にして、U-unitigsが、バラバラに散在するAlu配列と他の100〜400 bpの反復配列セグメントの93%以上を占めるようにした。

ユニティグ検出機(Unitigger)を作動させた結果、推測でヒトゲノムの73.6%をカバーし、正しくアセンブリされたサブコンティグのセットが得られた。ついでアセンブリ骨格形成機(Scaffolder)が、メイト対情報を用いて、これらをアセンブリ骨格にリンクさせる作業を行う。2つ以上のメイト対がある場合(所定のU-unitigsのペアがお互いに特定の距離と方向にあることを意味する)、メイト対の間違いが2%未満であると仮定すると、距離と方向性が誤りとなる確率は、さらに約1010分の1となる。したがって、少なくとも2個の2 kbpまたは10 kbpのメイト対で連結できるU-unitigsは、極めて確実に全部リンクさせることができ、中程度の大きさのアセンブリ骨格が作成できる。このアセンブリ骨格を、さらに50 kpbのメイト対とBAC末端配列を確認しながら、繰返してお互いに連結させる。このプロセスにより、サイズ的に100万bpの桁のアセンブリ骨格が作られた。これにはコンティグ間のギャップを含んでいるが、それは通常は反復配列に相当する部分で、まれに配列決定できない部分に相当する小さなギャップである。これらのアセンブリ骨格でゲノムの中の特有の塩基配列の大部分を再構成できる。

キイロショウジョウバエの場合のアセンブリでは、3段階の反復配列解析戦略(repeat resolution strategy)をとった。各段階は次第に難しくなり、そのためミスを犯しやすくなった。ヒトの場合のアセンブリでも、第1段階で「大岩」サブステージを継続使用することにした。この段階では、良好であるが決定的だとはいえない識別プログラムスコアを持つ全てのunitigを、アセンブリ骨格にあるギャップに配置する。これは、2対以上のメイト対があって、そのメイト対の読み取り配列のうち1つは既にアセンブリ骨格上にあるために、unitigを所定のギャップに曖昧さが全くなく配置させ得るという条件で行った。この方法でunitigを正しくないギャップに挿入してしまうのは、確率解析では10-7未満であると予測している。

我々は、ヒトアセンブリの次の段階「石」サブステージを改訂し、我々のこれまでの研究から示唆された機構に近づけるようにした(43)。各ギャップについて、アセンブリ骨格のコンティグの中にあり読み取り配列(R)の位置を示すはずのメイト対(M)を使って、ギャップの中に配置され得る全ての読み取り配列(R)を集める。セレラ社のメイト対情報は99%以上正しいので、こうすればこのセット中の、全てではないがほとんどの読み取り配列は、ギャップに収まるはずである。ある読み取り配列がそのギャップに収まらない時は、セットの中の残りの読み取り配列がそのギャップにうまく収まることはまれである。従って、我々は、単純に、この読み取り配列セットをギャップ内でアセンブルし、このアセンブリのやり方に矛盾するような読み取り配列は排除した。この作業は、キイロショウジョウバエのアセンブリのために行った作業よりも信頼性が高いことが判明した。ヒト第22番染色体について、シミュレートしたショットガン法データセットをアセンブリしたところ、全ての「石」が正しく位置された。

ギャップ解消の最終段階は、アセンブリしたBACデータ(ギャップをカバーしているはず)でギャップを埋めることである。これを、外部ギャップ「歩行(walking)」と呼ぶ。今回は、キイロショウジョウバエの配列決定で報告した非常に攻撃的な「丸小石」サブステージは適用しなかった。「丸小石」サブステージは、99.62%しか正しくない長い散在配列の反復再構築物(repeat reconstructions)を生み出すミスを犯したからである。ヒトゲノムのためには、99.99%以上の精度にならないことが判っている方法は用いない方が論理的に良い、と我々は決めた。その代価として、サイズが若干大きく、数も若干多いギャップができてしまったが。

アセンブリプロセスの最終段階およびそれに至るいくつかの中間過程では、各コンティグのコンセンサス配列を作製する。我々のアルゴリズムは、各塩基の評価に品質価値加重評価法(quality-value-weighted measures)を用いつつ、最大節減(maximum parsimony)原則で動いている。この効果は、各過程で、報告されるべき正しい塩基のベイズ推定(Bayesian estimate)に現れる。コンセンサス配列の作製には、セレラ社のデータがあればそれを利用する。セレラ社のデータが所定の領域をカバーしていないときは、BACデータ配列を用いる。

ヒトゲノムのWGAを達成できた鍵となったのは、オーバーラップ検出機と中央コンセンサス配列構築サブルーチン(central consensus sequence-constructing subroutines)を平行して動かすことであった。また、メモリーは重要な問題であった。キイロショウジョウバエのために構築したソフトウエアをそのまま適用すると、600ギガバイトのRAMを持つコンピュータが必要であった。オーバーラップ検出機にユニティグ検出機を加増することで、最大限の即時データ連絡により、同じ計算を28ギガバイトのRAMで済ますことができた。さらに、最初の3段階が加増的に処理できるため、データが配達されるたびに、計算のこの部分の状態を連続更新できた。また必要なときはいつでも、アセンブリ骨格構築と反復配列解析を完了するために7日間連続稼働させることができた。我々のアセンブリ作業のために使われた計算機インフラストラクチャーは、全体で、クラスターあたりで4ギガバイトのメモリーを備えた10台の4プロセッサー SMP(Compaq社、ES40、Regatta)と、64ギガバイトのメモリーを備えた1台の16-プロセッサーNUMAマシン(Compaq社、GS160、Wildfire)から構成される。アセンブラを1回作動させるための総計算時間は、おおよそ20,000CPU時間であった。

セレラ社データと、断片化bactigのデータをアセンブリすることにより、スパンが計2.848 Gbpで、2.586 Gbpの配列から構成される一組のアセンブリ骨格が作られた。チャフ、すなわちアセンブリに組み込まれなかった読み取り配列のセットは、1127万個(26%)であった。これはキイロショウジョウバエでの経験と矛盾しない。ゲノムの84%以上が100 kbp長以上のアセンブリ骨格でカバーされ、総計2.297 Gbpの配列のうち、平均91%を塩基配列が占め、9%がギャップであった。100 kbp以上のアセンブリ骨格1637個中に、計93,857個のギャップがあった。アセンブリ骨格サイズは平均1.5 Mbpであり、平均コンティグサイズは24.06 kbp、平均ギャップサイズは2.43 kbpで、それぞれのサイズ分布は基本的に指数関数的であった。全ギャップの50%以上が500 bpよりも短く、全ギャップの62%以上が1 kbpよりも短く、100 kbpを超える長さのギャップはなかった。同様に、配列の65%以上が30 kbp以上のコンティグの中にあり、31%以上が100 kbp以上のコンティグにあり、最大のコンティグは1.22 Mbpであった。表3に、このアセンブリ構造の統計を詳しくまとめ、区画化ショットガンアセンブリと直接比較した。

2.4 区画化ショットガンアセンブリ

WGAアプローチに加えて、部分的にアセンブリするアプローチを行った。これは、ゲノムを区画に分け、それぞれの区画を、個々にショットガン法でアセンブリする方法である。我々は、これが、染色体間でも多数ある重複の解決と、U-unitigsを計算するための統計処理の改良に有用であると考えた。区画化アセンブリプロセスでは、セレラ社の読み取り配列とbactigをゲノムの多数の大メガベース領域にクラスター化する過程と、その後にセレラ社データとbactigデータから得られた人造読み取り配列断片に対してWGAアセンブラーを走らせる過程を含んでいる。

CSAストラテジーの最初の段階は、セレラ社の読み取り配列を、特定のPFP BAC登録に対してBACコンティグに一致するものと、公表データと一致しないものに分けることであった。セレラ社の読み取り配列を適切に配置するために、このような相互一致があることは保証されなければならない。そのために、最初に全読み取り配列に対して、一般的反復配列エレメント(既にライブラリーがある)部分を目隠ししておいて、読み取り配列の目隠しされていない部分で40 bp以上が一致したときのみ、ヒットありとした。こうしてセレラ社の2727万個の読み取り配列中、2076万個がbactigと一致した。その他の62万個の読み取り配列には一致するものがなかったが、にもかかわらずそれらのメイト対がbactigと一致したため、bactigのBAC自体の領域に所属するものと判断された。残りの読み取り配列のうち、292万個の読み取り配列は完全にふるい落とされ、一致させることができなかった。しかし、他の297万個の読み取り配列は、GenBankデータセットには認められない計1.189 Gbpの目隠しされない配列を持っていた。セレラ社のデータは全体でゲノムの5.11倍に及ぶ重複があるため、我々は、240 Mbpのセレラ社特有の配列が、GenBankデータセットには含まれていないと推測している。

CSAプロセスの次の段階では、結合アセンブラが、ゲノムに関わる5倍分のセレラ社読み取り配列と、BACに入っているbactigを取り込み、組み合わせデータによりその場所特定のためのアセンブリを作成した。この高品質な再構成配列は一時的なものであり、次の段階でより信頼性の高い情報を提供するため、オーバーラップしたり隣接したりしている骨格配列のセットの中にそれらを埋めていくのに役立てようという単純なものである。概略を述べると、結合アセンブラは、最初に、一致するセレラ社読み取り配列のセットがあるかどうか調べ、ふるい落としで見逃された反復配列であることを示す過剰データ累積部位があるかどうか判断する。ある場合は、反復配列由来の読み取り配列(そのメイト対が本来あるべき位置にまだマッピングされていないもの)を除去する。次に、2つのbactigの間で相対的にみて同じ位置を一貫して示している全てのメイト対のセットを、1つのリンクとして束ね、その束の中のメイト対の数に応じて重みを与える。次に「がつがつ(greedy)」戦略で、この重みの順にメイト対の束を選び、bactigに順位をつける。適宜選抜されたメイト対の束は、造形用の2つのアセンブリ骨格を一緒に結びつけることができる。その骨格内のコンティグ間で大部分のリンクと一致している場合のみ、単一の骨格に組み上げるのである。一旦アセンブリ骨格が完全になれば、WGAアセンブラのところで前述した「石」戦略によりギャップを埋める。

フェーズ1と2のBACに関するGenBankのデータは、BAC1個あたり平均19.8個、平均8099 bpのbactigから構成されていた。結合アセンブラの適用により、個々のセレラ式BACアセンブリで平均1.83個の骨格にまとめられた(メディアン値:1個)。これは、平均サイズが18,973 bpの、平均8.57個のコンティグから構成されている。配列断片の順序と方向性の決定に加えて、この組み合わせ方式では、ギャップが57%少ない結果となった。フェーズ0のデータに関しては、平均のGenBank入力データは、平均784 bpの91.52個の読み取り配列から構成されていた。結合アセンブラの適用により、平均サイズが873 bpの、平均58.1個のコンティグから構成される平均54.8個のアセンブリ骨格が作られた。基本的には、一部少量のアセンブリが行われたが、典型的なフェーズ0のBACに代表される0.5〜1倍のデータセットに合致する十分な数のセレラ社データが正確にアセンブリされたでわけはなかった。結合アセンブラは、フェーズ3のデータにも適用され、SNP の同定、アセンブリの確認、セレラ社読み取り配列の位置決定が行われた。フェーズ0のデータは、BAC組み合わせ方式の全ゲノムショットガンデータセットとBACの1倍分の軽ショットガン塩基配列決定からでは、BACに拾われた領域の良いアセンブリ結果を生み出さないことを示唆している。少なくとも各BACについて3倍の軽ショットガン塩基配列決定が必要である。

GenBenkデータと適合しなかった589万個のセレラ社断片は、われわれの全ゲノムアセンブラでアセンブリした。アセンブリによって、計442 Mbpのスパンで、326 Mbpの配列から構成される、1セットのアセンブリ骨格が得られた。この骨格全体の20%以上が5 kbp以上の骨格配列で占められ、計302 Mbpの配列の中で、平均63%が配列そのもので27%がギャップであった。全ての5 kbp以上の骨格は、結合アセンブラで作成された全てのアセンブリ骨格とともに、次の、タイル積み段階(tiling phase)に持ち込まれた。

この段階での典型例では、ゲノムに関連した配列の少なくとも95%を構成するBAC領域すべてについて、1〜2個の骨格と、つながり不明のセレラ特有の骨格があった。ゲノム構成を明らかにするための次の段階は、ゲノム全体を通じてこれらのBAC骨格とセレラ特有骨格の順序づけとオーバーラップしている配列を、タイルを積み上げるようにして並びを決定するタイリングであった。このために、我々はセレラ社の50 kbp メイト対情報、BAC末端対情報(18)、配列標識部位(STS)マーカー(44)を用いて、広範囲にわたるガイド設定と染色体長分の分離を行った。比較的扱いやすい数の骨格であれは、このタイリングを完全に自動化した方法では行わず、最初のタイリングを良好な帰納法で計算し、その後人手(管理責任者)で、矛盾点や未結合点を解決する。このために、我々は、タイリングするオーバーラップのグラフとそれぞれの証拠を表示するグラフィカル・ユーザーインターフェースを開発した。管理責任者は、マッピングされたSTSデータ、配列オーバーラップのドットプロット、その選択を裏付けるメイト対の証拠のビジュアル表示が暗示するものを探査できた。このプロセスの結果、「コンポーネント」群が得られた。各コンポーネントは、管理責任者が認証したBAC骨格とセレラ特有骨格が積み上げられ並べられたセットである。このプロセスでは、推定スパンが2.922 Gbpの、3845個のコンポーネントが得られた。

最終的なCSAを作成するために、各コンポーネントを、WGAアルゴリズムを用いてアセンブリした。WGAプロセスで行ったのと同様に、アセンブラに、これとは独立にデータをアセンブリする自由度を持たせるために、bactigデータを2倍の人為的なショットガンデータに断片化した。bactigではなくこの人造読み取り配列を用いることで、アセンブリアルゴリズムは、もとのbactigのアセンブリの誤りを正し、PFP登録データの中のキメラ配列を除去できた。キメラ、すなわち、(ゲノムの他の部分から)混入した配列は、もともとそこにはなかったものだということで、コンポーネントの再アセンブリには組み込まなかった。これが効果を現し、CSAプロセスにおけるこれ以前の操作段階は、結局、ゲノムの大きな隣接区画に関係するセレラ社断片とPFPデータを一緒にするのに役立つだけとなった。そこでは我々は、WGAに用いたアセンブラを対象領域の当初段階からのアセンブリに用いた。

コンポーネントのWGAアセンブリにより、スパンが計2.906 Gbpで、2.654 Gbpの配列から構成される一セットの骨格が作られた。チャフ、すなわちアセンブリに組み込まれなかった読み取り配列のセットは、617万個(22%)であった。ゲノムの90%以上が100 kbp以上の骨格でカバーされ、配列は計2.492 Gbpで、平均92.2%が配列そのものであり、7.8%がギャップであった。100 kbp以上の骨格1,940個に属する107,199個のコンティグの中に、計105,264個のギャップがあった。骨格サイズは平均1.4 Mbpであり、平均コンティグサイズは23.24 kbp、平均ギャップサイズは2.0 kbpで、それぞれのサイズは指数関数的に分布している。またそのため、平均値は、大部分のデータの過小評価値となりやすい。図5に、様々なサイズ範囲の骨格中の塩基数ヒストグラムを示す。全ギャップの49%以上が500 bpよりも短く、62%以上が1 kbpよりも短く、全て100 kbpより短かったことを考慮していただきたい。同様に、配列の73%以上が30 kbp以上のコンティグの中にあり、49%以上が100 kbp以上のコンティグにあって、最大のコンティグは1.99 Mbpであった。表3に、WGAアセンブリと直接比較のため、このアセンブリの統計値の要約を示す。



図5   CSAのアセンブリ骨格サイズの分布。各サイズ範囲に属するアセンブリ骨格が全配列に占める割合(%)を示す。

表3 全ゲノムショットガンアセンブリと区画化ショットガンアセンブリのアセンブリ骨格の統計値


Scaffold size
All >30 kbp >100 kbp >500 kbp >1000 kbp

Compartmentalized shotgun assembly
No. of bp in scaffolds 2,905,568,203 2,748,892,430 2,700,489,906 2,489,357,260 2,248,689,128
(including intrascaffold gaps)
No. of bp in contigs 2,653,979,733 2,524,251,302 2,491,538,372 2,320,648,201 2,106,521,902
No. of scaffolds 53,591 2,845 1,935 1,060 721
No. of contigs 170,033 112,207 107,199 93,138 82,009
No. of gaps 116,442 109,362 105,264 92,078 81,288
No. of gaps <=1 kbp 72,091 69,175 67,289 59,915 53,354
Average scaffold size (bp) 54,217 966,219 1,395,602 2,348,450 3,118,848
Average contig size (bp) 15,609 22,496 23,242 24,916 25,686
Average intrascaffold gap size (bp) 2,161 2,054 1,985 1,832 1,749
Largest contig (bp) 1,988,321 1,988,321 1,988,321 1,988,321 1,988,321
% of total contigs 100 95 94 87 79
Whole-genome assembly
No. of bp in scaffolds (including intrascaffold gaps) 2,847,890,390 2,574,792,618 2,525,334,447 2,328,535,466 2,140,943,032
No. of bp in contigs 2,586,634,108 2,334,343,339 2,297,678,935 2,143,002,184 1,983,305,432
No. of scaffolds 118,968 2,507 1,637 818 554
No. of contigs 221,036 99,189 95,494 84,641 76,285
No. of gaps 102,068 96,682 93,857 83,823 75,731
No. of gaps <=1 kbp 62,356 60,343 59,156 54,079 49,592
Average scaffold size (bp) 23,938 1,027,041 1,542,660 2,846,620 3,864,518
Average contig size (bp) 11,702 23,534 24,061 25,319 25,999
Average intrascaffold gap size (bp) 2,560 2,487 2,426 2,213 2,082
Largest contig (bp) 1,224,073 1,224,073 1,224,073 1,224,073 1,224,073
% of total contigs 100 90 89 83 77


2.5. WGAとCSAの骨格比較

別々の計算処理プロセス(WGAとCSA)を経て、ヒトゲノムのアセンブリが2種類得られたので、これらの完全性、整合性、連結性を調べる別手段として、2つのアセンブリから得られた骨格を比較した。それぞれのアセンブリから、少なくとも1000個の断片を含む基準骨格を1セット(セレラ社の読み取り配列またはbactig断片)得た。これらは2218 個のWGA 骨格と1717 個のCSA 骨格で、各々計2.087 Gbpと2.474 Gbpであった。少なくとも20個のフラグメントか小さい方の骨格では少なくとも20%を占めている一方の基準骨格の配列を、他方のアセンブリから得られた全ての骨格の配列と比較した。各比較に関して、ミスマッチが2%以下で200 bp以上の一致があるものを全て表にした。

この表から、2つの方法で、各アセンブリの特有の配列量を推定した。第一の方法は、他方のアセンブリの一致部分とは異なる各アセンブリの塩基の数を測定することである。WGAの約82.5 Mbp(3.95%)がCSAにカバーされておらず、CSAの約204.5 Mbp(8.26%)がWGAにカバーされていなかった。この推定では、アセンブリの整合性や一致部分のユニークさは必要とされない。したがって、別の解析を行った。その解析では、一対の比較骨格間で1 kbp未満の一致配列 については、他の一致配列で矛盾のない順序と方向性が確認されない限り、除外した。これにより、一致配列数の測定ができた。より厳密なこの方法によると、WGAの1.982 Gbp(95.00%)がCSAにカバーされており、CSAの約2.169 Gbp(87.69%)がWGAにカバーされていた。

WGAをCSAと比較すると、骨格構造上の不一致度合いが評価できる。一方のアセンブリから得られた大型部分骨格で、他方のアセンブリ骨格群にあってはそのうちのたった一つとしか一致しない例で、それでも一致片によって示される配列のオーバーラップする長さが全長にわたらない例を探した。候補となる最初のセットは自動的に同定し、各候補セットを手動で調査した。このプロセスから、そのアセンブリがその場所に限定されるとは限らないと考えられる例を31個見いだした。これらの例はさらに調査し、どのアセンブリが間違っていて、なぜなのかを判断した。

さらに、順序と方向性について、局所的な不一致性を評価した。以下に述べる結果は、一方のアセンブリの1つのコンティグが、他方のアセンブリの2つ以上のコンティグに一致する場合(後者の順序と方向性が、対応する前者の場所と一致する限り)を除外している。これらの小さな配置ミスは塩基対で数百個程度のものであり、1 kbp以上のものはめったに関与していない。我々は、CSAアセンブリのうち、計295 kbp(0.012%)が、局所的にWGAアセンブリと矛盾し、WGAアセンブリのうち、2.108 Mbp(0.11%)が、CSAアセンブリと矛盾することを見いだした。

CSAアセンブリは骨格のカバー倍数の点では数パーセント良好であり、WGAよりもわずかに矛盾が少なかった。というのは、CSAは実際、メガベースサイズの課題のショットガンアセンブリを数千回行っているが、WGAはギガベースサイズの課題のショットガンアセンブリを1回行っているためである。課題サイズが2.5桁も大きいことを考えれば、この2者間の情報ロスは非常に少ない。2つの結果のうちCSAは論理的に出力しやすくベターで、その後の解析を開始しなければならない時点で入手できたことから、以降の解析は全て、このアセンブリに関して行った。

2.6. アセンブリ骨格のゲノムへのマッピング

ゲノムのアセンブリの最終段階は、染色体上にアセンブリ骨格を並べ、方向付けさせることである。最初に、CSAのコンポーネント中の順序に基づき、アセンブリ骨格をグループ化した。これらのグループ化したアセンブリ骨格を、アセンブリ骨格間に残存しているメイト対データを検討して並べ直した。次に、物理的マッピングデータを用いて、アセンブリ骨格グループを染色体上にマッピングした。この段階は、各アセンブリ骨格が複数のマーカーと重なるような、信頼性の高い高分解能マップ情報を得られるかどうかに左右される。入手可能なゲノム全体にわたるマップ情報は2つある。高密度STSマップと、ワシントン大学で作製されたBACクローンのフィンガープリントマップ(WashU BACマップ)である(45)。ゲノム全体のSTSマップのうち、GeneMap99(GM99)は最多のマーカーを有し、アセンブリ骨格をマッピングするには最も有用であった。2種類のマッピング法は、お互いに相補的なものである。フィンガープリントマップは、オーバーラップしたBACクローンの比較により作製されたため、部分的な情報が良好である。一方、GM99は、枠組みとなるマーカーが、よく実証された遺伝子マップに由来しているため、広範囲にわたる情報はより信頼性が高い。両タイプのマップを、座位アセンブリのために入力したコンポーネントの、人手による管理基準として用いたが、これらは、アセンブラにより生み出された配列順を規定するものではない。

フィンガープリントマップとGM99の、アセンブリ骨格のマッピング効率を調べるため、まず、大きなアセンブリ骨格で比較することでこれらのマップの信頼性を検討した。GM99上では、10個の大きなアセンブリ骨格(9 Mbp以上)上のSTSマーカーのうち1%のみが異なる染色体にマッピングされた。STSマーカーの2%が5フレームワーク結合分以上位置が違っていた。しかし、フィンガープリントマップでは、2%が異なる染色体上にあり、アセンブリ骨格配列におけるBAC位置の平均23.8%が、5 BAC以上、フィンガープリントマップの配置と一致しなかった。さらに不一致の根源を検討したところ、不一致のほとんどは10個のアセンブリ骨格のうち4個に由来することが明らかになった。これは、マップかアセンブリ骨格のいずれかの質にバラツキがあることを示している。4個のアセンブリ骨格とも、残り6個のアセンブリ骨格と同じようにアセンブリしており、クローン配列カバー倍数解析で評価すると、GM99との不一致率が同じように低いことが示された。このことから、これらの場合ではフィンガープリントマップの全体像は信頼性が低いと結論した。より小さいアセンブリ骨格では、GM99との不一致率は高かった(STSの4.21%が5フレームワーク結合分以上一致しなかった)が、フィンガープリントマップとの不一致率は低かった(BACの11%がフィンガープリントマップで5 BAC分以上一致しなかった)。この結果は、セレラ社のアセンブリ骨格構築では、小さいアセンブリ骨格よりも大きいアセンブリ骨格の方が長距離メイト対により、よりうまく合うという、クローン配列カバー倍数解析と一致している(46)

我々は、これらのマップのマーカー(BACまたはSTS)に基づき、セレラ社のアセンブリ骨格の並び順を2種類作製した。アセンブリ骨格の順序がGM99とWashU BACマップで一致する場合は、その順序が正しいことに確信を持ち、これらのアセンブリ骨格には、「アンカーアセンブリ骨格」と名付けた。両方のマップで全般的な不一致性が低く信頼できるアセンブリ骨格のみ、アンカーアセンブリ骨格であると判断した。GM99 側にあるアセンブリ骨格は、大枠順(framework order)を乱さない場合、順序を入れ替えてWashUの順序と一致させることが可能であった。各アセンブリ骨格の方向付けは、マッピングされた複数のマーカー(順序に矛盾がないもの)の存在により決定した。1つのマーカーしか持たないアセンブリ骨格は、方向を決めるには情報が不十分である。我々は、ゲノムの70.1%がアンカーアセンブリ骨格中にあることを見いだし(表4)、そのアンカーアセンブリ骨格は99%以上が方向付けされている。GM99はWashUマップよりも解像度が低いため、STSマーカーで一致するもののない多くのアセンブリ骨格を、アンカーアセンブリ骨格に関連づけて並べた。それらが、WashU BACマップ上の同じか隣接したBAC由来の配列を持っていたためである。また、WashU BACマップ上に時折みられるおかしな並べ方のため、WashUマップ上に「マッピングできない」と判断された多くのアセンブリ骨格を、アンカーアセンブリ骨格に関連づけてGM99で並べることができた。これらのアセンブリ骨格には、「整列化アセンブリ骨格」と名付けた。13.9%のアセンブリ骨格がこれらの補足的方法で並べることができた。このようにして、ゲノムの84.0%を確実に整然と並べた。

表4 アセンブリ骨格マッピングのまとめ。アセンブリ骨格は様々な信頼度でゲノムにマッピングされた(アンカーアセンブリ骨格が最も信頼性が高く、マッピングされないアセンブリ骨格が最も信頼性が低い)。アンカーアセンブリ骨格は、WashU MADCマップとGM99により、矛盾なく並べられた。整列化アセンブリ骨格は、WashU BACマップか、GM99か、コンポーネントタイリング路程のうち1つ以上に基づいて並べられた。拘束アセンブリ骨格は、少なくとも2つの社外マップの間では順序が矛盾するが、これらの位置は近隣のアンカーアセンブリ骨格または整列化アセンブリ骨格に隣接している。マッピングされないアセンブリ骨格は、所属染色体の記載がなされている程度である。アセンブリ骨格のサブカテゴリーは、それぞれのカテゴリーの下に示してある。


Mapped scaffold category Number Length (bp) % Total length

Anchored 1,526 1,860,676,676 70
Oriented 1,246 1,852,088,645 70
Unoriented 280 8,588,031 0.3
Ordered 2,001 369,235,857 14
Oriented 839 329,633,166 12
Unoriented 1,162 39,602,691 2
Bounded 38,241 368,753,463 14
Oriented 7,453 274,536,424 10
Unoriented 30,788 94,217,039 4
Unmapped 11,823 55,313,737 2
Known 281 2,505,844 0.1
chromosome
Unknown chromosome 11,542 52,807,893 2


次に、アンカーアセンブリ骨格の間に配置はできたものの並び順を決めることができなかった全てのアセンブリ骨格を、アンカーアセンブリ骨格間にある間隔体とし、アンカーアセンブリ骨格間に「拘束されている」と考えた。例えば、小さなアセンブリ骨格で、GeneMap の同じ箱(bin)に対してSTS 該当部位を持っているか、または同じBACに当たる部位があるものは、相対的に並べることができないが、他のアンカーアセンブリ骨格か整列化アセンブリ骨格に関連づけて、境界域に割り当てることができる。残りのアセンブリ骨格は、位置情報や、矛盾情報を持たず、一般的な染色体座位に割り当てることしかできなかった。以上のアプローチを用いて、ゲノムの約98%が、位置固定されるか、並べられるか、間隔拘束された。

最後に、染色体ごとにアセンブリ骨格を拡散展開することによって、染色体上に配置されたアセンブリ骨格の局在位置を確定した。ゲノムの2%にあたる、マッピングされなかった残りのアセンブリ骨格は、ゲノム全体を通じて均等に分配されると仮定した。マッピングされなかったアセンブリ骨格長の合計を、マッピングされたアセンブリ骨格数の合計で割ることにより、アセンブリ骨格間のギャップが1,483 bpであると推定した。このギャップを用いて、各染色体上の全アセンブリ骨格を分別し、染色体内にきちんと配置できないものをオフセットとして設けた。

アセンブリ骨格をマッピングする作業の間、我々は多くの問題に遭遇し、追加的な品質評価とバリデーション解析を行うはめになった。少なくとも978個(33,173個の3%)のBACがゲノム内の2ヶ所以上からの配列データを持つことが判明した(47)。これはアセンブリ戦略の項で述べたbactigキメラ解析と一致する。これらのBACは、CSAアセンブリの中でユニークな位置に割り当てることはできず、そのためアセンブリ骨格を並べるために利用することもできなかった。似たようなものでSTSも、ゲノム重複、繰返し配列、偽遺伝子があるために、アセンブリで常にユニークな位置に割り当てることができたわけではない。

完璧なオーバーラップを徹底的に検索するのには時間がかかるため、CSAでは、21,607個のアセンブリ骨格内ギャップを作り出してしまった。メイト対データが当該コンティグはオーバーラップすべきだと示唆しているのに、実際にはオーバーラップが検出されなかったギャップである。これらのギャップは長さ50 bpの固定長と定義されており、CSAアセンブリの計116,442個のギャップの18.6%を占めている。

アセンブリ骨格を並べる手段として、cDNAやESTデータで示されているエキソンの並び順は用いないことにした。このデータを利用しなかったのは、利用した場合、翻訳データにあわせるためにアセンブリ骨格を並べ直すことにより、アセンブリの特定の領域に偏りができ、アセンブリと遺伝子の定義プロセスのバリデーションがより困難になるためである。

2.7. アセンブリとバリデーション解析

我々は、完成度(ゲノムカバー度)と精度(順序と方向の構造的正確さとアセンブリのコンセンサス配列)の観点から、ゲノムのアセンブリを解析した。

完成度

完成度は、アセンブリ中の真正染色質の配列の割合として定義される。これは、真正染色質の配列決定が完了するまでは絶対的な確実性をもって知ることはできない。しかし、次のことに基づいて完成度を推測することができる:(i)アセンブリ骨格内のギャップの推定サイズ、(ii)、既に配列が発表されている第21、22番染色体におけるカバー度(48,49)、(iii)アセンブリに含まれる別個のランダム配列のセット(STSマーカー)の割合の解析。全ゲノムライブラリーには、ヘテロクロマチン配列が含まる。それをアセンブリする試みは行われたことがないが、キイロショウジョウバエでみられたように、ヘテロクロマチンの領域に埋め込まれた特有の配列の例がある可能性がある(50,51)

ヒト第21、22番染色体の配列決定は、高品質で完了し、発表されている(48,49)。この配列もアセンブラに入力したが、配列決定の完了した配列を切断してショットガンデータセットにし、構造的な多型やBACデータのアセンブリミスがあった場合、アセンブラが元の配列とは違うアセンブリをするようにした。特に、アセンブラは、コンポーネント長のスケール(一般的に数メガベースサイズ)で、反復配列を解析できなければならない。しかも、この2者比較は、アセンブラが反復配列を解析できる解像度レベルを明らかにする。特定の領域では、アセンブリ構造は、発表されているヒト第21、22番染色体のものとは異なる(後述参照)。セレラ社データに基づき、配列決定の“終わったはず”の配列を、この違った形で柔軟にアセンブリした結果、第21、22番染色体の配列よりも多いDNAセグメントアセンブリとなった。なぜ第21、22番染色体配列よりもセレラ社配列ではギャップが多いのか理由を検討し、それらはゲノムの他領域におけるギャップに特有のものであるかもしれないと期待した。セレラ社アセンブリでは、25個のアセンブリ骨格があり、それぞれ10 kbp以上の配列から構成され、全部集めると第21番染色体の94.3%をカバーしている。また、62個のアセンブリ骨格が第22番染色体の95.7%をカバーしている。これら2つの染色体のセレラ社アセンブリで残っているギャップ長は合計3.4 Mbpである。これらのギャップ配列は反復配列目隠し機(RepeatMasker)を使い、全ゲノムアセンブリに対して検索を行って解析した(52)。ギャップ配列の約50%が、反復配列目隠し機で同定される一般的反復配列で構成されていた。残りの半分以上は、コピー数の少ない反復配列であった。

完成度を評価するもっと包括的な方法は、アセンブリに含まれる全く独自の配列データセットの含量を測定することである。Genemap99由来の48,938個のSTSマーカー(51)をアセンブリ骨格と比較した。これらのマーカーはアセンブリプロセスでは利用されていなかったため、これらは完成度について真に独自性のある基準となった。ePCR(53)とBLAST(54)を用いて、アセンブリされたゲノムにSTSを配置した。我々は、マッピングされたゲノムに、44,524個(91%)のSTSを見いだした。その他に2648個のマーカー(5.4%)が、アセンブリされていないデータ、すなわち「チャフ」の検索で見いだされた。我々は、セレラ社配列にも、2000年9月時点のBACデータにも見いだされない1283個のSTSマーカー(2.6%)を同定した。このことから、これらのマーカーがヒト由来でない可能性が考えられる。もしそうならば、アセンブリされたセレラ社配列はヒトゲノムの93.4%を表していることになり、アセンブリされていないデータは5.5%で、カバー度合いは全部で計98.9%となる。同様に、CSAを36,678個のTNG放射線ハイブリッドマーカー(55a)に対して比較した。32,371個のマーカー(88%)が、マッピングされたCSA アセンブリ骨格に位置し、2055個のマーカー(5.6%)が残りの位置に発見された。この場合は、後者のゲノム全体調査を通じて、ゲノムの94%のカバー度となる。

精度

精度は、アセンブリの構造的・配列的正確さとして定義される。セレラ社データとGenBankデータの配列原材料は異なる人間に由来するため、配列決定の終了した他の配列と比較して、直接、アセンブリのコンセンサス配列のヌクレオチドレベルで正確さを調べることはできない。しかし6章に示すように、この作業は多型同定のため行われている。コンセンサス配列の精度は、読みとり配列の品質値から導き出される統計値をもとにすると、99.96%以上である。

アセンブリが構造的に正しいということは、メイト対解析で測定できる。正しいアセンブリでは、読み取り配列の全てのメイト対が、正しい距離と向きをもって、コンセンサス配列に配置されるはずである。ある一組の読みとり配列ペアを見て、それらが正しく方向付けされており、読みとり配列間の距離が、その読みとり配列がサンプリングされたライブラリのインサートサイズ分布の平均値(3標準偏差値以内であれば、その読みとり配列ペアは「正当(valid)」と呼ばれる。正しく方向付けされていないときは、その読みとり配列ペアは「方向ミス(misoriented)」と呼ばれ、読みとり配列間の距離が正しい範囲にないが正しく方向付けされている場合はその読みとり配列ペアは「分離ミス(misseparated)」と呼ばれる。アセンブラが用いる各ライブラリの平均(標準偏差はこのように決定される。これらを評価するため、第21番染色体の決定終了配列(48)にマッピングされた全ての読みとり配列を検討し、試験室内追跡エラーとキメラ化(ゲノムの2つの異なる部位が同一プラスミドにクローニングされている)の結果として不正メイト対がどのぐらいの量あったか、そして、正確な場合はインサートサイズの分布がどの程度きっちりしていたかを決定した(表5)。セレラ社の全ライブラリーの標準偏差は非常に小さく、少数の50 kbpライブラリーを除き、インサート長の15%以下であった。2 kbpおよび10 kbpのライブラリーは不正メイト対を2%未満しか含まなかったが、50 kbpのライブラリーは幾分多かった(約10%)。メイト対情報は完璧ではなかったにも関わらず、その精度は、特に数個のメイト対で方向性を確認するか否定するかする場合、対象アセンブリに関して正当・方向ミス・分離ミスの各メイト対量測定は、バリデーションにおいて信頼性の高い手段となると考えられる程であった。

表5 メイト対のバリデーション。セレラ社のフラグメント配列を、報告されている第21番染色体の配列にマッピングした。独自ににマッピングされた各メイト対を、方向付けと位置が正しいかどうか評価した(テストしたメイト対の数)。2つのメイトの相対方向や位置が間違っていた場合は、不正と判断した(不正メイト対の数)。BES:BAC末端配列

 


Library type Library no. Chromosome 21
Genome
Mean insert size (bp) SD (bp) SD/mean (%) No. of mate pairs tested No. of invalid mate pairs % invalid Mean insert size (bp) SD (bp) SD/ mean (%)

2 kbp 1 2,081 106 5.1 3,642  38 1.0 2,082 90 4.3
2 1,913 152 7.9 28,029 413 1.5 1,923 118 6.1
3 2,166 175 8.1 4,405  57 1.3 2,162 158 7.3
10 kbp 4 11,385 851 7.5 4,319  80 1.9 11,370 696 6.1
5 14,523 1,875 12.9 7,355 156 2.1 14,142 1,402 9.9
6 9,635 1,035 10.7 5,573 109 2.0 9,606 934 9.7
7 10,223 928 9.1 34,079 399 1.2 10,190 777 7.6
50 kbp 8 64,888 2,747 4.2 16   1 6.3 65,500 5,504 8.4
9 53,410 5,834 10.9 914 170 18.6 53,311 5,546 10.4
10 52,034 7,312 14.1 5,871 569 9.7 51,498 6,588 12.8
11 52,282 7,454 14.3 2,629 213 8.1 52,282 7,454 14.3
12 46,616 7,378 15.8 2,153 215 10.0 45,418 9,068 20.0
13 55,788 10,099 18.1 2,244 249 11.1 53,062 10,893 20.5
14 39,894 5,019 12.6 199   7 3.5 36,838 9,988 27.1
BES 15 48,931 9,813 20.1 144  10 6.9 47,845 4,774 10.0
16 48,130 4,232 8.8 195  14 7.2 47,924 4,581 9.6
17 106,027 27,778 26.2 330  16 4.8 152,000 26,600 17.5
18 160,575 54,973 34.2 155   8 5.2 161,750 27,000 16.7
19 164,155 19,453 11.9 642  44 6.9 176,500 19,500 11.05
Sum 102,894 2,768 2.7
(mean = 2.7)


ゲノムの配列クローンカバー倍数は39倍であった。これは、どの塩基対も平均39個のクローンに含まれる、言い換えれば、メイト対でつなぎ合わされた39個の読みとり配列のどこかにいることを意味している。配列クローンカバー倍数が低い領域、または不正メイト対の比率が高い領域は、アセンブリに問題がある可能性がある。正当メイト対を用いて、アセンブリの各塩基のカバー倍数を計算した(表6)。まとめると、30 kbp以上のアセンブリ骨格では、セレラ社アセンブリのうち配列クローンカバー倍数が3倍未満の領域にあるのは1%未満であった。このように、順序と方向を含めて、この測定をするだけでアセンブリの99%以上が強い支持を受けるのである。

表6 区画化ショットガンアセンブリ(CSA)とPFPアセンブリのゲノム全体のメイト対解析*


Genome library CSA
PFP
% valid % mis-oriented

T表6 区画化ショットガンアセンブリ(CSA)とPFPアセンブリのゲノム全体のメイト対解析.*


Genome library CSA
PFP
% valid %mis-oriented % mis-separateddagger % valid %mis-oriented % mis-separateddagger

2 kbp 98.5 0.6 1.0 95.7 2.0 2.3
10 kbp 96.7 1.0 2.3 81.9 9.6 8.6
50 kbp 93.9 4.5 1.5 64.2 22.3 13.5
BES 94.1 2.1 3.8 62.0 19.3 18.8
Mean 97.4 1.0 1.6 87.3 6.8 5.9

*Data for individual chromosomes can be found in Web fig. 3 onScience Online atwww.sciencemag.org/cgi/content/full/291/5507/1304/DC1.
dagger Mates are misseparated if their distance is >3 SD from the mean library size.


全ての方向ミスメイト対・分離ミスメイト対の場所と数を調べた。この解析をCSAアセンブリ(2000年10月1日の時点)に関して行ったのに加えて、2000年9月5日の時点(30, 55b)でのPFPアセンブリの調査も行った。後者の場合、セレラ社のメイト対をPFPアセンブリにマッピングしなければならなかった。高忠実度反復配列によるマッピングのエラーを避けるために、マッピングするメイト対は、両方の読みとり配列が6%未満の違いで一箇所だけで一致(match)しているもののみとした。閾値として5個以上の不正メイト対セットが存在すれば、そこは潜在的分断点(potential breakpoint)を示すものと設定した。そこでは2つのアセンブリの構成が異なっていた。発表されている配列とCSAの第21番染色体アセンブリとを図示して比較すれば(図6A)この方法論のバリデーションの一つとなる。図中の青いチェックマークは分断点を示す。双方の染色体配列に、同様な(少数の)分断点があった。例外はセレラ社のアセンブリの12セットのアセンブリ骨格で(212個の単一コンティグアセンブリ骨格にあって全部で染色体長の3%)、間違った位置にマッピングされていた。これらは小さすぎて信頼のおけるマッピングができなかったからである。図6と7と表6に、2つのアセンブリのメイト対と分断点の差を示す。両方のアセンブリとも、大きいインサートのライブラリー(50 kbpとBAC末端)の方が、小さいインサートのライブラリーより、方向ミス・分離ミスメイト対の割合が高かった(表6)。大きいインサートのライブラリーは、ゲノムのより大きな部分に及ぶという単純な理由から、不一致を検出しやすかった。第8番染色体の2つのアセンブリを図形的に比較(図6、BとC)したところ、セレラ社アセンブリよりもPFPアセンブリに分断点が多いことが示された。図7に、各染色体の両方のアセンブリの分断点地図(青いチェックマーク)を横に並べて示す。セレラ社アセンブリの順序と方向付けは、2つの配列決定終了染色体を除けば、分断点が実質的に少ない。図7には、両方のアセンブリの大きいギャップ(10kbp以上)も赤いチェックマークで示してある。CSAアセンブリでは、全てのギャップのサイズはメイト対データに基づいて推定されている。2つのアセンブリは、異なるヒトゲノムを用いているため、分断点は、構造多型により生じている可能性がある。これらは、双方のゲノムアセンブリの未完成状態をも反映している。

図6
図7




1 Celera Genomics, 45 West Gude Drive, Rockville, MD 20850, USA. 2 GenetixXpress, 78 Paci_c Road, Palm Beach, Sydney 2108, Australia. 3 Berkeley Drosophila Genome Project, University of California, Berkeley, CA 94720, USA. 4 Department of Biology, Penn State Uni-versity, 208 Mueller Lab, University Park, PA 16802, USA. 5 Department of Genetics, Case Western Reserve University School of Medicine, BRB-630, 10900 Euclid Avenue, Cleveland, OH 44106, USA. 6 Johns Hopkins University School of Medicine, Johns Hopkins Hospital, 600 North Wolfe Street, Blalock 1007, Baltimore, MD 21287_4922, USA. 7 Rockefeller University, 1230 York Avenue, New York, NY 10021_6399, USA. 8 New England BioLabs, 32 Tozer Road, Beverly, MA 01915, USA. 9 Division of Biology, 147-75, California Institute of Technology, 1200 East California Boulevard, Pasa- dena, CA 91125, USA. 10 Yale University School of Medicine, 333 Cedar Street, P.O. Box 208000, New Haven, CT 06520_8000, USA. 11 Applied Biosystems, 850 Lincoln Centre Drive, Foster City, CA 94404, USA. 12 The Institute for Genomic Research, 9712 Medical Center Drive, Rockville, MD 20850, USA. 13 Faculty of Life Sciences, Bar-Ilan University, Ramat-Gan, 52900 Israel. 14 Grup de Recerca en Informa `tica Me`dica, In-stitut Municipal d'Investigacio _ Me `dica, Universitat Pompeu Fabra, 08003-Barcelona, Catalonia, Spain.

連絡先:To whom correspondence should be addressed.
E- mail: humangenome@celera.com

Copyright © 2001 by The American Association for the Advancement of Science.