The function pldg(alpha) plots the gamma distribtion with shape parameter alpha; I have set the scale parameter equal to the shape parameter so the mean is always 1 and the graphs remain comparable as alpha changes.
> pldg function(al = 2) { xx <- seq(0, 5, 0.05) plot(xx, dgamma(xx * al, al) * al, type = "l") }
The function pldgm() applies pldg() to an array of alpha values.
> pldgm function(al = 2) { sapply(al, pldg) NULL }
Hence calling pldgm with an increasing sequence of values like
pldgm(seq(.1, 30, .2))
gives a "movie" showing the convergence from a J-shaped distribution to normality as the shape parameter increases.