RPy2を使ってPythonからRを呼び出してデータ解析

Pythonで書いたコードの中でRを呼び出す必要に迫られたので、備忘録を兼ねてやり方を記録しておきます。

この前から取り組んでいるPPiプロジェクトの一環で、タンパク質相互作用ネットワークを超幾何分析してスコアを記録しよう、みたいな話になった。

超幾何分析自体の話はbeyond the scopeすぎるので専門家に譲るとして、分析そのものはphyperというR言語の関数でちゃちゃっとできる。

phyper(q, m, n, k, lower.tail = FALSE, log.p = FALSE)

q:調べる蛋白質のモジュール内のリンク数

m:モジュール内のリンクの総数

n:モジュール外のリンクの総数

k:調べるタンパク質の全リンク数

マニュアル(英語)

それで、この前Pythonで生成したファイルにRを走らせるのが一番単純な手法かもしれないけど、面倒だしスマートじゃないので、やはりPythonコードのなかでRを動かしたいところ。

少し調べると、PythonとRの連携はRPy2というパッケージが人気の模様。 ちゃんとPython 3.xでも動作するということで、試してみる。

まずは$ sudo pip install rpy2 みたいな感じでローカルにRPy2をインストール。

それからPythonコードで読み込みます。

from rpy2.robjects import r

Rの関数の引数としてPython側の変数を使いたいものの、どうもうまくいかない。 調べてみると直接使わずにR側の新しい変数にポートして活用するらしい。

r.assign('RVariable', PyVariable)

という感じ。 これを引数に使う変数全部に対して行って、いよいよ関数にぶちこみます。

r('関数名(引数1,引数2,...)')

実際のコードでは、まずnamelist1wholelistからそれぞれモジュール内のリンクの総数moduleedgesとタンパク質相互採用ネットワーク全体のリンク数wholeedgesを計算し、wholeedgesからmoduleedgesを引くことでモジュール外のリンクの総数nonmoduleedgesを得ました。

それから、モジュール内のタンパク質の名簿であるpnamelistをループで回して、各タンパク質のモジュール内リンク数と総リンク数をそれぞれRにポートし、最後にそれをもとに超幾何分析を行ってdict型に収納、という手順を踏みました。

from rpy2.robjects import r
moduleedges = len(namelist1)
wholeedges = len(wholelist)
nonmoduleedges = wholeedges - moduleedges

from rpy2.robjects import r
#超幾何分析結果を収納するdictを用意
hypergeo = {}
#pythonの変数をRの変数にポート
r.assign('module', moduleedges)
r.assign('nonmodule', nonmoduleedges)

#モジュール内の各タンパクについて超幾何分析
for x in pnamelist:
    r.assign('intralink', intralink[x])
    r.assign('wholelink', wholelink[x])
    p = r('phyper(intralink, module, nonmodule, wholelink, lower.tail = FALSE, log.p = FALSE)')
#dictに分析結果を収納
    hypergeo[x] = p[0]

あとは前回のようにcsvに出力して終わり。

僕以外にも統計分析にはRを使うけどそれ以外は全部Python、みたいな人は多いらしく、今回使ったRPy2のようにPythonの中でネイティブなRを走らせる以外にも、Rの大人気グラフ作成パッケージggplot2をPythonで再現する試みとかいろいろあるみたいです。

やっぱりデータ分析分野はR言語で書かれた遺産が多いので、うまく使っていきたいところです。

今回参考にしたウェブサイト

バイオインフォマティクスによる遺伝子発現解析

R: The Hypergeometric Distribution

http://rpy.sourceforge.net/rpy2/doc-2.1/html/

http://www.okada.jp.org/RWiki/?Python%20%A4%C7%A1%A1R

生物系学部4年生が2013年に買って良かったもの

今年もいろいろ買いましたなあ。

Macbook Air 2013

ずっとVAIO使ってて、Macに若干苦手意識があったものの、今年からわりとプログラム書いたりし始めてVirtual BoxのUbuntuを立ち上げるのが面倒になってMacに書い変えました。最初はサブ機で使おうと思っていたのに、結局VAIOはほとんど放置。

f:id:pioneerboy:20131017225255j:plain

Macbook Air 13inch (Mid-2013)届きました! - 雨に煙る

Roost

これがなかったらMacbook Airの魅力も半減しているところ。どこへ行くにも持って行っています。

f:id:pioneerboy:20131026113655j:plain

最強ノートパソコンスタンド"Roost"が届いたよ - 雨に煙る

土屋鞄のショルダーバック

春先にバイトのボーナスで買って以来、かばんはこれひとつ。すこし味が出てきたかな? f:id:pioneerboy:20131231154134j:plain

iPad Air

論文読んだりRSSフィードチェックしたりマンガ読んだり、iPhoneから仕事を奪いまくっています。Twitterだけはお気に入りクライアントのThe WorldがiPad未対応なのでiPhoneメイン。

f:id:pioneerboy:20131116142920j:plain

iPad Airをお迎えしました - 雨に煙る

ハンドプレッソ

好きなときに電源いらずでエスプレッソを飲める、なかなか頼りになる手動エスプレッソメーカー。 やっぱり粉をタンピングするのが億劫になるのでカフェポッドと併用した方がいいですね。

f:id:pioneerboy:20131012165636j:plain

ハンドプレッソでカフェラテがぶがぶする - 雨に煙る

PCメガネ

JINSでしばらく前にフレームだけ買ってたヱヴァンゲリヲンコラボのアスカモデルに、あとからPCレンズを入れてもらいました。やや目の疲労が軽減された気がしますが、プラシーボかもしれません。

f:id:pioneerboy:20140105110356j:plain

「クリぼっち」のバカらしさを統計から攻めてみる

クリスマスといえば、この前Twitterに流れてきたこの画像を見て日本の統計リテラシーの低さを痛感しました。

f:id:pioneerboy:20140625115235j:plain

なにを問題としているか説明するために、まずは統計の手法について復習してみましょう。

日本に暮らしている人のなかでクリスマスをひとりきりで過ごす割合がどれくらいか調べたいと思った時、考えられる方法はざっくり分けて2つあります。

1つめは、日本に暮らしているひと全員に「あなたはクリスマスをひとりぼっちで過ごす予定ですか?^^」と聞いて回ることです。これを全数調査と呼びます。しかし、1.2億人ほどいる対象にそんなことをするのは非現実的です*1

そんなわけで、ほとんどの計量調査は2つめの手法、標本調査をつかって行われます。全数調査が不可能な場合、母集団(この場合は日本人全体?)の情報をうまく推定するため、構成因子(この場合は個人)を何人か抽出して統計分析を行います。

ここで気をつけなければいけないのが、「どれだけ調査結果を信用できるか?」という問題です。 全数調査の場合は全部の因子を総当りしているので、調査に使われた集団、つまり標本が母集団そのものと等しいわけですから、調査結果は母集団の性質そのものを表していると自信を持って良いでしょう。

しかし、標本調査の場合は、常に標本誤差を意識しなければなりません。これは、標本の集団が母集団の一部である関係上、母集団の真の性質と、標本から得られた性質が必ずしも一致しないことから起きる誤差のことです。たとえば、本来は日本人の中のたった数%しかいないキリスト教徒だけを標本として調査して、そのデータを元に「日本人の100%がキリスト教徒である」と言ってしまう危険性を意識しなければならないということです。

こういう事態を避けるために、標本調査を行う際は必ず標本の偏りを小さくする努力が行われます。それは、一般的に標本に使う因子を無差別に抽出することで成し遂げられます*2。先の例で言えば、1.2億人いる日本人をランダムに調査すれば、おそらくキリスト教徒は100回のうち数回しか出てこないので、そのデータを元により「真の性質」に近いデータを出すことが可能となります。身近なところだと、政党支持率などの世論調査で使われる、無作為に生成した電話番号に電話をかける手法ですね。

以上のことを踏まえて最初の画像を見なおしてみましょう。

f:id:pioneerboy:20140625115235j:plain

僕はこのようなデータをみたときに、「誰が調査したのか?」と「どうやって調査したのか?」の二点が一番気になります。標本調査の場合、それが「どれくらいこのデータを信用できるのか」に関わってくるからです。 この数字の調査元のネオマーケティング社は、インターネット上のアイリサーチというアンケートプラットフォーム上でこんなかんじで調査を行っている会社です。 つまり、アイリサーチのモニターとして登録してアンケートに答えた奇特な方々が調査対象なのですね。おやおや、標本としてめちゃくちゃ偏っていそうですね。

実際の調査レポートを見てみましょう。 こちらからダウンロードできます。https://www.i-research.jp/report_dl/dl68.html

f:id:pioneerboy:20131225004410p:plain

人数は500人、性別は男女およそ半々、年齢層は20代のみ、そして未婚。 未婚の20代に絞った調査とは、日本人全体を推定するにはずいぶん条件が厳しいですね。

結果はこちら。

f:id:pioneerboy:20131225004421p:plain

あれ?「クリスマスに一緒に過ごす異性のパートナーはいますか」?? この質問、クリスマスをひとりぼっちで過ごすのかではなく、一緒に過ごす異性のパートナーの有無を聞いていますね。 もはや調査内容そのものが違っていました。

そしてアイリサーチによる要約文がこちら。

f:id:pioneerboy:20131225004754p:plain

いや〜傑作ですね。クリスマスに異性のパートナーと過ごさない割合の項目のタイトルが「クリスマスに一人きり」ですよ。よっぽどの馬鹿か、よっぽどの悪意をもった馬鹿が書いたのでしょう。そして女性のほうが男性に比べて「クリスマスは異性と過ごす」と答えた割合が高かったというデータの結論(感想?)が「近年、‘草食系男子’という言葉が流行したように、恋愛に奥手な男性が増えているのかもしれません。 」ときました。このデータから非リア男子が増えていることがわかるんですか?^^

以上、当調査の問題をまとめると

  1. モニター登録者の20代未婚者というめちゃ偏った標本集団
  2. 少人数(男女約250人ずつ)
  3. 質問内容と結論の乖離

という感じですね。

いやはや、20代のネットユーザー対象の統計をさも日本人一般の性質と言わんばかりに放映したNHKの統計リテラシーの低さもさることながら、調査元のネオマーケティング社による「クリスマスを異性のパートナーと過ごさないならお前は一人ぼっちだ!」という超論理にも脱帽です。

よいこのみんな、まねしないでね。

*1:もっと予算がつきそうな調査なら話は別です、国勢調査とか

*2:もちろんランダマイズは標本のバイアスを避ける十分条件ではありません

Pythonでタンパク質相互作用ネットワークを分析する

金曜日頃からやや風邪気味だったので週末は自宅でひっそり過ごし、今日研究所バイトのために外出したら冬の到来を体感しました。寒い…。

冬生まれの自負と、汗っかきゆえの蒸し暑さ耐性の低さのために、いつも夏と冬だったら冬のほうが過ごしやすいと答えることにしているけど、こういう時はちょっと夏もいいかなあと思ったりする*1

研究所で鼻すすりながらPythonのお勉強をしつつコードを書いて、とりあえず欲しいデータを入手することに成功したので達成感に浸る。 備忘録を兼ねて、スパゲティージェノヴェーゼコードを投下します。

*1:そもそもその質問は夏場に発される場合が多い気がする

続きを読む
Related Posts Plugin for WordPress, Blogger...