2016年12月20日(火)

Arraystar Mouse LncRNA microarrayのGEOサーフィンリターンズ

投稿者:

「ネット」=「妖しいもの」というのが一昔前は定番だったと思いますが、今ではすっかり市民権を得て、新聞やテレビよりも信頼できる、なんて言っている人もいるようです。その真偽はともかくとして、ネット上のSNSは情報収集という観点からはなかなかバカにならないものがありまして、先日の当ブログのエントリーにも、早速rnacintosh LC475さんからツッコミがありました。

えーと。ありがたいツッコミなんですが、まず、ベンチ屋にはこの呪文が何を意味しているのかわかりません。こういう時はGoogle先生、ということで、まずはsedあたりをググってみますと、要はWordで言うところの検索&置換をするコマンドのようです。ということはこれは単に文字を変換しているだけ、、、あっ、そうか。プローブのリストを次世代シークエンサーで得られるリードのファイルのように扱って、それをマッピングしてやれば、UCSC Genome BrowserやIGVなどで見ることで、好きな遺伝子のところにデザインされてるプローブがすぐ分かる、ということですね。アッタマいいー!

csvファイルじゃなくてtxtファイルを読み込ませているのは、多分何か複雑な事情があるのだろうということで、そこは空気を読んでArraystar_probe_list.csvの拡張子だけ.txtに入れ替えたファイルを作り、この呪文を早速Rのコンソールに放り込んで、リターン!

> cat Arraystar_probe_list.txt | sed 's///>/g' | tr ',' '\n' > out.fa
Error: unexpected symbol in "cat Arraystar_probe_list.txt"

ん?何も起きぬ。。

そりゃそうですね。よく考えたら教えていただいたのはターミナルで使う呪文でした。これだから素人は、、、
気を取り直して、Macのターミナルを立ち上げて、同じ呪文を打ち込んで、ナッツリターン!

ん?変換されているけれども、何か変。。。out.faの中身を見てみると、

>AA014456_P1
GGGCGTTACCTGCCTTACACTTGAAGATCTTCCTCTGCATACAGTATTTTTTCCCTTCAC
AA017743_P1
CTTGAGTCCAGCTTTGAGGTTGGCAAATTGAACACAATACTTGGCTGTACCTTCAGCCCT
AA066037_P1
GTAATGGAGACTGCCAGGCTTTCGCCTCCAGACTCAAGACCCAAGACCTTTCAAAGAACC
..........

Fasta形式は遺伝子名をが>の後に来て、改行して配列、その繰り返し、というフォーマットのはずですが、なぜか一行目しか正しく変換されていません。うーむ。こまった。どうやら改行コードが正しく認識されていないようです。これが木立さんの言うところのフォーマット変換の罠というやつでしょうか。プロでも引っかかるわけですから、素人なら当然。

これまた改行コードについてGoogle先生に聞いてみると、大阪大学の山中卓さんという方のページに何やら役に立ちそうな情報があります。ちなみにこのへんの記事は「教授でもできる...」というカテゴリーにまとめられておりまして、おっさんでもこの程度はできるんだから若い学生さんたち頑張ってね、というメッセージがひしひし伝わってきます。レベルはだいぶ違いますが、当ブログの一連の「ベンチ屋がコマンド打ってみた」恥さらし企画に通じるところも、、、

ともあれ、ここに書いてある通りodなるコマンドで改行コードを調べてみると、

od -c Arraystar_probe_list.txt|less

0001120    _   P   1   ,   A   A   C   T   G   C   A   G   A   G   A   G
0001140    G   A   T   A   C   T   T   T   G   A   A   A   G   C   A   A
0001160    A   T   G   A   A   C   C   T   C   T   T   A   C   T   T   G
0001200    T   G   G   T   T   G   G   C   T   C   C   C   T   G   T   G
0001220   \r   A   A   1   1   7   5   4   4   _   P   1   ,   A   A   A

確かに、\rという改行コードになっています。unix系だと\nにしないといけないようです。
trなるコマンドを使えば改行コードは簡単に変換できるそうなので、

 tr \\r \\n  Arraystar_probe_list_unix.txt

として、unix風の改行コードにしたファイルを作っておいて、再び復活の呪文を入れて参勤交代リターンズ!

 cat Arraystar_probe_list_unix.txt | sed "s/^/>/g" | tr ',' '\n' > Arraystar_probe_list_unix.fa

Arraystar_probe_list_unix.faの中身を見てみると、、、

>AA014456_P1
GGGCGTTACCTGCCTTACACTTGAAGATCTTCCTCTGCATACAGTATTTTTTCCCTTCAC
>AA017743_P1
CTTGAGTCCAGCTTTGAGGTTGGCAAATTGAACACAATACTTGGCTGTACCTTCAGCCCT
>AA066037_P1
GTAATGGAGACTGCCAGGCTTTCGCCTCCAGACTCAAGACCCAAGACCTTTCAAAGAACC
....

無事、Bowtieが食ってくれるFasta形式のファイルを作ることができました。あとは、普通にマッピングして、sam>bam変換して、ソートして、インデックスつけるだけ。

bowtie2 -f -x genome -U Arraystar_probe_list_unix.fa -S Arraystar_probe_list_unix_bt2.sam
samtools view -bS Arraystar_probe_list_unix_bt2.sam > Arraystar_probe_list_unix_bt2.bam
samtools sort Arraystar_probe_list_unix_bt2.bam -o Arraystar_probe_list_unix_bt2.bam
samtools index Arraystar_probe_list_unix_bt2.bam

一応、mm10にマップして出来上がったbamとbaiファイルを下のリンクからダウンロードできるようにしておきました。IGVで自分の好きな遺伝子を入れれば、そこにマップされてるリードがすなわち目的のプローブということになります。これで、ここ数日の懸案だったArraystarのプローブ名教えてくれえ事案は、一件落着、と相成りました。rnacintosh LC475さんありがとうございました。

中川 真一

北海道大学 薬学研究院 教授
▶ プロフィールはこちら

コメントする

コメントはどなたでも歓迎します。メンバーの方はログインしてからコメントして下さい。

ブログアーカイブ

ログイン

サイト内検索