# High-resolution Population Maps Derived from Sentinel-1 and Sentinel-2

Nando Metzger  Rodrigo Caye Daudt  Devis Tuia  Konrad Schindler

###### Abstract

Detailed population maps play an important role in diverse fields ranging from
humanitarian action to urban planning. Generating such maps in a timely and
scalable manner presents a challenge, especially in data-scarce regions. To
address it we have developed Popcorn, a population mapping method whose only
inputs are free, globally available satellite images from Sentinel-1 and
Sentinel-2; and a small number of aggregate population counts over coarse
census districts for calibration. Despite the minimal data requirements our
approach surpasses the mapping accuracy of existing schemes, including several
that rely on building footprints derived from high-resolution imagery. E.g.,
we were able to produce population maps for Rwanda with 100 m GSD based on
less than 400 regional census counts. In Kigali, those maps reach an $R^{2}$
score of 66$\%$ w.r.t. a ground truth reference map, with an average error of
only $\pm$10 inhabitants/ha. Conveniently, Popcorn retrieves explicit maps of
built-up areas and of local building occupancy rates, making the mapping
process interpretable and offering additional insights, for instance about the
distribution of built-up, but unpopulated areas, e.g., industrial warehouses.
Moreover, we find that, once trained, the model can be applied repeatedly to
track population changes; and that it can be transferred to geographically
similar regions, e.g., from Uganda to Rwanda). With our work we aim to
democratize access to up-to-date and high-resolution population maps,
recognizing that some regions faced with particularly strong population
dynamics may lack the resources for costly micro-census campaigns.

###### keywords:

Population mapping , Deep learning , Weakly Supervised learning , Sentinel-1 ,
Sentinel-2

††journal: XXX (under review)

[label1]organization=Photogrammetry and Remote Sensing, ETH Zurich,
city=Zurich, country=Switzerland

[label2]organization=Environmental Computational Science and Earth Observation
Laboratory, EPFL,city=Sion, country=Switzerland

[label3]organization=Corresponding Author

##  1 Introduction

High-resolution population maps are an invaluable resource for multiple
domains, ranging from humanitarian aid to spatial analysis and planning. Fine-
grained information about the distribution of the resident population is
needed, for instance, to organize the distribution of food and medical
supplies, to optimize disaster response, to monitor migration patterns, and as
a basis for planning urban development. Arguably, there is a particularly
pressing need for such maps in developing countries that are often confronted
with rapid urbanization and internal migration, and associated challenges
regarding infrastructure, social development, and pressure on the environment.

Present methods for retrieving fine-grained, spatially explicit population
estimates often depend on extensive stacks of off-the-shelf geodata (Stevens
et al., 2015; Tu et al., 2022; Metzger et al., 2022). These datasets, in turn,
often depend on expensive base data such as very high-resolution (VHR)
satellite imagery. Moreover, they quickly become outdated in regions with
strong population dynamics and require frequent updating.

Especially high-resolution building polygons, available from organizations
like Google (Sirko et al., 2021, 2023), Ecopia (Dooley et al., 2020),
Microsoft (Microsoft, 2022), and DLR (Esch et al., 2022) offer undeniable
benefits for population mapping, but their sustained and comprehensive
maintenance presents a significant challenge. These datasets are typically
proprietary. It remains unclear how frequently they will release updates,
which regions will be covered or updated, and what their future cost may be.
There is a need for more flexible, scalable, and timely population mapping
approaches that do not depend on such high-quality, static, proprietary
datasets.

Bottom-up approaches use sparse micro-census counts and extrapolate from those
counts with the help of densely available auxiliary data to achieve the
required spatial coverage. Hillson et al. (2014) evaluate the uncertainty of
population estimates from satellite imagery and limited survey data, in a case
study for Bo City, Sierra Leone. Weber et al. (2018) combine micro-census
counts with high-resolution satellite imagery to create gridded population
estimates at 90 m GSD in northern Nigeria, as well as associated uncertainty
metrics via Monte Carlo simulation. Leasure et al. (2020) use a hierarchical
Bayesian framework to account for uncertainty in national population maps,
with a focus on sparse data. Similarly, Boo et al. (2022) apply micro-surveys
and Bayesian hierarchical models to yield detailed population estimates in the
Democratic Republic of Congo. A limitation of bottom-up approaches is the need
for a representative sample of micro-census counts, whose collection is a
major effort.

Top-down approaches redistribute coarse census counts to mapping units of
about 100$\times$100 to 1000$\times$1000 m2. This disaggregation process,
often termed _dasymetric mapping_ , can be understood as estimating _relative_
population counts: its goal is to determine what fraction of the overall
population is located within each pixel (respectively, mapping unit) in a
region. Those fractional weights are derived from various context variables
that, contrary to the population counts, are available at the target scale.
Current top-down models mainly use readily available geodata products,
including building polygons as mentioned above, OpenStreetMap (OSM-
Contributors, 2023), human settlement layers from land cover maps (Pesaresi
and Politis, 2022), settlement growth models (Nieves et al., 2020), and night
light composites (DAAC, 2018; NOAA’s National Centers for Environmental
Information, 2023). These covariates are mapped to disaggregation weights
either with tabular machine learning algorithms (e.g., random forests,
XGBoost), trained in a fully supervised fashion to reproduce more fine-grained
census units (Stevens et al., 2015; Tu et al., 2022; Sapena et al., 2022); or
with neural networks trained in a weakly supervised fashion. Specifically, our
previous work (Metzger et al., 2022) employs a guided super-resolution
technique to directly predict per-pixel counts that add up to the available
census units. One key takeaway from that work is that it is more effective to
predict local occupancy rates and then multiply them by the number of
buildings, rather than directly predict population counts.

When attempting to retrieve population maps from raw satellite images, rather
than from geodata layers with a more immediate relation to the number of
residents, the challenge becomes to extract features that are predictive of
population counts (or densities). The task is further complicated by the fact
that learning the feature extraction in a data-driven manner requires training
data, but ground truth population counts are only available either at few,
scattered locations or in aggregated form over large areas. Several studies
have explored ways to circumvent that bottleneck. Islam et al. (2017) proposed
to identify built-up regions with a Gaussian maximum-likelihood classifier
(a.k.a. quadratic discriminant analysis) and then use the resulting built-up
area maps to refine population density estimates. Along the same lines, Meta
and CIESIN (2022) detect built-up areas in VHR satellite images and generate a
binary built-up area layer with a 30$\times$30 m2 grid, which is subsequently
used to concentrate the per-region population in the built-up grid cells.
Similarly, Grippa et al. (2019) also utilize VHR imagery to generate a series
of land cover maps, which then serve as a basis for population disaggregation
in Dakar, Senegal. Also using high-resolution imagery, Jacobs et al. (2018)
combine Planet images (GSD 3 m) with fine-grained U.S. census blocks and train
a neural network to regress population density maps, using the aggregate
density per block as weak supervision signal. Hafner et al. (2023) aim to
estimate population growth rather than density, with a Siamese neural network
that takes as input a bi-temporal pair of Sentinel-2 images and regresses the
population change. A different, more interactive (and thus less scalable)
approach is described by Fibæk et al. (2022). Images from Sentinel-1 and
Sentinel-2 are segmented into four different classes of settlement structure
with a neural network. The class probabilities are then mapped to population
densities with the help of hand-crafted formulas. The neural network training
involves an active learning loop, during which the user must iteratively
select visually well-estimated training regions.

In summary, bottom-up methods in population mapping offer high accuracy but
lack scalability, making them unsuitable for large-scale applications.
Conversely, top-down models based on geodata are scalable, but limited by
sporadic updates of the geospatial data they rely on. Previous attempts to
incorporate satellite imagery into top-down models mitigate the dependence on
derived geodata products but do not scale well, since they require either VHR
imagery or manual interventions.

We argue that high-resolution population maps can be retrieved solely from
open-access Earth observation imagery. To that end, we introduce Popcorn
(“POPulation from COaRse census Numbers”), a model engineered to overcome the
above-mentioned limitations. Popcorn, and its ensemble variant Bag-of-Popcorn,
offer a solution for population mapping that

  1. 1.

is scalable: Our top-down approach can be trained with as few as 400
population counts over coarse census regions, making it suitable for country-
wide mapping. In addition, it generalizes beyond the training region rather
well.

  2. 2.

allows for frequent updates: The model is based exclusively on data from
Sentinel-1 and Sentinel-2, eliminating the reliance on commercial or one-off
data products. Besides the option of near real-time monitoring of population
dynamics, this ensures cost-effective and sustainable mapping.

  3. 3.

provides accurate estimates in data-scarce regions: Despite using satellite
imagery with moderate GSD and only weak supervision, the model reaches state-
of-the-art mapping accuracy in data-scarce environments.

We hope that the ability to produce population maps in a timely manner, from
publicly available data, and at low cost, will particularly benefit developing
countries with limited resources, as well as non-governmental organizations
operating in such countries.

We start by describing our data in Section 2. Section 3 presents the Popcorn
model in detail. Section 4 explains the experimental setup and Section 5
discusses the results. Sections 6 and 7 conclude the paper and give an outlook
on future research.

##  2 Data

###  2.1 Census counts and geographic boundaries

Our experiments use coarse-grained census data from Rwanda, Switzerland, and
Puerto Rico for training and fine-grained data for validation, summarized in
Table 1. The setup involves maintaining a consistent geographical area while
varying the level of detail (or granularity) of the census data. This approach
is a standard practice in the field of population disaggregation and super-
resolution studies. It allows us to effectively assess the performance of our
models in generating detailed population maps from coarser, less detailed
census data.

| #regions | avg. area | upscaling | difficulty
---|---|---|---|---
| $N$ | of regions | $S$ | $D=\frac{S}{N}$
Switzerland | 2318 | 18 km2 | $42^{2}$ | 0.8
Rwanda | 381 | 69 km2 | $83^{2}$ | 18
Puerto Rico | 945 | 14 km2 | $38^{2}$ | 1.5
Uganda | 1377 | 185 km2 | $136^{2}$ | 14

Table 1: Summary of training data and dataset difficulty. Upscaling from
census regions to a 1 ha raster is more difficult when the upscaling factor
$S$ is large, and when the available number of census regions $N$ is small.

#### Switzerland

Our training data originates from the most recent census data, broken down to
the level of Swiss municipalities. These are compatible with the 10-meter
pixel grid of Sentinel images, thanks to available high-resolution shape
files. While some previous studies (CIESIN, 2018; WorldPop, ; Schiavina et
al., 2019; Meta and CIESIN, 2022) have used older, but more numerous municipal
boundaries to maintain historical consistency, we opt for the most up-to-date,
post-merger boundaries. These offer a better match with the Sentinel imagery
and allow for population map aggregation at the images’ native resolution.

As reference data for evaluation, we leverage Switzerland’s register of
residence, a comprehensive national dataset with a grid resolution of 100
meters. This unique dataset contains 4.1 million individual population counts
and provides a robust benchmark to validate the satellite-based predictions of
our model. For completeness, we mention that the dataset does not show
population counts of one or two people: due to privacy concerns those are
rounded up to three. For our purposes this tiny bias is irrelevant.

#### Rwanda

For training, we utilize data based on Rwanda’s 2012 census and extrapolated
to 2020 with the help of a population growth model provided by the United
Nations, Department of Economic and Social Affairs, Population Division
(2022). Unfortunately, the boundaries for the 415 administrative regions
associated with the provided census are not accurate enough for our purposes.
We therefore match their population counts with high-resolution boundary
polygons obtained from OCHA ROSEA (2022), resulting in a new dataset that can
be meaningfully combined with Sentinel images, but has only 381 administrative
regions.

For evaluation we utilize detailed, publicly accessible census data for the
city of Kigali, made available by Hafner et al. (2023). That dataset offers
population counts on a 100 m grid, with a total of 72,000 grid cells to cover
Kigali’s 1,132,000 residents. It provides a high-resolution reference to
assess the performance of our model when faced with the conditions of a
developing country in the Global South. We note that the population counts in
the dataset are not integer numbers, and include a few implausible counts
below zero. We conjecture that this is probably the consequence of some
interpolation or re-gridding procedure during dataset creation. We have opted
to leave those values unchanged, avoid inadvertently introducing further
biases, and ensure our results remain comparable.

#### Puerto Rico

In Puerto Rico, our training data comes from the second-level administrative
tracts of the 2020 U.S. Census (U.S. Census Bureau, 2020) and contains 945
individual regions. The finest census level (so-called _blocks_) serves as a
reference for evaluation, with $\approx$41,000 blocks that have an average
area of about 400$\times$400 m2. While not as fine-grained as in Switzerland
and Kigali, the granularity is still sufficient for many purposes and can
serve as a reference to evaluate our image-based maps (after aggregating them
to census blocks).

####  2.1.1 Uganda

To evaluate the model’s ability to generalize across national borders, we
employ extrapolated census data from Uganda spanning from 2014 to 2020 from
United Nations, Department of Economic and Social Affairs, Population Division
(2022) together with the administrative boundaries provided by WorldPop . This
dataset comprises 1,377 samples and is rated at difficulty level 14, based on
the calculation scheme delineated in Table 1. For testing, we utilize the
fine-grained population counts from Kigali, as previously described.

###  2.2 Input Data

Figure 1: Schematic overview of our approach to population mapping from
Sentinel-1 and Sentinel-2 imagery. A pre-trained dual-stream (DS) building
detector estimates a per-pixel built-up score. Concurrently, a second,
trainable dual-stream block estimates occupancy rates. The population map is
derived as the per-pixel product of built-up score and occupancy. To supervise
the training of the occupancy branch, the predicted population counts are
aggregated within administrative regions and compared to the corresponding
census data.

#### Sentinel-1

Sentinel-1 is an active synthetic aperture radar (SAR) system deployed by the
European Space Agency (ESA). The mission consists of two satellites,
Sentinel-1A and Sentinel-1B, in phase-shifted sun-synchronous
orbits111Sentinel-1B failed in December 2021. The launch of its replacement
Sentinel-1C is planned for early 2024.. The constellation provides a revisit
time of six days, ensuring regular, frequent coverage for monitoring purposes.
For our purposes, we source data collected in interferometric wide swath mode
from Google Earth Engine. The data come in the form of ground-range corrected
back-scatter (log-)amplitudes for the VV and VH polarisations, resampled to 10
m ground sampling distance. We average images from the entire year 2020 into
seasonally averaged mosaics for the time frames as defined in Appendix A.

#### Sentinel-2

Sentinel-2 is an optical satellite mission developed by ESA. It also consists
of phase-shifted twin satellites, Sentinel-2A and Sentinel-2B, that operate in
sun-synchronous orbits and offer a revisit time of five days for regular,
consistent Earth observation and monitoring. Sentinel-2 captures high-
resolution multi-spectral imagery in 13 bands across the visible, near-
infrared, and shortwave infrared spectrum. The bands are particularly suited
for observing vegetation, land cover, and water bodies. For our purposes, we
source atmospherically corrected surface reflectance images (Level-2A product)
for the entire twelve months between March 2020 and February 2021 from Google
Earth Engine and produce four seasonal composites to minimize cloud coverage.
As commonly done for urban mapping (e.g., Pesaresi et al., 2016; Vigneshwaran
and Vasantha Kumar, 2018) we only use the four highest resolution (10 m GSD)
bands in red (B4), green (B3), blue (B2), and near-infrared (B8).

##  3 Methodology

The core of our method is a neural network model, termed Popcorn. That model
has two components: (1) a pre-trained, frozen built-up area extractor; and (2)
a building occupancy module that we train through weak supervision with coarse
census counts, as illustrated in Figure 1 and discussed in the following
subsections.

The model operates at the full Sentinel-1/-2 resolution, i.e., its output has
a nominal spatial resolution of 10 m. However, for the final product and
evaluation, we aggregate the raw output to a 1 ha (100$\times$100 m) grid.
Population counts at 10 m resolution do not have practical advantages, are
conceptually questionable, and are impossible to verify, as the inhabitants of
a single dwelling unit would in most cases be spread out over multiple pixels.

###  3.1 Building Extractor

To detect built-up areas we adapt the method developed by Hafner et al.
(2022). In their work, they introduce a classification-based neural network
architecture that separately extracts features from Sentinel-1 and Sentinel-2
with two parallel U-Net-like streams (called dual-stream architecture), which
are then concatenated and mapped to a built-up score, cf. Figure 2. The model
is initially trained with high-resolution building labels from urban areas in
North America and Australia, and made globally applicable in a self-supervised
fashion, exploiting the fact that predictions based on Sentinel-1 or
Sentinel-2 alone should be consistent also in the absence of reference labels.

We make the following modifications to the original architecture: We reduce
the channel depth in the U-Net feature extractors to 8 and 16 channels (as
opposed to the original 64 and 128 channels), see Figure 2. This turns out to
hardly affect performance (see Appendix B for a quantitative comparison), but
shrinks the total parameter count from 1.8 million to a mere 30’000. This
brings a drastic reduction in memory consumption which enables training with
much larger image patches. Moreover, it makes it possible to reuse the same
architecture for the occupancy branch (see below), which needs to be fine-
tuned with a much more limited dataset.

As the final layer of our building extraction branch, we employ a 1$\times$1
convolution followed by the sigmoid (a.k.a. logistic) activation function.
This means that the branch outputs a fractional _BuiltUp score_ , i.e., a
(pseudo-)probability that a given pixel lies in the built-up area, which can
be interpreted as a proxy for the local building density.

The training of our modified version of the building extractor closely follows
Hafner et al. (2022). We use the same training sites and target labels and
also maintain the same split into training, validation, and testing sites. In
other words, our model differs from the original one exclusively in terms of
architectural modifications.

###  3.2 Occupancy and Population Estimation

In earlier work (Metzger et al., 2022) we showed that, in a scenario where
per-pixel building counts have been observed, it is advantageous to leave
those counts unchanged and to only estimate a map of the per-building
occupancy rate $\text{OccRate}_{i}$, where $i$ denotes the pixel index. In
much the same way, also for Popcorn, we factor the population number
${\hat{p}_{i}}$ into the built-up score and an occupancy rate,

${\hat{p}_{i}}=\text{BuiltUp}_{i}\times\text{OccRate}_{i}\;.$ (1)

Even though both variables are derived from the same input images, this
strategy outperforms a direct estimation of the population counts, as we show
in Section 5.5. The reason for this is presumably the additional information
contributed by the built-up area labels that are used to train the building
detector, but for which there are no corresponding population counts to train
the direct predictor.

To estimate building occupancy we instantiate another dual-stream branch with
the same U-Net architecture as the building detector. The prediction head for
occupancy, after concatenating the Sentinel-1 and Sentinel-2 features,
consists of three 1$\times$1 convolution layers, each having a hidden
dimensionality of 64, followed by a ReLU activation function. We empirically
found that initializing the U-Net with the weights of the pretrained building
extractor improves the learning in data-scarce regions, see Section 5.5.

###  3.3 Weakly Supervised Training

At the core of our learning approach is a loss function that provides weak
supervision from coarse census counts. The function aggregates the fine-
grained population estimates from the model into the administrative regions
for which the census counts are known. A similar training strategy has also
been used by Jacobs et al. (2018); Metzger et al. (2022); Hafner et al.
(2022).

The disagreement between the aggregated population estimates and the ground
truth counts serves as the loss to be minimized. Like Metzger et al. (2022) we
use the log-$L1$ distance as a measure of disagreement:

$\mathcal{L}=\sum_{j\in\mathcal{N}}\left(\log(1+c_{j})-\log\left(1+\sum_{k\in
A_{j}}\hat{p}_{k}\right)\right)\;,$ (2)

where $\mathcal{N}$ denotes the set of administrative regions under
consideration. For each region $j$ in $\mathcal{N}$, $c_{j}$ is the true
census count, and the term $\sum_{k\in A_{j}}\hat{p}_{k}$ sums up the model
estimates $\hat{p}_{k}$ over all grid cells $k$ that fall into region $A_{j}$.
To stabilize the loss, the constant offset of $1$ bounds the logarithmic
counts from below in regions with very few inhabitants (or none at all).

Figure 2: Dual-stream (DS) architecture proposed by Hafner et al. (2022).

###  3.4 Bagging

As a standard way of making stochastically learned predictors more reliable,
we employ a model ensemble (Dietterich, 2000; Lakshminarayanan et al., 2017),
in the following referred to as a ”Bag-of-Popcorn”. Five instances of the
Popcorn network are trained independently, each time changing the random seed
that determines the values of all randomly initialized network weights as well
as the composition of the batches during stochastic gradient descent.
Additionally, each of the five models is applied to each of the four seasonal
image composites. The 20 resulting predictions are averaged to obtain our
final estimate.

###  3.5 Dasymetric Mapping

When estimating high-resolution population maps in the setting where census
data for the target region are available, we employ the dasymetric mapping
technique (McCleary Jr, 1969; Eicher and Brewer, 2001). I.e., in each census
region $A_{j}$ the estimated population numbers $\hat{p}_{i}$ of the pixels
are understood as _relative_ values w.r.t. the region aggregate, and re-scaled
to adjusted values $\hat{p}^{adj}_{i}$, such that their sum matches the actual
census count $c_{j}$ of the region:

$\hat{p}^{adj}_{i}=\frac{\hat{p}_{i}}{\sum_{k\in A_{j}}\hat{p}_{k}}\times
c_{j}\;,$ (3)

where the first term $\frac{\hat{p}_{i}}{\Sigma_{k\in A_{j}}\hat{p}_{k}}$
corresponds to the fraction of the region’s total population that lives within
pixel $i$. The post-processing step ensures that the estimated numbers add up
exactly to the census counts $c_{j}$. Unless those counts are further from the
truth than the satellite-based model estimates, the procedure can thus be
expected to improve the population maps, mitigating both over- and under-
estimation errors by the model.

We point out that in practice there can be situations where, indeed, the
available census counts are not more reliable than the model estimates, e.g.,
when the census is outdated or relies on projections rather than actual
counting. In that case, it may be better to refrain from dasymetric
calibration.

###  3.6 Implementation Details

Our model is implemented using PyTorch (Paszke et al., 2019). We maintain the
same hyper-parameter settings across all datasets and experiments, except for
the strength $\lambda_{\text{wd}}$ of the weight decay regularizer. The latter
is derived from the dataset difficulty scores $D$ (c.f. Table 1) via the
heuristic relation $5\lambda_{\text{wd}}=D\times 10^{-6}$. A regularizer on
the model outputs is also included to encourage sparsity via an additional
loss term defined as 0.01$\times$ the mean predicted outputs, similarly to
Jacobs et al. (2018).

To minimize the loss function we utilize the Adam optimizer (Kingma and Ba,
2015) with default parameters (i.e., $\beta_{1}=0.9$ and $\beta_{2}=0.999$).
The batch size is set to 2, such that the model can be trained on a single GPU
with 24 GB of onboard memory (in our case an Nvidia GeForce RTX™ 3090 Ti). The
base learning rate is progressively decayed by a factor of 0.75 after every
fifth epoch.

Considering the small batch size we opt to freeze the pretrained batch
normalization layers to avoid discrepancies between the outlier-sensitive
batch statistics while training and testing. The bias term of the final layer
in the occupancy branch is initialized with the country-wise disaggregation
factor, i.e. the average occupancy rate of the respective dataset.

The following standard data augmentation techniques are applied to increase
the variability of the training data: random adjustments of brightness and
linear contrast, random rotations by 0, 90, 180, or 270 degrees, and random
mirror reflection (“flipping”). Additionally, we account for seasonal
variations by randomly selecting images from all four seasons during training.

##  4 Experimental Setup

###  4.1 Evaluation Metrics

To evaluate the performance of our method we compare the predicted population
maps to the corresponding, fine-grained reference data in terms of three error
metrics: the coefficient of determination (R2), the Mean Absolute Error (MAE)
and the Root Mean Squared Error (RMSE).

##### Coefficient of Determination (R2)

The R2 score measures the proportion of variance explained by our predictions,
calculated as::

$R^{2}=1-\frac{\sum_{i=1}^{n}(p_{i}-\hat{p}_{i})^{2}}{\sum_{i=1}^{n}(p_{i}-\bar{p})^{2}}\;,$
(4)

with $p_{i}$ and $\hat{p}_{i}$ the true and predicted population counts for
the $i^{th}$ administrative cell, and $\bar{p}$ the mean of the true counts.

##### Mean Absolute Error (MAE)

The average absolute deviation between predicted and actual population counts,
in inhabitants per high-resolution unit, quantifies the expected prediction
error at a target unit in a robust fashion:

$\text{MAE}=\frac{1}{n}\sum_{i=1}^{n}|p_{i}-\hat{p}_{i}|$ (5)

##### Root Mean Squared Error (RMSE)

The mean _squared_ deviation between predicted and true counts, normalized
back to an intuitive scale of inhabitants per unit by taking the square root.
Due to the amplification of large deviations, the metric is strongly impacted
by rare, extreme prediction errors:

$\text{RMSE}=\sqrt{\frac{1}{n}\sum_{i=1}^{n}(p_{i}-\hat{p}_{i})^{2}}$ (6)

###  4.2 Baselines

To put the quality of the estimated, fine-grained population maps in context,
we benchmark our method against a broad range of competitors. The methods we
compare to fall into two categories: self-implemented baselines and
established third-party models.

  1. 1.

Self-implemented Baselines:

     * Plain Building Disaggregation: This elementary baseline simply apportionates the population within a census region among its constituting 1 ha pixels according to their building count. For Switzerland, the buildings are extracted from the national mapping agency’s Topographic Landscape Model (TLM) (Swisstopo, 2023), for Puerto Rico and Rwanda we use Google Open Buildings (Sirko et al., 2021).

     * BuiltUp Disaggregation: This baseline is closer to our work in that it does not require external building counts, but instead redistributes the aggregate numbers per census region to pixels based on the built-up scores derived from satellite imagery. We use the BuiltUp scores obtained with the extractor of Section 3.1

     * Bag-of-Popcorn+count: This variant combines the proposed Popcorn model with external building counts. The image-based built-up scores are replaced by high-resolution building count datasets based on TLM, respectively Google Open Buildings, whereas the per-pixel occupancy rates are predicted from satellite images as in Popcorn.

  2. 2.

Third-Party Models: For Switzerland and Rwanda, we compare to several other
population mapping schemes. A fair comparison for Puerto Rico is unfortunately
not possible, because those schemes use the highest-resolution census blocks
to generate their maps (whereas we intentionally build on a coarser census
level, so that the blocks can serve as ground truth for evaluation).

     * WorldPop (Stevens et al., 2015): WorldPop employs dasymetric mapping with a random forest model based on various off-the-shelf geo-datasets. There are two products, a basic one that distributes the population across all pixels in a census region, and one that additionally uses a built-up area map at the target resolution and constrains the disaggregation to the built-up pixels.

     * GPWv4 (CIESIN, 2018): The Gridded Population of the World Version 4 dataset offers global population estimates derived from the 2010 census round, and extrapolated for the years 2000, 2005, 2010, 2015, and 2020 using UN World Population Prospects. The population, according to the extrapolated numbers, is uniformly dispersed across each census region and stored in the form of counts on a 30 arc-second grid ($\approx 1$ km at the equator), taking into account the latitude-dependent area of the grid cells.

     * GHS-Pop (Carneiro Freire et al., 2016; Schiavina et al., 2019): Fuses sub-national census data with the GHS-Built layer (Pesaresi and Ehrlich, 2009) and a gridded built-up area map derived from Landsat.

     * HRPDM (Meta and CIESIN, 2022): The High Resolution Population Density Maps were produced by extracting a binary built-up layer from high-resolution satellite imagery and redistributing the available census counts to only the built-up pixels.

     * Pomelo (Metzger et al., 2022): This method predicts building occupancy from a collection of geodata layers (nightlights, terrain height, distance to roads and waterways, etc.) and combines it with building counts derived from high-resolution footprints (Google Open Buildings, Maxar Ecopia Maps). Optionally, it offers dasymetric disaggregation of available census counts.

##  5 Results and Discussion

In the following, we examine the prediction quality of the proposed Bag-of-
Popcorn model and compare it with the alternative population mapping
techniques described above. The evaluation encompasses three very different
regions with diverse geographical conditions, chosen also to have sufficiently
high-resolution reference data: Switzerland, Kigali (Rwanda), and Puerto Rico.
Quantitative results are shown in Tables 2, 3, and 4, respectively, where bold
numbers indicate the best-performing map per category. Moreover, we provide
visual comparisons and an exemplary use case for monitoring population
dynamics across time.

###  5.1 Evaluation of Population Maps

Table 2: Quantitative evaluation for Switzerland.

|  | R2 $\uparrow$ [%] | MAE $\downarrow$ | RMSE $\downarrow$
---|---|---|---|---
High Resolution | TLM Disaggregation | 40 | 1.60 | 10.2
WorldPop-Builtup | 37 | 2.11 | 10.5
Meta | 41 | 1.77 | 10.2
Pomelo | 53 | 1.45 | 9.1
Bag-of-Popcorn+count | 39 | 1.45 | 10.4
Medium Resolution | BuiltUp Disaggregation | 53 | 1.68 | 9.1
GPWv4 | 6 | 4.3 | 12.9
GHS-Pop | 45 | 1.88 | 9.8
WorldPop | 38 | 2.4 | 10.4
Bag-of-Popcorn | 59.9 | 1.35 | 8.4

Table 3: Quantitative evaluation for Kigali, Rwanda.

|  | R2 $\uparrow$ [%] | MAE $\downarrow$ | RMSE $\downarrow$
---|---|---|---|---
High Resolution | Google Buildings Disag. | 43 | 10.5 | 26.2
WorldPop-Builtup | 58 | 10.2 | 22.6
Meta | 42 | 12.2 | 26.5
Pomelo | 69 | 9.8 | 19.5
Bag-of-Popcorn+count | 61 | 10.2 | 21.6
Medium Resolution | BuiltUp Disaggregation | 34 | 12.0 | 28.1
GPWv4 | 12 | 20.4 | 32.6
GHS-Pop | 20 | 14.6 | 31.1
WorldPop | 33 | 16.7 | 28.4
Bag-of-Popcorn | 66 | 10.1 | 20.2

Table 4: Quantitative evaluation for Puerto Rico.

|  | R2 $\uparrow$ [%] | MAE $\downarrow$ | RMSE $\downarrow$
---|---|---|---|---
High Res. | Google Buildings Disag. | 74.6 | 26.8 | 65
Bag-of-Popcorn+count | 77.6 | 24.3 | 61
Medium Res. | BuiltUp Disaggregation | 74.0 | 28.7 | 65
GPWv4 | -77.0 | 67.3 | 171
Bag-of-Popcorn | 81.8 | 23.6 | 55

Switzerland |  | Rwanda |  | Puerto Rico |
---|---|---|---|---|---
|  |  |  |  |

Figure 3: Scatter plots for Switzerland, Rwanda, and Puerto Rico. Note the
logarithmic scale of the axes. Values close to zero (below 0.5) have been
grouped into a single bin.

Tables 2, 3, and 4 summarize the evaluation metrics for the three tested
regions. We note in particular that the Bag-of-Popcorn model, based
exclusively on medium-resolution Sentinel imagery, is competitive even with
methods that have access to additional, higher-resolution data products such
as building footprints extracted from VHR images or topographic – in most
cases even overtaking them. On the Swiss dataset it achieves the best results
for all three error metrics, i.e., the highest correlation with the ground
truth counts as well as the smallest (absolute and squared) deviations from
the true counts. On the Rwandan dataset, Bag-of-Popcorn surpasses all other
medium-resolution methods by large margins, boosting the $R^{2}$ value from
34% (for the strongest competitor) to an impressive 66%. The good mapping
quality underlines the robustness of our approach under conditions that are
challenging for map down-scaling and for machine learning: to map Kigali at 1
ha resolution one must disaggregate the census counts by a factor >6800, while
only 381 regions are available to train the method. Even among methods that
make use of high-resolution inputs, only Pomelo narrowly outperforms our
proposed method, but to do so requires up-to-date building footprints. The
mapping performance in Puerto Rico had to be evaluated at the level of census
blocks that are on average $\approx$400$\times$400 m2 in size, as no ha-level
reference data is available. The corresponding results are consistent with
those for Switzerland and Rwanda. Bag-of-Popcorn has a clear edge over simpler
disaggregation schemes driven only by building counts. Notably, dropping the
retrieval of built-up scores from Sentinel and instead using building counts
derived from Google Open Buildings (Bag-of-Popcorn+count) does improve the
mapping performance compared to simple disaggregation, i.e., the estimated
occupancy values bring an added value. Interestingly, our default Bag-of-
Popcorn is still better than this variant, in other words, the counts derived
from high-resolution building outlines do not bring an advantage over built-up
scores estimated from Sentinel images. Generally speaking, we find that
methods that allow for varying building occupancy work best, while methods
that assume a constant occupancy and disaggregate building information are
less accurate, but still fairly robust. Finally, the GPWv4 maps at 1 km
resolution perform worst, with very low (in Puerto Rico even strongly
negative) $R^{2}$ values. This provides evidence that down-scaling to ha-level
resolutions is indeed meaningful, and that readily available covariates like
satellite images or building counts do contain the necessary information to
resolve population maps to that resolution.

In Figure 3, we show scatter plots between Bag-of-Popcorn predictions and true
population counts for the three different geographic regions. While,
naturally, there are some large relative discrepancies (especially for low
counts), the values are clustered near the identity line across the entire
value range, and the distribution of the prediction errors is symmetric. I.e.,
most predictions deviate only little from the true values, and the overall fit
appears largely unbiased.

|  VHR |  Sentinel-2 |  Bag-of-Popcorn |  WorldPop Builtup. |  WorldPop |  Pomelo |  GHS-Pop
---|---|---|---|---|---|---|---
Mahama Refugee Camp, Rwanda |  |  |  |  |  |  |
Kingi, Rwanda |  |  |  |  |  |  |
Zurich, Switzerland |  |  |  |  |  |  |
|

Figure 4: Comparison of population maps for the Mahama refugee camp (Rwanda),
Kingi (Rwanda), and Zurich, Switzerland. The first column shows very high-
resolution images (Google Maps©, 2023), while the second column shows the
Sentinel-2 RGB composites (European Space Agency, 2020) for reference. Bag-of-
Popcorn maps were resampled to the same grid as WorldPop to ease visual
comparison.

Figure 4 provides a visual comparison, contrasting population maps created by
the Bag-of-Popcorn model with those of WorldPop, Pomelo, and GHS-Pop. For the
examples from Rwanda, no high-resolution ground truth is available to verify
the estimates. Bag-of-Popcorn’s map prediction for the Mahama Refugee Camp
indicates that the model has identified the human-made structures of the camp.
At the same time, it assigns them comparably lower densities than Pomelo and
GHS-Pop. By contrast, both WorldPop maps completely missed the camp, possibly
owing to the use of outdated base data. Turning to Kingi in Rwanda, the Bag-
of-Popcorn estimates are very similar to those of its most accurate competitor
Pomelo. In contrast, WorldPop seems to fall short when it comes to detecting
thin rows of buildings along the roads, a rather frequent building pattern in
low-density, rural regions that have basic road infrastructure. In Zurich,
Switzerland, Bag-of-Popcorn and Pomelo again predict the most credible overall
density pattern, although the Pomelo map contains one implausible spike. Those
two methods also stand out as the only ones that correctly handle zero-density
regions within the urban fabric such as the green fields in the image center.
The WorldPop maps exhibit excessive smoothing, possibly due to inaccuracies of
the underlying built-up area layer.

###  5.2 Further Observations

##### Built-up scores vs. building footprints

Somewhat surprisingly, our results show that, when it comes to population
disaggregation, high-resolution building counts do not necessarily provide
superior guidance compared to low-resolution built-up area scores. We
hypothesize that the continuous, soft built-up scores contain implicit
information about (aggregate) building size. We note that the effect is more
apparent in developed regions (Switzerland and Puerto Rico), possibly because
the built-up area detector was originally trained with data from the US and
Australia and is less familiar with the building structures of Rwanda.

##### Resolution and Synchronization

Another unexpected result was that the variant of Popcorn that has access to
high-resolution building counts tends to perform worse than its purely
Sentinel-based counterpart. On the one hand, this could be one more instance
of the issue discussed in the previous paragraph. But other factors might also
be at play. The high-resolution building footprints, derived from VHR
satellite data or airborne imagery, are not temporally aligned with the
Sentinel data used to derive occupancy, so in regions with substantial
construction (or demolition) activity the observed state might differ, leading
to inconsistencies between the two branches. The Google Open Buildings, in
particular, do not specify the period for which a building footprint is valid.

###  5.3 Visual Examples

VHR
---
Sentinel-2
Prediction
Occupancy
Built-up score

Figure 5: Comparison of population density estimates in the Goma-Gisenyi
border region: The vertical strip is the boundary, the left part lies in the
DRC, and the right in Rwanda. Sources: VHR (Google Maps©, 2023), Sentinel-2
(European Space Agency, 2020).

Figure 5 displays a part of the Goma and Gisenyi, two cities at the border
between the DRC and Rwanda respectively. The example illustrates the two
streams of the Bag-of-Popcorn model. In the upper part of the left side lies
the very densely populated Mapendo neighborhood. The built-up score is
saturated, still, the model manages to recognize the substantially higher
occupancy, and consequently output a higher population density. The figure
displays also VHR imagery to visually confirm the prediction. Arguably, the
explicit maps of occupancy and built-up area score also afford our model a
degree of interpretability, in the sense that physically meaningful
intermediate quantities can be separately inspected and their distributions
checked for plausibility.

|
---|---
|

Figure 6: Population density maps for four distinct locations in Puerto Rico,
overlaid on desaturated VHR images Google Maps© (2023) (used for visualization
only, not involved in map generation). Not only vegetation but also industrial
warehouses are recognized as uninhabited, whereas isolated buildings in the
forest are picked up.

In a further example from Puerto Rico, Figure 6, we showcase the ability of
the model to distinguish between high-density (top-left) and low-density (top-
right) residential zones, while correctly identifying industrial zones as
built-up, nevertheless unpopulated (bottom-left). Thanks to the dedicated,
separately trained built-up area extractor, the model is capable of detecting
isolated buildings that are likely inhabited (bottom-right). Together, the
four examples attest to a rather nuanced representation of population patterns
across different settlement structures, especially considering the fact that
the estimates are derived only from Sentinel imagery at 10 m GSD.

###  5.4 Generalization

We have evaluated how well the Bag-of-Popcorn model generalizes by training it
on data from Uganda and then applying it to Kigali, Rwanda. We present the
resulting performance in Table 5. As expected there is a noticeable drop in
performance, with an $R^{2}$ score of 44% compared to 66% for the model
trained for Rwanda. Still, even the Bag-of-Popcorn model trained on a
different, neighboring country surpasses all medium-resolution baselines that
are specialized to Rwanda, c.f. Table 3. In absolute terms, there is only a
moderate increase in mapping error from 10.1 to 11.6 people/ha in MAE, which
for many applications may still be acceptable.

Table 5: Transferring a Bag-of-Popcorn model to a new geographic location. Two model instances trained on Rwanda, respectively Uganda, are evaluated on the same target region, Kigali (Rwanda). Training set | R2 $\uparrow$ [%] | MAE $\downarrow$ | RMSE $\downarrow$
---|---|---|---
Rwanda | 66 | 10.1 | 20.2
Uganda | 44 | 11.6 | 26.1
Figure 7: Estimated population time-series for Bunia, DRC, generated by
independent, repeated mapping with the Bag-of-Popcorn model. The figure
illustrates the model’s stability over time and its potential for monitoring
population dynamics. Notably, it captures the progressive formation of the
Kigonze refugee camp north of the city center, but also the gradual
densification of the urban fringe west of the center.

To showcase the potential of our model for monitoring population dynamics, we
applied the instance trained in Rwanda to a time series of images showing
Bunia, a rapidly growing city in the neighboring DRC. Figure 7 chronicles the
estimated evolution of the city’s population from 2019 to 2023. The satellite-
based analysis confirms the expected strong growth of Bunia’s population
(predicted counts are posted above the maps). Conspicuously, the model can
track the locally concentrated growth of a refugee camp in the upper half of
the depicted region, and diffuse densification, especially in the southwestern
periphery.

###  5.5 Ablation Studies

#### Model Architecture

We conducted an ablation to illustrate the influence of the separately trained
building extractor, and of using it as initialization also for the occupancy
branch. The experiment is illustrated in Figure 6, and with the results
subsequently below. The analysis confirms that using separate modules for the
built-up score and the occupancy rate, rather than a monolithic population
estimator, offers noticeable advantages in the low data regime. The difference
is most pronounced in Rwanda, which features the smallest number of census
regions and the largest upscaling factor. In contrast, for the larger datasets
from Switzerland and Puerto Rico, the performance of both variants is
comparable.

Moreover, the study highlights the importance of initializing the occupancy
branch with the weights of the pre-trained building extractor, benefiting from
the larger underlying training dataset. Again the difference is largest in
Rwanda, whereas in Switzerland and Puerto Rico, the number of census samples
is sufficient to learn a passable occupancy predictor from scratch. Note that
starting from pre-trained weights not only enables the model to capitalize on
features learned from a larger training set but also promotes coherence
between the built-up area and occupancy branches, which may simplify learning
based on small gradient steps.

Table 6: Ablations of the model architecture. Cases A and C employ the dual-
branch structure with separate branches for the builtUp score and the
occupancy rate, c.f. Figure 1. Cases A and B initialize the trainable weights
with those of the pretrained building detector.

|  |  | Switzerland | Rwanda | Puerto Rico
---|---|---|---|---|---
| Weight Init. | OccRate | R2 $\uparrow$ [%] | MAE $\downarrow$ | RMSE $\downarrow$ | R2 $\uparrow$ [%] | MAE $\downarrow$ | RMSE $\downarrow$ | R2 $\uparrow$ [%] | MAE $\downarrow$ | RMSE $\downarrow$
Case A | ✓ | ✓ | 59.9 | 1.35 | 8.4 | 66 | 10.1 | 20.2 | 81.8 | 23.6 | 54.8
Case B | ✓ | ✗ | 59.2 | 1.34 | 8.5 | 54 | 11.6 | 23.7 | 82.0 | 23.3 | 54.4
Case C | ✗ | ✓ | 55.0 | 1.47 | 8.9 | 46 | 11.3 | 25.6 | 77.6 | 26.6 | 60.7
Case D | ✗ | ✗ | 52.7 | 1.48 | 9.1 | 48 | 12.5 | 25.2 | 74.6 | 28.9 | 64.6

#### Bagging

We evaluated the impact of model ensembling, transitioning from a singular
Popcorn model to a Bag-of-Popcorn configuration, using the Kigali dataset, as
detailed in Table 7. Five independent instances of the Popcorn model were
trained, moreover, we have four seasonal composites of Sentinel-1 and
Sentinel-1 images for every location. The base case is a prediction with a
single model and a single seasonal composite (including Sentinel-1 and
Sentinel-2). We repeat that experiment with each of the five model instances
and each of the four composites and report average performance metrics. In a
second scenario, a single model is applied to each of the four seasonal
composites. The predictions are averaged, corresponding to an ensemble with a
single model instance, but four different inputs (”test-time augmentation”).
Again we repeat the experiment with all five model instances and report
average metrics. Note the implicit assumption that population changes between
seasons of the same year can be neglected. The third scenario is the standard
ensemble setting where the same input data is fed to five different model
instances. Again we run the experiment for each of the four seasons and report
average metrics. Finally, we regard each of the 20 possible combinations (5
instances $\times$ 4 seasons) as a member of one large ensemble and average
their results.

Our findings indicate that each ensembling method markedly enhances the
model’s predictive capability. Test-time augmentation using different seasonal
composites contributed more significantly to performance improvement than the
random ensemble of model instances. However, the most robust performance was
observed in the “Large Bag” scenario, which integrated both multiple models
and seasonal inputs.

Table 7: Effect of ensembling evaluated on the Rwanda Dataset, evaluated on
Kigali.

|  | Ensemble over |
---|---|---|---
|  | members | seasons | #estimates | R2 $\uparrow$ [%] $\downarrow$ | MAE $\downarrow$ | RMSE $\downarrow$
Popcorn |  |  |  | 1 | 55 | 11.2 | 23.2
Bag-of-Popcorn | Small |  | ✓ | 4 | 65 | 10.3 | 20.7
Medium | ✓ |  | 5 | 60 | 10.9 | 21.9
Large | ✓ | ✓ | 20 | 66 | 10.1 | 20.2

#### Input Modalities

Finally, to assess the contributions of the two modalities – Sentinel-1 SAR
amplitudes in two polarization, respectively Sentinel-2 RGB-NIR spectra – we
train versions of the model that utilize only one modality. From the dual-
stream architecture depicted in Figure 2, we selectively remove the backbone
corresponding to either Sentinel-1 or Sentinel-2 and pass only the features
from the other, active backbone to the prediction heads for built-up score and
occupancy. That procedure disentangles the relative contributions of the two
modalities and makes it possible to quantitatively assess how much each of
them contributes to the overall model output.

Results are displayed in Table 8. The experiment demonstrates that the two
input data sources contribute complementary information, and leveraging both
together gives the best results, across all three test sites. On the contrary,
there is no clear trend that one or the other is a more important data source
for population mapping. In Switzerland, Sentinel-1 alone gives better results
than Sentinel-2 alone, albeit by a moderate margin. In Rwanda and Puerto Rico,
Sentinel-2 alone nearly matches the performance of the combined model, whereas
Sentinel-1 alone leads to a considerable performance drop.

Table 8: Contribution of Sentinel-1 SAR and Sentinel-2 optical inputs.

|  | Switzerland | Rwanda | Puerto Rico
---|---|---|---|---
S1 | S2 |  $R^{2}\uparrow$ [%] | MAE $\downarrow$ | RMSE $\downarrow$ |  $R^{2}\uparrow$ [%] | MAE $\downarrow$ | RMSE $\downarrow$ |  $R^{2}\uparrow$ [%] | MAE $\downarrow$ | RMSE $\downarrow$
✓ | ✓ | 59.9 | 1.35 | 8.4 | 66 | 10.1 | 20.2 | 81.8 | 23.6 | 54.8
✓ | ✗ | 54.9 | 1.50 | 8.9 | 48 | 12.3 | 25.1 | 69.0 | 35.7 | 71.4
✗ | ✓ | 51.7 | 1.52 | 9.2 | 66 | 10.5 | 20.4 | 80.9 | 24.1 | 56.1

##  6 Conclusion

We have introduced Bag-of-Popcorn, a neural network model capable of
estimating gridded population counts at 1 ha resolution from Sentinel-1 and
Sentinel-2 satellite imagery. The proposed model has been shown to work well
in a range of geographic regions, and greatly simplifies fine-grained
population mapping: it requires neither a ground-based micro-census that is
difficult to conduct at scale in large parts of the world; nor high-resolution
geodata products that are expensive and/or not guaranteed to be updated. The
proposed model can be trained solely based on coarse census data and requires
only a small amount of training data to yield quite convincing results (in our
experiments in Rwanda, only 400 regional census counts).

We have experimentally demonstrated that the proposed Bag-of-Popcorn model
compares favorably to existing population mapping tools. In particular, it
excels in locations that lack high-resolution building information. However,
we found that even in locations for which high-quality building polygons are
available the model often outperforms methods that use those polygons.

The Popcorn framework separates population mapping into two separate streams
that estimate built-up areas and building occupancy rates. This strategy makes
it possible to exploit training data for built-up area mapping, which is much
more plentiful than population counts. By separately training a built-up area
detector on publicly available data, which also can serve as initial values
for building occupancy estimation, we are able to compensate for the
coarseness of census data in many parts of the world.

##  7 Outlook

The Bag-of-Popcorn model can retrieve high-resolution population maps only
from free satellite images, and optionally coarse census counts for dasymetric
rescaling. Still, several limitations remain that should be addressed in
future work. A main bottleneck is scalability. Although we used comparatively
small countries for our study, we had to limit model complexity to process
them on a GPU with 24 GB of onboard memory. To enable population mapping for
larger countries or even entire continents it will be necessary to develop
leaner, more memory-efficient variants of our model, or to port it to a high-
performance computing system (or a combination of both).

Going beyond one-shot population mapping, our preliminary experiments suggest
that Bag-of-Popcorn yields consistent results when applied to imagery of the
same location captured at different times. This could constitute a compelling
alternative to conventional census projection models – an intriguing prospect
that merits deeper investigation. If image-based tracking of population
dynamics turns out to be viable, it would not only provide a means to bridge
the periods between census rounds in a more qualified manner, but could also
help to better understand demographic trends, including population growth and
migration patterns.

Taking the idea of tracking population dynamics to the extreme, the short
revisit times of contemporary Earth observation satellites raise the question
whether it is possible to continuously track population numbers in near-
realtime. This could have transformative implications but requires further
research to ensure consistent and comparable estimates over short time scales,
despite inevitable fluctuations of the observed reflectance and backscatter
values.

From an engineering point of view, it may be interesting to explore the link
to other down-scaling (a.k.a ”guided super-resolution”) techniques, for
instance, guided filtering (He et al., 2012) or image diffusion (Metzger et
al., 2023). Moreover, we speculate that techniques for regression in
unbalanced datasets (Yang et al., 2021) could be useful to more precisely
estimate very high population densities.

Finally, we are convinced that timely, spatially explicit, high-resolution
population maps have a lot of untapped application potential. We hope that
easy-to-use models based on reliable and widely accessible data sources, in
the spirit of our Bag-of-Popcorn, will make population mapping accessible to a
wider community of interested users and service providers in sectors like
urban planning, disaster relief, and public health.

##  8 Declaration of Competing Interest

The authors declare no competing interests.

##  9 Author contributions

KS, DT, RCD, and NM together developed the concept and designed the
methodology. NM curated the data, implemented the methodology, and conducted
the experiments and analysis. KS, DT, and RCD supervised the study. All
authors contributed to writing the manuscript, based on an initial draft by
NM.

## References

  * Boo et al. (2022) Boo, G., Darin, E., Leasure, D.R., Dooley, C.A., Chamberlain, H.R., Lázár, A.N., Tschirhart, K., Sinai, C., Hoff, N.A., Fuller, T., et al., 2022.  High-resolution population estimation using household survey data and building footprints.  Nature Communications 13, 1330.
  * Carneiro Freire et al. (2016) Carneiro Freire, S., Macmanus, K., Pesaresi, M., Doxsey-Whitfield, E., Mills, J., 2016.  Development of new open and free multi-temporal global population grids at 250 m resolution, in: AGILE Conference on Geographic Information Science.
  * CIESIN (2018) CIESIN, 2018.  Gridded population of the world, version 4 (GPWv4): Population density, revision 11.  doi:https://doi.org/10.7927/H49C6VHW. accessed 9/10/2023.
  * DAAC (2018) DAAC, O., 2018.  MODIS and VIIRS land products global subsetting and visualization tool.  Subset of MOD13Q1 product at various sites .
  * Dietterich (2000) Dietterich, T.G., 2000.  Ensemble methods in machine learning, in: International Workshop on Multiple Classifier Systems, pp. 1–15.
  * Dooley et al. (2020) Dooley, C., Tatem, A., Bondarenko, M., et al., 2020.  Gridded maps of building patterns throughout sub-Saharan Africa, version 1.0.  University of Southampton, UK.
  * Eicher and Brewer (2001) Eicher, C.L., Brewer, C.A., 2001.  Dasymetric mapping and areal interpolation: Implementation and evaluation.  Cartography and Geographic Information Science 28, 125–138.
  * Esch et al. (2022) Esch, T., Brzoska, E., Dech, S., Leutner, B., Palacios-Lopez, D., Metz-Marconcini, A., Marconcini, M., Roth, A., Zeidler, J., 2022.  World settlement footprint 3d-a first three-dimensional survey of the global building stock.  Remote Sensing of Environment 270, 112877.
  * European Space Agency (2020) European Space Agency, 2020.  Sentinel-2 composite imagery.  URL: https://dataspace.copernicus.eu/. imagery acquired between 2020-03-01 and 2021-03-01.
  * Fibæk et al. (2022) Fibæk, C.S., Keßler, C., Arsanjani, J.J., Trillo, M.L., 2022.  A deep learning method for creating globally applicable population estimates from Sentinel data.  Transactions in GIS 26, 3147–3175.
  * Google Maps© (2023) Google Maps©, 2023.  Satellite imagery.  URL: https://www.google.com/maps/. accessed: 2023-10-05.
  * Grippa et al. (2019) Grippa, T., Linard, C., Lennert, M., Georganos, S., Mboga, N., Vanhuysse, S., Gadiaga, A., Wolff, E., 2019.  Improving urban population distribution models with very-high resolution satellite information.  Data 4, 13.
  * Hafner et al. (2022) Hafner, S., Ban, Y., Nascetti, A., 2022.  Unsupervised domain adaptation for global urban extraction using Sentinel-1 SAR and Sentinel-2 MSI data.  Remote Sensing of Environment 280, 113192.
  * Hafner et al. (2023) Hafner, S., Georganos, S., Mugiraneza, T., Ban, Y., 2023.  Mapping urban population growth from Sentinel-2 MSI and census data using deep learning: A case study in Kigali, Rwanda, in: Joint Urban Remote Sensing Event (JURSE).
  * He et al. (2012) He, K., Sun, J., Tang, X., 2012.  Guided image filtering.  IEEE Transactions on Pattern Analysis and Machine Intelligence 35, 1397–1409.
  * Hillson et al. (2014) Hillson, R., Alejandre, J.D., Jacobsen, K.H., Ansumana, R., Bockarie, A.S., Bangura, U., Lamin, J.M., Malanoski, A.P., Stenger, D.A., 2014.  Methods for determining the uncertainty of population estimates derived from satellite imagery and limited survey data: a case study of Bo City, Sierra Leone.  PloS one 9, e112241.
  * Islam et al. (2017) Islam, M.B., Becker, M., Bargiel, D., Ahmed, K.R., Duzak, P., Emana, N.G., 2017.  Sentinel-2 satellite imagery based population estimation strategies at FabSpace 2.0 Lab Darmstadt, in: CLEF Working Notes.
  * Jacobs et al. (2018) Jacobs, N., Kraft, A., Rafique, M.U., Sharma, R.D., 2018.  A weakly supervised approach for estimating spatial density functions from high-resolution satellite imagery, in: SIGSPATIAL International Conference on Advances in Geographic Information Systems, pp. 33–42.
  * Kingma and Ba (2015) Kingma, D., Ba, J., 2015.  Adam: A method for stochastic optimization, in: International Conference on Learning Representations (ICLR).
  * Lakshminarayanan et al. (2017) Lakshminarayanan, B., Pritzel, A., Blundell, C., 2017.  Simple and scalable predictive uncertainty estimation using deep ensembles.  Advances in Neural Information Processing Systems (NeurIPS) 30.
  * Leasure et al. (2020) Leasure, D.R., Jochem, W.C., Weber, E.M., Seaman, V., Tatem, A.J., 2020.  National population mapping from sparse survey data: A hierarchical Bayesian modeling framework to account for uncertainty.  Proceedings of the National Academy of Sciences 117, 24173–24179.
  * McCleary Jr (1969) McCleary Jr, G.F., 1969.  The dasymetric method in thematic cartography.  The University of Wisconsin-Madison.
  * Meta and CIESIN (2022) Meta and CIESIN, 2022.  High resolution population density maps + demographic estimates.  URL: https://registry.opendata.aws/dataforgood-fb-hrsl.
  * Metzger et al. (2023) Metzger, N., Daudt, R.C., Schindler, K., 2023.  Guided depth super-resolution by deep anisotropic diffusion, in: IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR), pp. 18237–18246.
  * Metzger et al. (2022) Metzger, N., Vargas-Muñoz, J.E., Daudt, R.C., Kellenberger, B., Whelan, T.T.T., Ofli, F., Imran, M., Schindler, K., Tuia, D., 2022.  Fine-grained population mapping from coarse census counts and open geodata.  Scientific Reports 12, 20085.
  * Microsoft (2022) Microsoft, 2022.  Worldwide building footprints derived from satellite imagery.  Dataset.  URL: https://github.com/microsoft/GlobalMLBuildingFootprints. gitHub Repository.
  * Nieves et al. (2020) Nieves, J.J., Sorichetta, A., Linard, C., Bondarenko, M., Steele, J.E., Stevens, F.R., Gaughan, A.E., Carioli, A., Clarke, D.J., Esch, T., et al., 2020.  Annually modelling built-settlements between remotely-sensed observations using relative changes in subnational populations and lights at night.  Computers, Environment and Urban Systems 80, 101444.
  * NOAA’s National Centers for Environmental Information (2023) NOAA’s National Centers for Environmental Information, 2023.  Defense meteorological satellite program – operational linescan system.  URL: https://ngdc.noaa.gov/eog/dmsp.html.
  * OCHA ROSEA (2022) OCHA ROSEA, 2022.  Rwanda - Subnational Administrative Boundaries.  https://data.humdata.org/dataset/cod-ab-rwa.  Contains shapefiles for Rwanda. Published by OCHA’s Regional Office for Southern and Eastern Africa.
  * OSM-Contributors (2023) OSM-Contributors, 2023.  Openstreetmap.  https://www.openstreetmap.org/.
  * Paszke et al. (2019) Paszke, A., Gross, S., Massa, F., Lerer, A., Bradbury, J., Chanan, G., Killeen, T., Lin, Z., Gimelshein, N., Antiga, L., et al., 2019.  Pytorch: An imperative style, high-performance deep learning library.  Advances in Neural Information Processing Systems (NeurIPS) 32.
  * Pesaresi et al. (2016) Pesaresi, M., Corbane, C., Julea, A., Florczyk, A.J., Syrris, V., Soille, P., 2016.  Assessment of the added-value of Sentinel-2 for detecting built-up areas.  Remote Sensing 8.
  * Pesaresi and Ehrlich (2009) Pesaresi, M., Ehrlich, D., 2009.  A methodology to quantify built-up structures from optical VHR imagery, in: Global Mapping of Human Settlement. CRC Press, pp. 55–86.
  * Pesaresi and Politis (2022) Pesaresi, M., Politis, P., 2022.  GHS-BUILT-S R2022A – GHS built-up surface grid, derived from Sentinel-2 composite and Landsat, multitemporal (1975–2030).
  * Sapena et al. (2022) Sapena, M., Kühnl, M., Wurm, M., Patino, J.E., Duque, J.C., Taubenböck, H., 2022.  Empiric recommendations for population disaggregation under different data scenarios.  PloS one 17, e0274504.
  * Schiavina et al. (2019) Schiavina, M., Freire, S., MacManus, K., 2019.  GHS population grid multitemporal (1975-1990-2000-2015), R2019A URL: http://data.europa.eu/89h/0c6b9751-a71f-4062-830b-43c9f432370f.
  * Sirko et al. (2023) Sirko, W., Brempong, E.A., Marcos, J.T., Annkah, A., Korme, A., Hassen, M.A., Sapkota, K., Shekel, T., Diack, A., Nevo, S., et al., 2023.  High-resolution building and road detection from sentinel-2.  arXiv preprint arXiv:2310.11622 .
  * Sirko et al. (2021) Sirko, W., Kashubin, S., Ritter, M., Annkah, A., Bouchareb, Y., Dauphin, Y., Keysers, D., Neumann, M., Cisse, M., Quinn, J., 2021.  Continental-scale building detection from high resolution satellite imagery.  arXiv preprint arXiv:2107.12283 .
  * Stevens et al. (2015) Stevens, F.R., Gaughan, A.E., Linard, C., Tatem, A.J., 2015.  Disaggregating census data for population mapping using random forests with remotely-sensed and ancillary data.  PloS one 10, e0107042.
  * Swisstopo (2023) Swisstopo, 2023.  TLM – topographic landscape model.  URL: https://www.swisstopo.admin.ch/en/home.html.
  * Tu et al. (2022) Tu, W., Liu, Z., Du, Y., Yi, J., Liang, F., Wang, N., Qian, J., Huang, S., Wang, H., 2022.  An ensemble method to generate high-resolution gridded population data for China from digital footprint and ancillary geospatial data.  International Journal of Applied Earth Observation and Geoinformation 107, 102709.
  * United Nations, Department of Economic and Social Affairs, Population Division (2022) United Nations, Department of Economic and Social Affairs, Population Division, 2022.  World population prospects 2022: Methodology of the United Nations population estimates and projections.  UN DESA/POP/2022/TR/NO. 4.
  * U.S. Census Bureau (2020) U.S. Census Bureau, 2020.  2020 census data for puerto rico.  URL: https://www.census.gov/geographies/island-areas/puerto-rico.html. accessed: 2023-07-10.
  * Vigneshwaran and Vasantha Kumar (2018) Vigneshwaran, S., Vasantha Kumar, S., 2018.  Extraction of built-up area using high resolution Sentinel-2A and Google satellite imagery.  The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences 42, 165–169.
  * Weber et al. (2018) Weber, E.M., Seaman, V.Y., Stewart, R.N., Bird, T.J., Tatem, A.J., McKee, J.J., Bhaduri, B.L., Moehl, J.J., Reith, A.E., 2018.  Census-independent population mapping in northern Nigeria.  Remote Sensing of Environment 204, 786–798.
  * (46) WorldPop, .  WorldPop: Open spatial demographic data and research.  https://www.worldpop.org/.  Accessed: 16-10-2023.
  * Yang et al. (2021) Yang, Y., Zha, K., Chen, Y., Wang, H., Katabi, D., 2021.  Delving into deep imbalanced regression, in: International Conference on Machine Learning (ICML), pp. 11842–11851.

## Appendix

###  A Definition of Seasons

For our work, we define the temporal windows for the seasonal composites as
listed in Table 9.

| Start | End
---|---|---
Spring | 01-03-2020 | 31-05-2020
Summer | 01-06-2020 | 31-08-2020
Autumn | 01-09-2020 | 30-11-2020
Winter | 01-12-2020 | 28-02-2021
Table 9: Start and end dates of seasons for image compositing.

###  B Building Extraction Performance

Table 10 presents a comparative analysis of the built-up area detection
results with latent representations of varying size (i.e., channel depth per
U-net layer). The _original_ configuration is the one recommended in Hafner et
al. (2022), with 64 channels at half-resolution and 128 channels at quarter-
resolution. The _slim_ configuration is our down-sized variant with 8,
respectively channels. The higher capacity of the original model yields
substantially better built-up area detection in the “Source Domain”, where the
model’s training labels are located, i.e., North America and Australia. But
that advantage does not persist in the “Target Domain” covered only by a
consistency loss between Sentinel-1 and Sentinel-2 prediction on unlabeled
data, i.e., on all other continents. There, our slim version of the model
performs as well. It appears that the extra capacity allows for a more
specific (over-)fit to the local patterns, but the domain adaptation scheme is
not able to carry over that advantage to unseen places with different
settlement structures.

Table 10: Performance comparison of original and slim model configurations in
built-up area detection. Metrics include F1 score and Intersection over Union
(IoU) for both source and target domains.

|  $\\#$params | Source Domain | Target Domain
---|---|---|---
| F1 [%] $\uparrow$ | IoU [%] $\uparrow$ | F1 [%] $\uparrow$ | IoU [%] $\uparrow$
Original | 1.8M | 79 | 65 | 64 | 48
Slim | 30k | 72 | 57 | 65 | 48