A look at the Gamma distribution, possibly for rainfall modeling

Andy Long
Spring, 2018

Issues:

Rainfall data has a lot of zeros, which are troublesome for gamma distributions. On the other hand, zeros aren’t really zeros. That’s a curious notion, but actually rainfall data at zero just means “below the detection limit”. So it seems to me that we could find the detection limit dl, and then replace zeros with random uniform values in the interval (0,dl). Then do a fit, perhaps to the cumulative (cdf, which is better behaved than the pdf) to estimate  α and β, the shape and scale parameters.

Here’s the pdf of the gamma distribution, with its two parameters α (shape) and β (scale):

PDF[GammaDistribution[α,β],x]
DGam[x_]=Together[D[(((E^(-(x/β))) (x^(-1+α)) (β^(-α)))/(Gamma[α])),x]]
(* an example pdf: *)
Plot[PDF[GammaDistribution[3.5,2],x],{x,0,20},Filling->Axis]
(* an example cdf: *)
Plot[CDF[GammaDistribution[3.5,2],x],{x,0,20},Filling->Axis]

and here’s a smattering of different shapes and sizes of gammas:

Plot[Table[PDF[GammaDistribution[α,β],x],{α,.2+{1,4,6}},{β,2.3}]//Evaluate,{x,0,20},Filling->Axis,PlotRange-> All]
Plot[Table[PDF[GammaDistribution[α,β],x],{β,.2+{1,4,6}},{α,2.3}]//Evaluate,{x,0,20},Filling->Axis,PlotRange-> All]

The mean and variance are given this way in terms of α and β:

gammaMean=α*β
gammaVariance=α*β^2
soln=Solve[DGam[x]==0,x]
gammaMax=x/.First[First[soln]]

This animation may explain why its two parameters α and β are called the shape and scale parameters:

Note that I’m adjusting domain and range as one varies alpha and beta:

Manipulate[
Module[
{},
gammaMean=α*β;
gammaVariance=α*β^2;
gammaMax=(α-1)β;
Plot[CDF[GammaDistribution[α,β],x],{x,0,2*gammaMean},Filling->Axis,PlotRange->{{0,2*gammaMean},{0,1}},PlotLabel->"α shape, β scale"]
]
,{{α,2},0.1,10,.01}
,{{β,.5},0.1,10,.01}
]

Manipulate[
Module[
{},
gammaMean=α*β;
gammaVariance=α*β^2;
gammaMax=(α-1)β;
Plot[PDF[GammaDistribution[α,β],x],{x,0,2*gammaMean},Filling->Axis,PlotRange->{{0,2*gammaMean},{0,1.2*PDF[GammaDistribution[α,β],gammaMax]}},PlotLabel->"α shape, β scale"]
]
,{{α,2},0.1,10,.01}
,{{β,.5},0.1,10,.01}
]