Singapore Public Bus Communting Patterns Study

packages = c('rgdal', 'sf', 'spdep',  'tmap', 'tidyverse','scales')
for (p in packages){
if(!require(p, character.only = T)){
install.packages(p)
}
library(p,character.only = T)
}

Readme

This notebook has code folding on for brevity, to show all code chucks, toggle the option on the top right button. Alternatively, a toggle button is available on each code chuck to show the underlying code.

Abstract

With the incorporation of automated fare collection systems in Singapore’s public transport system and her push for Open Data initiative, Singapore is well-positioned to provide opportunities for researchers to analyse big data sets and draw meaningful insights. However, most of the efforts to date tend to be siloed within individual government agencies and fail to fully utilise the datasets available across multiple domains. In this study, we aim to identify the relationship between Singapore public bus commuting patterns from smart card data and urban characteristics from census data.

1. Introduction

2. Study Area and Data Preparation

2.1 Study Area

Singapore consist of 63 islands and ongoing land reclamation projects had increased singapore size from 581.5km² to 721.5km² (increase of 23%). Her population size in 2019 is 5.8 million with a population density of 8358 per Km².

Planning areas in Singapore are known as Development Guide Plan (DGP) areas (or Zones) and are the main urban planning and census divisions of Singapore delineated by Urban Redevelopment Authority (URA). The planning areas were introduced in early 1990s after the release of the 1991 Concept plan. The Statistics Department publish census data based on the planning areas since 2000.

There are a total of 55 DGP areas and they are are futher subdivided into subzones for statastical purposes. A total of 323 subzones exist.

Fig: Map of Singapore - Subzones

2.1 Singapore Bus Network

There are four public bus operators in Singapore SBS Transit, SMRT Buses, Tower Transit Singapore and Go-Ahead Singapore. In total the network covers 5040 bus stops including 5 bus stops in Malaysia.

busStopLcn <- st_read(dsn = "data/geospatial", 
                layer = "BusStop")
## Reading layer `BusStop' from data source `/Users/daniel/Google Drive (daniel.soh.2017@smu.edu.sg)/SMU | GD/1. Geospatial Analytics and Applications/Take Home Ex/Ex01/data/geospatial' using driver `ESRI Shapefile'
## Simple feature collection with 5040 features and 3 fields
## geometry type:  POINT
## dimension:      XY
## bbox:           xmin: 4427.938 ymin: 26482.1 xmax: 48282.5 ymax: 52983.82
## proj4string:    +proj=tmerc +lat_0=1.366666666666667 +lon_0=103.8333333333333 +k=1 +x_0=28001.642 +y_0=38744.572 +datum=WGS84 +units=m +no_defs
busStopLcn <- busStopLcn %>%select('BUS_STOP_N', 'geometry','LOC_DESC')

tmap_mode("view")

tm_shape(busStopLcn)+
tm_bubbles(col = "red",
           size = 0.001,
           border.col = "red",
           border.lwd = 1)

preserve72744598b0986fa4 Fig: Bus stop network

2.2 Data Preparation

2.2.1 Data Source

2.2.1.1 Passenger Volume by Bus Stop

Source: https://www.mytransport.sg/content/mytransport/home/dataMall/dynamic-data.html
Date Range: Jan 20

Import Dataset:

paxVolbyBusStop <- read_csv ("data/aspatial/passenger volume by busstop.csv")

## Parsed with column specification:
## cols(
##   YEAR_MONTH = col_character(),
##   DAY_TYPE = col_character(),
##   TIME_PER_HOUR = col_double(),
##   PT_TYPE = col_character(),
##   PT_CODE = col_character(),
##   TOTAL_TAP_IN_VOLUME = col_double(),
##   TOTAL_TAP_OUT_VOLUME = col_double()
## )

paxVolbyBusStop

2.2.1.2 Bus Stop Location

Source: https://www.mytransport.sg/content/mytransport/home/dataMall/static-data.html

busStopLcn <- st_read(dsn = "data/geospatial", 
                layer = "BusStop")

## Reading layer `BusStop' from data source `/Users/daniel/Google Drive (daniel.soh.2017@smu.edu.sg)/SMU | GD/1. Geospatial Analytics and Applications/Take Home Ex/Ex01/data/geospatial' using driver `ESRI Shapefile'
## Simple feature collection with 5040 features and 3 fields
## geometry type:  POINT
## dimension:      XY
## bbox:           xmin: 4427.938 ymin: 26482.1 xmax: 48282.5 ymax: 52983.82
## proj4string:    +proj=tmerc +lat_0=1.366666666666667 +lon_0=103.8333333333333 +k=1 +x_0=28001.642 +y_0=38744.572 +datum=WGS84 +units=m +no_defs

busStopLcn

2.2.1.3 Subzone Area

Source: https://data.gov.sg/dataset/master-plan-2014-subzone-boundary-web
Date Range:

subzone <- st_read(dsn = "data/geospatial", 
                layer = "MP14_SUBZONE_WEB_PL")

## Reading layer `MP14_SUBZONE_WEB_PL' from data source `/Users/daniel/Google Drive (daniel.soh.2017@smu.edu.sg)/SMU | GD/1. Geospatial Analytics and Applications/Take Home Ex/Ex01/data/geospatial' using driver `ESRI Shapefile'
## Simple feature collection with 323 features and 15 fields
## geometry type:  MULTIPOLYGON
## dimension:      XY
## bbox:           xmin: 2667.538 ymin: 15748.72 xmax: 56396.44 ymax: 50256.33
## proj4string:    +proj=tmerc +lat_0=1.366666666666667 +lon_0=103.8333333333333 +k=1 +x_0=28001.642 +y_0=38744.572 +datum=WGS84 +units=m +no_defs

subzone

2.2.1.4 Population

Source: https://data.gov.sg/dataset/singapore-residents-by-subzone-and-type-of-dwelling-2011-2019
Date Range: 2011-2019 (Select 2019 only)

popVol <-read_csv("data/aspatial/respopagesextod2011to2019.csv")

## Parsed with column specification:
## cols(
##   PA = col_character(),
##   SZ = col_character(),
##   AG = col_character(),
##   Sex = col_character(),
##   TOD = col_character(),
##   Pop = col_double(),
##   Time = col_double()
## )

popVol

2.2.2 Geospatial Data Wrangling

Population volume

Select only 2019 data, group population numbers by Subzone. Change subzone to uppercase for merging data.

popVol <-read_csv("data/aspatial/respopagesextod2011to2019.csv")

## Parsed with column specification:
## cols(
##   PA = col_character(),
##   SZ = col_character(),
##   AG = col_character(),
##   Sex = col_character(),
##   TOD = col_character(),
##   Pop = col_double(),
##   Time = col_double()
## )

popVol <- popVol %>% filter(Time == 2019) %>%
                        mutate(SZ = toupper(SZ)) %>%
                        group_by(SZ) %>% 
                        summarise(population = sum(Pop))

popVol

Remove bus stop roof numbers

busStopLcn <- busStopLcn %>%select('BUS_STOP_N', 'geometry','LOC_DESC')

busStopLcn

Encode Bus Stop Location with spaital data

busStopLcnSpatial <- left_join(paxVolbyBusStop, busStopLcn, by = c("PT_CODE" = "BUS_STOP_N"))

## Warning: Column `PT_CODE`/`BUS_STOP_N` joining character vector and factor,
## coercing into character vector

busStopLcnSpatial

Encode Bus stop location with subzone classification

sz_N_Geo <- subzone %>% select("SUBZONE_N","geometry")
busStopLcnSpatial <- st_as_sf(busStopLcnSpatial)
busStopSz <- st_join(busStopLcnSpatial, sz_N_Geo, join=st_within)

busStopSz

Encode population volume with Bus Stop Volume data

busStopSzPop <- left_join(busStopSz, popVol, by = c("SUBZONE_N" = "SZ"))

## Warning: Column `SUBZONE_N`/`SZ` joining factor and character vector, coercing
## into character vector

busStopSzPop

Select subzone with geometry

subzone_filter_geo <- subzone %>% select("SUBZONE_N","geometry")

subzone_filter_geo

Remove spot location of bus stop and input geometry of subzone

busStopSzPopRemoveGeo <- busStopSzPop %>% select("DAY_TYPE","TIME_PER_HOUR","PT_CODE","TOTAL_TAP_IN_VOLUME","TOTAL_TAP_OUT_VOLUME","LOC_DESC","SUBZONE_N","population")
busStopSzPoly <- st_join( subzone_filter_geo, busStopSzPopRemoveGeo, by = c("SUBZONE_N" = "SUBZONE_N"))

busStopSzPoly

Group by Subzone and get total popluation volume for each subzone

# get tap in vol for each subzone
tapInVolSubzone <- busStopSzPoly %>% 
                        group_by(SUBZONE_N.x) %>% 
                        summarise(tapInVolAll = sum(TOTAL_TAP_IN_VOLUME)) %>% 
                        select("SUBZONE_N.x","tapInVolAll")

tapInVolSubzone

Join population volume with Tap in volume

tapInVolSubzone_df <- tapInVolSubzone %>% st_drop_geometry()

LRTablePopVolTapVol <- merge(popVol, tapInVolSubzone_df, by.x = "SZ", by.y = "SUBZONE_N.x")

Construct sf dt for analysis

Subzone with tapInVol

SZ_TapInVol <- busStopSzPoly %>% 
                        group_by(SUBZONE_N.x) %>% 
                        summarise(tapInTotal = sum(TOTAL_TAP_IN_VOLUME))

add population

SZ_TapInVol_Pop <- left_join(SZ_TapInVol, popVol, by = c("SUBZONE_N.x" = "SZ"))

## Warning: Column `SUBZONE_N.x`/`SZ` joining factor and character vector, coercing
## into character vector

3. Geospatial Analysis

3.1 Compute Spatial Weight Matrix

Creating (QUEEN) contiguity based neighbours

SZ_TapInVol_Pop_q <- poly2nb(SZ_TapInVol_Pop, queen=TRUE)
summary(SZ_TapInVol_Pop_q)

## Neighbour list object:
## Number of regions: 323 
## Number of nonzero links: 1934 
## Percentage nonzero weights: 1.853751 
## Average number of links: 5.987616 
## 5 regions with no links:
## 176 208 227 255 258
## Link number distribution:
## 
##  0  1  2  3  4  5  6  7  8  9 10 11 12 14 17 
##  5  2  6 10 26 77 87 51 34 16  3  3  1  1  1 
## 2 least connected regions:
## 42 109 with 1 link
## 1 most connected region:
## 40 with 17 links

Creating (ROOK) contiguity based neighbours

SZ_TapInVol_Pop_r <- poly2nb(SZ_TapInVol_Pop, queen=FALSE)
summary(SZ_TapInVol_Pop_r)

## Neighbour list object:
## Number of regions: 323 
## Number of nonzero links: 1680 
## Percentage nonzero weights: 1.610291 
## Average number of links: 5.201238 
## 5 regions with no links:
## 176 208 227 255 258
## Link number distribution:
## 
##  0  1  2  3  4  5  6  7  8  9 10 11 13 14 
##  5  2  6 19 73 97 65 29 17  5  2  1  1  1 
## 2 least connected regions:
## 42 109 with 1 link
## 1 most connected region:
## 40 with 14 links

Plotting Queen and Rook contiguity based neighbours maps

par(mfrow=c(1,2))
plot(SZ_TapInVol_Pop$geometry, border="lightgrey", main="Queen Contiguity",)
plot(SZ_TapInVol_Pop_q, sf::st_centroid(SZ_TapInVol_Pop$geometry), pch = 19, cex = 0.3, add = TRUE, col= "red", main="Queen Contiguity")
plot(SZ_TapInVol_Pop$geometry, border="lightgrey", main="Rook Contiguity",)
plot(SZ_TapInVol_Pop_r, sf::st_centroid(SZ_TapInVol_Pop$geometry), pch = 19, cex = 0.3, add = TRUE, col = "red", main="Rook Contiguity")

Determine Cut off distance

coords <- sf::st_centroid(SZ_TapInVol_Pop$geometry)
k1 <- knn2nb(knearneigh(coords))
k1dists <- unlist(nbdists(k1, coords, longlat = TRUE))

## Warning in nbdists(k1, coords, longlat = TRUE): dnearneigh: longlat argument
## overrides object

summary(k1dists)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   182.5   618.1   904.8  1012.4  1180.0  5954.4

dist <- 5
wm_d62 <- dnearneigh(coords, 0, dist, longlat = TRUE)

## Warning in dnearneigh(coords, 0, dist, longlat = TRUE): dnearneigh: longlat
## argument overrides object

wm_d62

## Neighbour list object:
## Number of regions: 323 
## Number of nonzero links: 0 
## Percentage nonzero weights: 0 
## Average number of links: 0 
## 323 regions with no links:
## 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323

Plotting fixed distance weight matrix

plot(SZ_TapInVol_Pop$geometry, border="lightgrey")
plot(wm_d62, coords, add=TRUE)
plot(k1, coords, add=TRUE, col="red", length=0.08)

Computing adaptive distance weight matrix

wm_knn6 <- knn2nb(knearneigh(coords, k=6))
wm_knn6

## Neighbour list object:
## Number of regions: 323 
## Number of nonzero links: 1938 
## Percentage nonzero weights: 1.857585 
## Average number of links: 6 
## Non-symmetric neighbours list

Plot distance based neighbours

plot(SZ_TapInVol_Pop$geometry, border="lightgrey")
plot(wm_knn6, coords, pch = 19, cex = 0.6, add = TRUE, col = "red")

Weight based on IDW

dist <- nbdists(SZ_TapInVol_Pop_q, coords, longlat = TRUE)

## Warning in nbdists(SZ_TapInVol_Pop_q, coords, longlat = TRUE): dnearneigh:
## longlat argument overrides object

ids <- lapply(dist, function(x) 1/(x))

Row-standardised weights matrix

rswm_q <- nb2listw(SZ_TapInVol_Pop_q, style="W", zero.policy = TRUE)
print(rswm_q, zero.policy=TRUE)

## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 323 
## Number of nonzero links: 1934 
## Percentage nonzero weights: 1.853751 
## Average number of links: 5.987616 
## 5 regions with no links:
## 176 208 227 255 258
## 
## Weights style: W 
## Weights constants summary:
##     n     nn  S0       S1       S2
## W 318 101124 318 111.2781 1309.875

set.ZeroPolicyOption(TRUE)

## [1] FALSE

Measure of Global Spatial Autocorrelation

moran.test(SZ_TapInVol_Pop$population, listw=rswm_q, zero.policy = TRUE, na.action=na.omit)

## 
##  Moran I test under randomisation
## 
## data:  SZ_TapInVol_Pop$population  
## weights: rswm_q  n reduced by no-neighbour observations
##   
## 
## Moran I statistic standard deviate = 9.9171, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       0.318260394      -0.003154574       0.001050413

Computing Monte Carlo Moran’s I

set.seed(1234)
bperm= moran.mc(SZ_TapInVol_Pop$population, listw=rswm_q, nsim=999, zero.policy = TRUE, na.action=na.omit)
bperm

## 
##  Monte-Carlo simulation of Moran I
## 
## data:  SZ_TapInVol_Pop$population 
## weights: rswm_q  
## number of simulations + 1: 1000 
## 
## statistic = 0.31826, observed rank = 1000, p-value = 0.001
## alternative hypothesis: greater

summary(bperm$res[1:999])

##      Min.   1st Qu.    Median      Mean   3rd Qu.      Max. 
## -0.088914 -0.025625 -0.005633 -0.003581  0.014853  0.186210

hist(bperm$res, freq=TRUE, breaks=20, xlab="Simulated Moran's I")
abline(v=0, col="red") 

summary(bperm$res[1:999])

##      Min.   1st Qu.    Median      Mean   3rd Qu.      Max. 
## -0.088914 -0.025625 -0.005633 -0.003581  0.014853  0.186210

hist(bperm$res, freq=TRUE, breaks=20, xlab="Simulated Geary c")
abline(v=1, col="red") 

Visualising spatial lag value using Moran Scatterplot

MS <- moran.plot(SZ_TapInVol_Pop$population, rswm_q, zero.policy = TRUE, spChk=FALSE, labels=as.character(SZ_TapInVol_Pop$SUBZONE_N.x), xlab="GDPPC(RMB)", ylab="Spatially Lag GDPPC (RMB)")

3.2 Spatial Autocorrelation Test on residuual of Simple Regression Model

Construct simple linear regresion Model Tap-In vol (Y, target variable) with residential population (X, predictor variable)

y = mX + c

lm_TapInVol = lm(tapInVolAll ~ population, data = LRTablePopVolTapVol)

Plot:

ggplot(data = LRTablePopVolTapVol, aes(x = population, y = tapInVolAll)) + 
  geom_point(color='blue') +
  geom_smooth(method = "lm", se = FALSE) + 
  ggtitle("Tap In Volume (Million) vs Population (Thousands)") + 
  labs(y="Tap in Volume (Millions)", x = "Population (Thousands)") +
  scale_y_continuous(labels = unit_format(unit = "M", scale = 1e-6)) +
  scale_x_continuous(labels = unit_format(unit = "K", scale = 1e-4))

## `geom_smooth()` using formula 'y ~ x'

## Warning: Removed 18 rows containing non-finite values (stat_smooth).

## Warning: Removed 18 rows containing missing values (geom_point).

summary(lm_TapInVol)

## 
## Call:
## lm(formula = tapInVolAll ~ population, data = LRTablePopVolTapVol)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -863762 -123199  -63290   34580 1950650 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 1.245e+05  2.107e+04   5.907 9.34e-09 ***
## population  1.924e+01  9.370e-01  20.535  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 297800 on 303 degrees of freedom
##   (18 observations deleted due to missingness)
## Multiple R-squared:  0.5819, Adjusted R-squared:  0.5805 
## F-statistic: 421.7 on 1 and 303 DF,  p-value: < 2.2e-16

plot(lm_TapInVol, pch = 16, col = "blue") #Plot the results

abline(lm_TapInVol) #Add a regression line

#
plot(lm_TapInVol$residuals,pch = 16, col = "red")

# plot(LRTablePopVolTapVol, pch = 16, col = "blue")
# plot(lm_TapInVol, pch = 16, col = "blue") #Plot the results

3.3 Local Spatial Clustering

Cluster and Outlier Analysis

fips <- order(SZ_TapInVol_Pop$SUBZONE_N.x)
localMI <- localmoran(SZ_TapInVol_Pop$population, rswm_q)
head(localMI)

##              Ii        E.Ii    Var.Ii         Z.Ii Pr(z > 0)
## 1  5.918641e-03 -0.00310559 0.1588547  0.022641761 0.4909680
## 2  4.850453e-01 -0.00310559 0.4823921  0.702836048 0.2410790
## 3 -2.414453e-02 -0.00310559 0.1184125 -0.061139949 0.5243761
## 4  7.637794e-05 -0.00310559 0.1588547  0.007983546 0.4968151
## 5  3.593989e-02 -0.00310559 0.1049318  0.120536135 0.4520292
## 6 -7.518059e-02 -0.00310559 0.1588547 -0.180835891 0.5717518

Mapping the local Moran’s I

SZ_TapInVol_Pop.localMI <- cbind(SZ_TapInVol_Pop,localMI)
localMI.map <- tm_shape(SZ_TapInVol_Pop.localMI) +
  tm_fill(col = "Ii", 
          style = "pretty", 
          title = "local moran statistics") +
  tm_borders(alpha = 0.5)

Mapping local Moran’s I p-values

pvalue.map <- tm_shape(SZ_TapInVol_Pop.localMI) +
  tm_fill(col = "Pr.z...0.", 
          breaks=c(-Inf, 0.001, 0.01, 0.05, 0.1, Inf),
          palette="-Blues", 
          title = "local Moran's I p-values") +
  tm_borders(alpha = 0.5)

Plot

tmap_arrange(localMI.map, pvalue.map, asp=1, ncol=2)

## Warning: The shape SZ_TapInVol_Pop.localMI is invalid. See sf::st_is_valid

## Variable(s) "Ii" contains positive and negative values, so midpoint is set to 0. Set midpoint = NA to show the full spectrum of the color palette.

## Warning: The shape SZ_TapInVol_Pop.localMI is invalid. See sf::st_is_valid

preserveecc073060d9f6185

Creating a LISA Cluster Map

Plotting Moran scatterplot

nci <- moran.plot(SZ_TapInVol_Pop$population, rswm_q, labels=as.character(SZ_TapInVol_Pop$SUBZONE_N.x), xlab="GDPPC 2012", ylab="Spatially Lag GDPPC 2012")

#### Plotting Moran scatterplot with standardised variable

SZ_TapInVol_Pop$Z.population <- scale(SZ_TapInVol_Pop$population) %>% as.vector 

nci2 <- moran.plot(SZ_TapInVol_Pop$Z.population, rswm_q, labels=as.character(SZ_TapInVol_Pop$SUBZONE_N.x), xlab="z-GDPPC 2012", ylab="Spatially Lag z-GDPPC 2012")

### Prepare LISA Map

quadrant <- vector(mode="numeric",length=nrow(localMI))
DV <- SZ_TapInVol_Pop$population - mean(SZ_TapInVol_Pop$population)     
C_mI <- localMI[,1] - mean(localMI[,1])    
signif <- 0.05  
quadrant[DV >0 & C_mI>0] <- 4      
quadrant[DV <0 & C_mI<0] <- 1      
quadrant[DV <0 & C_mI>0] <- 2
quadrant[DV >0 & C_mI<0] <- 3
quadrant[localMI[,5]>signif] <- 0

Plot LISA Map

SZ_TapInVol_Pop.localMI$quadrant <- quadrant
colors <- c("#ffffff", "#2c7bb6", "#abd9e9", "#fdae61", "#d7191c")
clusters <- c("insignificant", "low-low", "low-high", "high-low", "high-high")

tm_shape(SZ_TapInVol_Pop.localMI) +
  tm_fill(col = "quadrant", style = "cat", palette = colors[c(sort(unique(quadrant)))+1], labels = clusters[c(sort(unique(quadrant)))+1])

## Warning: The shape SZ_TapInVol_Pop.localMI is invalid. See sf::st_is_valid

preserve2b2fed9ae6984cf5

Hot Spot and Cold Spot Area Analysis

knb <- knn2nb(knearneigh(coords, 5, longlat = TRUE), row.names=row.names(SZ_TapInVol_Pop$population))

## Warning in knearneigh(coords, 5, longlat = TRUE): dnearneigh: longlat argument
## overrides object

knb_lw <- nb2listw(knb, style = 'B')
summary(knb_lw)

## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 323 
## Number of nonzero links: 1615 
## Percentage nonzero weights: 1.547988 
## Average number of links: 5 
## Non-symmetric neighbours list
## Link number distribution:
## 
##   5 
## 323 
## 323 least connected regions:
## 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 with 5 links
## 323 most connected regions:
## 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 with 5 links
## 
## Weights style: B 
## Weights constants summary:
##     n     nn   S0   S1    S2
## B 323 104329 1615 2849 33228

Plot

plot(SZ_TapInVol_Pop$geometry, border="lightgrey")
plot(knb, coords, pch = 19, cex = 0.6, add = TRUE, col = "red")

# knb <- knn2nb(knearneigh(coords, k=k, longlat = TRUE), row.names=row.names(weekday$population))
# 
# dnb <- dnearneigh(coordinates(hunan), 0, 85, longlat = TRUE)
# dnb
# dnb_lw <- nb2listw(dnb, style = 'B')
# summary(dnb_lw)
fips <- order(SZ_TapInVol_Pop$SUBZONE_N.x)
gi.fixed <- localG(SZ_TapInVol_Pop$population, knb_lw)
gi.fixed

##   [1]  0.446286864  0.402501047 -0.324952677 -0.381309405 -0.279876292
##   [6] -0.666288247  4.418048297  1.279560147 -1.081179190  0.870516053
##  [11]  1.278953012 -1.551161739  2.910956708  3.574306980  4.416219428
##  [16]  2.752208443 -1.397017477  0.054419370 -1.570723446 -0.132436401
##  [21]  0.105513972 -1.535321692 -0.241460182  2.634572073  0.049386496
##  [26] -1.364085264  0.318870026 -1.509279623  0.976479856 -1.441800275
##  [31]  0.224138478  0.682618181  0.371987299  0.982444882 -0.359153136
##  [36] -0.105410783 -1.523963230 -1.266405584 -1.550465659  0.244753436
##  [41] -1.072903083 -1.072903083  0.095964981  0.254504693 -1.111996619
##  [46] -0.115980754  1.799909141 -1.298390642 -1.503977038  2.148983478
##  [51]  2.009260464  0.725658909 -1.523323926 -1.045065269 -1.365041495
##  [56]  0.112210401 -0.746460755 -0.650810589 -0.530220594 -1.568967909
##  [61] -0.514771571  4.128671525  0.910353055 -0.787251298 -1.490236861
##  [66]  0.619699660  0.402501047 -0.439122200 -1.488555924 -0.542001177
##  [71] -0.906987806  0.618432197 -1.055815455 -0.323092599  1.141034099
##  [76] -0.642676133 -1.133487197  1.197359379  3.243785843 -1.523990434
##  [81]  0.470024277  2.147004044 -0.198696943  0.942756588 -0.525432192
##  [86]  1.763918184 -1.388497166  0.463443288  0.965380284 -1.570723446
##  [91] -1.570723446 -0.457676008 -0.915678425 -0.394620879 -0.573890301
##  [96] -0.105365805  0.901957596  1.012378012  1.627207344  1.031029123
## [101]  1.567440334 -1.013277227 -0.885060538 -1.493479783  0.795091746
## [106] -1.570723446  0.031269586 -0.555069888 -1.570723446 -1.570723446
## [111] -0.072497325  0.920218967  0.346047339  0.216163249  0.065437789
## [116] -0.530286630 -1.228094430 -0.947131026 -0.302395819  1.404307728
## [121]  0.334490506  0.985289583  1.386668619  0.291574364  1.966446254
## [126] -0.967320023  1.073889896  1.762039682  0.355748654  1.838962419
## [131]  0.137414839  0.557514878 -1.365632429 -0.391787298 -1.368230074
## [136] -0.526146755 -1.419465473 -1.570723446 -0.432476787 -0.457088504
## [141]  0.289977420 -0.022088311  1.794642527  1.299982726  0.409765109
## [146]  0.103128098 -1.431985397  0.096622738 -1.048583327  2.805832303
## [151]  0.332253522  1.166410784 -0.323296306 -1.324195780 -1.294100846
## [156] -0.007542426 -1.570723446  0.555040144 -0.348428831  1.367686297
## [161]  2.939813623 -1.239428383 -0.356578628  1.948941467 -1.462552955
## [166] -1.160627378 -1.525506008 -0.077238182  0.626215213 -1.378880705
## [171] -0.377060847  0.092469954 -0.946536441 -1.149038156  0.156679272
## [176] -1.466895925  1.658130693  1.483727390 -0.493755093 -1.299512085
## [181] -1.180116428 -1.395127183 -0.049675333  0.287387925 -1.465278168
## [186] -0.808397165  4.321598236 -0.125219712  1.599681472  1.808185247
## [191] -0.497535876 -1.420376476  1.134309521 -0.603672907 -0.141464881
## [196] -1.217822437 -0.466544664  1.371649895 -0.939983792 -1.100579468
## [201] -1.533606361 -1.570723446 -1.570723446 -1.368558086  0.029336110
## [206] -1.564202877  0.249016885 -0.073751280  1.434506486  2.769005948
## [211]  0.975020080 -0.594761688 -1.530095286 -0.640902778 -0.719538402
## [216] -0.344217968  2.412254869 -0.834810501 -1.365827105  1.762039682
## [221] -1.570723446  1.640575261 -1.448726355  0.982651603  0.268829383
## [226]  0.149370439 -1.362065239 -0.048731495  0.618260713  0.356772935
## [231] -0.349882430  2.028065735 -0.922573567  4.216718288  0.276102325
## [236]  1.781745682 -1.007049370  1.550371981  2.407325221  1.793388572
## [241] -0.973907106  0.033518637  0.776666326  2.052314802  1.489429740
## [246]  0.588691087 -1.570723446  1.105602225  7.110190358 -0.516899181
## [251]  0.693669532 -0.636526543  0.455953011 -1.354229960 -1.504013010
## [256]  0.421556833 -1.531349241 -1.383131692 -1.176939044  0.139215056
## [261] -0.250156990 -0.092879541 -0.949798112  0.568479488  3.560509583
## [266]  5.478506586  4.670983628 -1.290645045 -0.263234996 -0.810324785
## [271] -1.275976117 -1.020198605 -1.404352496  1.354308654 -0.701684489
## [276] -0.677583110 -0.272222799  2.154527777 -1.570723446  0.624952766
## [281] -0.407941195 -0.270439831  0.436825467  0.522904802  0.133276917
## [286] -0.523168959  1.099674660  2.042242243 -1.570723446 -1.570723446
## [291] -1.570723446 -1.570723446 -1.570723446  1.056062693 -0.711683573
## [296] -0.864434979  0.064113693 -0.314240676  0.888569067 -1.474413530
## [301]  1.717248979  3.863315041  0.137787318  2.714341448  1.196127668
## [306]  1.290610319  1.322152066  2.676618507  0.503858235 -0.613603310
## [311]  4.515270787  0.339773451 -0.146229913 -0.381430962 -0.400532103
## [316] -0.486742669  3.335792432  1.667408106  2.087137798 -0.256567890
## [321] -0.294779715  1.194184922  1.045473784
## attr(,"gstari")
## [1] FALSE
## attr(,"call")
## localG(x = SZ_TapInVol_Pop$population, listw = knb_lw)
## attr(,"class")
## [1] "localG"

SZ_TapInVol_Pop.gi <- cbind(SZ_TapInVol_Pop, as.matrix(gi.fixed))
names(SZ_TapInVol_Pop.gi)[6] <- "gstat"

Adaptative Distance

fips <- order(SZ_TapInVol_Pop$SUBZONE_N.x)
gi.adaptive <- localG(SZ_TapInVol_Pop$population, knb_lw)
SZ_TapInVol_Pop.gi <- cbind(SZ_TapInVol_Pop, as.matrix(gi.adaptive))
names(SZ_TapInVol_Pop.gi)[6] <- "gstat_adaptive"

3.4.1 Analysis

4. Conclusion

5. References

https://en.wikipedia.org/wiki/Singapore https://en.wikipedia.org/wiki/Planning_Areas_of_Singapore https://www.worldometers.info/world-population/singapore-population/