Robust and Conjugate Spatio-Temporal Gaussian Processes

We thank the reviewer for their careful consideration of our paper and helpful feedback. We address below the valuable suggestions made to further improve our work:

Scaling is shown to be linear in the number of time steps (but not in the number of spatial locations).

We agree with this concern. In the camera-ready version, we will replace the claim with “An alternative approach which has linear-in-time cost is to use a state-space representation.”

It would be useful (...) if the authors ran an experiment in which the centering function was adaptive, but all other parameters were kept as in RCGP to show that this does in fact reduce sensitivity.

We agree that isolating the effect of the adaptive centering function would strengthen our claim regarding reduced sensitivity. We will add an experiment comparing ST-RCGP and RCGP where all parameters are kept as in RCGP apart from the centering function. We further develop it by also altering the shrinking function. This analysis can be found here: https://anonymous.4open.science/r/ST-RCGP-21DD/tests/RCGP-vs-ST-RCGP-centering-function-sensitivity-analysis.ipynb. We welcome any further feedback.

Consider including coverage plots for at least one other dataset, for example the temperature dataset considered.

We agree that having more datasets for which we demonstrate the (frequentist) coverage of ST-RCGP would help further strengthen our claim. We will include in the paper a coverage analysis of the ST-RCGP for the temperature dataset. This analysis can be found at https://anonymous.4open.science/r/ST-RCGP-21DD/experiments/weather-forecasting/model-fitting.ipynb (at the end of the notebook).

Why is decontaminated data used for selection of hyperparameters? Wouldn't it be more realistic to use the contaminated data?

Fitting the RCGP on contaminated data leads to poor estimates due to overfitting outliers during hyperparameter optimisation. Using decontaminated data improves RCGP results. This choice thus strengthens our claim: even with this advantage, RCGP underperforms compared to ST-RCGP.

I'm not convinced this is a meaningful timing comparison.

We don’t include optimisation time because it can be easily skewed by tempering with the optimisation procedure, especially for STGPs, which rely on gradient-based optimisation without clear stopping criteria, number of optimisation steps, learning rate, and even optimisers. Either way, the two optimisation objectives (for STGP and ST-RCGP) are extremely similar, and their computational cost are nearly identical.

Focusing on inference time is meaningful in online settings where one-step inference matters. For example, real-time applications like stock price estimation may preclude methods that iterate at each time step (e.g. variational STGP) due to time constraints, favoring closed-form solutions such as the Kalman filter, the robust filter from Durán-Martín et al. (2024), or the ST-RCGP.

Appendix C.10 are the parameter presented the parameter selected after optimization, or the initialized values of parameters?

The parameters presented are the ones post-optimization.

There should be a more clear description of the contribution relative to the prior work Duran-Martin et al, 2024. In particular, what is new?

The key differences are:

1. We deal with STGP problems such as hyperparameter optimisation and smoothing.

2. We carefully specify the weight function (see answer to vhF8 on IMQ) and its parameters to improve robustness (see Weight Function paragraph on line 237), whereas Duran-Martin et al. (2024) provide few justifications for the centering and shrinking functions and learning rate, which are parameters intrinsically connected to robustness and crucial to improving performance.

3. Duran-Martin et al (2024) use a weighted log-likelihood, while we use the weighted score-matching loss from Altamirano et al. (2024), leading to a different update step.

We will include in the paper, at line 226, “However, these methods do not deal with STGP problems such as hyperparameter optimisation and smoothing and do not investigate downweighting optimally.”

I don't understand the caption of figure 2.” “Comment on appendix A about abuse of notation.

We thank the reviewer for bringing this to our attention. We will remove this sentence in the caption in the camera-ready version, and make clear that we are abusing notation.

How does the point in time at which outliers occur affect the sensitivity of your method?.

An empirical analysis we conducted, available at https://anonymous.4open.science/r/ST-RCGP-21DD/tests/outliers-early-vs-late.ipynb, demonstrates that the method is not affected by the impact of outlier position, which we expect to be because of smoothing.

Finally, we want to thank the reviewer again for the careful review and consideration given to our paper. We hope this rebuttal addresses any remaining questions and concerns.

We thank the reviewer for their thoughtful review and positive remarks on our empirical evaluation and the ST-RCGP’s novel combination of robustness and efficiency. We comment on the feedback below:

Since the data is synthetic, I have some concerns about how meaningful the results really are.

Robustness studies often inject outliers due to the lack of datasets with labelled outliers and the ability to control testing across varied outlier types (as was done in our work). However, we will add an experiment in the paper with the well-log dataset from [8, 11], which consists of 4,050 nuclear magnetic resonance measurements recorded while drilling a well. Outliers were identified scientifically and correspond to short-term events in geological history. The results can be found here: https://anonymous.4open.science/r/ST-RCGP-21DD/experiments/well-log-data/NMR-data-fit.pdf with the code in the same folder.

The statement of Proposition 3.2 does not seem to be very strong.

Whilst we empathise with this comment, the study of the influence function (IF) dates from the seminal work of Hampel in 1968 [1] and has been the gold standard in robust statistics ever since [2,3,4]. The typical approach to show robustness is to bound the influence function. While one could work out explicitly the bounds, there is no clear gain in doing so as the values of the KL divergence are not very interpretable. ​​To prove the robustness of an algorithm, bounding the IF is then sufficient and remains one of the strongest (and in many cases the only) available criteria within the robust Bayesian literature [see e.g. 5,6,7].

Can you comment on any alternatives to the IMQ kernel for assigning the weights?

We choose the IMQ for several reasons. To elaborate on those, in the final version, we change lines 237-243 with: “The statements of Propositions 3.1 and 3.2 impose some constraints on the choice of weight function. In particular, it should be strictly positive and differentiable over its domain to ensure quantities in Proposition 3.1 are well-defined, and it must have a bounded supremum to satisfy Proposition 3.2. Moreover, to avoid attributing weight to arbitrarily large $y\in \mathbb{R}$, we require $\lim_{|y| \rightarrow \infty}w(\mathbf{x}, y) = 0$ (decaying property). These requirements make the IMQ an appropriate choice of weight function, in addition to the fact that it has been well-studied and recommended in prior literature [6, 9, 10,11]. We further justify our choice on the basis that the IMQ hyper-parameters γ,$\beta$ and c are related to concepts of robustness—relations we exploit to specify the IMQ. In particular, to select γ, β and c, we follow four guiding principles; we want to: …”

We acknowledge the need for a sensitivity analysis to understand how the shape of the weight function affects ST-RCGP's posterior estimates. In the camera-ready version, we will examine how varying the IMQ exponent—currently $\alpha:=-1/2$—impacts results. We conduct this analysis because Proposition 3.2 shows robustness requires weights to decay faster than $1/\sqrt{|y|},$, i.e., $\alpha<-1/4$; but, overly fast decay can reduce statistical efficiency (overly robust)—highlighting a tradeoff worth exploring. The result of this analysis can be found at https://anonymous.4open.science/r/ST-RCGP-21DD/tests/IMQ-sensitivity/IMQ-exponent-testing.pdf, with the code in the same folder.

Thank you again for the thorough review. We hope this rebuttal addresses any remaining concerns and contributes to a stronger evaluation of our paper.

[1] Hampel, Frank Rudolf. Contributions to the theory of robust estimation. 1968.

[2] Huber, Peter J., and Elvezio M. Ronchetti. Robust statistics. 2011.

[3] Maronna, Ricardo A., et al. Robust statistics: theory and methods (with R). 2019.

[4] Hampel, Frank et al. Robust statistics: The approach based on influence functions. 1987.

[5] Ghosh, A. and Basu, A. Robust Bayes estimation using the density power divergence. Annals of the Institute of Statistical Mathematics, 68(2), 2016.

[6] Matsubara, T., Knoblauch, J., Briol, F.-X., and Oates, C. J. Robust generalised Bayesian inference for intractable likelihoods. Journal of the Royal Statistical Society: Series B (Statistical Methodology). 2022.

[7] Duran-Martin, Gerardo, et al. "Outlier-robust Kalman filtering through generalised Bayes." 2024.

[8] Altamirano, Matias, François-Xavier Briol, and Jeremias Knoblauch. "Robust and scalable Bayesian online changepoint detection." International Conference on Machine Learning. PMLR, 2023.

[9] Chen, Wilson Ye, et al. "Stein point markov chain monte carlo." International Conference on Machine Learning. PMLR, 2019.

[10] Riabiz, Marina, et al. "Optimal thinning of MCMC output." Journal of the Royal Statistical Society Series B: Statistical Methodology 84.4 (2022).

[11] Ruanaidh, J. J. O. and Fitzgerald, W. J. Numerical Bayesian methods applied to signal processing. Springer Science & Business Media, 1996.

We thank the reviewer’s feedback and are glad our theoretical results and experimentation on our proposed method were appreciated. It has been mentioned that:

“The baseline models used for comparison vary across experiments due to data characteristics. It would be beneficial for the authors to include an additional state-of-the-art model beyond RCGP to better assess the overall effectiveness of the proposed method.”

We agree that it would be beneficial to illustrate the fit of additional methods apart from RCGP. For this reason, we will include in the paper the following plot: https://anonymous.4open.science/r/ST-RCGP-21DD/experiments/financial-applications/high-frequeny-data/HFT-data-induced-crash-fit.pdf, which was produced by code available here: https://anonymous.4open.science/r/ST-RCGP-21DD/experiments/financial-applications/high-frequeny-data/HFT-index-futures-speed-comp.ipynb. This plot compares ST-RCGP against some of the state-of-the-art methods offered by the “BayesNewton” package and is an extension of Table 2.

Again, we are grateful for the consideration given to our paper and hope our answer satisfies the concern raised.

We appreciate the thoughtful feedback and the expertise on the matter. We are glad you found the paper clearly written and easy to follow. We comment on the feedback below:

1) “provides inferences that are comparable to state-of-the-art non-Gaussian STGPs in the presence of outliers, but at a fraction of the cost" (...) I find this claim a bit too general.

We acknowledge that our claim is too broad and will adjust it as follows: “The ST-RCGP provides inferences comparable to state-of-the-art non-Gaussian STGPs with a computational cost similar to classical STGPs.” Additionally, we believe further accelerating the ST-RCGP is an exciting direction for future work, so we will replace lines 422-426 (right side) with: “Our method builds on classical STGPs, which remain vulnerable to scaling in the spatial dimension and are not optimal if parallel computing is available. These issues were addressed by Hamelijnck et al. (2021) through variational approximations and parallel-scan algorithms (see also [1]), and could be adapted to ST-RCGPs, providing the same benefit.”

2) Regarding the EWR: it is stated that "since the STGP is not robust to outliers, the closer EWR is to one, the less robust our posterior is to outliers, and vice-versa." I don't think this statement is strictly accurate

We agree and will replace the highlighted sentence with: “Therefore, since the STGP is not robust to outliers, EWRs that are near one are not necessarily optimal, since an EWR of 1 implies a solution—the vanilla STGP—which is not robust."

3) I think it could avoid possible confusion to more clearly state that Proposition 3.1 assumes the modified Fisher divergence.

Thank you for raising this point. We will ensure to make this clearer in the camera-ready version by replacing the current text in Proposition 3.1 with: “Then, we obtain a (generalised) posterior (…) with the score-matching loss $\mathcal{L}$ given by X for Y” where X will be replaced by the loss as provided in lines 674-675, and Y the mathematical objects used in the loss, which are defined on line 677.

4) There is a small body of work using 3σ-rejection-esque filters for Markovian GPs that are not cited.

​​Thank you for bringing this to our attention. We will include the following statement after the last sentence of the paragraph on line 066: “We also highlight a small body of work that uses outlier-rejection Kalman filters with STGPs to improve robustness in inference tasks such as prediction and anomaly detection. While these methods offer robustness at lower computational cost compared to non-conjugate STGPs, they are generally considered less expressive and not as strongly supported by theoretical foundations.”

5) minor editorial comments

We agree with these editorial comments and will make the necessary changes in the final version.

6) It is written "assuming a stationary and separable kernel (...)" but stationarity is not sufficient for representation as an LTI SDE (e.g., the squared exponential kernel).

Thank you for highlighting this. We will update the final version accordingly.

7) Are there obvious technical barriers in also combining these works, i.e., having computation-aware ST-RCGPs?

We believe there should be no fundamental technical barriers to combining these approaches because the filtering update in ST-RCGP closely mirrors that of the standard STGP. Given this structural similarity, the main arguments presented by Pförtner, Marvin, et al (2024) [2] should remain applicable.

8) How does the smoothing work?

The generalised Bayes approach used on the filtering distribution produces a Gaussian posterior, and thus, the resulting smoothing distribution remains Gaussian. Therefore, as for classical STGPs, we can apply RTS smoothing using the GB filtering posterior. To clarify ambiguities about smoothing, the above explanation will be added in the paragraph on lines 224-235 (Left side, Methodology section).

Finally, we want to thank the reviewer again for the careful consideration given to our paper. We hope we have responded to any remaining questions or concerns in this rebuttal.

References

[1] Särkkä, Simo, and Ángel F. García-Fernández. "Temporal parallelisation of Bayesian smoothers." IEEE Transactions on Automatic Control 66.1 (2020): 299-306.

[2] Pförtner, Marvin, et al. "Computation-aware Kalman filtering and smoothing." arXiv preprint arXiv:2405.08971 (2024).

[3] Bock, Christian, et al. "Online time series anomaly detection with state space Gaussian processes." arXiv preprint arXiv:2201.06763 (2022).

[4] Waxman, Daniel, and Petar M. Djurić. "A Gaussian Process-based Streaming Algorithm for Prediction of Time Series With Regimes and Outliers." 2024 27th International Conference on Information Fusion (FUSION). IEEE, 2024.