gradient function that treats nan values like edges

Question

I am looking for a function like numpy.gradient, which has a built-in treatment for NaN values. It should recognize where a NaN value occurs and compute one-sided derivatives at those points.

For example if f[x-epx], f[x], and f[x+epx] are all not NaN, I would like it to calculate a derivative like:

grad_f[x] = (f[x+epx] - f[x-epx])/(2*eps)

however, if f[x-eps] is Nan then I would like to obtain the following derivative:

grad_f[x] = (f[x+eps] - f[x])/eps

I am looking for a function that does this also in higher dimensions. Do you know any such function? Thanks a lot in advance.