# Non-stationary Random Fourier Features

Random Fourier Features (RFFs) were introduced by A. Rahimi and B. Recht in their 2007 publication Random Features for Large-Scale Kernel Machines. RFFs are based on Bochner's theorem, which applies to stationary kernels on Euclidean spaces i.e. any positive semi-definite function $$k$$ satisfying $$k(x, y) = \kappa(||x - y||), ~\forall~x, y~ \in \mathbb{R}^d.$$ Bochner's theorem states that such a function can be written as a Fourier transform : $$k(x, y) = \int_{\mathbb{R}^d} \exp{\left[i \omega^T (x - y)\right]}\mathbb{P}(d\omega),$$ where $$\mathbb{P}$$ is a probability measure. This is an expectation with respect to $$\mathbb{P}$$ and we can therefore approximate $$k$$ via Monte-Carlo based on a sample drawn from $$\mathbb{P}$$.

It turns out that Bochner's theorem can be extended to non-stationary kernels on Euclidean spaces as well. This is known as Yaglom's theorem or the Yaglom-Genton theorem. A non-stationary kernel can therefore also be approximated via a Monte Carlo scheme. This is what is done in a 2017 paper by J-F. Ton et al., Spatial Mapping with Gaussian Processes and Nonstationary Fourier Features. I do not understand the first steps derivation of the Monte Carlo scheme in that paper.

For some context, Yaglom's theorem says that a non-stationary kernel $$k$$ on $$\mathbb{R}^d$$ can be written as $$k(x, y) = \int_{\mathbb{R}^d \times \mathbb{R}^d} \exp{\left[i \left( \omega_1 ^T x - \omega_2 ^T y \right) \right]} \mu(d\omega_1, d\omega_2),$$ where $$\mu$$ is the Lebesgue-Stieltjes measure associated to some $$f$$ which is a positive semi-definite function with bounded variation. This is the same as $$k(x, y) = \int_{\mathbb{R}^d \times \mathbb{R}^d} \exp{\left[i \left( \omega_1 ^T x - \omega_2 ^T y \right) \right]} f(\omega_1, \omega_2)d\omega_1 d\omega_2,$$ where $$f$$ is as above. The authors first claim that $$f$$ can be written as $$f(\omega_1, \omega_2) = g(\omega_1, \omega_2) + g(\omega_2, \omega_1) + g(\omega_1, \omega_1) + g(\omega_2, \omega_2),$$ where $$g$$ is some density. It is not clear to me why there should exist a density $$g$$ satisfying that.

Further, the authors then plug in that form for $$f$$ in the integral. That should lead to $$k(x, y) = \frac{1}{4} \int_{\mathbb{R}^d \times \mathbb{R}^d} \exp{\left[i \left( \omega_1 ^T x - \omega_2 ^T y \right) \right]} \left( g(\omega_1, \omega_2) + g(\omega_2, \omega_1) + g(\omega_1, \omega_1) + g(\omega_2, \omega_2) \right) d\omega_1 d\omega_2.$$ The authors obtain instead $$k(x, y) = \frac{1}{4} \int_{\mathbb{R}^d \times \mathbb{R}^d} \left(\exp{\left[i \left( \omega_1 ^T x - \omega_2 ^T y \right) \right]} + \exp{\left[i \left( \omega_2 ^T x - \omega_1 ^T y \right) \right]} + \exp{\left[i \omega_1 ^T \left(x - y \right) \right]} + \exp{\left[i \omega_2 ^T \left(x - y \right) \right]} \right) \mu(d\omega_1, d\omega_2).$$ I do not understand how they get to that expression. Could someone unblock me ?

The authors did not mean that all $$f$$'s may necessarily be written as this sum of four terms in $$g$$. They chose to place themselves in the settings where $$f$$ may be written this way (in so-called "symmetrised" form) because in that case the imaginary part disappears automatically when deriving the RFF scheme.
In the stationary case, we do not need to "symmetrise" : in fact what is typically done when deriving the RFF scheme is to use the fact that as the kernel is real-valued, $$k(x, y) = \text{Re}\left[k(x, y)\right]$$ (i.e. its real part). This is what was done in the original 2007 paper on RFFs by A. Rahimi and B. Recht, by the way.
In the nonstationary case, this "trick" can be used but you end up with a form of $$k$$ from which you cannot write the kernel Gram matrix in factorised form...which is a problem since what makes RFFs work is that the kernel Gram matrix can be written in factorised form $$\text{K}_{\text{X}} = \Psi \Psi^{\text{T}}$$.
On the second part of my question, the problem is mostly a matter of notations : what is meant by $$g(\omega_1, \omega_1)$$ is $$g(\omega_1, \omega_2) \mathbb{1}\{\omega_1 = \omega_2\}$$. Likewise for $$g(\omega_2, \omega_2)$$.
Further : in a 2015 paper by Samo and Roberts, section 2.3 specifically, it is shown how to derive a RFF scheme for some nonstationary kernels without having to assume $$f$$ can be written as the sum of the four terms in $$g$$. What the results in this Samo & Roberts paper imply is that not all nonstationary kernels lend themselves to RFFs. This is a clear difference with the stationary case : all stationary kernels may be approximated by RFFs.