Mixture models

In the context of malaria and other infectious diseases, mixture models are developed under the assumption that the population of interest is indeed a mixture of latent seropositive and seronegative populations [Reference Sepúlveda4, Reference Cook20]. More formally, let Y _i denote the log-transformed antibody measurement for the i-th individual. Let S ⁺ and S ⁻ be a shorthand notation for ‘seropositive’ and ‘seronegative’ classifications, respectively. Assuming independent and identically distributed realisations for a sample of n individuals, we write the density function of Y _i as

(1)$$f\lpar y_i\rpar = \prod\limits_{i = 1}^n \lsqb \lpar 1-p\rpar f_{S^-}\lpar y_i\semicolon \;\mu _{S^-}\comma \;\sigma _{S^-}^2 \rpar + pf_{S^ + }\lpar y_i\semicolon \;\mu _{S^ + }\comma \;\sigma _{S^ + }^2 \rpar \rsqb \comma \;$$

where $f_{S^ + }$ is a univariate log-Gaussian distribution with mean $\mu _{S^ + }$ and variance $\sigma _{S^ + }^2$ for the S ⁺ population, and analogously for S ⁻; finally, pis the probability of being S ⁺.

Let C _i and $C_i^\ast$ denote the random variables representing classification based on the mixture model and true classification of the i-th individual, respectively. One approach is to define a seropositivity threshold, usually $\mu _{S^-} + 3\sigma _{S^-}$, above which C _i is S ⁺, and S ⁻ if below [Reference Sepúlveda4, Reference Bollaerts15, Reference Hens16, Reference Bosomprah19, Reference Cook21]. An alternative, more elaborate, approach is to first calculate the probability of belonging to group $C_i^\ast$, conditional on the antibody measurement Y _i = y _i, i.e.

(2)$$\eqalign{P\lpar C_i^\ast = S^ + \vert y_i\rpar & = \displaystyle{{\,pf_{S^ + }\lpar y_i\semicolon \;\theta _{S^ + }\rpar } \over {\lpar 1-p\rpar f_{S^-}\lpar y_i\semicolon \;\theta _{S^-}\rpar + pf_{S^ + }\lpar y_i\semicolon \;\theta _{S^ + }\rpar }} \cr P\lpar C_i^ \ast = S^-\vert y_i\rpar & = 1-P\lpar C_i^\ast = S^ + \vert Y_i = y_i\rpar.} $$

Based on two probability thresholds, c ⁻ and c ⁺, the classification C _i is

(3)$$C_i = \left\{{\matrix{ {S^-} \hskip-19pt\hfill {{\rm if\ }P\lpar C_i^\ast{ = } S^-\vert Y_i = y_i\rpar \le c^-} \hfill \cr I \hfill {{\hskip11pt \rm if \ }c^- \lt P\lpar C_i^\ast{ = } S^-\vert Y_i = y_i\rpar \lt c^ + } \hfill \cr {S^ + } \hfill {{\hskip-21pt \rm if\ }P\lpar C_i^\ast{ = } S^ + \vert Y_i = y_i\rpar \ge c^ + } \hfill \cr } } \right.\comma \;$$

where I is an additional classification label introduced to denote inconclusive cases. In serology analysis, a common approach is to exclude these cases, depending on the type of disease, and report the proportion of inconclusive cases [Reference Bollaerts15, Reference Hens16, Reference Del Fava22].

In malaria serology, most studies favour the first threshold-based approach that does not introduce the classification for inconclusive cases [Reference Yman17–Reference Bosomprah19, Reference von Fricken23–Reference Ashton25]. This is likely due to the nature of antibody responses to malaria infections which result in a large proportion of ‘inconclusive’ cases, as reported by Sepúlveda et al. [Reference Sepúlveda4].

However, both of these threshold-based approaches are prone to misclassification, which can create bias in estimating epidemiological parameters [Reference Sepúlveda4, Reference Bollaerts15, Reference Hens16]. Furthermore, current applications of mixture models in malaria serology analysis do not take into account the age-dependence of antibody levels, and assume that the mixing of S ⁺ and S ⁻ is the same across all ages, which may further exacerbate the issue of misclassification.

The two component mixture Gaussian models also do not account for antibody boosting upon re-exposure to malaria parasites. Sepúlveda et al. [Reference Sepúlveda4] present an extension to the traditional mixture model where more components are added in order to account for this boosting effect. These components can be interpreted as varying degrees of malaria exposure; unexposed, once exposed, twice exposed, etc. Assuming a known number of components, say K, the sampling distribution is given by

(4)$$f\lpar y_i\rpar = \prod\limits_{i = 1}^n {\left[{\sum\limits_{k = 1}^K {\,p_kf_k\lpar y_i\semicolon \;\theta_k\rpar } } \right].} $$

The number of components K is then treated as an additional parameter to estimate using the profile likelihood. However, the interpretation of the components of the model is problematic due to ambiguity about classification rules, particularly when component means are close together. This approach also further compounds the problem of inconclusive cases as they occur across multiple components.

Fig. 1. An illustration of the mixture model showing the bi-modal distributions for the S ⁻ (red) and S ⁺ (blue) populations. The dotted line in (a) shows the seropositivity threshold $\mu _{S^-} + 3\sigma _{S^-}$, above which individuals are classified as S ⁺. The grey rectangle in (b) shows the inconclusive cases as defined by equation (3). In this case, the probability thresholds c ⁻ and c ⁺ have been set to 90%. Individuals below this grey region are classified as S ⁻, while individuals above this region are classified as S ⁺. These data are taken from the Pf AMA1 analysis in section ‘Analysis of malaria serology data from Western Kenya’.

Reversible catalytic models

Following the dichotomisation of the continuous antibody measurements through the application of a mixture model, the resulting S ⁺ and S ⁻ outcomes are modelled using an RCM. A common assumption of the RCM is that individuals are born S ⁻ and, after becoming S ⁺ upon exposure to malaria, can revert to S ⁻ in the absence of exposure. This mechanistic approach is illustrated in Figure 2a. Since antibody data are assumed to represent the cumulative exposure of individuals during their lifespan, the age of individual prior to the sample collection is used as a proxy for historical time.

Fig. 2. (a) This figure is a representation of the reversible catalytic model (RCM) where individuals transition between seronegative (S ⁻) and seropositive (S ⁺) states through the SCR, λ(a)and the SRR, ω. (b) This figure is a representation of the superinfection model (SIM) where individuals can have their antibodies ‘boosted’ through increasing seropositive (S ^+…) states depending on the cumulative exposure to malaria parasites.

Let λ(a) denote the seroconversion rate for an individual at age a and ω the seroreversion rate. According to the RCM, the temporal dynamics that regulate the proportion of S ⁺ individuals of age a, hence p(a), are expressed by the following differential equation

(5)$$\displaystyle{{{\rm d}p} \over {{\rm d}a}} = \lambda \lpar a \rpar \lpar {1-p\lpar a \rpar } \rpar -\omega p\lpar a \rpar.$$

In the above equation, λ(a) is a measure of the underlying transmission intensity which is associated with the gold standard indicator of transmission, the EIR [Reference Corran8], while ω is typically fixed and assumed to be constant [Reference Sepúlveda4]. However, some authors Bosomprah [Reference Bosomprah19] and Akpogheneta et al. [Reference Akpogheneta26] suggest that ω may be age-dependent. Sepúlveda et al. [Reference Sepúlveda4] argue that the malaria serology data often carry little information in the estimation of ω, a problem which will persist also in our novel modelling framework. Hence, throughout this paper, we shall make the working assumption of a constant ω. Note that the reciprocals of λ and ω estimates, i.e. 1/λ and 1/ω, indicate the estimated number of years within which seroconversion and seroreversion would occur, respectively.

Three transmission profiles have so far been proposed to model the seroconversion rate λ(a). The simplest assumes a constant transmission, hence λ(a) = λ for all a. In this case, the differential equation in (5) gives the following solution

(6)$$p\lpar a \rpar = \displaystyle{\lambda \over {\lambda + \omega }}\lpar 1-{\rm e}^{-\lpar {\lambda + \omega } \rpar a}\rpar.$$

In the equation above, the proportion of S ⁺ at older ages reaches a maximum value of about λ/(λ + ω). In other words, in a cohort of an initially malaria-naive population, p(a) will ultimately reach a plateau at which the number of individuals seroconverting is the same as the number of individuals seroreverting [Reference Sepúlveda4, Reference Corran8]. However, these assumptions may be too stringent as they ignore changes in transmission that may be due, for example, to the introduction of control interventions [Reference Sepúlveda4, Reference Cook20, Reference Cook21].

To tackle this issue, one approach is to assume a transmission profile with a sharp drop in transmission at the time of intervention. In this model, two transmission rates are estimated: λ ₁ and λ ₂ which represent the transmission rates before and after the drop, respectively. An alternative approach to account for control interventions is to assume a linear reduction in the seroconversion rate λ(a), rather than a step-change as we have just illustrated. However, in this case, the differential equation in (5) cannot be solved analytically and numerical procedures must instead be used.

In the study by Yman et al. [Reference Yman17], the two transmission profiles that do not assume a constant λ(a) provide a better fit to the data. However, assumption of a step-change or linear drop in λ(a) may be inappropriate in the presence of major or prolonged malaria outbreaks within the historical time-frame considered. In general, the validity of any of these profiles is dependent on a variety of factors, including intervention history, climate and vector characteristics. More recently, Varela et al. [Reference Varela27] propose a model where the number of times that λ changed in the past, which is also estimated from the data.

Where seropositivity is defined using the traditional two-component Gaussian mixture model, there is still the issue of how to account for antibody boosting due to repeated exposure to malaria parasites. Bosomprah [Reference Bosomprah19] suggests an extension to the RCM, which involves creating more seropositive classes in a superinfection model (SIM), similar to the multi-component mixture model described by Sepúlveda et al. [Reference Sepúlveda4]. In this framework, a seronegative individual can transition to the first seropositive class, S ⁺, upon first exposure, and subsequently to a higher seropositive class S ⁺⁺ upon re-exposure, and so on, as illustrated in Figure 2b. The SIM also faces challenges with interpretation of results where initial exposure and boosting between the multiple seropositive classes may be conflated [Reference Pothin3, Reference Sepúlveda4].

Antibody acquisition models

An alternative modelling approach to estimate transmission intensity is to use AAMs [Reference Yman17, Reference Weber18]. Unlike RCMs, AMMs use the full antibody measurements without requiring any dichotomisation of the data. More specifically, AAMs are used to estimate the boosting rate, i.e. the rate at which antibodies are acquired, a marker for transmission intensity [Reference Sepúlveda4, Reference Yman17, Reference Weber18, Reference White28]. Let μ(a) denote the average antibody level in the general population of individuals of age a. Assuming that following exposure to parasites, μ(a) is boosted at a rate γ(a) and assuming a constant decay rate r, we can express this mechanism through the following differential equation

(7)$$\displaystyle{{{\rm d}\mu } \over {{\rm d}a}} = \gamma \lpar a \rpar -r\mu \lpar a \rpar.$$

We can then use the above equation to infer changes in average antibody levels as a function of age a. Finally, in order to fit (7) using likelihood-based methods of inference, the antibody levels of individuals at age a are assumed to follow a log-Gaussian distribution with mean μ(a) and variance σ ² [Reference Sepúlveda4, Reference Yman17, Reference Weber18].

Similar to the way the way seroconversion rates have been modelled in RCMs (section ‘Reversible catalytic models’), previous studies have considered three transmission profiles for the specification of γ(a). The simplest approach assumes that γ(a) = γ is constant which leads to the following solution of (7)

(8)$$\mu \lpar a \rpar = \displaystyle{\gamma \over r}\lpar 1-{\rm e}^{{-}ra}\rpar.$$

Similarly to RCMs, extensions of the AAM assumes either a step-change or linear reduction in the boosting rate γ; see Sepúlveda et al. [Reference Sepúlveda4], Yman et al. [Reference Yman17] and Weber et al. [Reference Weber18] for more details.

Direct comparison of γ and λ from the AAM and RCM, respectively, may not be possible as these estimate different serological indicators. However, Yman et al. [Reference Yman17] find that the AAMs provide a more consistent fit to age-dependent antibody data compared to RCM fit to age-dependent seroprevalence data. Additionally, AAMs provide better precision in parameter estimation and appear to be more robust to sample size reduction. It has been found that AMMs often provide a good fit to serological data in high to moderate transmission settings, where a large proportion of individuals may be seropositive [Reference Yman17], or where an antigen is highly immunogenic, leading to high seropositivity to its antibody in the population [Reference Weber18].

Alternative empirical approaches to model age-dependency

When the interest is in describing the effect of age on the distribution of antibody data, an empirical, rather than mechanistic approach, may provide a better statistical solution. Additionally, the empirical approach outlined in this section can be used to validate the unified mechanistic model by assessing the discrepancy between the age distributions generated by the two modelling approaches.

To this end, we modify the framework introduced in the previous section by replacing the modelling of mixing probability based on RCMs, and the mean level of antibodies based on AAMs, with their empirical counterparts. More specifically, we model the age-dependency in λ(a) and p(a) using a log-linear and logit-linear regression as

(13)$$\eqalign{ \log \left\{{\displaystyle{{\,p\lpar a\rpar } \over {1-p\lpar a\rpar }}} \right\} = \alpha _1 + f_1\lpar a\rpar \cr \qquad \!\log \lcub \mu \lpar a\rpar \rcub = \alpha _2 + f_2\lpar a\rpar \comma \;} $$

where f ₁(a) and f ₂(a) are functions that can be specified with the aid of simple graphical tools, such as scatter plots. The resulting structure of the empirical model is summarised in Figure 3b, and we give examples of this in the application of section ‘Analysis of malaria serology data from Western Kenya’.