Numerical solution for inverse of a function

[Torsten Ponik]

Dear all,

I need to find the inverse of a slightly modified Gompertz function. To my understanding, there is no analytical solution to this problem. Is there a way in Stata to numerically approximate the inverse of such a function?

Function:

y = n + b1*exp(-exp(-(b2*-m)(x-b3)))

(Note: I only estimate b1, b2, b3 and set n and m manually)
(Note 2: it only has to work for x > 0 and y <=1)

[Joseph Coveney]

Maybe something along these lines?

version 18.0

clear *

mata:
mata set matastrict on

struct Stuff {
    real scalar m, n, x, y
}

void function evaluator(real scalar todo, real rowvector P,
    struct Stuff scalar o, real scalar v, ///
    real rowvector G, real matrix H) {

    pragma unused todo; pragma unused G; pragma unused H

    v = (o.y - o.n + P[1] * exp( -exp( -(P[2] * - o.m) * (o.x - P[3]))) )^2
}

o = Stuff(1)
o.m = 1; o.n = 1; o.x = 1; o.y = 1

s = optimize_init()
optimize_init_which(s, "min")
optimize_init_technique(s, "nm")
optimize_init_nmsimplexdeltas(s, 1e-5)
optimize_init_evaluator(s, &evaluator())
optimize_init_argument(s, 1, o)
optimize_init_params(s, (1, 1, 1))

optimize(s)

end

exit

Above is intended to be just illustrative of a possible approach. And it's prototype code; for production work, you would encapsulate things in a class with a public method, or write a front function that you can call and have it handle the nitty-gritty work behind the scenes.

Also, if you've got better ideas for starting values, then you could probably get better (faster) performance with the default Newton-Raphson algorithm.

[Joseph Coveney]

Got a little too fast there. Better that

v = (o.y - o.n + P[1] * exp( -exp( -(P[2] * - o.m) * (o.x - P[3]))) )^2

be written

v = (o.n + P[1] * exp( -exp( -(P[2] * - o.m) * (o.x - P[3]))) - o.y )^2

[Torsten Ponik]

Thanks a lot Joseph!