Methodological guidance for predictor variable selection for adolescent smoking outcomes in Global Youth Tobacco Survey using R and Python

Background

Adolescent smoking remains a major public health concern worldwide, with particularly serious implications in low- and middle-income countries where prevention and cessation resources are often limited [1, 2]. Most adult smokers begin using tobacco during adolescence, making this life stage critical for intervention. The Global Youth Tobacco Survey (GYTS), led by the World Health Organization and the US Centers for Disease Control and Prevention, was designed to address this need by providing standardized, comparable data on youth tobacco use across more than 140 countries [3–5]. Globally, analyses of GYTS data have documented wide variation in smoking prevalence and strong associations with peer smoking, exposure to tobacco marketing, and awareness of tobacco products [5–7]. However, despite the richness of these data, published studies often apply simplified analytical approaches, most commonly descriptive statistics or logistic regression models that use survey weights only, while ignoring clustering and stratification. This practice can lead to biased standard errors, overly narrow confidence intervals, and results that are difficult to compare across countries and studies [8].

Methodologically, GYTS data present several recurring challenges that are not always handled consistently in the literature. The survey uses a two-stage clustered sampling design, requiring analysts to account for weights, strata, and primary sampling units to obtain valid population level estimates [7]. In multi country analyses, additional complexity arises from differences in school systems, questionnaire wording, and variable coding, all of which require careful harmonization [9]. At the same time, the survey includes a large number of potential predictors, creating a tension between building parsimonious models and retaining variables that are theoretically important [10]. Many studies rely on automated variable selection based purely on statistical significance, which can result in the exclusion of key demographic variables such as age and sex, increasing the risk of residual confounding [11]. Although SAS and SPSS dominate the published GYTS literature [12], open source tools such as R and Python offer greater transparency and reproducibility, yet practical guidance on how to use them correctly for complex survey analysis remains fragmented [7].

This study responds to these gaps by presenting a clear and reproducible analytical framework for explanatory modeling of adolescent smoking behavior using GYTS data in both R and Python [12, 13]. The objectives are fivefold. First, we standardize data preprocessing and harmonization procedures for single country and pooled multi country analyses. Second, we demonstrate comparable workflows for survey-weighted logistic regression in R and Python, while recognizing that only the R implementation fully accounts for weights, stratification, and clustering. Third, we apply constrained stepwise variable selection methods that protect a priori demographic variables regardless of marginal statistical significance. Fourth, we assess cross platform reproducibility by comparing results obtained from identical model specifications in R and Python. Finally, we provide explicit justification for predictor selection, linking variables to established behavioral theory and evidence from systematic reviews, rather than relying solely on data driven selection.

The initial predictor variable selection in this framework is guided by Social Cognitive Theory that explains health behaviors as the result of dynamic interactions between personal, behavioral, and environmental factors in adolescence smoking [14, 15]. In the context of adolescent smoking, personal factors include age, sex, and perceptions about tobacco, behavioral factors capture prior experimentation and intentions, and environmental factors reflect peer influence, parental smoking, media exposure, and school-based education. Social Cognitive Theory is particularly well suited to adolescent tobacco use because adolescence is a period of intense social learning, identity formation, and sensitivity to norms [16, 17]. Young people are strongly influenced by what they observe among peers and role models, while repeated exposure to tobacco marketing and media shapes awareness and perceived acceptability of smoking. By organizing predictors within this framework, the analysis moves beyond purely statistical associations [18] and provides a coherent explanation of how social environments, individual characteristics, and observed behaviors jointly shape smoking risk among adolescents [19, 20].

Methods

Data sources and study design

This study used data from the Global Youth Tobacco Survey, a school-based survey designed to produce nationally representative estimates of tobacco related behaviors among adolescents [21–27]. GYTS applies a two-stage cluster sampling design, with schools selected as primary sampling units and classrooms randomly selected within schools [21–27]. All analyses accounted for survey weights, stratification, and clustering to ensure valid population level inference [7]. In this study, we considered two GYTS data analysis scenarios as outlined below:

Study 1: single country analysis, Zambia 2021

The first analysis used the Zambia GYTS 2021 dataset, which initially included 6,499 school-going adolescents [28]. After excluding participants with missing outcome or predictor information required for the complete-case model, the final analytical sample comprised 2,959 adolescents aged 11 to 17 years (see Figure 1). The 2021 GYTS for Zambia was amongst the first surveys to collect data on shisha and electronic cigarettes in the Southern Africa Development Community (SADC) region. The survey instrument consisted of 56 standardized questions covering cigarette smoking, use of other tobacco products, exposure to second-hand smoke, access to tobacco, media and marketing exposure, cessation attempts, and tobacco related knowledge. Zambia used the standard GYTS core questionnaire with limited country specific adaptations, allowing direct use of core variables for most constructs. This analysis was used to demonstrate the full analytical workflow in a context where survey design variables and questionnaire structure were fully aligned.

Study 2: multi country analysis, 2015 to 2019

The second analysis pooled GYTS data from four Sub Saharan African countries collected between 2015 and 2019. These included Ghana in 2017 with 5,430 participants [29], Mauritius in 2016 with 4,123 participants [30], Seychelles in 2015 with 2,453 participants [31], and Togo in 2019 with 3,908 participants [32]. The pooled sample comprised 15,914 adolescents aged 11 to 17 years. These countries were selected to reflect diversity in educational systems and tobacco control contexts while maintaining sufficient overlap in standardized GYTS core variables. Pooling the data allowed assessment of cross-country consistency in predictors of adolescent smoking and enabled robust cross platform validation. Missing data were assessed before model development and were generally low at the individual variable level. Across the pooled dataset of 15,914 adolescents, 2,554 participants (16.1%) had at least one missing value among variables considered for analysis, with country-specific proportions ranging from 11.5% in Togo to 21.7% in Seychelles. Consequently, 13,360 participants (83.9%) had complete information and were retained in the complete-case analytical dataset (see Figure 2).

Predictor variable selection framework

The primary outcome was current cigarette smoking, defined using the standard GYTS item asking: “During the past 30 days, on how many days did you smoke cigarettes?” This corresponds to GYTS core question C7, or the country-specific equivalent where local questionnaire numbering differed. Responses indicating smoking on 0 days were coded as non-current smoking (0), while responses indicating smoking on 1 or more days in the past 30 days were coded as current smoking (1). Participants with missing, refused, or non-classifiable responses to this item were excluded from the analytical model requiring the outcome. For pooled analyses, country-specific variable names and response categories were harmonized to this common binary definition before model fitting. The framework governing the predictor variable selection is presented in Figure 3. Predictor selection followed a structured framework grounded in epidemiology and behavioral science [33]. Four guiding principles informed this process. First, theoretical coherence was ensured by mapping all predictors to the domains of Social Cognitive Theory. Second, empirical support was required, with variables drawn from established GYTS methodology and systematic review evidence [26]. Third, modifiability was considered to include both risk and protective factors relevant for intervention. Fourth, cross cultural validity was prioritized by focusing on standardized GYTS core questions to enable meaningful multi country comparisons.

Demographic core

A core set of demographic variables was defined a priori and retained in all models regardless of statistical significance [34]. These included age, sex, grade or educational level, and geographic region or country. These variables are well established confounders in adolescent tobacco research and are explicitly recommended for adjustment in GYTS analyses. Their definitions, theoretical justification, empirical support, and expected associations are summarized in Table 1. Excluding these variables could bias estimates of modifiable risk factors, particularly when automated variable selection procedures are used. Age was modeled to allow non-linear effects, reflecting evidence that susceptibility peaks in early to mid-adolescence.

Table 1. A priori variables: definition, justification, and expected associations.

Individual level predictors

Individual level predictors captured intrapersonal characteristics related to susceptibility and behavioral intention [35, 36]. These variables correspond to the personal factors' domain of Social Cognitive Theory and are described in Table 2. They included prior experimentation with cigarettes, smokeless tobacco, and shisha, intention to smoke in the next 12 months, and a proxy measure of socioeconomic status. Prior tobacco use reflects reduced self-efficacy to resist and increased familiarity with tobacco products, while intention represents the most proximal predictor of future smoking behavior.

Table 2. Individual-level predictors with associated theoretical and empirical justification.

Social factors

Social environment predictors represented interpersonal influences from peers, family, school, and the wider community [37–39]. These variables, summarized in Table 3, included having close friends who smoke, parental or guardian smoking, exposure to smoking at school, and exposure to smoking in public places. Social Cognitive Theory highlights observational learning as a key mechanism through which adolescents adopt behaviors observed in significant others. Consistent with the literature, peer smoking emerged as the strongest and most consistent social predictor of adolescent smoking [38].

Table 3. Social environment predictors with associated theoretical and empirical justification.

Media and marketing exposure

Media and marketing variables captured environmental influences shaping awareness, social norms, and outcome expectations [40–42]. As shown in Table 4, these included awareness of tobacco products, exposure to tobacco advertising on television and radio, ownership of tobacco branded items, and exposure to anti-tobacco media messages. Awareness was treated as a prerequisite for experimentation, while advertising and branding reflect industry driven normalization of tobacco use [43]. Anti-tobacco media exposure was included as a protective factor that may counterbalance pro tobacco influences.

Table 4. Media and marketing predictors with theoretical and empirical justification.

Policy and school based protection

Structural and policy related variables reflected protective environments aligned with the WHO MPOWER framework [44–46]. These predictors, described in Table 5, included being taught about the dangers of smoking in school and the presence of smoke free policies in enclosed and public places. From a behavioral perspective, these factors represent environmental restructuring that reduces exposure, denormalizes smoking, and supports non-smoking norms. These variables are directly relevant for policy and prevention planning [47].

Table 5. Policy and education predictors with theoretical and empirical justification.

Variables excluded and rationale

Several variables were considered but excluded from the primary analysis due to limited comparability or methodological concerns. These exclusions and their rationale are summarized in Table 6. Mental health [48] indicators and sensation seeking measures were excluded because they are not consistently captured across GYTS waves. Internet and social media exposure variables were excluded due to inconsistent availability in earlier surveys. Product specific variables such as brand preferences were excluded to avoid selection bias, as they apply only to current users [49].

Table 6. Excluded variables and rationale.

Data harmonization procedures

Single country variable mapping

For the Zambia dataset, standardized GYTS core variables were mapped to country specific adaptations using a systematic extraction protocol [9]. Where country specific variables were missing, fallback rules were applied to preserve sample size while maintaining consistency. Age categories were harmonized into single year values to allow flexible modeling and comparability with pooled analyses.

Multi country educational system harmonization

Educational systems varied substantially across the four countries included in the pooled analysis. To enable valid cross-country comparison [50, 51], grade levels were harmonized to a unified UK based year system (see Supplementary material S1). This approach aligns educational stage with age and is widely used in global health research. The harmonization protocol for each country is presented in Table 7, along with justification based on age of entry and progression patterns.

Table 7. Educational system harmonization protocol.

Why harmonize to UK system?

• International comparability: UK educational stages widely understood in global health research.
• Age-grade alignment: UK Years roughly correspond to age (Year 7 = 11–12 years, Year 11 = 15–16 years).
• WHO recommendation: GYTS targeting ages 13–15 corresponds to UK Years 8–10.
• Previous studies: Published multi-country GYTS analyses used similar harmonization [52].

Missing data management

Missing data were assessed before model development for the outcome, a priori demographic covariates, and all candidate predictors. Item-level missingness was summarized separately by country and for the pooled dataset. For ease of interpretation, missingness is presented graphically in Supplementary Figure S1 rather than as a large table. The figure reports the percentage missing for each outcome and predictor variable across the Zambia 2021 dataset, each pooled-country dataset, and the combined pooled dataset. Complete-case analysis was used for the primary models because the purpose of the study was to demonstrate a transparent, reproducible modeling workflow using harmonized variables. However, because complete-case analysis can reduce precision and may introduce selection bias when missingness is associated with exposure or outcome status, the implications of missing data are considered in the limitations.

Software environment

Analyses were conducted using R version 4.2+ with the survey package version 4.5 and Python version 3.10+ with statsmodels version 0.13.1+. All scripts were run using a fixed random seed of 123. The full computational environment is documented using an renv.lock file for R and an environment.yml file for Python. These files, together with the analytical scripts, variable harmonization files, model output tables, renv.lock, and environment.yml, are archived in a public version-controlled GitHub repository and in a DOI-linked Zenodo archive. The GitHub repository is available at https://github.com/andreasngambi-star/GYTS_analysis_R_Python.git, and the persistent Zenodo DOI is 10.5281/zenodo.21254547. Exact package versions and installation instructions are provided in the Supplementary material S4.

Results

Analytical sample derivation

For Zambia 2021, 6,499 adolescents participated in the survey. Of these, 225 participants were excluded because the current smoking outcome could not be classified, and 3,315 were excluded because they had missing information for one or more predictors included in the complete-case modeling framework. The final analytical sample for the Zambia model therefore comprised 2,959 adolescents (see Figure 1). For the pooled multi-country analysis, the initial dataset included 15,914 adolescents from Ghana 2017, Mauritius 2016, Seychelles 2015, and Togo 2019. A total of 2,554 participants had missing information for at least one outcome or predictor variable required for the complete-case model. The final pooled analytical sample therefore comprised 13,360 adolescents. The full participant flow is shown in Figure 2.

Analytical framework implementation

To demonstrate the transferability of the proposed workflow, two related but distinct modeling examples were used. The Zambia 2021 dataset was used to demonstrate the single-country workflow for current cigarette smoking, while the pooled Ghana, Mauritius, Seychelles, and Togo dataset was used to demonstrate the same framework for current e-cigarette use.

Survey design specification

All analyses explicitly accounted for the full complex survey design, including stratification, clustering at the school level, and sampling weights. Applying weights alone was considered insufficient, as ignoring clustering leads to underestimated standard errors and overly narrow confidence intervals. Design based variance estimation using Taylor series linearization was therefore used to obtain valid standard errors and confidence intervals. This approach accounts for the design effect arising from school-based clustering, which is rarely negligible in adolescent surveys.

Why weights alone are insufficient

The majority of published GYTS analyses either apply only sampling weights or ignore weighting entirely [7]. Because GYTS uses stratified two-stage cluster sampling, valid inference requires accounting for both weights and design structure (strata and clusters) [8]. Using weights without design leads to underestimated standard errors and artificially narrow confidence intervals, inflating Type I error.

The design effect quantifies variance inflation from clustering: it increases as average cluster size and intra-class correlation rise. Ignoring clustering implicitly assumes zero intra-class correlation, which is rarely true in school-based surveys [53]. We therefore specify the full complex survey design and use design-based variance estimation (Taylor series linearization) to obtain correct standard errors and confidence intervals (see Supplementary material S2a).

Weighted baseline characteristics

Baseline characteristics of the analytical samples were summarized using survey weights. For each categorical variable, we report unweighted numbers, unweighted percentages, and weighted percentages. These descriptive analyses were conducted separately for the Zambia 2021 sample (Table 8) and the pooled multi-country sample (Table 9) to describe the distribution of participants retained in the analytical datasets after exclusions.

Table 8. Descriptive characteristics of adolescents in the full Zambia 2021 GYTS sample, by region, n = 6,499.

Values are n (%) for categorical variables. Percentages are column percentages.

Table 9. Descriptive characteristics of adolescents included in the pooled multi-country analysis, Ghana, Mauritius, Seychelles, and Togo, 2015–2019, n = 15,914.

Statistical analysis

Constrained stepwise logistic regression

Survey weighted logistic regression was used to model current smoking behavior. A constrained stepwise variable selection approach was applied to balance parsimony with theoretical integrity. A priori demographic variables were protected from removal throughout the selection process [10, 18], while candidate predictors from individual, social, media, and policy domains were evaluated based on improvements in model fit using the Akaike Information Criterion. This approach avoids excluding key confounders while still allowing data driven refinement of the final model. The objective of this variable selection process was not automated discovery, but disciplined reduction under theory-informed constraints. Candidate predictors were pre-specified based on Social Cognitive Theory and existing empirical evidence, and selection procedures were applied solely to reduce redundancy and overfitting while preserving conceptual coherence.

Alternative variable selection methods, including Least Absolute Shrinkage and Selection Operator (LASSO) and Elastic Net regression, were considered during framework development. These approaches are valuable for high-dimensional prediction problems and situations involving large numbers of correlated predictors. However, the primary objective of this study was explanatory modeling rather than prediction. Because GYTS analyses are intended to generate interpretable estimates for policy and public health decision-making, we prioritized a theory-driven framework in which key demographic and contextual variables were retained regardless of statistical significance. Furthermore, implementation of design-based penalized regression remains limited in complex survey settings. Consequently, a constrained stepwise approach based on Akaike Information Criterion was selected to balance model parsimony, interpretability, and compatibility with survey-weighted regression.

To assess the potential for model overfitting, events per variable were calculated as the number of current cigarette-smoking events divided by the number of estimated model parameters. In the final current cigarette-smoking model, there were 1,479 current smokers and 33 estimated parameters, giving an EPV of 44.82. This value exceeds the commonly used minimum threshold of 10 events per variable, suggesting adequate event support for the final explanatory model. Because the analysis was explanatory rather than intended for individual risk prediction, EPV was interpreted as a model-stability diagnostic rather than as evidence of predictive validity.

The multivariable weighted logistic regression model identified several factors associated with current cigarette smoking among adolescents (Table 10). Region, age, sex, and grade were retained a priori, while other candidate predictors were selected using stepwise selection from the full model. In the final model, current cigarette smoking was positively associated with ever trying cigarette smoking, intention to use tobacco in the next 5 years, ownership of an item with a tobacco logo, and seeing teachers smoking at school. Protective associations were observed for being taught about the dangers of tobacco use, believing that long-term tobacco use leads to cancer, and reporting ever use of shisha. The multivariable weighted logistic regression results for current e-cigarette use among adolescents in the pooled sample from Ghana, Mauritius, Seychelles, and Togo are presented in Table 11.

Table 10. Multivariable weighted logistic regression model of factors associated with current cigarette smoking among adolescents in Zambia.

OR, odds ratio; CI, confidence interval.

The model was adjusted using stepwise selection starting from the full model, with region, age, sex, and grade retained a priori. Statistically significant associations at P < 0.05 are shown in bold. The events-per-variable assessment indicated adequate model support: 1,479 current smokers divided by 33 estimated parameters, giving an EPV of 44.82. Model diagnostics for the current cigarette-smoking model were: AIC = 2812.3, BIC = 3010.1, C-statistic = 0.856 (95% CI: 0.841 to 0.870), Brier score = 0.094, and Hosmer-Lemeshow P = 0.104.

Table 11. Multivariable weighted logistic regression model of factors associated with current e-cigarette use among adolescents (n = 15,914) in Ghana, Mauritius, Seychelles, and Togo.

The asterisk (*) identifies the demographic variables—age, sex, and school grade—that were retained a priori in the multivariable model regardless of statistical significance; OR, odds ratio; CI, confidence interval.

Values are adjusted odds ratios from a weighted multivariable logistic regression model of current e-cigarette use. The model was fitted using stepwise selection starting from the full model, with region, age, sex, and grade retained a priori. Statistically significant associations at P < 0.05 are shown in bold. Model diagnostics for the e-cigarette-use model were: C-statistic = 0.847 (95% CI: 0.834 to 0.859), Brier score = 0.050, and Hosmer-Lemeshow χ²(8) = 117.57, P < 0.0001. The Hosmer-Lemeshow result indicates calibration misfit; therefore, this model should be interpreted as an explanatory and methodological demonstration rather than as an individual-level prediction model.

Cross platform implementation and validation

Parallel analytical workflows were implemented in R (see Supplementary material S3a) and Python (see Supplementary material S3b). R used native complex survey methods, while Python applied survey weighted generalized linear models with normalized weights. Identical model specifications were fitted in both platforms, and point estimates were compared to assess reproducibility. Cross-platform validation results are presented in Table 12 and demonstrate near-identical adjusted odds ratios for key predictors. Although Python does not natively account for survey stratification and clustering, the present comparison focused on reproducibility of point estimates. Formal variance estimation and statistical inference were based on the R implementation, which fully incorporated the complex survey design. Therefore, validation was restricted to consistency of point estimates rather than equivalence of statistical inference. Python's statsmodels.GLM with freq_weights implements survey-weighted estimation but does not natively support stratification or clustering adjustments. For full complex survey support in Python, use specialized packages (samplics, surveystat) or interface with R via rpy2. Our validation shows equivalent point estimates but potentially different standard errors.

Table 12. Cross-platform validation of adjusted odds ratios.

The objective of the cross-platform comparison was not to establish statistical equivalence between R and Python, but rather to evaluate whether identical model specifications produced similar substantive conclusions regarding predictor importance and direction of association. Variance estimation remained the responsibility of the R implementation, which fully accounted for the complex survey design.

Results of cross platform analysis

Cross platform validation was carried out using the pooled multi country dataset comprising 15,914 adolescents. When comparable model specifications, convergence criteria, and random seeds were applied, the R and Python implementations produced almost identical effect estimates. Differences in adjusted odds ratios were consistently below 0.01, indicating strong consistency in point estimation across platforms. Although Python does not natively account for survey stratification and clustering, the present comparison focused primarily on reproducibility of point estimates. The adjusted odds ratios obtained from R and Python were virtually identical, indicating consistency in substantive conclusions across platforms. Formal variance estimation and inference were based on the R implementation, which fully incorporated the complex survey design. This validation step strengthens confidence in the reproducibility of point estimates and supports the use of either platform for exploratory modeling, provided that formal variance estimation and inference continue to rely on software capable of incorporating the full complex survey design.

Key methodological findings

The methodological contribution of this study lies in formalizing a constrained, theory-informed variable selection strategy for complex survey data that is reproducible across software platforms and transferable to other multi-stage population surveys. Across both the single-country and pooled multi-country analyses, the framework showed theory-constrained consistency, with a priori demographic variables retained by design and several modifiable predictors selected under the same modeling rules. Because empirical selection frequencies, bootstrap resampling, and cross-validation-based stability analyses were not performed, these findings should not be interpreted as formal evidence of algorithmic stability. Instead, the observed consistency reflects the use of a constrained modeling framework, harmonized predictor definitions, and pre-specified theory-informed covariates.

The results were also strikingly consistent across analytical platforms. When identical constraints and harmonized inputs were applied, the R and Python implementations produced comparable point estimates and retained the same predictors. This cross-platform agreement indicates that the framework is driven by methodology rather than software choice, and that substantive conclusions remain intact across computational environments, provided that complex survey design features are correctly specified.

Crucially, the constrained selection strategy ensured that key confounders particularly demographic and design-related variables were protected throughout the modeling process, regardless of statistical significance. This led to models that were both interpretable and theoretically grounded, reducing redundancy among correlated predictors without sacrificing conceptual coherence. In the final current cigarette-smoking model, a priori demographic variables were retained by design alongside selected modifiable predictors. Model diagnostics were reported to improve transparency, but were interpreted descriptively because the objective was explanatory modeling rather than individual-level risk prediction. For the current cigarette-smoking model, discrimination was good, with a C-statistic of 0.856 (95% CI: 0.841 to 0.870). The Brier score was 0.094, and the Hosmer-Lemeshow test showed no strong evidence of poor calibration (P = 0.104). For the secondary current e-cigarette-use model, discrimination was also good, with a C-statistic of 0.847 (95% CI: 0.834 to 0.859), and the Brier score was 0.050. However, the Hosmer-Lemeshow test indicated calibration misfit for the e-cigarette model, χ²(8) = 117.57, P < 0.0001. These diagnostics support transparent model reporting but should not be interpreted as evidence that either model is intended for individual-level prediction.

Together, these findings highlight the value of disciplined, theory-informed constraint when modeling adolescent smoking behavior using complex survey data and point to the broader applicability of this approach for population-based health research. As a sensitivity check, selected predictors were compared under AIC- and BIC-based selection. Predictors retained under both criteria were interpreted as more robustly supported, while predictors selected only under AIC were interpreted as more model-dependent (see Supplementary Tables S3, S4).

Discussion

Methodological innovations

This study contributes a practical and generalizable solution to a long-standing methodological challenge in complex survey analysis: how to reduce high-dimensional predictor sets without sacrificing theoretical coherence, confounding control, or design-based validity. This framework introduces several important methodological advances for analysing GYTS data in a way that is both rigorous and practical. First, the use of constrained stepwise selection allows key demographic variables to be retained in the model while still permitting data driven selection of other predictors. This approach strikes a balance between theory and parsimony, and directly addresses well known weaknesses of fully unconstrained stepwise methods that may drop variables of clear public health importance [54]. By protecting a priori factors such as age, sex, and grade, the resulting models remain interpretable and policy relevant [10]. Although modern regularization approaches such as LASSO and Elastic Net can improve prediction accuracy, they do not naturally preserve theoretically important confounders and currently have limited support within design-based survey analysis frameworks. Future work should explore survey-weighted penalized regression as an extension of the proposed methodology.

A second innovation relates to cross platform reproducibility [55]. Providing parallel implementations in both R and Python enables researchers to replicate findings across analytical environments and select software based on institutional capacity rather than technical constraints. This is particularly important in global health settings where access to software, training, and long-term support varies widely. The availability of parallel code also strengthens transparency and reproducibility, which are increasingly expected in multi country epidemiological research.

The framework also introduces an educational harmonization protocol that supports valid comparisons across countries with different schooling systems [50]. Mapping grades to comparable educational stages reduces misclassification and improves the interpretability of cross-country estimates. This step is often overlooked or inconsistently applied in previous multi country GYTS analyses, yet it is critical when comparing adolescent behaviors across diverse educational contexts.

Software comparison and recommendations

Table 13 shows the comparison between R and Python highlights clear strengths and trade offs that are relevant for GYTS analysis [56]. R offers strong native support for complex survey designs, including full specification of weights, strata, and clusters, making it the preferred option for producing official prevalence estimates and survey weighted descriptives [57]. It also supports constrained stepwise selection through established tools, which aligns well with the proposed modeling framework (see Table 13). Python, on the other hand, excels in machine learning integration and advanced predictive modeling, supported by a rapidly growing data science ecosystem [58]. While its support for complex survey designs remains limited, it is well suited for extending analyses into machine learning applications once core survey-based estimates have been established. Overall, R is recommended for primary GYTS analysis and official reporting, while Python is best positioned as a complementary tool for exploratory modeling and advanced analytics, depending on the research objective and available expertise.

Table 13. Platform comparison for Global Youth Tobacco Survey analysis.

✓ indicates full, native, or strong support; χ indicates limited support or a feature that requires a manual/custom implementation; and Δ indicates growing or emerging support.

Predictor variable selection, lessons for future studies

Several clear lessons emerge from the way predictors behaved across countries and platforms. Age should always be treated as a categorical variable rather than a continuous one, as risk does not increase in a straight line and tends to peak in early to mid-adolescence, particularly among those aged 13 to 14 years [34]. Core demographic factors such as age, sex, grade, and region should be protected as a priori variables to reduce confounding and improve comparability across settings. Awareness related variables, especially awareness of tobacco products, consistently show the strongest associations and should be central to any GYTS smoking model. Social context also matters deeply, with peer and parental smoking repeatedly emerging as dominant predictors, underscoring the need to capture the social environment in a comprehensive way. Policy related variables should be retained even when they are not statistically significant, as they reflect modifiable levers for public health action [46]. In contrast, product specific variables such as brand or flavor preferences are best excluded from initiation models, as they introduce selection bias and blur causal interpretation.

Strengths and limitations

This study has several notable strengths. First, it addresses a major methodological gap in the GYTS literature by providing a transparent and reproducible framework for analysing complex survey data. Second, predictor selection was guided by Social Cognitive Theory and empirical evidence rather than relying solely on automated statistical procedures. Third, the framework explicitly incorporates survey weights, stratification, and clustering, thereby improving the validity of population-level inference. Fourth, all analytical workflows are implemented using open-source software and accompanied by reproducible code. Finally, the framework was evaluated using both single-country and pooled multi-country datasets and is readily adaptable to other large-scale health surveys that employ multi-stage sampling designs.

This study has several important limitations that should be considered when interpreting the findings. In Python, survey weighted analysis relies on frequency weights and does not fully account for stratification and clustering, which may lead to small differences in standard errors compared with full design-based approaches [59]. Stepwise selection, while practical and widely used in epidemiological modeling, can produce sample-dependent variable selection and may increase the risk of overfitting. We mitigated these concerns through theory-based predictor specification, protection of a priori confounders, AIC-based selection criteria, comparison with BIC-selected predictors where available, and transparent reporting of all candidate variables. The cross-sectional nature of GYTS data also limits causal inference, with temporal ordering inferred based on theory rather than direct observation. The reliance on complete-case analysis is an important limitation. Although some individual variables had low item-level missingness, the cumulative effect of requiring complete information across all variables led to exclusion of a substantial number of participants, particularly in the Zambia analysis. If missingness was related to smoking status, social exposures, or school attendance patterns, complete-case estimates may be biased. Future applications of the framework should consider multiple imputation or inverse-probability weighting as sensitivity analyses when missingness is substantial. Finally, tobacco use is self-reported and may be underestimated due to social desirability, although the anonymous design of GYTS helps reduce this bias [60].

Implications for tobacco control

Despite these limitations, the findings have clear implications for tobacco control policy and programming. The consistent importance of awareness, peer influence, and marketing exposure points to priority intervention areas. Strengthening counter marketing campaigns can help reduce the appeal of tobacco products among adolescents [61]. Enforcing bans on smoking in public places and school environments remains essential, both to reduce exposure and to shift social norms. School based education programs continue to show protective effects and should be expanded and better integrated into curricula [62]. Interventions that directly address peer norms, including social network-based approaches, may be especially effective during adolescence when social influence is strongest.

Replicability to other complex surveys

The proposed framework is highly transferable to other large-scale surveillance systems that use stratified multi stage cluster sampling with survey weights. Surveys such as GATS [63], WHO STEPS [64], and GSHS [65] share core design features with GYTS, making direct adaptation feasible. With appropriate mapping of strata, primary sampling units, and weights, the same principles of harmonization, constrained variable selection, and design-based analysis can be applied. While R remains the strongest option for full survey specification, Python can complement these analyses where advanced modeling is required. Across all platforms and surveys, careful attention to design effects is essential, as ignoring clustering and stratification leads to underestimated uncertainty and misleading confidence intervals.

Conclusion

This methodological framework provides researchers with adaptable, reproducible workflows for GYTS smoking behavior analysis that address critical gaps in current practice. The R implementation serves as the reference standard, offering robust, well-validated native support for complex, two-stage survey designs and official variance estimation consistent with GYTS methodological requirements. The Python implementation is provided to enhance reproducibility, transparency, and extensibility, particularly for custom algorithms and machine-learning–based explorations, but is not intended as the primary platform for variance estimation. Both platforms produce highly comparable point estimates when properly configured, although R remains the preferred platform for valid variance estimation because it fully supports complex survey designs. The comprehensive predictor variable justification framework, grounded in Social Cognitive Theory and supported by systematic review evidence, ensures that analyses are theoretically coherent, empirically supported, and policy relevant. By standardizing variable selection and harmonization protocols, we enable valid cross-country comparisons and meta-analytic synthesis, strengthening the evidence base for adolescent tobacco control globally.

Generative AI statement

The author(s) declared that Generative AI was not used in the creation of this manuscript.

Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.

Publisher's note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

The GitHub repository is available at: https://github.com/andreasngambi-star/GYTS_analysis_R_Python.git, and the persistent Zenodo DOI is 10.5281/zenodo.21254547.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fams.2026.1895038/full#supplementary-material

Supplementary material Supplementary material at the source

Supplementary material Supplementary material at the source

Footnotes

Edited by: Ashis SenGupta, Augusta University, United States

Reviewed by: Shiva Shankar Bugude, Qassim University, Saudi Arabia

Abdulkadir Keskin, Istanbul Medeniyet University, Türkiye

Footnotes

AIC, Akaike Information Criterion; aOR, Adjusted Odds Ratio; CDC, Centers for Disease Control and Prevention; CI, Confidence Interval; CVD, Cardiovascular Disease; FCTC, Framework Convention on Tobacco Control; GATS, Global Adult Tobacco Survey; GLM, Generalized Linear Model; GSHS, Global School-based Student Health Survey; GYTS, Global Youth Tobacco Survey; ICC, Intra-class Correlation Coefficient; LMICs, Low- and Middle-Income Countries; MPOWER, Monitor, Protect, Offer, Warn, Enforce, Raise (WHO tobacco control framework); NCD, Noncommunicable Disease; OR, Odds Ratio; PSU, Primary Sampling Unit; R, R Statistical Software; SADC, Southern African Development Community; SES, Socioeconomic Status; SPSS, Statistical Package for the Social Sciences; SAS, Statistical Analysis System; STEPS, STEPwise Approach to NCD Risk Factor Surveillance; UK, United Kingdom; WHO, World Health Organization.