Implement pre-built BG/NBD/BetaGeoFitter model

[larryshamalama]

It would be good to add a BetaGeoFitter function that returns a ContNonContract with some default priors. A signature that resembles what is provided in the lifetimes package would be a good idea. Something along the lines of the following snipet.

def BetaGeoFitter(name, a, b, r, alpha, T, T0, *, observed, **kwargs):
    p = pm.Beta(f"{name}_beta", a, b, size=size, shape=shape)
    lam = pm.Gamma(f"{name}_gamma", r, 1/alpha, size=size, shape=shape)
    return ContNonContract(name, lam, p, T, T0, size=size, shape=shape, **kwargs)

We should also add some useful summary stats / plots. If they are not specific to the BG/NBD the better!