Recovering noise-free variables in binary trials I have a collection of $N$ different objects. An expert "looks" at the object $O_i$ and produces a prediction: binary label $y_i\in\{0,1\}$. The predicted binary label $y_i$ is observed (data) and is known to be noisy. Also, objects are of different complexity: an expert with poor "skill" will likely mislabel a work item.
Now assume that there are $K$ experts (each with a different "skill") with a consensus on what it means to label "1" or "0". The $k$th expert produces an output label $y^k_i$ (data) for object $O_i$. Therefore, each object has $K$ observations.
Additionally, each object has observed real-valued features vector $x_i$ (data). Note that the feature vectors are not used by experts; features are only noisy observations that are believed to be useful for prediction given no experts!
The main goal is to estimate noise-free/hidden $\hat y_i$ label for each object. It is expected that feature vectors and multiple experts are complimentary and helpful to discover the hidden latent labels.
For example, there are 5 objects with hidden true labels: [1 0 0 1 1]. An expert's intrinsic (hidden) flipping probability is $r=0.1$. This means that the observed outputs for this experts is a vector of length 5 where each label is flipped independently with a probability of 10%.
It is unclear to me how to build a model that attempts to:


*

*uncover hidden labels $\hat y_i$

*infer the $k$th worker "skill" (flipping probability?)


I started by fitting a logistic regression model but realized that I am fitting different quantity.
My package of choice is PyMC3 but I am open to other probabilistic programming packages such as Stan, JAGS.

Toy (made up) data:
Suppose there are 3 objects and 4 experts.
The feature vectors are given (4 dimensions):
$x_1 = [0.3, 0.4, 0.01, 0.9]$
$x_2 = [0.1, 0.9, 0.05, 0.6]$
$x_3 = [0.9, 0.1, 0.15, 0.3]$
One can train a binary classifier that takes object vectors and predicts binary label. The problem is that the feature vectors are known to contain fairly little information about the prediction task, thus the hope on crowd of experts.
The experts (independently) make their label predictions on the three objects (observed):
Expert 1: $y^1 = [1,1,1]$
Expert 2: $y^2 = [0,1,0]$
Expert 3: $y^3 = [1,1,0]$
Expert 4: $y^4 = [1,0,0]$
One can simply use majority vote to approximate the true label. The problem is that the expert predictions are noisy and experts have different expertise.
The data are object feature vectors $x_i$ and expert label vectors $y^k$.
The hidden/latent factors are: noise-free object labels $\hat y$ and expert flipping probabilities $r^k$.
 A: After few weeks, I believe that I found the model that addresses my problem pretty accurately. The Dawid-Skene model, originally presented in EM version in 1979, models worker abilities as confusion matrices. The probabilistic programming version of the model, presented in the Stan manual, mixes and converges quickly as the discrete random variables for categories are summed out.
There are other models that take into account work item difficulties that are matched against worker skills. Some also attempt to model communities in the population of annotators for improved work performance estimates. I decided to start with a simpler model.
/*
 * Dawid and Skene. 1979.
 *
 * Extended with fixed priors.
 *
 * Discrete parameters for category of items marginalized out.
 *
 *  p(y | pi, theta)

 *  = PROD_i SUM_k p(k | pi) 
 *                 * PROD_{n: ii[n]=i} p(y[n] | theta[jj[n], k])
 * 
 *  log p(y | pi, theta)
 *
 *  = SUM_i log SUM_k p(k | pi) 
 *                    * PROD_{n: ii[n]=i} p(y[n] | theta[jj[n], k])
 *
 *  = SUM_i log SUM_k exp(log ( p(k | pi) 
 *                              * PROD_{n: ii[n]=i} p(y[n] | theta[jj[n], k])))
 *
 *  = SUM_i log SUM_k exp(log p(k | pi)
 *                        + SUM_{n: ii[n]=i} log p(y[n] | theta[jj[n], k]))
 *
 * = SUM_i LOG_SUM_EXP_k log p(k | pi)
 *                        + SUM_{n: ii[n]=i} log p(y[n] | theta[jj[n], k])
 */

data {
  int<lower=2> K;              // number of categories
  int<lower=1> I;              // number of items
  int<lower=1> J;              // number of coders
  int<lower=1> N;              // number of observations
  int<lower=1,upper=J> jj[N];  // coder for observation n
  int<lower=1,upper=I> ii[N];  // item for observation n
  int<lower=1,upper=K> y[N];   // label for observation n
  vector<lower=0>[K] alpha;    // prior for prevalence (positive)
  vector<lower=0>[K] beta[K];  // prior for coder responses (positive)
}
parameters {
  simplex[K] pi;               // prevalence of categories
  simplex[K] theta[J,K];       // response of anotator j to category k
}
transformed parameters {
  real cat_log[I,K];           // vector of log probs (up to const) for item i
  for (i in 1:I)
    for (k in 1:K)
      cat_log[i,k] <- log(pi[k]);  // efficiency: compute log(pi[k]) once
  for (n in 1:N)
    for (k in 1:K)
      cat_log[ii[n],k] <- cat_log[ii[n],k] + log(theta[jj[n],k,y[n]]);
}
model {
  pi ~ dirichlet(alpha);
  for (j in 1:J)
    for (k in 1:K)
      theta[j,k] ~ dirichlet(beta[k]);

  for (i in 1:I)
    increment_log_prob(log_sum_exp(cat_log[i]));
}

The model naturally handles missing data (not all workers annotate every work item) and also provides probabilistic estimates for parameters of interest.
There are still a lot of experiments to do with this model but it appears that the results are pretty meaningful. At least on my crowdsourced data.
