JP6624736B2 - インフルエンザウイルス解析装置、インフルエンザウイルス解析方法及びインフルエンザウイルス解析プログラム - Google Patents

Description

本発明は、次世代シーケンサで解読したインフルエンザウイルスの大量なデータ解析を行うことができるインフルエンザウイルス解析装置、インフルエンザウイルス解析方法、及びインフルエンザウイルス解析プログラムに関する。

H5N1亜型高病原性鳥インフルエンザウイルス（HPAIV）は、２００３年以降世界各地で家禽での発生を引き起こし、アジアを中心に常在化の様相を呈している。２０１４年には、H5N8亜型HPAIがユーラシア大陸だけでなくアメリカ大陸までに侵入し、HPAIVの世界的蔓延による食糧の供給にも影響が及ぶことが懸念される。我が国においても、HPAIVの侵入を数度経験しており、特に２０１０年−２０１１年にかけては家禽及び野鳥での全国的な発生が認められ、２０１４年にもH5N8亜型HPAIVの侵入が認められた。これらHPAIVの防疫のために、常に監視していなければならない状況である。

国際獣疫事務所において、高病原性鳥インフルエンザウイルス（HPAIV）は、A型インフルエンザウイルスで鶏に病原性を示すもの又はHA亜型がH5またはH7亜型で赤血球凝集素タンパクの開裂部位のアミノ酸配列が既知のHPAIVと類似しているものと定義されている。我が国の法律においては、家禽での症状の有無にかかわらず、家禽に感染したA型インフルエンザウイルスのHA亜型がH5またはH7型であることが確認されれば摘発淘汰の対象となる。従って、A型インフルエンザウイルスの亜型の決定を正確かつ迅速に決定することが大変重要である。

また、インフルエンザウイルスは家畜衛生分野だけでなく、公衆衛生学的にも大変重要な疾病である。人による季節性インフルエンザのサーベイランスにより得られたA型及びB型インフルエンザウイルスは、季節性インフルエンザ予防・治療のためのワクチン選定、薬剤耐性解析等に用いられるため、大量のウイルスデータの解析を迅速に行うことが必須となる。さらに、人でのインフルエンザパンデミックの発生時には、そのウイルスの亜型同定は、迅速なパンデミックワクチン開発に重要である。

最近の技術の進歩により、次世代シーケンサにより遺伝子を解読した後、例えば、特許文献１〜３に示されるような、次世代シーケンサに付帯するソフトウエア又は特殊な解析ソフトを用いてデータを解析している。これにより、次世代シーケンサによりインフルエンザウイルスの遺伝子を大量に解析することが可能となった。

特開２００５−１３５０５３号公報 特開２００６−０３９８６７号公報 特開２００９−１３１２４２号公報

しかし、インフルエンザウイルスの解析に上記従来技術を適用しても、解析にはゲノムが８つの分節に分かれているという特徴を持つインフルエンザウイルスに特異的なストラテジーを要し、それぞれの遺伝子分節ごとにトライアンドエラーを繰り返す必要があり、手間と時間を要するため、緊急な場合には、インフルエンザウイルスの解析を迅速に行うことができずに、役に立たない場合がある。

また、上記従来技術では、大量の解読されたデータを取り扱う為、インフルエンザウイルスの解析を迅速に処理できないという問題もある。

本発明は、上述した課題を解決するために創案されたものであり、インフルエンザウイルスの解析を正確かつ迅速に行うことができるインフルエンザウイルス解析装置、インフルエンザウイルス解析方法及びインフルエンザウイルス解析プログラムを提供することを目的とする。

上記目的を達成するために、本発明のインフルエンザウイルス解析装置は、参照用インフルエンザウイルスデータが入力されるリファレンスデータ入力部と、前記リファレンスデータ入力部に入力された参照用インフルエンザウイルスデータの塩基配列を複数の分節に分類して分節毎にリファレンス配列を集積する分節分割部と、検体の塩基配列が入力される検体データ入力部と、前記分節毎に集積されたリファレンス配列に前記検体の塩基配列をマッピングするマッピング処理部と、前記マッピング処理部によるマッピング結果からリファレンス配列の塩基長に対するマッピング済領域の塩基長の割合であるカバー率を算出するカバー率算出部と、前記算出されたカバー率の最大の値又は上位３位までの高カバー率を有するリファレンス配列のマッピング結果に基づいて変異を反映したコンセンサス配列を分節毎に作成するコンセンサス配列作成部と、分節毎に前記コンセンサス配列と前記集積されたリファレンス配列の相同性検索を実行し所定の類似度以上のリファレンス配列のうち最大の塩基長を有する１つの改定リファレンス配列を決定する相同性検索処理部とを備え、前記マッピング処理部において改定リファレンス配列に前記検体の塩基配列をマッピングしてマッピング結果に基づいて変異を反映した改定コンセンサス配列を作成することを主要な特徴とする。

また、本発明のインフルエンザウイルス解析方法は、参照用インフルエンザウイルスデータの塩基配列を複数の分節に分類して分節毎にリファレンス配列を集積する分節分割工程と、前記分節毎に集積されたリファレンス配列に検体の塩基配列をマッピングする第１のマッピング処理工程と、前記第１のマッピング処理工程によるマッピング結果からリファレンス配列の塩基長に対するマッピング済領域の塩基長の割合であるカバー率を算出するカバー率算出工程と、前記算出されたカバー率の最大の値又は上位３位までの高カバー率を有するリファレンス配列を分節毎に抽出する抽出工程と、前記最大のカバー率又は上位３位までの高カバー率を有するリファレンス配列のマッピング結果に基づいて変異を反映したコンセンサス配列を分節毎に作成するコンセンサス配列作成工程と、分節毎に前記コンセンサス配列と前記集積されたリファレンス配列の相同性検索を実行し所定の類似度以上のリファレンス配列のうち最大の塩基長を有する１つの改定リファレンス配列を決定する相同性検索処理工程と、前記改定リファレンス配列に前記検体の塩基配列をマッピングしてマッピング結果に基づいて変異を反映した改定コンセンサス配列を作成する第２のマッピング処理工程とを備えたことを主要な特徴とする。

また、本発明のインフルエンザウイルス解析プログラムは、上記インフルエンザウイルス解析装置としてコンピュータを機能させるためのコンピュータプログラムである。

本発明によれば、２度のマッピング処理を行うことにより、インフルエンザウイルスの解析を迅速かつ正確に行うことができる。また、インフルエンザウイルスゲノムの全長配列を決定することができる。

インフルエンザウイルス解析装置のブロック構成を示す図である。 FASTA形式ファイルのフォーマットを示す図である。 公開データより取得したFASTA形式ファイルの一例を示す図である。 参照用インフルエンザウイルスデータを分節毎に分割する方法を示すフローチャート図である。 分節に分類する場合の判定に用いられる正規表現を示す図である。 分節分割処理において、NAの分節に属すると判断される場合、いずれの分節にも属さないと判断される場合の具体例を示す。 NA分節として集積された複数のリファレンス配列が記録されているファイルの一例を示す図である。 FASTQ形式ファイルのフォーマットとクオリティ値をチェックする状態を示す図である。 NA分節における相同性検索の結果例を示す。 HA分節に対応する改定コンセンサス配列とha分節に属するリファレンス配列との相同性検索を行った結果例を示す。 図１０の正規表現適応部分に基づき、型判定を行う方法を示す。 インフルエンザウイルス解析方法のフローチャートを示す図である。

以下、図面を参照して本発明の一実施形態を説明する。

図１は、本発明のインフルエンザウイルス解析装置１００の構成例を示すブロック図で
ある。

インフルエンザウイルス解析装置は、検体データ入力部１、トリミング処理部１２、マッピング処理部２、カバー率算出部３、コンセンサス配列作成部４、相同性検索処理部５、型判定部６、登録用データ作成部７、リファレンスデータ入力部１１、分節分割部１３、制御部１０等から構成される。検体データ入力部１には、シーケンサ２１で読み取られた検体の塩基配列データが入力される。一方、あらかじめ、公開データベース２２等から参照用のインフルエンザウイルスの塩基配列データをリファレンスデータ入力部１１で取得しておき、これを分節分割部１３で、所定の分類方法により分割を行って集積し、インフルエンザウイルスの分節毎に集積されたリファレンス配列にしておく。また、分節分割部１３で分節毎に集積されたリファレンス配列に対し、相同性分析部１３１で分節毎にリファレンス配列相互の類似度をすべて算出し、類似性の高いリファレンス配列を除外し、一定の類似性以下になるように調整する。分節毎にその後、検体データをリファレンス配列にマッピングする等の処理を行いながら、検体のインフルエンザウイルスの解析を進めることになる。

なお、制御部１０は、検体データ入力部１、トリミング処理部１２、マッピング処理部２、カバー率算出部３、コンセンサス配列作成部４、相同性検索処理部５、型判定部６、登録用データ作成部７、リファレンスデータ入力部１１、分節分割部１３の各部を制御し、表示端末部２３のデータ入力を受け付けて制御や演算等を行ったり、表示端末部２３へのデータ出力や表示制御等を行う。

図１２は、インフルエンザウイルス解析方法を示すフローチャートであり、図１〜図１２を用いて以下具体的にインフルエンザウイルス解析処理について説明する。

あらかじめ、インフルエンザウイルス公開データを用いてインフルエンザウイルスをインデックス化する。インデックス化は、インフルエンザウイルスデータを８個の各分節に分割し、各分節に塩基配列データを作成することにより行う。これらの８分節に分割されたリファレンス配列は、分節毎に集積される。

元となる公開データは、例えば、NCBI（http://www.ncbi.nlm.nih.gov/)やGISAID（http://platform.gisaid.org/epi3/frontend#60d323）などのインフルエンザウイルス情報をもつ公開データを利用する。リファレンス配列作成の概要は、まず、上記等の公開サイトからインフルエンザウイルスFASTA形式ファイルをリファレンスデータ入力部１１にダウンロード等の方法で取得する。取得したファイルを特定のルールに従い、８種類の分節とその他のファイルという形でそれぞれに分割して集約する。次に、８種類の分節配列データを解析に用いることができる特定のフォーマットに変換する。

FASTA形式ファイルとは、良く知られているように、塩基／アミノ酸の配列を記述したテキストファイルであり、行頭が「>」にて始まる配列名と、配列名行の改行後に続く塩基又はアミノ酸行から構成される。

具体的にインフルエンザウイルスのインデックス化について述べる。まず、公開サイト等から取得したFASTA形式のインフルエンザウイルスデータは、図２のフォーマットで記載されている。図２のように、「>」記号に続けて記載されている文字列が当該配列の配列名となり、次行から当該配列名に対応する塩基配列が記載される。

FASTA形式のインフルエンザウイルスデータの具体的な一例を図３に示す。これは、ftp://ftp.ncbi.nlm.nih.gov/genomes/INFLUENZA/influenza.fna.gzの参照配列データベースのリンクよりダウンロードしたファイルを解凍したものである。ファイルの拡張子がfnaとなっているのは、FASTA形式のファイルであることを示し、gzはアーカイブの形式を示す。なお、すべての塩基データを掲載することはできないので、一部は省略している。

図４は、参照用インフルエンザウイルスデータを分節毎に分割する方法を示すフローチャート図である。図３のinfluenza.fna.gzファイルをリファレンスデータ入力部１１にシーケンシャルに読み込み（Ｓ１）、「>」記号に続けて記載されている配列名を取得する。この配列名を利用して、分節分割部１３で８個の分節塩基配列と、その他の塩基配列とに分割する。

Ｓ２に示されるように、ファイルヘッダを参照する。すなわち、配列名の部分を参照し
、PB1、PB2、PA、HA、NP、NA、MP、NSの８個の分節と、分節が識別できない塩基配列とに
分類する。

配列名の文字列に、正規表現１〜８のいずれかが存在するかを判定していき、その正規表現が含まれる場合は、識別コードを付与する。正規表現１〜８は、それぞれ、図５（ａ）〜（ｈ）に対応するもので、正規表現１から順に、PB1、PB2、PA、HA、NP、NA、MP、NS各分節に対応している。図５に示された正規表現は、／で挟まれた文字列が１個のキーワードとなっている。

したがって、正規表現１が含まれる場合とは、図５（ａ）示すように、[ |(|,]pb1[ |)|,]、basic polymerase 1、polymerase basic 1、polymerase basic protein 1、segment 2、segment: 2の６種類のキーワードのいずれかが配列名に含まれる場合を示す。正規表現２〜８についても同様であり、図５（ｂ）〜（ｈ）に記載されている各キーワードが含まれる場合は、それぞれに対応する分節に属すると判定する。

図４のＳ３では、最初に正規表現１が配列名に含まれるか判断し、含まれる場合（YES）には、PB1分節の識別コードを当該塩基配列に付与した（Ｓ４）後、ステップＳ５に進み、含まれない場合（NO）は、そのまま次のステップＳ５に進む。

Ｓ５では、正規表現２が配列名に含まれるか否かを判断し、含まれる場合（YES）には、pb2分節の識別コードを当該塩基配列に付与した（Ｓ6）後、ステップＳ７に進み、含まれない場合（NO）は、そのまま次のステップＳ７に進む。このように、Ｓ７では正規表現３が配列名に含まれるか否かを、Ｓ９では正規表現４が配列名に含まれるか否かを、Ｓ１１では正規表現５が配列名に含まれるか否かを、Ｓ１３では正規表現６が配列名に含まれるか否かを、Ｓ１５では正規表現７が配列名に含まれるか否かを、Ｓ１７では正規表現８が配列名に含まれるか否かを次々に判断する。

そして、各々のステップにおいて当該正規表現が含まれる場合（YES）には、Ｓ８ではPA分節の識別コードを、Ｓ１０ではHA分節の識別コードを、Ｓ１２ではNP分節の識別コードを、Ｓ１４ではNA分節の識別コードを、Ｓ１６ではMP分節の識別コードを、Ｓ１８ではNS分節の識別コードを、当該塩基配列に付与する。

１つの塩基配列に対して、Ｓ３〜Ｓ１８の処理が終了すると、Ｓ１９で、識別コードが割り当てられており、かつ、それが単一の識別コードであるか否かを判断する。YESの場合は、Ｓ２０に示すように、識別コードを参照して、当該塩基配列を割り当てられた識別コードに対応する分節の記憶領域に格納する。NOの場合は、Ｓ２１のように、識別できない塩基配列として、８種類の分節とは異なるその他の領域に格納する。以上のようにして、取得した参照用インフルエンザウイルスデータのすべての塩基配列の分類が終了したか否かを判断し（Ｓ２２）、終了していなければＳ２からの処理に戻り、終了していれば、分節分割に関する処理を終了する。これにより、分節分割部１３には、分節毎にリファレンス配列データが集積されて格納された、分節毎リファレンス配列データ記憶部１３Ａが形成される。

図６は、図４の分節分割処理において、NAの分節に属すると判断される場合、いずれの分節にも属さないと判断される場合の具体例を示す。図６（ａ）に示される配列名には、neuraminidaseのキーワードが存在し、それ以外のキーワードがないので、NAの分節に分類される例である。図６（ｂ）は、図５（ａ）〜（ｈ）のいずれの正規表現も存在しないので、その他に分類される例である。図６（ｃ）は正規表現のキーワードが、segment7というMP分節に関するキーワード、nonstructural proteinというNS分節に関するキーワードが存在するため、２種類以上の分節に属することになり、単一の分節に決定することができず識別できないので、その他に分類される例である。

このように、いずれの分節にも属さない、又は、２つ以上の分節に属する場合は、その他の識別できない塩基配列に分類し、解析に用いる分節からは除く。

図７は、インフルエンザウイルスの参照配列データベースのリンクよりダウンロードしたファイルから、８種類の分節に分割したもののうち、NA分節として集積された複数のリファレンス配列が記録されているファイルの一例を示す。この例では、４つの塩基配列が記録されている。

以上により、８個の分節FASTA形式ファイルと、１個のいずれの分節にも属さないその他のFASTA形式ファイルとの合計９個のファイルが作成されるので、分節分割部１３の分節毎リファレンス配列データ記憶部１３Ａでこれらのデータを保持しておく。次に、Ｓ２３に示すように、Ｓ２０で作成されたファイルを用い、各分節内においてリファレンス配列として集積された塩基配列間相互のすべての相同性分析を相同性分析部１３１で行い、類似度が高い塩基配列を除外する。これにより、１つの分節内における塩基配列相互の類似度が一定のレベル（80％から90％）以下になるように調整される。これは、分節内における２つの塩基配列のすべての組み合わせについて相同性を計算した場合、そのいずれの数値も一定のレベル以下になることを意味する。各分節内において、前記のように塩基配列相互の類似度が一定のレベル（80％から90％）以下になる調整されたリファレンス配列は、分節毎リファレンス配列データ記憶部１３Ｂで保持しておく。この作業には一般的に配列間の相同性が指定されたパーセンテージになるように配列の取捨選択を行うフリーツールが利用される。一例として、CD-HITを用いた。CD-HITでは、例えば、類似度が高い塩基配列１、塩基配列２、塩基配列３の各塩基配列相互の類似度が一定のレベルを超えている（類似度が高い）場合は、この３つの塩基配列のうち塩基長の最も長い配列を残し、他の２つの塩基配列は捨てるという処理を行う。このように、塩基長が長い配列を残すことにより、公開データベース２２に登録されているインフルエンザウイルスには配列の読み取りを行った登録者がウイルスゲノムの全長を完全にシーケンスしていなくても登録した質の低い塩基配列が含まれているため、これらの配列を除去することができる。また、類似度の算出方法は、一般によく知られた方法が使用されており、二つの塩基配列間で塩基の並びを少しずつずらして照合し、もっとも塩基配列間の一致率が高くなった数値を使用している。なお、上記の類似度を調整するレベル（閾値）を下げすぎると、抽出されるリファレンス配列の数がかなり減少し、自然界にこれまで想定していなかったウイルスの配列が出現した場合、リファレンス配列と大幅に配列が異なるために型を判定できない危険性がある。

一方で、シーケンサ２１で検体のシーケンシングを実施する。例えば、シーケンサ２１は、次世代シーケンサと呼ばれる装置を用い、シーケンシングは、イルミナ社MiSeqもしくはその他のDNAシーケンス装置のプロトコールに従って行われる。シーケンシングした結果をFASTQ形式のデータフォーマットで出力し、検体データ入力部１に入力することで、ＳＴ１の検体データの取り込みが行われる。ＳＴ１へ取り込まれるシーケンシング結果のデータフォーマットは、FASTQ形式にかかわらず塩基配列のみ又は塩基配列と各塩基のクオリティスコアが含まれたデータ形式のいずれも対応可能である。これは、後述するマッピングソフトでマッピング可能なデータフォーマットとするためである。

FASTQ形式ファイルとは、良く知られているように、塩基配列と各塩基のクオリティスコアを記述したテキストファイルであり、行頭が「@」にて始まる４行に対して、１リード分の塩基配列情報から構成される。具体例を図８（ａ）に示す。各行は、以下のように構成されている。

１行目は＠から始まるリード識別名、２行目はリードの塩基、３行目は＋から始まるリード機器別名（リード機器別名自体は省略可能）、４行目は塩基のクオリティスコア、となっている。なお、なお、具体例では、４行目のクオリティスコアは、ASCIIコードで表示されたフレッドスコアとなっており、２行目と同じ文字数になる。

トリミング処理部１２では、検体のFASTQ形式のデータをトリミングして分析に用いる。検体のデータをすべて使用しても信頼性に欠けるため、一定のクオリティ値を持つデータのみを抽出して信頼性を向上させる処理である（ＳＴ２）。

トリミング処理部１２では、ユーザが設定した塩基配列幅の平均のクオリティ値（ＱＶ）が指定の値以下になったところまでを切り取る。まず、検体データのリードの先頭付近の塩基はクオリティ値が低いため先頭から１塩基ずつＱＶをチェックし行き、所定の値以上のＱＶが出現した時点から設定された塩基配列幅の平均ＱＶを算出する。このとき、先頭の塩基から所定の値以上のＱＶが見つかった塩基の１つ手前の塩基までを削除する。そして、設定した塩基配列幅を１塩基ずつずらし、設定された塩基配列幅の平均ＱＶを算出する。平均ＱＶが所定のレベルを下回ったときに、その時点で塩基配列幅の最後尾の塩基から以降の塩基配列を削除する。

このようにして、トリミングされたリードの塩基長が一定の値を超えていれば、そのトリミング後の塩基配列を残し、一定の値以下であれば、トリミング後の塩基配列は削除され、使用されない。図８（ｂ）は、リード塩基に対して、ASCIIコードで表示されたフレッドスコアから１０進のクオリティ値に変換したものを示している。また、一例として、平均ＱＶを検査するための塩基配列幅は２０、トリミング後のリードの塩基長の最小閾値は１２０、クオリティ値のカットオフ値は７に設定されているものとする。

先頭部分の塩基のクオリティ値は、７より小さいので、この部分は削除する。先頭の塩基から８番目の塩基は、７以上になったので、８番目の塩基（図で始点と記載されている箇所）から順に取り出せる範囲を決定する。８番目の塩基から２０塩基分の平均ＱＶを算出する。この平均ＱＶは、単純平均なので、図のＷで示された、８番目の塩基から２７番目の塩基までについて、クオリティ値を加算して２０で割れば算出される。

そして、平均ＱＶを算出するための始点の塩基を１塩基ずつずらして、その塩基から２０塩基分の平均ＱＶを逐次算出する。この処理を平均ＱＶが７より小さくなるまで続ける。平均ＱＶが７より小さくなった時点で、この処理を中止し、このときの最後尾の塩基から後の塩基配列はすべて削除する。このようにして、トリミングされて残ったリードの塩基長が１２０を超えていれば、トリミング後の塩基配列を残し、１２０以下であれば、トリミング後の塩基配列は削除する。

トリミング処理部１２では、上記のクオリティ及び塩基長の判定処理（トリミング処理）を検体のすべてのFASTQデータに対して行い（ＳＴ２）、抽出された塩基配列をFASTA形式でマッピング処理部２に出力する。

次に、マッピング処理部２における図１２のＳＴ３の処理では、２種類のリファレンス配列データと２つの処理ルートの組み合わせにより処理が行われる。２種類のリファレンス配列データは、図１２に示される分節毎リファレンス配列データ記憶部１３Ａ及び１３Ｂに記憶されているデータである。第１の処理ルートは、図１２のＰ１で示される処理ルートであり、第２の処理ルートはＰ２で示される処理ルートである。インフルエンザウイルス解析装置１００では、分節毎リファレンス配列データ記憶部１３Ａ又は１３Ｂのどちらのデータを用いるかを選択できるようになっており、この選択により、マッピング処理部２は、分節毎リファレンス配列データの読み出しルートをＲ１とするか、あるいはＲ２とするかを切り替える。また、Ｐ１とＰ２の２つの処理ルートは、インフルエンザウイルス解析装置１００において、どちらかを選択できるように構成されている。後述するが、Ｐ１の処理ルートでは、１種類のみのインフルエンザウイルス型を分析出力することができ、Ｐ２の処理ルートでは、最大３種類までのインフルエンザウイルス型をすべて分析出力することができる。このように、処理方法としては、読み出しルートＲ１対して処理ルートＰ１又はＰ２のいずれかを組み合わせることができ、同様に、読み出しルートＲ２に対して処理ルートＰ１又はＰ２のいずれかを組み合わせることができるので、合計４種類の処理方法を選択することができる。Ｐ１の処理ルートは、分節毎リファレンス配列データ記憶部１３Ａ又は１３Ｂから読みだされた８分節に分割されたリファレンス配列データに、トリミング処理部１２でトリミング処理された検体データをマッピングする（ＳＴ３）。マッピングは、各分節毎のリファレンス配列データに対して行われる。マッピングには一例としてフリーツールであるBowtie2を使用した。これにより、各分節のリファレンス配列と検体データとが一致する部分が探索され、検体データが、各分節のリファレンス配列の一致した領域に貼り付けられる。一方、第２のルートは、図１２のＰ２で示される処理ルートである。これは、分節毎リファレンス配列データ記憶部１３Ａ又は１３Ｂから読みだされた８分節に分割されたリファレンス配列データに、トリミング処理部１２でトリミング処理された検体データをＰ１の処理と同様なフリーツールを用いてマッピングする（ＳＴ３）。

各分節のリファレンス配列に対してすべての検体データのマッピングが行われた後、そのマッピング結果は、良く使用されるＳＡＭ形式等で出力される。次に、カバー率算出部３では、マッピング結果に対してカバー率を算出する。カバー率算出は、以下のように行われる。
カバー率＝マッピング済領域の塩基長／各分節毎のリファレンス配列の塩基長・・（１）

カバー率算出部３では、マッピング済領域を算出するため、例えば、ＳＡＭ形式ファイルをＢＡＭ形式ファイルに変換し、これをソート済みBAM形式ファイルとし、次にPileup形式ファイルに変換する。Pileup形式ファイルは、以下の項目が含まれている。マッピングされたリファレンス位置、マッピングされたリファレンス側塩基、リファレンス側塩基にマッピングされた塩基数、マッピングされた検体側塩基等である。

（１）特定の分節において、１つのリファレンス配列の総塩基数を、マッピングを実行した際のFASTAファイルから取得する。（２）Pileupファイル中に記録されているリファレンス側塩基にマッピングされた塩基数を取得する。（３）リファレンス配列に対するカバー率を上記（１）式により算出する。（４）上記（１）〜（３）を、pileupファイル中の特定の分節のリファレンス配列すべてに対して行う。Ｐ１の処理ルートでは、（５）上記（４）の結果で最大のカバー率を有するリファレンス配列を取得する。（６）上記（１）〜（５）の処理を各分節毎に行う。これにより、一番カバー率の高いリファレンス配列が分節毎に抽出される（ＳＴ４）。したがって、合計８個のリファレンス配列が抽出される。
一方、Ｐ２の処理ルートでは、（７）上記（４）の結果で、カバー率の高い順に並べた場合、１位〜３位（上位３位）までのカバー率を有するリファレンス配列を取得する。（８）上記（１）〜（４）、（７）の処理を各分節毎に行う。これにより、上位３位までのリファレンス配列が分節毎に抽出される（ＳＴ４１）。したがって、最大で、２４個（８分節×カバー率上位３個）のリファレンス配列が抽出される。

次に、Ｐ１の処理ルートでは、コンセンサス配列作成部４で、最大カバー率を有するリファレンス配列を基に分節毎にコンセンサス配列を作成する（ＳＴ５）。マッピング済みファイルより作成されたpileupファイルをもとに、最大カバー率を有するリファレンス配列とマッピングパターンの違いを検出（変異検出）し、変異を反映したコンセンサス配列を作成する。この場合、最大カバー率を有するリファレンス配列とマッピングパターンが相違する塩基位置においては、対応するリード側の塩基の出現頻度がオペレーターの指定した最小頻度と最大頻度の間に収まればIUPAC表記に従いコンセンサス配列へ反映し、最大頻度以上であればリード側の塩基がコンセンサス配列に反映される。最小頻度未満であればリファレンス配列側の塩基がコンセンサス配列に反映される。このようにして、コンセンサス配列も、各分節毎に作成され、合計８個が形成される。

一方、Ｐ２の処理ルートでは、コンセンサス配列作成部４で、上位３位までのカバー率を有するリファレンス配列を基に分節毎にコンセンサス配列を作成する（ＳＴ５１）。マッピング済みファイルより作成されたpileupファイルをもとに、カバー率の高い上位3件のリファレンス配列とマッピングパターンの違いを検出（変異検出）し、変異を反映したコンセンサス配列を作成する。この場合、カバー率の高い上位3件のリファレンス配列とマッピングパターンが相違する塩基位置においては、対応するリード側の塩基の出現頻度がオペレーターの指定した最小頻度と最大頻度の間に収まればIUPAC表記に従いコンセンサス配列へ反映し、最大頻度以上であればリード側の塩基がコンセンサス配列に反映される。最小頻度未満であればリファレンス配列側の塩基を反映させるようにする。このようにして、コンセンサス配列も、各分節毎において各リファレンス配列毎に作成され、最大で２４個が形成される。

次に、以下のＳＴ６〜ＳＴ１０の処理については、Ｐ１の処理ルートについてもＰ２の処理ルートについても同様である。作成されたコンセンサス配列は、FASTA形式で相同性検索処理部５に出力される。相同性検索処理部５では、コンセンサス配列をクエリーとして、分節毎リファレンス配列データ記憶部１３Ａから読みだされた対応する各８分節のリファレンス配列に対し相同性検索を実行する（ＳＴ６）。相同性検索は、例えば、コンセンサス配列がPB1分節の最大カバー率を有するリファレンス配列又はカバー率上位３位までのリファレンス配列から作成されている場合には、PB1分節に対応したコンセンサス配列であるため、PB1分節のリファレンス配列のみを用いる。

同様に、PB2分節に対応したコンセンサス配列には、PB2分節のリファレンス配列を、PA分節に対応したコンセンサス配列には、PA分節のリファレンス配列のみを、HA分節に対応したコンセンサス配列には、HA分節のリファレンス配列のみを、HA分節に対応したコンセンサス配列には、HA分節のリファレンス配列のみを、NP分節に対応したコンセンサス配列には、NP分節のリファレンス配列のみを、NA分節に対応したコンセンサス配列には、NA分節のリファレンス配列のみを、MP分節に対応したコンセンサス配列には、MP分節のリファレンス配列のみを用いる。

このように、コンセンサス配列を作成する基になったリファレンス配列が属する分節と同じ分節に分類されているリファレンス配列のみと当該コンセンサス配列とで相同性検索を実行する。

相同性検索には、一例としてNCBI-BLASTNを利用した。BLAST結果はTAB区切りフォーマットの特定出力項目を指定および閾値を設定する。これにより、所定の類似度以上のリファレンス配列を表示することができる。BLAST結果より、所定の類似度以上（相同性が高く）であり、最大塩基長を有するリファレンス配列をFASTA形式にて取得し、改定リファレンス配列とする。これを分節毎に行う（ＳＴ７）。これにより、最大カバー率を有するリファレンス配列を用いた場合には、合計８個の改定リファレンス配列が作成され、カバー率上位３位までのリファレンス配列を用いた場合には最大で２４個の改定リファレンス配列が作成される。

図９は、NA分節におけるBLASTによる相同性検索の結果例を示す。図９では、各配列の１行目に記載された100.00、99.79等の数値が類似度（相同性の高さ）を表し、次に記載されている、1425の数値が塩基長を示している。

図９の例では、コンセンサス配列と類似度の高い順に上位５つのNA分節におけるリファレンス配列が示される。この場合は、一番最初の配列が改定リファレンス配列として決定される。なぜなら、この場合５つのリファレンス配列の塩基長が同じなので、類似度が最も高い唯一のリファレンス配列であるこの配列が選ばれるからである。仮に、５つのリファレンス配列の塩基長が同一ではない場合には、塩基長が最長のリファレンス配列が（たとえ類似度が５つのうちで最も高くなくとも）改定リファレンス配列として決定される。仮に、５つのリファレンス配列の塩基長が同一ではなく、かつ、塩基長が最長のリファレンス配列が複数ある場合には、その複数のうち類似度が最も高いリファレンス配列が改定リファレンス配列として決定される。改定リファレンス配列は、分節毎に決定されるので、合計８個作成される。

次に、上記で得られた８個又は最大２４個の改定リファレンス配列に対して検体データ入力部１で得られた検体データをマッピング処理部２で再度マッピング処理を行う（ＳＴ８）。そして、改定リファレンス配列のマッピングパターンから２回目のコンセンサス配列を作成する。２回目のコンセンサス配列（改定コンセンサス配列）をコンセンサス配列作成部４で作成する（ＳＴ９）。

上記のように、２度のマッピングを行うことにより、検体のインフルエンザウイルスゲノム８分節分の全長配列を決定することができる。１度のマッピングでは、リファレンス配列の数が多いため、リードがまばらに貼り付けられる可能性が高く、かつ、ある領域に対してリードが貼り付けられる数が必然的に減ってしまうため、全長配列は得られない。
そこで、本実施形態では、２度目のマッピングの際に使用されるリファレンス配列は、１度目のマッピングで得られたコンセンサス配列に対して相同性が高いリファレンス配列を分節毎に１つに決定しているので、その改定リファレンス配列は最小限の数となり、上記のような問題はなくなり、全長配列を得ることができる。また、改定リファレンス配列を相同性の高さだけで決定すると、塩基配列長が短いリファレンス配列の方が、相同性は高く算出されてしまい、インフルエンザウイルスゲノムの全長をカバーできない可能性がある。本実施形態では、相同性が高く、最大塩基長を有するリファレンス配列を改定リファレンス配列としており、このように、塩基配列の長さについても条件を設けているので、全長をカバーした改定リファレンス配列となり、全長配列を決定できるという効果がある。

さて、ＳＴ９では、VarScanというフリーツールを使用した。この過程はリファレンス配列とマッピングパターンの違い(変異)を反映したデータを出力させる機能を持った何れかのプログラムを用いて実施する。一例としてVarScanに改定リファレンス配列に貼りついたリードのパターンを読み込ませ、その改定リファレンス配列とマッピングパターンの違い(変異)を反映したデータを出力させる。このVarScanから出力されたデータに基づき、ユーザが指定したアレル頻度（ある場所における変異を支持するリードの数と支持しないリードの数の比率）の閾値を超えた変異についてSNV (Single Nucleotide Variant)、Insertion、Deletion 塩基を反映した改定コンセンサス配列を作成する。またマッピングに利用した改定リファレンス配列より3’および5’側が短い場合は塩基「N」を補完する。上記のように、IUPAC形式塩基へ変換され、分節毎に改定コンセンサス配列が作成される。最大カバー率を有するリファレンス配列を用いた場合には、改定コンセンサス配列は８個作成され、カバー率上位３位までのリファレンス配列を用いた場合には最大で２４個の改定コンセンサス配列が作成される。

上記各分節に対応した改定コンセンサス配列をFASTA形式にてファイルへ出力する。次に、登録用データ作成部７でFASTA形式からGISAID登録形式に変換される。

８分節における改定コンセンサス配列の作成はＡ型インフルエンザに由来するリファレンス配列とＢ型インフルエンザに由来するリファレンス配列それぞれについて同じ方法で別々に実施され、改定コンセンサス配列が作成されたインフルエンザ型がＡ型であったかＢ型であったかが判別される。また、カバー率上位３件までのリファレンス配列から作成された改定コンセンサス配列の場合は、改定コンセンサス配列毎に実施される。

改定コンセンサス配列を作成する基となった改定リファレンス配列が属している分節がHA又はNAの場合は、相同性検索処理部５で、改定コンセンサス配列をクエリーにして、分節毎リファレンス配列データ記憶部１３Ａに保持されている分節毎のリファレンス配列データとの相同性検索を実行し、型判定部６でインフルエンザウイルスの型判定を行う（ＳＴ１０）。相同性検索は、改定コンセンサス配列を作成する基となった改定リファレンス配列が属している分節がHAの場合は、HA分節に属するリファレンス配列のみを用いて行う。コンセンサス配列を作成する基となった改定リファレンス配列が属している分節がNAの場合は、NA分節に属するリファレンス配列のみを用いて行う。相同性検索には、NCBI-BLASTNを利用した。

例えば、インフルエンザウイルスの型判定は以下のように行われる。最終ウイルスの型判定を行うため、HA又はNAの分節に対応している改定コンセンサス配列をクエリーに、BLASTNを実行した。図１０は、HA分節に対応した改定コンセンサス配列とHA分節に対応したリファレンス配列との相同性検索を行った結果例を示す。本実施例では、相同性（類似度）が高い順に５位までを表示させるようにしている。

図１１（ａ）は、図１０の正規表現適応部分を示す。すなわち、図１０の表示で、A virus〜complete cdsまでの項目を示している。ここで、注目されるのは、図１１において、H5N8等と表現されている箇所である。図１１（ａ）は、その部分を実線の枠で囲んで示している。実線の枠で囲まれた部分を用いて型判定が行われる。

図１１（ｂ）は、図１１（ａ）の実線の枠で囲まれた部分を取り出して表示したものである。ＨＡ分節に対応している改定コンセンサス配列を用いている場合は、ＨＡ遺伝子であることはわかるが、亜型はわからないので、その亜型を判定する。図の破線の枠で囲まれた部分を見ると、Ｈの文字に続いて記載されている数値の部分が全て同値の５であるため、当該検体データはＨ５亜型と決定される。

このように、亜型判定に際し、ＨＡ亜型の判定においては、相同性検索結果において、正規表現のＨに続く数値がすべて一致した場合にＨＡ亜型の判定が当該数値により行われる。同様に、ＮＡ亜型を判定するＮＡ分節についても、相同性検索結果において、正規表現のＮに続く数値がすべて一致した場合にＮＡ亜型の判定が当該数値により行われる。

図１１（ｃ）は、ＨＡ分節の相同性検索を行った結果の正規表現適応部分を表示しているが、図１１（ａ）とは異なる例を示す。この例では、図１１（ａ）で実線の枠で囲った領域に相当する部分を比較すると、１つの結果が異なる。Ｈの文字に続いて記載されている数値の部分において上から３番目の数値は３で、他はすべて５であるので、Ｈ５とＨ３が混在しているため、当該サンプルに対して亜型判別は無効となる。すなわち、亜型判定は行われない。また、ＨＡ分節は、ＨＡ亜型の判定に用いられるものなので、Ｎに続く数値の部分は全て判定対象とならない。例えば、上から５番目だけＮ１で、後の４つの結果はＮ８となっているが、これはＨＡ亜型の判定に影響を与えない。

また、上記の型判定の結果は、登録用データ作成部７でGISAID登録形式に変換される。

以上のように、最終的な結果として出力される全長配列により、Ａ型インフルエンザかＢ型インフルエンザかの判定が可能で、Ａ型インフルエンザについてはHA、NAの亜型を迅速にかつ正確に判定できる。 HA遺伝子の配列の解明により、H5またはH7亜型HAタンパク質の開裂部位のアミノ酸配列や、HA, NA及び他の6つの遺伝子分節がコードするアミノ酸配列の特徴や、ノンコーディングリージョンでの塩基置換も同定することができる。また、カバー率上位３位までのリファレンス配列を用いた場合は、上述した処理を行うことにより、３種類までの異なった型のインフルエンザが混入していた場合も判定できる。

すなわち、上記処理ルートＰ１の方法では、最大カバー率のリファレンス配列を用いているため、体内でもっとも多く存在しているインフルエンザウイルス株のみが検出され、ウイルス量が少ないもう一方のインフルエンザウイルス株が検出されない。検出されないウイルス株が毒性の強い亜型であった場合、感染力の強いウイルスの感染を見逃すことになる。処理ルートＰ２の方法では、カバー率上位３位までのリファレンス配列を用いるため、複数のインフルエンザウイルス株がひとつの個体に重複感染している場合でも量が少ないほうのウイルスを見逃すことがなくなる。また、特に読み出しルートＲ２と処理ルートＰ２を組み合わせることにより、マッピングに用いるリファレンス配列相互の類似度が一定のレベル以下に調整されることで、存在量が１位のインフルエンザウイルスと２位のインフルエンザウイルスを見分けることが可能になる。類似度に上限がない場合、例えばＡ型インフルエンザのHAで一番存在量が多かった亜型がHA5で２番目がHA3だとすると、HA5のリファレンス配列間の類似度が高いためマッピング時に複数のHA5リファレンス配列にリードが張り付き、結果的に存在量1位から3位まですべてHA5に由来するリファレンス配列を基にした型の判定とコンセンサス配列の作成が行われてしまい、HA3が存在量の下位に沈んでしまうため検出できない。亜型内のリファレンス配列の類似度に制限を設けることで、HA5ウイルスに由来するリードはリファレンス配列のHA5の中でもっとも相同性の高いリファレンス配列にマッピングされ、ほかのHA5リファレンス配列に対しては相同性が、例えば８０％以下と低くなるためマッピングされない。この結果、検体中に僅かに含まれていたHA3亜型ウイルス由来のリードがHA3にマッピングされたことが上位3位までの検出結果に含まれることになり、３種類までの異なった型のインフルエンザウイルスを判定できることに繋がる。

なお、インフルエンザウイルス解析装置としてコンピュータを機能させるためのインフルエンザウイルス解析プログラムは、半導体メモリ、磁気ディスク、光ディスク、光磁気ディスク、磁気テープなどのコンピュータ読み取り可能な記録媒体に記録でき、また、インターネットなどの通信網を介して伝送させて、広く流通させることができる。

１ 検体データ入力部
２ マッピング処理部
３ カバー率算出部
４ コンセンサス配列作成部
５ 相同性検索処理部
６ 型判定部
７ 登録用データ作成部
１０ 制御部
１１ リファレンスデータ入力部
１２ トリミング処理部
１３ 分節分割部
１３Ａ 分節毎リファレンス配列データ記憶部
１３Ｂ 分節毎リファレンス配列データ記憶部
２１ シーケンサ
２２ 公開データベース
２３ 表示端末部
１００ インフルエンザウイルス解析装置
１３１ 相同性分析部

Claims (6)

1. 参照用インフルエンザウイルスデータが入力されるリファレンスデータ入力部と、
   前記リファレンスデータ入力部に入力された参照用インフルエンザウイルスデータの塩基配列を複数の分節に分類して分節毎にリファレンス配列を集積する分節分割部と、
   検体の塩基配列が入力される検体データ入力部と、
   前記分節毎に集積されたリファレンス配列に前記検体の塩基配列をマッピングするマッ
   ピング処理部と、
   前記マッピング処理部によるマッピング結果からリファレンス配列の塩基長に対するマ
   ッピング済領域の塩基長の割合であるカバー率を算出するカバー率算出部と、
   前記算出されたカバー率の最大の値又は上位３位までの高カバー率を有するリファレンス配列のマッピング結果に基づいて変異を反映したコンセンサス配列を分節毎に作成するコンセンサス配列作成部と、
   分節毎に前記コンセンサス配列と前記集積されたリファレンス配列の相同性検索を実行し所定の類似度以上のリファレンス配列のうち最大の塩基長を有する１つの改定リファレンス配列を決定する相同性検索処理部とを備え、
   前記マッピング処理部において改定リファレンス配列に前記検体の塩基配列をマッピングしてマッピング結果に基づいて変異を反映した改定コンセンサス配列を作成することを特徴とするインフルエンザウイルス解析装置。
2. 前記分節分割部では、前記参照用インフルエンザウイルスデータの塩基配列をヘッダーに所定のキーワードが含まれるか否かにより、８種類の分節といずれの分節にも分類できない塩基配列とに分類することを特徴とする請求項１に記載のインフルエンザウイルス解析装置。
3. 前記分節分割部に集積されたリファレンス配列間の相同性分析を分節毎に行い、一定の類似度以下のリファレンス配列に調整する相同性分析部を備え、前記マッピング処理部によりリファレンス配列に前記検体の塩基配列をマッピングする際には、前記一定の類似度以下に調整されたリファレンス配列を用いることを特徴とする請求項１又は２に記載のインフルエンザウイルス解析装置。
4. 参照用インフルエンザウイルスデータの塩基配列を複数の分節に分類して分節毎にリファレンス配列を集積する分節分割工程と、
   前記分節毎に集積されたリファレンス配列に検体の塩基配列をマッピングする第１のマッピング処理工程と、
   前記第１のマッピング処理工程によるマッピング結果からリファレンス配列の塩基長に対するマッピング済領域の塩基長の割合であるカバー率を算出するカバー率算出工程と、
   前記算出されたカバー率の最大の値又は上位３位までの高カバー率を有するリファレンス配列を分節毎に抽出する抽出工程と、
   前記最大のカバー率又は上位３位までの高カバー率を有するリファレンス配列のマッピング結果に基づいて変異を反映したコンセンサス配列を分節毎に作成するコンセンサス配列作成工程と、
   分節毎に前記コンセンサス配列と前記集積されたリファレンス配列の相同性検索を実行し所定の類似度以上のリファレンス配列のうち最大の塩基長を有する１つの改定リファレンス配列を決定する相同性検索処理工程と、
   前記改定リファレンス配列に前記検体の塩基配列をマッピングしてマッピング結果に基づいて変異を反映した改定コンセンサス配列を作成する第２のマッピング処理工程とを備えたことを特徴とするインフルエンザウイルス解析方法。
5. 前記分節分割工程で集積されたリファレンス配列間の相同性分析を分節毎に行い、一定の類似度以下のリファレンス配列に調整する相同性分析工程を備え、前記第１のマッピング処理工程では、前記一定の類似度以下に調整されたリファレンス配列を用いることを特徴とする請求項４に記載のインフルエンザウイルス解析方法。
6. 請求項１に記載のインフルエンザウイルス解析装置としてコンピュータを機能させるた
   めのインフルエンザウイルス解析プログラム。

Images