actually, here’s an attempt…
def piecewise_constant(x, beta_region, knots):
# Start with the left-most region
y = beta_region[0] * T(x, -np.inf, knots[0])
# Add constants for each middle region
for i in range(len(knots) - 1):
y += beta_region[i + 1] * T(x, knots[i], knots[i + 1])
# Right-most region
y += beta_region[-1] * T(x, knots[-1], None)
No attempt has been made to check the priors or make them robust to other datasets. But this should help hopefully?
Sorry it’s been a long day. You’ll want to put priors over the break points. Probably need to enforce ordering. Can take a look in the morning if it’s not clear how to do that