Getting a likelihood based on the root of an equation

Modeling
7

Sorry, framing it as integrating out a parameter was confusing.

You said you can solve the problem with ABC. How does that go? I assume it’s

  1. Draw A_1,A_2,A_3 from the prior
  2. Solve for t
  3. Accept the draw if \left|t-t^{\star}\right|<\epsilon, otherwise try again

Stan can sample from the same distribution like so:

  1. The prior is the prior for parameters A_1,A_2
  2. The likelihood is the probability of generating an accepted draw, conditional on A_1,A_2.

At step (2), the draw is accepted if A_3 is such that solving the equation yields a root in a narrow interval around the observed t^\star. This in turn means A_3 falls in a narrow interval around A_3^\star (the value of A_3 that implies t=t^\star). The equation is

1=A_{3}tf\left(t,A_{1},A_{2}\right)

(where the complicated parts are hidden inside the function f) and while solving for t requires iterative methods, solving for A_3 when t is known is straightforward

A_{3}^\star=\frac{1}{t^\star f\left(t^\star,A_{1},A_{2}\right)}

(I don’t know why I brought up implicit differentiation previously; you can get \frac{\partial A_3}{\partial t} from the above directly.)
Since we’re dealing with a super narrow interval we can approximate

\left|t-t^{\star}\right|<\epsilon\Leftrightarrow\left|A_{3}-A_{3}^{\star}\right|<\left|\frac{\partial A_{3}}{\partial t}\right|\epsilon

and the probability of hitting that interval

\mathbb{P}\left(\left|A_{3}-A_{3}^{\star}\right|<\left|\frac{\partial A_{3}}{\partial t}\right|\epsilon\right)\approx 2\epsilon \left|\frac{\partial A_{3}}{\partial t}\right| p\left(A_{3}^{\star}\right)

where p\left( A_3\right) is the probability density function for the (prior) distribution of A_3 and \epsilon is a constant (so it doesn’t matter for the posterior distribution).

The above really just amounts to computing the distribution of t — you don’t need a analytic solution for the equation, the density function is the density function of A_3 multiplied by the Jacobian \left|\frac{\partial A_{3}}{\partial t}\right|.

I think I did the derivative wrong; terms always seem to go missing when I apply the chain rule more than once…

Here’s a second attempt (using the identity 1-inv_logit(x)==inv_logit(-x))

for (i in 1:N) {
  real rt_ = rt[i];
  real A1rt2 = A_1[i]*rt_*rt_;
  real invlgA1 = inv_logit(A1rt2);
  real d = inv_logit(-A_2*rt_*invlgA1)*invlgA1*(1+2*A1rt2*inv_logit(-A1rt2));
  target += log(A_3[i]) + log(1/rt_ + A_2*d);
}