2-Dimensional Model with Clusters

library(bayesplot)
library(coda)
library(ggplot2)
library(ggrepel)
library(tidyverse)

Load 2-dimensional DP prior model MCMC object

load("Out/mcmc_list2dDP.Rdata")

thetas <- grep("^theta", colnames(mcmc_list2dDP[[1]]), value = TRUE)
gammas <- grep("^alpha", colnames(mcmc_list2dDP[[1]]), value = TRUE)
clusters <- grep("^cluster", colnames(mcmc_list2dDP[[1]]), value = TRUE)

Summary of MCMC chain output

View a summary of model developed here: https://rpubs.com/steveneheil/pwc-MCMC-CDAT

load("Out/summary2dDP.Rdata")
summary2dDP
##               Mean    SD Naive SE Time-series SE   2.5%    25%    50%    75%
## theta1.1121  1.194 0.406    0.014          0.016  0.450  0.912  1.178  1.469
## theta1.1154 -0.012 0.393    0.013          0.014 -0.733 -0.284 -0.029  0.246
## theta1.1222  0.636 0.353    0.012          0.014 -0.036  0.387  0.638  0.874
## theta1.1637  0.608 0.357    0.012          0.013 -0.059  0.370  0.599  0.835
## theta1.1995  0.753 0.414    0.014          0.015  0.040  0.461  0.715  1.007
## theta1.3369  1.579 0.457    0.015          0.016  0.713  1.267  1.570  1.877
## theta1.4004 -0.769 0.430    0.014          0.016 -1.613 -1.058 -0.749 -0.479
## theta1.4209 -1.285 0.419    0.014          0.014 -2.168 -1.560 -1.255 -1.011
## theta1.4268 -0.770 0.440    0.015          0.016 -1.737 -1.047 -0.747 -0.471
## theta1.4841 -0.884 0.368    0.012          0.013 -1.653 -1.125 -0.866 -0.637
## theta1.5074  1.156 0.400    0.013          0.015  0.386  0.881  1.126  1.435
## theta1.5495 -2.339 0.466    0.016          0.016 -3.262 -2.654 -2.323 -2.018
## theta1.5710 -0.864 0.396    0.013          0.013 -1.664 -1.123 -0.863 -0.593
## theta1.8294 -0.184 0.384    0.013          0.014 -0.927 -0.449 -0.165  0.082
## theta1.8980  0.651 0.374    0.012          0.012 -0.120  0.401  0.651  0.914
## theta1.9899  0.471 0.380    0.013          0.013 -0.290  0.222  0.476  0.733
## theta2.1121 -1.101 0.361    0.012          0.014 -1.849 -1.326 -1.090 -0.869
## theta2.1154 -1.554 0.348    0.012          0.012 -2.217 -1.790 -1.557 -1.334
## theta2.1222 -0.193 0.318    0.011          0.011 -0.865 -0.394 -0.184  0.014
## theta2.1637 -0.515 0.328    0.011          0.012 -1.177 -0.734 -0.497 -0.289
## theta2.1995 -1.192 0.355    0.012          0.013 -1.891 -1.424 -1.189 -0.944
## theta2.3369  1.315 0.410    0.014          0.015  0.499  1.026  1.325  1.597
## theta2.4004  1.280 0.336    0.011          0.012  0.658  1.061  1.270  1.488
## theta2.4209  0.762 0.346    0.012          0.012  0.160  0.512  0.730  0.995
## theta2.4268  1.083 0.323    0.011          0.012  0.480  0.869  1.069  1.295
## theta2.4841 -0.339 0.338    0.011          0.013 -0.993 -0.576 -0.331 -0.105
## theta2.5074 -1.071 0.353    0.012          0.014 -1.732 -1.317 -1.062 -0.822
## theta2.5495  0.518 0.342    0.011          0.011  0.025  0.246  0.470  0.747
## theta2.5710 -0.929 0.355    0.012          0.013 -1.615 -1.167 -0.946 -0.679
## theta2.8294  0.737 0.327    0.011          0.010  0.110  0.512  0.730  0.959
## theta2.8980  0.904 0.340    0.011          0.012  0.254  0.658  0.895  1.136
## theta2.9899  0.614 0.332    0.011          0.012 -0.048  0.396  0.621  0.846
## gamma.1      0.854 0.133    0.004          0.006  0.592  0.762  0.857  0.949
## gamma.25     0.838 0.144    0.005          0.006  0.542  0.743  0.844  0.940
## gamma.26     0.854 0.133    0.004          0.006  0.600  0.761  0.855  0.949
## gamma.27     0.854 0.133    0.004          0.006  0.593  0.763  0.859  0.949
## gamma.28     0.378 0.213    0.007          0.008  0.023  0.196  0.382  0.540
## gamma.29     0.856 0.133    0.004          0.006  0.605  0.763  0.860  0.950
## gamma.30     0.290 0.204    0.007          0.008  0.014  0.108  0.249  0.462
## gamma.31     0.129 0.102    0.003          0.003  0.006  0.046  0.106  0.186
## gamma.34     0.857 0.135    0.004          0.006  0.593  0.765  0.861  0.951
## gamma.35     1.420 0.111    0.004          0.004  1.153  1.360  1.442  1.506
## gamma.36     1.418 0.111    0.004          0.004  1.150  1.354  1.439  1.501
## gamma.37     1.417 0.112    0.004          0.004  1.148  1.354  1.439  1.502
## gamma.38     1.424 0.107    0.004          0.004  1.159  1.363  1.445  1.507
## cluster.1    2.388 1.180    0.039          0.114  1.000  1.000  2.000  4.000
## cluster.25   2.401 1.171    0.039          0.101  1.000  1.000  2.000  4.000
## cluster.26   2.381 1.185    0.039          0.114  1.000  1.000  2.000  4.000
## cluster.27   2.396 1.180    0.039          0.111  1.000  1.000  2.000  4.000
## cluster.28   2.532 1.063    0.035          0.056  1.000  2.000  3.000  3.000
## cluster.29   2.376 1.171    0.039          0.108  1.000  1.000  2.000  3.250
## cluster.30   2.559 1.058    0.035          0.066  1.000  2.000  3.000  3.000
## cluster.31   2.634 1.071    0.036          0.071  1.000  2.000  3.000  4.000
## cluster.34   2.386 1.184    0.039          0.108  1.000  1.000  2.000  4.000
## cluster.35   2.371 1.136    0.038          0.106  1.000  1.000  2.000  3.000
## cluster.36   2.374 1.132    0.038          0.089  1.000  1.000  2.000  3.000
## cluster.37   2.373 1.135    0.038          0.109  1.000  1.000  2.000  3.000
## cluster.38   2.399 1.139    0.038          0.090  1.000  1.000  2.000  3.000
## n.clusters   3.790 0.408    0.014          0.015  3.000  4.000  4.000  4.000
## alpha        1.625 0.979    0.033          0.037  0.328  0.915  1.430  2.071
##              97.5% Model
## theta1.1121  1.970  2dDP
## theta1.1154  0.777  2dDP
## theta1.1222  1.283  2dDP
## theta1.1637  1.361  2dDP
## theta1.1995  1.669  2dDP
## theta1.3369  2.484  2dDP
## theta1.4004  0.012  2dDP
## theta1.4209 -0.495  2dDP
## theta1.4268  0.052  2dDP
## theta1.4841 -0.179  2dDP
## theta1.5074  1.956  2dDP
## theta1.5495 -1.467  2dDP
## theta1.5710 -0.086  2dDP
## theta1.8294  0.547  2dDP
## theta1.8980  1.345  2dDP
## theta1.9899  1.153  2dDP
## theta2.1121 -0.456  2dDP
## theta2.1154 -0.901  2dDP
## theta2.1222  0.417  2dDP
## theta2.1637  0.109  2dDP
## theta2.1995 -0.516  2dDP
## theta2.3369  2.128  2dDP
## theta2.4004  1.987  2dDP
## theta2.4209  1.489  2dDP
## theta2.4268  1.752  2dDP
## theta2.4841  0.342  2dDP
## theta2.5074 -0.420  2dDP
## theta2.5495  1.316  2dDP
## theta2.5710 -0.245  2dDP
## theta2.8294  1.410  2dDP
## theta2.8980  1.564  2dDP
## theta2.9899  1.222  2dDP
## gamma.1      1.099  2dDP
## gamma.25     1.099  2dDP
## gamma.26     1.099  2dDP
## gamma.27     1.099  2dDP
## gamma.28     0.775  2dDP
## gamma.29     1.100  2dDP
## gamma.30     0.686  2dDP
## gamma.31     0.368  2dDP
## gamma.34     1.100  2dDP
## gamma.35     1.564  2dDP
## gamma.36     1.563  2dDP
## gamma.37     1.563  2dDP
## gamma.38     1.564  2dDP
## cluster.1    4.000  2dDP
## cluster.25   4.000  2dDP
## cluster.26   4.000  2dDP
## cluster.27   4.000  2dDP
## cluster.28   4.000  2dDP
## cluster.29   4.000  2dDP
## cluster.30   4.000  2dDP
## cluster.31   4.000  2dDP
## cluster.34   4.000  2dDP
## cluster.35   4.000  2dDP
## cluster.36   4.000  2dDP
## cluster.37   4.000  2dDP
## cluster.38   4.000  2dDP
## n.clusters   4.000  2dDP
## alpha        4.047  2dDP

Plot 2dDP Model posterior parameters from samples

Trace Plots

First traceplots of theta parameter iterations. This bayesplot function has a default fixed scale for x axis, which enables better comparison between parameters than the traceplot function in coda. So fixed scales work best with same parameters.

Start with thetas

# Trace plot only for theta parameters
mcmc_trace(mcmc_list2dDP, pars = thetas, facet_args = list(ncol = 2, scales = "fixed"))

Next the gammas with their own fixed y axis.

# Trace plot only for gamma
mcmc_trace(mcmc_list2dDP, pars = gammas, facet_args = list(ncol = 2, scales = "fixed"))

Start with finally the clusters with their own fixed y axis.

# Trace plot only for cluster parameters
mcmc_trace(mcmc_list2dDP, pars = clusters, facet_args = list(ncol = 2, scales = "fixed"))

Next density plots of all posterior parameter distributions

Density Plots

densplot(mcmc_list2dDP, right = FALSE) 

# Right = FALSE helps makes each cluster integer fall neatly to the right of its own bin in histogram

Autocorrelation Plots

Next the autocorrelation plots

autocorr.plot(mcmc_list2dDP, auto.layout = FALSE, lag.max=20)

Create a tidy data table

summary2dDP <- summary2dDP %>%
  mutate(
    Variable = str_extract(rownames(.), "^[^\\.]+"),
    ID = str_extract(rownames(.), "(?<=\\.)\\d+")
  ) %>%
  `rownames<-`(NULL) %>%
  select(Model, Variable, ID, everything() )

summary2dDPWideThetas <- summary2dDP %>%
  filter(Variable == "theta1" | Variable == "theta2") %>%
  select(ID, Variable, '50%', '2.5%', '97.5%') %>%
  pivot_wider(names_from = Variable, values_from = c('50%', '2.5%', '97.5%'))

Plot the posterior means of items

ggplot(summary2dDPWideThetas, aes(x = `50%_theta1`, y = `50%_theta2`, label = ID)) +
  # Add horizontal error bars (theta1 CI)
  geom_errorbarh(aes(xmin = `2.5%_theta1`, xmax = `97.5%_theta1`), height = 0.05, color = "gray60", linewidth = 0.3) +
  # Add vertical error bars (theta2 CI)
  geom_errorbar(aes(ymin = `2.5%_theta2`, ymax = `97.5%_theta2`), width = 0.05, color = "gray60", linewidth = 0.3) +
  # Add point and label
  geom_point() +
  geom_text(vjust = -0.5, size = 3) +
  labs(
    title = "Posterior Medians of and 95% Credible Intervals",
    x = "Theta 1",
    y = "Theta 2"
  ) +
  theme_minimal()

Plot the gamma median for each rater.

gammaAngles <- summary2dDP %>%
  filter(Variable == "gamma") %>%
  select(ID, gamma = `50%`) %>%
  mutate(
    x = cos(gamma) * 2,  # x component
    y = sin(gamma)* 2   # y component
  )

ggplot(gammaAngles, aes(x = 0, y = 0, xend = x, yend = y)) +
  geom_segment(arrow = arrow(length = unit(0.15, "cm")),
               color = "black", linewidth = 0.7) +
  coord_fixed() +  # 1:1 aspect ratio
  labs(
    title = "Rater Preference Vectors (Gamma Angles)",
    x = "Theta 1",
    y = "Theta 2"
  ) +
  # Smart label positioning
  geom_text_repel(
    data = gammaAngles,
    aes(x = x, y = y, label = ID),
    size = 3,
    segment.color = NA,
    nudge_x = 0.1, nudge_y = 0.1,
    max.overlaps = Inf) +
  
  theme_minimal()