2.1 Study Design and Population

The MCS is a nationally representative longitudinal survey of around 19,000 children born in the UK between 2000 and 2002 [18]. To date, seven MCS waves have been conducted at ages 9 months (wave 1 in 2001), 3 years (wave 2 in 2004), 5 years (wave 3 in 2006), 7 years (wave 4 in 2008), 11 years (wave 5 in 2012), 14 years (wave 6 in 2015), and 17 years (wave 7 in 2018). Trained interviewers used standardized methodologies to measure child height and weight (to the nearest 0.1 cm and nearest 0.1 kg, respectively) and data on sociodemographic and family characteristics [19]. Our assessment of transitions between bodyweight categories spanned children aged 3 (when height was first measured) to 17 years. The datasets used for this study are publicly available and details for accessing them are provided in the electronic supplementary material (“Data availability” section).

The study population for our primary analysis included singletons with at least two body mass index (BMI) measurements across MCS waves and complete data on relevant variables (described below). BMI was calculated as kg/m2 and then converted into standard deviation (SD) scores (BMI z-scores) using the British Growth Reference (UK90) for age-adjusted and sex-specific categories [20]. Conversions to BMI z-scores were implemented using Stata 17 statistical software [21, 22]. Children in our study population were then grouped into bodyweight categories for each BMI record based on the UK cut-offs for population monitoring: less than or equal to the 2nd centile (underweight), greater than the 2nd and less than the 85th centile (healthy weight), greater than or equal to the 85th and less than the 95th centile (overweight), and greater than or equal to the 95th centile (obesity) [23]. In sensitivity analyses, we also applied the World Health Organization (WHO) references and classifications for children aged 60 months or younger (underweight [or thinness]: BMI z-score < −2; healthy weight: BMI z-score ≥ −2 and ≤ 2; overweight: BMI z-score > 2 and ≤ 3; and obesity: BMI z-score > 3) and children aged 61 months or older (underweight [or thinness]: BMI z-score < −2; healthy weight: BMI z-score ≥ −2 and ≤ 1; overweight: BMI z-score > 1 and ≤ 2; and obesity: BMI z-score > 2) [24]. This study followed the Strengthening the Reporting of Observational Studies in Epidemiology (STROBE) reporting guideline for cohort studies [26].

2.2 Covariates

The choice of individual-level child, mother, and sociodemographic factors to include in our analysis was informed by the relevant literature [8,9,10, 12, 19, 27]. Child and mother factors were as follows: child's sex (male or female); ethnicity (white and non-white); gestational age at birth (preterm [< 37 weeks], early term [≥ 37 and < 39 weeks], full term [≥ 39 and < 41 weeks], late term [≥ 41 and < 42 weeks], and post-term [≥ 42 weeks]) [28]; mother’s age at birth of child (12–19, 20–29, and 30 plus years); mother’s BMI category during pregnancy (underweight, healthy weight, and overweight/obesity); mother’s frequency of alcohol consumption during pregnancy (monthly or more frequently, less than once a month, and never); and mode of delivery (normal, assisted, planned cesarean section, and emergency cesarean section). Sociodemographic factors were as follows: mother’s highest academic qualification (first/higher degree, diplomas in higher education, A/AS/S levels, O level/General Certificate of Secondary Education (GCSE) grades A–C, GCSE grades D–G, other academic qualifications, and none of these qualifications); and family income categorized as Organisation for Economic Co-operation and Development (OECD)-weighted quintiles [29].

2.3 Statistical Analysis2.3.1 Model Derivation and Validation

Our evaluation was based on a multistate transition modeling approach within the maximum-likelihood estimation framework specified by the MSM function [15] in the statistical software R [30]. Our multistate model aims to estimate transitions between bodyweight categories. It is defined by a continuous-time, finite-state stochastic process, with the Markovian assumption that subsequent movement to another state depends only on the current state and not on past states [14, 31]. The principal outputs from our Markov multistate transition model were transition hazard rates (also known as transition intensities), which represent the instantaneous risk of moving from one state to another, such as from healthy weight to overweight. Transition probabilities could then be derived for cycle lengths of interest (such as weeks/months/years) by taking the matrix exponential of the transition intensity matrix [15].

Our models specified four predefined BMI categories: underweight, healthy weight, overweight, and obesity. We allowed for only clinically plausible transition intensities between adjacent waves of data. For example, even if a child had a healthy weight status in wave 2 and obesity in wave 3, an instantaneous progression was only possible between adjacent states. Therefore, though unobserved, the child must have progressed to an overweight status first. This assumption is consistent with existing studies on childhood BMI transitions [9, 10]. Consequently, our four-state model consisted of six possible transitions: underweight to healthy weight; healthy weight to underweight; healthy weight to overweight; overweight to healthy weight; overweight to obesity; and obesity to overweight. The model structure is illustrated and further described in the electronic supplementary material (“Multistate model structure” section).

We estimated transition intensities and annual transition probabilities for a primary analysis (our base case) and two scenario analyses. The primary model reflected our hypothesis that transition probabilities are not time-homogeneous across all waves, i.e., they depend significantly on a child’s age in the year of observation. For the primary analysis, we extracted a dataset of children with measures of BMI for each wave. We then extracted five corresponding pairwise datasets from adjacent sequential waves as follows: waves 2 and 3 (ages 3 and 5), waves 3 and 4 (ages 5 and 7), waves 4 and 5 (ages 7 and 11), waves 5 and 6 (ages 11 and 14), and waves 6 and 7 (ages 14 and 17). Repeated measures of BMI were required to construct each pairwise dataset such that the same children in the earlier wave of a pair were followed into the later wave of the pair (see the electronic supplementary material: “Dataset extraction: primary analysis” section). A time-homogeneous process was assumed only within these five adjacent pairs of MCS waves with piecewise linkage between successive pairwise datasets. In other words, the transition hazards estimated were constant within pairwise datasets but varied between pairwise datasets (see the electronic supplementary material, “Dataset extraction: primary analysis” section, for details).

Scenarios 1 and 2 assumed that transition probabilities were time-homogeneous across all waves, i.e., throughout childhood. We extracted a single dataset for scenario 1 consisting of individuals with at least one repeated measure of BMI between wave 2 (age 3) and any other wave (see the electronic supplementary material: “Dataset extraction: scenario analyses” section). In scenario 2, we assumed that the transition probabilities between waves 2 and 3 were identical to those between successive pairs of adjacent waves. Therefore, for this scenario, we only extracted a dataset of children with repeated BMI measurements between waves 2 and 3 (ages 3 and 5).

Finally, we derived a validation dataset, which we used to test the performance of the three models by comparing observed and predicted (model) estimates of prevalence. We identified 5486 children with complete BMI measurements from waves 2 to 7, forming the basis of our validation dataset. Subsequently, we normalized and scaled the survey weights attached to each child a hundred times, replicating an additional 218,366 children. We merged the replicated datasets with the original validation dataset, which resulted in 221,436 children. Of these, 0.78%, 72.70%, 15.90%, and 10.63% fell in the underweight, healthy weight, overweight, and obesity categories at baseline (wave 2/age 3), respectively. These baseline percentages were then inputted as baseline cohorts into Markov traces of our derived base-case and scenario analyses to simulate the projected yearly prevalence for each bodyweight category. The validation process is further described in the electronic supplementary material (“Model validation” section), and estimates of weighted and unweighted prevalence by wave using both the UK cut-offs for population monitoring and the WHO cut-offs are reported in the electronic supplementary material (“Observed prevalence” section).

2.3.2 Integrating Complex Survey Design Within Models

The sample design of the MCS data aimed to ensure an adequate representation of the UK childhood population as well as guarantee sufficient sampling for critical subgroups [32]. Consequently, the survey was geographically clustered and disproportionately stratified with oversampling of children in disadvantaged socioeconomic circumstances and, in England, children from minority ethnic backgrounds. Despite these measures, the MCS surveys reported complex response patterns characterized by higher non-response rates for ethnic minorities and families in deprived areas and inter-wave attrition [18, 33]. The MCS database constructors provide survey weights that account for the clustered sample design, the unequal probability of being sampled, survey nonresponse, and adjust for inter-wave attrition [18]. To avoid bias in the estimation of point estimates (transition hazards and probabilities), we incorporated these survey weights into a validated adaption of the MSM R package [34].

Confidence intervals derived without adequate handling of variance will be artificially narrow [34]. Variance estimation for parameter point estimates posed an additional layer of complexity to our estimation strategy because the survey design variables provided by the MCS database constructors are suited to the linearization/Taylor series estimation method [35], which does not lend itself to the MSM function. The computationally intensive replication method we adopted is an alternative to linearization [36]. We derived 398 jackknife (jk-n) replicate weights (matching the number of clusters in the MCS as stipulated in the methods for variance calculation by Valliant and Dever [37]) in Stata [21] for each modeling dataset and then used the replicated weights to calculate estimates of variance [37,38,39]. We accelerated the time-consuming variance estimation process using parallel computing approaches [40]. Technical details of the steps taken in our estimation of replicate weights are described in the electronic supplementary material (“Incorporating complex survey design” section).

2.3.3 Robustness Check for Complete Case Analysis

Similar to other longitudinal studies of childhood excess weight using the MCS, where maximum likelihood estimation formed the basis of the statistical method [41, 42] or where other methods such as regression analyses were applied [27, 43], our analyses were based on complete-case analysis, where we included only children with complete data on the child, mother, and sociodemographic variables, outlined earlier. We tested the robustness of the derived transition hazards from our primary analysis by estimating and comparing transition hazards before and after excluding children with incomplete data on explanatory variables. Our sensitivity analysis entailed the estimation of unadjusted transition hazards before excluding children with incomplete data on the child, mother, and sociodemographic variables.