INSDC Feature Table 地獄めぐり

バイオインフォマティシャンたるもの一度は扱ったことがあるであろうファイルの一つに GenBank と EMBL 形式ファイルがあるかと思います。その中には配列の様々な特徴 (Feature) を記した Feature Table という部分があります。この部分のフォーマットは International Nucleotide Sequence Database Collaboration (INSDC) という NCBI GenBank、EMBL ENA、NIG DDBJ の三者によって成り立っている国際組織によって定義されています (The DDBJ/ENA/GenBank Feature Table Definition)。さて、最近とあるソフトウェアを開発するにあたってこの Feature Table のパーサを書くことになったんですが、あまりにも地獄すぎたのでこの機会に大放出します。

1. Location 地獄

それぞれの Feature には配列に対応する位置情報、INSDCのドキュメントで言うところの Location が与えられています。多くの場合 Location は start..end という形でレンジで与えられるのですが、ドキュメントを見れば分かる通りいくつかパターンがあります。数字1つで位置塩基（残基）を指していたり <start..>end というように Feature の全長がカバーされていないという注釈が入っていたりします。このくらいまでは多少なりともプログラミングができればある程度処理できるでしょう。

さて、問題はここからです。配列の 3' 側で配列に欠けがある場合は start..>end と記述することが仕様で要求されているわけですが、一部古いエントリでは start..end> となっている場合があるらしいです（ソース発見できず）。当然そういう例は簡単には発見できないんですがこれを許容するデータを突っ込んでくることは当然想定しておく必要が（古いデータ・ソフトウェアを使うのがごく普通にありえる世界なので）あるわけです。この時点で大抵の人は血反吐を吐きながら例外処理を加えることになります。

しかし Location 地獄の本当の恐ろしさに比べればこの程度はジャブみたいなものです。そう、Locationの中には関数のよう記述できる Locationも存在するのです。このようなLocationは3つ存在します。それが complement、join、orderです。

まず complement はそのまんまで与えられた領域の相補鎖を指していて complement(location) と書きます。次に join と order は前者が与えられた複数領域を結合するのに対して後者はただ複数領域を指すだけです。これらはそれぞれ join(location1, location2, ...)、order(location1, location2, ...) という具合で記述します。

さて、聡明な方はもうお気づきかと思いますがこれらに渡すLocationは当然前述のような単純なLocationでもいい上に関数型Locationでもいいわけです。するとネストされたLocationをパースするということが必要になってきます。一応規約には

Note : location operator "complement" can be used in combination with either "join" or "order" within the same location; combinations of "join" and "order" within the same location (nested operators) are illegal.

と書かれていますが、今後の例外なども踏まえるとこれを守らないケースがあると思っておいたほうがいいということは容易に想像できます。このネストを真面目に解決しようとすると常人ならば死にます。パーサコンビネータなど使おうねとなってくるので並大抵のバイオインフォマティシャンであればここで心が折れるでしょう。

2. Qualifier Value 地獄

さぁなんとか Location を正しくハンドリングできたぞ！ともなれば次は Qualifier のパースです。Qualifier なくして Feature なし。Qualifier は名前と値のセットのことを言い、Feature の様々な属性を記述するために使われます。Qualifier の書式は一見とてもシンプルで、以下の3通りがあります。

1. /name="value"
2. /name=value
3. /name

「なんだ！簡単じゃないか！」そう思ったそこのあなた。簡単だったら地獄なわけがないじゃないですか。もちろんここにも地獄が潜んでいます。

まず 1. のケースについてですがここでは便宜的に「文字列値」と呼ぶことにします。文字列値の最大の地獄はこれが複数行にまたいで存在しうるということです。バイオインフォマティクスに限らず多くのインフォマティクスの領域でファイル内の行の幅が80文字を超えないようにするというのはよくある習慣ですね。これがここにも適用されます。そのため、「改行した後の空白の数は直前の行の空白の数と一致しなければフォーマット不備なのか？」という疑問が自然と浮かび上がるわけですが、規約にはそれについて何も書かれていません。どうしろと？しかもGenBankフォーマットの場合（EMBLでもあるかもしれませんが）基本的にほとんどのフィールドが複数行をまたぐ場合に改行とインデントを取り除いても許されるんですが E. coli K-12株 (NC_000913) の COMMENT フィールドはなぜかインデントの深さがちょっと深い行が二行あり、この二行は列挙なのでインデントの深さをを無視して改行を取り除くと明らかに文章がおかしくなるというトラップが潜んでいるのでシンプルに「まぁええやろｗ」としてしまうわけには行きません。

さて、腹を決めていざパーサを書くかとなったところでプログラマであればふと疑問が浮かんでくるでしょう。テキスト内に出てくるダブルクォートはどうするんだ？多くの場合、テキスト内にダブルクォートが登場する場合はエスケープをします（\" と記述する）。もちろんここでもエスケープをするのですが、そのエスケープの記述方法は規約でダブルクォート2つ ("") で記述すると明記してあります。いや、なんでだよ。意味がわからないよ。なんで色気を出して独自エスケープ記法を編み出してるんだよ。キレた私はそっとこの仕様を見なかったことにしました（今書いてるソフトウェアではそのうち対応します......）。というかこれ、規約読まずにパーサ書いた人絶対これに起因するバグ仕込んでるでしょ。世の中に存在するパーサの何割にこの致命的な地雷が埋まっているのか考えたくもないですね。おーやだやだ。

ここまでの流れを汲んで文字列値は「ダブルクォートの間のインデントだけ取り除いた文字列として扱う」という取り決めをしてしまえばだいたい解決します。ところがどっこい translation という Qualifier が存在してコイツが悪さをします。この Qualifier はつまるところ CDS が翻訳された後のタンパク質の配列を指しているわけですが当然これも改行されているので配列として扱いたければ改行を落とす必要があります。勘弁してくれ。

残りの2つのケースについてですが、 2. のケースを「リテラル値」、3. のケースを「トグル値」とそれぞれここでは呼びますがこの2つに関しては文字列値のような改行の不安がないため意外とすんなりとパーサは書けます。そう、「パーサは」ね......

3. Qualifier Name 地獄

Qualifier Valueをうまくパースできるようになった皆さん、おめでとうございます。では次なる地獄は Qualifier Name です。「え、こっちが先じゃないの？」と思われる方もいらっしゃるかもしれませんが、最初にパーサを書くときはnameはそれとしてパースしてValueが3種類のどれかを判断してよしなに処理するというやり方になるかと思います。まぁ確かにこれでもいいんですがフォーマッタを書くときにどのQualifierがどういう値を持っていたかってわからないんですよね。聡明な読者の方は「リテラル値は基本的に数値なのでは？それなら数値かどうかさえわかれば難しい話ではないだろう」と思われるかもしれません。が、そんなわけないじゃん。数値でも unknown という値が許されていたりします。やめろ。というわけでフォーマットするときにどの Qualifier Name がどの Qualifier Value の書式を取るのかを判断する必要があります。「そんなの、規約の中で定義されてる Qualifier Name を元に分岐すればいいじゃん」と言われればそうなんですが、めんどくせぇぞ？なにせ Qualifier は100種類もありますからね。がんばってくださいね。

あ、あと言い忘れてたんですが、Feature Table Definitionには定義されてない Qualifier Name が普通にGenBankのレコードの中に登場したりするぞ！（例：ヒトのインシュリン には calculated_mol_wt や site_type などが出てくる）がんばってね！！！

4. 重複地獄

晴れてパーサが書けたぞ！となってテストを走らせると入力と出力が一致しない...... なんでだ？と思い一致していない場所を見てみると元のファイルでは2つあった同名の Qualifier が出力では1つになっていたというケースがありえます。例えば db_xref フィールドは1つの Feature に複数付与されていることがあります。そう、Qualifier Listの名の通り Qualifier の集合はディクショナリではなく Key-Value Pair のリストなのです！アホか。

何が面倒かと言うと「論理的には一意な Qualifier と差別化する必要があるか？」という問いに答える必要があります。規格はこれについて何も教えてはくれません。世のプログラマの良識に判断が委ねられているわけですね。

5. Mandatory/Optional 地獄

ここまでくればあとは Feature Key + Location 行と Qualifier List のパーサを組み合わせておしまい！と言いたいところですがそれは叶いません。なぜなら厳密に規格に沿ったパーサを書く場合、Feature Key ごとに許可された Qualifier が制限されているため、このバリデーションを加える必要があります。Feature Keyは 53種類あるので、それぞれについて100種類の Qualifierのどれがついていていいかを判断する必要があります。普通に考えれば許可された Qualifier Name 以外は弾くようにすればいいんですが、上述の通り規格に記されていない Qualifier Name がついていることが想定されるため、WhitelistではなくBlacklistにする必要があります。もうやだ。

6. フォーマッタ地獄

様々な妥協と血と涙とアルコールと罵倒を経てようやく完成したパーサ。あなたの多大なる時間的コストはしかし完璧なFeature Tableパーサという形で実を結びました。全ては報われ例外処理だらけのスパゲティコードは成仏し涅槃への道が開けたと思ったか。最後にして究極の地獄、フォーマッタ地獄があなたを待ち受けています。今までよりひどい地獄があるのか！？そう思われるかもしれませんね。あるんだなぁ、これが。

読者の多くにとってはフォーマッタがそんなに地獄であるように感じないかもしれません。私もそうでしたが、一番の地獄の始まりはここからです。なぜなら、これまでに出てきた地獄全てに触れている最大にして最悪の地獄がフォーマッタ地獄だからです。

まず Qualifier Value 地獄。パーサを書いたときに改行を落としましたか？ではまたこの改行をフォーマットする際に戻さなくてはなりませんね。このときにたいがいのプログラマは思うでしょう。「80文字より手前の位置にある空白を改行で置換すればいい、ない場合はより後ろにある空白を探して最初の空白を改行で置換しよう」と。Feature Table Definitionではこういった安易な判断をする愚かで矮小なプログラマを危惧してきちんとその例外となる Qualifier の例を掲載しています。その1つがこちらです：

Qualifier       /experiment=
Definition      a brief description of the nature of the experimental 
                evidence that supports the feature identification or assignment.
Value format    "[CATEGORY:]text"
                where CATEGORY is one of the following:
                "COORDINATES" support for the annotated coordinates
                "DESCRIPTION" support for a broad concept of function such as that
                based on phenotype, genetic approach, biochemical function, pathway
                information, etc.
                "EXISTENCE" support for the known or inferred existence of the product
                where text is free text (see examples)
Example         /experiment="5' RACE"
                /experiment="Northern blot [DOI: 12.3456/FT.789.1.234-567.2010]"
                /experiment="heterologous expression system of Xenopus laevis
                oocytes [PMID: 12345678, 10101010, 987654]"
                /experiment="COORDINATES: 5' and 3' RACE"
Comment         detailed experimental details should not be included, and would
                normally be found in the cited publications; PMID, DOI and any 
                experimental database ID is allowed to be used in /experiment
                qualifier; value "experimental evidence, no additional details
                recorded" was used to  replace instances of /evidence=EXPERIMENTAL in
                December 2005

ポイントはこの中の Example の部分です。「これのどこが変なんだ？」と思われるかもしれませんが、例えば

                /experiment="Northern blot [DOI: 12.3456/FT.789.1.234-567.2010]"

これは DOI: の後の空白で本来改行するべき長さですが例では改行されていません。「いやいや、ただここで改行入ってないだけだからテストケースとしては改行を入れればいいじゃん」と思いました？

                /experiment="heterologous expression system of Xenopus laevis
                oocytes [PMID: 12345678, 10101010, 987654]"

こっちはばっちり改行が入ってるんですよね。しかもこれ、本来は Xenopus の直後で改行しているべきなんですがおそらく気を利かせた筆者が「80文字超えちゃうけど生物種名は一行になってて欲しいよね〜」とあえてこのようにしたんだろうと思います。ふざけろ。生物種名を判定するだけでも一苦労ですが、超えた文字数が何文字以上の場合に同じ行に残すかなどの判断は完全に人間の好みです。こんなの規格じゃねえよ。

全てを悟った私は改行を残し、適切に改行を挿入するのはユーザの責任と割り切ることにしましたとさ。

まだまだ続く地獄めぐり

他にも私が妥協したことによって訪れずに済んだ地獄がいくつか存在します（Qualifier Value が取れる値の語彙が決められている予約語地獄、値を更に複雑な構造としてパースできる型地獄など）。もっと深く規格を読み込んでいけばまだまだ発見されていない地獄も数多くあるのではないかと思います。と言いつつも歴史的経緯から仕方がない部分も当然ありますし、毎年開催されているBioHackathonなどでどうにかしようと皆さんががんばっていることは当然知っています。その上で、一人でも多くのバイオインフォマティシャンがこの数々の地獄の存在を知ることで苦しみを事前に避けることができればいいなと思いこんなエントリを書いてみました。あなたが日々何気なく使っているパーサに感謝の念を払いつつ、どこにバグが潜んでいるかわからないということを常に念頭に入れておきましょう。

最後に1つアドバイス、不必要にGenBankのパーサは書くな。
（※ 私はGoにパーサがないので仕方がなく書いています）