Matplotlib: Plot path integral of curve under bivariate gaussian

Question

I am trying to plot the path (or line) integral under a (negative log) bivariate gaussian. I am aiming to produce something that looks like the different stages of this neat little gif I found a wikipedia.

So far, I've produced a (negative log) bivariate normal in 3-d:

import matplotlib
import random
import numpy as np
import seaborn as sns
from mpl_toolkits.mplot3d import Axes3D
from matplotlib.mlab import bivariate_normal

sns.set(style="white", palette="muted", color_codes=True)

%matplotlib inline

def pathIntegral():

    # Bivariate gaussian 

    a = np.linspace(-15, 15, 100)
    b = a
    X,Y = np.meshgrid(a, b)
    Z = bivariate_normal(X, Y)
    surprise_Z = -np.log(Z)

    # Figure

    fig = plt.figure(figsize=(10, 7))    

    ax = fig.gca(projection='3d')    
    surf = ax.plot_surface(X, Y, surprise_Z, rstride=1, cstride=1, 
        cmap = plt.cm.gist_heat_r, alpha=0.1, linewidth=0.1)

    ax.scatter(-5, -5, 0, c='k', marker='o')
    ax.scatter(5, 5, 0, c='k', marker='o')

Which produces this:

---

Questions:

1) How do i produce a sufficiently smooth path from two points? E.g. from the two black points that I've drawn in with scatter. I have found this post, but for the life of me can't seem to project that curve on my own x,y axis!

2) Having done so, how would I go about drawing a "curtain" over that path up to the (negative log) bivariate gaussian?

3) An extra plot 2-d (like one of the last ones in the .gif) of the scalar value (negative log gaussian) and euclidiance distance (stright line on x,y) is a bonus question.

Thanks in advance!

Answer 1

import matplotlib.pyplot as plt
import numpy as np
import seaborn as sns
from mpl_toolkits.mplot3d import Axes3D
from matplotlib.mlab import bivariate_normal
import scipy.interpolate


def func(X, Y):
    # Bivariate gaussian
    return -np.log(bivariate_normal(X, Y))

sns.set(style="white", palette="muted", color_codes=True)


def create_smoth_path(points):
    data = np.array(points).T

    tck, u = scipy.interpolate.splprep(data, s=0)
    unew = np.arange(0, 1.01, 0.01)
    out = scipy.interpolate.splev(unew, tck)
    return out


def pathIntegral(path, func_):
    x = path[0]
    y = path[1]
    z = func_(x, y)
    z1 = np.zeros((2, x.shape[0]))
    z1[1] = z

    x1 = np.tile(x, (2, 1))
    y1 = np.tile(y, (2, 1))

    start = min(x.min(), y.min())-2
    stop = max(x.max(), y.max())+2
    a = np.arange(start, stop, 0.5)
    b = a
    X, Y = np.meshgrid(a, b)
    Z = func_(X, Y)

    ax = plt.subplot2grid((3,1), (0,0), rowspan=2, projection='3d')
    ax.plot_surface(X, Y, Z, rstride=1, cstride=1, cmap=plt.cm.gist_heat_r, alpha=0.2, linewidth=0.1)

    ax.plot_surface(x1, y1, z1, rstride=1, cstride=1, color='g')
    ax.scatter(x[0], y[0], 0, c='k', marker='o')
    ax.scatter(x[-1], y[-1], 0, c='k', marker='o')
    ax.plot(x1, y1, color='b')
    dx = np.diff(x)
    dy = np.diff(y)
    dr = np.sqrt(dx**2 + dy**2)
    r = np.concatenate((np.zeros(1), np.cumsum(dr)))
    ax = plt.subplot2grid((3,1), (2,0))
    ax.set_xlim([np.min(r), np.max(r)])
    plt.plot(r, z, 'r')
    plt.xlabel('r')
    plt.ylabel('z(x, y)')
    plt.show()

# path
points = [(-5, -6), (-1, -3), (1, -2), (2, 1), (0, 0), (2, 2), (4, 3), (4, 3.5), (5, 5)]

r = create_smoth_path(points)
pathIntegral(r, func)

Output: