Solución:
De acuerdo con la documentación, desea usar el parámetro de escala (theta), pero como está definiendo beta, que es el inverso de theta, entonces pasa la escala con el valor de 1 / beta, que en su ejemplo sería 1 / 3 o 0,33333.
Por lo tanto, intente:
y1 = stats.gamma.pdf(x, a=29, scale=0.33333)
Como respondió @Hielke, por lo que se explica en la documentación de scipy.stats 1.4.1, parece que el parámetro escalar es igual a beta. De hecho, la función desarrollada originalmente es:
gamma.pdf(x, a) = x^(a-1) * exp(-x) / gamma(a)
Si uno reemplaza X por una combinación de los dos parámetros opcionales loc y escala como :
x = (y - loc) / scale
Uno debería tener:
gamma.pdf(x, a) = (y - loc)^(a-1) * exp( -(y - loc)/scale ) / (scale^(a-1) * gamma(a))
Si lo tomas loc = 0 entonces reconoció la expresión de la distribución Gamma como se define normalmente. Multiplica por el inverso de escala y puedes concluir que escala = beta en esta función y loc es un desplazamiento.
En realidad he intentado detallar la explicación de la documentación:
Específicamente, gamma.pdf (x, a, loc, scale) es idénticamente equivalente a gamma.pdf (y, a) / scale con y = (x – loc) / scale.