# Ordered vector uniform in angle

I need to sample angle uniformly which is simple enough:

``````parameters {

vector<lower=0, upper=1>[N] angle_raw;

}

transformed parameters {

vector[N] angle;

angle = acos(1-angle_raw); // will be uniformly distributed between 0 and pi/2

}

``````

However, if I need the angle to be ordered, I lose the ability to place constraints. So, following the manual:

``````parameters {

ordered angle_rawer[N];

}

transformed parameters {

vector<lower=0, upper=1>[N] angle_raw;
vector[N] angle;

angle_raw = inv_logit(angle_rawer);

angle = acos(1-angle_raw); // will be uniformly distributed between 0 and pi/2

}

model {

angle_raw ~ uniform(0,1);

target += log(angle_raw) + log(1 - angle_raw); // jacobian adjustment.

}

``````

Is this the correct way to handle this?

Sorry, no idea what the correct way is, just noticed the following line which might be a problem:

2 Likes

Oops. That was just a typo in the example! Thanks!

Ok, after hunting a bit (and understanding), I think I might have a solution. While I’m not sure if it works in the bigger model I’m building, I have tested a bit and I think it is the right path.

I tested this with a little python code to make sure the general idea worked:

``````def simplex_sample(n):
"""
sample from an n-dim unit simplex
"""
sample = np.random.uniform(0,1,size=n-1)
sample = np.append(sample,1.)
sample.sort()
sample = np.append(0,sample)

diffs = np.diff(sample)

return diffs

def ordered_interval(n, a = 0, b=1.):
"""
sample an ordered (distinguishable) interval
"""

simplex = simplex_sample(n+1)

return a+ np.cumsum(simplex)[:n]*(b-a)
``````

Which yields:

``````
sample = ordered_interval(n=10000, a = 0, b= 1)

len(np.unique(sample)) # it equals 10000

``````

And then I transform this to variable uniform in angle:

``````

``````

Thus, the method should produce an ordered angle uniform in angle.

The stan code should be:

``````parameters {

simplex angle_rawer[N+1];

}

transformed parameters {

postive_ordered[N] angle_raw;
vector[N] angle;

angle_raw = head( cumulative_sum(angle_rawer), N );

angle = acos(1-angle_raw); // will be uniformly distributed between 0 and pi/2 in angle

}

model {

angle_raw ~ uniform(0,1);

}

``````

And no jacobian adjustments needed (I believe! )

Update: This works great in the bigger model.