That implementation looks very useful, thanks, and will certainly be a lot easier to tweak than the C code referenced above.
The problem I face is that I am trying to find a universal method of initialising the MCMC runs for around 1,200 different demand patterns (time series) as I do not have the resource to tweak each one. For instance, I have one simple case where there appear to be local minima at 3 and 5 knots. Sometimes the MCMC will have a modal number of knots of 3, sometimes 5, and sometimes it is bimodal (3 and 5). I can see from the BIC values that usually 3 knots gives a lower minimum BIC - and hence a ‘better model’ by some measure - than 5 knots. But with a run in which the modal number of knots is 5, the lowest BIC is given by a model with 8 knots.
Any help with this would be very welcome.