vectorized PL

anno-difficulty/paper/pvalues.py

@@ -19,16 +19,18 @@
 data = { 'I': I, 'J': J, 'N': N,
             'item': item, 'rater': rater, 'rating': rating }
 init = { 'prev': 0.2,
+         'sigma': 1,             
          'sens': np.full(J, 2),
          'spec': np.full(J, 2),
          'diff': np.full(I, 0.0) }
 model = csp.CmdStanModel(stan_file = '../stan/revised-1pl.stan')
 sample = model.sample(data = data, show_console = True, refresh = 10,
-                          iter_warmup=200, iter_sampling=200,
-                          parallel_chains = 4,
-                          chains = 2, inits = init,
-                          seed = 92584)
+                          iter_warmup=1000, iter_sampling=1000,
+                          parallel_chains = 4, adapt_delta=0.9,
+                          chains = 4, inits = init,
+                          seed = 837689)
 
+# ACTUAL W DENTISTRY: [1880,1065,404,247,173,100]
 
 # ORIGINAL
 # y = df.pivot(index='item', columns='coder', values='response').to_numpy()

anno-difficulty/stan/revised-1pl.stan

@@ -30,13 +30,16 @@ parameters {
   real<lower=0, upper=1> prev;
   vector[J] spec;
   vector<lower=-spec>[J] sens;
-  vector[I] diff;
+  real<lower=0> sigma;
+  vector<multiplier=sigma>[I] diff;
+  // vector[I] diff;
 }
 model {
   prev ~ uniform(0, 1);
   sens ~ normal(0, 3);
   spec ~ normal(0, 3);
-  diff ~ normal(0, 1);
+  sigma ~ lognormal(0, 1);
+  diff ~ normal(0, 1.5);
   vector[I] lp_pos = rep_vector(log(prev), I);
   vector[I] lp_neg = rep_vector(log1m(prev), I);
   for (n in 1:N) {

sushi-rating/plackett-luce-vectorized.stan

@@ -0,0 +1,22 @@
+functions {
+  real plackett_luce_lpmf(array[] int y, vector beta) {
+    vector[size(y)] beta_y = beta[y];
+    // return sum(log(beta_y ./ cumulative_sum(beta_y)));
+    return sum(log(beta_y) - log(cumulative_sum(beta_y)));
+  }
+}
+data {
+  int<lower=1> I;                       // # items
+  int<lower=1> K;                       // # items ranked per rater
+  int<lower=1> R;                       // # raters
+  array[R, K] int<lower=1, upper=I> y;  // rankings (y[r, 1] < y[r, 2] < ...)
+}
+parameters {
+  simplex[I] beta;                     // item quality
+}
+model {
+  vector[I] log_beta = log(beta);
+  for (r in 1:R) {
+    y[r] ~ plackett_luce(beta);
+  }
+}