区間に関するクエリ

TL;DR

Genomicなデータを扱っていると、区間に対するクエリを扱いたいことがよくある。

例えば、遺伝子領域と重複するピークやSNVを探したり、Open Chromatin Regionにあるヒストンマークを探したりするような利用方法が考えられる。

こういったクエリは、ナイーブに投げると になってしまう。 現実問題として、この計算量はあまりよろしくない。

解決策

CLI

CLIであるbedtoolsなどを使えばこういったクエリをある程度効率的に処理することができる。基本的には くらいになる。速度面で問題なければこういったツールを活用するのが良いと考えられる。

• bedtools
• bedops
• bedtk

テスト済みのツールが使えるので定型的な処理をする分には安心感がある。また、bedなどのパーサーを書かなくていいので嬉しい。 一番メジャーかつできる処理が多いのはbedtools、高速なのはbedtkかなという感じを受ける。

自作

CLIで解決できない問題の場合は自分でコードを書くことになる。その場合に使えるようなデータ構造は当然ながら古くから研究されている。計算量のオーダーとしてはそこまで変化はないっぽいが、ベンチマークとか見るとだいぶ違うので最適化はされていっているらしい。

基本的なアイデアとしては、Interval TreeやR-treeなどが使われている。

このあたりの基本的なデータ構造の実装としては、色々あるがrust-bioのデータ構造だったり、python実装だったりがある。 ALV木として実装されているので、自分で実装するよりこういったものを利用すると苦労は少なそう。

最近開発されているより高速・省メモリなデータ構造としては、以下のような物がある。ctrangesはbedtkで使用されているデータ構造。

Name	Language	Github
Augmented Interval List (AIList)	C, Python	https://github.com/databio/AIList
cgranges	C, C++	https://github.com/lh3/cgranges
Cache Oblivious Interval Trees (COITree)	Rust	https://github.com/dcjones/coitrees

Githubで公開されているベンチマーク的にはCOITreeが一番早い。個人的にはCほとんどわからないので普通にCOITree使いそう。