Plotting polynomials roots NOTE: See follow up to this question here

I created a simple python script to plot quadratic, cubic and quartic polynomials with integer coefficients between -4 and 4. It uses numpy to find the roots for the polynomials and matplotlib for the actual plotting of the points. Because I want to plot all possible polynomial with integer coefficients between -4 and 4 I figured I needed to do a recursive function with a loop. However if I'm not mistaken this will give a time complexity of \$O(n!)\$.

import matplotlib

# Use Qt4Agg because otherwise matplotlib crashes on my computer running Arch Linux
matplotlib.use('Qt4Agg')
import matplotlib.pyplot as plt
import matplotlib.cm as cm
import numpy as np

min_degree = 2
max_degree = 3
colours = []

def polynomial_scatter(ax, min_x, max_x, min_degree, max_degree, cur_degree, coeff):
    global colours

    if cur_degree &lt;= max_degree:
        for c in range(min_x, max_x):
            new_coeff = coeff + [c]
            polynomial_scatter(ax, min_x, max_x, min_degree, max_degree, cur_degree + 1, new_coeff)
            if cur_degree &gt;= min_degree:
                roots = np.roots(new_coeff)

                # check if no possible solution exists
                if len(roots) == 0:
                    continue

                # the real part of the roots will be our x values, the imaginary parts will be our y values
                x, y = zip(*[(t.real, t.imag) for t in roots])

                # scatter the roots for this polynomial with a colour corresponding to the degree of the polynomial
                ax.scatter(x, y, c=colours[cur_degree - max_degree], alpha=.5, s=1)

if __name__ == "__main__":
    f, ax = plt.subplots()
    colours = cm.rainbow(np.linspace(0, 1, max_degree - min_degree + 1))
    polynomial_scatter(ax, -4, 4, min_degree, max_degree, 0, [])
    plt.savefig("plot.svg", format="svg")
    plt.show()


Here are some pictures showing the result:




if I set max_degree to 5 it takes between 7-15 hours on my computer which is very long. Perhaps there is a way to implement multithreading.

How can I improve the performance of the polynomial_scatter algorithm?
I know that you subsequently posted c++ code for this and that post has answers, one of which is accepted. 

I did notice that the code above references the global variable colours. This may be seen as bad programming style because it can cause side effects that may be challenging to detect. Alternatively, function parameters could be used, or else a  class with a class/data/instance member.