Comparison between the Gaussian kernel matrix and the matrix computed from the dot products of the vectors of CRandomFourierDotFeatures.

random_fourier_kernel_comparison.cpp

#include <shogun/base/init.h>
#include <shogun/features/RandomFourierDotFeatures.h>
#include <shogun/kernel/GaussianKernel.h>
#include <shogun/kernel/normalizer/IdentityKernelNormalizer.h>

#include <stdio.h>

using namespace shogun;

int main(int argv, char** argc)
{
	init_shogun_with_defaults();
	
	int32_t dims[] = {100, 200, 300, 400, 500, 650};
	for (index_t d=0; d<6; d++)
	{
		int32_t num_dim = dims[d];
		SG_SPRINT("Starting experiment for number of dimensions = %d\n", num_dim);
		int32_t num_vecs = 10;
		SGMatrix<float64_t> mat(num_dim, num_vecs);
		for (index_t i=0; i<num_vecs; i++)
		{
			for (index_t j=0; j<num_dim; j++)
			{
				mat(j,i) = CMath::random(0,1) + 0.5;
			}
		}

		SGVector<float64_t> params(1);
		params[0] = num_dim - 20;
		SG_SPRINT("    Using kernel_width = %f\n", params[0]);

		CDenseFeatures<float64_t>* dense_feats = new CDenseFeatures<float64_t>(mat);
		SG_REF(dense_feats);

		CIdentityKernelNormalizer* normalizer = new CIdentityKernelNormalizer();
		SG_REF(normalizer);

		CGaussianKernel* kernel = new CGaussianKernel(dense_feats, dense_feats, params[0]);
		kernel->set_normalizer(normalizer);
		
		SGMatrix<float64_t> kernel_mat = kernel->get_kernel_matrix();
		
		int D[] = {50, 100, 200, 300, 400, 500, 1000, 2000};
		for (index_t i=0; i<8; i++)
		{
			CRandomFourierDotFeatures<float64_t>* rand_feats =
					new CRandomFourierDotFeatures<float64_t>(dense_feats, D[i], KernelName::GAUSSIAN, params); 
			SGMatrix<float64_t> rand_mat(num_vecs, num_vecs);
			for (index_t j=0; j<num_vecs; j++)
			{
				for (index_t k=0; k<num_vecs; k++)
					rand_mat(j,k) = rand_feats->dot(j, rand_feats, k);
			}
			
			SG_UNREF(rand_feats);
			float64_t max = -1;
			for (index_t j=0; j<num_vecs; j++)
			{
				for (index_t k=0; k<num_vecs; k++)
				{
					rand_mat(j,k) = CMath::abs(kernel_mat(j,k)-rand_mat(j,k));
					if (rand_mat(j,k) >= max)
						max = rand_mat(j,k);
				}
			}

			SG_SPRINT("\tMax diff for D = %d : %f\n", D[i], max);
		}
		
		SG_UNREF(kernel);
		SG_UNREF(dense_feats);
		SG_UNREF(normalizer);
	}
	exit_shogun();
}