Venue: ICLR 2022 Conference

Reviews:

final_score: rating: 8: accept, good paper
review_content: summary_of_the_paper: The paper describes a method that combines variational Gaussian processes with a neural network that produces inducing points and associated variational parameters depending on the evaluation location.  The approach reduces training time compared to other sparse GP approaches, and has excellent performance in regression experiments (and good performance in classification experiments).  The method also does a good job at respecting the (exact) GP prior structure on the latent functions.
main_review: I liked the paper a lot once I understood what was going on.  But even with some familiarity with variational GPs, it took me a while to process this.  Most of my comments are individually minor, but they share a common theme that they are not easy to process without either reading out of order or knowing something about variational GPs already.  I think the idea here is quite interesting, and hope thinking about these points will make the paper more accessible!

Detailed comments:

- Top of p2: "This reduces the training cost to O(M^3)" -- perhaps say "O(M^3) per iteration"?

- The kernel is introduced as k(x,x') = E[f(x) f(x')], but perhaps it would be better to write E[(f(x)-m(x)) (f(x')-m(x'))] since the mean function is only "often set to zero" at this point?

- On page 2, we have "a popular covariance function k(.,.) is the squared exponential", but a Matern 3/2 kernel is used later in the experiments.  Perhaps it would be better to replace "squared exponential" with "Matern family" or the like?

- "In Titsias (2009), they do optimize q(u) in closed-form" -> "In Titsias (2009), they optimize q(u) in closed form"

- "The ELBO can be expressed as a sum..." -- there is no prior explicit reference to the evidence lower bound, though of course the bound itself has already been introduced.  It's probably worth a word or two here.

- "To achieve this, we consider a meta-point xtilde" -- I still don't completely understand the role of the meta-points (or what is meta about them)!  Some intuition would be useful.  Also, it is unfortunate that the meta-points are quickly given the same notation as the data points.  It seems like they end up coinciding in the algorithm as implemented, but are conceptually distinct from the data points?

- In Algorithm 1, I think standard fonts (vs the math spacing) would be better for ELBO, KL_div, and Log_lk.

- In related work: "it is introduced a new interpretation" -> "the authors introduce a new interpretation"; "that allow to consider" -> "that allow one to consider"; "it is also described a mechanism" -> "a mechanism is also described".

- As a general comment: the approximation for the pointwise posterior variance is clear, and this is the most important part of the uncertainty in many cases.  But it's unclear to me how one would compute the approximate posterior covariance between two points, which is also sometimes useful.
summary_of_the_review: Though some points of the presentation could be improved, I think this is an interesting idea illustrated by good experience.  I would like to see it appear in the conference.
review_title: null

final_score: rating: 3: reject, not good enough
review_content: summary_of_the_paper: The manuscript proposes a variant of a sparse variational approximation in order to render Gaussian process inference scalable. The central idea of the work is to use variationa inference and have a neural network predict the parameters of the variational distribution directly from the input data rather than optimising them using ELBO. Experiments on several substantial regression and classification datasets conclude the paper. The available code was not checked for this review.
main_review: 1) Strengths
 - a) The manuscript is rather well written and simple to understand.
 - b) The building blocks are established and well-understood.
 - c) The evaluation is done on relevant public datasets.
 - d) Code is available.

2)Concerns
 - a) It remains open which properties of the GP framework used in the motivation for the work remain intact after all.
 - b) The empirical evaluation is lacking at least three simple baselines.
     i) A neural network predicting the class probability rather than the variational parameters with sigmoid activation in the last layer. Along the same line, a regression network predicting mean and standard deviation (or log standard deviation). Both networks should be trained with maximum likelihood.
     ii) A sparse GP using FITC approximation.
     iii) A neural network predicting the inducing points rather than the variational parameters which is then used in combination with FITC or similar.

3) Remarks
 - a) The Abstract states that "for some learning tasks a large number of inducing points may be required" and that the proposed method is going to address this point. However, the conclusion on the "drastically reduced" number of inducing points is a little unfair as one could argue that all the complexity is shifted into the network for inducing point prediction. One should rather worry about number of parameters instead.
 - b) I have the impression that a direct neural network prediction of the class probability or the mean/std (possible mean/log(std)) would do a similar job as IDSGP and that the training would be much faster.
summary_of_the_review: The proposed algorithm is certainly an option that has not been explored before but the practical value is at least questionnable due to the missing baseline experiments.
review_title: null

final_score: rating: 3: reject, not good enough
review_content: summary_of_the_paper: Paper considers the sparse approximation Gaussian processes (GPs) which tries to scale GPs for large data sets. The authors in this paper propose a novel method (IDSGP) to improve the computational cost of sparse GPs. With the proposed method, the number of the inducing points is reduced drastically and only a small number of inducing points are estimated for each input point. The inducing points locations, as well as the parameters of the variational posterior distribution, are the outputs of a Deep neural network (DNN). 

First, they sample a mini-batch of the training data. For each data point in the mini-batch, their solution produces a small set of inducing points (and associated variational approximation) that are specific to each input data point.

It is an extension for Sparse within a Sparse GP (SWSGP) method Tran et al. (2020). For M inducing points, the SWSGP operates on a subset of H inducing points for each input point, with $H \ll M$, while maintaining a sparse approximation with M inducing variables. Here, instead of H nearest inducing points, the model estimates M inducing points and variational distribution parameters in a DNN. The claim is the proposed IDSGP method can provide competitive predictions while it improves the time complexity of the state-of-the-art sparse GPs baselines.
main_review: The main contribution of the paper is using a deep neural network (DNN) to find inducing points' locations and variational distribution parameters. In general form, it can reduce the computational cost of the original methods (SWSGP and SVGP) by training the variational parameters as the output layer of a DNN. But such a reduction in computational cost happens when it uses a small number of inducing points for each data point and small network architecture. However, there are concerns about the paper:

1- Novelty and Contribution: The paper has been inspired by the SWSGP method Tran et al. (2020). SWSGP samples a mini-batch from the input data set $\mathcal{D}$ and for each data point in the mini-batch, it finds H nearest inducing points. Instead of finding the nearest inducing points, IDSGP uses a deep Neural network and amortized variational inference Kingma and Welling (2014) and Shu et al. (2018) (available method) to find a small number of inducing points for each data point in mini-batch. The other parts of the algorithm in IDSGP are the same as the proposed algorithm in Tran et al. (2020). Indeed, using DNN to estimate variational parameters is not a novel approach. For instance, this network architecture is used in the encoder part of the conventional variational auto-encoder (VAE) (and also variational inference with the normalizing flow) where the parameters of the variational distribution (and flow parameters in normalizing flow) are estimated via DNN.  
    
 2-  Limitations of the Method: The limitations of the method are not properly discussed. The objective function has severely pathological properties which have not been described in the paper. For instance, consider the $n=M=1$ case, where n is the batch size, and a stationary kernel with unit variance. The optimized bound is exactly equal to a sum of GPs with one point each i.e. $\sum_i N(y_i, | 0, 1 + \sigma^2)$, and this does not depend on the prior at all. In fact,  as all priors are equivalent to the white noise prior. The SWSGP method does indeed have a similar problem, considering the $n=M=H=1$, but in the paper, it is explicitly described how the effective prior is modified. 
    
3- Fall back to the Prior:  In regions of the space without training data the extrapolation is determined by the DNN and there is no guarantee the predictions fall back to the prior. This is actually apparent in Figure 2 where VSGP is guaranteed to give the prior when there is no data or inducing points, whereas IDSGP visibly extrapolates in a way that doesn't seem sensible.
    
4- Performance: It is not clear why the IDSGP method can provide better prediction quality than SWSGP or SVGP while it only uses a small number of inducing points. Why is the DNN is better than the nearest point sampling in SWSGP, given that they maximize the same lower bound? SWSGP uses some nearest inducing points for each data point while IDSGP uses a small set of inducing points that have been estimated for each data point. Why does the second one have better performance? Since H-inducing points in SWSGP are selected based on input points, SWSGP is also an input-dependent SGP method. There is no reason to accept the inducing point set (in IDSGP) can lead to better predictive posterior and its improvement in prediction quality is ambiguous. For instance, in Appendix D.1, when the number of nearest inducing points (H) increases, SWSGP presents better results and its uncertainty is close to the standard GP. 
   
 5- Computational cost: The paper claims that since it uses a fewer number of inducing points, computational cost can be improved. However, it is not always true and the paper has not presented a related analysis about that. For instance, let consider the toy example in Section 5.1 and Figure 2. In this small data set, SWSGP is a little faster than IDSGP based on the computational cost discussions in both papers and experiment settings (i.e. the number of the layers and hidden units in DNN, number of inducing points in IDSGP, and the value of H and M in SWSGP). Besides, for a non-diagonal variational covariance matrix S, the output layer of DNN has $2M + \frac{M^2+M}{2}$ variables (M for inducing points, M for mean and $\frac{M^2+M}{2}$ for Cholesky matrix). When it needs more inducing points for each input point, the computational cost of DNN increases because of the dimension of the output layer and also network parameters (number of layers and hidden units).      


summary_of_the_review: In general, the paper is well-written. However, the proposed model has not been explained appropriately in Section 3. 

Section 3.2 is the main contribution of the paper which proposes the computational graph that estimates inducing points and variational parameters. Here it is necessary that the authors address the major issues of their method: why maximizing the lower bond (Equation (4)) in the proposed DNN leads to better results than the available baselines (e.g. SWSGP and SVGP which maximize the same objective function)?  Why does a small number of inducing points provide a better posterior than other baselines which find more inducing points (for instance SWSGP provides more inducing points but uses only H nearest inducing points for each input point)? Indeed, it is not clear how the batch permutation-invariance is achieved.  Something like a sum over the batch dimension is clearly not a good idea because the encoder is sensitive to the ordering of the batch. However, the paper does not mention how this issue is considered.

Besides, here are some omissions of related work that are not much discussed but could be mentioned, e.g. 
https://arxiv.org/pdf/1303.0383.pdf, 
https://mlg.eng.cam.ac.uk/zoubin/papers/aistats07localGP.pdf, 
https://www.dbs.ifi.lmu.de/~tresp/papers/bcm6.pdf,
http://proceedings.mlr.press/v84/salimbeni18a.html.
review_title: null

final_score: rating: 6: marginally above the acceptance threshold
review_content: summary_of_the_paper: This paper proposes to combine sparse GPs with a neural network architecture to compute the
inducing points locations associated to each input point. In particular, the paper employs a neural network to carry out amortized VI to compute the parameters of the approximate variational distribution
q modeling the posterior distribution associated to the values of the outputted inducing points. The inducing points are given by a mapping from the inputs provided by a neural network. The paper demonstrates improvement on training and prediction time while the proposed method is able to perform similar or better than the standard sparse GP approaches.

main_review: The paper has made a trade-off between training speed (pros) and throwing away the readily available gradient from the ELBO (cons).

Strengths:
* Optimizing the ELBO for Sparse GP using deep neural network is interesting for improving the running time (computational cost)

Weaknesses:
* Training a deep neural network may require a large amount of data as opposed to the standard sparse GP.

* It could be a drawback from the design in high-dimensional input setting. This can be inefficient to learn a mapping from one input (in high-dimensional) to multiple output (again in high-dimensional) of the inducing points. Learning such mapping may require a significantly large amount of training data.

* It is unclear about how to make use of the ELBO gradient in learning the neural network. On the one hand, existing variational inference techniques for sparse GP make use of such gradient when taking partial derivative with respect to the variational parameters. On the other hand, with deep neural network, one can utilize the automatic differentiation to perform back-propagation to learn the network (as a black-box). Note that optimizing the ELBO as a black-box (with deep neural network) using automatic differentiation may be inefficient if we have the gradient signal readily available from the ELBO.

Novelty:
* The paper follows Tran et al 2020 to derive the lower bound on the log-marginal likelihood.

* Note that existing works has considered using deep neural network to optimize the ELBO, (but not in the setting of sparse GP), such as [1]

Writing and presentation:
* The paper is well written and nicely presented.
* The related work is great.

Suggestions:
* One of the key advantages of GP against other supervised learning approaches is the uncertainty quantification. Therefore, additionally to presenting the classification and regression performance on UCI datasets, it is recommended to study/analyse/compare the uncertainty estimation. Without such uncertainty estimation, people can simply perform standard deep learning method to learn these classification and regression.


Minor points:
* in the last paragraph of the introduction, the paper claims to have different sets of inducing points associated to each input location. However, given the network, it seems that each input will only produce a single set of inducing point, not a set of including points?

* what is the dimension of Cholesky factor of the covariance matrix S=LL^T and the number of inducing points per data point Z(x_i)

* In Fig 1, there is a typo at the first column, last row which should be h_P^(1), instead of h_P^(L)


[1] Blundell, C., Cornebise, J., Kavukcuoglu, K., & Wierstra, D. (2015, June). Weight uncertainty in neural network. In International Conference on Machine Learning (pp. 1613-1622). PMLR.


summary_of_the_review: The paper proposes to improve the training time and flexibility of sparse GP approximation using deep neural network. While it makes sense to optimize the ELBO using deep neural network, the proposed framework may have certain limitation, such as learning the mapping from one (high-dim) input to multiple (high-dim) output can be inefficient, especially when we have limited data. Another suggestion is to demonstrate the uncertainty estimation of the proposed approach.

The paper is well written and presented. The related work section is well covered.

The paper is currently in the border-line.
review_title: null

final_score: rating: 5: marginally below the acceptance threshold
review_content: summary_of_the_paper: The paper proposes a variant of the method in Tran, 2021 for approximate Gaussian process regression. The major changes with respect to the classic Titsias, 2009 paper on variational Gaussian process regression are:
* The inducing inputs are data-dependent, in the sense that $Z=Z(x_i)$. 
* The inducing inputs, posterior mean, and posterior covariance are learned via amortized variational inference through learning a neural network.

The authors show some promising empirical evidence on commonly used regression and classification benchmarks.

main_review: 
## Significant Comments and Questions
*  Suppose I am interested in making a prediction of the joint distribution at a pair of inputs $x,x'$. How would I do this with you method? If I understand correctly, I would feed both $x$ and $x'$ into a neural net, which would output two distinct sets of inducing inputs, as well as distributions over them. Can you clarify how these can be combined to approximate, for example, the posterior covariance between $x$ and $x'$?
* Relatedly, it seems to me as though the posterior predictive distribution $f$ changes between train and test time. In particular, the posterior distribution depends on $p(\tilde{x})$ which differs at train and test time. Do you foresee this leading to any issues?
* Matthews et al 2015 established that maximizing the lower bound in Titsias, 2009 is equivalent to minimizing a KL-divergence between stochastic processes. This is essential to justifying maximization of the lower bound with respect to $Z$ as a form of variational inference. Given the augmented model, it is no longer clear to me that a similar argument applies. I am curious to hear if the authors have any thoughts as to whether it does or does not. If not, does the method risk over-fitting by optimizing $Z$? Relatedly, can the quality of inference (as measured by KL-divergence to the posterior) be shown to be monotonically increasing in $Z$?
* The authors should consider adding the closed form solution for approximate Gaussian process regression in experiments (Titsias, 2009). While more computationally expensive per iteration that the method of Hensman et al 2013, using the analytic form of the optimal posterior often improves inference quality signficicantly, and convergence of the bound proposed in Hensman et al 2013 can be slow.

## Minor Comments and Questions
* The authors repeatedly refer to the "training cost" when they mean the "cost per iteration of training". The training cost should incorporate how long it takes for the training of the method to converge.
* At several points you state that at least several hundred inducing points are needed for complicated problems. It seems as though several thousands of inducing points might be a more appropriate statement. Generally using several hundred inducing points would not represent a major computational obstacle. The lower bounds on the number of inducing points needed  for the KL-divergence between approximate posterior and prior to not be large in Burt et al, 2020 might also be of interest.
* On page two, $k(x,x') = \mathbb{E}[f(x)f(x')]$ is only correct if the mean is $0$. While you state that you generally assume this, it should be more clear this is needed in this definition.
* The last paragraph in section 2.1 is copied almost word for word from an earlier paragraph.
* In experiments, when results are compared, are model parameters selected via ELBO maximization (approximate maximum marginal likelihood)?

## References
``On Sparse variational methods and the Kullback-Leibler divergence between stochastic processes'' Matthews, Hensman, Turner, Ghahramani. 2015.

``Convergence of Sparse Variational Inference in Gaussian Processes
Regression'', Burt, van der Wilk, Rasmussen.
summary_of_the_review: ## Summary of Review
The method is interesting and the empirical results seem promising. I do have some concerns that the posterior distribution defined over $f$ is not consistent. Perhaps this is not essential if only marginal predictions are needed but it does seem to be a significant drawback in cases when covariances or samples from the approximate posterior are needed.  This also makes me slightly concerned that the objective function is not quite as well founded as that of Titsias, 2009 (see the discussion in Matthews et al 2015). I think some of these critiques apply equally to Tran et al, and I am curious to hear the authors thoughts on these limitations. 

The writing is generally reasonably clear, though aspects could be improved. While I have rated the paper as marginally below the acceptance threshold, I think it is quite close to the threshold.
review_title: null

Meta-review: Three reviewers recommend Reject. Two reviewers recommend Accept although do not champion the paper. I believe the paper develops an interesting idea to better estimate the location of the inducing inputs in sparse GP models. However, I still think the paper would benefit from another careful revision and therefore I do not recommend Acceptance at this stage. I agree with reviewers that 1) currently the method is unable to estimate the covariance between two data points. This is important in applications of GPs for uncertainty quantification such as Bayesian optimisation. For example, including a BayesOpt example would clearly strengthen the paper. 2) empirical evaluation lacks simple baselines, e.g. Titsias (2009). The authors claim that Titsias (2009) does not scale and that's why they don't care for it. Even if this is true, including an example that helps to better compare against this method at a different scale might strengthen the model proposed here.
Decision: Reject