Strange behaviour of QR re-parameterized splines

Modeling
1

As I understand, the QR re-parameterization is useful in cases with correlated predictors, of which I believe splines would be a good example, so I coded up a spline fitting model with QR re-parameterization of the splines. It seems to get higher ESS/time as compared to the non-qr version, but when I plot the posterior samples of the function being fit, there’s some strange jagged portions at the left-hand side that don’t manifest when I do the non-qr version:

Rplot02

Below is the Stan code; any ideas what’s going awry?

data {
    // n_y: number of observations in y
    int n_y ;
    // y: vector of observations for y
    vector[n_y] y ;
    // x: unique values of x
    vector[n_y] x ;
    //num_splines: number of splines
    int num_splines ;
    //splines: matrix of splines
    matrix[n_y,num_splines] splines ;
}
transformed data{
    matrix[num_splines,n_y] Q_ast;
    Q_ast = transpose( qr_Q(splines) [, 1:num_splines] * sqrt(n_y - 1.0) ) ;
}
parameters {
    // noise: measurement noise
    real<lower=0> noise ;
    // volatility: wiggliness of the functions
    real<lower=0> volatility ;
    //intercept: intercept of the function
    real intercept ;
    //slope: slope of the function
    real slope ;
    // splines_helper: helper variable for splines
    row_vector[num_splines] splines_helper ;
}
transformed parameters{
    // f: latent functions
    vector[n_y] f ;
    f = intercept + slope*x + //linear part
        to_vector((volatility*cumulative_sum(splines_helper))*Q_ast) ; //wiggly part
}
model {
    intercept ~ normal(0, 1) ;
    slope ~ normal(0, 1) ;
    volatility ~ normal(0, 1) ;
    noise ~ normal(0, 1) ;
    splines_helper ~ normal(0,1);
    y ~ normal(f,noise);
}
2

My guess is you are giving the splines too much flexibility, which when combined with the QR reparameterization means that it can actually get into and out of the region of the parameter space where it overfits.

1 Like