December 14, 2014

Rstanでブログが検索結果に出ない話

で表示されないけど、 でなら表示される。 Rstan単独で検索してブログが表示される人もいるので、検索エンジンに嫌われてる可能性が高くて涙目です。 ちなみに、「stan ワイブル」で検索すれば、上位にヒットしました。

December 10, 2014

Rstanで重回帰分析

ワイブルばかりだったので、メジャーな重回帰分析をしてみた。Rstanにハマりつつある。 linear_regression.stan data { int<lower=1> J; // number of data int<lower=1> K; // number of covariate real y[J]; matrix[J,K] x; } parameters { real<lower=0> sd0; vector[K] beta; real beta0; } model { for(j in 1:J){ increment_log_prob(normal_log(y[j], x[j] * beta + beta0, sd0)); } } linear_regression.R library(rstan) n <- 100 x1 <- rnorm(n, mean=1, sd=0.1) x2 <- rnorm(n, mean=1, sd=0.1) beta1 <- 1 beta2 <- 2 beta0 <- 3 y <- beta1 * x1 + beta2 * x2 + beta0 data. Read more

December 4, 2014

Rstanでワイブル回帰の推定(説明変数複数バージョン)

matrixやvectorの使い方を覚えた。 weibull_regression.stan data { int<lower=1> J; // number of data int<lower=1> K; // number of covariate int status[J]; real<lower=0> tt[J]; matrix[J,K] x; } parameters { real<lower=1, upper=3> m; vector[K] beta; real beta0; } transformed parameters { real<lower=0> eta[J]; for (j in 1:J) eta[j] <- exp(x[j]*beta + beta0); } model { for(j in 1:J){ if(status[j] == 1){ increment_log_prob(weibull_log(tt[j], m, eta[j])); } if(status[j] == 0){ increment_log_prob(weibull_ccdf_log(tt[j], m, eta[j])); } } } weibull_regression. Read more

December 4, 2014

Rstanでワイブル回帰の推定

weibull_regression.stan data { int<lower=1> J; // number of data int status[J]; real<lower=0> tt[J]; real x[J]; } parameters { real<lower=1, upper=3> m; real beta; real beta0; } transformed parameters { real<lower=0> eta[J]; for (j in 1:J) eta[j] <- exp(beta * x[j] + beta0); } model { for(j in 1:J){ if(status[j] == 1){ increment_log_prob(weibull_log(tt[j], m, eta[j])); } if(status[j] == 0){ increment_log_prob(weibull_ccdf_log(tt[j], m, eta[j])); } } } mを1から3ぐらいにしておかないと、推定値が遥か彼方に飛ばされることがある。ワイブル分布は厄介やね。 weibull_regression.R library(rstan) set.seed(0) n <- 1000 # 台数 m <- 2 # 形状パラメータ beta <- -1 # 係数 beta0 <- 10 # 切片 x <- rnorm(n) data. Read more

December 4, 2014

Rstanでワイブル分布のパラメータ推定(打切りのあるデータで)

weibull.stan data { int<lower=1> J; // number of data int status[J]; real<lower=0> tt[J]; } parameters { real<lower=1> m; real<lower=1> eta; } model { for(j in 1:J){ if(status[j] == 1){ increment_log_prob(weibull_log(tt[j], m, eta)); } if(status[j] == 0){ increment_log_prob(weibull_ccdf_log(tt[j], m, eta)); } } } weibull.R サンプルデータの生成 > set.seed(0) > data.df <- data.frame(tt = rweibull(1000, shape = 2, scale = 100), + status = rep(1,1000)) > data.df[data.df$tt > 30,]$status <- 0 > data. Read more

December 3, 2014

Rstanでワイブル分布のパラメータ推定

Rstanでワイブル分布のパラメータ推定をしてみる。打切りを考慮するにはどうしたら良いのだろうか? Stan Modeling Language User’s Guide and Reference Manual, v2.5.0の83ページ辺りを読めば出来そう。 weibull.stan data { int<lower=0> J; // number of data real<lower=0> tt[J]; // time } parameters { real<lower=0> m; real<lower=0> eta; } model { for(j in 1:J) tt[j] ~ weibull(m, eta); } weibull.R library(rstan) data.ls <- list(J = 1000, tt = rweibull(1000, shape = 2, scale = 100)) fit <- stan(file = 'weibull.stan', data = data.ls, iter = 1000, chains = 4) 結果 > fit Inference for Stan model: weibull. Read more

© gepuro 2013