ようやく修論の日本語版が一通り出来上がったので、ふたたびアセンブリ沼に帰ってきました。ショットガンシーケンスからのゲノム作成はちょっと心が折れたので（参照：ぐだぐだde novoゲノムアセンブリ日記 - 研究メモ）、一旦ミトゲノムでも作ってみますか～という感じ。
参照したのは以下のサイトです。これらを見ればほぼできます。
使ったプログラムはこちらです。
本記事では動物ミトゲノムのアセンブリについて書きますが、植物のミトゲノムや葉緑体など、ほかのオルガネラのアセンブリにも使えるみたいです。
README の下のほうに何のアセンブリができるか書いてあるはずです。
GetOrganelle/README.md at master · Kinggerm/GetOrganelle · GitHub
以下に、自分がやったことを書いていきます。
諸々のインストール
導入はほぼ上坂先生のページそのままです。
１．プログラムのインストール
※新しい仮想環境にプログラムをインストールする場合です
※新しい環境にインストールしない場合は一番下の一文のみでOK
conda create -n getorganelle -y
conda activate getorganelle
conda install -c bioconda -c biojj getorganelle -y
２．アセンブリに用いるデータベースのインストール
git clone https://github.com/Kinggerm/GetOrganelle.git
ちなみにこのプログラムで動かすスクリプトは二種類
①get_organelle_config.py：データベースを読み込むためのスクリプト
②get_organelle_from_reads.py：実際にアセンブリを行うためのスクリプト
３．対象種に合わせたデータベースを読み込む
使えるデータベースはembplant_pt, embplant_mt, embplant_nr, fungus_mt, animal_mt, and/or other_pt の６種類
今回は動物なので
get_organelle_config.py -a animal_mt
４．プログラムの実行
時間がかかるのでnohup で実行したほうが良いです。
nohup get_organelle_from_reads.py -1 Hoge_R1.fastq.gz -2 Hoge_R2.fastq.gz -o Hoge -R 10 -k 14 -F animal_mt -s reference_mitchondriongenome.fasta
-1, -2 #インプットファイル（生データ）の読み込み
-o #アウトプットファイルの名前
-R #公式サイトにはMaximum extension rounds と書いてあるのですが、何だかわかりません。README で動物のミトコンゲノムでは10 に指定されているので今回は10 にしました。
-k #最適k 値
-F #データベースの設定
-s #リファレンス指定。README では指定されていないんですけど、わたしのときはリファレンス指定しないでやったらプログラムが止まっちゃいました。
止まった…
ERROR: No animal_mt seed reads found!
2021-01-13 15:48:38,235 - ERROR: Please check your raw data or change your animal_mt seed!
しかし、生成されたSeed フォルダにはちゃんとリファレンスのミトゲノムが入っていました。
どうして…
強いてリファレンスにおかしいところがあるとすれば、サンプル名がスペースで区切られていることですが、そんなのが関係あるんだろうか…
とりあえずリファレンスを治してもう一周してみます。
わたしはどうやらアセンブリの神に愛されていないらしいです。
20210115
強いてリファレンスにおかしいところがあるとすれば、サンプル名がスペースで区切られていることですが、そんなのが関係あるんだろうか…>関係なかった。
やっぱり同じエラー吐いて止まりました。
あと私が変えられる数値はR とk くらいです。
R のことはちょっとよくわからないので、k をいじってみることにしました。
今回与えた14 という数字は、自分が持っているショットガンシーケンスデータをk-mer ごとに頻度分布をつくって比較し、一番収束しそうなものを使いましたが、ダメだったので、kmergenie というソフトで算出した最適k 値を使ってみました。
さて、どうなるかな…
20210116
はい、止まった。まあたしかにseed 配列に問題があるって言われてるのにk-mer いじっても仕方ないんですわ。はぁ…。
よく考えてみたらリファレンスのファイル名がめちゃめちゃ長くて、もしかして…と思いました。とりあえず、短くして回してみます。
修論には間に合わなくても論文投稿には間に合わせたいな…ほんとはゲノムのショットガンシーケンスのほうもどうにかしたいのですけど…
【追記】
2017/01/27
結局上手くいかず、ソフトウェアの開発者のJianjun 先生にログを送ったところ、データの中にターゲットリードがない/制限されているようだ、という返信が。
ひとまず、example でも実行して様子を見てみることにしました。
スパコンとの相性が悪いのかな…。
