Order statistics in PyMC3

No worries. I got this in this way.

import numpy as np
import pymc3 as pm
import theano.tensor as tt

# data
K = 5 # number of items being ranked
J = 20 # number of raters
yreal = np.argsort(np.random.randn(K, 1), axis=0)
print(yreal)
y = np.argsort(yreal + np.random.randn(K, J), axis=0)
print(y)
# transformed data{
y_argsort = np.argsort(y, axis=0)

oorder=pm.distributions.transforms.Ordered()
with pm.Model():
    mu_hat = pm.Normal('mu_hat', 0, 1, shape=K-1)
    # set first value to 0 to avoid unidentified model
    mu = tt.concatenate([[0.], mu_hat])
    #sd = pm.HalfCauchy('sigma', 1.)
#     latent = pm.Normal('latent',
#                        mu=mu[y_argsort],
#                        sd=1., # using sd does not work yet
#                        transform=Ordered2D(), 
#                        shape=(K,J),
#                        testval=np.repeat(np.arange(K)[:,None], J, axis=1))
#                         # There are some problems with Ordered 
#                         # right now, you need to give a testval
    latent_probs=[]
    for j in np.arange(J):
        _latent=pm.Normal('latent%i' %j, mu[y_argsort[:,j]],1,shape=K,transform=oorder,testval=np.arange(K))
        latent_probs.append(_latent)
    #_latent=pm.Normal('latent%i' %j, mu[y_argsort[:,j]],1,shape=,transform=oorder,testval=np.arange(K))
    latent=pm.Deterministic("latent",tt.stack(latent_probs,axis=1))
    trace = pm.sample(init="adapt_diag")