In [1]:
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline

In [2]:
def gaussian2D(x, y, mean_xy, cov):
    mean_xy = np.array(mean_xy)[:,np.newaxis]
    xy = np.vstack((x, y)) - mean_xy
    Cinv = np.linalg.inv(cov)
    
    arg = np.einsum('ni,ij,jn->n', xy.T, Cinv, xy)
    return np.exp(-0.5 * arg)

In [3]:
x = np.arange(0., 10., 0.1)
y = np.arange(0., 10., 0.1)
x,y = np.meshgrid(x,y)
mean_xy = np.array([5.,4.])
cov = np.array([[1.,1.],
                [1.,2.]])**2.
  
f = gaussian2D(x.ravel(), y.ravel(), mean_xy, cov)

plt.imshow(f.reshape(100,100))


Out[3]:
<matplotlib.image.AxesImage at 0x10e2c2a58>