Support vector based approach to density estimation

Question

I'm interested in the possiblity of using support vector machines for density estimation. The main source I found for this is a paper by Mukherjee & Vapnik 'Support vector method for multivariate density estimation' (Link) which claims that this approach performs better than for example the Parzen-Window approach. The scikit-learn documentation mentions the possibility of using SVMs for this purpose, however there's no clear description of it to be found.

The paper by Mukherjee & Vapnik presents a different optimization problem that is to be solved using the support vector approach for this, while the scikit-learn documentation implies the use of one-class SVMs for this.

Do these two approaches amount to different results? Would an estimation as presented in the Mukherjee & Vapnik paper amount to the same as a one-class-SVM density estimation as shown in the scikit-learn example, if done for both classes in a binary problem?

I have failed to find a reference implementation of the approach by Mukherjee & Vapnik. Is anyone aware of an implementation that I could have a look at?

Can anyone provide further material on this matter?

Answer 1

I think the two approaches are distinct. The paper by Mukherjee & Vapnik describes a method that is akin to a SVR of the empirical distribution function. They also have adaptive kernel width. Conversely one-class SVMs tries to build "a function that takes the value +1 in a small region capturing most of the data points, and -1 elsewhere" (see Schölkopf et al.).

I am not aware of implementations of the former method and I have never seen it applied. Regarding one-class SVMs, they are widely used.

If you want to use them for density estimation, my guess would be that you should be able to use the values returned by decision_function of a one-class SVM to compute the pdf, after adding back the bias $\rho$ (which is effectively used in a one-class SVM as threshold for outlier detection) and scaling the result.

Looking at sklearn/svm/src/libsvm/svm.cpp, we can see that (as expected) decision_function returns:

$$ h(\mathbf{x}) = \left[ \sum_{s \in \mathcal S} \alpha_s K(\mathbf{x}^{(s)},\mathbf{x}) \right] - \rho $$ where $\mathcal S$ is the set of support vectors and $K()$ is $$ K(\mathbf{x},\mathbf{z}) = \exp\left(-\gamma \, || \mathbf{x} - \mathbf{z} ||_2^2 \right). $$

We need to add back $\rho$ so that what we obtain is an unnormalised mixture of $d$-dimensional multivariate Gaussians with covariance matrix $(2\gamma)^{-1} \mathbf{I}$.

Additionally, we need to normalise the result so that it integrates to one. Each Gaussian in the summation above misses a normalisation factor: $$ \left(\frac{\gamma}{\pi}\right)^{d/2} $$ which needs to be added back. Also, the coefficients $\mathbf{\alpha}$ must add to one.

Overall, I think the density estimation can be computed as: $$ p(\mathbf{x}) = \left(\frac{\gamma}{\pi}\right)^{d/2} \, |\mathbf{\alpha}|_1^{-1} \, (h(\mathbf{x}) + \rho) $$

Something like (not tested):

n_dim = 2
gamma = ...
clf = svm.OneClassSVM(nu=..., kernel="rbf", gamma=gamma)
clf.fit(X_train)
xx, yy = np.meshgrid(np.linspace(...), np.linspace(...))
Z = clf.decision_function(np.c_[xx.ravel(), yy.ravel()])
Z = (Z - clf.intercept_[0]) * (gamma/np.pi)**(n_dim/2.) / np.sum(clf.dual_coef_)
Z.reshape(xx.shape)
plt.contourf(xx, yy, Z)

(Note: it looks like scikit-learn stores $-\rho$, rather than $\rho$, in clf.intercept_[0].)