Neural Splines: Fitting 3D Surfaces with Infinitely-Wide Neural Networks

Francis Williams, Matthew Trager, Joan Bruna, Denis Zorin

Introduction

Estimating a 3D surface from a scattered point cloud is a classical and important problem in computer vision and computer graphics. In this task, the input is a set of 3D points sampled from an unknown surface and, possibly, normals of the surface at those points. The goal is to estimate a representation of the complete surface from which the input samples were obtained, for example, a polygonal mesh or an implicit function. This problem is challenging in practice: It is inherently ill-posed, since an infinite number of surfaces can interpolate the data. Furthermore, the input 3D points are often incomplete and noisy, as they are acquired from range sensors such as LIDAR, structured light, and laser scans. Ideally, the recovered surface should not interpolate noise but preserve key features and surface details.

Many early surface reconstruction techniques consider a kernel formulation of the surface reconstruction problem, using translation-invariant kernels such as biharmonic RBFs , Gaussian kernels , or compactly supported RBFs . Currently, the most widely used method for surface reconstruction is Screened Poisson Surface Reconstruction , which solves a variant of the Poisson equation to find an implicit function whose zero-level set produces an approximating surface. We show in this paper that this method can also be viewed as a kernel method for a particular choice of kernel.

More recently, many papers have used neural networks to represent an implicit function or a local chart in a manifold atlas as a means of reconstructing a surface . These methods can be integrated into a data-driven learning pipeline or directly applied in the so called “overfitting” regime, where a massively overparameterized (\ie, more parameters than input points) neural network is fitted to a single input point cloud as a functional representation for a surface. Empirical evidence has shown that these methods enjoy some form of “implicit regularization” that biases the recovered surface towards smoothness. Moreover, employing early stopping in gradient descent can prevent these approaches from interpolating noise.

Under certain parameter initializations, infinitely-wide neural networks de facto behave as kernel machines , defining Reproducing Kernel Hilbert Spaces (RKHS), whose kernel is obtained by linearizing the neural network mapping around its initialization. While the kernel regime simplifies non-linear neural network learning into a convex program and provides a simple explanation for the successful optimization of overparameterized models, it cannot explain the good generalization properties observed for high-dimensional problems due to the inability of the RKHS to approximate non-smooth functions . However, the situation for low-dimensional problems such as surface reconstruction is entirely different, and RKHS can provide powerful characterizations of regularity. In this context, shows that in the univariate case the RKHS norm associated with a wide ReLU network is a weighted curvature, and leads to cubic spline interpolants. In higher dimensions, similar (albeit more complex) characterizations of the RKHS norm exist (see also Proposition 2). In order to assess the benefits of neural networks on such low-dimensional problems, it is thus important to first understand their linearized counterparts, given by their associated kernel machines.

In this work, we demonstrate that in fact kernels arising from shallow ReLU networks are extremely competitive for 3D surface reconstruction, achieving state-of-the-art results: outperforming classical methods as well as non-linear methods based on far more complex neural network optimization.

Kernels provide many advantages over neural networks in our context: (i) They are well understood theoretically. (ii) Kernel regression boils down to solving a linear system, and avoids gradient descent that suffers from slow convergence. (iii) Kernel-based interpolants are represented using a number of parameters that is linear in the size of the input, whereas overparameterized neural networks require many more parameters than points. (iv) The inductive bias of kernel methods can be characterized explicitly via the RKHS norm (Section 3.3). (v) Kernel methods lend to scalable and efficient implementations (Section 3.5). We provide explicit expressions for two kinds of infinite-width ReLU kernels, their derivatives, and their corresponding RKHS norms. We further argue that these kernels can be viewed as a multidimensional generalization of cubic spline interpolation in 1D. Moreover, we show that Poisson Surface Reconstruction can itself be viewed as a kernel method and give an expression for its RKHS norm, suggesting that kernels are a broad framework which enable the rigorous understanding both traditional and modern surface reconstruction techniques.

Methods for 3D surface reconstruction can mostly be divided by their choice of surface representation: These are (i) the zero level set of a volumetric scalar function , (ii) the fixed point of a projection operator onto locally fitted patches , (iii) a mesh connecting the input points , (iv) a collection of local parametric maps , or (v) the union of parametric shapes . For a recent up-to-date survey of surface reconstruction techniques, we refer the reader to .

The “random feature” kernels used this paper arise from training the top-layer weights in two-layer networks of infinite width and were described in . More recently, a lot of work has focused on a different kernel, known as the “neural tangent kernel” , that linearly approximates the training of both layers of a neural network. We chose to use random feature kernels since, when the input dimension is small, typical initializations in neural networks lead to mainly training the top layer weights. More broadly, the function spaces associated with shallow neural networks were studied in .

Neural Spline Formulation

We begin by introducing the Neural Spline family of functions through the lens of finite-width shallow neural networks in Section 2.1, and turn to the infinite-width limit formulation in Section 2.2. While our approach is simple and our derivations straightforward, we differ slightly from a standard kernel regression setting since we employ a multi-dimensional kernel that includes the gradient of the fitted function. This formulation allows us to fit points and normals simultaneously.

We first assume that the model f(x;θ)f(x;\theta) is a two-layer ReLU neural network with mm neurons, but we keep the bottom layer weights fixed from initialization:

In typical initialization schemes for neural networks (\eg, Kaiming He initalization ), the weights of each layer are initialized with a variance that is inversely proportional to the number of inputs of that layer. Our choice of fixing (ai,bi)(a_{i},b_{i}) is motivated by the fact that using a standard two-layer network, for a given target model and as m→∞m\rightarrow\infty, only the top layer weights tend to vary throughout training .

This minimizer is given explicitly by θ∗=W†δ\theta^{*}=W^{\dagger}\delta where W†W^{\dagger} denotes the Moore-Penrose pseudo-inverse of WW and

where for compactness we used wi=(ai,bi)w_{i}=(a_{i},b_{i}), φwi(x)=[aix+bi]+\varphi_{w_{i}}(x)=[a_{i}x+b_{i}]_{+} and ∇xφwi(x)=1[aix+bi]+ai\nabla_{x}\varphi_{w_{i}}(x)={\bf 1}[a_{i}x+b_{i}]_{+}a_{i}.

2 Infinite-Width Kernel

As the number of neurons mm tends to infinity, the kernel (4) converges to K∞(x,x′)K_{\infty}(x,x^{\prime}) defined by

The uniform distribution (8) corresponds to the default initialization of linear layers in PyTorch and, as we argue in Section 3.1, it leads to a direct generalization of cubic spline interpolation. However, the Gaussian initialization actually leads to simpler analytical expressions for K∞K_{\infty} and produces almost the same results. We refer to Appendix B for a discussion and comparison of the two distributions.

Discussion

The expression assuming −k≤x′<x≤k-k\leq x^{\prime}<x\leq k is obtained by swapping xx and x′x^{\prime}. For fixed xx, the map Kspline(x,⋅)K_{\rm spline}(x,\cdot) is piecewise cubic and twice continuously differentiable (C2C^{2}). This implies that kernel regression with (10) yields cubic spline interpolation. Applying the Neural Spline objective (1) with derivative constraints at the samples in d=1d=1 also yields a piecewice cubic interpolant, although this curve is in general only C1C^{1}. For other distributions of aa and bb, the kernel is no longer cubic, but the norm in the RKHS is a weighed norm of curvature (see for details). In this sense, our approach with the initialization (8) can be viewed as a multi-dimensional version of spline interpolation.

2 Regularization and Robustness to Noise

To deal with noisy data, we can optionally add a simple regularizer term to our formulation that corresponds to penalizing the RKHS norm of the interpolant (“kernel ridge regression”). Concretely, we replace (6) with

for j=1,…,sj=1,\ldots,s. The regularizer term affects the spectrum of the Gram matrix K∞(xj,xi)K_{\infty}(x_{j},x_{i}) by smoothing its smallest eigenvalues, with a similar effect to early stopping in gradient descent (see \eg ). Figure 2 shows examples of applying this regularizer on 2D and 3D problems.

3 Inductive Bias of Interpolants

Our interpolant f∗f^{*} belongs to the Hilbert space given by

where τ(w)\tau(w) is a measure over the weights w=(a,b)w=(a,b) (\eg, (8) or (9)), and c(a,b)c(a,b) can be viewed as a continuous analog of the outer-layer weights cic_{i} of the finite network given by (2). The inductive bias of Neural Splines is thus determined by the RKHS norm in (12) since our method outputs the interpolating function which minimizes that norm. As noted in , if f(x)=∫c(w)φw(x)dτ(w)f(x)=\int c(w)\varphi_{w}(x)d\tau(w) and dτ(w)=dτa(a)dτb(b)d\tau(w)=d\tau_{a}(a)d\tau_{b}(b) then, by differentiating twice, we note that the Laplacian of ff is given by

By comparing (12) and (13), we see that the RKHS norm and the Laplacian of ff are closely related. More precisely, the Laplacian Δf(x)\Delta f(x) is the Dual Radon Transform of c(a,b)c(a,b). Under certain assumptions, (13) can be inverted, yielding an explicit expression for cc in terms of Δf\Delta f. Intuitively, this shows that bounding the RKHS norm imposes a constraint on the Laplacian of ff, and thus encourages ff to be smooth. We report the following statement from and we refer to Appendix D for more details.

Let f(x)=∫c(a,b)[ax+b]+dτ(a,b)f(x)=\int c(a,b)[ax+b]_{+}d\tau(a,b), and the constant γd=12(2π)d−1\gamma_{d}=\frac{1}{2(2\pi)^{d-1}}. If we assume that c(a,b)=c(−a,−b)c(a,b)=c(-a,-b) holds, then

where R{f}(a,b)\mathcal{R}\{f\}(a,b) is the Radon Transform of ff.

4 Poisson Surface Reconstruction as a Kernel

We cast Screened Poisson Surface Reconstruction in kernel form to facilitate comparisons. In its simplest form, Poisson reconstruction, extracts the level set of a smoothed indicator function determined as the solution of

The vector field VV is obtained by interpolating the normals using a fixed-grid spline basis and barycentric coordinates of the sample points with respect to the grid cell containing it. This is equivalent to using a non-translation invariant non-symmetric locally-supported kernel KB(z,x)K_{B}(z,x):

where KB(z,x)=∑jB1(x−cj)Bn(z−cj)K_{B}(z,x)=\sum_{j}B_{1}(x-c_{j})B_{n}(z-c_{j}).

To study the qualitative properties of this kernel, we replace KB(z,x)K_{B}(z,x) with a radial kernel Bn1(∣z−x∣)B^{1}_{n}(|z-x|) (see Appendix C) which has qualitatively similar behavior (see Figure 3). Since both Bn1B^{1}_{n} and the Laplace kernel 1∣x−z∣\frac{1}{|x-z|} are radial functions, their convolution is also radial, yielding a translation-invariant radial approximation KPRapproxK_{\text{PR}}^{\text{approx}} of KPRK_{\text{PR}}:

The RKHS norm of the corresponding to the approximate Poisson kernel KPRapproxK_{{\rm PR}}^{{\rm approx}} is

where F[⋅]\mathcal{F}[\cdot] is the Fourier transform.

We discuss the kernel formulation of Poisson Reconstruction in more detail in Appendix C.

5 Fast and Scalable Implementation

We provide a fast and scalable implementation of Neural Spline kernels based on FALKON , a recently-proposed solver for kernel-ridge-regression which runs in parallel on the GPU. While naïve kernel ridge-regression with NN points requires solving and storing an N×NN\times N dense linear system, FALKON uses conjugate gradient descent requiring only O(N)\mathcal{O}(N) storage, and N\sqrt{N} convergence (though in practice we find that even for very large inputs, we converge in fewer than 10 iterations). To speed up convergence, FALKON can optionally store an M×MM\times M preconditioner matrix in CPU memory (where M≪NM\ll N, see paragrph below). To maximize performance and reduce memory overhead, we rely on KeOps to evaluate kernel matrix-vector products symbolically on the GPU, which means our implementation uses only a small constant amount GPU memory and can be readily used on commodity hardware. Section 4.4 compares the performance of our implementation against other state of the art surface reconstruction techniques. We note that in principle, low-dimensional kernel methods can be accelerated using fast multipole-based approaches (in particular, in the context of 3D surface reconstruction this was used in ); this yields optimal O(N)O(N) time complexity, for dense matrix-vector multiplication.

Full kernel ridge regression predicts a function which is supported on every input point xix_{i} as in (6), requiring NN coefficients to store the resulting function. We rely on Nyström sampling to instead produce a kernel function which is supported on a small MM-sized subset of the input points (while still minimizing a loss on all the points). This is equivalent to approximating the kernel matrix with a low rank linear system. To choose Nyström samples, we leverage the geometric nature of our problem and select these by downsampling the input point cloud to have a blue-noise distribution using Bridson’s algorithm . We demonstrate the effect of varying the number of Nyström samples qualitatively in Figure 4.

Experiments and Results

We now demonstrate the effectiveness of Neural Splines on the task of surface reconstruction. For all the experiments in this section, we used the analytical form of the kernel (5) with Gaussian initialization (8). Appendix A.7 compares the uniform (8) and Gaussian (9) kernels, showing almost no measurable difference in the reconstructions produced by either. We compare the empirical and analytical kernels in detail in Appendix A.6. An implementation of Neural Splines is available at https://github.com/fwilliams/neural-splines.

We performed a quantitative evaluation on a subset of the Shapenet dataset to demonstrate that the inductive bias of Neural Splines is particularly effective for reconstructing surfaces from sparse points. We chose 1024 random points and normals sampled from the surface of 20 shapes per category across 13 categories (totalling 260 shapes). Using this dataset, we compared our method against Implicit Geometric Regularization (IGR) , SIREN , Fourier Feature Networks , Biharmonic RBF (Biharmonic) , SVM surface modelling (SVR) , and Screened Poisson Surface Reconstruction (Poisson) . The first three techniques are modern neural network based methods, while the latter three techniques are classical methods based on kernels or solving a PDE. As criteria for the benchmark, we consider the Intersection over Union (IoU) and Chamfer Distance between the reconstructed shapes and the ground truth shapes. The former metric captures the accuracy of the predicted occupancy function, while the latter metric captures the accuracy of the predicted surface. Under both metrics, our method outperforms all other methods by a large margin. Table 1 reports quantitative results for the experiment and Figure 5 shows visual results on a few models. We report per-category results in Appendix A.4 and show many more figures in Appendix A.3.

Several of the above methods have parameters which can be tuned to increase performance. To ensure a fair comparison, we ran parameters sweeps for each Shapenet model for these methods, reporting the maximimum of each metric under consideration (see Appendix A.1 for a detailed description). For our method we did not tune parameters. We used no regularization and 1024 as Nyström samples for all models.

2 Surface Reconstruction Benchmark

We evaluated our method on the Surface Reconstruction Benchmark which consists of simulated noisy range scans (points and normals) taken from 5 shapes with challenging properties such as complex topologies, sharp features, and small surface details. We evaluate our method against Deep Geometric Prior (DGP) , Implicit Geometric Regularization (IGR) , SIREN , and Fourier Feature Networks (FFN) . We remark that DGP establishes itself as superior to a dozen other classical methods on this benchmark and IGR furhter outperforms DGP. As in and , Table 2 reports the Hausdorff (dHd_{H}) and Chamfer (dCd_{C}) distances between the reconstruction and ground-truth. We also report the one sided Hausdorff (dH⃗d_{\vec{H}}) and Chamfer (dC⃗d_{\vec{C}}) distances between the scan and the reconstruction, which measures how much the reconstructions overfits noise in the input. Our reconstructions are quantitatively closer to ground truth on all but one model. Figure 6 shows visual examples of a few models from the benchmark. All models are shown in Appendix A.3. As with the Shapenet benchmark, we did parameter sweeps for those methods which have tunable parameters, choosing the best model for each metric under consideration. See Appendix A.2 for details.

3 Large Scale Reconstruction of Full Scenes

Figure 7 shows a full scene consisting of 9 million points reconstructed using our method. This is the same scene used in , and contains many thin features that are difficult to reconstruct (\eg, curtains, plant, and lampshade). To generate the input point cloud for this experiment, we densely sampled a mesh extracted from a 204832048^{3} occupancy grid of the scene. To perform the reconstruction, we subdivided the space into 8x8x8 cells containing between 10k and 500k samples and reconstructed each cell interdependently using up to 15k Nyström samples. The whole process takes 1.5 hours on a machine with an NVIDIA 1080Ti GPU.

4 Timing and Performance

Table 3 compares the average running time and GPU memory usage of our method and others when reconstructing point clouds from the Surface Reconstruction Benchmark described in Section 4.2. To ensure a fair comparison we ran all neural network based methods for 5000 iterations. We do not report exact CPU memory usage for the methods since it is hard to measure but we remark that, by observation, system memory usage never exceeded 3GiB for any of the methods. Appendix A.5 reports the running times and memory usages for individual models in the benchmark.

5 Our Method versus Neural Networks

Implicit Geometric Regularization demonstrates that ReLU networks have a natural inductive bias making them good at reconstructing shapes. However, methods based on ReLU networks suffer from slow convergence (see Figure 8). SIREN drastically improves convergence speed by replacing ReLU with sinusoidal activations and a clever initialization, however, SIREN occasionally underfits on sparse inputs and requires an additional loss computed on points in the volume around the shape to prevent reconstruction artefacts. This loss increase runtimes, and suggest that SIREN’s inductive bias may not be ideal for sparse reconstruction tasks. By projecting input points onto random Fourier features, Fourier Feature Networks (FFNs) present a principled approach rooted in Kernel methods to control the inductive bias of the solution as well as convergence speed. While FFNs are capable of producing high quality solutions, they are sensitive to the choice of feature distribution when data is sparse (see Figure 9) and thus require tuning to work well. In contrast, our technique converges in seconds (Section 4.4), has a well suited inductive bias for shape representation (Section 3.3), and requires minimal parameter tuning.

Conclusion and Future Work

We have shown that Neural Spline kernels arising from infinitely wide shallow ReLU networks are very effective tools for 3D surface reconstruction, outperforming state-of-the-art methods while being computationally efficient and conceptually simple. In a sense, our work bridges the gap between traditional reconstruction methods and modern neural network based methods by leveraging the deep connection between neural networks and kernels.

We remark that our kernel formulation is fully differentiable. In the future, we hope to integrate Neural Splines into deep learning pipelines and apply them to other 3D tasks such as shape completion and sparse reconstruction. On the theory side, we would like to investigate and compare the approximation properties of different kernels (in particular those arising from infinite width sinusoidal networks) in the context of 3D reconstruction.

This work is partially supported by the Alfred P. Sloan Foundation, NSF RI-1816753, NSF CAREER CIF 1845360, NSF CHS-1901091, and Samsung Electronics.

References

Appendix A Additional Experiments

Many of the methods we compared against on Shapenet have tunable parameters which can drastically alter the quality of reconstructed outputs. To ensure a proper comparison, we ran sweeps over these parameters where appropriate choosing the best reconstruction for each model under both metrics (Chamfer and IoU). For our method we did no parameter sweeps, using no regularization and 1024 Nyström samples for each model in the dataset. We describe the exerimental methodology for each method in the benchmark below.

We trained each model for 5k iterations with Adam and a learning rate of 0.001 using the same parameters and architecture as proposed in the original paper. We included the normals in the loss with the parameter τ\tau set to 1. The Eikonal regularization term λ\lambda was set to 0.1. While IGR can slightly improve by using a very large number of iterations. doing so is prohibitively slow over many models. Figure 8 motivates our choice of iterations, demonstrating only a slight improvements between 5k and 100k iterations (the latter which required 2 hours of fitting on a NVIDIA-1080-Ti GPU).

Screened Poisson Surface Reconstruction [31]

We considered every possible combination of the following parameters: the octree depth in ,thenumberofpointsperleafin, the number of points per leaf in, the point weight (which controls the degree to which the method interpolates the input) in [4.0,100.0,1000.0][4.0,100.0,1000.0]. Since all the shapes in the benchmark are watertight meshes, we used Dirichlet boundary constraints for the reconstruction.

SIREN [40]

We trained used a SIREN network with 5 hidden layers each containing 256 neurons for 5000 iterations, using a learning rate of 1e-4 with Adam. We used the same loss for shapes as in the original SIREN paper, sampling an additional NN (where NN is the number of input points) points in an axis aligned bounding box whose diagonal is 10% larger than the object’s bounding box. We minimize the same loss for shapes as the SIREN paper (see section 4.2 in the original SIREN paper).

Fourier Feature Networks [42]

We used an 8-layer ReLU MLP with 256 Fourier features sampled from a Gaussian distribution. This is the same architecture as the shape representation experiment in the original paper. For each model in the benchmark, we did a parameter sweep on the variance σ\sigma of the Gaussian distribution, considering σ∈{0.1,0.25,0.5.0.6,0.7,0.8,0.9,1.0,1.25,1.5,3.0}\sigma\in\{0.1,0.25,0.5.0.6,0.7,0.8,0.9,1.0,1.25,1.5,3.0\}. The range of parameters was chosen by empirical verification on 3 models from the airplanes, benches, and cars categories.

SVR [22]

As in the original paper, we use a Gaussian kernel to perform support vector regression. To generate occupancy samples, we augmented the input points with an ”inside“ and ”outside“ point by perturbing them by ±ϵ\pm\epsilon along the normal at that point. We used ϵ=0.01\epsilon=0.01 forall the models. For each model we did a joint parameter sweep over the regularization parameter C∈{1.0,0.1,0.01,0.001,0.0001}C\in\{1.0,0.1,0.01,0.001,0.0001\} and the variance parameter σ∈{0.002.0.001,0.0004,0.0002,0.0001}\sigma\in\{0.002.0.001,0.0004,0.0002,0.0001\}. The range of parameters was chosen by empirical verification on 3 models from the airplanes, benches, and cars categories.

Biharmonic RBF [9]

To generate occupancy samples, we augmented the input points with an ”inside“ and ”outside“ point by perturbing them by ± ϵ\pm\ \epsilon along the normal at that point. We used ϵ=0.01\epsilon=0.01 for all the models. The biharmonic function is very simple ϕ(r)=r\phi(r)=r where r=∥xi−xj∥r=\|x_{i}-x_{j}\| and does require tuning parameters.

A.2 Detailed Description of Surface Reconstrucion Benchmark Experiment

A.3 Additional Figures

Figures 10 and 11 show at least one model reconstructed from each ShepeNet category using our method, Implicit Geometric Regularization (IGR) , SIREN , Fourier Feature Networks (FFN) , Screened Poisson Surface Reconstruction , Biharmonic RBF , and Support Vector Regression (SVR) . Figure 12 shows the reconstructions of all the models from the Surface Reconstruction Benchmark using our method, IGR , SIREN , and Fourier Feature Networks .

A.4 Quantitative Results Per ShapeNet Class

Tables 4 show the per ShapeNet category IoU and Chamfer distance statistics for the benchmark described in Section 4.1.

A.5 Per Model Performance Numbers

Table 5 shows the runtime and GPU usage required to reconstruct each model in the Surface Reconstruction Benchmark. For our model, we used 15k Nyström samples and a regularization of 1e-11. We do not report CPU memory usage since it is hard to profile exactly, however we observed that none of the methods used more than 4GiB of CPU memory. All timings were done on a machine with a single NVIDIA-V100 GPU with 16GiB of VRAM, 32GiB of CPU RAM, and and an 8 core Intel Xeon processor.

A.6 Empirical versus Analytical Kernel

Figure 14 compares results using the empirical kernel with mm neurons and using the analytical kernel. Figure 13 shows the convergence of the empirical Kernel to the analytic one as the number mm of neurons grows.

A.7 Quantitative Comparison between Gaussian and Uniform Kernels

Table 6 shows a quantitative comparison between Neural Splines using the Gaussian initialization (9) and the Uniform initialization (8) on the benchmark described in Section (4.1). The results in both cases are very close to each other in both Chamfer distance and in IoU.

Appendix B Derivation of the Infinite Width Kernels

We derive an explicit expression for the kernel K∞K_{\infty} in (5) in the case of uniform initialization (8). We first prove the following Lemma which we will use in our calculations.

where Fg=∫02πg(ψ)F(∥x∥cos⁡(ψ),∥x′∥cos⁡(ψ−α))dψF_{g}=\int_{0}^{2\pi}g(\psi)\mathcal{F}(\|x\|\cos(\psi),\|x^{\prime}\|\cos(\psi-\alpha))d\psi (for g=1,cos⁡,cos⁡2g=1,\cos,\cos^{2}, etc.), Q∈SO(d)Q\in SO(d) is such that Qx=(0,…,∥x∥,0)TQx=(0,\ldots,\|x\|,0)^{T},Qx′=(0,…,∥x′∥cos⁡(α),∥x′∥sin⁡(α))TQx^{\prime}=(0,\ldots,\|x^{\prime}\|\cos(\alpha),\|x^{\prime}\|\sin(\alpha))^{T} and

We now consider the integral ∫aaTF(aTx,aTx′)dΩ\int aa^{T}\mathcal{F}(a^{T}x,a^{T}x^{\prime})d\Omega. For any two indices i≤j≤d−2i\leq j\leq d-2, we have

This is now a product of d−1d-1 one-dimensional integrals. Since ∫0πsin⁡s=0\int_{0}^{\pi}\sin^{s}=0 if ss is odd, we have that the integral (18) vanishes if i≠ji\neq j. If instead i=ji=j, we use the fact that

This proves the diagonal part in our expression for ∫aaTF(aTx,aTx′)dΩ\int aa^{T}\mathcal{F}(a^{T}x,a^{T}x^{\prime})d\Omega. All remaining terms as well as the two integrals ∫aF(aTx,aTx′)dΩ\int a\mathcal{F}(a^{T}x,a^{T}x^{\prime})d\Omega and ∫F(aTx,aTx′)dΩ\int\mathcal{F}(a^{T}x,a^{T}x^{\prime})d\Omega follow from very similar (and slightly simpler) calculations. ∎

We now apply Lemma 5 to compute the kernel K∞K_{\infty} with the uniform initialization (8).

where α=arccos⁡(x⋅x′∥x∥∥x′∥)\alpha=\arccos\left(\frac{x\cdot x^{\prime}}{\|x\|\|x^{\prime}\|}\right), \tau=\arctan\bigg(\frac{||x||-||x^{\prime}||\cos(\alpha)}{||x^{\prime}||\sin(\alpha)}\bigg{missing}), Q∈SO(d)Q\in SO(d) is such that Qx=(0,…,∥x∥,0)TQx=(0,\ldots,\|x\|,0)^{T},Qx′=(0,…,∥x′∥cos⁡(α),∥x′∥sin⁡(α))TQx^{\prime}=(0,\ldots,\|x^{\prime}\|\cos(\alpha),\|x^{\prime}\|\sin(\alpha))^{T} and

The idea is to compute the integral with respect to the bias term bb and then split the result into homogeneous expressions where Lemma 5 can be applied. In particular, assuming that −k≤s,t≤k-k\leq s,t\leq k:

B.2 Gaussian Initialization

The Gaussian initialization (9) yields the following simpler formula for K∞K_{\infty}. The first term is well known and is derived in . The second and third terms are easily derived by taking derivatives of the first term.

If a∼N(0,Idd−1)a\sim\mathcal{N}(0,Id_{d-1}) and b∼N(0,1)b\sim\mathcal{N}(0,1), then

Appendix C Poisson Surface Reconstruction Kernel

In its simplest form, Poisson reconstruction of a surface , extracts the level set of a smoothed indicator function determined as the solution of

The vector field VV is obtained by interpolating the normals using a fixed-grid spline basis and barycentric coordinates of the sample points with respect to the grid cell containing it. This is equivalent to using a non-translation invariant non-symmetric locally-supported kernel KB(z,x)K_{B}(z,x):

Then KB(z,x)=∑jB1,3(x−cj)Bn,3(z−cj)K_{B}(z,x)=\sum_{j}B_{1,3}(x-c_{j})B_{n,3}(z-c_{j}). where only 8 terms corresponding to the vertices cjc_{j} of the grid cube containing xx are nonzero. This yields the following expression for the kernel corresponding to Poisson reconstruction,

, the convolution of the Laplacian kernel 1/∣x−z∣1/|x-z| and the gradient of KBK_{B}. Using the identity ∇(f∗g)=(∇f∗g)\nabla(f*g)=(\nabla f*g), we can write this as the gradient of KPR(x,x′)gK_{\text{PR}}(x,x^{\prime})_{g}, defined as

To make it easier to understand the qualitative behavior of the kernel, replacing KB(z,x)K_{B}(z,x) with a radial kernel Bn1(∣z−x∣)B^{1}_{n}(|z-x|), with qualitatively similar behavior (see Figure 3) yields a translation-invariant radial approximation KPRapproxK_{\text{PR}}^{\text{approx}} of the kernel KPRK_{\text{PR}}, as the convolution of two radial kernels is a radial function.

As both Bn1B^{1}_{n} and the Laplace kernel are radial functions, their convolution is also radial. It can be expressed in a more explicit form using the relation between Fourier and Hankel transforms for radial functions. For n=3n=3, the Hankel transform is related to Fourier transform by

The Hankel transform is an involution, so the relationship for the inverse Fourier transform is similar. Writing g∗h=F−1[F[g]F[h]]g*h=\mathcal{F}^{-1}[\mathcal{F}[g]\mathcal{F}[h]], we obtain the expression for the radial convolution in terms of one-dimensional integrals,

where we use Hr[1/r]=1/s\mathcal{H}_{r}[1/r]=1/s. and Kgapprox(x)K^{\text{approx}}_{g}(x) is just the gradient of this, i.e., a derivative times ∣x∣/x|x|/x.

The RKHS norm for the space corresponding to this kernel is given by

with F[Kapprox]\mathcal{F}[K^{\text{approx}}] obtained using the Hankel transforms as above.

Appendix D RKHS Norm of the Neural Spline Kernel

We now discuss how c(a,b)c(a,b) in (12) is related to the Laplacian of the fucntion. If we make the mild assumption that our functions contain a linear and bias term (Lemma 8), then c(a,b)c(a,b) is the Radon Transform of the laplacian of the function. Thus, the least norm minimizers of the least squares problem (1) are related to the laplacian of the function and the RKHS norm corresponds to the integral of the laplacian over hyperplanes in the domain. In our experiments, we added an option to include the linear and bias terms to the solution. They appear to have no effect on the final reconstruction. The derivation below is borrowed from .

Then, flim(x)f_{\text{lim}}(x) can always be rewritten as

We can split the integral in flimf_{\text{lim}} into even and odd parts:

where c+c^{+} and c−c^{-} are the even and odd parts of cc respectively. Observing that [t]++[−t]+=∣t∣[t]_{+}+[-t]_{+}=|t| and [t]+−[−t]+=t[t]_{+}-[-t]_{+}=t, we have that

Using Lemma 8, we will consider without loss of generality, neural networks of the form (25) with even measures c(a,b)c(a,b). We now give a few useful definitions and lemmas.

where ds(x)ds(x) is a measure on the (d−1)(d-1)-hyperplane aTx=ba^{T}x=b. Intuitively the Radon transform represents a function in terms of its integrals along all possible hyperplanes.

Since the hyperplane aTx=ba^{T}x=b is the same as the hyperplane −aTx=−b-a^{T}x=-b, the Radon transform is an even function. i.e. R{f}(a,b)=R{f}(−a,−b)\mathcal{R}\{f\}(a,b)=\mathcal{R}\{f\}(-a,-b).

The Radon Transform satisfies the intertwining property. i.e. for any positive integer ss

where R{f}(a,b)\mathcal{R}\{f\}(a,b) is the Radon Transform of ff. In particular, for d=3d=3,

The Laplacian of flimf_{\text{lim}} in is (25)

which is precisely the Dual Radon Transform of c(a,b)c(a,b). Since cc is even, and assuming it decays rapidly with bb, we can invert it using Lemma 12 yielding

The RKHS norm of the function flimf_{\text{lim}} is