Open athos opened 2 years ago
ありがとうございます。ただ、このエラーはこれまで生じなかったので、、例えばこのエラーを起こすような共有できるテストデータなどあれば嬉しいですがいかがでしょうか?
エラーを再現するためのサンプルデータを用意しました。
こちらのzipファイルに以下のファイルを同梱しています:
seq1.fq.gz
, seq2.fq.gz
realigned.fixed.sam
seq1.fq.gz
およびseq2.fq.gz
は以下の融合遺伝子を検出するために人工的に作ったhg38用の配列です:
chr19:42295064
(strand: +
), 3'側 chr4:190174745
(strand: +
)GGGTGGAGGCAGTG
realigned.fixed.sam
は上記seq1.fq.gz
およびseq2.fq.gz
をbwaでアライメントしたうえで、出力されたSAMをfusionfusionが期待するSTARのSAM形式に調整した結果のファイルです。手許の環境では、このSAMファイルを使って以下のようにエラーが再現することを確認しました:
$ fusionfusion --reference_genome <reference_file> --star realigned.fixed.sam --out <out_dir> --genome_id hg38 --no_blat --grc
Traceback (most recent call last):
File "/usr/local/bin/fusionfusion", line 11, in <module>
load_entry_point('fusionfusion==0.5.0', 'console_scripts', 'fusionfusion')()
File "build/bdist.linux-x86_64/egg/fusionfusion/__init__.py", line 10, in main
File "build/bdist.linux-x86_64/egg/fusionfusion/run.py", line 149, in fusionfusion_main
File "build/bdist.linux-x86_64/egg/fusionfusion/parseJunctionInfo.py", line 596, in parseJuncInfo_STAR
File "build/bdist.linux-x86_64/egg/fusionfusion/parseJunctionInfo.py", line 494, in getFusInfo_STAR
NameError: global name 'read' is not defined
$
修正を最新のdevelブランチにリベースしました。
エラーが再現したのはbwaでアライメントしたうえで、出力されたSAMをfusionfusionが期待するSTARのSAM形式に調整した結果のファイル、とのことですが、fastqを直接STARでアライメントした結果のファイルをfusionfusionに入力していただけないでしょうか。 よろしくお願いいたします。
コメントありがとうございます。
こちらの件について、https://github.com/Genomon-Project/fusionfusion/pull/17#issuecomment-1028692226 で共有致しましたサンプルデータをこちらでSTARでマッピングしてみましたが、Chimeric.out.sam
に期待するようなCICとDUX4に跨がるアライメントの組が出力されない(そのため、insertionも期待したように検出されない)ことが分かりました。
STARのオプションの調整で意図したアライメントの組が得られるようになるか、また別のfusionでよりエラーを再現しやすいものがあるか等、もう少しこちらで調べてみたいと思います。
こちらの件について、改めて60件ほどinsertionを含むfusionの配列を準備して検証してみましたが、STARによって直接アライメントすることでエラーを再現できるサンプルは見つけられませんでした。
しかし、エラーが発生する箇所はfusionfusionがinsertionを検出した場合には必ず実行される処理で、サンプルが異なることで、あるいは将来的なSTARやfusionfusionの更新によってinsertionが検出可能になった場合、必ず未定義変数エラーを引き起こすことになるため、修正しておくのが望ましいではないかと考えています。
なお、このエラーは以下のようにPythonの静的解析ツール(pylint等)でも検出することができます:
$ pip install pylint
$ pylint fusionfusion | grep -E 'E[[:digit:]]+'
...
fusionfusion/parseJunctionInfo.py:494:30: E0602: Undefined variable 'read' (undefined-variable)
fusionfusion/parseJunctionInfo.py:566:30: E0602: Undefined variable 'read' (undefined-variable)
...
$
16 の修正版のプルリクエストになります。
基本的には #16 と同様に、insertionを持つ融合遺伝子を検出した場合に発生するエラーを修正する変更です。 前回の修正では、
read.seq
という未定義変数の参照を避け、F[9]
(SAMファイルのSEQカラムに当たる)からinsertionの配列を抽出するようにしてエラーを修正しようしていました。しかし、前回の修正ではアライメントの並びによっては正しいアライメントからinsertion配列を抽出できない場合がありました。具体的には、同一QNAMEを持つプライマリ、サプリメンタリ、およびそれらのペアからなるアライメントの3つ組が、入力SAMファイル中で、
もしくは、
の順に並んでいる場合、ペアのアライメントから誤ってinsertion配列を抽出してしまう問題がありました。
このプルリクエストでは、insertion配列を必ずプライマリのアライメントから抽出するように変更することで、このアライメントの順序依存性による問題を解決します。