# Read in the data, create a data frame called “Crabs”
rm(list=ls()); set.seed(20210405);
filename <- "http://www.stat.ufl.edu/~aa/cat/data/Crabs.dat"
Crabs <- read.table(filename, header=T)
dim(Crabs); names(Crabs);
# Here we consider a logistic regression with response
# variable y = 1 if satellite(s), 0 if none
y <- Crabs$y; x <- Crabs$width; rm(Crabs); n <- length(y);
# Fit Bayesian logistic regression of y = "presence of
# satellite" on x = width
library(rstan)
stan_model <- "
data{
int
int
vector [n] x;
}
parameters{
real alpha;
real beta;
}
transformed parameters{
real EL_50;
EL_50 = -alpha / beta;
}
model{
# alpha ~ normal(0, 100);
# beta ~ normal(0, 100);
y ~ bernoulli_logit(alpha + beta*x) ;
}
”
data <- list(n=n, y=y, x=x) fit <- stan(model_code=stan_model, data=data, iter=2000, chains=5) sims <- extract(fit); names(sims); Results <- summary(fit)$summary; Results <- Results[1:3, ]; round(Results, 2);