Derivative of Log-CDF

Questions
1

Hi all,

I’m writing a model which uses a Gamma distributed RV. I compute the cdf of the RV, and use this to perform later computation. When I want to perform sampling, I get errors which state that the derivative of the logcdf function is not available.

Is there a way I can implement this? The derivative of the CDF is available in closed form (since it is pdf of the Gamma distribution, so I thought this should be easy to implement?

Thanks!

2

If the CDF is implemented using theano, the derivative is available automatically. Also, logcdf is available in Gamma distribution so you should consider using it instead: https://github.com/pymc-devs/pymc3/blob/67155fa4d000e73f99d761970ece73c2236a97d5/pymc3/distributions/continuous.py#L2474

3

Hi,

Sorry for the slow response! In fact, I’m using logcdf.

Here is the code I use:

self.GI_mean = pm.Normal('GI_mean', gi_mean_mean, gi_mean_sd)
self.GI_sd = pm.Normal('GI_sd', gi_sd_mean, gi_sd_sd)

gi_beta = self.GI_mean / self.GI_sd ** 2
gi_alpha = self.GI_mean ** 2 / self.GI_sd ** 2

GI_dist = pm.Gamma.dist(alpha=gi_alpha, beta=gi_beta)

bins = np.zeros(gi_truncation + 1)
bins[1:] = np.arange(gi_truncation)
bins[2:] += 0.5
bins[:2] += 1e-5

cdf_vals = T.exp(GI_dist.logcdf(bins))
pmf = cdf_vals[1:] - cdf_vals[:-1]
GI_rev = T.repeat(T.reshape(pmf[::-1] / T.sum(pmf), (1, 1, gi_truncation)), 2, axis=0)

This throws the following error:

MethodNotDefined: ('grad', <class 'theano.scalar.basic_scipy.GammaInc'>, 'GammaInc')

Which is the reason for this question - this error makes me think that the derivative of the logcdf here is not setup. Maybe I’ve done something else wrong instead?

Thanks!

4

Oh I see, :frowning: Could you raise an issue on theano-pymc https://github.com/pymc-devs/Theano-PyMC?

1 Like
5

Cool, sure - I just did.

2 Likes
6

I ran into the same issue. I needed to use logcdf to implement a model with left censored data. A convenience function for something like gamma_lcdf hasn’t been implemented yet, so I wanted to use logcdf rather than write my own, but encountered this problem.

Edit: Adding a little more detail.

I wanted to implement a censored model with pm.Potential, as is generally recommended. I was looking to implement something convenient like:

left_cens = pm.Potential('left_censored', pm.Gamma.dist(alpha, beta).logcdf(cens))

But encountered the same error of gradient being undefined. I believe gradient is only being used for NUTS sampling. If I set step=pm.Metropolis() like I might have to do if I were dealing with a likelihood that is not differentiable, then I can sample and nothing breaks. I’d just prefer to use NUTS while sampling for obvious reasons.

1 Like
7

I remember I asked the same question few years ago. Would it be possible to remove logcdf for gamma if it is simply not implemented?

8

The logcdf of the gamma function is correctly implemented and there is no reason why it should be removed.

What is not implemented is the gradient of the logcdf as it depends on the gradient of the incomplete gamma function, which has a complicated form.