Subscribed unsubscribe Subscribe Subscribe

CPSAM.org

computer, programming, statistics and more!

牛丼問題の簡素化

require(rstan)

oneGroup <- "
  data {
    int n;
    vector[n] dat;
  }

  parameters {
    real mu;
    real<lower=0> sigma;
  }

  model {
    dat ~ normal(mu, sigma);
  }
"

set.seed(1)
x <- c(76.5,83.9,87.9,70.8,84.6,85.1,79.6,79.8,79.7,78.0)
don <- list(dat=x, n=length(x))

oneGroupMean <- stan(
  model_code = oneGroup,
  data = don,
  iter = 51000,
  warmup = 1000,
  thin = 1,
  chains = 4
)

]

summary(oneGroupMean)$summary

m <- summary(oneGroupMean)$summary[1,1]
s <- summary(oneGroupMean)$summary[2,1]
x.seq <- seq(60,100,1)

plot(x.seq, dnorm(x.seq, m, s))
par(new=T)
hist(x, xlim=c(60,100))
Remove all ads