Numpy power devuelve un valor negativo.

Quiero trazar la distribución de Poisson y obtener probabilidades negativas para lambda >= 9.

Este código genera gráficos para diferentes valores de lambda:

import numpy as np
import matplotlib.pyplot as plt
from scipy.special import factorial

for lambda_val in range(1, 12, 2):

plt.figure()

k = np.arange(0, 20)
y = np.power(lambda_val, k)*np.exp(-lambda_val)/factorial(k)

plt.bar(k, y)
plt.title('lambda = ' + str(lambda_val))
plt.xlabel('k')
plt.ylabel('probabilidad')
plt.ylim([-0.1, 0.4])
plt.grid()
plt.show()

Por favor, vea estos dos gráficos:

Lambda = 5 se ve bien en mi opinión.
Imagen del gráfico lambda=5

Lambda = 9 no se ve bien.
Imagen del gráfico lambda=9

Estoy bastante seguro de que tiene algo que ver con np.power porque
np.power(11, 9) me da: -1937019605, mientras que
11**9 me da: 2357947691 (lo mismo en WolframAlpha).

Pero si evito np.power y uso
y = (lambda_val**k)*math.exp(-lambda_val)/factorial(k)
para calcular la probabilidad, obtengo valores negativos también. Estoy totalmente confundido. ¿Alguien puede explicarme el efecto o qué estoy haciendo mal? Gracias de antemano. 🙂