STAR-RSEM による発現量推定

UB3/informatics/bioinformatics/ngs_expression_matrix

このページの最終更新日: 2023/11/23

  1. STAR-RSEMによる発現量推定

広告

STAR-RSEMによる発現量推定

2023 年 10 月、Mac でトライ。brew install star で install されたが、いざ実行してみると star: Bad Option: --runMode という エラー。これは、brew では違う star というプログラムを インストールしているということで、やはりソースファイルからコンパイルする必要がある。

github に従い進める。gcc のバージョンがマニュアルと異なっている。以下を実行。

make STARforMacStatic CXX=/usr/local/Cellar/gcc/13.2.0/bin/g++-13

このエラーが出てしまう。brew install gcc@8 でマニュアルにあるバージョンの gcc をインストールし、実行しても同じエラー。

Assertion failed: (resultIndex < sectData.atoms.size()), function findAtom, file Relocations.cpp, line 1336. collect2: error: ld returned 1 exit status

XCode エラーらしいので、Apple Store から 15.0.1 をインストールするが、状況は変わらず。ネットで XCode 15.1 のベータ版ではこのエラーが解決されているという書き込みがあったので、このバージョンをダウンロード。.xip という見慣れない拡張子。これは圧縮ファイルのようで、解凍すると XCode-beta.app になる。このアプリを立ち上げると、関連コンポーネントのインストール画面が現れる。

Command_Line_Tools_for_Xcode_15.1_beta.dmg というのもダウンロード。これは普通のdmg なので、指示に従ってインストール。これで make がちゃんと動くようになり、gcc 13 でもコンパイルできた。コンパイル後は source に STAR というファイルが現れるので、それをパスの通っている /usr/local/bin に移して完了。

star に必要なのは、ゲノムの fasta file とアノテーションの GTF file。まずは このページ に従い、インデックスを作成。

ゲノムファイルがあるディレクトリを指定した上で、.fa ファイルと .gtf ファイルを両方指定する。

続いてマッピング。

設定したオプションは以下の通り。

  • --outFilterMismatchNmax 999: ミスマッチの数。large number switches off this filter と書いてあるが、このオプションがオフになるとどうなるのか?
  • --outFilterMultimapNmax 20: 一つのリードがマップされる最大数。もしこの数以上にマップされた場合、そのリードは "unmapped" に分類される。

このあたりで、やはり大きなメモリが必要なこと、Windows では動かないことなどがネックになってくる (参考)。できれば、16 GB RAM ぐらいの Laptop で走らせたい。ということで徐々に Rsubread へ移行。STAR は使わなくなった。


STAR-RSEMによる発現量推定 の通りに進めてみる。2021 年 12 月、Ubuntu 18.4 を使用。

まずは fastq-dump で SRR065504, SRR065495, SRR065509, SRR065524 の 4 つをダウンロード。O/N でやったので何時間かかったかわからないが、これには時間がかかる。

続いてゲノム配列をダウンロード。新しいアセンブリが出ているかもしれないが、とりあえずは UCSC hg38 で練習ということで 常染色体と性染色体のみのゲノム配列ファイル genome.fa を作成する に従い hg38.fa をダウンロード。

さらに git clone を使って DROMPA3 をダウンロードし、splitmultifasta.pl を実行しようとすると Can't locate Path/Class.pm in @INC というエラー。

これは fasta を分割するステップなので、seqkit でも同じことができる。seqkit stat で配列が 455 個ということがわかったので、seqkit split -p 455 として hg38.fa を 455 個のファイルに分割。無事に分割できるのだが、splitmultifasta.pl はファイル名を配列名と同じにしてくれるのに対し、seqkit split だと hg38.part_001.fa のようになってしまう。ここから、どれが chr1 か探すのはちょっと手間がかかるので、エラーを解決することに。

You may need to install the Path::Class module というメッセージがあったので、これに従い cpan Path::Class を実行。これで splitmultifasta.pl が動くようになる。常染色体と性染色体のみを含む genome.fa を作成。ミトコンドリアは入れなかった。

Gene annotation データを用意する (gtf形式) に従い、Ensembl からアノテーションデータをダウンロード。chr を追加し、コメントアウトの行を除いた Homo_sapiens.GRCh38.94.addchr.gtf を作成。これを使用する。

RSEM, STAR のインストールは順調に進み、一時的にパスを通す。ミトコンドリア DNA の名前も置換して Homo_sapiens.GRCh38.94.addchr.new.gtf を作成。


広告

References

  1. リードカウントのマトリクスを作成する. Link: Last access 2023/10/15.

コメント欄

サーバー移転のため、コメント欄は一時閉鎖中です。サイドバーから「管理人への質問」へどうぞ。