3 — Methodology
Thesis · APA section
Thesis · Section 3 of 7 · Paper landing
3.1 Research design
Longitudinal ecological panel design; unit of analysis is the district-month (77 districts × 108 months, January 2015 – December 2023; N = 7,327 valid observations with non-missing, non-zero case data) (Bao et al., 2023; Hashizume et al., 2010). The analytical pipeline follows six sequential stages: (1) data acquisition and integration; (2) harmonization, quality assurance, and cleaning; (3) feature engineering; (4) exploratory analysis; (5) ML model development, training, and hyperparameter optimization; (6) model evaluation, feature-importance analysis, and result interpretation.
3.2 Study area
3.2.1 Geographic and administrative context
Nepal covers ~147,181 km² across 77 administrative districts spanning seven provinces. Districts range from ~241 km² (Bhaktapur) to ~3,422 km² (Dolpa), and from ~83,000 (Manang) to ~1,744,000 (Kathmandu) inhabitants (Central Bureau of Statistics, 2021).
3.2.2 Ecological zones and climatic characteristics
Three ecological zones span Nepal’s altitudinal gradient (60 m – 8,848 m asl):
- Terai (60–300 m asl; ~17 % land, > 50 % population): humid subtropical; 25–30 °C; 1,400–2,500 mm annual precipitation.
- Middle Hills (300–3,000 m asl; ~68 % land): temperate to subtropical; Kathmandu Valley ~18 °C, 1,400 mm.
- High Mountain (> 3,000 m asl; ~15 % land): alpine to arctic; sparse populations.
The southwest monsoon (June–September) delivers ~75–80 % of annual precipitation.
3.2.3 Health infrastructure and WASH coverage
Three-tier system: district hospitals, primary health care centers, village-level health posts — all reporting through HMIS. Provincial WASH coverage ranges from 78 % (Bagmati) to 43 % (Karnali); improved sanitation access remains < 75 % in several Terai districts (Central Bureau of Statistics and UNICEF, 2019).
3.3 Data sources
Four primary institutional and open-access sources documented in detail in §6 Data:
| Dataset | Variable | Resolution | Provider | Reference |
|---|---|---|---|---|
| HMIS | Monthly outpatient typhoid case counts | District-month | MoHP / DoHS | (Ministry of Health and Population, Government of Nepal, 2023) |
| DRR Portal | Flood event counts | District-month | NDRRMA | (NPDRR Nepal, 2025) |
| CHIRPS v2.0 | Precipitation | 0.05° (~5 km) | Climate Hazards Group | (Funk et al., 2015) |
| ERA5-Land | 2 m air temperature, relative humidity | 0.1° (~9 km) | ECMWF · Copernicus | (Muñoz-Sabater et al., 2021) |
| District shapefiles | Administrative boundaries | District polygons | CBS (Central Bureau of Statistics, 2021) | (Central Bureau of Statistics, 2021) |
3.4 Data cleaning and integration
- Calendar conversion — Bikram Sambat → Gregorian via 12-month mapping table (max alignment error ~15 days per record).
- District-code reconciliation — CBS 3-digit codes as primary spatial merge key; 12 name-spelling variants resolved.
- Numeric cleaning — comma-formatted HMIS counts parsed (847 records corrected); extreme values audited.
- Spatial aggregation — CHIRPS and ERA5-Land grids aggregated to district polygons by area-weighted zonal statistics.
- Missing data — 2,463 district-months excluded; no outcome imputation. Final analytical sample: N = 7,327.
3.5 Feature engineering
Eleven features engineered into five functional categories:
| Feature | Category | Derivation |
|---|---|---|
log_cases |
Outcome transform | log(cases + 1) |
cases_lag1 |
Autoregressive | cases value from prior calendar month, same district |
precip_lag1 |
Lagged climate | CHIRPS precipitation shifted 1 month |
temp_mean_lag1 |
Lagged climate | ERA5-Land temperature shifted 1 month |
humidity_lag1 |
Lagged climate | ERA5-Land humidity shifted 1 month |
flood_lag1 |
Lagged climate | DRR flood count shifted 1 month |
precip_roll3 |
Cumulative climate | rolling mean over t, t−1, t−2 |
temp_roll3 |
Cumulative climate | rolling mean over t, t−1, t−2 |
monsoon |
Seasonal indicator | 1 if month ∈ {Jun, Jul, Aug, Sep} |
month_sin, month_cos |
Cyclical encoding | sin/cos(2π × month / 12) |
All lag and rolling features computed within each district to prevent cross-district leakage.
3.6 Train-test split and validation design
Strictly chronological split — last 12 calendar months (January 2022 – December 2023) reserved as held-out test set; preceding ~84 months for training. Within training, TimeSeriesSplit (k = 5) for XGBoost hyperparameter selection.
3.7 Machine-learning model specifications
3.7.1 Random Forest
scikit-learn 1.3 (Pedregosa et al., 2011); final hyperparameters from 5-fold CV grid search:
n_estimators = 500;max_depth = None;min_samples_split = 5;min_samples_leaf = 2;max_features = 'sqrt'; seed = 42.
3.7.2 XGBoost
xgboost 1.7 (Chen & Guestrin, 2016); Bayesian optimization with 50 evaluation trials:
n_estimators = 400;max_depth = 6;learning_rate = 0.05;subsample = 0.80;colsample_bytree = 0.80;min_child_weight = 3;reg_alpha = 0.10;reg_lambda = 1.00;gamma = 0.10; seed = 42.
3.7.3 Multilayer Perceptron
scikit-learn MLPRegressor (Pedregosa et al., 2011):
- Hidden layers
(128, 64, 32); activation ReLU; α = 0.001 L2 regularization. - Adam optimizer @
learning_rate_init = 0.001; max 500 epochs; early stopping triggered by validation-loss change < 1 × 10⁻⁴ for 20 consecutive iterations. - Features standardized using
StandardScalerfitted only on training set.
3.7.4 Weighted Ensemble
Combines out-of-sample predictions of three base models using non-negative weights optimized to minimize RMSE on a held-out validation subset of training data, via scipy.optimize.minimize with SLSQP subject to weights ≥ 0 and Σweights = 1.
Optimized weights: XGBoost = 0.45, MLP = 0.35, RF = 0.20.
\[\hat{y}_{\text{ensemble}} = 0.20 \times \hat{y}_{\text{RF}} + 0.45 \times \hat{y}_{\text{XGB}} + 0.35 \times \hat{y}_{\text{MLP}}\]
3.8 Model evaluation metrics
- RMSE = \(\sqrt{\sum (y_i - \hat{y}_i)^2 / n}\) — cases per district per month.
- MAE = \(\sum |y_i - \hat{y}_i| / n\) — scale-interpretable.
- R² = \(1 - SS_{\text{res}} / SS_{\text{tot}}\) — variance explained.
- MAPE = \(100\% \times \sum |y_i - \hat{y}_i| / y_i / n\) — percentage error.
RMSE and MAE computed on back-transformed case-count predictions; R² computed on log-scale predictions.
3.9 Hyperparameter search summary
| Model | Hyperparameter | Search space | Optimized | Method |
|---|---|---|---|---|
| RF | n_estimators | 100, 200, 300, 500 | 500 | Grid |
| RF | max_depth | 5, 10, 20, None | None | Grid |
| XGB | n_estimators | 100–600 | 400 | Bayesian |
| XGB | max_depth | 3–10 | 6 | Bayesian |
| XGB | learning_rate | 0.01–0.30 | 0.05 | Bayesian |
| XGB | reg_alpha (L1) | 0.00–1.00 | 0.10 | Bayesian |
| XGB | reg_lambda (L2) | 0.50–2.00 | 1.00 | Bayesian |
| MLP | hidden_layer_sizes | multiple configs | (128, 64, 32) | Manual |
| MLP | alpha (L2) | 0.0001–0.01 | 0.001 | Grid |
| Ensemble | RF / XGB / MLP weights | 0–1, Σ=1 | 0.20 / 0.45 / 0.35 | SLSQP |
3.10 Ethical considerations and data limitations
The study used only aggregate, de-identified, district-level surveillance and climate data; no individual patient records were accessed. All HMIS data obtained through official institutional channels with appropriate data-access authorization. Five categories of limitation acknowledged:
- HMIS under-reporting likely most severe in Karnali and Sudurpashchim.
- Diagnostic misclassification without lab confirmation (Parry et al., 2002).
- Spatial aggregation averages over intra-district heterogeneity.
- Flood-event underregistration for small-scale events.
- Ecological fallacy — findings describe population-level associations.
← Back: 2 — Literature Review · → Next: 4 — Results