Chapter 4 Practice Problems - Alissa Rivero
2024-09-19
Chapter 4 Questions
______________________________________________________________________________________
Easy Question 1:
The likelihood is yi ∼ Normal(μ, σ)
Easy Question 2:
2 parameters.
Easy Question 3
Easy Question 4
μi = α + βxi
Easy Question 5:
3 parameters
______________________________________________________________________________________
Medium Question 1
samp_mu <- rnorm(10000, 0, 10)
samp_sig <- runif(10000, 0, 10)
pri_y <- rnorm(10000, samp_mu, samp_sig)
dens(pri_y)
Medium Question 2
form <- alist(
y~dnorm(mu, sigma),
mu~dnorm(0, 10),
sigma~dunif(0,10))Medium Question 4
ht ∼ Normal(μ, σ),
μi = α + β * xi,
α ∼ Normal(150, 25),
β ∼ Normal(3,1),
σ ∼ Uniform(0,50),
α ∼ Normal(150, 20)
set.seed(42)
dens(rnorm(10000, 150, 25))
N <- 1000
a <- rnorm(10000, 150, 25)
b <- rnorm(N, 3, 1)
sigma <- runif(N, 0, 50)
# Checking on height of dist & plot
height_0 <- rnorm(N, a + b * 0, sigma)
dens(height_0)
Medium Question 8
library(splines)
data(cherry_blossoms)
d <- cherry_blossoms
d2 <- d[complete.cases(d$temp),]
numknots<-15
knots<-quantile(d2$year,probs=seq(0,1,length.out=numknots))
B <- bs(d2$year, knots=knots[-c(1,numknots)], degree=3, intercept=TRUE)
m4.7 <- quap(alist(T ~ dnorm(mu, sigma), mu <- a + B%*%w, a ~ dnorm(6, 10), w ~ dnorm(0, 1), sigma ~ dexp(1)), data=list(T=d2$temp, B=B), start=list(w=rep(0, ncol(B))))
post <- extract.samples(m4.7)
w <- apply(post$w, 2, mean)
plot(NULL, xlim=range(d2$year), ylim=c(-2,2), xlab="year", ylab="basis * weight")
for (i in 1:ncol(B)) lines(d2$year, w[i]*B[,i])
mu <- link(m4.7)
mu_PI <- apply(mu, 2, PI, 0.97)
plot(d2$year, d2$temp, col=col.alpha(rangi2,0.3), pch=16)
shade(mu_PI, d2$year, col=col.alpha("black",0.5))
numknots<-30
knots<-quantile(d2$year,probs=seq(0,1,length.out=numknots))
B <- bs(d2$year, knots=knots[-c(1,numknots)], degree=3, intercept=TRUE)
m4.7 <- quap(
alist(T ~ dnorm(mu, sigma), mu <- a + B%*%w, a ~ dnorm(6, 10),
w ~ dnorm(0, 1), sigma ~ dexp(1)), data=list(T=d2$temp, B=B),
start=list(w=rep(0, ncol(B)))
)
post <- extract.samples(m4.7)
w <- apply(post$w, 2, mean)
plot(NULL,
xlim=range(d2$year),
ylim=c(-2,2),
xlab="year",
ylab="basis * weight")
for (i in 1:ncol(B)) lines(d2$year, w[i]*B[,i])
mu <- link(m4.7)
mu_PI <- apply(mu, 2, PI, 0.97)
plot(d2$year, d2$temp, col=col.alpha(rangi2,0.3), pch=16)
shade(mu_PI, d2$year, col=col.alpha("black",0.5))
numknots<-30
knots<-quantile(d2$year,probs=seq(0,1,length.out=numknots))
B <- bs(d2$year, knots=knots[-c(1,numknots)], degree=3, intercept=TRUE)
m4.7 <- quap(alist(
T ~ dnorm(mu, sigma),
mu <- a + B%*%w,
a~dnorm(6,10),
w ~ dnorm(0, 10),
sigma ~ dexp(1)
),
data=list(T=d2$temp, B=B),
start=list(w=rep(0, ncol(B)))
)
post <- extract.samples(m4.7)
w <- apply(post$w, 2, mean)
plot(NULL, xlim=range(d2$year), ylim=c(-2,2), xlab="year", ylab="basis * weight")
for (i in 1:ncol(B))
lines(d2$year, w[i]*B[,i])
Hard Question 1
library(tibble)
library(dplyr)
library(tidyr)
library(glue)
library(tidybayes)## Warning: package 'tidybayes' was built under R version 4.4.1##
## Attaching package: 'tidybayes'## The following object is masked from 'package:tune':
##
## parameters## The following object is masked from 'package:dials':
##
## parameterslibrary(knitr)
library(kableExtra)##
## Attaching package: 'kableExtra'## The following object is masked from 'package:dplyr':
##
## group_rowslibrary(brms)## Warning: package 'brms' was built under R version 4.4.1## Loading required package: Rcpp##
## Attaching package: 'Rcpp'## The following object is masked from 'package:rsample':
##
## populate## Loading 'brms' package (version 2.22.0). Useful instructions
## can be found by typing help('brms'). A more detailed introduction
## to the package is available through vignette('brms_overview').##
## Attaching package: 'brms'## The following objects are masked from 'package:tidybayes':
##
## dstudent_t, pstudent_t, qstudent_t, rstudent_t## The following object is masked from 'package:lme4':
##
## ngrps## The following object is masked from 'package:psych':
##
## cs## The following object is masked from 'package:dials':
##
## mixture## The following objects are masked from 'package:rethinking':
##
## LOO, stancode, WAIC## The following object is masked from 'package:stats':
##
## arlibrary(here)## here() starts at /Users/alissarvr/Evolutionary Statisticsdata(Howell1)
how_dat <- Howell1 %>%
filter(age >= 18) %>%
mutate(weight_c = weight - mean(weight))
medium <- quap(alist(height ~ dnorm(mu, sigma),
mu <- a + b * (weight_c),
a ~ dnorm(178, 20),
b ~ dlnorm(0, 1),
sigma ~ dunif(0, 50)),
data = how_dat)
round(vcov(medium), 3)## a b sigma
## a 0.073 0.000 0.000
## b 0.000 0.002 0.000
## sigma 0.000 0.000 0.037model <- brm(height ~ 1 + weight_c, data = how_dat, family = gaussian,
prior = c(prior(normal(178, 20), class = Intercept),
prior(lognormal(0, 1), class = b, lb = 0),
prior(uniform(0, 50), class = sigma)),
iter = 28000, warmup = 27000, chains = 4, cores = 4, seed = 1234)## Warning: It appears as if you have specified an upper bounded prior on a parameter that has no natural upper bound.
## If this is really what you want, please specify argument 'ub' of 'set_prior' appropriately.
## Warning occurred for prior
## <lower=0> sigma ~ uniform(0, 50)## Compiling Stan program...## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'## Trying to compile a simple C file## Running /Library/Frameworks/R.framework/Resources/bin/R CMD SHLIB foo.c
## using C compiler: ‘Apple clang version 16.0.0 (clang-1600.0.26.3)’
## using SDK: ‘’
## clang -arch arm64 -I"/Library/Frameworks/R.framework/Resources/include" -DNDEBUG -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/Rcpp/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/unsupported" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/BH/include" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/src/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppParallel/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/rstan/include" -DEIGEN_NO_DEBUG -DBOOST_DISABLE_ASSERTS -DBOOST_PENDING_INTEGER_LOG2_HPP -DSTAN_THREADS -DUSE_STANC3 -DSTRICT_R_HEADERS -DBOOST_PHOENIX_NO_VARIADIC_EXPRESSION -D_HAS_AUTO_PTR_ETC=0 -include '/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp' -D_REENTRANT -DRCPP_PARALLEL_USE_TBB=1 -I/opt/R/arm64/include -fPIC -falign-functions=64 -Wall -g -O2 -c foo.c -o foo.o
## In file included from <built-in>:1:
## In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp:22:
## In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Dense:1:
## In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Core:19:
## /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/src/Core/util/Macros.h:679:10: fatal error: 'cmath' file not found
## 679 | #include <cmath>
## | ^~~~~~~
## 1 error generated.
## make: *** [foo.o] Error 1## Start sampling## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'tibble(individual = 1:5,
weight = c(46.95, 43.72, 64.78, 32.59, 54.63)) |>
mutate(weight_c = weight - mean(how_dat$weight)) |>
add_predicted_draws(model) |>
group_by(individual, weight) |>
mean_qi(.prediction, .width = 0.89) |>
mutate(range = glue("[{sprintf('%0.1f', .lower)}--{sprintf('%0.1f', .upper)}]"),
.prediction = sprintf("%0.1f", .prediction)) |>
select(individual, weight, exp = .prediction, range) |>
kbl(align = "c", booktabs = TRUE,
col.names = c("Individual", "Weight", "Expected Height", "89% Interval"))| Individual | Weight | Expected Height | 89% Interval |
|---|---|---|---|
| 1 | 46.95 | 156.4 | [148.4–164.4] |
| 2 | 43.72 | 153.6 | [145.8–161.7] |
| 3 | 64.78 | 172.7 | [164.7–180.8] |
| 4 | 32.59 | 143.3 | [135.0–151.5] |
| 5 | 54.63 | 163.3 | [154.9–171.5] |
Hard Question 8
cb_dat <- cherry_blossoms |>
drop_na(doy)
knots_15 <- quantile(cb_dat$year, probs = seq(0, 1, length.out = 15))
B_15 <- bs(cb_dat$year, knots = knots_15[-c(1, 15)], degree = 3, intercept = TRUE)
cb_dat_15 <- cb_dat |>
mutate(B = B_15)
b4h8 <- brm(doy ~ 0 + B, data = cb_dat_15, family = gaussian,
prior = c(prior(normal(0, 10), class = b),
prior(exponential(1), class = sigma)),
iter = 2000, warmup = 1000, chains = 4, cores = 4, seed = 1234)## Compiling Stan program...## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'## Trying to compile a simple C file## Running /Library/Frameworks/R.framework/Resources/bin/R CMD SHLIB foo.c
## using C compiler: ‘Apple clang version 16.0.0 (clang-1600.0.26.3)’
## using SDK: ‘’
## clang -arch arm64 -I"/Library/Frameworks/R.framework/Resources/include" -DNDEBUG -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/Rcpp/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/unsupported" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/BH/include" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/src/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppParallel/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/rstan/include" -DEIGEN_NO_DEBUG -DBOOST_DISABLE_ASSERTS -DBOOST_PENDING_INTEGER_LOG2_HPP -DSTAN_THREADS -DUSE_STANC3 -DSTRICT_R_HEADERS -DBOOST_PHOENIX_NO_VARIADIC_EXPRESSION -D_HAS_AUTO_PTR_ETC=0 -include '/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp' -D_REENTRANT -DRCPP_PARALLEL_USE_TBB=1 -I/opt/R/arm64/include -fPIC -falign-functions=64 -Wall -g -O2 -c foo.c -o foo.o
## In file included from <built-in>:1:
## In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp:22:
## In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Dense:1:
## In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Core:19:
## /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/src/Core/util/Macros.h:679:10: fatal error: 'cmath' file not found
## 679 | #include <cmath>
## | ^~~~~~~
## 1 error generated.
## make: *** [foo.o] Error 1## Start sampling## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'epred_draws(b4h8, cb_dat_15) |>
summarize(mean_qi(.epred, .width = 0.89), .groups = "drop") |>
ggplot(aes(x = year, y = doy)) +
geom_point(alpha = 0.2) +
geom_hline(yintercept = mean(cb_dat$doy), linetype = "dashed") +
geom_lineribbon(aes(y = y, ymin = ymin, ymax = ymax),
alpha = 0.8, fill = "#009FB7") +
labs(x = "Year", y = "Day in Year")
b4h8_2 <- brm(doy ~ 0 + B, data = cb_dat_15, family = gaussian,
prior = c(prior(normal(100, 10), class = b),
prior(exponential(1), class = sigma)),
iter = 2000, warmup = 1000, chains = 4, cores = 4, seed = 1234)## Compiling Stan program...## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'## Trying to compile a simple C file## Running /Library/Frameworks/R.framework/Resources/bin/R CMD SHLIB foo.c
## using C compiler: ‘Apple clang version 16.0.0 (clang-1600.0.26.3)’
## using SDK: ‘’
## clang -arch arm64 -I"/Library/Frameworks/R.framework/Resources/include" -DNDEBUG -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/Rcpp/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/unsupported" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/BH/include" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/src/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppParallel/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/rstan/include" -DEIGEN_NO_DEBUG -DBOOST_DISABLE_ASSERTS -DBOOST_PENDING_INTEGER_LOG2_HPP -DSTAN_THREADS -DUSE_STANC3 -DSTRICT_R_HEADERS -DBOOST_PHOENIX_NO_VARIADIC_EXPRESSION -D_HAS_AUTO_PTR_ETC=0 -include '/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp' -D_REENTRANT -DRCPP_PARALLEL_USE_TBB=1 -I/opt/R/arm64/include -fPIC -falign-functions=64 -Wall -g -O2 -c foo.c -o foo.o
## In file included from <built-in>:1:
## In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp:22:
## In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Dense:1:
## In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Core:19:
## /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/src/Core/util/Macros.h:679:10: fatal error: 'cmath' file not found
## 679 | #include <cmath>
## | ^~~~~~~
## 1 error generated.
## make: *** [foo.o] Error 1## Start sampling## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'## Found more than one class "stanfit" in cache; using the first, from namespace 'rethinking'## Also defined by 'rstan'epred_draws(b4h8_2, cb_dat_15) |>
summarize(mean_qi(.epred, .width = 0.89), .groups = "drop") |>
ggplot(aes(x = year, y = doy)) +
geom_point(alpha = 0.2) +
geom_hline(yintercept = mean(cb_dat$doy), linetype = "dashed") +
geom_lineribbon(aes(y = y, ymin = ymin, ymax = ymax),
alpha = 0.8, fill = "#009FB7") +
labs(x = "Year", y = "Day in Year")