JAGSでマルコフ連鎖モンテカルロ法

連休は旅行いったりとか子供と遊んでたんだけどそれ以外はJAGSいじってた。

ProductName マルコフ連鎖モンテカルロ法 (統計ライブラリー)
豊田 秀樹
朝倉書店 / ¥ 4,410 ()
在庫あり。

5章をWikiに書いていくことにした。

マルコフ連鎖モンテカルロ法をJAGSで

分散分析のあたりを出来るようにしておくとよいと思った。

  • 5.1 ロジスティック回帰モデル 
  • 5.2 メタ分析 
  • 5.3 多項ロジットモデル 
  • 5.4 対数線形モデル 
  • 5.5 ポアソン回帰 
  • 5.6 2値データに対する回帰分析 
  • 5.7 トービット回帰モデル 
  • 5.8 変曲点のある回帰分析 
  • 5.9 生存時間分析(ワイブル回帰) 
  • 5.10 生存時間分析(コックス回帰) 
  • 5.11 時系列モデル 
  • 5.12 分散分析 
  • 5.13 分散成分分析 
  • 5.14 分散分析(枝分かれ配置) 
  • 5.15 一般化可能性理論 
  • 5.16 反復測定データの分散分析 
  • 5.17 階層線形モデル 
  • 5.18 項目反応理論(2母数2値モデル) 
  • 5.19 項目反応理論(段階反応モデル) 
  • 5.20 項目反応理論(名義反応モデル) 
  • 5.21 項目反応理論(部分採点・評定尺度モデル)
  • 5.22 項目反応理論(連続反応モデル)
  • 5.23 多次元IRT
  • 5.24 項目反応理論(混合名義反応モデル) 
  • 5.25 項目反応理論における特異項目機能(DIF) の分析 
  • 5.26 正規混合モデル 
  • 5.27 潜在クラス分析 
  • 5.28 成長曲線モデル 
  • 5.29 非線形成長曲線モデル 
  • 5.30 因子分析 
  • 5.31 多母集団分析 
  • 5.32 非線形SEM 

あとこの本だけだと難しすぎるのでPRMLの下巻が必須(多分)。

ProductName パターン認識と機械学習 下 - ベイズ理論による統計的予測
C. M. ビショップ
シュプリンガー・ジャパン株式会社 / ¥ 8,190 ()
在庫あり。

hikiにコードを貼付けるためのちょっとしたコマンドのメモ

markdownだとspace4つだがhikiの場合はspace一つなのでrubyで

#!/usr/bin/env ruby
File::open(ARGV[0]) {|f| f.each {|line| print " " + line}}

コードを貼付けるためにスペース入れたりするのは結構めんどうなのではてなの記法はなかなか便利(ただ忘れやすい)

wikiの記法をつくる記法(metaな記法)があるとよいんではなかろうかとか思う。wiki記法ってDSLみたいだもんね。

Jericho / Weekend Players

ソウルフル、エモーショナル。静かにしぼりあげる感じの歌声が好き。

ProductName Pursuit of Happiness
Weekend Players
Multiply / ¥ 1,355 (2003-03-11)
在庫あり。

レビューにあったけど、EBTGよりはポジティブな感じ(比較したことなかったが。)


最近iTunesの曲にレートをふるブームが来ている。レートのついてない曲をランダムで250曲ほど表示するスマートリストを使ってるが、ショートカットでレーティングする方法がわからんので使いづらい。

amazon認証プロキシ

amazonの認証に対応するためにhikiのバージョンを0.8.8.1にあげたついでにamazonプラグインを覗いてみたらプロキシで対応していた。

drkcoreのほうはNet::Amazon使っているので、ローカルでSignature組み立ててるんだけど、wikiりたいだけのユーザーがわざわざsecret_keyとか取得すんのもあほくさいからwikiみたいのはプロキシーでいいのかと。

ProductName パターン、Wiki、XP ~時を超えた創造の原則 (WEB+DB PRESS plusシリーズ)
江渡 浩一郎
技術評論社 / ¥ 2,394 ()
在庫あり。

Re:やさしい気持ち

これはどうやろか。言葉が増えすぎて情景が失われてる感じ。

ProductName Re:やさしい気持ち
HALCALI
ERJ(SME)(M) / ¥ 1,020 (2009-02-11)
在庫あり。

リミックスすんだったら、もっと切り刻んだほうが。あと元の曲の声いれんだったら、voいらん気がした。remixに期待するのはメタな部分の抽出か、別の解釈なんだなぁと改めて。

JAGSでメタアナリシス

マルコフ連鎖モンテカルロ法の5-2をJAGSで

ProductName マルコフ連鎖モンテカルロ法 (統計ライブラリー)
豊田 秀樹
朝倉書店 / ¥ 4,410 ()
在庫あり。

モデル

model
{
  for( i in 1 : Num ) {
    rc[i] ~ dbin(pc[i], nc[i]);
    rt[i] ~ dbin(pt[i], nt[i]);
    logit(pc[i]) <- mu[i];
    logit(pt[i]) <- mu[i] + delta[i];
    mu[i] ~ dnorm(0.0,1.0E-5);
    delta[i] ~ dnorm(d, tau);
  }
  d ~ dnorm(0.0,1.0E-6);
  tau ~ dgamma(0.001,0.001);
  delta.new ~ dnorm(d, tau);
  sigma <- 1 / sqrt(tau);
}

Rのコード

library(rjags)

data <- list(rt = c(49, 44, 102, 32, 85, 246, 1570),
     nt = c(615, 758, 832, 317, 810, 2267, 8587),
     rc = c(67, 64, 126, 38, 52, 219, 1720),
     nc = c(624, 771, 850, 309, 406, 2257, 8600),
    Num = 7)

inits <- list(d = 0, delta.new = 0, tau=1,
     mu = c(0, 0, 0, 0, 0, 0, 0),
  delta = c(0, 0, 0, 0, 0, 0, 0))

m <- jags.model(
  file = "model.txt",
  data = data,  
  inits = list(inits),
  n.chains = 1
)


update(m,5000)

x <- coda.samples(
    m,
    c("d","sigma"),
    thin = 1, n.iter = 15000
)

summary(x)

plot(x)

結果

> summary(x)

Iterations = 6001:21000
Thinning interval = 1 
Number of chains = 1 
Sample size per chain = 15000 

1. Empirical mean and standard deviation for each variable,
   plus standard error of the mean:

         Mean      SD  Naive SE Time-series SE
d     -0.1400 0.08028 0.0006555       0.001885
sigma  0.1285 0.08188 0.0006686       0.002196

2. Quantiles for each variable:

          2.5%      25%     50%     75%    97.5%
d     -0.31570 -0.18244 -0.1338 -0.0908 0.003366
sigma  0.02855  0.06916  0.1110  0.1672 0.331475

mcmc_5_2

J / RATN

「あともう一回だけ」が超切ない。でも前向き。

ProductName J
RATN
disque corde / ¥ 2,415 (2005-08-20)
在庫あり。

Dummy / Portishead

Wandering Star 聴くとやる気がでてくる。

ProductName Dummy
Portishead
Go! Discs/London / ¥ 1,357 (1994-10-17)
在庫あり。

Grott / DUB ARCHANOID TRIM

スピリチュアルダブ!

ProductName Grott
DUB ARCHANOID TRIM(D.A.T.)
インディーズ・メーカー / ¥ 2,730 (2002-10-25)


これはやばかった。

最近、itunesでレーティング付けてるので、良さげなアルバムを掘り起こしてる。

JAGSでロジスティック回帰モデル

macでwinbugsを頑張って動かすのもアレなので、潔くJAGSで。

マルコフ連鎖モンテカルロ法の5-1

ProductName マルコフ連鎖モンテカルロ法 (統計ライブラリー)
豊田 秀樹
朝倉書店 / ¥ 4,410 ()
在庫あり。

生態学のデータ解析 - R で JAGSによるとjagsの場合にはBUGS言語の文(式)の終わりにはセミコロンが必要。

model
{
  for(i in 1:N){
    r[i] ~ dbin(p[i],n[i]);
    e[i] ~ dnorm(0.0,tau);
    logit(p[i]) <- beta0+beta1*x1[i]+beta21*x21[i]+beta22*x22[i]
    +beta23*x23[i]+beta24*x24[i]+beta25*x25[i]
    +beta26*x26[i]+beta27*x27[i]+e[i];
  }
  beta0 ~ dnorm(0.0,1.0E-6);
  beta1 ~ dnorm(0.0,1.0E-6);
  beta21 ~ dnorm(0.0,1.0E-6);
  beta22 ~ dnorm(0.0,1.0E-6);
  beta23 ~ dnorm(0.0,1.0E-6);
  beta24 <- 0;
  beta25 ~ dnorm(0.0,1.0E-6);
  beta26 ~ dnorm(0.0,1.0E-6);
  beta27 ~ dnorm(0.0,1.0E-6);
  tau ~ dgamma(0.001,0.001);
  sigma<-1/sqrt(tau);
}

Rのコード

library(rjags)

data <- list(
r = c(1,13,44,155,92,100,18,0,0,2,1,0,2,0),
n = c(13,70,115,301,153,141,26,6,16,25,32,8,9,4),
x1 = c(0,0,0,0,0,0,0,1,1,1,1,1,1,1),
x21 = c(1,0,0,0,0,0,0,1,0,0,0,0,0,0),
x22 = c(0,1,0,0,0,0,0,0,1,0,0,0,0,0),
x23 = c(0,0,1,0,0,0,0,0,0,1,0,0,0,0),
x24 = c(0,0,0,1,0,0,0,0,0,0,1,0,0,0),
x25 = c(0,0,0,0,1,0,0,0,0,0,0,1,0,0),
x26 = c(0,0,0,0,0,1,0,0,0,0,0,0,1,0),
x27 = c(0,0,0,0,0,0,1,0,0,0,0,0,0,1),
N = 14)

inits <- list(beta0 = 0, beta1 = 0, beta21 = 0,  beta22 = 0, beta23 = 0, beta24 = 0, beta25 = 0, beta26 = 0,  beta27 = 0, tau = 1)

m <- jags.model(
  file = "model.txt",
  data = data,  
  inits = list(inits),
  n.chains = 1
)

update(m,10000)

x <- coda.samples(
    m,
    c("beta0","beta1","beta21","beta22","beta23","beta24","beta25","beta26","beta27", "sigma"),
    thin = 1, n.iter = 40000
)

summary(x)

plot(x)

サマリーの表示

> summary(x)

Iterations = 11001:51000
Thinning interval = 1 
Number of chains = 1 
Sample size per chain = 40000 

1. Empirical mean and standard deviation for each variable,
   plus standard error of the mean:

           Mean     SD Naive SE Time-series SE
beta0  -0.02576 0.4188 0.002094        0.03425
beta1  -3.11738 0.6059 0.003030        0.01815
beta21 -3.03168 1.4556 0.007278        0.03269
beta22 -1.52877 0.6757 0.003379        0.02993
beta23 -0.28549 0.6355 0.003177        0.04374
beta24  0.00000 0.0000 0.000000        0.00000
beta25  0.35897 0.5926 0.002963        0.02957
beta26  1.05996 0.6389 0.003194        0.04301
beta27  0.70287 0.7677 0.003838        0.03026
sigma   0.34110 0.3980 0.001990        0.03152

2. Quantiles for each variable:

           2.5%      25%      50%      75%   97.5%
beta0  -1.11529 -0.15412  0.02137  0.16865  0.7406
beta1  -4.48731 -3.46881 -3.06188 -2.70570 -2.0882
beta21 -6.39669 -3.83976 -2.86884 -2.05688 -0.6382
beta22 -2.90198 -1.87085 -1.54216 -1.21238 -0.0738
beta23 -1.27301 -0.63488 -0.39022 -0.08097  1.4261
beta24  0.00000  0.00000  0.00000  0.00000  0.0000
beta25 -0.87361  0.10421  0.34604  0.59078  1.7201
beta26  0.08425  0.72545  0.96155  1.25666  2.7809
beta27 -0.80720  0.33750  0.71085  1.10396  2.1284
sigma   0.02816  0.08266  0.19387  0.45606  1.3892

5章は事例がたくさん載っているので、この勢いで一通りやっておこう。

  • 5.1 ロジスティック回帰モデル 
  • 5.2 メタ分析 
  • 5.3 多項ロジットモデル 
  • 5.4 対数線形モデル 
  • 5.5 ポアソン回帰 
  • 5.6 2値データに対する回帰分析 
  • 5.7 トービット回帰モデル 
  • 5.8 変曲点のある回帰分析 
  • 5.9 生存時間分析(ワイブル回帰) 
  • 5.10 生存時間分析(コックス回帰) 
  • 5.11 時系列モデル 
  • 5.12 分散分析 
  • 5.13 分散成分分析 
  • 5.14 分散分析(枝分かれ配置) 
  • 5.15 一般化可能性理論 
  • 5.16 反復測定データの分散分析 
  • 5.17 階層線形モデル 
  • 5.18 項目反応理論(2母数2値モデル) 
  • 5.19 項目反応理論(段階反応モデル) 
  • 5.20 項目反応理論(名義反応モデル) 
  • 5.21 項目反応理論(部分採点・評定尺度モデル)
  • 5.22 項目反応理論(連続反応モデル)
  • 5.23 多次元IRT
  • 5.24 項目反応理論(混合名義反応モデル) 
  • 5.25 項目反応理論における特異項目機能(DIF) の分析 
  • 5.26 正規混合モデル 
  • 5.27 潜在クラス分析 
  • 5.28 成長曲線モデル 
  • 5.29 非線形成長曲線モデル 
  • 5.30 因子分析 
  • 5.31 多母集団分析 
  • 5.32 非線形SEM