Calculating Kernel Matrix Using Numpy Methods

I have a data of shape d X N (each column is a vector of features) I have this code for calculating the kernel matrix: def kernel(x1, x2): return x1.T @ x2 data = np.array([[1,2

Solution 1:

When tested for bigger dimensions of your problem, and random numbers to ensure the robustness, for instance with dimensions (100,200), there are several ways:

import numpy as np

defkernel(x1, x2):
    return x1.T @ x2

defkernel_kenny(a):
    result = []
    for i inrange(a.shape[1]):
      current_result = []
      for j inrange(a.shape[1]):
        x1 = a[:, i]
        x2 = a[:, j]
        current_result.append(kernel(x1, x2))
    
      result.append(current_result)

    return np.array(result)

a = np.random.random((100,200))

res1 = kernel_kenny(a)

# perhaps einsum signature might help you to understand the calculations
res2 = np.einsum('ji,jk->ik', a, a, optimize=True)
# or the following if you want to explicitly specify the transpose# res2 = np.einsum('ij,jk->ik', a.T, a, optimize=True)# or simply ...
res3 = a.T @ a

Hera are the sanity checks:

np.allclose(res1,res2)
>>>True

np.allclose(res1,res3)
>>>True

and timings:

%timeit kernel_kenny(a)
>>>83.2 ms ± 425 µs per loop (mean ± std. dev. of7 runs, 10 loops each)

%timeit np.einsum('ji,jk->ik', a, a, optimize=True)
>>>325 µs ± 4.15 µs per loop (mean ± std. dev. of7 runs, 1000 loops each)

%timeit a.T @ a
>>>82 µs ± 9.39 µs per loop (mean ± std. dev. of7 runs, 10000 loops each)