There are some problems about the inverse Gaussian process. When I estimate the parameters of an oil leakage process as an inverse Gaussian process, that is, the leakage increment obeys the inverse Gaussian distribution,
with pm.Model() as model:
w = pm.Uniform('w', lower=0, upper=100)
k = pm.Uniform('k', lower=0, upper=100)
q = pm.Uniform('q', lower=0, upper=20)
lamb = pm.Uniform("lamb", lower=0, upper=1000)
deltaLam = pow(time_data[2], q) - pow(time_data[1], q)
mu = pm.TruncatedNormal("mu", mu=w, tau=k ** 2,lower=0)
mu_delta = mu * deltaLam
la = mu ** 3 * deltaLam ** 2 / lamb
Y_obs = pm.Wald(
"Y_obs", mu=mu_delta, lam=la, observed=data)
start = pm.find_MAP(fmin=optimize.fmin_powell)
step = pm.Metropolis()
trace = pm.sample(2000, step=step, chains=3, tune=1000, cores=1)
pm.traceplot(trace)
plt.show()
for RV in model.basic_RVs:
print(RV.name, RV.logp(model.test_point))
print(model.logp(model.test_point))