3

I have a 3D polygon plot and want to smooth the plot on the y axis (i.e. I want it to look like 'slices of a surface plot').

Consider this MWE (taken from here):

from mpl_toolkits.mplot3d import Axes3D
from matplotlib.collections import PolyCollection
import matplotlib.pyplot as plt
from matplotlib import colors as mcolors
import numpy as np
from scipy.stats import norm

fig = plt.figure()
ax = fig.gca(projection='3d')

xs = np.arange(-10, 10, 2)
verts = []
zs = [0.0, 1.0, 2.0, 3.0]

for z in zs:
    ys = np.random.rand(len(xs))
    ys[0], ys[-1] = 0, 0
    verts.append(list(zip(xs, ys)))

poly = PolyCollection(verts, facecolors=[mcolors.to_rgba('r', alpha=0.6),
                                         mcolors.to_rgba('g', alpha=0.6), 
                                         mcolors.to_rgba('b', alpha=0.6), 
                                         mcolors.to_rgba('y', alpha=0.6)])
poly.set_alpha(0.7)
ax.add_collection3d(poly, zs=zs, zdir='y')
ax.set_xlabel('X')
ax.set_xlim3d(-10, 10)
ax.set_ylabel('Y')
ax.set_ylim3d(-1, 4)
ax.set_zlabel('Z')
ax.set_zlim3d(0, 1)
plt.show()

Now, I want to replace the four plots with normal distributions (to ideally form continuous lines).

I have created the distributions here:

def get_xs(lwr_bound = -4, upr_bound = 4, n = 80):
    """ generates the x space betwee lwr_bound and upr_bound so that it has n intermediary steps """
    xs = np.arange(lwr_bound, upr_bound, (upr_bound - lwr_bound) / n) # x space -- number of points on l/r dimension
    return(xs)

xs = get_xs()

dists = [1, 2, 3, 4]

def get_distribution_params(list_):
    """ generates the distribution parameters (mu and sigma) for len(list_) distributions"""
    mus = []
    sigmas = []
    for i in range(len(dists)):
        mus.append(round((i + 1) + 0.1 * np.random.randint(0,10), 3))
        sigmas.append(round((i + 1) * .01 * np.random.randint(0,10), 3))
    return mus, sigmas

mus, sigmas = get_distribution_params(dists)

def get_distributions(list_, xs, mus, sigmas):
    """ generates len(list_) normal distributions, with different mu and sigma values """
    distributions = [] # distributions

    for i in range(len(list_)):
        x_ = xs
        z_ = norm.pdf(xs, loc = mus[i], scale = sigmas[0])
        distributions.append(list(zip(x_, z_)))
        #print(x_[60], z_[60])

    return distributions

distributions = get_distributions(list_ = dists, xs = xs, mus = mus, sigmas = sigmas)

But adding them to the code (with poly = PolyCollection(distributions, ...) and ax.add_collection3d(poly, zs=distributions, zdir='z') throws a ValueError (ValueError: input operand has more dimensions than allowed by the axis remapping) I cannot resolve.

Ivo
  • 3,890
  • 5
  • 22
  • 53
  • Please add your complete non-working example with imports, norm is missing and topics = dists? – Joe Jan 09 '20 at 19:22

2 Answers2

2

The error is caused by passing distributions to zs where zs expects that when verts in PolyCollection has shape MxNx2 the object passed to zs has shape M. So when it reaches this check

cpdef ndarray broadcast_to(ndarray array, shape):
    # ...
    if array.ndim < len(shape):
        raise ValueError(
            'input operand has more dimensions than allowed by the axis '
            'remapping')
    # ...

in the underlying numpy code, it fails. I believe this occurs because the number of dimensions expected (array.ndim) is less than the number of dimensions of zs (len(shape)). It is expecting an array of shape (4,) but receives an array of shape (4, 80, 2).

This error can be resolved by using an array of the correct shape - e.g. zs from the original example or dists from your code. Using zs=dists and adjusting the axis limits to [0,5] for x, y, and z gives

enter image description here

This looks a bit odd for two reasons:

  1. There is a typo in z_ = norm.pdf(xs, loc = mus[i], scale = sigmas[0]) which gives all the distributions the same sigma, it should be z_ = norm.pdf(xs, loc = mus[i], scale = sigmas[i])
  2. The viewing geometry: the distributions have the positive xz plane as their base, this is also the plane we are looking through.

Changing the viewing geometry via ax.view_init will yield a clearer plot:

enter image description here


Edit

Here is the complete code which generates the plot shown,

from mpl_toolkits.mplot3d import Axes3D
from matplotlib.collections import PolyCollection
import matplotlib.pyplot as plt
import matplotlib.colors as mcolors
import numpy as np
from scipy.stats import norm

np.random.seed(8)
def get_xs(lwr_bound = -4, upr_bound = 4, n = 80):
    return np.arange(lwr_bound, upr_bound, (upr_bound - lwr_bound) / n)

def get_distribution_params(list_):
    mus = [round((i+1) + 0.1 * np.random.randint(0,10), 3) for i in range(len(dists))]
    sigmas = [round((i+1) * .01 * np.random.randint(0,10), 3) for i in range(len(dists))]
    return mus, sigmas

def get_distributions(list_, xs, mus, sigmas):
    return [list(zip(xs, norm.pdf(xs, loc=mus[i], scale=sigmas[i] if sigmas[i] != 0.0 
            else 0.1))) for i in range(len(list_))]

dists = [1, 2, 3, 4]
xs = get_xs()
mus, sigmas = get_distribution_params(dists)
distributions = get_distributions(dists, xs, mus, sigmas)

fc = [mcolors.to_rgba('r', alpha=0.6), mcolors.to_rgba('g', alpha=0.6), 
      mcolors.to_rgba('b', alpha=0.6), mcolors.to_rgba('y', alpha=0.6)]

poly = PolyCollection(distributions, fc=fc)
fig = plt.figure()
ax = fig.gca(projection='3d')
ax.add_collection3d(poly, zs=np.array(dists).astype(float), zdir='z')
ax.view_init(azim=115)
ax.set_zlim([0, 5])
ax.set_ylim([0, 5])
ax.set_xlim([0, 5])

I based it off the code you provide in the question, but made some modifications for brevity and to be more consistent with the usual styling.


Note   -   The example code you have given will fail depending on the np.random.seed(), in order to ensure it works I have added a check in the call to norm.pdf which ensures the scale is non-zero: scale = sigma[i] if sigma[i] != 0.0 else 0.1.

William Miller
  • 9,839
  • 3
  • 25
  • 46
  • Hi William, Thanks a lot for your response! I couldn't quite reproduce your code (i.e. understand which snippets to include where) without further errors making me think I got it wrong. Could you clarify how the code you're suggesting looks? :) Many thanks! – Ivo Jan 21 '20 at 10:10
  • 1
    @Ivo I would be happy to, but it may take me a couple days to get around to it.... feel free to ping me if I haven’t updated my answer in a few days – William Miller Jan 21 '20 at 10:11
  • have you found the time yet to look at this? :) – Ivo Jan 30 '20 at 16:35
  • @Ivo I haven't quite yet - I got kinda busy, but I plan to soon. Thanks for pinging me – William Miller Jan 30 '20 at 18:10
  • @Ivo I've added the complete code I use to generate the plot shown - hopefully that will help clarify things – William Miller Feb 02 '20 at 00:25
1

Using ax.add_collection3d(poly, zs=dists, zdir='z') instead of ax.add_collection3d(poly, zs=distributions, zdir='z') should fix the issue.


Additionally, you might want to replace

def get_xs(lwr_bound = -4, upr_bound = 4, n = 80):
    """ generates the x space betwee lwr_bound and upr_bound so that it has n intermediary steps """
    xs = np.arange(lwr_bound, upr_bound, (upr_bound - lwr_bound) / n) # x space -- number of points on l/r dimension
    return(xs)

xs = get_xs()

by

xs = np.linspace(-4, 4, 80)

Also, I believe the scale = sigmas[0] should actually be scale = sigmas[i] in the line

z_ = norm.pdf(xs, loc = mus[i], scale = sigmas[0])

Finally, I believe you should adjust the xlim, ylim and zlim appropriatly, as you swapped the y and z dimensions of the plot and changed its scales when comparing to the reference code.

William Miller
  • 9,839
  • 3
  • 25
  • 46
Maxpxt
  • 189
  • 4