多重遺伝子の分子進化や、複数遺伝子座を使った系統解析で、個々の遺伝子のゲノム上の位置を示した図を作りたい時がありますが、そういう時に使える簡単な方法を調べたところ、以下のサイトがありました。
https://stackoverflow.com/questions/33727432/how-to-plot-positions-along-a-chromosome-graphic
Rのbarplot関数を使う方法が簡単で、これは自分のデータでもできました。ggplot2を使う方法は試していませんが、こちらもできそうです。遺伝子名と一緒に示したい場合は、こちらが良い感じです。
ちなみに自分のデータでやってみたら、こんな感じになりました（ドジョウのゲノム上にある1,798個の遺伝子の位置を横線で表示）。
気がついたらもう梅雨になっていました。毎日何かと忙しく過ごしているうちに、日々が過ぎていきます。息子の成長と、自分が少し歳を取ったことを感じます。(2023.6.11)
伊勢田の弁天様にお参りをしています。
3月に入ってだいぶ春らしくなってきました。昔からこの時期は花粉症が酷かったのですが、数年前に舌下免疫療法を試したところ、昨年くらいから症状がかなり改善して、だいぶ過ごしやすくなりました。ありがたいことです。(2023.3.7)
BRAKER 3.0.0をインストールしたので、それぞれの遺伝子予測方法の結果を比較してみた。
GeneMark-EXおよびAUGUSTUSのtrainingに使用する配列情報により、3つの方法が提案されている。
いくつかの異なる組織由来のRNA-seqデータを使用 これらをHisat2でリファレンスゲノムにマッピングし、BAM形式で読み込ませた。
DRR328565, SRR5997795, SRR7586852, SRR7586854,SRR14160170, SRR14160176
遺伝子構造予測のtraining（Braker2, Braker3）に使用するタンパク質のアミノ酸配列は、OrthoDBのVertebrataをダウンロードして使用した（proteins.fasta）。
詳細は以下のサイトから。
braker.pl --threads 20 --species=Mang_china \ --genome=JALDXJ000000000_Mang_assembly.sm.rename.fasta \ --bam=Mang_china_RNAseq/DRR328565.bam,\ Mang_china_RNAseq/SRR5997795.bam,\ Mang_china_RNAseq/SRR7586852.bam,\ Mang_china_RNAseq/SRR7586854.bam,\ Mang_china_RNAseq/SRR14160170.bam,\ Mang_china_RNAseq/SRR14160176.bam
braker.pl --threads 20 --species=Mang_china \ --genome=JALDXJ000000000_Mang_assembly.sm.rename.fasta \ --prot_seq=proteins.fasta
braker.pl --threads 20 --species=Mang_china \ --genome=JALDXJ000000000_Mang_assembly.sm.rename.fasta \ --prot_seq=proteins.fasta \ --bam=Mang_china_RNAseq/DRR328565.bam,\ Mang_china_RNAseq/SRR5997795.bam,\ Mang_china_RNAseq/SRR7586852.bam,\ Mang_china_RNAseq/SRR7586854.bam,\ Mang_china_RNAseq/SRR14160170.bam,\ Mang_china_RNAseq/SRR14160176.bam
| algorithm | N_genes | N_transcripts |
|---|---|---|
| Braker1 | 45159 | 54210 |
| Braker2 | 48408 | 52575 |
| Braker3 | 30882 | 35050 |
Braker1, Braker2と比較して、Braker3ではgenes, transcriptsの数がかなり少なかった。単一の遺伝子を複数の異なる遺伝子として予測してしまうケースが減少しているのかもしれない。
| Complete (single + duplicated) | Complete (Single-copy) | Complete (Duplicated) | Fragmented | Missing | |
|---|---|---|---|---|---|
| Braker1 | 93.9 | 72.6 | 21.3 | 0.7 | 5.4 |
| Braker2 | 90.9 | 82.2 | 8.7 | 3.2 | 5.9 |
| Braker3 | 93.8 | 70.4 | 23.4 | 2.0 | 4.2 |
Braker1とBraker3のcompletenessの数値は大きく違わなかったが、Braker1の方がDuplicated BUSCOsの割合はやや少なかった。Braker2は、Braker1, Braker3と比較すると、completenessがやや低かったが、duplicated BUSCOsの割合はBraker1, Braker3よりもかなり低かった。また、missingの割合はBraker3がもっとも低かった。
今回の解析で気付いたことの1つだが、出力されるGTFファイルの形式が、Braker3のみ少し違っていた。Braker3では、geneおよびtranscriptを示す行がなく、各transcriptが"###"で区切られて出力される形式になっていた（以下の例）。
# Braker3 ### CM050755.1 AUGUSTUS start_codon 220137 220139 . + 0 transcript_id "g5.t1"; gene_id "g5"; supported "False"; CM050755.1 AUGUSTUS CDS 220137 224024 0.72 + 0 transcript_id "g5.t1"; gene_id "g5"; cds_type "single"; CM050755.1 AUGUSTUS stop_codon 224022 224024 . + 0 transcript_id "g5.t1"; gene_id "g5"; supported "True"; ### CM050755.1 AUGUSTUS start_codon 239666 239668 . + 0 transcript_id "g6.t1"; gene_id "g6"; supported "False"; CM050755.1 AUGUSTUS CDS 239666 239817 0.92 + 0 transcript_id "g6.t1"; gene_id "g6"; cds_type "initial"; CM050755.1 AUGUSTUS CDS 240465 240543 0.91 + 1 transcript_id "g6.t1"; gene_id "g6"; cds_type "terminal"; CM050755.1 AUGUSTUS intron 239818 240464 0.92 + 0 transcript_id "g6.t1"; gene_id "g6"; supported "True"; CM050755.1 AUGUSTUS stop_codon 240541 240543 . + 0 transcript_id "g6.t1"; gene_id "g6"; supported "False"; # Braker1 & Braker2 CM050755.1 AUGUSTUS gene 34204 47808 . + . g1 CM050755.1 AUGUSTUS transcript 34204 47808 0.99 + . g1.t1 CM050755.1 AUGUSTUS start_codon 34204 34206 . + 0 transcript_id "g1.t1"; gene_id "g1"; CM050755.1 AUGUSTUS CDS 34204 34495 1 + 0 transcript_id "g1.t1"; gene_id "g1"; CM050755.1 AUGUSTUS exon 34204 34495 . + . transcript_id "g1.t1"; gene_id "g1"; CM050755.1 AUGUSTUS intron 34496 47032 1 + . transcript_id "g1.t1"; gene_id "g1"; CM050755.1 AUGUSTUS CDS 47033 47140 1 + 2 transcript_id "g1.t1"; gene_id "g1"; CM050755.1 AUGUSTUS exon 47033 47140 . + . transcript_id "g1.t1"; gene_id "g1"; CM050755.1 AUGUSTUS intron 47141 47278 1 + . transcript_id "g1.t1"; gene_id "g1"; CM050755.1 AUGUSTUS CDS 47279 47416 1 + 2 transcript_id "g1.t1"; gene_id "g1"; CM050755.1 AUGUSTUS exon 47279 47416 . + . transcript_id "g1.t1"; gene_id "g1"; CM050755.1 AUGUSTUS intron 47417 47587 1 + . transcript_id "g1.t1"; gene_id "g1"; CM050755.1 AUGUSTUS CDS 47588 47808 0.99 + 2 transcript_id "g1.t1"; gene_id "g1"; CM050755.1 AUGUSTUS exon 47588 47808 . + . transcript_id "g1.t1"; gene_id "g1"; CM050755.1 AUGUSTUS stop_codon 47806 47808 . + 0 transcript_id "g1.t1"; gene_id "g1";
2023年3月3日、遺伝子予測プログラムBrakerの改良版であるBraker3が公開された。Brakerは、RNA-seqなどのtranscriptsをhintとしてゲノム中の遺伝子予測を行うBraker1と、近縁種などのタンパク質アミノ酸配列をhintとして遺伝子予測を行うBraker2があるが、Braker3はそれらの両方を用いて遺伝子予測を行うとのこと。
（2023.3.8現在）今回、Ubuntu20.04の環境にBraker3を導入したので、インストールについての備忘録的なメモを書いておく。実際の解析結果については、結果が得られたら改めてまとめていくつもり。
OS: Ubuntu 20.04 LTS
CPU: Intel Xeon(R) CPU ES-2620 @ 2.10 GHz x 24
RAM: 256 GB
Install disk: 2 TB HDD
Brakerは依存関係が多いので、一般的にはAnacondaを使ってインストールする（たとえばこのサイト）のが手間が少ないと思われるが、今回はBraker3を使用するため、Anacondaを使わず、依存関係にある全てのプログラムを手動でインストールすることにした（すぐに対応すると思いますが・・・）。
Brakerを動かすにあたって必要なプログラムの一覧、またPerlのモジュールなどは、Brakerのサイトに記載されているので、まずはこれらを一つずつインストールした。
（230318追記）Anaconda3 (miniconda3)で仮想環境を作って、Perl自体も含めてモジュールをconda installするのが良い（参考URL）
（230327追記）↑この場合、GeneMark-ES/ET/EPのPerl and Python scriptのヘッダ部分（使用するPerl等のPATHを指定する部分）をAnacondaに対応して書き換える必要がある こちらのURLを参考にした
# "my_path_to_AUGUSTUS"は自身のaugustusのプログラムがあるところのPATH export AUGUSTUS_CONFIG_PATH=/my_path_to_AUGUSTUS/augustus/config/
braker.plのあるディレクトリにパスを通す。
今回は自分ではやらなかったけど、これが一番簡単&確実なような気がしてきた・・・
hub.docker.com
