AoY Hybrid Model
AoY-Hybrid.RmdBecause of the utility of the AoY as a predictor of parasite densities, we developed hybrid variables describing the change in the mean AoY.
\[ \frac{\partial y}{\partial a} + \frac{\partial y}{\partial \alpha}= - (F(r,m)+h) y, \]
with the boundary condition \[y_\tau(a,0)=h_\tau(a).\]
We write it in this form, in part, because we do not know how to write the function \(F(r,m).\) If we assumed \(F(r,m) = r\), then we would be assuming that every infection that cleared reset the AoY to zero. If some individuals remained infected, however, the AoY would not reset to zero. It would, instead, change to the age of the next youngest.
In the next section, we derive an equation for \(dy/da,\) and then we find a functional form that works pretty well.
Derivation for \(dy/da\)
We start by multiplying by \(\alpha\) and integrating with respect to \(\alpha\):
\[ \int_0^a \alpha \frac{\partial y(\alpha, a)}{\partial a} d \alpha + \int_0^a \alpha \frac{\partial y(\alpha, a)}{\partial \alpha} d\alpha = - \left(F(r,m) + h \right) \int_0^a \alpha y(\alpha, a) d\alpha \]
but
\[\int_0^a \alpha y(\alpha, a) d\alpha = \left<Y \right> p\]
so
\[ \frac{d \left< Y \right> p}{da} + \int_0^a \alpha dy = - p \left(F(r,m) + h \right) \left< Y \right> \]
Using the chain rule:
\[ \frac{\partial \left<A \right> p}{\partial a}= p \frac{\partial \left<A \right>}{\partial a} + \left<A \right> \frac{\partial p}{\partial a} \]
and substituting:
\[ p \frac{d \left< Y \right>}{da} + \left<Y \right> \frac{dp}{da} + \int_0^a \alpha dy = - p \left(F(r,m) + h \right) \left<Y \right> \]
but \(dp/da = h(1-p) - R(m) p\) so
\[ p \frac{d \left< Y \right>}{da} = - \left[h (1-p) - R(m) p \right] \left<Y \right> - p \left(F(r,m) + h \right) \left<Y \right> - \int_0^a \alpha dy \]
and now lots of things cancel, so:
\[ p \frac{d \left< Y \right>}{da} = - h \left<Y \right> -(F(r,m)-R(m)) \left<Y \right> - \int_0^a \alpha dy \]
For the last term, we integrate by parts:
\[ \int_0^a \alpha \partial z = \alpha y(a,\alpha) \big|_{\alpha = a} - \alpha y(a,\alpha) \big|_{\alpha = 0} - \int_0^a y d\alpha \]
which simplifies to: (The first terms are wrong, I think is wrong)
\[ \int_0^a \alpha \partial y = a h_d(0)e^{-ra} \]
\[ \frac{dy}{da} = 1 - \frac{h}{p} y - \left( F(r,m) - R(m) \right) y \]
\(F(r,m)\)
Our intuition was that the changes in the AoY could be worked out by treating the MoI cases. We arrived at a function:
\[F(r, \zeta) = R(m) \left( 1 + \frac{\zeta^{\sqrt{2}}}{4} + \sum_{k \geq 3} \frac{1}{k} \frac{\zeta^3}{k!} \right). \]
\[F(r, \zeta) = r \frac{\mbox{Pr}(\zeta=1; m(t))}{\mbox{Pr}(\zeta>0; m(t))} \left( 1 + \frac{\zeta^{\sqrt{2}}}{4} + \sum_{k \geq 3} \frac{1}{k} \frac{\zeta^3}{k!} \right). \]
Performance did not seem to improve by adding terms for \(k>5.\)
\(\zeta=1\) If MoI=1, then the term should clear. This cancels out the other \(R(m)\) term.
We reasoned that if the MoI were 2 or greater, then a change would occur with probability \(1/\zeta\), increasing the AoY to something close to \(y\).
\(\zeta=2\): Somehow, we made a mistake when we wrote down the first term, but we didn’t notice until much later. After fiddling around with the exponent, we found that \(\zeta^{\sqrt{2}}\) worked better than \(\zeta\) (where we started, by mistake), and curiously, better than \(\zeta^2\)
\(\zeta>2\): These seemed to work reasonably well.
Ultimately, the justification for using this formula is that it seems to work quite well.