Random sampling of sorted numbers within range

[junpenglao]

You might find some previous discussion helpful here:

Sampling uniformly in a triangular support

There are a couple of ways to do it: You can use the ordered transformation here. You can add jacobian correction for the volume change in your case above as well to make it proper. You can jointly sample from a bivariate standard Gaussian then rotating them by 45 degrees. The rotated sample is transformed into rates that lie in the unit square. See here. You can add a potential to restrict the order, see here cell [4] for an example

TL;dr:

tr = pm.transforms
Order = tr.Ordered()
Logodd = tr.LogOdds()
chain_tran = tr.Chain([Logodd, Order])

n = 5
with pm.Model() as m0:
    x = pm.Uniform("x", 0.0, 1.0, shape=n, transform=chain_tran,
                   testval=np.linspace(.1, .9, n))  # test value need to be ordered
    trace = pm.sample()