How to sort the two change points in a longitudinal model?

Dear Stan community,

I have a two change points binary longitudinal model and I would like to make sure the change point 2 (cp[2]) is later than change point 1 (cp[1]). Does anyone one know how to do this?

The model is:

p[i] <-  c0 + c1 * (X[i]-cp[1]) + c2 * (X[i]-cp[1]) * step(X[i]-cp[1]) + c3 * (X[i]-cp[2]) * step(X[i]-cp[2]) ;

I have set the prior of ‘cptemp’ follows the uniform distribution for some period of time and the ‘cp’ variable might link to the ‘cptemp’ such as using cp[1] = min(cptemp) and cp[2]=max(cptemp) or cp[1:2] ← sort_asc(cptemp), but I do not know where to put them and I got some errors for them.

    for (k in 1:2){
	 cptemp[k] ~ uniform(0.07,21.45);		

Or should I set some rule of the ‘cp’ variable ?

The full code is below:

// Binary long model for PA+
  data {
    // Define variables in data
    // Number of observations (an integer)
    int N;
    // Number of individual
    int M;

    // Covariates
    int X[N];
    int num[N];
    // Outcome
    int Y[N];
  parameters {
    // Define parameters to estimate
    real c0;
    real c1;
    real c2;
    real c3;
    real cp[2];
    real cptemp[2];
    real s[N];
    real u[N];
    real ga;
    real sig[2];
  transformed parameters  {
    real p[N];
    for(i in 1:N){ 
       p[i] <-  c0 + c1 * (X[i]-cp[1]) + c2 * (X[i]-cp[1]) * step(X[i]-cp[1]) + c3 * (X[i]-cp[2]) * step(X[i]-cp[2]) + ga*s[num[i]] + u[num[i]];
   model {
    // Prior part of Bayesian inference
    c0 ~ normal(0,10000);
    c1 ~ normal(0,10000);
    c2 ~ normal(0,10000);	
    c3 ~ normal(0,10000);	
    sig[1] ~ uniform(0,100);
    sig[2] ~ uniform(0,100);
    for(i in 1:M){ 
         u[i] ~ normal(0,sig[1]);
         s[i] ~ normal(0,sig[2]);
    ga ~ normal(0,10000);

    for (k in 1:2){
	 cptemp[k] ~ uniform(0.07,21.45);		

   Y ~ bernoulli_logit(p);

Could you please help? Thanks a lot!

You could consider defining cp as an ordered vector?

What about defining change point 2 as change point 1 plus some positively-constrained value, i.e.

cp[2] <- cp[1] + cp_diff

And depending on your knowledge of the change points, you could use a log-transformation to force cp_diff to be positive, or just a prior with a zero lower bound.