Add sample plot to Interpolated Distribution docs

[canyon289]

Description of your problem

Interpolated Docs are missing sample plot. One should be added
https://docs.pymc.io/api/distributions/continuous.html#pymc3.distributions.continuous.Interpolated

Please provide any additional information below.
See example from Normal plot for how how this can be done
https://github.com/pymc-devs/pymc3/blob/master/pymc3/distributions/continuous.py#L433

[Dpananos]

This looks fun, I'll take it on. Looks like this is meant to interpolate a set of points on a grid and constrict a distribution out of those points. The other plots show distributions under various parameterizations, which I don't think this is amenable to.

Would it be sufficient to construct a lattice, sample from it, and then plot a histogram and lattice points? Here is an example of what I'm thinking

[canyon289]

Sounds good to me :)

[Rishabh261998]

Hi, can I take this issue? Thanks

[mjhajharia]

so, im a little confused about something, the input into the interpolated function - lattice points and their pdf, are those supposed to be sampled from a pymc3 distribution(intentionally picking them from uneven intervals), or they can be any set of monotonically increasing random points. also i'm still trying to make sense of what would pdf at a point mean, as in wouldnt probability distribution at that point be a better measure. i might be absolutely wrong, but i guess probability density function is sort of measured around a point right, so maybe distribution would make more sense, its just a lil ambiguous in my head. also is it any different from the scipy interpolation thing, or is it the same thing.

[mjhajharia]

so, im a little confused about something, the input into the interpolated function - lattice points and their pdf, are those supposed to be sampled from a pymc3 distribution(intentionally picking them from uneven intervals), or they can be any set of monotonically increasing random points. also i'm still trying to make sense of what would pdf at a point mean, as in wouldnt probability distribution at that point be a better measure. i might be absolutely wrong, but i guess probability density function is sort of measured around a point right, so maybe distribution would make more sense, its just a lil ambiguous in my head. also is it any different from the scipy interpolation thing, or is it the same thing.

okay got the density thing, that was my bad. still confused about some of it

[canyon289]

@Rishabh261998 Since you were first you're welcome to work on this, although it looks like @almostmeenal is interested as well.

If you're still interested could you comment in the next 3 days indicating your interest

[mjhajharia]

I couldn't figure out a way to plot points which aren't a part of a pymc3 distribution in arviz, so i just did matplotlib to get an idea, on a small sample(100) , it was terrible lol, so I'm still wondering if it makes sense to sample a pymc3 distribution randomly as input for the interpolation or not, because it needs a static input.

so, im a little confused about something, the input into the interpolated function - lattice points and their pdf, are those supposed to be sampled from a pymc3 distribution(intentionally picking them from uneven intervals), or they can be any set of monotonically increasing random points. also i'm still trying to make sense of what would pdf at a point mean, as in wouldnt probability distribution at that point be a better measure. i might be absolutely wrong, but i guess probability density function is sort of measured around a point right, so maybe distribution would make more sense, its just a lil ambiguous in my head. also is it any different from the scipy interpolation thing, or is it the same thing.

okay got the density thing, that was my bad. still confused about some of it

[OriolAbril]

I think this example plot needs to be significantly more complicated than the plots for other distributions.

The plots for all the other distributions show the pdf (probability density function) for several values of the input parameters. The plot for the normal distribution for example shows the pdf of 4 different normal distributions, each with different mean and standard deviation.

In the case of the Interpolated distribution, we don't have the pdf, so there is no way to reproduce these kind of plots, we have to generate a completely different plot customized to the Interpolated case. The main use of the interpolated distribution is to be able to use distributions whose pdf is not known in PyMC3 models. Instead of knowing the pdf like we do for the other cases, here we only have available some evaluations of the pdf on a grid (not necessarily uniform). With this info we can create an Interpolated distribution that will behave like this unknown distribution of which we only know the pdf at some given points.

In order to showcase this conceptual difference, I think we can follow these steps:

1. Define a mixture pdf from two normals, which will act as our unknown pdf
2. Evaluate this pdf on a uneven grid, this is the info that we will use to try and reproduce the original pdf with Interpolated
3. Create the Interpolated distribution
4. Plot something

Note that the Interpolated distribution defines the logp and the random method that should both approximate the original and unknown pdf, so I think it is a good idea to have the following elements in the plot:

1. Original (unkown to Interpolate pdf) plotted on a fine grid (uniform grid too to make things easy) -> lineplot
2. Known (to Interpolate) points of the pdf evaluated on the nonuniform grid (rough grid) -> scatterplot
3. Interpolated distribution plotted on the fine grid -> lineplot
4. Histogram of ~100/500 samples from the interpolated distribution -> histogram

The example image provided above by Dpananos has elements 2 and 4 only. We could also use ArviZ for 4 and plot a KDE instead of a histogram, but I think it is an overkill and its better to keep the example with minimal requirements.

Implementation proposals

# preprocessing: definitions taken from https://github.com/pymc-devs/pymc3/blob/master/pymc3/distributions/continuous.py#L445
import matplotlib.pyplot as plt
import numpy as np
import scipy.stats as st
plt.style.use('seaborn-darkgrid')
x = np.linspace(-5, 5, 1000)

# step 1: define mixture
# use scipy.stats (like in other pdf plots) to define this unkown pdf
# I propose this two normals of different loc and scale
original_pdf = lambda x: st.norm.pdf(x, -2, 1) / 2 + st.norm.pdf(x, 3, .6) / 2

# step 2: evaluation of rough grid
# as the normals have different scale, we "need" more points on >0 because changes are steeper
# so we use a non uniform grid
xgrid = np.concatenate((np.linspace(-5 ,0, 10, endpoint=False), np.linspace(0, 5, 20)))
# evaluate pdf on xgrid

# step 3: create Interpolated
# note this is outside a model, so you'll need to use pm.Interpolated.dist(xgrid, pdf_points), not pm.Interpolated(...)
# evaluate the interp_dist at x array, you can use `.logp(x).eval()` to get the log, then exponentiate before plotting
# generate random samples from interp_dist

# step 4: plot
# add the 4 plots, with a legend, x/ylabel...

[mjhajharia]

In the case of the Interpolated distribution, we don't have the pdf, so there is no way to reproduce these kind of plots, we have to generate a completely different plot customized to the Interpolated case. The main use of the interpolated distribution is to be able to use distributions whose pdf is not known in PyMC3 models. Instead of knowing the pdf like we do for the other cases, here we only have available some evaluations of the pdf on a grid (not necessarily uniform). With this info we can create an Interpolated distribution that will behave like this unknown distribution of which we only know the pdf at some given points.

okay yeah that makes sense, so essentially its line fitting, except the linear limitation is removed, I'll see if I can get it done today, also I guess it's cool if I play around with the input as long as the final figure makes sense? I mean I was trying random stuff and often the interpolation wasn't that good for some randomized points, so yeah.

[OriolAbril]

I guess it's cool if I play around with the input as long as the final figure makes sense?

Play with it as much as you want, I added this toy example because I figured it would help in making the comment and idea on how to move forward clear, not as the final example to be used in the docs

[mjhajharia]

hi, does this work? I'm not sure if i did it right, also should I adhere to the color theme in previous plots?

here's the code:

x = np.linspace(gamma.ppf(0.01, a),gamma.ppf(0.99, a), 10000)
rv = gamma(a)
dictionary = dict(zip(x, rv.pdf(x)))
sampled_dict = {k: v for k, v in random.sample(dictionary.items(),50)}
input_dict = OrderedDict(sorted(sampled_dict.items()))
points = np.array(list(input_dict.keys()))
pdf = np.array(list(input_dict.values()))
interpolated = pm.Interpolated.dist(points, pdf)

plt.style.use('seaborn-darkgrid')
fig, ax = plt.subplots(1, 1)
ax.plot(x, rv.pdf(x), color='deepskyblue', linestyle = '--', lw=4, label='Original Gamma pdf',alpha=0.8)
ax.plot(points, pdf, color='black', marker='o', label='Lattice Points',alpha=0.5,linestyle='')
ax.plot(interpolated.x_points, interpolated.pdf_points,color='forestgreen',lw=4,label='Interpolated pdf',alpha=0.8,)
r = interpolated.random(size=1000)
ax.hist(r, density=True, histtype='stepfilled', alpha=0.4,bins = 25,align ='right',color='grey')
ax.legend(loc='best', frameon=False)
plt.show()

[mjhajharia]

like ofc the plot isnt very nice and I think playing with other distributions would be better, I'll fine tune the plot later, but as a rough idea of what we're trying to do does this make sense

[OriolAbril]

This looks good, I will add some comments here, but I'd encourage you to already submit a PR with these changes as it will be better to review and comment on specific parts of the code. You can add [WIP] or use github to mark it as draft and work at your own pace until it's finished and can be merged.

x = np.linspace(gamma.ppf(0.01, a),gamma.ppf(0.99, a), 10000)
rv = gamma(a)

This is great, using the ppf is clever to get the right limits for the linspace. I only have two very minor comments. First, using 10000 points is probably a bit too much, all the other plots use 1000 and do already look perfectly smooth. Second, I'd define rv before x and use rv.ppf instead of generating a new gamma with parameter a every time.

dictionary = dict(zip(x, rv.pdf(x)))
sampled_dict = {k: v for k, v in random.sample(dictionary.items(),50)}
input_dict = OrderedDict(sorted(sampled_dict.items()))
points = np.array(list(input_dict.keys()))
pdf = np.array(list(input_dict.values()))

I strongly discourage using a dictionary for this. All the code above can be reduced to two steps: defining the rough grid points and evaluating the pdf on these points to create pdf array, like so:

points = np.random.choice(x, size=50, replace=False)
pdf = rv.pdf(points)

Moreover, I'd personally not use a random grid, here is why. It is true that Interpolated accepts nonuniform grids, but it is also doing a linear interpolation between the points of the grid, with random points, these can be extremely close (as it happens a few times by seeing that there are black dots even with the alpha set) and therefore, interpolating between them doesn't really make sense. If this were a real case and we had no other option than to use a random grid, then I'd say go for it, but as this is the official docummentation and there is only one example, I think that the example should illustrate best practices. Otherwise, I fear users will use random grids as defaults because it's what they have seen in the docs, even when a uniform grid could be used and do a much better job.

ax.plot(interpolated.x_points, interpolated.pdf_points,color='forestgreen',lw=4,label='Interpolated pdf',alpha=0.8,)

The values used for plotting should not be x_points and pdf_points, we have already plotted them in the previous line. here we should be evaluating the pdf of the interpolated distribution on all the points in x. As I commented in the sample code above, you can do this with logpdf(x).eval():

ax.plot(x, np.exp(interpolated.logpdf(x).eval()), ...)

As for general aesthetics of the plot, I think it would be great to try to keep the style similar to the other plots of the section. The black dots and grey histogram are already great as is, for the other two lines, I'd use C0 and C1 as colors, and I would not set the lw manually so it falls back to the style default. The dashed linestyle, alpha and so on are already great, no need to change them