Mapping Brunei: An Introduction to Spatial Data Analysis with {bruneimap} in R

Brunei R User Group Meetup

Alvin Bong

February 28, 2025

Preliminaries

Welcome to the First Brunei R User Group Meetup of 2025!

Celebrating our 1st Anniversary! 🎉

“The RUGS mission is to facilitate the person-to-person exchange of knowledge in small group settings on a global scale.” — R Consortium

About Us

• A community of UBD-ians and R enthusiasts
• Building an R user community in Brunei
• Advocating for Open Source collaboration

{width=“5% fig-align=”right”}

Prelim (Cont.)

Past Events: Intro to R Workshop, R>aya Meetup & Sharing Sessions

Goal: In conjunction with National Day, our team present you {bruneimap} package, containing spatial data on Brunei (boundaries & amenities)

Plan for today

📍 Part I: Refresher on spatial data & analysis

🗺️ Part II: Demonstration of {bruneimap}datasets & features

🎉 BONUS: Refreshment & networking

I) Introduction

Types of spatial data

• Points
• Lines
• Polygons
• Raster

Location of Masjids in Brunei

Brunei Highways

Population by kampong in Brunei

An important spatial data is raster or grid data. These are Pixelated (or gridded) data where each pixel is associated with a specific geographical location E.g. satellite imagery, digital elevation models, etc.

Used alot in geoscience!

Not discussing these today…

Spatial data in R

• Simple features refers to a formal standard (ISO 19125-1:2004) that describes how objects in the real world can be represented in computers, with emphasis on the spatial geometry of these objects.
• The {sf} package (Pebesma and Bivand 2023) provides support for handling spatial data in R.

R has well-supported classes for storing spatial data (sp) and interfacing to the above mentioned environments (rgdal, rgeos), but has so far lacked a complete implementation of simple features, making conversions at times convoluted, inefficient or incomplete. The package sf tries to fill this gap, and aims at succeeding sp in the long term.

The standard is widely implemented in spatial databases (such as PostGIS), commercial GIS (e.g., ESRI ArcGIS) and forms the vector data basis for libraries such as GDAL. A subset of simple features forms the GeoJSON standard.

Spatial data in R (cont.)

pop_sf |>
  select(kampong, mukim, district, population, household, geometry)

Simple feature collection with 438 features and 5 fields
Geometry type: GEOMETRY
Dimension:     XY
Bounding box:  xmin: 114.0759 ymin: 4.002635 xmax: 115.3603 ymax: 5.045878
Geodetic CRS:  WGS 84
# A tibble: 438 × 6
   kampong         mukim district population household                  geometry
   <chr>           <chr> <chr>         <dbl>     <dbl>             <POLYGON [°]>
 1 Kg. Biang       Muki… Temburo…         75        13 ((115.1546 4.66665, 115.…
 2 Kg. Amo         Muki… Temburo…        394        83 ((115.1626 4.664215, 115…
 3 Kg. Sibut       Muki… Temburo…        192        37 ((115.1263 4.625564, 115…
 4 Kg. Sumbiling … Muki… Temburo…         91        23 ((115.1301 4.596809, 115…
 5 Kg. Batang Duri Muki… Temburo…        108        23 ((115.1923 4.580851, 115…
 6 Kg. Sumbiling … Muki… Temburo…        143        23 ((115.1427 4.593578, 115…
 7 Kg. Selapon     Muki… Temburo…        199        38 ((115.2182 4.697223, 115…
 8 Kg. Lepong      Muki… Temburo…        123        26 ((115.088 4.60364, 115.0…
 9 Kg. Semabat     Muki… Temburo…         95        26 ((115.1047 4.578915, 115…
10 Kg. Temada      Muki… Temburo…         90        23 ((115.1473 4.527759, 115…
# ℹ 428 more rows

• E.g. Data set with (\(n=438\)) features containing spatial attributes (geometry) and feature-related variables (population, household, occ_liv_q).

• Notice the use of {dplyr} verbs! Link: Tidyverse methods for sf objects

Spatial analysis

To assess – or at least, to account – for spatial patterns in the data.

• Crucial to acknowledge the inherent non-randomness and interdependence within spatial data, increasing the validity and accuracy.

• Numerous motivations for doing so including

  • Agriculture (Besag and Higdon 1999)
  • Ecology (Brewer and Nolan 2007)
  • Education (Wall 2004)
  • Epidemiology (Xie et al. 2017; Li et al. 2011)
  • Image analysis (Gavin and Jennison 1997)

• Spatial statistics are indispensable for capturing the complexity of spatial data, leading to more accurate, reliable, and nuanced insights and decisions.

• Melanie Wall: To provide an example of this type of ecological spatial regression using these models we consider state level summary data related to the SAT college entrance exam for the year 1999.

• Spatial GAM to include X, Y coordinates for better disease mapping.

• Athma exacerbation and proximity of residence to major roads in Detroit, Michigan.

Spatial analysis (cont.)

Many areas of analysis:

• spatial autocorrelation
• spatial interpolation
• disease (risk) modelling
• mobility networks

Statistical tests: Moran I (global, local), Getis Ord, Geary’s C

Will not go too theoretical, focus more on applications in Section II

Types of spatial patterns

• Zero autocorrelation
• Positive autocorrelation
• Negative autocorrelation

• Values are randomly distributed in space, independently of each other (and particularly its neighbours)

• Lack of a spatial pattern.

• Values are non-randomly distributed in space.

• Large values are “attracting” other large values, and vice versa.

• A clustering effect is present.

• Values are non-randomly distributed in space.

• Large values are “repelling” other large values, and vice versa.

• A dispersion effect is present.

How to meaure this pattern?

One way is to use Global Moran’s I test (1948), defined as: \[ I = \frac{M}{\sum_{i}{\sum_{j}{w_{ij}}}} \frac{\sum_{i}{\sum_{j}{w_{ij}(y_i - \bar{y})(y_j - \bar{y})}}}{\sum_{i}{(y_i - \bar{y})^2}} \in [-1, 1] \] where \(w_{ij}\) is the “spatial weight” between locations \(i\) and \(j\).

Interpretation

• Value \(\approx +1\): +ve autocorrelation

• Value \(\approx -1\): -ve autocorrelation

• Values \(\approx 0\): randomness

Hypothesis test

• \(H_0: I = 0\) (no spatial autocorrelation),

• \(H_1: I \neq 0\) (presence of spatial autocorrelation).

–>

II) The {bruneimap} package

Spatial datasets

• Brunei
• District
• Mukim
• Kampong
• Forest Reserves

# install.packages("bruneimap")       # library(bruneimap)
brn_sf

Simple feature collection with 17 features and 1 field
Geometry type: GEOMETRY
Dimension:     XY
Bounding box:  xmin: 114.0759 ymin: 4.002644 xmax: 115.3647 ymax: 5.045878
Geodetic CRS:  WGS 84
# A tibble: 17 × 2
   name                                                                 geometry
   <chr>                                                          <GEOMETRY [°]>
 1 Mainland              MULTIPOLYGON (((115.1704 4.868704, 115.1704 4.868901, …
 2 Pulau Bedukang        POLYGON ((115.0622 4.982716, 115.0623 4.982782, 115.06…
 3 Pulau Muara Besar     POLYGON ((115.1302 4.996881, 115.1301 4.996724, 115.13…
 4 Pulau Si Mangga Besar POLYGON ((115.0451 4.981587, 115.0452 4.981653, 115.04…
 5 Pulau Suhung Damit    POLYGON ((115.0415 4.999762, 115.0416 4.999756, 115.04…
 6 Pulau ...             POLYGON ((115.0435 4.997858, 115.0436 4.997846, 115.04…
 7 Pulau ...             POLYGON ((115.0419 4.995619, 115.042 4.995607, 115.042…
 8 Pulau Chermin         POLYGON ((115.025 4.931068, 115.0251 4.931383, 115.025…
 9 Pulau Kaingaran       POLYGON ((115.029 4.950193, 115.0291 4.950195, 115.029…
10 Pulau Baru-Baru       POLYGON ((115.0413 4.908514, 115.0413 4.908514, 115.04…
11 Pulau Berbunut        POLYGON ((115.0469 4.896099, 115.0469 4.896099, 115.04…
12 Pulau Selirong        MULTIPOLYGON (((115.1292 4.861138, 115.1295 4.86126, 1…
13 Pulau Siarau          POLYGON ((115.0518 4.81406, 115.0518 4.814058, 115.051…
14 Pulau Sibungor        POLYGON ((114.9557 4.866652, 114.9558 4.866528, 114.95…
15 Pulau Tarap           POLYGON ((115.0733 4.836794, 115.0735 4.836788, 115.07…
16 Pulau Selanjak        POLYGON ((115.0662 4.829475, 115.0662 4.829385, 115.06…
17 Pulau Pepatan         POLYGON ((115.0443 4.910401, 115.0441 4.910399, 115.04…

• Notice geometry type is Geometry due to “multipolygon”
• Essentially still areal data

dis_sf

Simple feature collection with 4 features and 6 fields
Geometry type: POLYGON
Dimension:     XY
Bounding box:  xmin: 114.0759 ymin: 4.002644 xmax: 115.3647 ymax: 5.045878
Geodetic CRS:  WGS 84
# A tibble: 4 × 7
     id district                           geometry     X     Y perimeter   area
* <dbl> <chr>                         <POLYGON [°]> <dbl> <dbl>       [m]  [m^2]
1     1 Belait       ((114.5232 4.708587, 114.5235…  115.  4.36   320382. 2.81e9
2     2 Tutong       ((114.7945 4.924489, 114.7946…  115.  4.62   205103. 1.22e9
3     3 Temburong    ((115.1705 4.86852, 115.1705 …  115.  4.59   263505. 1.26e9
4     8 Brunei-Muara ((115.0594 5.044369, 115.0594…  115.  4.90   198744. 5.46e8

mkm_sf

Simple feature collection with 39 features and 7 fields
Geometry type: GEOMETRY
Dimension:     XY
Bounding box:  xmin: 114.0759 ymin: 4.002644 xmax: 115.3647 ymax: 5.045878
Geodetic CRS:  WGS 84
# A tibble: 39 × 8
      id mukim   district                  geometry     X     Y perimeter   area
 * <dbl> <chr>   <chr>                <POLYGON [°]> <dbl> <dbl>       [m]  [m^2]
 1     1 Mukim … Belait   ((114.8557 4.279265, 114…  115.  4.15   123890. 6.10e8
 2     2 Mukim … Belait   ((114.4861 4.164966, 114…  115.  4.32   140813. 8.19e8
 3     3 Mukim … Belait   ((114.421 4.558823, 114.…  114.  4.45    79425. 3.21e8
 4     4 Mukim … Belait   ((114.2525 4.596173, 114…  114.  4.56    57841. 8.83e7
 5     5 Mukim … Belait   ((114.554 4.646316, 114.…  115.  4.52    81222. 2.99e8
 6     6 Mukim … Belait   ((114.4084 4.648287, 114…  114.  4.59    56930. 1.65e8
 7     7 Mukim … Temburo… ((115.1572 4.669601, 115…  115.  4.49   152257. 5.69e8
 8     8 Mukim … Belait   ((114.5232 4.708587, 114…  114.  4.63    53176. 1.45e8
 9     9 Mukim … Temburo… ((115.1016 4.761736, 115…  115.  4.66    96041. 2.29e8
10    10 Mukim … Tutong   ((114.6139 4.781696, 114…  115.  4.72    44979. 7.99e7
# ℹ 29 more rows

• Consists 39 mukims
• Related variables: centre x, y, perimeter, area, etc

kpg_sf

Simple feature collection with 438 features and 8 fields
Geometry type: GEOMETRY
Dimension:     XY
Bounding box:  xmin: 114.0759 ymin: 4.002635 xmax: 115.3603 ymax: 5.045878
Geodetic CRS:  WGS 84
# A tibble: 438 × 9
      id kampong  mukim district                  geometry     X     Y perimeter
 * <dbl> <chr>    <chr> <chr>                <POLYGON [°]> <dbl> <dbl>       [m]
 1     1 Kg. Bia… Muki… Temburo… ((115.1546 4.66665, 115.…  115.  4.67    21056.
 2     2 Kg. Amo  Muki… Temburo… ((115.1626 4.664215, 115…  115.  4.63    38684.
 3     3 Kg. Sib… Muki… Temburo… ((115.1263 4.625564, 115…  115.  4.61    16291.
 4     4 Kg. Sum… Muki… Temburo… ((115.1301 4.596809, 115…  115.  4.60    12994.
 5     5 Kg. Bat… Muki… Temburo… ((115.1923 4.580851, 115…  115.  4.57    33791.
 6     6 Kg. Sum… Muki… Temburo… ((115.1427 4.593578, 115…  115.  4.59    13796.
 7     7 Hutan S… Muki… Temburo… ((115.2727 4.566723, 115…  115.  4.44   114731.
 8     8 Kg. Sel… Muki… Temburo… ((115.2182 4.697223, 115…  115.  4.63    57100.
 9     9 Kg. Lep… Muki… Temburo… ((115.088 4.60364, 115.0…  115.  4.60    16377.
10    10 Kg. Sem… Muki… Temburo… ((115.1047 4.578915, 115…  115.  4.57    18679.
# ℹ 428 more rows
# ℹ 1 more variable: area [m^2]

• Standardised to “Kg. XXX”

fr_sf

Simple feature collection with 20 features and 5 fields
Geometry type: POLYGON
Dimension:     XY
Bounding box:  xmin: 114.1894 ymin: 4.002635 xmax: 115.3602 ymax: 4.998087
Geodetic CRS:  WGS 84
# A tibble: 20 × 6
      id name                  category   area    prop                  geometry
   <dbl> <chr>                 <chr>     [m^2]   <dbl>             <POLYGON [°]>
 1    18 Batu Apoi FR. 1st Ex… Proposed 2.15e8 3.68e-2 ((115.2432 4.741331, 115…
 2    20 Forest Plantation Ar… Plantat… 5.35e8 9.15e-2 ((114.4972 4.609854, 114…
 3    19 Biang FR.             Gazetted 2.82e7 4.83e-3 ((115.0791 4.702919, 115…
 4     1 Belait Peat Swamp FR. Proposed 6.91e8 1.18e-1 ((114.1894 4.54777, 114.…
 5     2 Badas FR.             Gazetted 1.42e6 2.43e-4 ((114.4237 4.5686, 114.4…
 6     3 Anduki FR.            Gazetted 1.17e7 2.01e-3 ((114.44 4.645843, 114.4…
 7     4 Andulau FR.           Gazetted 1.35e8 2.31e-2 ((114.5371 4.661383, 114…
 8     5 Arboretum FR.         Gazetted 1.54e6 2.63e-4 ((114.4909 4.672378, 114…
 9     6 Labi FR.              Gazetted 1.22e9 2.09e-1 ((114.4915 4.551124, 114…
10     9 Ladan FR.             Gazetted 2.82e8 4.83e-2 ((114.7375 4.703889, 114…
11     8 Ladan FR. Ext.        Proposed 2.06e7 3.53e-3 ((114.7806 4.409001, 114…
12     7 Labi FR. Ext.         Proposed 3.78e7 6.48e-3 ((114.738 4.390782, 114.…
13    10 Bukit Shahbandar FR.  Gazetted 1.93e7 3.30e-3 ((114.851 4.944191, 114.…
14    11 Berakas FR.           Gazetted 4.46e6 7.64e-4 ((114.9109 4.990735, 114…
15    12 Subok FR.             Gazetted 1.52e6 2.60e-4 ((114.9632 4.885511, 114…
16    13 Selirong FR.          Gazetted 2.40e7 4.11e-3 ((115.1378 4.903892, 115…
17    14 Labu FR.              Gazetted 1.43e8 2.45e-2 ((115.1497 4.872417, 115…
18    15 Peradayan FR.         Gazetted 9.89e6 1.69e-3 ((115.1735 4.756919, 115…
19    16 Batu Apoi FR.         Gazetted 4.68e8 8.01e-2 ((115.2747 4.566675, 115…
20    17 Batu Apoi FR. 2nd Ex… Proposed 2.78e7 4.76e-3 ((115.0817 4.514593, 115…

EDA (Biggest area)

• District
• Mukim
• Kampong
• Forest Reserves

dis_sf %>%
  select(district, perimeter, area) %>%
  arrange(desc(area))

Simple feature collection with 4 features and 3 fields
Geometry type: POLYGON
Dimension:     XY
Bounding box:  xmin: 114.0759 ymin: 4.002644 xmax: 115.3647 ymax: 5.045878
Geodetic CRS:  WGS 84
# A tibble: 4 × 4
  district     perimeter        area                                    geometry
  <chr>              [m]       [m^2]                               <POLYGON [°]>
1 Belait         320382. 2813605686. ((114.5232 4.708587, 114.5235 4.708474, 11…
2 Temburong      263505. 1259343398. ((115.1705 4.86852, 115.1705 4.868362, 115…
3 Tutong         205103. 1222213870. ((114.7945 4.924489, 114.7946 4.924442, 11…
4 Brunei-Muara   198744.  546395598. ((115.0594 5.044369, 115.0594 5.044021, 11…

mkm_sf %>%
  select(mukim, district, area) %>%
  arrange(desc(area))

Simple feature collection with 39 features and 3 fields
Geometry type: GEOMETRY
Dimension:     XY
Bounding box:  xmin: 114.0759 ymin: 4.002644 xmax: 115.3647 ymax: 5.045878
Geodetic CRS:  WGS 84
# A tibble: 39 × 4
   mukim             district        area                               geometry
   <chr>             <chr>          [m^2]                          <POLYGON [°]>
 1 Mukim Sukang      Belait    818701681. ((114.4861 4.164966, 114.4888 4.29140…
 2 Mukim Melilas     Belait    610141519. ((114.8557 4.279265, 114.8558 4.27923…
 3 Mukim Amo         Temburong 569489572. ((115.1572 4.669601, 115.1572 4.66960…
 4 Mukim Rambai      Tutong    556510693. ((114.6646 4.65737, 114.6654 4.648215…
 5 Mukim Labi        Belait    365579560. ((114.4861 4.164966, 114.486 4.165036…
 6 Mukim Kuala Balai Belait    321489759. ((114.421 4.558823, 114.4216 4.558641…
 7 Mukim Bukit Sawat Belait    299075426. ((114.554 4.646316, 114.5544 4.646191…
 8 Mukim Labu        Temburong 244496279. ((115.168 4.873203, 115.1682 4.873133…
 9 Mukim Batu Apoi   Temburong 229337396. ((115.1016 4.761736, 115.1056 4.74522…
10 Mukim Seria       Belait    164891541. ((114.4084 4.648287, 114.4086 4.64782…
# ℹ 29 more rows

kpg_sf %>%
  select(kampong, mukim, district, area) %>%
  arrange(desc(area))

Simple feature collection with 438 features and 4 fields
Geometry type: GEOMETRY
Dimension:     XY
Bounding box:  xmin: 114.0759 ymin: 4.002635 xmax: 115.3603 ymax: 5.045878
Geodetic CRS:  WGS 84
# A tibble: 438 × 5
   kampong                  mukim      district   area                  geometry
   <chr>                    <chr>      <chr>     [m^2]             <POLYGON [°]>
 1 Kg. Bukit Tuding         Mukim Mel… Belait   5.63e8 ((114.8566 4.278318, 114…
 2 Hutan Simpan Batu Apoi   Mukim Amo  Temburo… 4.23e8 ((115.2727 4.566723, 115…
 3 Kg. Dungun               Mukim Suk… Belait   3.56e8 ((114.8444 4.315076, 114…
 4 Hutan Simpan Kuala Balai Mukim Kua… Belait   2.70e8 ((114.4216 4.557172, 114…
 5 Hutan Simpan Bukit Ladan Mukim Ram… Tutong   2.47e8 ((114.7998 4.656117, 114…
 6 Kg. Sukang               Mukim Suk… Belait   1.90e8 ((114.6014 4.375045, 114…
 7 Kg. Selapon              Mukim Bat… Temburo… 1.46e8 ((115.2182 4.697223, 115…
 8 Hutan Simpan Labu        Mukim Labu Temburo… 1.37e8 ((115.1705 4.86852, 115.…
 9 Kg. Kukub                Mukim Suk… Belait   1.30e8 ((114.7368 4.396723, 114…
10 Badas                    Mukim Ser… Belait   1.11e8 ((114.4287 4.608144, 114…
# ℹ 428 more rows

• On mukim & kampong level, not just those in Belait are largest in size

fr_sf %>%
  select(name, area) %>%
  arrange(desc(area))

Simple feature collection with 20 features and 2 fields
Geometry type: POLYGON
Dimension:     XY
Bounding box:  xmin: 114.1894 ymin: 4.002635 xmax: 115.3602 ymax: 4.998087
Geodetic CRS:  WGS 84
# A tibble: 20 × 3
   name                          area                                   geometry
   <chr>                        [m^2]                              <POLYGON [°]>
 1 Labi FR.               1223008338. ((114.4915 4.551124, 114.4968 4.549981, 1…
 2 Belait Peat Swamp FR.   691025583. ((114.1894 4.54777, 114.1936 4.5507, 114.…
 3 Forest Plantation Area  534737587. ((114.4972 4.609854, 114.5959 4.607429, 1…
 4 Batu Apoi FR.           467847280. ((115.2747 4.566675, 115.2745 4.566074, 1…
 5 Ladan FR.               282206090. ((114.7375 4.703889, 114.7384 4.7038, 114…
 6 Batu Apoi FR. 1st Ext.  214762907. ((115.2432 4.741331, 115.2427 4.739851, 1…
 7 Labu FR.                143194103. ((115.1497 4.872417, 115.1504 4.872185, 1…
 8 Andulau FR.             135181331. ((114.5371 4.661383, 114.5416 4.661389, 1…
 9 Labi FR. Ext.            37826575. ((114.738 4.390782, 114.7378 4.394926, 11…
10 Biang FR.                28216005. ((115.0791 4.702919, 115.1072 4.701784, 1…
11 Batu Apoi FR. 2nd Ext.   27828335. ((115.0817 4.514593, 115.0865 4.515219, 1…
12 Selirong FR.             24006733. ((115.1378 4.903892, 115.1381 4.903722, 1…
13 Ladan FR. Ext.           20605750. ((114.7806 4.409001, 114.7818 4.352634, 1…
14 Bukit Shahbandar FR.     19293602. ((114.851 4.944191, 114.8629 4.944299, 11…
15 Anduki FR.               11728221. ((114.44 4.645843, 114.4288 4.643316, 114…
16 Peradayan FR.             9885611. ((115.1735 4.756919, 115.1736 4.757048, 1…
17 Berakas FR.               4463776. ((114.9109 4.990735, 114.9232 4.99749, 11…
18 Arboretum FR.             1535690. ((114.4909 4.672378, 114.4912 4.672647, 1…
19 Subok FR.                 1519648. ((114.9632 4.885511, 114.9633 4.884182, 1…
20 Badas FR.                 1419018. ((114.4237 4.5686, 114.4213 4.562911, 114…

Areal Plots

• Brunei
• District
• Mukim
• Kampong
• Forest Reserves

mapview(brn_sf, zcol="name", layer.name = "Country Boundary")

mapview(dis_sf, zcol="district", layer.name = "District")

mapview(mkm_sf, zcol="mukim", layer.name = "Mukim")

mapview(kpg_sf, zcol="kampong", layer.name = "Kampong")

mapview(fr_sf, zcol="name", layer.name = "Forest Reserves")

• 70% of Brunei are forests but not all designated as forest reserves

Study Variable Datasets: Census

• Census
• Pop. by district
• Pop. by mkm
• Pop. by kpg

census2021

# A tibble: 365 × 11
      id kampong      mukim district population pop_male pop_female pop_bruneian
   <dbl> <chr>        <chr> <chr>         <dbl>    <dbl>      <dbl>        <dbl>
 1     1 Kg. Biang    Muki… Temburo…         75       46         29           37
 2     2 Kg. Amo      Muki… Temburo…        394      218        176          280
 3     3 Kg. Sibut    Muki… Temburo…        192       98         94          174
 4     4 Kg. Sumbili… Muki… Temburo…         91       48         43           55
 5     5 Kg. Batang … Muki… Temburo…        108       60         48           57
 6     6 Kg. Sumbili… Muki… Temburo…        143       68         75           64
 7     8 Kg. Selapon  Muki… Temburo…        199      115         84          114
 8     9 Kg. Lepong   Muki… Temburo…        123       65         58           88
 9    10 Kg. Semabat  Muki… Temburo…         95       52         43           63
10    12 Kg. Temada   Muki… Temburo…         90       46         44           35
# ℹ 355 more rows
# ℹ 3 more variables: pop_pr <dbl>, household <dbl>, occ_liv_q <dbl>

census2021 %>%
  group_by(district) %>%
  summarise(pop = sum(population)) %>%
  arrange(desc(pop))

# A tibble: 4 × 2
  district        pop
  <chr>         <dbl>
1 Brunei Muara 320323
2 Belait        65531
3 Tutong        47210
4 Temburong      9444

census2021 %>%
  group_by(mukim, district) %>%
  summarise(pop = sum(population)) %>%
  arrange(desc(pop))

# A tibble: 39 × 3
# Groups:   mukim [39]
   mukim              district       pop
   <chr>              <chr>        <dbl>
 1 Mukim Sengkurong   Brunei Muara 40972
 2 Mukim Mentiri      Brunei Muara 39324
 3 Mukim Berakas B    Brunei Muara 39284
 4 Mukim Gadong B     Brunei Muara 38067
 5 Mukim Gadong A     Brunei Muara 35424
 6 Mukim Kuala Belait Belait       28793
 7 Mukim Berakas A    Brunei Muara 28311
 8 Mukim Kilanas      Brunei Muara 24981
 9 Mukim Serasa       Brunei Muara 18569
10 Mukim Seria        Belait       18313
# ℹ 29 more rows

census2021 %>%
  group_by(kampong, mukim, district) %>%
  summarise(pop = sum(population)) %>%
  arrange(desc(pop))

# A tibble: 364 × 4
# Groups:   kampong, mukim [364]
   kampong             mukim            district       pop
   <chr>               <chr>            <chr>        <dbl>
 1 Kg. Panchor Mentiri Mukim Mentiri    Brunei Muara 13358
 2 Kg. Tanah Jambu     Mukim Mentiri    Brunei Muara 11695
 3 Kg. Panaga          Mukim Seria      Belait       10301
 4 Kg. Bukit Beruang   Mukim Telisai    Tutong        9835
 5 Kg. Meragang        Mukim Serasa     Brunei Muara  9190
 6 Kg. Mata-Mata       Mukim Gadong B   Brunei Muara  7159
 7 Kg. Beribi          Mukim Gadong B   Brunei Muara  6490
 8 Kg. Kilanas         Mukim Kilanas    Brunei Muara  6357
 9 Kg. Sungai Akar     Mukim Berakas B  Brunei Muara  6129
10 Kg. Mulaut          Mukim Sengkurong Brunei Muara  5981
# ℹ 354 more rows

• On kampong-level, Belait & Tutong population stands out!

Census EDA (Cont.)

census2021 %>%
  summarise(sum(population), sum(pop_male), sum(pop_female))

# A tibble: 1 × 3
  `sum(population)` `sum(pop_male)` `sum(pop_female)`
              <dbl>           <dbl>             <dbl>
1            442508          233152            209356

census2021 %>%
  summarise(bruneian_pct = sum(pop_bruneian)/sum(population)*100,
            pr_pct = sum(pop_pr)/sum(population)*100)

# A tibble: 1 × 2
  bruneian_pct pr_pct
         <dbl>  <dbl>
1         75.8   5.85

• Other available datasets on masjid, schools, healthcare.

• Link: Bruneimap Reference

• If need additional datasets, can use {osmdata} package (Open Street Map)

Visualisation of Census

mapview(pop_sf, zcol = "population", label = "kampong", layer.name = "Pop. by kampong")

Example Analysis on School Data

Topic Question

• Are schools in Brunei evenly distributed?

• Any spatial patterns?

• Do school locations effectively serve the population? (are school clusters aligned with high-population areas?)

1) Are schools in Brunei evenly distributed?

• Raw Data
• By Sector
• By Cluster

Simple feature collection with 252 features and 13 fields
Geometry type: POINT
Dimension:     XY
Bounding box:  xmin: 114.1788 ymin: 4.258506 xmax: 115.198 ymax: 5.036068
Geodetic CRS:  WGS 84
# A tibble: 252 × 14
     No. School   Sector Cluster Education.Level            geometry    id
 * <dbl> <chr>    <chr>  <chr>   <chr>                   <POINT [°]> <dbl>
 1     1 Sekolah… MOE    Cluste… Pre-Primary & …  (114.8973 4.87885)    45
 2     2 Sekolah… MOE    Cluste… Pre-Primary & … (114.9132 4.892325)    66
 3     3 Sekolah… MOE    Cluste… Pre-Primary & … (114.9005 4.898354)    87
 4     4 Sekolah… MOE    Cluste… Pre-Primary & … (114.8835 4.921296)   251
 5     5 Sekolah… MOE    Cluste… Pre-Primary & … (114.8788 4.901515)   245
 6     6 Sekolah… MOE    Cluste… Pre-Primary & … (114.8848 4.902309)    87
 7     7 Sekolah… MOE    Cluste… Pre-Primary & … (114.9555 4.872499)   113
 8     8 Sekolah… MOE    Cluste… Pre-Primary & … (114.9358 4.896647)   116
 9     9 Sekolah… MOE    Cluste… Pre-Primary & …  (114.8862 4.92387)   954
10    10 Maktab … MOE    Cluste… Secondary       (114.9382 4.899599)   116
# ℹ 242 more rows
# ℹ 7 more variables: kampong <chr>, mukim <chr>, district <chr>, X <dbl>,
#   Y <dbl>, perimeter [m], area [m^2]

ggplot() +
  geom_sf(data = kpg_sf, fill = NA, col = "grey", lwd = 0.3) +
  geom_sf(data = mkm_sf, fill = NA, col = "black", lwd = 0.5) +
  geom_sf(data = dis_sf, fill = NA, col = "black", lwd = 1) +
  geom_sf(data = sch_sf, 
          aes(col = Sector),
          size = 2) +
  scale_color_brewer(palette = "Set2") +
  theme_void()

• MOE schools are more evenly distributed
• Private schools more clustered in city centers
• Limited MORA schools

ggplot() +
  geom_sf(data = kpg_sf, fill = NA, col = "grey", lwd = 0.3) +
  geom_sf(data = mkm_sf, fill = NA, col = "black", lwd = 0.5) +
  geom_sf(data = dis_sf, fill = NA, col = "black", lwd = 1) +
  geom_sf(data = sch_sf %>% filter(Cluster != "NA"), 
          aes(col = Cluster),
          size = 2) +
  scale_color_brewer(palette = "Set2") +
  theme_void()

• MOE Clusters seems to be clustered spatially, by districts/mukims

School by mukim

mkm_sch <- mkm_sf %>%
  left_join(
    sch_sf %>%
      as_tibble() %>%
      count(mukim,, name = "count_of_schools"),
    by = "mukim"
  ) %>%
  select(mukim, count_of_schools)

mapview(mkm_sch,
        zcol = "count_of_schools",
        layer.name = "School Count",
        label = "mukim")

• Most schools clustered around shoreline

• Seems to be spatially autocorrelated, lets explore more!

2) Spatial autocorrelation on school

mkm_sch %>%
  mutate(rand_count = sample(count_of_schools)) %>%
  mapview(zcol = "rand_count", layer.name = "School count")

• Hypothetical map if not clustered (random shuffle)

Spatial autocorrelation on school (cont.)

To verify, we apply Global Moran I Test :

library(sfdep)
mkm_sch$count_of_schools[is.na(mkm_sch$count_of_schools)] <- 0
nb <- st_contiguity(mkm_sch)
wt <- st_weights(nb)
global_moran_test(mkm_sch$count_of_schools, nb, wt)

    Moran I test under randomisation

data:  x  
weights: listw    

Moran I statistic standard deviate = 4.4379, p-value = 4.542e-06
alternative hypothesis: greater
sample estimates:
Moran I statistic       Expectation          Variance 
       0.45711421       -0.02631579        0.01186631 

• I value of 0.457 > 0, implying Positive Autocorrelation
• statistically significant results low p-value of \(\mathbf{4.54 \times 10^{-6}} < 0.001\)

3) School clusters aligned with high-population areas?

mapview(pop_sf %>%
          group_by(mukim) %>%
          summarise(pop = sum(population, na.rm = TRUE)),
        zcol = "pop",
        layer.name = " Population"
) |
  mapview(mkm_sch,
          zcol = "count_of_schools",
          layer.name = "School Count",
          label = "mukim")

• Generally yes except Brunei-Muara where,
• some high-population mukims do not have many schools

Implications

Q1: Are schools in Brunei evenly distributed?

A: No, mostly clustered around shoreline. MOE schools more evenly distributed. Private schools mostly in city.

Q2: Any spatial patterns?

A: Yes. Positive Autocorrelation. There are clusters, mukims with high school count are near each other.

Q3: Do school locations effectively serve the population?

A: Most school clusters aligned with high-population areas. Except Brunei-Muara, where some high-population mukims do not have many schools.

Key Takeway

• Basic recap of spatial data & spatial analysis
• Introduce & visualise datasets on {bruneimap}
• Sample analysis on sch_sf

Explore Further: Can you apply similar analysis to masjid locations in Brunei? How do their distributions align with population needs? 🕌

Next Steps

• Engage with the Brunei R User Group
• Join our future workshops and meetups
• Follow us @brunei.r.user for updates! 🚀

Thanks!

https://bruneir.github.io/brm-bruneimap

References

Besag, Julian, and David Higdon. 1999. “Bayesian Analysis of Agricultural Field Experiments.” Journal of the Royal Statistical Society Series B: Statistical Methodology 61 (4): 691–746.

Brewer, Mark J, and Andrew J Nolan. 2007. “Variable Smoothing in Bayesian Intrinsic Autoregressions.” Environmetrics: The Official Journal of the International Environmetrics Society 18 (8): 841–57.

Gavin, John, and Christopher Jennison. 1997. “A Subpixel Image Restoration Algorithm.” Journal of Computational and Graphical Statistics 6 (2): 182–201.

Li, Shi, Stuart Batterman, Elizabeth Wasilevich, Huda Elasaad, Robert Wahl, and Bhramar Mukherjee. 2011. “Asthma Exacerbation and Proximity of Residence to Major Roads: A Population-Based Matched Case-Control Study Among the Pediatric Medicaid Population in Detroit, Michigan.” Environmental Health 10 (1): 1–10.

Moran, Patrick AP. 1948. “The Interpretation of Statistical Maps.” Journal of the Royal Statistical Society. Series B (Methodological) 10 (2): 243–51.

Pebesma, Edzer, and Roger Bivand. 2023. Spatial Data Science: With applications in R. Chapman and Hall/CRC. https://doi.org/10.1201/9780429459016.

Wall, Melanie M. 2004. “A Close Look at the Spatial Structure Implied by the CAR and SAR Models.” Journal of Statistical Planning and Inference 121 (2): 311–24.

Xie, Sherrie, Rebecca Greenblatt, Michael Z Levy, and Blanca E Himes. 2017. “Enhancing Electronic Health Record Data with Geospatial Information.” AMIA Summits on Translational Science Proceedings 2017: 123.