So-net無料ブログ作成
検索選択
統計学もどき ブログトップ

RでJAGSとやらを使うつもり [統計学もどき]

以前に検索で見つかったこちらのPDFから
http://www.rri.kyoto-u.ac.jp/NSRG/FSEB/abstracts2015-8.pdf

村瀬香・名古屋市立大学大学院システム自然科学研究科

階層ベイズと統計モデリングを通じて、原発事故が野生動物とヒトに与えた影響を考える

という発表が気になりまして、これらを勉強するために

データ解析のための統計モデリング入門――一般化線形モデル・階層ベイズモデル・MCMC (確率と情報の科学)   久保 拓弥 (著)

というのが定評があり、幸い県立図書館に置いてましたのでリクエストして借りました。
なかなか簡単に読み進めませんが、階層ベイズモデルを使うにはRだけではだめで、WinBUGSというのを使わないといけないと書いてありました。一応インストだけしてみるかとチャレンジしたところ、いまどきパッチをあてるなどと言うことを要求したり、Win7で使うにはいろいろややこしいことを要求したのであきらめかけたところ、村瀬香氏の論文を探すと#原子力発電_原爆の子さんブログに翻訳があり、ここを見ると分析に使ったソフトが書いていました。

http://besobernow-yuima.blogspot.jp/2015/04/httpt_3.html
#Nature 誌サイエンティフィック・リポーツ「オオタカの繁殖に対する福島第一原発事故の影響」

「われわれは、R 3.0.2上でJAGS 3.3.0とrjags 3–12を用いてMCMCサンプリングを実施した。サンプリング過程で、最初の10,000回試行はしくじりとして破棄され、MCMC連鎖ごとに100回試行をサンプリングとして、200,000回試行を実施した。5連鎖に対して、同じ手順が適用され、10,000件のサンプルが生成された。われわれは、収束兆候を求め、また値が1.0に非常に近接することを確かめるために潜在スケール減少係数R30の多変数版を計算した 30。」

JAGSというのがWinBUGSの代わりで、rjagsがRからJAGSを動かすものとわかりました。
久保 拓弥氏のHPにちゃんとインストの仕方が出ておりまして

生態学のデータ解析 - R で JAGS (rjags 編)
http://hosho.ees.hokudai.ac.jp/~kubo/ce/RtoJags.html

一応 RやEZRから JAGSが使えるようになったはずです。サンプルを試してみたけどうまくいかずあきらめかけていたところ、こちらにあったサンプルは

JAGSを使ってギブスサンプリングを試してみた (1)
http://www.singularpoint.org/blog/r/mcmc-jags-1/
JAGSサンプル.png
意味はわかっていませんが、動いたみたい。

Rでポアソン分布の共役事前分布とやらを試したつもり


nice!(1)  コメント(24)  トラックバック(0) 
共通テーマ:blog

エクセルでやったらどえらい簡単にできたが間違ってないか? [統計学もどき]

EZRの記事がコメントが増えたのと、今月も1つくらいは書いておかねばと思いましたので、新たに設けたカテゴリーで行ってみたいと思います。

アワプラさんが貴重なまとめをしていただいています。

甲状腺がん悪性・悪性疑い152人〜福島県民健康調査
http://www.ourplanet-tv.org/?q=node/2004

検討委員会での記者会見にて、基本調査の①-8ページの福島県住民の外部被ばく線量分布と、甲状腺検査(本格検査)の②-6ページの甲状腺がん症例の外部被ばく線量分布が大きく異なり、甲状腺症例の外部被ばく線量が高い方に分布しているという指摘があった。このオッズ比に関して、性別を調整して表2として示した。この調整オッズ比とその95%信頼区間の推定結果は、EpiInfo 7TMのStatCalcの最尤推定値(Mid-P)を用いた。95%信頼区間の下の値が1より大きいと、いわゆる「統計的有意差がある」ということになる。
  その結果、明瞭な量反応関係が観察された。つまり、外部被ばく線量に関連する何らかの要因(恐らく内部被ばく)が甲状腺がんの発症に関与していることが、より明瞭に推察できることになった。1mSv未満もまた、被ばくしていると考えられるので、示されたオッズ比は過小評価されている。先行検査では、平成23年度、24年度、25年度の別に、甲状腺がん症例の線量分布が発表されておらず、平均有病期間が平成26年度としてほぼそろっている本格検査とは異なり、明瞭ではなかった。基本調査の住民の性・年齢・市町村別線量分布、先行検査および本格検査における甲状腺がん症例の性・年齢・市町村別線量分布の発表が待たれるところである。
2015113014.jpg 
これだけでも十分なのですが、ちと統計学の相関関係が気になってあれこれ探して、Rで試すのはなんだかまだ私には無理そうだったので、エクセルのはないかと見たところ

散布図と相関係数(相関分析)
http://hitorimarketing.net/tools/correlation-analysis.html

これに倣ってやってみたところ(発症者数は100万人当たりに換算してます。)

相関係数.png

案外 楽にピアソンの積率相関係数(単相関係数)とやらが表示されました。簡単すぎたのでこれでいいのかと思ってしまいましたが、点が三つしかない散布図でも、正の相関になっていますので、多分これでいいのでしょう。Rに慣れていけば、P値なんかもわかるので、相関係数の正確さもわかってくるでしょう。

毎月3日のアベ政治を忘れておった。流行語大賞ノミネートまでいったんでした。
sho_fj.jpg
nice!(9)  コメント(25)  トラックバック(0) 
共通テーマ:blog

いくらか使えたEZR [統計学もどき]

あんまり奥村教授のところのお世話にはなりたくなかったのですが、EZRがなんぼか使えたので記念撮影しました。

EGR-X2検定.png

2×2の表,オッズ比,相対危険度
http://oku.edu.mie-u.ac.jp/~okumura/stat/2by2.html
----------------------------------------------------------------------------------------
カイ2乗検定による方法
これはオッズ比などを求めるものではありませんが,行・列の独立性の検定をします。デフォルトでは連続性の補正をします。

> chisq.test(x) # 補正あり

Pearson's Chi-squared test with Yates' continuity correction

data: x
X-squared = 3.4808, df = 1, p-value = 0.06208

> chisq.test(x, correct=FALSE) # 補正なし

Pearson's Chi-squared test

data: x
X-squared = 4.8577, df = 1, p-value = 0.02752
補正あり・なしで p 値はそれぞれ 0.06208,0.02752 です。
----------------------------------------------------------------------------------------
EZRの上のフォームに1行づつ入力して実行ボタンを押すと、ちゃんと同じ数値が出ましたね。
ちょうど「統計のはなし」で読んでるX(カイ)2乗検定の部分で
X-squared  カイ2乗値が 4.8577、 df  自由度が1  p-value P値が 0.02752
ちゅうとこなんでしょうが、まだ説明はできません。

 これだけではあまりにも芸がないので、数値や言葉をこの前の茨城のケースに変えてみました。
 まあ単純なオッズ比はエクセルと同じなのは当然ですが、フィッシャーテストとやらのオッズ比も近い数字が出ました。カイ2乗値がやたらでかい。
EZR-X2テスト.png

このままの数値で放置しておくわけにいかないので、追記します。
先月ツイッターで見つけた『生活の党』を支援する道民の会(浜菊会)さんのブログで津田先生のオッズ比についてツイッターをまとめたものです。
http://blog.goo.ne.jp/hamagikukai/e/df8da73b1f964e21eb0cf501d70256e1
こちらの方にもあります。
http://blog.goo.ne.jp/hamagikukai/e/c5f58cfc1c97a05bb73f02cb922e1928
「②がんセンター数値の年少罹患率、男0.06、女性0.12を基本とし、スクリーニング効果100倍を見積もって、18/十万とする。因みにこれは全年齢の女性(男性より何倍か好発傾向)における年齢調節罹患率よりも高い。よって、対照群は、甲状腺癌アリ18、ナシ99982、となる。」

これを使わせていただいて、暴露なしを「甲状腺癌アリ18、ナシ99982」として再びEZRを試すと
EZR-X2テスト2.png
カイ2乗検定、Fisherの正確検定の数値もまともになりました。
スクリーニング効果を認めるわけではないのですが、暴露なしデータがはっきりしませんので、事故前全国平均を使わずに緩めに数値を見ても、事故後の甲状腺がんの発病は危険水域です。
あ これは茨城のデータです。
nice!(6)  コメント(32)  トラックバック(0) 
共通テーマ:blog

少しもEasyでないEZR [統計学もどき]

いろいろやってる場合ではないといいながら、エクセルばかりになるといけんのでRも試したいが、RStudioが英語につき日本語で使いやすいのはないかと探してみたらEZR(Easy R)なるものがあるというのでちょいとインストしてみました。

なんと自治医科大学のさいたま医療センターなんですねえ。
http://www.jichi.ac.jp/saitama-sct/SaitamaHP.files/statmed.html

勝手が違うのか、弁当屋のサンプルを試してみようとしてもでてこん、というかデモの起動の仕方がわからないという、初歩の初歩状態です。

EZR.png

一応起動はしますが、簡単なマニュアルはあるけど、これで使いこなすようになるのは、私には無理なので、またぞろ本を物色する必要があるのかな。

初心者でもすぐにできるフリー統計ソフトEZR(Easy R)で誰でも簡単統計解析
http://www.nankodo.co.jp/g/g9784524261581/
EZRでやさしく学ぶ統計学 改訂2版
http://www.chugaiigaku.jp/item/detail.php?id=1660

以上2冊が出てますが、結構高い。
nice!(3)  コメント(22)  トラックバック(0) 
共通テーマ:blog

RStudioに挑戦したつもり [統計学もどき]

今月は時間がないといながら、仕事の合間に気分転換です。
挑戦したと言えるほどのことではなく、「とある弁当屋の統計技師」のサンプル表示がうまくいっただけのことですが、一応動いたのでこの先いろいろ練習はできそうな気配です。
画像は、ロジスティック回帰分析の解説のある第6章のサンプルを表示したものです。
Rstadio.png
途中なにやらエラーは出るのですが、まだ意味不明につき、一応サンプルが最後まで表示されたのでよしとしましょう。
こちらに説明があります。
とある弁当屋の統計技師 第6章 ロジスティック
http://rmecab.jp/ranko/chap6.html

まあ先は長いな。
nice!(7)  コメント(24)  トラックバック(0) 
共通テーマ:blog
統計学もどき ブログトップ