# NOT RUN {
# define a variable greta array, and another that is calculated from it
# then calculate what value y would take for different values of x
x <- normal(0, 1, dim = 3)
a <- lognormal(0, 1)
y <- sum(x ^ 2) + a
calculate(y, list(x = c(0.1, 0.2, 0.3), a = 2))
# if the greta array only depends on data,
# you can pass an empty list to values (this is the default)
x <- ones(3, 3)
y <- sum(x)
calculate(y)
# define a model
alpha <- normal(0, 1)
beta <- normal(0, 1)
sigma <- lognormal(1, 0.1)
mu <- alpha + iris$Petal.Length * beta
distribution(iris$Petal.Width) <- normal(mu, sigma)
m <- model(alpha, beta, sigma)
# calculate intermediate greta arrays, given some parameter values
calculate(mu[1:5], list(alpha = 1, beta = 2, sigma = 0.5))
calculate(mu[1:5], list(alpha = -1, beta = 0.2, sigma = 0.5))
# fit the model then calculate samples at a new greta array
draws <- mcmc(m, n_samples = 500)
petal_length_plot <- seq(min(iris$Petal.Length),
max(iris$Petal.Length),
length.out = 100)
mu_plot <- alpha + petal_length_plot * beta
mu_plot_draws <- calculate(mu_plot, draws)
# plot the draws
mu_est <- colMeans(mu_plot_draws[[1]])
plot(mu_est ~ petal_length_plot, type = "n",
ylim = range(mu_plot_draws[[1]]))
apply(mu_plot_draws[[1]], 1, lines,
x = petal_length_plot, col = grey(0.8))
lines(mu_est ~ petal_length_plot, lwd = 2)
# }
# NOT RUN {
# }
Run the code above in your browser using DataLab