Dealing with missing data and custom distribution

This is the minimal working example I can get to work.

    ...: with pm.Model() as model: 
    ...:  
    ...:     ## Independent priors 
    ...:     alpha = pm.Normal('User', mu = 0, sigma = 1) 
    ...:     gamma = pm.Normal('Task', mu = 0, sigma = 1) 
    ...:  
    ...:     ## Log-Likelihood 
    ...:     def logp(d): 
    ...:         return pm.Normal.dist(alpha, pm.math.exp(gamma)).logp(d) 
    ...:  
    ...:     ll = pm.DensityDist('ll', logp, observed = observations)          
    ...:     trace = pm.sample() 

Having the dictionary in observed seems to fail.
But removing it is not enough for your example. I still think it has something to do with your logp returning an invalid value or gradient at some point. The sampler always seems to fail after 38%… or something. We can also look at where it fails to have a better idea.