えーと。ありがたいツッコミなんですが、まず、ベンチ屋にはこの呪文が何を意味しているのかわかりません。こういう時は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 \\nArraystar_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さんありがとうございました。
2020