100日後にベイズ統計ができるゆーみん

とある大学院でベイズをマスターしようと目論む博士課程

100日後にベイズ統計ができるゆーみん

とある大学院でベイズをマスターしようと目論む博士課程

最近の記事

【2回目】100日後にベイズ統計ができるゆーみん【28日目】

こんにちは! 今日はポアソン回帰の続きをやっていきます。 5.4.3 Stanで実装 ↓モデル式の実装例 λ[n] = exp(b1+b2A[n]+b3Score[n]) n=1,...,N M[n] ~ Poisson(λ[n]) n=1,...,N data { int N; int<lower=0, upper=1> A[N]; real<lower=0, upper=1> Score[N]; int<lower=0> M[

    • 【2回目】100日後にベイズ統計ができるゆーみん【27日目】

      おはようございます! もはや朝活になりつつある100日ベイズです。 次はポアソン回帰です! データは二項ロジスティック回帰で使ったデータと同じです。 左から、学生のID、バイトが好きかどうか、試験の成績、総授業回数、実際に授業に出席した回数です。 今回は、総授業回数Mが応答変数となり、今回の解析の目的は 「説明変数AとScoreが総授業回数Mにどれほど影響しているか知りたい」です。 総授業回数のような離散的で負の値を取らないカウントデータに対しては、ポアソン分布を用

      • 【2回目】100日後にベイズ統計ができるゆーみん【25日目】

        おはようございます! 今日も進めて参ります。 5.3.7 図によるモデルのチェック ・確率と実測値のプロット 出席確率の高さが実際の出席に繋がっているのかを確認しましょう! 横軸にq[i]の事後分布の中央値、縦軸にY[i]をとる散布図を描きます。 library(ggplot2)set.seed(123)load('output/result-model5-5.RData')ms <- rstan::extract(fit)logistic <- function(x

        • 【2回目】100日後にベイズ統計ができるゆーみん【24日目】

          おはようございます! 今日も少しでも進めていきます。 5.3.5 モデル式の記述 q[i]=inv_logit(b1 + b2A[i] + b3Score[i] + b4Wether[i]) i=1,...,I Y[i] ~ Bernoulli (q[i]) i=1,...,I ここでは添字が変更されています。それは添字が学生ごとに対応していたのに対して出欠ごとに対応しているからです。 このモデル式では出欠データごとに独立に確率分布が生成されているため

        【2回目】100日後にベイズ統計ができるゆーみん【28日目】

          【2回目】100日後にベイズ統計ができるゆーみん【23日目】

          今日は5章ロジスティック回帰です! 今回のデータはいままでのデータに天気の情報が入ります. Wether:天気(A:晴れ,B:曇り,C:雨) 5.3.1 解析の目的 ここでは3つの説明変数(バイトが好きかどうか,学問への興味,天気)で出席確率がどれほど予測できるか それぞれの説明変数が出席確率にどれほど影響しているのかを知りたい の2つを目的とする. 今回は出席確率を考えるので,ロジスティック回帰を行う. 5.3.2 データ分布の確認 このデータの場合,A列

          【2回目】100日後にベイズ統計ができるゆーみん【23日目】

          【2回目】100日後にベイズ統計ができるゆーみん【22日目】

          おはようございます! 昨日は少し体調が悪く、おやすみしてしまいましたが 今日は頑張って進めていきますよ〜 一昨日はStanで実装までしました。 5.2.6 推定結果の解釈 得られた事後平均の値をモデル式に代入してみます。 q[n] = inv_logit(0.09-0.62A[n] + 1.90Score[n] / 200) n=1,...,N Y[n] ~ Binomial(M[n], q[n]) n=1,...,N ここでinv_lo

          【2回目】100日後にベイズ統計ができるゆーみん【22日目】

          【2回目】100日後にベイズ統計ができるゆーみん【20日目】

          おはようございます! 今日も少しずつ進めていきます。 5.2.5 Stanで実装 data { int N; int<lower=0, upper=1> A[N]; real<lower=0, upper=1> Score[N]; int<lower=0> M[N]; int<lower=0> Y[N];}parameters { real b1; real b2; real b3;}transformed parameters { real q[N]; for (n

          【2回目】100日後にベイズ統計ができるゆーみん【20日目】

          【2回目】100日後にベイズ統計ができるゆーみん【19日目】

          おはようございます! 今日もちょっとだけ進めていきますよ〜 5.2.3 メカニズムの想像 今回のデータは重回帰だとちょっとまずいのです。 なぜならば出席率の実測値がほぼ1であるデータがいくつかあるため、説明変数の組み合わせによっては出席率の予測値が1を超えてしまいます。 この場合は説明変数の線形結合(b1+b2A+b3Score)を、(0,1) の範囲に変換します。方法はロジスティック関数1/{1+exp(-x)}を使って変換します。Stanではinv_logit関数

          【2回目】100日後にベイズ統計ができるゆーみん【19日目】

          【2回目】100日後にベイズ統計ができるゆーみん【18日目】

          おはようございます! 一日遅れで朝やるのが癖になりつつあります。 やらないよりはいいかな() 次は、二項ロジスティック回帰です! 今回はこのようなデータを扱う。 PersonID:学生のI D A:バイトが好きか、そうでないか Score:学問への興味の強さを数値化したもの M:3ヶ月における履修登録した科目の総授業回数 Y:そのうち実際に出席した回数 5.2.1 解析の目的 解析の目的は前回と同じで 「2つの説明変数AとScoreで応答変数がどれほど予測できる

          【2回目】100日後にベイズ統計ができるゆーみん【18日目】

          【2回目】100日後にベイズ統計ができるゆーみん【17日目】

          おはようございます... 色々あってボロボロですが、続けていきたいと思います。 重回帰の続きです。 5.1.8 図によるモデルのチェック ・実測値と予測値のプロット この図からは多くの点がy=xの直線の近くにあって、80%区間がy=xの直線を含んでいるので、2つの説明変数から応答変数を十分に予測できそうだと判断した。 またデータ点全体のばらつきに系統だった偏りがないので今回の解析には問題ないと判断する。 ・推定されたノイズの分布 その他の影響(ノイズ)ε[n]はNo

          【2回目】100日後にベイズ統計ができるゆーみん【17日目】

          【2回目】100日後にベイズ統計ができるゆーみん【16日目】

          おはようございます! もういっつも寝落ちしていて、睡眠は万全なゆーみんです。 こんな自分も愛していきたいと思います(反省しろ)。 では、重回帰の続きをやっていきます。 昨日はStanコードを書いたので、Rコードで実行していきます。 library(rstan)d <- read.csv(file='input/data-attendance-1.txt')data <- list(N=nrow(d), A=d$A, Score=d$Score/200, Y=d$Y)fit

          【2回目】100日後にベイズ統計ができるゆーみん【16日目】

          【2回目】100日後にベイズ統計ができるゆーみん【15日目】

          こんばんは! 重回帰の続きをやっていきます〜 5.1.1 解析の目的 目的は 「2つの説明変数AとScoreで応答変数Yがどれほど予測できるのか知りたい」 「それぞれの説明変数が出席率にどれほど影響しているのかを知りたい」 です。 5.1.2 データの分布の確認 データが与えられたら 分布 を確認していきます。 載ってあったコードをちょっと噛み砕いた(私が分かるように)ものも併記しておきます。 A:アルバイトが好きかどうかを表す。0が好きでないで、1が好き。 Sco

          【2回目】100日後にベイズ統計ができるゆーみん【15日目】

          【2回目】100日後にベイズ統計ができるゆーみん【14日目】

          おはようございます... また寝落ちていたので朝に勉強しています、ゆーみんです。 情けないです。 本当にちょこっとだけ進めます。 4章を無事終えたので、次は5章の基本的な回帰とモデルのチェックに進みます。 ここでは統計モデリングでよく使用される重回帰、二項ロジスティック回帰、ロジスティック回帰、ポアソン分布が紹介されています。 5.1 重回帰 4章では年齢という一つの説明変数から年収を予測していました。 今回は復習の説明変数があるパターンを学びます。 A:アルバイト

          【2回目】100日後にベイズ統計ができるゆーみん【14日目】

          【2回目】100日後にベイズ統計ができるゆーみん【13日目】

          こんばんは! 今日もへろへろのゆーみんです。 今日は復習回です。 はじめからデータの読み込みを行い、 Stanを実行し、 MCMCサンプルを取り出して要約統計量を算出し、 基本年収のベイズ信頼区間と年収のベイズ予測区間のグラフを描いていきます! 全体的なコードは、サポートページに上がっています↓ library(rstan)library(ggplot2)d <- read.csv(file='input/data-salary.txt')X_new <- 23:60d

          【2回目】100日後にベイズ統計ができるゆーみん【13日目】

          【2回目】100日後にベイズ統計ができるゆーみん【12日目】

          こんばんは!今日は早めに帰れたので今日中にアップできそうです。 4.4.12 transformed parameters ブロックとgenerated quantitiesブロック 処理はR側でするのではなく、 できるだけStan側で実行するのが望ましいそうな。 特にtransformed parameters ブロックとgenerated quantitiesブロックが利用できるそうです。 data { int N; real X[N]; real Y[N];}pa

          【2回目】100日後にベイズ統計ができるゆーみん【12日目】

          【2回目】100日後にベイズ統計ができるゆーみん【11日目】

          はい!ゆーみんです! 無事、一日飛ばしました() やらかしました。 そして今も朝というやらかしぶりですが、 なんとか進めて参ります。 おそらく7月いっぱいはこんな感じです(諦めるな) 4.4.10 並列計算の実行方法 MCMCのアルゴリズムはchainごとに独立しているので、 各chainを並列に計算しておいてまとめることができる。 並列化するとなんやかんやでデバックが難しくなるので、 並列化なしでうまくいくぜ〜ってことを確かめた上でやるのがいいらしい。 ちなみにや

          【2回目】100日後にベイズ統計ができるゆーみん【11日目】