Rebase Inducing Points Using QR

[akleeman]

Prereq: #475

The existing rebase_inducing_points method relied heavily on the fit_from_prediction method. It was computing a prediction of the new inducing points, then calling fit_from_prediction on that result,

fit_from_prediction(new_inducing_points, fit_model.predict(new_inducing_points).joint());

This involved a few questionable steps. For example, if the prediction is not full rank (which could happen when you have numerically identical inducing points) then this step:

    const Eigen::SerializableLDLT C_ldlt(prediction.covariance);
    const Eigen::MatrixXd sigma_inv_sqrt = C_ldlt.sqrt_solve(K_zz);

could become unstable since C_ldlt may have near zero diagonals and we then invert them in C_ldlt.sqrt_solve.

This PR switches to a new approach which is more directly based on the QR decomposition. While there is still a step involving the cholesky of a potentially singular matrix, the subsequent decomposition is not used in any inverses, so the consequences should be much smaller.

A description of the algorithm was added to the documentation (screenshot below)

[akleeman]

I tried pulling this into our models and it actually significantly changed performance 🤔. Going to need to root cause why, in the meantime I won't merge this (but don't let that stop you from taking a look with an eye out for any mistakes I may have made).

[peddie]

I don't see any obvious performance missteps here. But Householder QR is always twice as expensive as an equivalently-sized Cholesky decomposition.

I am glad to see the permutation matrix cleanups going in -- might it make sense to pull those out as a separate first PR? I would happily merge that stuff, and it would also leave fewer changes to examine for performance.

[akleeman]

@peddie #475

[akleeman]

I don't see any obvious performance missteps here. But Householder QR is always twice as expensive as an equivalently-sized Cholesky decomposition.

Both this and the existing require LDLT and QR decompositions. But it is true that this version will likely be slower since it requires two LDLTs and a larger QR.

That said, this approach opens the opportunity to do both a rebase and update in one single step (not implemented yet) at which point the cost might even out?