Estimating age-time-dependent malaria force of infection accounting for unobserved heterogeneity

Summary

Despite well-recognized heterogeneity in malaria transmission, key parameters such as the force of infection (FOI) are generally estimated ignoring the intrinsic variability in individual infection risks. Given the potential impact of heterogeneity on the estimation of the FOI, we estimate this quantity accounting for both observed and unobserved heterogeneity. We used cohort data of children aged 0·5–10 years evaluated for the presence of malaria parasites at three sites in Uganda. Assuming a Susceptible–Infected–Susceptible model, we show how the FOI relates to the point prevalence, enabling the estimation of the FOI by modelling the prevalence using a generalized linear mixed model. We derive bounds for varying parasite clearance distributions. The resulting FOI varies significantly with age and is estimated to be highest among children aged 5–10 years in areas of high and medium malaria transmission and highest in children aged below 1 year in a low transmission setting. Heterogeneity is greater between than within households and it increases with decreasing risk of malaria infection. This suggests that next to the individual's age, heterogeneity in malaria FOI may be attributed to household conditions. When estimating the FOI, accounting for both observed and unobserved heterogeneity in malaria acquisition is important for refining malaria spread models.

INTRODUCTION

Estimating the burden of malaria and evaluating the impact of control strategies, requires reliable estimates of transmission intensities [1]. Measures of malaria transmission intensity include the entomological inoculation rate (EIR), parasite prevalence and force of infection (FOI) [1–6]. The EIR is defined as the number of infectious bites per person per unit time [2, 7], whereas the FOI is defined as the number of infections per person per unit time [4] or the per capita rate at which a susceptible individual acquires infection [8, 9]. The malaria FOI counts all incident (i.e., new) human malaria infections in a specified time interval regardless of clinical symptoms and recurrent infections [4]. The EIR and FOI are related but differ; the EIR considers the number of infective bites delivered by the mosquito vector, whereas the FOI focuses on the infections acquired by the human host. In theory, there should be close relationship between the EIR and the FOI, especially in children with less developed immunity. In practice, however, there is a discrepancy between the two because not every infectious bite results in an infection due to various factors [10]. The efficiency of transmission can be estimated by taking the ratio of the two measures, i.e., the ratio of the EIR to the FOI, the number of infectious bites required to cause an infection [10]. A smaller ratio of the EIR to the FOI implies higher transmission efficiency. Most studies have shown that malaria transmission is highly inefficient [4]; whereas more recently malaria FOI has been estimated from serological data [1, 11] by detecting past exposure to malaria infection, here we focus on estimating malaria FOI from parasitaemia data [12–14].

Despite well-recognized heterogeneity in malaria transmission [15, 16], the FOI is often estimated ignoring intrinsic variability in the individual risk of malaria infection. Heterogeneity in malaria infection arises due to variability in risk factors, including environmental, vector and host-related factors [17]. Taking these sources of heterogeneity into account [15, 17] in population-based epidemiological studies has been shown to be important [8].

Ronald Ross first published a mathematical model for malaria transmission in 1908 [16, 18]. This model was only firmly established in 1950 by the work of George Macdonald who used Ross's idea [16]. The ‘Ross–Macdonald’ model describes a simplified set of concepts that serves as a basis for studying mosquito-borne pathogen transmission [16]. Using this concept, mathematical methods to estimate the FOI in relation to the EIR have been proposed by, e.g., Smith et al. [3, 4], Keeling and Rohani [19] and Aguas et al. [20]. Some of the parameters involved in these models are often unknown and should be estimated from data [21]. A solution proposed by Ross in 1916 is to iterate between two modelling frameworks, that is, mathematical and statistical models [21, 22]. The major difference in these two is that the mathematical models (priori) are based on differential equations describing the biological mechanism and causal pathway of transmission, whereas the statistical models (posteriori) start by the statistical analysis of observations and work backwards to the underlying cause [21]. These two frameworks complement each other and, here, we provide an explicit link between them.

In this paper, we use the well-known generalized linear mixed model (GLMM) framework (see, e.g., [19]) to estimate the point prevalence accounting for both observed and unobserved heterogeneity and show how the FOI can be obtained from the point prevalence based on a mathematical Susceptible–Infected–Susceptible (SIS) model. We derive an expression and easy-to-calculate bounds of the FOI for varying parasite clearance distributions. Our results can be used to refine mathematical malaria transmission models.

METHODS

Source of data

The results in this paper are based on cohort data from children aged 0·5–10 years in three regions in Uganda: Nagongera sub-county, Tororo district; Kihihi sub-county, Kanungu district; and Walukuba sub-county, Jinja district. The data were collected as part of the Program for Resistance, Immunology, Surveillance and Modelling of malaria (PRISM) study. The study regions are characterized by distinct transmission intensities. The EIR was previously estimated to be 310, 32 and 2·8 infectious bites per unit year, respectively, for Nagongera, Kihihi and Walukuba [6]. The study participants were recruited from 300 randomly selected households (100 per region) located within the catchment areas. Data were routinely collected every 3 months (routine visits) and for non-routine clinical (symptomatic) visits. Individuals were tested for the presence of Plasmodium parasites using microscopy from August 2011 to August 2014 (3 years). All symptomatic malaria infections were treated with artemether–lumefantrine (AL) anti-malarial medications. More detailed information regarding the study design can be found in Kamya et al. [6]. Given that for clinical visits the sampling process is outcome dependent (see the ‘Discussion’ section), the analysis here is restricted to the planned routine visits yielding unbiased estimates (simulation study, not shown).

The SIS model, point prevalence and FOI

A simplified version of malaria transmission can be described using the so-called SIS compartmental transmission model. This mathematical model classifies the population into two compartments, i.e., the susceptible (S) and the infected (I) class, which can be graphically depicted as shown in Figure 1.

Fig. 1. A schematic diagram of the SIS compartmental model illustrating the simplified dynamics for malaria transmission.

Here, the rate λ(t) at which individuals leave the susceptible state S at time t and flow to the infected state I, as they are infected with malaria parasites, is referred to as the FOI. Furthermore, γ represents the time-invariant clearance rate at which individuals regain susceptibility after clearing malaria parasites from their blood. Let s(t) denote the proportion of susceptible individuals in the population and i(t) the proportion of infected individuals at calendar time t, i.e., the (point) prevalence, then the following set of ODEs (ordinary differential equations) describes transitions in the compartmental SIS model:

(1) $$\left\{ {\matrix{ {{s}^{\prime}(t)\; = \; - \lambda (t)s(t)\; + \; \gamma i(t),{\rm \; \;}} \cr {{i}^{\prime}(t)\; = \; \lambda (t)s(t)\; - \; \gamma i(t).\; \; {\rm \;}} \cr}} \right.$$

As individuals are either susceptible to infection or malaria infected (at least in the aforementioned simplified SIS model), we have s(t) = 1 − i(t). Substituting this expression for s(t) into (1) yields:

(2) $${\rm \lambda} ({\rm t}) = \displaystyle{{i(t)\gamma \; + \; {i}^{\prime}(t)} \over {1{\rm \;} - {\rm \;} i(t)}},$$

where i′(t) is the derivative of the point prevalence with respect to t. The FOI λ(t) can thus be estimated using an estimate for the prevalence i(t) and the clearance rate γ.

Relaxing the assumption of an exponentially distributed parasite clearance distribution in the SIS model can be done by dividing the I compartment into J sub-compartments, such that infected individuals move from the first sub-compartment I ₁ to the second I ₂, and later to the J ^th sub-compartment I _J during the different phases of clearing malaria parasites. Using identical rates γ for the transitions between these sub-compartments and for moving from I _J back to the S compartment results in an Erlang distribution with shape parameter J and rate γ for the time spent in all of the sub-compartments [23]. It is easily shown that equation (2) yields an upper bound for the FOI when compared with the aforementioned Erlang clearance distribution (see the Appendix). A lower bound is readily obtained by taking γ = 0 in equation (2) (SI model – see the Appendix). The FOI is thus bounded by [λ _L (t), λ _U (t)] = [i′(t)/(1 − i(t)), ((i′t) + γi(t))/(1 − i(t)))].

Estimates for both the exponential assumption (upper bound) as well as the lower bound are presented in this paper. In order to estimate the prevalence π(t) ≡ i(t), we use a GLMM to account for individual- and household-specific clustering. This will enable us to explicitly model the observed and unobserved heterogeneity in the acquisition of malaria infection.

Generalized linear mixed model

GLMMs extend the well-known generalized linear models by explicitly taking into account (multiple levels of) clustering of observations [24].

Let Y _ijk denote the binary response variable indicating parasitaemia in the blood (1 if parasites are present – malaria infected; and 0 if not – malaria uninfected) for the i ^th individual nested in the J ^th household at the k ^th visit. Similarly, let X _ijk be a (p + 1) × 1 vector containing covariate information on p independent variables, and Z _ijk be a q × 1 vector of information associated with q random effects. Given the subject-specific random effects b _ij and the covariate information X _ijk , the random variables Y _ijk |X _ijk are assumed to be conditionally independent with conditional mean π(X _ijk | b _ij ) = E(Y _ijk |X _ijk , b _ij ) = P(Y _ijk  = 1|X _ijk , b _ij ). The GLMM relates the conditional mean to the covariates X _ijk and Z _ijk as follows:

(3) $$\eqalign{g\left[ {\pi (X_{ijk}{\rm \vert} {\bi b}_{ij})} \right] & = g\left[ {P(Y_{ijk} = 1{\rm \vert} X_{ijk},\; \;{\bi b}_{ij})} \right] \cr & = \; X_{ijk}^T {\bi \beta} + Z_{ijk}^T {\bi b}_{ij}.}$$

Here, g is a monotonic link function (e.g., logit, cloglog and log); $\eta (X_{ijk}\; \vert\; {\bi b}_{ij})=X_{ijk}^T {\bi \beta} + Z_{ijk}^T {\bi b}_{ij}$ is the linear predictor with β a vector of unknown regression parameters for the fixed effects; ${\bi b}_{ij} \sim N(0,{\bf D})$ a vector of subject-specific random effects for subject i in household j for which elements are assumed to be mutually independent; and ${\bf D}$ a q × q variance–covariance matrix [25].

Using equations (2) and (3), the FOI can be obtained using different link functions. Table 1 presents the prevalence and FOI when selecting either the logit, cloglog or log-link function in the GLMM.

Table 1. General structures for the FOI according to different link functions in a GLMM framework

η refers to the linear predictor η(X _ijk | b _ij ) and η′ represents the derivative of the linear predictor with respect to the predictor of interest.

Flexible parametric modelling

In a parametric framework such as the GLMM, fractional polynomials provide a very flexible modelling tool for the linear predictor η(X _ijk | b _ij ) [21, 26, 27]. In this paper, a GLMM using a fractional polynomial of degree one with regard to age, with power p selected from a grid (−3, −2, −1, −0·5, 0, 0·5, 1, 2, 3) using Akaike's information criterion (AIC), is used [28]. More precisely, we use

(4) $$\eqalign{\eta (X_{ijk}\; \vert \; {\bi b}_{ij}) & = \eta (a_{ijk}\;, {\rm \;} l_{ij} \vert {\bi b}_{ij}) \cr & = \beta _0 + \beta _1\; {\rm age}_{ijk}^p + \beta _2l_{ij} + b_{0i\left( j \right)} \cr & \quad + b_{1i(\,j)}{\rm age}_{ijk}^p,}$$

where b _0i(j) is the nested random intercept and b _1i(j) is the nested random slope for age. Nesting is done to explicitly acknowledge that individuals make up households. Furthermore, shifted year of birth: l _ij , defined as the child's birth year minus the birth year of the oldest child in the cohort (i.e., baseline year 2001), is used in the model to account for the (calendar) time effect since [calendar time] = [birth year] + [age]. The linear predictor (4) can be further extended to include additional covariates.

Age-time-dependent FOI

In equation (3), the conditional mean π(X _ijk | b _ij ) is the point prevalence conditional on the random and fixed effects. In this paper, we use the logit-link function, which enables easy calculation of the ICC (intra-cluster correlation coefficient) through an approximation indicating how much the elements within a cluster are correlated [24, 29, 30].

The age-time-dependent FOI, conditional on random effects, is estimated by plugging in the parameter estimates obtained from the final fit in equation (2). More specifically, using a logit-link, the conditional age-time-dependent FOI is estimated as follows:

(5) $$\eqalign{{\hat \lambda}_{l_{ij}}(a_{ijk} \vert {\bi b}_{ij}) & = \hat \gamma e^{\hat \eta (a_{ijk},\;l_{ij} \vert {\bi b}_{ij})} \cr & \quad + \hat \eta ^{\prime}(a_{ijk},\;l_{ij} \vert {\bi b}_{ij})\hat \pi _{l_{ij}}(a_{ijk},\;l_{ij} \vert {\bi b}_{ij})}$$

where $\hat \gamma $ is an estimate for the clearance rate and $\hat \pi _{l_{ij}}(a_{ijk},\;l_{ij}\vert{\bi b}_{ij})$ is the estimated age- and time-dependent conditional prevalence. For the lower boundary of FOI, $\hat \gamma e^{\hat \eta (a_{ijk},\;l_{ij}\vert{\bi b}_{ij})}$ is omitted in equation (5). In the above expression, an estimate for the clearance rate γ is required. Previously, Bekessy et al. [12] estimated annual clearance rates of 1·643, 0·584 and 0·986 years⁻¹ for children aged <1, 1–4 and 5–8 years, respectively. Later, Singer et al. [14] estimated these rates as 1·917, 1·425 and 2·364 years⁻¹ for ages <1, 1–4 and 5–8 years, respectively. Sama et al. [13] estimated a constant annual clearance rate of 1·825 years⁻¹ by assuming an exponential distribution for infection duration or parasite clearance. Most recently, Bretscher et al. [31] studied the parametric distributions of the infection durations using Ghanaian data, and concluded based on AIC that a Weibull distribution gave a better fit to the data followed by a gamma distribution, while an exponential one was performing worst. In this paper, we use both exponential and Erlang clearance distributions to derive estimates for the malaria FOI obtained based on the aforementioned clearance rates as distributional parameters.

Often, an investigator may wish to observe population averaged estimates. Under the random effects framework, this can be achieved by taking the expectation of the conditional estimates (e.g., the FOI in (5)) resulting into unconditional or marginal estimates. Using the logit-link function, the unconditional (population) FOI is given by

(6) $$\eqalign{\lambda _{l_{ij}}(a_{ijk}) & = E(\lambda _{l_{ij}}(a_{ijk} \vert {\bi b}_{ij})) \cr & = E(\gamma e^{\eta (a_{ijk} , l_{ij}\vert {\bi b}_{ij})} \cr & \quad + {\eta} ^{\prime}(a_{ijk},l_{ij} \vert {\bi b}_{ij}){^\ast}\pi _{l_{ij}}(a_{ijk},l_{ij} \vert {\bi b}_{ij})).}$$

Calculation of the marginalized FOI in (6), requires integrating out the random effects, b _ij , over their fitted distribution. This can be done using numerical integration techniques or based on numerical averaging [24].

Model selection

Model building was done using both AIC [32] and a likelihood ratio test for the random effects based on the appropriate mixture of chi-square distributions [33]. Backward model building was performed starting with the random effects and then the fixed effects. The covariates considered in the model-building process included study site, age, time since enrolment, shifted birth year (i.e., shifted birth year = birth year–birth year of the oldest child), previous use of AL treatment and the infectious status at the previous visit. The covariates, ‘time since enrolment’ and ‘shifted birth year’ were generated to represent the calendar time, albeit we preferred the latter one since participants were not enrolled at the same time point.

RESULTS

Of 989 children, recruited between August 2011 and August 2014, 334 (33·8%), 355 (35·9%) and 300 (30·3%) were from Nagongera, Kihihi and Walukuba, respectively. The baseline parasite prevalence among children aged below 5 years was 38·2%, 12·8% and 9·5% for Nagongera, Kihihi and Walukuba, respectively. The monthly parasite prevalence was higher in Nagongera (range: 26·7–68·4%) followed by Kihihi (range: 7·0–68·0%) and lastly by Walukuba (range: 0–42·9%). Other summary statistics are presented in Table 2. In general, the prevalence was higher among older children (5–10 years).

Table 2. Recruited number of children, baseline and monthly parasite prevalence, by study site and age group

^a Parasite prevalence.

The parasite prevalence increases with age particularly for children <3 years of age and after 7 years of age a decrease is observed (Fig. 2, panel A). The prevalence increases with calendar time in Kihihi with increasing variability, while it decreases in Walukuba, and slightly increases in Nagongera (Fig. 2, panel B). These observations suggest a difference in malaria infection risk between the three study sites. Also, the infection risk seems to vary with age and calendar time and it tends to take different trends between sites indicating a possibility for a site-time interaction effect. The relationship with age seems to be non-linear. These observed effects were taken into consideration when building the GLMM.

Fig. 2. Proportion of children infected with malaria parasites (parasitaemia) in a cohort followed for 3 years, by study site (Nagongera, Kihihi and Walukuba) in Uganda based on data from August 2011 to August 2014 with the size of the dots proportional to the number of observations. (A) Observed parasitaemia varying with age; (B) observed parasitaemia varying with calendar time.

The mean structure in our model consists of a fractional polynomial of age with power −1 (selected based on AIC) and the following covariates (based on significance testing at 5% significance level): shifted year of birth; infection status at previous visit and AL use; and study site. Goodness-of-fit of the final model was assessed using the ratio of the generalized Chi-square statistic to its degrees of freedom. A value of 0·74 was obtained, which is fairly close to 1, indicating that the variability in these data seems to be adequately modelled and little residual overdispersion remains present [34].

The parameter estimates, standard errors, and corresponding test results of the final GLMM fit are shown in Table 3. More details about the candidate models can be found in the Appendix (Tables A1 and A2) together with the fitted conditional and marginal prevalences for the different AL use categories (Fig. A2). The results in Table 3 show an overall significant effect of age and shifted year of birth; the effect of age and shifted year of birth is non-significant and borderline significant, respectively, for Walukuba, whereas the effect of age is significant for Kihihi and Nagongera. Shifted year of birth is significant for Kihihi and non-significant for Nangongera. There is significant heterogeneity in the rate of acquiring malaria infection between households (Walukuba: variance = 2·80; Kihihi: variance = 1·16; Nagongera: variance = 0·21) and between household members (variance = 0·24). The intra-household correlation coefficients are 0·44, 0·25 and 0·06 for Walukuba, Kihihi and Nagongera, indicating moderate, low and very low correlation within households, respectively. The intra-individual correlation coefficients are 0·04, 0·05 and 0·06 for Walukuba, Kihihi and Nagongera, respectively, indicating very low correlation in all sites.

Table 3. Estimates of the fitted GLMM using a fractional polynomial of degree 1 for age and a logit-link function

^a Birth year – min (birth year).

^b Note that the OR here should be interpreted at the Age⁻¹ level.

Based on the final model fit and using equations (5) and (6) both the conditional (given the random effects) and marginal (population averaged) FOIs can be calculated provided that γ can be estimated. However, estimating γ from the same data is not possible due to an identifiability problem: two or more distinct values of γ give rise to the same (log)likelihood (see Fig. A1 in the Appendix). Therefore, we use γ equal to the annual clearance rates given by Bekessy et al. [12] as 1·643, 0·584 and 0·986 years⁻¹ for children aged <1, 1–4 and 5–10 years, respectively, to calculate the conditional and marginal FOIs. We further conduct a sensitivity analysis by considering different clearance rates ranging from 0 to 3 motivated by the ranges estimated by Bekessy et al. [12], Singer et al. [14], Sama et al. [13] and Bretscher et al. [31] (see Fig. 5, top row). As discussed before, we also provide lower bounds for the FOI.

Figure 3 shows estimates for the marginal FOI together with the corresponding lower bound estimates. We focused on children who were born in the baseline year for graphical reasons. Similar plots were obtained (not shown) for other birth years. Estimates for the lower boundary of the FOI were higher in Nagongera followed by Kihihi and Walukuba. For Nagongera and Walukuba, the lower bound for the FOI was highest for children aged below 1 year and least in those aged 5–10 years, yet. In Kihihi, it is highest among those aged 1–4 years.

Fig. 3. The lower bound (green) for the marginal annual FOI and the difference between upper and lower bound (yellow) with full bar showing the upper bound for the FOI, by study site, age group (A, <1 year; B, 1–4 years; and C, 5–10 years) and the infection status at the previous visit and past use of AL (negative and no AL in the past (left column), negative and AL in the past (second left column), symptomatic (second right column) and asymptomatic (right column)) for children assumed to be born in the baseline year (2001). Top row, Nagongera; middle row, Kihihi; bottom row, Walukuba.

Figure 3 further shows that in Nagongera and Kihihi, the estimates for the marginal FOI were highest among children aged 5–10 years; yet in Walukuba it was highest among those aged below 1 year. The values for the marginal FOI obtained using the upper boundary estimator, stratified by site, age group and the previous infection status and use of AL are given in Table A3 in the Appendix. At the extreme, the previously symptomatic children acquire up to four infections per year in Nagongera, and eight infections per year both in Kihihi and Walukuba. Overall, the FOI is highest among the asymptomatic children and smallest among previously symptomatic children across all age groups and sites (Fig. 3 and Table 3A). Although Figure 3 clearly shows the impact of different distributional assumptions with regard to the clearance time, the lower and upper bound estimates do not fully capture uncertainty around the point estimates. In Table A4 of the Appendix, we show the 95% confidence bounds for the age- and time-dependent FOI.

Figure 4 (top row) shows the predicted conditional FOIs for 50 randomly selected individual profiles at each of the three sites based on the lower boundary estimator for the FOI. For graphical purposes, we focused on subjects who were symptomatic at the previous visit and who were born in the baseline year. However, similar plots are obtained for other levels of the infection status at the previous visit and for different birth years. Figure 4 (bottom row) shows the predicted marginal FOIs again based on the lower boundary estimator, by age (continuous scale) and infection status at the previous visit and past AL use. In general, the lower boundary estimator indicates that younger children have the greatest FOI. In all sites, individuals that were asymptomatic at the previous visit have the highest FOI, regardless of age. The depicted conditional FOI curves show that individuals have different profiles, indicating substantial unobserved heterogeneity. The increasing trend in the FOI from 6 months of age is likely attributed to loss of maternal immunity in infants [35].

Fig. 4. Top row: Individual-specific evolutions for the conditional annual FOI obtained using the lower boundary estimator, by study site for children assumed to be symptomatic at the previous visit and who were born in the baseline year (2001). Bottom row: The marginal annual FOI, obtained using the lower boundary estimator, by study site and the infection status at the previous visit and past use of AL (negative and no AL in the past (solid lines), negative and AL in the past (dotted lines), symptomatic (dash-dotted lines) and asymptomatic (long-dashed lines)). Left column, Nagongera; middle column, Kihihi; right column, Walukuba.

Figure 5 (top row) shows the marginal FOIs for different clearance rates ranging from 0 up to 3 years⁻¹ (y-axis). For graphical purposes, and without loss of generality, we again focused on subjects who were symptomatic at the previous visit and who were born in the baseline year. The colour gradient from green (dark) to brown (light) in Figure 5 (top row) corresponds to an increasing FOI. The figure indicates that in Nagongera and Kihihi, children who are below 1 year of age have a lower FOI (green colour) regardless of the presumed clearance rate. Also, in Nagongera and Kihihi, the risk for malaria infection increases with increasing clearance rate, except for the younger children <1–2 years. In Walukuba, the FOI increases with increasing clearance rate regardless of age.

Fig. 5. Top row: The marginal annual FOI (contour lines) considering different values for the clearance rate ranging from 0 to 3 years⁻¹ by study site for individuals assumed to be symptomatic at the previous visit and were born in the baseline year. Bottom row: The marginal annual FOI, obtained using the upper boundary estimator, for individuals assumed to be symptomatic at the previous visit, by study site, birth year (2001, 2004, 2007 and 2010) and by age group (A: <1 year, B: 1–4 years, and C: 5–10 years). Left panel, Nagongera; middle panel, Kihihi; right panel, Walukuba.

Figure 5 (bottom row) shows how the FOI varied with age groups (A ⩽ 1 year, B = 1–4 years, C = 5–10 years) and calendar time among subjects assumed to be symptomatic at the previous visit. In Kihihi, the risk of acquiring a new malaria infection is slightly higher for children born in 2010 compared with those born in earlier years across age groups but not for Nagongera and Walukuba. This would be expected since children born at a later year are younger than those born at an earlier year, and hence are at a higher risk of infection.

DISCUSSION

In this paper, we use data from a cohort study to estimate the malaria FOI among Ugandan children while accounting for observed and unobserved heterogeneities. The results clearly demonstrate the existence of heterogeneity in the acquisition of malaria infections, which is greater between households than between household members. These observations emphasize the claim by White et al. [17] that heterogeneity in malaria infection can arise due to several unobserved factors, including environmental, vector and host-related factors. This implies that estimating the malaria transmission parameters assuming homogeneity in the acquisition of infection may yield misleading results.

The findings were based on the use of a readily available statistical method, the GLMM, which takes into account heterogeneity between individuals and households in the acquisition of malaria infection. In particular, a fractional polynomial of age of degree 1 and power −1, adjusted for the calendar time, by means of the so-called ‘shifted birth year’ (i.e., shifted birth year = birth year–birth year of the oldest child), and other covariates, was considered. The fractional polynomial was chosen because it provides a very flexible modelling tool while retaining the strength of a parametric function. The random slope effects for the fractional polynomial function of age resulted in negative estimates for the FOI, which are biologically implausible and therefore the random slopes were dropped. This could be perceived as a drawback of using the GLMM in combination with fractional polynomials and a more mechanistic approach in which heterogeneity is taken into account at different levels could prove valuable here (further research). When allowing for serial correlation in the model through the specification of an AR(1) correlation structure, the model failed to converge, indicating that too little information was available in the PRISM data to accommodate serial correlation, at least when assuming that the AR(1) assumption is appropriate. An in-depth investigation thereof is an interesting topic for further research.

Based on the SIS model, we derived an expression relating the FOI to the prevalence for infectious diseases such as malaria where we cannot assume lifelong immunity. This expression is an extension of the one proposed by Hens et al. [21] for a so-called SIR model assuming lifelong immunity after recovery, an assumption, which is untenable for malaria. A compartmental model, which can account for temporally recovery due to prior use of treatment (induced immunity) or due to previous exposure to infection (acquired immunity), that is, Susceptible–Infected–Recovered(Treatment)–Susceptible (SIR(T)S), would potentially offer a better alternative compared with the more restrictive SIS model. However, an SIR(T)S model does not yield a closed-form expression for the point prevalence, and hence, for the FOI. Nevertheless, the derivations are approximately valid for an SIR(T)S model with short recovery duration (derivations not included here). Consequently, we focused on the SIS model, albeit that we adjusted for the previous infection status and treatment in our model. The standard SIS compartmental model assumes that the clearance rate is exponentially distributed. We derived two estimators for the FOI, which provide a lower and upper boundary for the FOI based on different Erlang distributions for the clearance rate. The lower boundary approximately holds for a scenario in which the clearance rate is small compared with the FOI. Although mathematical models encompassing more complicated and more realistic transmission dynamics for malaria could be considered, we defer their treatment to future research in which we will combine Nonlinear Mixed Model (NLMM) methodology and numerical approaches for the estimation of the model parameters in the presence of unobserved heterogeneity.

The temporal inhomogeneity observed in the data is not in contradiction with the SIS model we used. Heterogeneity, age and temporal aspects are addressed in the GLMM, through the specification of random effects as well as age- and calendar time variables; whereas derivations from the SIS model under endemic equilibrium enable the estimation of the age- and time-dependent FOI from the estimated age- and time-dependent parasite prevalence. Furthermore, estimation of the reproduction number can be done when focusing on the underlying mechanistic modelling of the FOI. However, we deem this to be beyond the scope of this specific manuscript. Seasonality is not explicitly modelled here; however, inclusion of a covariate describing the amount of rainfall, due to the absence of a clear distinction between the different seasons, and based on additional information (not part of the PRISM data) would be an interesting topic for further research.

When the clearance rate is considered negligible, the rate at which children get infected is highest among those between 1 and 2 years. When the clearance rate is non-negligible, the infection rate is higher among children older than 5 years in areas with high and medium transmission (e.g., Nagongera and Kihihi) and higher in children below 1 year in areas with low transmission (e.g., Walukuba). In Kihihi, the FOI was least for children aged <1 year and it is observed to increase as children grow up from 6 months to 1 year. This could be explained by the fact that children lose maternal immunity in their first year of life [35], which puts them at an increased risk of malaria infection. The higher FOI among children aged 5 years and older could be explained by the fact that these children are often asymptomatic malaria cases and are rarely treated, which makes them reservoirs for infections. This finding conquers with the work by Walldorf et al. [36] who reported that children aged 6–15 years were at higher risk of (asymptomatic) infection compared with the younger ones. They concluded that older children represent an underappreciated reservoir of malaria infection and have less exposure to antimalarial interventions.

A higher risk was seen among children in Nagongera compared with those in Kihihi and Walukuba with no significant difference between the latter two sites. This could be explained by the fact that Nagongera is a predominantly rural area with many semi-structured houses and many mosquitoes compared with Walukuba or Kihihi as was noted by Kilama et al. [5]. Our results also demonstrated the importance of prior treatment in lowering infection risk due to the post-treatment prophylactic effect of longer acting anti-malarials, such as AL. For example, children who were previously treated with AL (the symptomatic malaria cases) had a lower risk of getting re-infected compared with those who were asymptomatic or negative at the previous visit.

This study has two major limitations. First, the analysis was based on results of parasite prevalence determined by microscopy, which is less sensitive than molecular methods such as PCR (polymerase chain reaction) or LAMP (loop-mediated isothermal amplification method) [37, 38]. Thus, sub-microscopic infections would not have been detected. This could have resulted into lower estimates of the FOI. In addition, genotyping was not performed to distinguish new and recurrent infections. As a result, the FOI among individuals who were asymptomatic at the previous visit could have been overestimated. Secondly, the unscheduled clinical visits by the symptomatic individuals were triggered by the study outcome (i.e., parasitaemia). This creates a dependency between the observation-time and outcome processes. This dependence, if not accounted for, has a potential to introduce bias in the model estimates and hence in the estimation of the FOI. This bias was avoided by dropping clinical visits and by using only routine data, although the infection status and use of treatment during clinical visits was accounted for in the model. This implies that the analysis used less data than was actually available. The latter limitation will be dealt with in future research by modelling both the outcome and the observation-time processes concurrently using a joint model [39, 40].

To conclude, we have used longitudinal data from a cohort of Ugandan children to estimate the malaria FOI accounting for both observed and unobserved heterogeneity. First, we show how the FOI relates to parasite prevalence assuming an SIS compartmental model and giving both lower and upper boundaries thereof by relaxing the exponential assumption with regard to the parasite clearance distribution. We estimated the parasite prevalence using a GLMM, whose estimates were used to obtain an estimate for the FOI. The malaria FOI was highest among children aged 1–2 years based on the lower boundary estimator, and it was higher among children older than 5 years in areas of high and medium transmission based on the upper boundary estimator. In a low transmission setting, the FOI was highest in children aged below 1 year regardless of the boundary estimator for the FOI. The FOI varied between study sites highest in Nagongera and least in Walukuba. Heterogeneity increases with decreasing FOI and is greater between households than household members. We recommend that estimating the malaria FOI should be done accounting for both observed and unobserved heterogeneity to enable refining existing mathematical models in which the FOI may be unknown.