We thank the reviewers for their thoughtful comments and suggestions. We agree with most of the remarks and will incorporate them in the final version of the paper. In the following, we will address some of the raised questions.$> R1: the analysis is restricted to the squared exponential kernel
While in the paper we provided analysis for the widely employed squared exponential kernel, similar results can be obtained for other differentiable kernels. For example, the analysis can straightforwardly be extended to linear, periodic, rational quadratic, Matérn kernels, and any sums and products of these kernels.
> R1: the derived error bound cannot be used as it makes the quadrature too complex (Sec 5.3) […] Is the bound too conservative?
For short time horizons, the quadratures required for the error bound can be computed with little computational effort. For long term horizons, e.g., 100s or 1000 discrete time steps as in our examples, more computational effort is necessary. Our approach provides strong deterministic stability guarantees and, thus, has to consider the worst case scenario using conservative error bounds. In practice, accurate stability region estimates can be obtained quickly with coarse quadrature rules, while more computational effort allows to obtain stability guarantees. To our best knowledge, we present the first approach that is able to provide such guarantees.
> R2: What was missing was a time budget for the algorithms. Can they run online […]?
> Can I use them to prove stability for an offline identified system on my laptop?
For the implementation, there is a trade-off between memory usage and computation time. While the algorithm runs on a current laptop, it is significantly faster when more memory space is available. In this case computation of the stability region takes a few (e.g., 3-4) minutes. We will add details on computation time for all experiments in the final version.
> R1: It is not clear if the paper requires any assumptions on the control policy.
> R3: Sections 3 and 4 mention the constraints required on controllers […] only in the following passage […] I would like to have seen more elaboration upon this point throughout the paper: why are these conditions sufficient? Are they necessary?
We agree that the assumptions on the control policy should be more prominent and will highlight them in the final version. As stated in the paper (line 189 f.), the controller must be a differentiable function of the current state. Differentiability is necessary to prove Lemma 1 and to ensure that the GP predictive distribution depends smoothly on the state, enabling numerical quadrature. Dependency on the current state only is required for the (inductive) invariant set analysis. No further assumptions on the controller are required.
> R3: My major concern about the approach is its practical applicability in establishing the stability of real systems […] Intuitively, the grid should fail to scale to high dimension […] This is cursorily addressed by the authors in the following passage […] This last sentence needs far greater support.
Any currently-thinkable deterministic approach will at some point be crushed by scalability problems. Our current implementation enables first stability guarantees, even for reasonably small real world systems. 
Any quadrature can directly be incorporated in our approach. Only basic full product grid quadrature approaches are limited to few dimensions. Sophisticated quadratures exist, that scale polynomially with the problem dimension and work for up to 20 dimensions (which allow to handle a substantial amount of control problems), for example:
- Sparse grid quadrature [1,2,3,4]: smart combination of 1D rules to a multivariate rule of desired exactness, avoids redundancy of full product grid.
- Non-tensor rules [5,6,7]: computed directly and, thus, by construction more efficient than full product rules. An example is Stroud quadrature [5], which scales quadratically with the number of dimensions.
We currently run experiments with these approaches and will add the results in the final version.
In addition to sophisticated quadrature rules, dimensional reduction techniques can be employed. We already included a dimensionality reduction approach that exploits the often lower dimensional manifold of the system dynamics. We will provide details of the approach and evaluation in the final version of the paper.
[1] Estimation with Numerical Integration on Sparse Grids; Heiss, Winschel
[2] Smolyak cubature of given polynomial degree with few nodes for increasing dimension; Petras
[3] Sparse grids; Bungartz, Griebel
[4] High dimensional Integration of Smooth Functions over Cubes; Novak, Ritter
[5] Numerical Evaluation of Multiple Integrals; Hammer, Stroud
[6] A numerical algorithm for the construction of efficient quadrature rules in two and higher dimensions; Xiao, Gimbutas
[7] Extensions of Gauss Quadrature Via Linear Programming; Ryu, Boyd