if (FALSE) {
#----------------------------------------------------------------------------
# Blimp Example 4.3: Linear Regression
#..........
# Trace Plots
# Example 1a: Default setting, specifying name of the folder
blimp.plot("Posterior_Ex4.3")
# Example 1b: Default setting, specifying the posterior file
blimp.plot("Posterior_Ex4.3/posterior.csv")
# Example 1c: Print parameters 2, 3, 4, and 5
blimp.plot("Posterior_Ex4.3", param = 2:5)
# Example 1e: Arrange panels in three columns
blimp.plot("Posterior_Ex4.3", ncol = 3)
# Example 1f: Specify "Pastel 1" palette for the hcl.colors function
blimp.plot("Posterior_Ex4.3", palette = "Pastel 1")
#..........
# Posterior Distribution Plots
# Example 2a: Default setting, i.e., posterior median and equal-tailed interval
blimp.plot("Posterior_Ex4.3", plot = "post")
# Example 2b: Display posterior mean and maximum a posteriori
blimp.plot("Posterior_Ex4.3", plot = "post", point = c("m", "map"))
# Example 2c: Display maximum a posteriori and highest density interval
blimp.plot("Posterior_Ex4.3", plot = "post", point = "map", ci = "hdi")
# Example 2d: Do not display any point estimates and credible interval
blimp.plot("Posterior_Ex4.3", plot = "post", point = "none", ci = "none")
# Example 2d: Do not display histograms
blimp.plot("Posterior_Ex4.3", plot = "post", hist = FALSE)
#..........
# Save Plots
# Example 3a: Save all plots in pdf format
blimp.plot("Posterior_Ex4.3", saveplot = "all")
# Example 3b: Save all plots in png format with 300 dpi
blimp.plot("Posterior_Ex4.3", saveplot = "all", file = "Blimp_Plot.png", dpi = 300)
# Example 3a: Save posterior distribution plot, specify width and height of the plot
blimp.plot("Posterior_Ex4.3", plot = "none", saveplot = "post",
width = 7.5, height = 7)
#----------------------------------------------------------------------------
# Plot from misty.object
# Create misty.object
object <- blimp.plot("Posterior_Ex4.3", plot = "none")
# Trace plot
blimp.plot(object, plot = "trace")
# Posterior distribution plot
blimp.plot(object, plot = "post")
#----------------------------------------------------------------------------
# Create Plots Manually
# Load ggplot2 package
library(ggplot2)
# Create misty object
object <- blimp.plot("Posterior_Ex4.3", plot = "none")
#..........
# Example 4: Trace Plots
# Extract data
data.trace <- object$data$trace
# Plot
ggplot(data.trace, aes(x = iter, y = value, color = chain)) +
annotate("rect", xmin = 0, xmax = 1000, ymin = -Inf, ymax = Inf,
alpha = 0.4, fill = "gray85") +
geom_line() +
facet_wrap(~ param, ncol = 2, scales = "free") +
scale_x_continuous(name = "", expand = c(0.02, 0)) +
scale_y_continuous(name = "", expand = c(0.02, 0)) +
scale_colour_manual(name = "Chain",
values = hcl.colors(n = 2, palette = "Set 2")) +
theme_bw() +
guides(color = guide_legend(nrow = 1, byrow = TRUE)) +
theme(plot.margin = margin(c(4, 15, -10, 0)),
legend.position = "bottom",
legend.title = element_text(size = 10),
legend.text = element_text(size = 10),
legend.box.margin = margin(c(-16, 6, 6, 6)),
legend.background = element_rect(fill = "transparent"))
#..........
# Example 5: Posterior Distribution Plots
# Extract data
data.post <- object$data$post
# Plot
ggplot(data.post, aes(x = value)) +
geom_histogram(aes(y = after_stat(density)), color = "black", alpha = 0.4,
fill = "gray85") +
geom_density(color = "#0072B2") +
geom_vline(data = data.frame(param = levels(data.post$param),
stat = tapply(data.post$value, data.post$param, median)),
aes(xintercept = stat, color = "Median"), linewidth = 0.6) +
geom_vline(data = data.frame(param = levels(data.post$param),
low = tapply(data.post$value, data.post$param,
function(y) quantile(y, probs = 0.025))),
aes(xintercept = low), linetype = "dashed", linewidth = 0.6) +
geom_vline(data = data.frame(param = levels(data.post$param),
upp = tapply(data.post$value, data.post$param,
function(y) quantile(y, probs = 0.975))),
aes(xintercept = upp), linetype = "dashed", linewidth = 0.6) +
facet_wrap(~ param, ncol = 2, scales = "free") +
scale_x_continuous(name = "", expand = c(0.02, 0)) +
scale_y_continuous(name = "Probability Density, f(x)",
expand = expansion(mult = c(0L, 0.05))) +
scale_color_manual(name = "Point Estimate", values = c(Median = "#D55E00")) +
labs(caption = "95% Equal-Tailed Interval") +
theme_bw() +
theme(plot.margin = margin(c(4, 15, -8, 4)),
plot.caption = element_text(hjust = 0.5, vjust = 7),
legend.position = "bottom",
legend.title = element_text(size = 10),
legend.text = element_text(size = 10),
legend.box.margin = margin(c(-30, 6, 6, 6)),
legend.background = element_rect(fill = "transparent"))
}
Run the code above in your browser using DataLab