1. Background Information
This analysis aims to investigate if the distribution of Airbnb listings are affected by location factors such as near to existing hotels, MRT services and tourist attractions and analyse the impact of COVID-19 on Airbnb business in Singapore.
2. Dataset
- Aspatial dataset:
- Airbnb listings as at June 2019 downloaded from eLearn
- Airbnb listings as at June 2021, specifically 29 June, 2021.csv from Inside Airbnb under show archived data
- Hotels
- Tourist attractions
- Geospatial dataset:
- MRTLRT services
- CoastalOutline: a polygon feature data showing the national boundary of Singapore. It is provided by SLA and is in ESRI shapefile format.
- MP14_SUBZONE_WEB_PL: a polygon feature data providing information of URA 2014 Master Plan Planning Subzone boundary data. It is in ESRI shapefile format. This data set was also downloaded from Data.gov.sg.
Note*:
- There were some issues with a few of the datasets.
- For Airbnb listings as at June 2019, this dataset is provided by Prof. Kam as the dataset no longer exist on Inside Airbnb website.
- For Hotels, Tourist attractions and MRT services,
- The original plan was to download it from SLA OneMap Service by using onemapsgapi.
- However, as there were some issues with creating a onemap account due to the maintenance issue, Prof. Kam has kindly provided us with the data in CSV format.
3. Install and Load R packages
This code chunk performs 3 tasks:
- A list called packages will be created and will consists of all the R packages required to accomplish this hands-on exercise.
- Check if R packages on package have been installed in R and if not, they will be installed.
- After all the R packages have been installed, they will be loaded.
packages <- c('sf', 'tidyverse', 'tmap', 'rgdal', 'maptools', 'raster','spatstat', 'kableExtra', 'devtools')
for(p in packages){
if(!require(p, character.only = T)){
install.packages(p)
}
library(p, character.only = T)
}
devtools::install_github("gadenbuie/xaringanExtra")
library(xaringanExtra)
More on the packages used:
sf: used for importing, managing, and processing geospatial data
- specifically vector-based geospatial data
tidyverse: used for importing, wrangling and visualising data. It consists of a family of R packages, such as:
- readr for importing csv data,
- readxl for importing Excel worksheet,
- tidyr for manipulating data,
- dplyr for transforming data, and
- ggplot2 for visualising data
tmap: provides functions for plotting cartographic quality static point patterns maps or interactive maps by using leaflet API.
maptools: Manipulate geographic data.
- Will be used to convert spatial objects into ppp format of spatstat
raster: reads, writes, manipulates, analyses and model of gridded spatial data
spatstat: contains useful functions for point pattern analysis.
- Will be used to perform 1st- and 2nd-order spatial point patterns analysis and derive kernel density estimation (KDE) layer.
gridExtra: For arranging KDE maps when plotted into grid object
kableExtra: Designed to extend the basic functionality of tables. In this exercise, I will be using it to construct complex tables and customize styles using a readable syntax.
devtools: used for installing any R packages which is not available in RCRAN. In this exercise, I will be installing using devtools to install the package xaringanExtra which is still under development stage.
xaringanExtra: is an enhancement of xaringan package. As it is still under development stage, we can still install the current version using install_github function of devtools. This package will be used to add Panelsets to contain both the r code chunk and results whereever applicable.
4. Importing and Wrangling Geospatial Data
- The code chunks below uses st_read() of sf package to import the geospatial data sets
- The imported shapefile will be a simple feature object of sf.
- st_crs() of sf package will also be used to check the CRS.
4.1 Import and check CRS
4.1.1 CoastalOutline
Code Chunk
sg_sf <- st_read(dsn = "data/geospatial", layer="CostalOutline")
Reading layer `CostalOutline' from data source
`C:\aisyahajit2018\IS415\IS415_blog\_posts\2021-10-01-take-home-exercise-2\data\geospatial'
using driver `ESRI Shapefile'
Simple feature collection with 60 features and 4 fields
Geometry type: POLYGON
Dimension: XY
Bounding box: xmin: 2663.926 ymin: 16357.98 xmax: 56047.79 ymax: 50244.03
Projected CRS: SVY21st_crs
st_crs(sg_sf)
Coordinate Reference System:
User input: SVY21
wkt:
PROJCRS["SVY21",
BASEGEOGCRS["SVY21[WGS84]",
DATUM["World Geodetic System 1984",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]],
ID["EPSG",6326]],
PRIMEM["Greenwich",0,
ANGLEUNIT["Degree",0.0174532925199433]]],
CONVERSION["unnamed",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",1.36666666666667,
ANGLEUNIT["Degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",103.833333333333,
ANGLEUNIT["Degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",1,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",28001.642,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",38744.572,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["(E)",east,
ORDER[1],
LENGTHUNIT["metre",1,
ID["EPSG",9001]]],
AXIS["(N)",north,
ORDER[2],
LENGTHUNIT["metre",1,
ID["EPSG",9001]]]]From the results above, we can see that:
- sg has 60 features and 4 fields.
- The geometry type is polygon thus, when we convert it to a Spatial* class object later, we should keep in mind that the suitable class would be SpatialPolygons.
- The projected CRS for sg is SVY21.
- However, the EPSG code shown is 9001 which is wrong since the correct EPSG code for SVY21 should be 3414.
4.1.2 MPSZ
Code Chunk
mpsz_sf <- st_read(dsn = "data/geospatial", layer="MP14_SUBZONE_WEB_PL")
Reading layer `MP14_SUBZONE_WEB_PL' from data source
`C:\aisyahajit2018\IS415\IS415_blog\_posts\2021-10-01-take-home-exercise-2\data\geospatial'
using driver `ESRI Shapefile'
Simple feature collection with 323 features and 15 fields
Geometry type: MULTIPOLYGON
Dimension: XY
Bounding box: xmin: 2667.538 ymin: 15748.72 xmax: 56396.44 ymax: 50256.33
Projected CRS: SVY21st_crs
st_crs(mpsz_sf)
Coordinate Reference System:
User input: SVY21
wkt:
PROJCRS["SVY21",
BASEGEOGCRS["SVY21[WGS84]",
DATUM["World Geodetic System 1984",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]],
ID["EPSG",6326]],
PRIMEM["Greenwich",0,
ANGLEUNIT["Degree",0.0174532925199433]]],
CONVERSION["unnamed",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",1.36666666666667,
ANGLEUNIT["Degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",103.833333333333,
ANGLEUNIT["Degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",1,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",28001.642,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",38744.572,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["(E)",east,
ORDER[1],
LENGTHUNIT["metre",1,
ID["EPSG",9001]]],
AXIS["(N)",north,
ORDER[2],
LENGTHUNIT["metre",1,
ID["EPSG",9001]]]]From the results above, we can see that:
- mpsz has 323 features and 15 fields.
- The geometry type is multipolygon. Similarly, when we convert it to a Spatial* class object later, we should keep in mind that the suitable class would be SpatialPolygons.
- The projected CRS for mpsz is also SVY21 but the EPSG code shown is 9001 which is wrong since the correct EPSG code for SVY21 should be 3414.
4.1.3 MRT
Code Chunk
mrtlrt_sf <- st_read(dsn = "data/geospatial", layer="MRTLRTStnPtt")
Reading layer `MRTLRTStnPtt' from data source
`C:\aisyahajit2018\IS415\IS415_blog\_posts\2021-10-01-take-home-exercise-2\data\geospatial'
using driver `ESRI Shapefile'
Simple feature collection with 185 features and 3 fields
Geometry type: POINT
Dimension: XY
Bounding box: xmin: 6138.311 ymin: 27555.06 xmax: 45254.86 ymax: 47854.2
Projected CRS: SVY21st_crs
st_crs(mrtlrt_sf)
Coordinate Reference System:
User input: SVY21
wkt:
PROJCRS["SVY21",
BASEGEOGCRS["SVY21[WGS84]",
DATUM["World Geodetic System 1984",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]],
ID["EPSG",6326]],
PRIMEM["Greenwich",0,
ANGLEUNIT["Degree",0.0174532925199433]]],
CONVERSION["unnamed",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",1.36666666666667,
ANGLEUNIT["Degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",103.833333333333,
ANGLEUNIT["Degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",1,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",28001.642,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",38744.572,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["(E)",east,
ORDER[1],
LENGTHUNIT["metre",1,
ID["EPSG",9001]]],
AXIS["(N)",north,
ORDER[2],
LENGTHUNIT["metre",1,
ID["EPSG",9001]]]]From the results above, we can see that:
- mrtlrt has 185 features and 3 fields.
- The geometry type is point. So, Swhen we convert it to a Spatial* class object later, we should keep in mind that the suitable class would be SpatialPoints.
- Similar to sg and mpsz, the projected CRS for mrtlrt is also SVY21 but the EPSG code shown wrong.
Since all the projected CRS for the above sf dataframes is SVY21, we need to assign them the correct EPSG code which is 3414 in the next section.
4.2 Assign EPSG code to sf dataframes and check
4.2.1 Assign EPSG code to sg_sf
Code Chunk
sg_sf <- st_set_crs(sg_sf, 3414)
st_crs
st_crs(sg_sf)
Coordinate Reference System:
User input: EPSG:3414
wkt:
PROJCRS["SVY21 / Singapore TM",
BASEGEOGCRS["SVY21",
DATUM["SVY21",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4757]],
CONVERSION["Singapore Transverse Mercator",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",1.36666666666667,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",103.833333333333,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",1,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",28001.642,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",38744.572,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["northing (N)",north,
ORDER[1],
LENGTHUNIT["metre",1]],
AXIS["easting (E)",east,
ORDER[2],
LENGTHUNIT["metre",1]],
USAGE[
SCOPE["Cadastre, engineering survey, topographic mapping."],
AREA["Singapore - onshore and offshore."],
BBOX[1.13,103.59,1.47,104.07]],
ID["EPSG",3414]]4.2.2 Assign EPSG code to mpsz_sf
Code Chunk
mpsz_sf <- st_set_crs(mpsz_sf, 3414)
st_crs
st_crs(mpsz_sf)
Coordinate Reference System:
User input: EPSG:3414
wkt:
PROJCRS["SVY21 / Singapore TM",
BASEGEOGCRS["SVY21",
DATUM["SVY21",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4757]],
CONVERSION["Singapore Transverse Mercator",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",1.36666666666667,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",103.833333333333,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",1,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",28001.642,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",38744.572,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["northing (N)",north,
ORDER[1],
LENGTHUNIT["metre",1]],
AXIS["easting (E)",east,
ORDER[2],
LENGTHUNIT["metre",1]],
USAGE[
SCOPE["Cadastre, engineering survey, topographic mapping."],
AREA["Singapore - onshore and offshore."],
BBOX[1.13,103.59,1.47,104.07]],
ID["EPSG",3414]]4.2.3 Assign EPSG code to mrtlrt_sf
Code Chunk
mrtlrt_sf <- st_set_crs(mrtlrt_sf, 3414)
st_crs
st_crs(mrtlrt_sf)
Coordinate Reference System:
User input: EPSG:3414
wkt:
PROJCRS["SVY21 / Singapore TM",
BASEGEOGCRS["SVY21",
DATUM["SVY21",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4757]],
CONVERSION["Singapore Transverse Mercator",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",1.36666666666667,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",103.833333333333,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",1,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",28001.642,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",38744.572,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["northing (N)",north,
ORDER[1],
LENGTHUNIT["metre",1]],
AXIS["easting (E)",east,
ORDER[2],
LENGTHUNIT["metre",1]],
USAGE[
SCOPE["Cadastre, engineering survey, topographic mapping."],
AREA["Singapore - onshore and offshore."],
BBOX[1.13,103.59,1.47,104.07]],
ID["EPSG",3414]]The CRS has now been correctly assigned for all 3 geospatial data.
4.3 Check if geometries are valid
- In the following code chunks, st_is_valid() is used to check for any invalid geometries.
- We have to check this so that we make sure that each polygon doesn’t break any topological rules and ideally, all of the geometries should be valid.
[1] 1[1] 9[1] 0From the results above, there are:
- 1 invalid geometries for sg_sf,
- 9 invalid geometries for mpsz_sf,
- 0 invalid geometries for mrtlrt_sf
4.4 Handle the invalid geometries
- The following code chunk uses st_make_valid() function to make the geometries valid for sg_sf and mpsz_sf
sg_sf <- st_make_valid(sg_sf)
mpsz_sf <- st_make_valid(mpsz_sf)
Based on the results above, there are no longer any invalid geometries.
4.5 Mapping the geospatial layers
tm_shape(sg_sf) +
tm_polygons() +
tm_shape(mpsz_sf) +
tm_borders(alpha = 0.5) +
tm_shape(mrtlrt_sf) +
tm_dots(col="red", size=0.1)
- From the above map, we can see how the sg and mpsz differs.
- The CoastalOutline, sg, contains polygons outlining the coast of mainland and offshore islands in Singapore, indicated by the border.
- Whereas, the mpsz shows polygons of a subzone boundary, indicated by the polygons shaded in Grey.
- Our MrtLrt services are also widely scattered across Singapore, with a relatively higher number of services found in the Central, West and Northeast region.
5. Importing and Wrangling Aspatial Data
5.1 Import and inspect aspatial data
- Here, we will be importing our aspatial data using read_csv() of readr package
- Then, we will take a look at our aspatial data using glimpse() of dplyr package
5.1.1 Airbnb listings as at June 2019
Code Chunk
abb_jun19 <- read_csv("data/aspatial/30062019.csv")
glimpse
glimpse(abb_jun19)
Rows: 8,293
Columns: 16
$ id <dbl> 49091, 50646, 56334, 71609, 7~
$ name <chr> "COZICOMFORT LONG TERM STAY R~
$ host_id <dbl> 266763, 227796, 266763, 36704~
$ host_name <chr> "Francesca", "Sujatha", "Fran~
$ neighbourhood_group <chr> "North Region", "Central Regi~
$ neighbourhood <chr> "Woodlands", "Bukit Timah", "~
$ latitude <dbl> 1.44255, 1.33235, 1.44246, 1.~
$ longitude <dbl> 103.7958, 103.7852, 103.7967,~
$ room_type <chr> "Private room", "Private room~
$ price <dbl> 81, 80, 68, 200, 92, 102, 203~
$ minimum_nights <dbl> 180, 90, 6, 1, 1, 1, 1, 7, 30~
$ number_of_reviews <dbl> 1, 18, 20, 12, 20, 35, 23, 15~
$ last_review <date> 2013-10-21, 2014-12-26, 2015~
$ reviews_per_month <dbl> 0.01, 0.28, 0.21, 0.13, 0.21,~
$ calculated_host_listings_count <dbl> 2, 1, 2, 9, 9, 9, 9, 1, 4, 4,~
$ availability_365 <dbl> 365, 365, 365, 353, 353, 348,~5.1.2 Airbnb listings as at June 2021
Code Chunk
abb_jun21 <- read_csv("data/aspatial/listings.csv")
glimpse
glimpse(abb_jun21)
Rows: 4,238
Columns: 16
$ id <dbl> 49091, 50646, 56334, 71609, 7~
$ name <chr> "COZICOMFORT LONG TERM STAY R~
$ host_id <dbl> 266763, 227796, 266763, 36704~
$ host_name <chr> "Francesca", "Sujatha", "Fran~
$ neighbourhood_group <chr> "North Region", "Central Regi~
$ neighbourhood <chr> "Woodlands", "Bukit Timah", "~
$ latitude <dbl> 1.44129, 1.33432, 1.44041, 1.~
$ longitude <dbl> 103.7957, 103.7852, 103.7967,~
$ room_type <chr> "Private room", "Private room~
$ price <dbl> 81, 80, 67, 177, 81, 81, 52, ~
$ minimum_nights <dbl> 180, 90, 6, 90, 90, 90, 14, 1~
$ number_of_reviews <dbl> 1, 18, 20, 20, 24, 48, 20, 13~
$ last_review <date> 2013-10-21, 2014-07-08, 2015~
$ reviews_per_month <dbl> 0.01, 0.22, 0.16, 0.29, 0.34,~
$ calculated_host_listings_count <dbl> 2, 1, 2, 4, 4, 4, 50, 50, 7, ~
$ availability_365 <dbl> 365, 365, 365, 365, 365, 365,~5.1.3 Hotels
Code Chunk
hotels <- read.csv("data/aspatial/hotels.csv")
glimpse
glimpse(hotels)
Rows: 422
Columns: 9
$ NAME <chr> "Jayleen Clarke Quay Hotel", "JEN Singapor~
$ ADDRESSPOSTALCODE <int> 59390, 238858, 249716, 399041, 238485, 399~
$ ADDRESSSTREETNAME <chr> "25 New Bridge Road", "277 Orchard Road , ~
$ HYPERLINK <chr> "jayleenclarkequay@gmail.com", "singaporeo~
$ TOTALROOMS <int> 20, 499, 565, 42, 81, 33, 17, 634, 56, 451~
$ KEEPERNAME <chr> "James Federick Chong Kah Yean", "Kuok Oon~
$ Lat <dbl> 1.288715, 1.300510, 1.304271, 1.311586, 1.~
$ Lng <dbl> 103.8475, 103.8392, 103.8239, 103.8775, 10~
$ ICON_NAME <chr> "hotel.gif", "hotel.gif", "hotel.gif", "ho~5.1.3 Tourism
Code Chunk
tourism <- read.csv("data/aspatial/tourism.csv")
glimpse
glimpse(tourism)
Rows: 107
Columns: 17
$ NAME <chr> "Chinatown Heritage Centre, Singapore", "T~
$ DESCRIPTION <chr> "Experience how Singaporeâ\200\231s early ~
$ ADDRESSSTREETNAME <chr> "48 Pagoda Street", "158 Telok Ayer Street~
$ HYPERLINK <chr> "http://www.singaporechinatown.com.sg/", "~
$ PHOTOURL <chr> "www.yoursingapore.com/content/dam/desktop~
$ URL_PATH <chr> "www.yoursingapore.com/en/see-do-singapore~
$ IMAGE_ALT_TEXT <chr> "Learn more about local Chinese culture at~
$ PHOTOCREDITS <chr> "Joel Chua DY", "Joel Chua DY", "Joel Chua~
$ LASTMODIFIED <chr> "2015-11-02T10:16:52.847+08:00", "2015-11-~
$ LATITUDE <dbl> 1.283510, 1.280940, 1.310070, 1.277219, 1.~
$ LONGTITUDE <dbl> 103.8444, 103.8476, 103.8994, 103.8373, 10~
$ META_DESCRIPTION <chr> "At the Chinatown Heritage Centre, experie~
$ OPENING_HOURS <chr> "Daily, 9am â\200“ 8pm,Last entry at 7pm.*~
$ Lat <dbl> 1.283510, 1.280940, 1.310070, 1.277219, 1.~
$ Lng <dbl> 103.8444, 103.8476, 103.8994, 103.8373, 10~
$ ICON_NAME <chr> "tourist_spot.gif", "tourist_spot.gif", "t~
$ ADDRESSPOSTALCODE <int> NA, NA, NA, 0, NA, NA, NA, NA, NA, NA, NA,~From the above results, we can see that for:
- abb_jun19:
- The data contains 8,293 rows and 16 columns.
- The coordinates are in
latitudeandlongitudecolumns.
- abb_jun21:
- The data contains 4,238 rows and 16 columns.
- The coordinates are also in
latitudeandlongitudecolumns.
- hotels:
- The data contains 422 rows and 9 columns.
- The coordinates columns are in
LatandLng.
- tourism:
- The data contains 107 rows and 17 columns.
- Only tourism in particular has 2 geographic coordinates which are:
LONGTITUDEandLATITUDELatandLng
- Hence, since the coordinate columns are Latitude & Longtitude / Lat & Lng, their projected CRS will be WGS84.
- In the later section, we will need to:
- Assign them the respective EPSG code which is 4326
- Transform it to fit with our geospatial data, which is using SVY21 with EPSG code 3414.
5.2 Check for NA values
NOTE:
- When checking for NA values across the whole data set, we should be careful and not be too quick to omit them.
- It may give NA values under fields we will not be using, such as
last_review.- This might be because the listing is new and does not have reviews, hence there is no date for last_review.
- Since we are more concerned with variables such as
latitude,longitudeand particularlyroom_typefor Section B, we shall not remove the NA values if there are no NA values for these specific columns.
5.2.1 Check for NA values in abb_jun19
id name
0 2
host_id host_name
0 82
neighbourhood_group neighbourhood
0 0
latitude longitude
0 0
room_type price
0 0
minimum_nights number_of_reviews
0 0
last_review reviews_per_month
3137 3137
calculated_host_listings_count availability_365
0 0 From the above results, we can see that there are:
- 3137 rows under the columns
last_reviewandreviews_per_monthwith NA values. - 82 rows under
host_namecolumn with NA values. - 2 rows under
namecolumn with NA values. - As there are no NA values for the columns we are interested in, we shall not remove them.
5.2.2 Check for NA values in abb_jun21
id name
0 0
host_id host_name
0 8
neighbourhood_group neighbourhood
0 0
latitude longitude
0 0
room_type price
0 0
minimum_nights number_of_reviews
0 0
last_review reviews_per_month
1759 1759
calculated_host_listings_count availability_365
0 0 From the above results, we can see that there are:
- 1759 rows under the columns
last_reviewandreviews_per_monthwith NA values. - 8 rows under
host_namecolumn with NA values. - Similar to June 2019’s listings, since there are no NA values for the columns we are interested in, we shall not remove them.
5.2.3 Check for NA values in hotels
NAME ADDRESSPOSTALCODE ADDRESSSTREETNAME
0 0 0
HYPERLINK TOTALROOMS KEEPERNAME
129 0 0
Lat Lng ICON_NAME
0 0 0 - We can see from the results above that hotels does not contain any NA values.
5.2.4 Check for NA values in tourism
NAME DESCRIPTION ADDRESSSTREETNAME
0 0 2
HYPERLINK PHOTOURL URL_PATH
15 12 0
IMAGE_ALT_TEXT PHOTOCREDITS LASTMODIFIED
17 62 0
LATITUDE LONGTITUDE META_DESCRIPTION
1 1 0
OPENING_HOURS Lat Lng
29 0 0
ICON_NAME ADDRESSPOSTALCODE
0 96 - We can see from the results above that the columns,
LONGTITUDEandLATITUDE, contains NA values but theLatandLngcolumns does not contain any NA values. - We might also think that since
LatandLngcolumns does not contain any NA values, we should use these columns as our coordinates instead. - However, we should not be too quick to decide on this and instead, look at the data closely.
5.2.4.1 Inspect NA values in tourism coordinates
Code Chunk
tourism_na <- tourism %>%
filter_at(vars(LONGTITUDE, LATITUDE), any_vars(is.na(.)))
Head
kable(head(tourism_na)) %>%
kable_styling(latex_options="striped", full_width = F) %>%
scroll_box(width = "100%", height = "400px")
| NAME | DESCRIPTION | ADDRESSSTREETNAME | HYPERLINK | PHOTOURL | URL_PATH | IMAGE_ALT_TEXT | PHOTOCREDITS | LASTMODIFIED | LATITUDE | LONGTITUDE | META_DESCRIPTION | OPENING_HOURS | Lat | Lng | ICON_NAME | ADDRESSPOSTALCODE |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Cruises from Singapore | Watch the Singapore skyline disappear from view as you sail away into the horizon. As the gateway to Southeast Asia and the world, Singapore is the choice homeport of numerous cruise lines. | NA | NA | www.yoursingapore.com/content/dam/desktop/global/see-do-singapore/beyond-singapore/beyond-singapore-carousel01-rect.jpg | www.yoursingapore.com/en/see-do-singapore/beyond-singapore.html | NA | NA | 2015-11-03T14:54:39.639+08:00 | NA | NA | Take these luxurious cruise liners from Singapore and visit lush tropical rainforests, colourful open air markets, mystifying ancient cities and white sand beaches. | NA | 0 | 0 | tourist_spot.gif | 0 |
- From the results above, we can see that the tourism activity that has NA values for the
LONGTITUDEandLATITUDEcolumns is Cruises from Singapore. - It does not belong to a specific location on the map since it might be outside of Singapore’s geographical boundaries.
- Additionally, although the values for the
LatandLngcolumns are 0, having 0,0 coordinates can mean different things. - Thus, we should remove this row from our dataframe.
5.2.4.2 Remove NA values
tourism <- tourism %>%
filter_at(vars(LONGTITUDE, LATITUDE), any_vars(!is.na(.)))
5.2.4.3 Check if NA value is removed
tourism_na <- tourism %>%
filter_at(vars(LONGTITUDE, LATITUDE), any_vars(is.na(.)))
nrow(tourism_na)
[1] 0We can see from the above results that the row with NA value in the coordinates columns has been removed.
5.3 Convert aspatial data frame into sf objects, assign and transform CRS
- As mentioned earlier, the coordinate columns are in decimal degrees which means that they are using WGS84.
- We will need to assign them the respective epsg code 4326 before before transforming it to 3414 - the epsg code for SVY21.
abb_jun19_sf <- st_as_sf(abb_jun19,
coords = c("longitude",
"latitude"),
crs=4326) %>%
st_transform(crs = 3414)
abb_jun21_sf <- st_as_sf(abb_jun21,
coords = c("longitude",
"latitude"),
crs=4326) %>%
st_transform(crs = 3414)
hotels_sf <- st_as_sf(hotels,
coords = c("Lng",
"Lat"),
crs=4326) %>%
st_transform(crs = 3414)
tourism_sf <- st_as_sf(tourism,
coords = c("Lng",
"Lat"),
crs=4326) %>%
st_transform(crs = 3414)
5.4. Check the respective CRS
abb_jun19_sf crs
st_crs(abb_jun19_sf)
Coordinate Reference System:
User input: EPSG:3414
wkt:
PROJCRS["SVY21 / Singapore TM",
BASEGEOGCRS["SVY21",
DATUM["SVY21",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4757]],
CONVERSION["Singapore Transverse Mercator",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",1.36666666666667,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",103.833333333333,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",1,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",28001.642,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",38744.572,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["northing (N)",north,
ORDER[1],
LENGTHUNIT["metre",1]],
AXIS["easting (E)",east,
ORDER[2],
LENGTHUNIT["metre",1]],
USAGE[
SCOPE["Cadastre, engineering survey, topographic mapping."],
AREA["Singapore - onshore and offshore."],
BBOX[1.13,103.59,1.47,104.07]],
ID["EPSG",3414]]abb_jun21_sf crs
st_crs(abb_jun21_sf)
Coordinate Reference System:
User input: EPSG:3414
wkt:
PROJCRS["SVY21 / Singapore TM",
BASEGEOGCRS["SVY21",
DATUM["SVY21",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4757]],
CONVERSION["Singapore Transverse Mercator",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",1.36666666666667,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",103.833333333333,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",1,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",28001.642,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",38744.572,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["northing (N)",north,
ORDER[1],
LENGTHUNIT["metre",1]],
AXIS["easting (E)",east,
ORDER[2],
LENGTHUNIT["metre",1]],
USAGE[
SCOPE["Cadastre, engineering survey, topographic mapping."],
AREA["Singapore - onshore and offshore."],
BBOX[1.13,103.59,1.47,104.07]],
ID["EPSG",3414]]hotels_sf crs
st_crs(hotels_sf)
Coordinate Reference System:
User input: EPSG:3414
wkt:
PROJCRS["SVY21 / Singapore TM",
BASEGEOGCRS["SVY21",
DATUM["SVY21",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4757]],
CONVERSION["Singapore Transverse Mercator",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",1.36666666666667,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",103.833333333333,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",1,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",28001.642,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",38744.572,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["northing (N)",north,
ORDER[1],
LENGTHUNIT["metre",1]],
AXIS["easting (E)",east,
ORDER[2],
LENGTHUNIT["metre",1]],
USAGE[
SCOPE["Cadastre, engineering survey, topographic mapping."],
AREA["Singapore - onshore and offshore."],
BBOX[1.13,103.59,1.47,104.07]],
ID["EPSG",3414]]tourism_sf crs
st_crs(tourism_sf)
Coordinate Reference System:
User input: EPSG:3414
wkt:
PROJCRS["SVY21 / Singapore TM",
BASEGEOGCRS["SVY21",
DATUM["SVY21",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4757]],
CONVERSION["Singapore Transverse Mercator",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",1.36666666666667,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",103.833333333333,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",1,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",28001.642,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",38744.572,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["northing (N)",north,
ORDER[1],
LENGTHUNIT["metre",1]],
AXIS["easting (E)",east,
ORDER[2],
LENGTHUNIT["metre",1]],
USAGE[
SCOPE["Cadastre, engineering survey, topographic mapping."],
AREA["Singapore - onshore and offshore."],
BBOX[1.13,103.59,1.47,104.07]],
ID["EPSG",3414]]5.5. Plotting it altogether
In the code chunks below, mapping functions of tmap package is used.
5.5.1 Plot map without abb_jun21
- Here, we will plot SG Airbnb Jun 2019 listings, hotels, tourist attractions and MRT services to ensure that they are being assigned the right CRS.
tm_shape(sg_sf) +
tm_polygons() +
tm_shape(mpsz_sf) +
tm_borders(alpha = 0.5) +
tmap_options(check.and.fix = TRUE) +
tm_shape(mrtlrt_sf) +
tm_dots(col ="red", size = 0.02) +
tm_shape(abb_jun19_sf) +
tm_dots(col ="blue", size = 0.02) +
tm_shape(hotels_sf) +
tm_dots(col ="orange", size = 0.02) +
tm_shape(tourism_sf) +
tm_dots(col ="purple", size = 0.02) +
tm_layout(title= 'Map of Airbnb 2019 listings with other location factors',
title.position = c('right', 'bottom'))
5.5.2 Plot map without abb_jun19
- Here, we will plot SG Airbnb Jun 2021 listings, hotels, tourist attractions and MRT services to ensure that they are being assigned the right CRS.
tm_shape(sg_sf) +
tm_polygons() +
tm_shape(mpsz_sf) +
tm_borders(alpha = 0.5) +
tmap_options(check.and.fix = TRUE) +
tm_shape(mrtlrt_sf) +
tm_dots(col ="red", size = 0.02) +
tm_shape(abb_jun21_sf) +
tm_dots(col ="green", size = 0.02) +
tm_shape(hotels_sf) +
tm_dots(col ="orange", size = 0.02) +
tm_shape(tourism_sf) +
tm_dots(col ="purple", size = 0.02) +
tm_layout(title= 'Map of Airbnb 2021 listings with other location factors',
title.position = c('right', 'bottom'))
6. Spatial Data Wrangling
In this section, we will be converting our simple feature data frames to ppp objects for our analysis later.
6.1 Convert from sf dataframes to sp’s Spatial* Classes
- The first step is to convert from sf dataframes to sp’s Spatial* Classes using as_Spatial() of sf.
abb_jun19 <- as_Spatial(abb_jun19_sf)
abb_jun21 <- as_Spatial(abb_jun21_sf)
hotels <- as_Spatial(hotels_sf)
tourism <- as_Spatial(tourism_sf)
sg <- as_Spatial(sg_sf)
mpsz <- as_Spatial(mpsz_sf)
mrtlrt <- as_Spatial(mrtlrt_sf)
6.2 Convert from Spatial* class to generic sp format
- The spatstat package requires the analytical data in ppp object form
- As there is no direct way to convert a Spatial* class to ppp object, we have to convert the Spatial classes* into Spatial object first.
- The function as() here will force an object to belong to a class.
- Here, we will use different arguments like SpatialPoints and SpatialPolysons for the given input Class.
- Note: After this step, the Data dataframe will be dropped.
abb_jun19_sp <- as(abb_jun19, "SpatialPoints")
abb_jun21_sp <- as(abb_jun21, "SpatialPoints")
hotels_sp <- as(hotels, "SpatialPoints")
tourism_sp <- as(tourism, "SpatialPoints")
mrtlrt_sp <- as(mrtlrt, "SpatialPoints")
sg_sp <- as(sg, "SpatialPolygons")
mpsz_sp <- as(mpsz, "SpatialPolygons")
6.3 Convert from sp to spatstat ppp format
- Now, we will use the as() function again but with “ppp” as the argument for Class to convert the spatial data into spatstat’s ppp object format
- Note:
- In ppp objects, the projection information and bbox no longer exist.
- There is also a data called window in the ppp object which we can use to see our study area.
- Later on, we will need to define the study area for our owin object.
abb_jun19_ppp <- as(abb_jun19_sp, "ppp")
abb_jun21_ppp <- as(abb_jun21_sp, "ppp")
hotels_ppp <- as(hotels_sp, "ppp")
tourism_ppp <- as(tourism_sp, "ppp")
mrtlrt_ppp <- as(mrtlrt_sp, "ppp")
6.4 Summary statistics of newly created ppp objects
- Here, we can take a quick look at the summary statistics of the newly created ppp object by using the code chunk below.
abb_jun19_ppp summary
summary(abb_jun19_ppp)
Planar point pattern: 8293 points
Average intensity 9.345289e-06 points per square unit
*Pattern contains duplicated points*
Coordinates are given to 3 decimal places
i.e. rounded to the nearest multiple of 0.001 units
Window: rectangle = [7215.57, 44098.31] x [25166.35, 49226.35] units
(36880 x 24060 units)
Window area = 887399000 square unitsabb_jun21_ppp summary
summary(abb_jun21_ppp)
Planar point pattern: 4238 points
Average intensity 5.114513e-06 points per square unit
*Pattern contains duplicated points*
Coordinates are given to 3 decimal places
i.e. rounded to the nearest multiple of 0.001 units
Window: rectangle = [7406.99, 43337.89] x [25330, 48391.55] units
(35930 x 23060 units)
Window area = 828622000 square unitshotels_ppp summary
summary(hotels_ppp)
Planar point pattern: 422 points
Average intensity 5.58414e-07 points per square unit
*Pattern contains duplicated points*
Coordinates are given to 3 decimal places
i.e. rounded to the nearest multiple of 0.001 units
Window: rectangle = [5939.24, 45334.18] x [25379.44, 44562.4] units
(39390 x 19180 units)
Window area = 755712000 square unitstourism_ppp summary
summary(tourism_ppp)
Planar point pattern: 106 points
Average intensity 1.328016e-07 points per square unit
*Pattern contains duplicated points*
Coordinates are given to 3 decimal places
i.e. rounded to the nearest multiple of 0.001 units
Window: rectangle = [11380.23, 43659.54] x [22869.34, 47596.73] units
(32280 x 24730 units)
Window area = 798183000 square unitsmrtlrt_ppp summary
summary(mrtlrt_ppp)
Planar point pattern: 185 points
Average intensity 2.32988e-07 points per square unit
Coordinates are given to 3 decimal places
i.e. rounded to the nearest multiple of 0.001 units
Window: rectangle = [6138.31, 45254.86] x [27555.06, 47854.2] units
(39120 x 20300 units)
Window area = 794032000 square unitsFrom the results above:
- We can see that there are some warning message about duplicates. The ppp objects that contains duplicates are:
- Airbnb Jun 2019 listings,
- Airbnb Jun 2021 listings,
- hotels
- tourism
- In spatial point patterns analysis an issue of significant is the presence of duplicates.
- The statistical methodology used for spatial point patterns processes is based largely on the assumption that process are simple, that is, that the points cannot be coincident.
6.4 Handling duplicates
6.4.1 Check for number of duplicated points
- The function multiplicity() helps us count the number of duplicates for each point in a spatial point pattern.
- If we want to know how many locations have more than one point event, we can use the code chunk below:
sum(multiplicity(abb_jun19_ppp) > 1)
[1] 6sum(multiplicity(abb_jun21_ppp) > 1)
[1] 224sum(multiplicity(hotels_ppp) > 1)
[1] 10sum(multiplicity(tourism_ppp) > 1)
[1] 7sum(multiplicity(mrtlrt_ppp) > 1)
[1] 06.4.2 Handling duplicates with jittering
- Out of the 3 known ways to overcome duplicates, we will use the jittering approach.
- This approach will add a small perturbation to the duplicate points so that they do not occupy the exact same space.
- The code chunk below implements the jittering approach.
- rjitter() of spatstat package is used
abb_jun19_ppp_jit <- rjitter(abb_jun19_ppp,
retry=TRUE,
nsim=1,
drop=TRUE)
abb_jun21_ppp_jit <- rjitter(abb_jun21_ppp,
retry=TRUE,
nsim=1,
drop=TRUE)
hotels_ppp_jit <- rjitter(hotels_ppp,
retry=TRUE,
nsim=1,
drop=TRUE)
tourism_ppp_jit <- rjitter(tourism_ppp,
retry=TRUE,
nsim=1,
drop=TRUE)
6.4.3 Check if duplicates still exist
- Here, instead of checking the summary, we can use the functions any() and duplicated() to ensure that there are no duplicated points after the jittering approach
- If the function returns FALSE, it will help us confirm that there are indeed no duplicated points present.
any(duplicated(abb_jun19_ppp_jit))
[1] FALSEany(duplicated(abb_jun21_ppp_jit))
[1] FALSEany(duplicated(hotels_ppp_jit))
[1] FALSEany(duplicated(tourism_ppp_jit))
[1] FALSEFrom the above results, we can confirm that there are no longer duplicates after performing the jittering approach.
6.5 Create owin objects
- When analysing spatial point patterns, it is a good practice to confine the spatial point patterns analysis with a geographical area such as the Singapore boundary.
- The presence of water catchment areas may skew the calculations, so we need to limit by planning subzone area later on.
- In spatstat package, this polygonal region is known as owin, an object specially designed to represent this polygonal region.
- We can create this object by converting sg SpatialPolygon object into owin object of spatstat
- The ouput object can be displayed by using plot() function
sg_owin <- as(sg_sp, "owin")
plot(sg_owin)
6.6 Combine point events object and owin object
The code chunk below combines the point objects of abb_jun19, abb_jun21, hotels, tourism and mrtlrt with the owin object into one ppp object class.
abb_jun19SG_ppp = abb_jun19_ppp_jit[sg_owin]
abb_jun21SG_ppp = abb_jun21_ppp_jit[sg_owin]
hotelsSG_ppp = hotels_ppp_jit[sg_owin]
tourismSG_ppp = tourism_ppp_jit[sg_owin]
mrtlrtSG_ppp = mrtlrt_ppp[sg_owin]
6.6.1 Plot ppp object within Singapore boundary
The code chunk below plots the ppp object class using plot function of base R package.
par(mfrow=c(3,2), mar=c(1,1,1,1))
plot(abb_jun19SG_ppp)
plot(abb_jun21SG_ppp)
plot(hotelsSG_ppp)
plot(tourismSG_ppp)
plot(mrtlrtSG_ppp)
6.7 Rescale to kilometers
- As we will be computing the kernel density estimation (KDE) later on, we need to convert our point pattern data into kilometres to be able to obtain any valuable insights.
- The default unit of measurement of SVY21 is in metres.
- Hence, we will use rescale() of spatstat package which helps to converts a point pattern dataset to another unit of length.
- The arguments for the rescale(X, s, unitname) function are:
- X: our point pattern data (e.g. abb_jun19SG_ppp, hotelsSG_ppp)
- s: the new unit
- unitname: the name of new unit (e.g. inches or in this case, we specifically using km)
abb_jun19SG_ppp.km <- rescale(abb_jun19SG_ppp, 1000, "km")
abb_jun21SG_ppp.km <- rescale(abb_jun21SG_ppp, 1000, "km")
hotelsSG_ppp.km <- rescale(hotelsSG_ppp, 1000, "km")
tourismSG_ppp.km <- rescale(tourismSG_ppp, 1000, "km")
mrtlrtSG_ppp.km <- rescale(mrtlrtSG_ppp, 1000, "km")
7. Section A: Airbnb Distribution in 2019
This section aims to investigate if the distribution of Airbnb listings at June 2019 are affected by location factors such as near to existing hotels, MRT services and tourist attractions. For our investigation, this section will consist of 2 parts:
- Exploratory Spatial Data Analysis
- Where we will need to derive the Kernel Density Estimation maps with the objective of understanding and potentially predicting event patterns
- Second-Order Spatial Point Patterns Analysis
- Where we perform this analysis to understand when data points start to form clusters or when data points have spatial patterns that are statistically significant.
7.1 Exploratory Spatial Data Analysis
7.1.1 KDE maps
- density() of spatstat package is used to plot the re-scaled data set and plot the KDE map
- The arguments of this function includes:
- Bandwidth Selection
- sigma=bw.diggle helps to compute automatically.
- Other possible methods are “bw.CvL()”, etc. They are standalone functions, so we can extract them by calling this function externally with the appropriate arguments.
- Kernel Type - the smoothing method.
- Instead of Gaussian, there are also other methods like: Epanechnikov, Quartic and Dic. For this analysis, there are 2 methods that can be applied.
- People tend to rent an Airbnb when it is nearer to the transportation services and tourist attractions. Hence, this resembles Quartic.
- People also tend to rent an Airbnb that are better in terms of price, past reviews and other factors. Hence, they are willing to book an Airbnb that are slightly further from places of interest. This would resemble Gaussian.
- Based on the above explanations, I would be using Quartic for this analysis since it is more relevant.
- Instead of Gaussian, there are also other methods like: Epanechnikov, Quartic and Dic. For this analysis, there are 2 methods that can be applied.
- Edge : Refers to the intensity estimate corrected for edge effect bias, the default is FALSE
- Bandwidth Selection
7.1.1.1 Explore different bandwidth methods
To define the bandwidth of the density estimation, we can explore the different bandwidth methods such as:
- bw.diggle(), bw.ppl(), bw.scott() and bw.CvL.
- To help us determine and select the appropriate bandwidth, I have created a function to run these 4 bandwidth methods easily for all our ppp objects.
- The function below performs 3 steps:
- With the inputted ppp object, run the 4 different automatic bandwidth methods
- Add the output and respective bandwidth names to a list
- Print the output list
abb_jun19 bw
bw_methods(abb_jun19SG_ppp.km)
$CvL
sigma
5.700728
$scott
sigma.x sigma.y
1.0748842 0.7425303
$ppl
sigma
0.2254825
$diggle
sigma
0.0419688 hotels bw
bw_methods(hotelsSG_ppp.km)
$CvL
sigma
9.092461
$scott
sigma.x sigma.y
1.2702643 0.7999742
$ppl
sigma
0.9148847
$diggle
sigma
0.08289152 tourism bw
bw_methods(tourismSG_ppp.km)
$CvL
sigma
7.186458
$scott
sigma.x sigma.y
2.234335 1.889744
$ppl
sigma
1.801735
$diggle
sigma
0.2901203 mrtlrt bw
bw_methods(mrtlrtSG_ppp.km)
$CvL
sigma
5.893095
$scott
sigma.x sigma.y
3.081230 2.104299
$ppl
sigma
1.211568
$diggle
sigma
0.6382647 Based on the above results:
- For all of the spatial points, the bw.diggle() automatic bandwidth method gives us a narrower bandwith.
- With further research, Baddeley et (2016) suggested using bw.ppl() algorithm because in their experience, it tends to produce more appropriate values when the pattern consists predominantly of tight clusters.
- However, Baddeley et also mentioned that if we want to detect a single tight cluster in the midst of random noise, then the bw.diggle() method seems to work best.
First, we have to understand what is the meaning of tight clusters
- According to a study, data points in a tight cluster are close to each other. Therefore, they tend to share the same kernel centers and have similar local densities.
- Based on our previous plots, it is observed that all of our spatial points (Airbnb 2019, Hotels and Tourism, MrtLrt) have significant clusters in the Central region.
- Hence, we will plot the kde map using bw.diggle() automatic bandwidth method.
7.1.1.2 Plot KDE Maps of Airbnb 2019 listings, hotels, MRT services, and tourist attractions
- The following code chunk will help us to plot 2 maps (KDE with bw.diggle() and Point Map) for each of our spatial point datasets for a clearer view.
7.1.1.2.1 Airbnb June 2019 listings
kde_diggle(abb_jun19SG_ppp.km)
7.1.1.2.2 Hotels
kde_diggle(hotelsSG_ppp.km)
7.1.1.2.3 Tourism
kde_diggle(tourismSG_ppp.km)
7.1.1.2.4 MrtLrt
kde_diggle(mrtlrtSG_ppp.km)
From the above maps,
- We can see that the South-east area of the Central region have higher density of Airbnb listings in 2019, indicated by the pink-orange spots on the KDE map.
- Similarly, for Hotels and Tourism, the South-east area of the Central region have higher density of Airbnb listings in 2019, indicated by the pink-orange spots on the KDE map.
- For MrtLrt, there are higher density areas are found in the Central, North and North-east region indicated by the brighter orange-yellow spots.
- Although there are black patches on the point map, it is better and wiser for us to use the KDE map since it has different colour scales to represent the density of the supposed cluster.
- The darker the purple shade, the less dense it is while the lighter the orange to yellow shade, the more dense the area.
- From the legend as well, we can see clearly that there are about 4 colours (Purple, Pink, Orange and Yellow) that we can identify and use to describe the level of density.
- However, with the point map, it is difficult for us to say which area is denser than the other area due to the single black colour used.
7.1.2 Display kernel density maps on openstreetmap of Singapore
7.1.2.1 Compute KDE maps with metres
- In this section, we do not have to rescale our point pattern into kilometres before computing our values as OpenStreetMap’s default units are in metres.
- Also, we are using bw.diggle as our bandwidth method.
kde_abb_jun19SG_ppp.m <- density(abb_jun19SG_ppp, sigma=bw.diggle, edge=TRUE,kernel="quartic")
kde_hotelsSG_ppp.m <- density(hotelsSG_ppp, sigma=bw.diggle, edge=TRUE, kernel="quartic")
kde_tourismSG_ppp.m <- density(tourismSG_ppp, sigma=bw.diggle, edge=TRUE, kernel="quartic")
kde_mrtlrtSG_ppp.m <- density(mrtlrtSG_ppp, sigma=bw.diggle, edge=TRUE, kernel="quartic")
7.1.2.2 Convert KDE output into grid objects
- The result will be the same but, we just convert it so that it is suitable for mapping purposes.
- The following code chunk uses as.SpatialGridDataFrame.im() to convert KDE output into grid objects.
gridded_kde_abb_jun19SG_dg <- as.SpatialGridDataFrame.im(kde_abb_jun19SG_ppp.m)
gridded_kde_hotels_dg <- as.SpatialGridDataFrame.im(kde_hotelsSG_ppp.m)
gridded_kde_tourism_dg <- as.SpatialGridDataFrame.im(kde_tourismSG_ppp.m)
gridded_kde_mrtlrt_dg <- as.SpatialGridDataFrame.im(kde_mrtlrtSG_ppp.m)
7.1.2.3 Convert gridded output into raster
- The following code chunk uses raster() of raster package to convert the gridded KDE objects to RasterLayer object.
kde_abb_jun19SG_dg_raster <- raster(gridded_kde_abb_jun19SG_dg)
kde_abb_jun19SG_dg_raster
class : RasterLayer
dimensions : 128, 128, 16384 (nrow, ncol, ncell)
resolution : 417.0614, 264.7348 (x, y)
extent : 2663.926, 56047.79, 16357.98, 50244.03 (xmin, xmax, ymin, ymax)
crs : NA
source : memory
names : v
values : -2.289334e-19, 0.001766134 (min, max)kde_hotels_dg_raster <- raster(gridded_kde_hotels_dg)
kde_hotels_dg_raster
class : RasterLayer
dimensions : 128, 128, 16384 (nrow, ncol, ncell)
resolution : 417.0614, 264.7348 (x, y)
extent : 2663.926, 56047.79, 16357.98, 50244.03 (xmin, xmax, ymin, ymax)
crs : NA
source : memory
names : v
values : -2.992963e-20, 0.000190199 (min, max)kde_tourism_dg_raster <- raster(gridded_kde_tourism_dg)
kde_tourism_dg_raster
class : RasterLayer
dimensions : 128, 128, 16384 (nrow, ncol, ncell)
resolution : 417.0614, 264.7348 (x, y)
extent : 2663.926, 56047.79, 16357.98, 50244.03 (xmin, xmax, ymin, ymax)
crs : NA
source : memory
names : v
values : -4.173818e-21, 1.775046e-05 (min, max)kde_mrtlrt_dg_raster <- raster(gridded_kde_mrtlrt_dg)
kde_mrtlrt_dg_raster
class : RasterLayer
dimensions : 128, 128, 16384 (nrow, ncol, ncell)
resolution : 417.0614, 264.7348 (x, y)
extent : 2663.926, 56047.79, 16357.98, 50244.03 (xmin, xmax, ymin, ymax)
crs : NA
source : memory
names : v
values : -1.538645e-21, 3.368543e-06 (min, max)From the results above:
- The CRS for all Raster objects are NA. We will need to assign the relevant projection systems before plotting a map with it.
7.1.2.4 Assign projection systems
- The following code chunk assigns the EPSG code of 3414 to the Raster objects
projection(kde_abb_jun19SG_dg_raster) <- CRS("+init=EPSG:3414")
kde_abb_jun19SG_dg_raster
class : RasterLayer
dimensions : 128, 128, 16384 (nrow, ncol, ncell)
resolution : 417.0614, 264.7348 (x, y)
extent : 2663.926, 56047.79, 16357.98, 50244.03 (xmin, xmax, ymin, ymax)
crs : +proj=tmerc +lat_0=1.36666666666667 +lon_0=103.833333333333 +k=1 +x_0=28001.642 +y_0=38744.572 +ellps=WGS84 +units=m +no_defs
source : memory
names : v
values : -2.289334e-19, 0.001766134 (min, max)projection(kde_hotels_dg_raster) <- CRS("+init=EPSG:3414")
kde_hotels_dg_raster
class : RasterLayer
dimensions : 128, 128, 16384 (nrow, ncol, ncell)
resolution : 417.0614, 264.7348 (x, y)
extent : 2663.926, 56047.79, 16357.98, 50244.03 (xmin, xmax, ymin, ymax)
crs : +proj=tmerc +lat_0=1.36666666666667 +lon_0=103.833333333333 +k=1 +x_0=28001.642 +y_0=38744.572 +ellps=WGS84 +units=m +no_defs
source : memory
names : v
values : -2.992963e-20, 0.000190199 (min, max)projection(kde_tourism_dg_raster) <- CRS("+init=EPSG:3414")
kde_tourism_dg_raster
class : RasterLayer
dimensions : 128, 128, 16384 (nrow, ncol, ncell)
resolution : 417.0614, 264.7348 (x, y)
extent : 2663.926, 56047.79, 16357.98, 50244.03 (xmin, xmax, ymin, ymax)
crs : +proj=tmerc +lat_0=1.36666666666667 +lon_0=103.833333333333 +k=1 +x_0=28001.642 +y_0=38744.572 +ellps=WGS84 +units=m +no_defs
source : memory
names : v
values : -4.173818e-21, 1.775046e-05 (min, max)projection(kde_mrtlrt_dg_raster) <- CRS("+init=EPSG:3414")
kde_mrtlrt_dg_raster
class : RasterLayer
dimensions : 128, 128, 16384 (nrow, ncol, ncell)
resolution : 417.0614, 264.7348 (x, y)
extent : 2663.926, 56047.79, 16357.98, 50244.03 (xmin, xmax, ymin, ymax)
crs : +proj=tmerc +lat_0=1.36666666666667 +lon_0=103.833333333333 +k=1 +x_0=28001.642 +y_0=38744.572 +ellps=WGS84 +units=m +no_defs
source : memory
names : v
values : -1.538645e-21, 3.368543e-06 (min, max)We can see from the results above that the CRS for all Raster objects are assigned the relevant projection systems.
7.1.2.5 Visualise KDE Raster objects in OpenStreetMap
- To help facilitate the plotting of the various OpenStreetMaps for investigation, we will create a function.
- The function will allow us to plot the Airbnb 2019 listing with the other location factors to see if there are any insights that we can get out of the comparison.
- The function makes use of the tmap package to display the raster object in OpenStreetMap of Singapore.
- We also make use of tm_basemap() to specify the visualization to be in OpenStreetMap
- Note, we should avoid using contrasting colours like Blues and Purples in the same map as it would be harder to tell the data apart.
plot_oms <- function(raster1, raster2, rt_legend1, rt_legend2) {
tmap_mode('view')
tm_basemap("OpenStreetMap")+
tm_shape(raster1) +
tm_raster("v", title = rt_legend1, alpha= 0.7, palette = "Reds" ) +
tm_shape(raster2) +
tm_raster("v", title = rt_legend2, alpha= 0.7, palette = "Blues")
}
7.1.2.5.1 Plot Airbnb Jun2019 with Hotels
plot_oms(kde_hotels_dg_raster, kde_abb_jun19SG_dg_raster, "Hotels", "Airbnb Jun2019")
tmap_mode('plot')
Spatial patterns observations:
- Airbnb 2019 listings are represented in Blue shades while Hotels are represented in shades of Reds.
- Similar to our previous observations from the KDE maps, for Hotels, there are higher density areas are found in the Central region.
- We can also see that the Airbnb listings here are cluttered in this Central region, indicated by the darker blue shade.
7.1.2.5.2 Plot Airbnb Jun2019 with Tourism
plot_oms(kde_tourism_dg_raster, kde_abb_jun19SG_dg_raster, "Tourism", "Airbnb Jun2019")
tmap_mode('plot')
Spatial patterns observations:
- Airbnb 2019 listings are represented in Blue shades while Tourism are represented in shades of Reds.
- Similar to our previous observations from the KDE maps, for Tourism, there are higher density areas are found in the Central region.
- We can also see that the Airbnb listings here are cluttered in this Central region, indicated by the darker blue shade.
7.1.2.5.3 Plot Airbnb Jun2019 with MrtLrt
plot_oms(kde_mrtlrt_dg_raster, kde_abb_jun19SG_dg_raster, "MrtLrt", "Airbnb Jun2019")
tmap_mode('plot')
Spatial patterns observations:
- Airbnb 2019 listings are represented in Blue shades while MrtLrt are represented in shades of Reds.
- Similar to our previous observations from the KDE maps, for MrtLrt, there are higher density areas are found in the Central, North-east and North-west region.
- Out of the 3 most obvious clusters in the Central, Northeast and Northwest areas, the Central region can be seen as the region with most dense MrtLrt Servies, highlighted by the darker red shade.
- We can also see that the Airbnb listings here are cluttered in this Central region, indicated by the darker blue shade.
Advantage of Kernel Density Map over Point Map:
- For point maps, while some clustering can be identified, it is difficult to pin a number to the clusters.
- When there are large numbers of points, it may get difficult to count and calculate the actual figures.
- The size of the points also has to be carefully selected to show the distribution.
- The Kernel Density Map however, offers an advantage over Point Maps as it can show the intensity of clustering in areas that have undergone first order spatial point pattern analysis.
- This intensity, represented by the shades of colours, make it easier for users to understand that the brighter the colour, the denser the area is while the darker the colour, the lesser points / the less dense that area is, therefore highlighting the more important insights.
- It has the ability to offer information on the estimate of listings in the areas as well.
7.1.3 Compare Spatial Point Patterns Using KDE
From the previous KDE Maps, we can observe that there are specific areas in the Central Region, such as Geylang, Kallang, Downtown Core and Rochor which have a higher density. Hence, we will look into these areas as intensive computational power is required to conduct a Second Order Analysis on all areas of Singapore.
7.1.3.1 Extract, Convert and Create Owin for Study Area
The following code chunk will:
- Extract the target area (Geylang, Kallang, Downtown and Rochor)
- Convert Spatial Polygons Data Frame into generic sp format
- Create owin object for each target area
- The function also requires 2 arguments, sz (in this case our mpsz) and the specific planning area we are interested in.
get_study_area <- function(plsz, area) {
plsz[plsz$PLN_AREA_N == area,] %>%
as('Spatial') %>%
as('SpatialPolygons') %>%
as('owin')
}
gl_owin <- get_study_area(mpsz_sf, "GEYLANG")
kl_owin <- get_study_area(mpsz_sf, "KALLANG")
dc_owin <- get_study_area(mpsz_sf, "DOWNTOWN CORE")
rc_owin <- get_study_area(mpsz_sf, "ROCHOR")
7.1.3.2 Combining Spatial points and the study area
abbjun19_gl_ppp = abb_jun19_ppp_jit[gl_owin]
hotels_gl_ppp = hotels_ppp_jit[gl_owin]
tourism_gl_ppp = tourism_ppp_jit[gl_owin]
mrtlrt_gl_ppp = mrtlrt_ppp[gl_owin]
abbjun19_kl_ppp = abb_jun19_ppp_jit[kl_owin]
hotels_kl_ppp = hotels_ppp_jit[kl_owin]
tourism_kl_ppp = tourism_ppp_jit[kl_owin]
mrtlrt_kl_ppp = mrtlrt_ppp[kl_owin]
abbjun19_dc_ppp = abb_jun19_ppp_jit[dc_owin]
hotels_dc_ppp = hotels_ppp_jit[dc_owin]
tourism_dc_ppp = tourism_ppp_jit[dc_owin]
mrtlrt_dc_ppp = mrtlrt_ppp[dc_owin]
abbjun19_rc_ppp = abb_jun19_ppp_jit[rc_owin]
hotels_rc_ppp = hotels_ppp_jit[rc_owin]
tourism_rc_ppp = tourism_ppp_jit[rc_owin]
mrtlrt_rc_ppp = mrtlrt_ppp[rc_owin]
7.1.3.6 Rescale from metres to kilometres
abbjun19_gl_ppp.km = rescale(abbjun19_gl_ppp, 1000, "km")
hotels_gl_ppp.km = rescale(hotels_gl_ppp, 1000, "km")
tourism_gl_ppp.km = rescale(tourism_gl_ppp, 1000, "km")
mrtlrt_gl_ppp.km = rescale(mrtlrt_gl_ppp, 1000, "km")
abbjun19_kl_ppp.km = rescale(abbjun19_kl_ppp, 1000, "km")
hotels_kl_ppp.km = rescale(hotels_kl_ppp, 1000, "km")
tourism_kl_ppp.km = rescale(tourism_kl_ppp, 1000, "km")
mrtlrt_kl_ppp.km = rescale(mrtlrt_kl_ppp, 1000, "km")
abbjun19_dc_ppp.km = rescale(abbjun19_dc_ppp, 1000, "km")
hotels_dc_ppp.km = rescale(hotels_dc_ppp, 1000, "km")
tourism_dc_ppp.km = rescale(tourism_dc_ppp, 1000, "km")
mrtlrt_dc_ppp.km = rescale(mrtlrt_dc_ppp, 1000, "km")
abbjun19_rc_ppp.km = rescale(abbjun19_rc_ppp, 1000, "km")
hotels_rc_ppp.km = rescale(hotels_rc_ppp, 1000, "km")
tourism_rc_ppp.km = rescale(tourism_rc_ppp, 1000, "km")
mrtlrt_rc_ppp.km = rescale(mrtlrt_rc_ppp, 1000, "km")
7.1.3.7 Plot KDE Maps of the study area
layout(matrix(c(1,2,3,4), 2, 2, byrow = TRUE), widths=c(5,5), heights=c(5,5))
par(mfrow=c(2,2))
plot(density(abbjun19_gl_ppp.km, sigma=bw.diggle, edge=TRUE, kernel="quartic"), main="Geylang")
plot(density(abbjun19_kl_ppp.km, sigma=bw.diggle, edge=TRUE, kernel="quartic"), main="Kallang")
plot(density(abbjun19_dc_ppp.km, sigma=bw.diggle, edge=TRUE, kernel="quartic"), main="Downtown Core")
plot(density(abbjun19_rc_ppp.km, sigma=bw.diggle, edge=TRUE, kernel="quartic"), main="Rochor")
From the KDE plots above, it is easy to think that Rohor is the most dense area. However, we should not be too quick to judge as the intensity level of Rochor is lower than Downtown core and Geylang.
7.1.3.8 Plot Airbnb, Hotels, Mrt, Tourist Attraction in each area
7.1.3.8.1 Geylang
par(mfrow=c(2,2), mar=c(1,1,1,1))
plot(abbjun19_gl_ppp.km, main = "Airbnb")
plot(hotels_gl_ppp.km, main = "Hotels")
plot(mrtlrt_gl_ppp.km, main = "Mrt")
plot(tourism_gl_ppp.km, main = "Tourist Attractions")
From the plots above:
- Despite having a few Hotels around Geylang and only 1 tourist attraction there, there are still many more Airbnbs as compared to Hotels.
- This might be because of the convenience of transport, indicated by the 10 MRT stations in Geylang area.
7.1.3.8.2 Kallang
par(mfrow=c(2,2), mar=c(1,1,1,1))
plot(abbjun19_kl_ppp.km, main = "Airbnb")
plot(hotels_kl_ppp.km, main = "Hotels")
plot(mrtlrt_kl_ppp.km, main = "Mrt")
plot(tourism_kl_ppp.km, main = "Tourist Attractions")
From the plots above:
- Similarly, for Kallang area, despite having a few Hotels and only 1 tourist attraction there, there are still many more Airbnbs as compared to Hotels.
- This might be because of the convenience of transport, indicated by the 8 MRT stations in Kallang area.
7.1.3.8.3 Downtown Core
par(mfrow=c(2,2), mar=c(1,1,1,1))
plot(abbjun19_dc_ppp.km, main = "Airbnb")
plot(hotels_dc_ppp.km, main = "Hotels")
plot(mrtlrt_dc_ppp.km, main = "Mrt")
plot(tourism_dc_ppp.km, main = "Tourist Attractions")
From the plots above:
- For Downtown core, there are still many Airbnbs despite the number of hotels around it.
- However, apart from the access to the multiple MRTs, there are also quite a number of tourist attractions located in this area.
7.1.3.8.4 Rochor
par(mfrow=c(2,2), mar=c(1,1,1,1))
plot(abbjun19_rc_ppp.km, main = "Airbnb")
plot(hotels_rc_ppp.km, main = "Hotels")
plot(mrtlrt_rc_ppp.km, main = "Mrt")
plot(tourism_rc_ppp.km, main = "Tourist Attractions")
From the plots above:
- For Rochor, the Airbnbs are more scattered across the area, aimilar to Hotels.
- There are quite a few Mrt services and Tourist attractions in that area as well.
Overall, below are the observations about the target areas and location factors:
- There seems to be some relationship between the Airbnbs and Hotels as both the number of Airbnbs and Hotels are quite high. Hence, having a Hotel nearby might mean that there is a possibility of an Airbnb being located in that area.
- However, the accessibility of MRT is one of the possible factor that might affect the number of Airbnbs in the area.
- Only for Downtown Core area and Rochor, there might be some kind of relationship between the Airbnbs and Tourist Attractions.
7.2 Second-order Spatial Point Patterns Analysis
To draw statistical conclusions on whether the distribution of Airbnb listings are affected by location factors such as hotels, tourist attractions and MRT services, we need to perform Second-order Spatial Point Analysis.
- For this section, we will be focusing on 2 areas, Kallang and Downtown Core
- The presence of water catchment areas may skew the calculations, so for this exercise, we limit by planning subzone area.
- Here, we look at Downtown Core since there were many Airbnbs and tourist attractions located in that area.
- As Kallang and Downtown Core differ in terms of the number of MrtLrt in its area, we will only focus on these 2 areas for our analysis since Geylang and Kallang have similar observations.
- Furthermore, since we are interested in investigating whether the distribution of Airbnb listings are affected by location factors, for e.g. whether the locations of MRT stations affect the distribution of Airbnb listing, we need to perform Multi-type Point Patterns Analysis.
- For this analysis, we will be using Cross L-Function.
- Although there are many functions, the most commonly used function is the L-function which is a transformed version of Ripley’s K-function. L-function is good because there is a red diagonal line in Ripley’s K-function which can be further transformed by deducting the distance for easier interpretation.
7.2.1 Combine ppp objects
kl_ppp <- superimpose(kl_abbjun2019_mark=abbjun19_kl_ppp,
kl_hotels_mark=hotels_kl_ppp,
kl_tourism_mark=tourism_kl_ppp,
kl_mrtlrt_mark=mrtlrt_kl_ppp)
dc_ppp <- superimpose(dc_abbjun2019_mark=abbjun19_dc_ppp,
dc_hotels_mark=hotels_dc_ppp,
dc_tourism_mark=tourism_dc_ppp,
dc_mrtlrt_mark=mrtlrt_dc_ppp)
7.2.2 Multi-type Point Patterns Analysis: Cross L-Function
lcross <- function(data, i, j, title) {
compute <- Lcross(data, i, j, correction='border')
plot(compute, . - r ~ r, xlab="d", ylab="L(d)-r", main=title)
}
7.2.2.1 Geylang
layout(matrix(c(1,2,3, 0), 2, 2, byrow = TRUE), widths=c(5,5), heights=c(5,5))
par(mfrow=c(2,2))
lcross(data = kl_ppp, i="kl_abbjun2019_mark", j="kl_hotels_mark", title="Airbnb 2019 & Hotels Lcross")
lcross(data = kl_ppp, i="kl_abbjun2019_mark", j="kl_tourism_mark", title="Airbnb 2019 & Tourism Lcross")
lcross(data = kl_ppp, i="kl_abbjun2019_mark", j="kl_mrtlrt_mark", title="Airbnb 2019 & MrtLrt Lcross")
7.2.2.2 Downtown Core
layout(matrix(c(1,2,3, 0), 2, 2, byrow = TRUE), widths=c(5,5), heights=c(5,5))
par(mfrow=c(2,2))
lcross(data = dc_ppp, i="dc_abbjun2019_mark", j="dc_hotels_mark", title="Airbnb 2019 & Hotels Lcross")
lcross(data = dc_ppp, i="dc_abbjun2019_mark", j="dc_tourism_mark", title="Airbnb 2019 & Tourism Lcross")
lcross(data = dc_ppp, i="dc_abbjun2019_mark", j="dc_mrtlrt_mark", title="Airbnb 2019 & MrtLrt Lcross")
Based on the plots of both Geylang and Downtown Core, - It shows that there is a sign that the marked spatial point events are not independent spatially. - However, a hypothesis test is required to confirm the observation statistically.
7.2.3 Performing CSR testing on the Cross L-Function
The hypothesis and test are as follows:
- Ho = The distribution of Airbnbs in 2019 and location factors such as hotels, tourist attractions, MrtLrt services are spatially independent.
- H1 = The distribution of Airbnbs in 2019 and location factors such as hotels, tourist attractions, MrtLrt services are NOT spatially independent.
The null hypothesis will be rejected if p-value is smaller than alpha value of 0.1 (i.e. at 90% confidence interval) (This applies to all points events above.)
Note: I chose 90% confidence interval for this analysis as higher intervals would take a longer time to run.
7.2.3.1 Monte Carlo test with Cross L-Function
- envelope() of spatstat package is used to perform CSR test
- As I chose 90% confidence interval, the number of simulations that will be run is 99+1 = 100.
- In the function below, we should indicate 99 as nsim argument since the first simulation is 0 and not 1.
- A global set.seed() is also used to ensure that the same output is computed each time the code is being ran.
- Note: set.seed(1234) is used for this analysis and is indicated at the top of this r markdown file together with the r setup settings.
# Geylang
abb2019_hotels_Lcross_kl.csr <- envelope(kl_ppp, Lcross,i="kl_abbjun2019_mark", j="kl_hotels_mark", correction='border', nsim=99)
Generating 99 simulations of CSR ...
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.
Done.abb2019_tourism_Lcross_kl.csr <- envelope(kl_ppp, Lcross,i="kl_abbjun2019_mark", j="kl_tourism_mark",correction='border', nsim=99)
Generating 99 simulations of CSR ...
1, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
2, 3, 4, 5, 6, 7, 8, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
9, 10, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
11, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
12, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
13, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
14, 15, 16, 17, 18, 19, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
20, 21, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
22, 23, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
24, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
25, 26, 27, 28, 29, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
30, 31, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
32, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
33, 34, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
35,
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
36, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
37, 38, 39, 40, 41, 42, 43, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
44, 45, 46, 47, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
48, 49, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
50, 51, 52, 53, 54, 55, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
56, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
57, 58, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
59, 60, 61, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
62, 63, 64, 65, 66, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
67, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
68, 69, 70,
71, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
72, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
73, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
74, 75, 76, 77, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
78, 79, 80, 81, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
82, 83, 84, 85, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
86, 87, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
88, 89, 90, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
91, 92, 93, 94, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
95, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
96, 97, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
98, Error in Kcross(X, i, j, ..., correction = correction) :
No points have mark j = kl_tourism_mark
[retrying]
99.
Done.abb2019_mrtlrt_Lcross_kl.csr <- envelope(kl_ppp, Lcross,i="kl_abbjun2019_mark", j="kl_mrtlrt_mark",correction='border', nsim=99)
Generating 99 simulations of CSR ...
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.
Done.# Downtown Core
abb2019_hotels_Lcross_dc.csr <- envelope(dc_ppp, Lcross,i="dc_abbjun2019_mark", j="dc_hotels_mark",correction='border', nsim=99)
Generating 99 simulations of CSR ...
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.
Done.abb2019_tourism_Lcross_dc.csr <- envelope(dc_ppp, Lcross,i="dc_abbjun2019_mark", j="dc_tourism_mark",correction='border', nsim=99)
Generating 99 simulations of CSR ...
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.
Done.abb2019_mrtlrt_Lcross_dc.csr <- envelope(dc_ppp, Lcross,i="dc_abbjun2019_mark", j="dc_mrtlrt_mark",correction='border', nsim=99)
Generating 99 simulations of CSR ...
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.
Done.7.2.4 Plot Cross L-Function
- plot() is used to plot the simulated Cross L-Function
7.2.4.1 Plot Cross L-Function for Geylang
layout(matrix(c(1,2,3, 0), 2, 2, byrow = TRUE), widths=c(5,5), heights=c(5,5))
par(mfrow=c(2,2))
plot(abb2019_hotels_Lcross_kl.csr, . -r ~ r, xlab="distance(m)")
plot(abb2019_tourism_Lcross_kl.csr, . -r ~ r, xlab="distance(m)")
plot(abb2019_mrtlrt_Lcross_kl.csr, . -r ~ r, xlab="distance(m)")
Recall how we should read the plot: - Red line represents the theoretical curve - Black line represents the observed curve - Grey shaded area represents the envelope derived based on the simulation.
For Geylang Airbnbs & Hotels:
- For distance < 350, the black line(observed) is above the red line(theo) and it is beyond the upper confidence envelope, we reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Hotels in Geylang are not spatially independent for this distance.
- For 350 < distance < 380, the black line(observed) is lower than the red line(theo) and is within lower confidence envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Hotels in Geylang are spatially independent for this distance.
- For 380 < distance < 400, the black line(observed) is lower than the red line(theo) and is beyond the lower confidence envelope, we reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Hotels in Geylang are not spatially independent for this distance.
For Geylang Airbnbs & Tourist Attractions:
- For distance < 50, the black line(observed) is above than the red line(theo) and is within the confidence envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Tourist Attractions in Geylang are spatially independent for this distance.
- For 50 < distance < 100, the black line(observed) is above than the red line(theo) and is beyond the lower confidence envelope, we reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Tourist Attractions in Geylang are not spatially independent for this distance.
- For 100 < distance < 400, the black line(observed) is within the confidence envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Tourist Attractions in Geylang are spatially independent for this distance.
For Geylang Airbnbs & Mrtlrt:
- For 0 < distance < 400, the black line(observed) is within the confidence envelope, hence we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Mrtlrt in Geylang are spatially independent for this distance.
7.2.4.2 Plot Cross L-Function for Downtown Core
layout(matrix(c(1,2,3, 0), 2, 2, byrow = TRUE), widths=c(5,5), heights=c(5,5))
par(mfrow=c(2,2))
plot(abb2019_hotels_Lcross_dc.csr, . -r ~ r, xlab="distance(m)")
plot(abb2019_tourism_Lcross_dc.csr, . -r ~ r, xlab="distance(m)")
plot(abb2019_mrtlrt_Lcross_dc.csr, . -r ~ r, xlab="distance(m)")
For Downtown Core Airbnbs & Hotels:
- For distance < 10, the black line(observed) is above red line(theo) and it is within the upper confidence envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Hotels in Downtown Core are spatially independent for this distance.
- For 10 < distance < 150, the black line(observed) is above the red line(theo) and is beyond the upper confidence envelope, we reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Hotels in Downtown Core are not spatially independent for this distance.
- For 150 < distance < 180, the first half of the black line(observed) above the red line(theo) while the second half of the black line(observed) is below the red line(theo). This blac line is the within the upper and lower confidence envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Hotels in Downtown Core are spatially independent for this distance.
- For 180 < distance < 300, the black line(observed) is below the red line(theo) and is beyond the lower confidence envelope, we reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Hotels in Downtown Core are not spatially independent for this distance.
For Downtown Core Airbnbs & Tourist Attractions:
- For distance < 50, the black line(observed) is below red line(theo) and it is within the lower confidence envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Tourist Attractions in Downtown Core are spatially independent for this distance.
- For 50 < distance < 320, the black line(observed) is below the red line(theo) and is beyond the lower confidence envelope, we reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Tourist Attractions in Downtown Core are not spatially independent for this distance.
For Downtown Core Airbnbs & MrtLrt:
- For distance < 50, the black line(observed) is below above the red line(theo) and is within the lower confidence envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Mrtlrt in Downtown Core are spatially independent for this distance.
- For 50 < distance < 80, the black line(observed) is above the red line(theo) and is within the upper confidence envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Mrtlrt in Downtown Core are spatially independent for this distance.
- For 80 < distance < 220, the black line(observed) is above the red line(theo) and is beyond the upper confidence envelope, we reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Mrtlrt in Downtown Core are not spatially independent for this distance.
- For 220 < distance < 280, the black line(observed) is above the red line(theo) and is within the upper confidence envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Mrtlrt in Downtown Core are spatially independent for this distance.
- For 280 < distance < 320, the black line(observed) is above the red line(theo) and is beyond the upper confidence envelope, we reject the null hypothesis at 90% confidence interval. We conclude that Airbnbs and Mrtlrt in Downtown Core are not spatially independent for this distance.
To conclude this section,
- There are some sort of evidence that the location factors such as Hotels, Tourist attractions and MrtLrt services do affect the distribution of Airbnbs in 2019.
- However, this is only at certain distances.
- Here, we can clearly see that for Downtown Core, the Airbnbs in this area are more likely to be affected by all 3 factors depending on the distance.
- However, for Geylang, the distribution might be affected by only 2 of the factors which is Hotels and Tourist attraction.
- Hence, we have statistical conclusions that the 2 points are not spatially independent of the location factors.
8. Section B: Impact of COVID-19
- In Section B, we will only be focusing on the room_type field. Hence, we can drop the other variables that we are not interested in.
- To ensure a fair comparison, we should also drop room types that do not exist in both 2019 and 2021 listings.
- Out of all the results, we can be fully 90% confident that Geylang Airbnbs & Mrtlrt are spatially independent.
8.1 Keep only required columns for analysis
keep <- c("room_type","latitude", "longtitude") # list of col names
abb_jun19 <- abb_jun19[,(names(abb_jun19) %in% keep)] #remove columns
glimpse(abb_jun19)
Formal class 'SpatialPointsDataFrame' [package "sp"] with 5 slots
..@ data :'data.frame': 8293 obs. of 1 variable:
.. ..$ room_type: chr [1:8293] "Private room" "Private room" "Private room" "Private room" ...
..@ coords.nrs : num(0)
..@ coords : num [1:8293, 1:2] 23825 22646 23922 41778 42057 ...
.. ..- attr(*, "dimnames")=List of 2
..@ bbox : num [1:2, 1:2] 7216 25166 44098 49226
.. ..- attr(*, "dimnames")=List of 2
..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slotkeep <- c("room_type","latitude", "longtitude") # list of col names
abb_jun21 <- abb_jun21[,(names(abb_jun21) %in% keep)] #remove columns
glimpse(abb_jun21)
Formal class 'SpatialPointsDataFrame' [package "sp"] with 5 slots
..@ data :'data.frame': 4238 obs. of 1 variable:
.. ..$ room_type: chr [1:4238] "Private room" "Private room" "Private room" "Private room" ...
..@ coords.nrs : num(0)
..@ coords : num [1:4238, 1:2] 23815 22646 23924 41973 42052 ...
.. ..- attr(*, "dimnames")=List of 2
..@ bbox : num [1:2, 1:2] 7407 25330 43338 48392
.. ..- attr(*, "dimnames")=List of 2
..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slot8.2 Determine the unique room type
The following code chunks will help us check the unique room types using unique feature of base R:
unique(abb_jun19$room_type)
[1] "Private room" "Entire home/apt" "Shared room" unique(abb_jun21$room_type)
[1] "Private room" "Entire home/apt" "Shared room"
[4] "Hotel room" From the results above,
- There are 3 room types in Airbnb 2019: Private room, Entire home/apt and Shared room
- There are 4 room types in Airbnb 2021: Private room, Entire home/apt, Shared room and Hotel room.
- With further research: the definition of these room types are:
- Entire home/apt: Guests have the whole place to themselves which usually includes a bedroom, a bathroom, and a kitchen.
- Private room: In contrast to Entire Home, guests have their own private room for sleeping.
- Shared room: Guests sleep in a common bedrrom that could be shared with others.
- Hotel room: Hotels listed on Airbnb
- To ensure a fair comparison, we will exclude Hotel room from Airbnb 2021 in the following code chunk:
abb_jun21 <- abb_jun21[!abb_jun21$room_type == "Hotel room", ]
unique(abb_jun21$room_type)
[1] "Private room" "Entire home/apt" "Shared room" 8.3 Check marked field data type and convert to factor
- str() of Base R will be used to check the data type of both abb_jun19 and abb_jun21 SpatialPointsDataFrame.
str(abb_jun19)
Formal class 'SpatialPointsDataFrame' [package "sp"] with 5 slots
..@ data :'data.frame': 8293 obs. of 1 variable:
.. ..$ room_type: chr [1:8293] "Private room" "Private room" "Private room" "Private room" ...
..@ coords.nrs : num(0)
..@ coords : num [1:8293, 1:2] 23825 22646 23922 41778 42057 ...
.. ..- attr(*, "dimnames")=List of 2
.. .. ..$ : NULL
.. .. ..$ : chr [1:2] "coords.x1" "coords.x2"
..@ bbox : num [1:2, 1:2] 7216 25166 44098 49226
.. ..- attr(*, "dimnames")=List of 2
.. .. ..$ : chr [1:2] "coords.x1" "coords.x2"
.. .. ..$ : chr [1:2] "min" "max"
..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slot
.. .. ..@ projargs: chr "+proj=tmerc +lat_0=1.36666666666667 +lon_0=103.833333333333 +k=1 +x_0=28001.642 +y_0=38744.572 +ellps=WGS84 +to"| __truncated__
.. .. ..$ comment: chr "PROJCRS[\"SVY21 / Singapore TM\",\n BASEGEOGCRS[\"SVY21\",\n DATUM[\"SVY21\",\n ELLIPSOID["| __truncated__str(abb_jun21)
Formal class 'SpatialPointsDataFrame' [package "sp"] with 5 slots
..@ data :'data.frame': 4040 obs. of 1 variable:
.. ..$ room_type: chr [1:4040] "Private room" "Private room" "Private room" "Private room" ...
..@ coords.nrs : num(0)
..@ coords : num [1:4040, 1:2] 23815 22646 23924 41973 42052 ...
.. ..- attr(*, "dimnames")=List of 2
.. .. ..$ : NULL
.. .. ..$ : chr [1:2] "coords.x1" "coords.x2"
..@ bbox : num [1:2, 1:2] 7407 25330 43338 48392
.. ..- attr(*, "dimnames")=List of 2
.. .. ..$ : chr [1:2] "coords.x1" "coords.x2"
.. .. ..$ : chr [1:2] "min" "max"
..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slot
.. .. ..@ projargs: chr "+proj=tmerc +lat_0=1.36666666666667 +lon_0=103.833333333333 +k=1 +x_0=28001.642 +y_0=38744.572 +ellps=WGS84 +to"| __truncated__
.. .. ..$ comment: chr "PROJCRS[\"SVY21 / Singapore TM\",\n BASEGEOGCRS[\"SVY21\",\n DATUM[\"SVY21\",\n ELLIPSOID["| __truncated__The output above shows that room_type field is in character data type and not in factor data type as required by spatstat package. Hence, the code chunk below will be used to convert room_type field to factor data type.
8.4 Exploratory Spatial Data Analysis
8.4.1 Standardise colour palettes
The following code chunks will help us assign specific colours to the unique room types for us to visualise and differentiate the room types easily.
roomtype_colours = c(`Entire home/apt`='lightgreen', `Private room`='plum', `Shared room`='coral')
8.4.2 Plot Airbnb listings Maps by room_type
In Singapore, there are generally 5 regions demarcated by the Urban Redevelopment Authority of Singapore to aid in its planning efforts. The regions are identified as:
- Central Region
- East Region
- North Region
- North-East Region
- West Region
8.4.2.1 Plot Airbnb 2019 listings Map by room_type
tmap_mode("view")
tm_basemap("OpenStreetMap") +
tm_shape(abb_jun19) +
tm_dots(col = "room_type",
size = 0.02,
alpha = 0.8,
palette=roomtype_colours) +
tm_view(text.size.variable = 1, view.legend.position = c("right", "bottom"))
tmap_mode('plot')
tm_shape(mpsz) +
tm_borders(alpha = 0.5) +
tm_shape(abb_jun19) +
tm_dots(col = 'room_type',
size = 0.5,
palette=roomtype_colours) +
tm_facets(by="room_type")+
tm_scale_bar() +
tm_layout(asp=1, outer.margins = 0)
Spatial Pattern Observations from Airbnb listings 2019:
At a quick glance, we can roughly see that the concentration of Private rooms (plum) is the highest while the concentration of Shared rooms (coral) is the lowest.
Entire homes/apt (lightgreen) are more commonly found throughout the Central, North-East and East regions.
Private rooms (plum) are found throughout the 5 regions - North, North-East, West, East and Central regions.
Shared rooms (coral) are sparsely scattered throughout 5 regions - North, South, East, West and Central regions.
A higher concentration of Entire homes/apts (lightgreen) rooms can be found in the Central region while higher concentrations of Private rooms (plum) are found in the Central, North-East and East regions.
For Shared rooms (coral), there is a mix of locations where some Shared rooms are at random locations while others are more concentrated at the Central region.
From the map, it seems that Private rooms (plum) and Entire homes/apts (lightgreen) display a higher chance of clustering while Shared rooms (coral) display randomness.
Just by looking at the exploratory map, it is difficult to draw meaningful observations of clustering, dispersion or random that are true to the behavior of spatial points.
8.4.2.2 Plot Airbnb 2021 listings Map by room_type
tmap_mode("view")
tm_basemap("OpenStreetMap") +
tm_shape(abb_jun21) +
tm_dots(col = "room_type",
size = 0.02,
alpha = 0.8,
palette=roomtype_colours) +
tm_view(text.size.variable = 1, view.legend.position = c("right", "bottom"))
tmap_mode('plot')
tm_shape(mpsz) +
tm_borders(alpha = 0.5) +
tm_shape(abb_jun21) +
tm_dots(col = 'room_type',
size = 0.5,
palette=roomtype_colours) +
tm_facets(by="room_type")+
tm_scale_bar() +
tm_layout(asp=1, outer.margins = 0)
Spatial Pattern Observations from Airbnb listings 2021:
Compared to the patterns observed in 2019 listings, there are few different observations that can be made:
- At a quick glance, we can roughly see that the concentration of Private rooms (plum) is the highest while the concentration of Shared rooms (coral) is the lowest.
- Entire homes/apt (lightgreen) are more commonly found throughout the Central and East regions.
- Private rooms (plum) are found throughout the 5 regions - North, North-East, West, East and Central regions.
- Shared rooms (coral) are sparsely scattered throughout the 5 regions - North, North-East, West, East and Central regions.
- A higher concentration of Entire homes/apts (lightgreen) rooms can be found in the Central region while higher concentrations of Private rooms (plum) are found in the Central, North-East and East regions.
- For Shared rooms (coral), there is a mix of locations where some Shared rooms are at random locations while others are more concentrated at the Central region.
Comparing Spatial Pattern Observations of 2019 and 2021:
- The concentrations of airbnb listings of all the room_types except hotel rooms of 2021 are much lesser compared the concentrations of rooms in 2019.
- Specifically, Entire homes/apt rooms (lightgreen) listings in 2021 are much less concentrated in the Central region as compared to 2019.
- Similarly, Private rooms (plum) listings in 2021 are much less concentrated in the North and Central region as compared to 2019.
- Underlying reason here might be because of the ongoing pandemic that we are experiencing right now.
8.5 Spatial Data Wrangling
8.5.1 Converting the SpatialPointsDataFrame into ppp format
8.5.1.1 Airbnb 2019 listings
abb_jun19_ppp <- as(abb_jun19, "ppp")
plot(abb_jun19_ppp, which.marks = "room_type")
Figure above reveals that there are 3 sub-rooms in the marks list. They are: Entire home/apt, Private room and Shared room.
8.5.1.2 Airbnb 2021 listings
abb_jun21_ppp <- as(abb_jun21, "ppp")
plot(abb_jun21_ppp, which.marks = "room_type")
Figure above reveals that there are 3 sub-rooms in the marks list. They are: Entire home/apt, Private room and Shared room.
8.5.2 Summary of ppp objects
To examine the summary statistics of this spatial object, summary() of Base R will be used as shown in the code chunk below:
8.5.2.1 Summary of abb_jun19_ppp object
summary(abb_jun19_ppp)
Marked planar point pattern: 8293 points
Average intensity 9.345289e-06 points per square unit
*Pattern contains duplicated points*
Coordinates are given to 3 decimal places
i.e. rounded to the nearest multiple of 0.001 units
Multitype:
frequency proportion intensity
Entire home/apt 4264 0.51416860 4.805054e-06
Private room 3582 0.43193050 4.036516e-06
Shared room 447 0.05390088 5.037193e-07
Window: rectangle = [7215.57, 44098.31] x [25166.35, 49226.35] units
(36880 x 24060 units)
Window area = 887399000 square unitsThe report above reveals that for AirBnb 2019 listings, Entire home/apt is the largest room type with Airbnb listings in Singapore with a market share of 51%. This is followed by Private room and Shared room.
8.5.2.2 Summary of abb_jun21_ppp object
summary(abb_jun21_ppp)
Marked planar point pattern: 4040 points
Average intensity 4.875562e-06 points per square unit
*Pattern contains duplicated points*
Coordinates are given to 3 decimal places
i.e. rounded to the nearest multiple of 0.001 units
Multitype:
frequency proportion intensity
Entire home/apt 1817 0.44975250 2.192796e-06
Private room 2050 0.50742570 2.473986e-06
Shared room 173 0.04282178 2.087803e-07
Window: rectangle = [7406.99, 43337.89] x [25330, 48391.55] units
(35930 x 23060 units)
Window area = 828622000 square unitsThe report above reveals that for AirBnb 2021 listings, Private room is the largest room type with Airbnb listings in Singapore with a market share of 50%. This is followed by Entire home/apt and Shared room.
It is also important to note that both the Airbnb 2019 and 2021 spatial point objects contains duplicated points. The quality of our analysis will be compromised if we failed to resolve this data issue. Hence, we will be resolving this issue in the following code chunk.
8.5.3 Resolve duplicate issue with Jittering method
abb_jun19_ppp_jit <- rjitter(abb_jun19_ppp, retry=TRUE, nsim=1, drop=TRUE)
any(duplicated(abb_jun19_ppp_jit))
[1] FALSEabb_jun21_ppp_jit <- rjitter(abb_jun21_ppp, retry=TRUE, nsim=1, drop=TRUE)
any(duplicated(abb_jun21_ppp_jit))
[1] FALSE8.5.4 Creating owin
- The presence of water catchment areas may skew the calculations, so for this analysis, we limit by planning subzone area.
- Here, we look at Downtown Core since there were many Airbnbs and tourist attractions located in that area.
8.5.4.1 Extracting study area
dc = mpsz[mpsz@data$PLN_AREA_N == "DOWNTOWN CORE",]
plot(dc, main = "Downtown Core")
8.5.4.2 Converting the spatial point data frame into generic sp format
- Here we will convert these SpatialPolygonsDataFrame layers into generic spatialpolygons layers by using as.SpatialPolygons.tess(x) of maptools package.
Code Chunk
dc_sp = as(dc, "SpatialPolygons")
str
str(dc_sp)
Formal class 'SpatialPolygons' [package "sp"] with 4 slots
..@ polygons :List of 12
.. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
.. .. .. ..@ Polygons :List of 1
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 29145 28467
.. .. .. .. .. .. ..@ area : num 103238
.. .. .. .. .. .. ..@ hole : logi FALSE
.. .. .. .. .. .. ..@ ringDir: int 1
.. .. .. .. .. .. ..@ coords : num [1:86, 1:2] 29201 29194 29194 29193 29191 ...
.. .. .. ..@ plotOrder: int 1
.. .. .. ..@ labpt : num [1:2] 29145 28467
.. .. .. ..@ ID : chr "22"
.. .. .. ..@ area : num 103238
.. .. .. ..$ comment: chr "0"
.. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
.. .. .. ..@ Polygons :List of 1
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 29392 29097
.. .. .. .. .. .. ..@ area : num 63665
.. .. .. .. .. .. ..@ hole : logi FALSE
.. .. .. .. .. .. ..@ ringDir: int 1
.. .. .. .. .. .. ..@ coords : num [1:116, 1:2] 29342 29343 29353 29356 29359 ...
.. .. .. ..@ plotOrder: int 1
.. .. .. ..@ labpt : num [1:2] 29392 29097
.. .. .. ..@ ID : chr "26"
.. .. .. ..@ area : num 63665
.. .. .. ..$ comment: chr "0"
.. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
.. .. .. ..@ Polygons :List of 1
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 29730 29011
.. .. .. .. .. .. ..@ area : num 196620
.. .. .. .. .. .. ..@ hole : logi FALSE
.. .. .. .. .. .. ..@ ringDir: int 1
.. .. .. .. .. .. ..@ coords : num [1:85, 1:2] 29808 29866 29805 29747 29733 ...
.. .. .. ..@ plotOrder: int 1
.. .. .. ..@ labpt : num [1:2] 29730 29011
.. .. .. ..@ ID : chr "27"
.. .. .. ..@ area : num 196620
.. .. .. ..$ comment: chr "0"
.. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
.. .. .. ..@ Polygons :List of 1
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 30126 28683
.. .. .. .. .. .. ..@ area : num 1070723
.. .. .. .. .. .. ..@ hole : logi FALSE
.. .. .. .. .. .. ..@ ringDir: int 1
.. .. .. .. .. .. ..@ coords : num [1:197, 1:2] 30437 30510 30511 30488 30487 ...
.. .. .. ..@ plotOrder: int 1
.. .. .. ..@ labpt : num [1:2] 30126 28683
.. .. .. ..@ ID : chr "31"
.. .. .. ..@ area : num 1070723
.. .. .. ..$ comment: chr "0"
.. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
.. .. .. ..@ Polygons :List of 1
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 29707 29745
.. .. .. .. .. .. ..@ area : num 39438
.. .. .. .. .. .. ..@ hole : logi FALSE
.. .. .. .. .. .. ..@ ringDir: int 1
.. .. .. .. .. .. ..@ coords : num [1:36, 1:2] 29814 29807 29786 29760 29752 ...
.. .. .. ..@ plotOrder: int 1
.. .. .. ..@ labpt : num [1:2] 29707 29745
.. .. .. ..@ ID : chr "37"
.. .. .. ..@ area : num 39438
.. .. .. ..$ comment: chr "0"
.. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
.. .. .. ..@ Polygons :List of 1
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 29969 29573
.. .. .. .. .. .. ..@ area : num 188767
.. .. .. .. .. .. ..@ hole : logi FALSE
.. .. .. .. .. .. ..@ ringDir: int 1
.. .. .. .. .. .. ..@ coords : num [1:97, 1:2] 30138 30138 30160 30140 30147 ...
.. .. .. ..@ plotOrder: int 1
.. .. .. ..@ labpt : num [1:2] 29969 29573
.. .. .. ..@ ID : chr "38"
.. .. .. ..@ area : num 188767
.. .. .. ..$ comment: chr "0"
.. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
.. .. .. ..@ Polygons :List of 1
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 30806 29530
.. .. .. .. .. .. ..@ area : num 521201
.. .. .. .. .. .. ..@ hole : logi FALSE
.. .. .. .. .. .. ..@ ringDir: int 1
.. .. .. .. .. .. ..@ coords : num [1:32, 1:2] 30845 30794 30649 30649 30600 ...
.. .. .. ..@ plotOrder: int 1
.. .. .. ..@ labpt : num [1:2] 30806 29530
.. .. .. ..@ ID : chr "41"
.. .. .. ..@ area : num 521201
.. .. .. ..$ comment: chr "0"
.. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
.. .. .. ..@ Polygons :List of 1
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 30380 29776
.. .. .. .. .. .. ..@ area : num 261844
.. .. .. .. .. .. ..@ hole : logi FALSE
.. .. .. .. .. .. ..@ ringDir: int 1
.. .. .. .. .. .. ..@ coords : num [1:55, 1:2] 30437 30304 30215 30214 30211 ...
.. .. .. ..@ plotOrder: int 1
.. .. .. ..@ labpt : num [1:2] 30380 29776
.. .. .. ..@ ID : chr "43"
.. .. .. ..@ area : num 261844
.. .. .. ..$ comment: chr "0"
.. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
.. .. .. ..@ Polygons :List of 1
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 30938 30472
.. .. .. .. .. .. ..@ area : num 886955
.. .. .. .. .. .. ..@ hole : logi FALSE
.. .. .. .. .. .. ..@ ringDir: int 1
.. .. .. .. .. .. ..@ coords : num [1:299, 1:2] 31353 31355 31357 31359 31364 ...
.. .. .. ..@ plotOrder: int 1
.. .. .. ..@ labpt : num [1:2] 30938 30472
.. .. .. ..@ ID : chr "45"
.. .. .. ..@ area : num 886955
.. .. .. ..$ comment: chr "0"
.. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
.. .. .. ..@ Polygons :List of 1
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 30209 30537
.. .. .. .. .. .. ..@ area : num 710569
.. .. .. .. .. .. ..@ hole : logi FALSE
.. .. .. .. .. .. ..@ ringDir: int 1
.. .. .. .. .. .. ..@ coords : num [1:170, 1:2] 30671 30645 30512 30444 30404 ...
.. .. .. ..@ plotOrder: int 1
.. .. .. ..@ labpt : num [1:2] 30209 30537
.. .. .. ..@ ID : chr "48"
.. .. .. ..@ area : num 710569
.. .. .. ..$ comment: chr "0"
.. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
.. .. .. ..@ Polygons :List of 1
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 30528 31278
.. .. .. .. .. .. ..@ area : num 280175
.. .. .. .. .. .. ..@ hole : logi FALSE
.. .. .. .. .. .. ..@ ringDir: int 1
.. .. .. .. .. .. ..@ coords : num [1:131, 1:2] 30574 30568 30562 30527 30515 ...
.. .. .. ..@ plotOrder: int 1
.. .. .. ..@ labpt : num [1:2] 30528 31278
.. .. .. ..@ ID : chr "58"
.. .. .. ..@ area : num 280175
.. .. .. ..$ comment: chr "0"
.. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
.. .. .. ..@ Polygons :List of 1
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 29398 28706
.. .. .. .. .. .. ..@ area : num 145519
.. .. .. .. .. .. ..@ hole : logi FALSE
.. .. .. .. .. .. ..@ ringDir: int 1
.. .. .. .. .. .. ..@ coords : num [1:96, 1:2] 29612 29584 29562 29545 29519 ...
.. .. .. ..@ plotOrder: int 1
.. .. .. ..@ labpt : num [1:2] 29398 28706
.. .. .. ..@ ID : chr "62"
.. .. .. ..@ area : num 145519
.. .. .. ..$ comment: chr "0"
..@ plotOrder : int [1:12] 4 9 10 7 11 8 3 6 12 1 ...
..@ bbox : num [1:2, 1:2] 28896 27914 31528 31714
.. ..- attr(*, "dimnames")=List of 2
.. .. ..$ : chr [1:2] "x" "y"
.. .. ..$ : chr [1:2] "min" "max"
..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slot
.. .. ..@ projargs: chr "+proj=tmerc +lat_0=1.36666666666667 +lon_0=103.833333333333 +k=1 +x_0=28001.642 +y_0=38744.572 +ellps=WGS84 +to"| __truncated__
.. .. ..$ comment: chr "PROJCRS[\"SVY21 / Singapore TM\",\n BASEGEOGCRS[\"SVY21\",\n DATUM[\"SVY21\",\n ELLIPSOID["| __truncated__8.5.4.3 Creating owin object
- We will convert these SpatialPolygons objects into owin objects that is required by spatstat.
Code Chunk
dc_owin = as(dc_sp, "owin")
str
str(dc_owin)
List of 5
$ type : chr "polygonal"
$ xrange: num [1:2] 28896 31528
$ yrange: num [1:2] 27914 31714
$ bdry :List of 12
..$ :List of 2
.. ..$ x: num [1:85] 29201 29192 29189 29184 29169 ...
.. ..$ y: num [1:85] 28595 28595 28595 28595 28595 ...
..$ :List of 2
.. ..$ x: num [1:115] 29342 29341 29339 29339 29338 ...
.. ..$ y: num [1:115] 29190 29197 29202 29203 29211 ...
..$ :List of 2
.. ..$ x: num [1:84] 29808 29764 29764 29766 29783 ...
.. ..$ y: num [1:84] 28757 28784 28785 28788 28817 ...
..$ :List of 2
.. ..$ x: num [1:196] 30437 30304 30215 30214 30211 ...
.. ..$ y: num [1:196] 29517 29416 29420 29397 29397 ...
..$ :List of 2
.. ..$ x: num [1:35] 29814 29815 29817 29820 29821 ...
.. ..$ y: num [1:35] 29617 29619 29626 29636 29640 ...
..$ :List of 2
.. ..$ x: num [1:96] 30138 30119 30115 30113 30095 ...
.. ..$ y: num [1:96] 29843 29881 29878 29877 29867 ...
..$ :List of 2
.. ..$ x: num [1:31] 30845 30895 30947 30968 30994 ...
.. ..$ y: num [1:31] 29022 29099 29177 29209 29257 ...
..$ :List of 2
.. ..$ x: num [1:54] 30437 30444 30445 30646 30770 ...
.. ..$ y: num [1:54] 29517 29529 29530 29842 30034 ...
..$ :List of 2
.. ..$ x: num [1:298] 31353 31352 31350 31347 31346 ...
.. ..$ y: num [1:298] 30695 30695 30695 30695 30695 ...
..$ :List of 2
.. ..$ x: num [1:169] 30671 30679 30747 30787 30819 ...
.. ..$ y: num [1:169] 30868 30878 30965 31009 31042 ...
..$ :List of 2
.. ..$ x: num [1:130] 30574 30580 30591 30601 30612 ...
.. ..$ y: num [1:130] 31006 31014 31030 31042 31056 ...
..$ :List of 2
.. ..$ x: num [1:95] 29612 29630 29636 29650 29633 ...
.. ..$ y: num [1:95] 28665 28693 28704 28723 28733 ...
$ units :List of 3
..$ singular : chr "unit"
..$ plural : chr "units"
..$ multiplier: num 1
..- attr(*, "class")= chr "unitname"
- attr(*, "class")= chr "owin"8.5.4.4 Combining Airbnb points and the study area
abbjun19_dc_ppp = abb_jun19_ppp_jit[dc_owin]
summary(abbjun19_dc_ppp)
Marked planar point pattern: 467 points
Average intensity 0.0001045044 points per square unit
Coordinates are given to 3 decimal places
i.e. rounded to the nearest multiple of 0.001 units
Multitype:
frequency proportion intensity
Entire home/apt 345 0.73875800 7.720343e-05
Private room 109 0.23340470 2.439181e-05
Shared room 13 0.02783726 2.909115e-06
Window: polygonal boundary
12 separate polygons (no holes)
vertices area relative.area
polygon 1 85 103238.0 0.02310
polygon 2 115 63665.0 0.01420
polygon 3 84 196620.0 0.04400
polygon 4 196 1070720.0 0.24000
polygon 5 35 39437.9 0.00883
polygon 6 96 188767.0 0.04220
polygon 7 31 521201.0 0.11700
polygon 8 54 261844.0 0.05860
polygon 9 298 886955.0 0.19800
polygon 10 169 710569.0 0.15900
polygon 11 130 280175.0 0.06270
polygon 12 95 145519.0 0.03260
enclosing rectangle: [28896.262, 31528.047] x [27914.19, 31714.21]
units
(2632 x 3800 units)
Window area = 4468710 square units
Fraction of frame area: 0.447abbjun21_dc_ppp = abb_jun21_ppp_jit[dc_owin]
summary(abbjun21_dc_ppp)
Marked planar point pattern: 243 points
Average intensity 5.437807e-05 points per square unit
Coordinates are given to 3 decimal places
i.e. rounded to the nearest multiple of 0.001 units
Multitype:
frequency proportion intensity
Entire home/apt 148 0.60905350 3.311915e-05
Private room 92 0.37860080 2.058758e-05
Shared room 3 0.01234568 6.713342e-07
Window: polygonal boundary
12 separate polygons (no holes)
vertices area relative.area
polygon 1 85 103238.0 0.02310
polygon 2 115 63665.0 0.01420
polygon 3 84 196620.0 0.04400
polygon 4 196 1070720.0 0.24000
polygon 5 35 39437.9 0.00883
polygon 6 96 188767.0 0.04220
polygon 7 31 521201.0 0.11700
polygon 8 54 261844.0 0.05860
polygon 9 298 886955.0 0.19800
polygon 10 169 710569.0 0.15900
polygon 11 130 280175.0 0.06270
polygon 12 95 145519.0 0.03260
enclosing rectangle: [28896.262, 31528.047] x [27914.19, 31714.21]
units
(2632 x 3800 units)
Window area = 4468710 square units
Fraction of frame area: 0.4478.5.4.5 Plotting Airbnb points and the study area
8.6 First Order Spatial Point Pattern Analysis
8.6.1 Plot KDE Maps
- In the code chunk below, density() of spatstat package is used to compute the kernel density objects.
- Then, plot() is used to plot the output kernel density objects derived.
- Instead of writing them in two seperate lines, the code chunk below shows how they can be combined into one single line code chunk.
plot(density((split(rescale(abbjun19_dc_ppp, 1000))), sigma=bw.diggle, edge=TRUE, kernel="quartic"),main="KDE Maps of Airbnb 2019 Room types")
plot(density((split(rescale(abbjun21_dc_ppp, 1000))), sigma=bw.diggle, edge=TRUE, kernel="quartic"),main="KDE Maps of Airbnb 2021 Room types")
From the KDE Maps above,
- For 2019 Airbnb listings,
- Entire home/apt rooms are more concentrated at the South-west of Downtown Core area. The intensity however is the highest out of the other 2 types of rooms and it is close to 2500 and 3000.
- Private homes have 3 concentrated areas in the North, Southwest and Central of Downtown Core area. The Central area however, is the most dense with a brighter yellow shade.
- Shared rooms are not that concentrated in this area however, as we can see from the darker blue shade.
- For 2021 Airbnb listings,
- Similar observations can be made for 2021 listings for Entire home/apt rooms as it is also more concentrated at the South-west of Downtown Core area. The intensity however has dropped to 1200 from the previous 3000.
- Here, the most obvious change can be seen from the KDE Map of Private homes as the previous 3 concentrated areas are longer found in 2021 KDE Map for Private homes. However, we must also note that the intensity scale did change quite largely from 200 to 2000. This might mean that because of COVID 19, there are much more Private homes listed as Airbnb however, it is too early to confirm this.
- Shared rooms in 2021 can be found concentrated in North and Central West of Downtown Core. However, the intensity scale dropped drastically from 200 to 14.
8.6.2 Reveal the density of Airbnbs by room types
intensity(rescale(abbjun19_dc_ppp, 1000))
Entire home/apt Private room Shared room
77.203430 24.391808 2.909115 intensity(rescale(abbjun21_dc_ppp, 1000))
Entire home/apt Private room Shared room
33.1191527 20.5875814 0.6713342 The output reveals that in Downtown Core area, both the Entire home/apt rooms for Airbnb 2019 and 2021 listings has the highest density of 77.43 units per km square and 32.67 units per km square respectively. This is followed by Private Rooms for both 2019 and 2021, with 25.06 and 21.26 units per km square and Shared Rooms for both 2019 and 2021, with 2.69 and 0.67 units per km square.
8.6.5 Plot KDE Maps on openstreetmap of Singapore
8.6.5.1 Extract each room type from SpatialPointsDataFrame
abb2019_pr <- abb_jun19[abb_jun19@data$room_type == "Private room",]
abb2019_er <- abb_jun19[abb_jun19@data$room_type == "Entire home/apt",]
abb2019_sr <- abb_jun19[abb_jun19@data$room_type == "Private room",]
abb2021_pr <- abb_jun21[abb_jun21@data$room_type == "Private room",]
abb2021_er <- abb_jun21[abb_jun21@data$room_type == "Entire home/apt",]
abb2021_sr <- abb_jun21[abb_jun21@data$room_type == "Private room",]
8.6.5.2 Converting the Spatial* class into generic sp format
abb2019_pr_sp <- as(abb2019_pr, "SpatialPoints")
abb2019_er_sp <- as(abb2019_er, "SpatialPoints")
abb2019_sr_sp <- as(abb2019_sr, "SpatialPoints")
abb2021_pr_sp <- as(abb2021_pr, "SpatialPoints")
abb2021_er_sp <- as(abb2021_er, "SpatialPoints")
abb2021_sr_sp <- as(abb2021_sr, "SpatialPoints")
8.6.3.3 Converting the generic sp format into spatstat’s ppp format
abb2019_pr_ppp <- as(abb2019_pr_sp, "ppp")
abb2019_er_ppp <- as(abb2019_er_sp, "ppp")
abb2019_sr_ppp <- as(abb2019_sr_sp, "ppp")
abb2021_pr_ppp <- as(abb2021_pr_sp, "ppp")
abb2021_er_ppp <- as(abb2021_er_sp, "ppp")
abb2021_sr_ppp <- as(abb2021_sr_sp, "ppp")
8.6.3.4 Check for duplicates in each ppp’s objects
any(duplicated(abb2019_pr_ppp))
[1] FALSEany(duplicated(abb2019_er_ppp))
[1] TRUEany(duplicated(abb2019_sr_ppp))
[1] FALSEany(duplicated(abb2021_pr_ppp))
[1] TRUEany(duplicated(abb2021_er_ppp))
[1] TRUEany(duplicated(abb2021_sr_ppp))
[1] TRUEAll ppp object except abb2019_pr_ppp and abb2019_sr_ppp have duplicates.
8.6.3.5 Handle duplicates and check again
abb2019_pr_ppp_jit <- rjitter(abb2019_pr_ppp, retry=TRUE, nsim=1, drop=TRUE)
abb2019_er_ppp_jit <- rjitter(abb2019_er_ppp, retry=TRUE, nsim=1, drop=TRUE)
abb2019_sr_ppp_jit <- rjitter(abb2019_sr_ppp, retry=TRUE, nsim=1, drop=TRUE)
abb2021_pr_ppp_jit <- rjitter(abb2021_pr_ppp, retry=TRUE, nsim=1, drop=TRUE)
abb2021_er_ppp_jit <- rjitter(abb2021_er_ppp, retry=TRUE, nsim=1, drop=TRUE)
abb2021_sr_ppp_jit <- rjitter(abb2021_sr_ppp, retry=TRUE, nsim=1, drop=TRUE)
any(duplicated(abb2019_pr_ppp_jit))
[1] FALSEany(duplicated(abb2019_er_ppp_jit))
[1] FALSEany(duplicated(abb2019_sr_ppp_jit))
[1] FALSEany(duplicated(abb2021_pr_ppp_jit))
[1] FALSEany(duplicated(abb2021_er_ppp_jit))
[1] FALSEany(duplicated(abb2021_sr_ppp_jit))
[1] FALSE8.6.3.6 Combining point events object and owin object
abb2019_pr_ppp = abb2019_pr_ppp_jit[sg_owin]
abb2019_er_ppp = abb2019_er_ppp_jit[sg_owin]
abb2019_sr_ppp = abb2019_sr_ppp_jit[sg_owin]
abb2021_pr_ppp = abb2021_pr_ppp_jit[sg_owin]
abb2021_er_ppp = abb2021_er_ppp_jit[sg_owin]
abb2021_sr_ppp = abb2021_sr_ppp_jit[sg_owin]
8.6.3.7 Compute KDE
kde_abb2019_pr_bw <- density(abb2019_pr_ppp, sigma=bw.diggle, edge=TRUE, kernel="quartic")
kde_abb2019_er_bw <- density(abb2019_er_ppp, sigma=bw.diggle, edge=TRUE, kernel="quartic")
kde_abb2019_sr_bw <- density(abb2019_sr_ppp, sigma=bw.diggle, edge=TRUE, kernel="quartic")
kde_abb2021_pr_bw <- density(abb2021_pr_ppp, sigma=bw.diggle, edge=TRUE, kernel="quartic")
kde_abb2021_er_bw <- density(abb2021_er_ppp, sigma=bw.diggle, edge=TRUE, kernel="quartic")
kde_abb2021_sr_bw <- density(abb2021_sr_ppp, sigma=bw.diggle, edge=TRUE, kernel="quartic")
8.6.3.8 Convert KDE output into grid object
gridded_kde_abb2019_pr_bw <- as.SpatialGridDataFrame.im(kde_abb2019_pr_bw)
gridded_kde_abb2019_er_bw <- as.SpatialGridDataFrame.im(kde_abb2019_er_bw)
gridded_kde_abb2019_sr_bw <- as.SpatialGridDataFrame.im(kde_abb2019_sr_bw)
gridded_kde_abb2021_pr_bw <- as.SpatialGridDataFrame.im(kde_abb2021_pr_bw)
gridded_kde_abb2021_er_bw <- as.SpatialGridDataFrame.im(kde_abb2021_er_bw)
gridded_kde_abb2021_sr_bw <- as.SpatialGridDataFrame.im(kde_abb2021_sr_bw)
8.6.3.9 Convert grid object into raster and assign project system
kde_abb2019_pr_bw_raster <- raster(gridded_kde_abb2019_pr_bw)
kde_abb2019_er_bw_raster <- raster(gridded_kde_abb2019_er_bw)
kde_abb2019_sr_bw_raster <- raster(gridded_kde_abb2019_sr_bw)
projection(kde_abb2019_pr_bw_raster) <- CRS("+init=EPSG:3414")
projection(kde_abb2019_er_bw_raster) <- CRS("+init=EPSG:3414")
projection(kde_abb2019_sr_bw_raster) <- CRS("+init=EPSG:3414")
kde_abb2021_pr_bw_raster <- raster(gridded_kde_abb2021_pr_bw)
kde_abb2021_er_bw_raster <- raster(gridded_kde_abb2021_er_bw)
kde_abb2021_sr_bw_raster <- raster(gridded_kde_abb2021_sr_bw)
projection(kde_abb2021_pr_bw_raster) <- CRS("+init=EPSG:3414")
projection(kde_abb2021_er_bw_raster) <- CRS("+init=EPSG:3414")
projection(kde_abb2021_sr_bw_raster) <- CRS("+init=EPSG:3414")
Check that Raster has been assigned correct CRS:
kde_abb2019_pr_bw_raster
class : RasterLayer
dimensions : 128, 128, 16384 (nrow, ncol, ncell)
resolution : 417.0614, 264.7348 (x, y)
extent : 2663.926, 56047.79, 16357.98, 50244.03 (xmin, xmax, ymin, ymax)
crs : +proj=tmerc +lat_0=1.36666666666667 +lon_0=103.833333333333 +k=1 +x_0=28001.642 +y_0=38744.572 +ellps=WGS84 +units=m +no_defs
source : memory
names : v
values : -1.197734e-19, 0.0008313512 (min, max)kde_abb2019_er_bw_raster
class : RasterLayer
dimensions : 128, 128, 16384 (nrow, ncol, ncell)
resolution : 417.0614, 264.7348 (x, y)
extent : 2663.926, 56047.79, 16357.98, 50244.03 (xmin, xmax, ymin, ymax)
crs : +proj=tmerc +lat_0=1.36666666666667 +lon_0=103.833333333333 +k=1 +x_0=28001.642 +y_0=38744.572 +ellps=WGS84 +units=m +no_defs
source : memory
names : v
values : -1.74862e-19, 0.001358564 (min, max)kde_abb2019_sr_bw_raster
class : RasterLayer
dimensions : 128, 128, 16384 (nrow, ncol, ncell)
resolution : 417.0614, 264.7348 (x, y)
extent : 2663.926, 56047.79, 16357.98, 50244.03 (xmin, xmax, ymin, ymax)
crs : +proj=tmerc +lat_0=1.36666666666667 +lon_0=103.833333333333 +k=1 +x_0=28001.642 +y_0=38744.572 +ellps=WGS84 +units=m +no_defs
source : memory
names : v
values : -1.006599e-19, 0.0008043117 (min, max)kde_abb2021_pr_bw_raster
class : RasterLayer
dimensions : 128, 128, 16384 (nrow, ncol, ncell)
resolution : 417.0614, 264.7348 (x, y)
extent : 2663.926, 56047.79, 16357.98, 50244.03 (xmin, xmax, ymin, ymax)
crs : +proj=tmerc +lat_0=1.36666666666667 +lon_0=103.833333333333 +k=1 +x_0=28001.642 +y_0=38744.572 +ellps=WGS84 +units=m +no_defs
source : memory
names : v
values : -1.458615e-19, 0.0007245677 (min, max)kde_abb2021_er_bw_raster
class : RasterLayer
dimensions : 128, 128, 16384 (nrow, ncol, ncell)
resolution : 417.0614, 264.7348 (x, y)
extent : 2663.926, 56047.79, 16357.98, 50244.03 (xmin, xmax, ymin, ymax)
crs : +proj=tmerc +lat_0=1.36666666666667 +lon_0=103.833333333333 +k=1 +x_0=28001.642 +y_0=38744.572 +ellps=WGS84 +units=m +no_defs
source : memory
names : v
values : -7.325024e-20, 0.0003351125 (min, max)kde_abb2021_sr_bw_raster
class : RasterLayer
dimensions : 128, 128, 16384 (nrow, ncol, ncell)
resolution : 417.0614, 264.7348 (x, y)
extent : 2663.926, 56047.79, 16357.98, 50244.03 (xmin, xmax, ymin, ymax)
crs : +proj=tmerc +lat_0=1.36666666666667 +lon_0=103.833333333333 +k=1 +x_0=28001.642 +y_0=38744.572 +ellps=WGS84 +units=m +no_defs
source : memory
names : v
values : -1.488199e-19, 0.0007517389 (min, max)8.6.3.10 Visualise output in openstreetmap
8.6.3.10.1 Private Room Airbnb 2019 and 2021
plot_oms(kde_abb2019_pr_bw_raster, kde_abb2021_pr_bw_raster, "Private Room 2019", "Private Room 2021")
tmap_mode('plot')
Spatial patterns observed:
- The concentration of Private Rooms for both years are not high.
- There are more private rooms in more areas in 2019 as compared to 2021 as we can see more red shades than blue shades.
- The areas shaded in purple (combination of red and blue), are those areas where there are significant concentrations of 2019 and 2021 Private Rooms listings. These areas are found obvious near Kallang and People’s Park.
- Other not so obvious areas are found in Little India, Selegie, Golden Mile and near Pasir Panjang where 2019 listings can be found.
8.6.3.10.2 Entire home/apt Airbnb 2019 and 2021
plot_oms(kde_abb2019_er_bw_raster, kde_abb2021_er_bw_raster, "Entire home/apt 2019", "Entire home/apt 2021")
tmap_mode('plot')
Spatial patterns observed:
- Here, we can see that there are more listings for Entire home/apt scattered across the Central region as compared to the previous map about Airbnb private rooms.
- Entire home/apt for 2021 listings have more listings in other areas as compared to 2019.
- The concentrations for both 2019 and 2021 Entire home/apt listings can be found near areas like Geylang, Kallang and Tanjong Pagar.
8.6.3.10.3 Shared Room Airbnb 2019 and 2021
plot_oms(kde_abb2019_sr_bw_raster, kde_abb2021_sr_bw_raster, "Shared Room 2019", "Shared Room 2021")
tmap_mode('plot')
Spatial patterns observed:
- Airbnb shared rooms for 2019 and 2021 listings are more concentrated within areas in Southeast area of the Central region such as People’s Park, Little India and Geylang.
- In 2019, the shared rooms are more concentrated in People Park and Kallang area as compared to Little India which is highlighted by a lighter shade of red.
- The shared rooms in 2021 are more concentrated at areas like Little India and People Park.
- This might signify the impact of COVID 19, where it shows that lesser Airbnb listings for shared room in Kallang but more in Little India in 2021.
8.7 Second-order Spatial Point Patterns Analysis
8.7.1 Combine ppp objects with owin object
abb2019_pr_ppp_dc = abb2019_pr_ppp_jit[dc_owin]
abb2019_er_ppp_dc = abb2019_er_ppp_jit[dc_owin]
abb2019_sr_ppp_dc = abb2019_sr_ppp_jit[dc_owin]
abb2021_pr_ppp_dc = abb2021_pr_ppp_jit[dc_owin]
abb2021_er_ppp_dc = abb2021_er_ppp_jit[dc_owin]
abb2021_sr_ppp_dc = abb2021_sr_ppp_jit[dc_owin]
8.7.2 Combine ppp objects
dc_ppp_s2 <- superimpose(abb2019_pr_dc_mark=abb2019_pr_ppp_dc,
abb2019_er_dc_mark=abb2019_er_ppp_dc,
abb2019_sr_dc_mark=abb2019_sr_ppp_dc,
abb2021_pr_dc_mark=abb2021_pr_ppp_dc,
abb2021_er_dc_mark=abb2021_er_ppp_dc,
abb2021_sr_dc_mark=abb2021_sr_ppp_dc)
8.7.3 Multi-type Point Patterns Analysis: Cross L-Function
layout(matrix(c(1,2,3,4,5,6), 3, 2, byrow = TRUE), widths=c(7,7), heights=c(7,7))
par(mfrow=c(3,2))
lcross(data = dc_ppp_s2, i="abb2019_pr_dc_mark", j="abb2019_er_dc_mark", title="PR 2019 & ER 2019")
lcross(data = dc_ppp_s2, i="abb2019_pr_dc_mark", j="abb2019_sr_dc_mark", title="PR 2019 & SR 2019")
lcross(data = dc_ppp_s2, i="abb2019_sr_dc_mark", j="abb2019_er_dc_mark", title="SR 2019 & ER 2019")
lcross(data = dc_ppp_s2, i="abb2021_pr_dc_mark", j="abb2021_er_dc_mark", title="PR 2021 & ER 2021")
lcross(data = dc_ppp_s2, i="abb2021_pr_dc_mark", j="abb2021_sr_dc_mark", title="PR 2021 & SR 2021")
lcross(data = dc_ppp_s2, i="abb2021_sr_dc_mark", j="abb2021_er_dc_mark", title="SR 2021 & ER 2021")
Based on the plots above,
- It shows that there is a sign that the marked spatial point events are not independent spatially.
- However, a hypothesis test is required to confirm the observation statistically.
8.7.4 Performing CSR testing on the Cross L-Function
- The hypothesis and test are as follows:
- Ho = The distribution of the room_type_a listings and room_type_b listings are spatially independent in Downtown Core area are spatially independent.
- H1 = The distribution of the room_type_a listings and room_type_b listings are spatially independent in Downtown Core area are NOT spatially independent.
- The null hypothesis will be rejected if p-value is smaller than alpha value of 0.1 (i.e. at 90% confident interval).
- The above hypothesis applies for the following pairs of points (room_type_a & room_type_b):
- Airbnb 2019 Private rooms & Airbnb 2019 Entire home/apt rooms
- Airbnb 2019 Private rooms & Airbnb 2019 Shared rooms
- Airbnb 2019 Shared rooms & Airbnb 2019 Entire home/apt rooms
- Airbnb 2021 Private rooms & Airbnb 2021 Entire home/apt rooms
- Airbnb 2021 Private rooms & Airbnb 2021 Shared rooms
- Airbnb 2021 Shared rooms & Airbnb 2021 Entire home/apt rooms
8.7.4.1 Monte Carlo test with Cross L-Function
- envelope() of spatstat package is used to perform CSR test
- As I chose 90% confidence interval, the number of simulations that will be run is 99+1 = 100.
- In the function below, we should indicate 99 as nsim argument since the first simulation is 0 and not 1.
- A global set.seed() is also used to ensure that the same output is computed each time the code is being ran.
- Note: set.seed(1234) is used for this analysis and is indicated at the top of this r markdown file together with the r setup settings.
abb2019pr_abb2019er_Lcross_dc.csr <- envelope(dc_ppp_s2, Lcross,
i="abb2019_pr_dc_mark", j="abb2019_er_dc_mark",
correction='border', nsim=99)
Generating 99 simulations of CSR ...
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.
Done.abb2019pr_abb2019sr_Lcross_dc.csr <- envelope(dc_ppp_s2, Lcross,
i="abb2019_pr_dc_mark", j="abb2019_sr_dc_mark",
correction='border', nsim=99)
Generating 99 simulations of CSR ...
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.
Done.abb2019sr_abb2019er_Lcross_dc.csr <- envelope(dc_ppp_s2, Lcross,
i="abb2019_sr_dc_mark", j="abb2019_er_dc_mark",
correction='border', nsim=99)
Generating 99 simulations of CSR ...
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.
Done.abb2021pr_abb2021er_Lcross_dc.csr <- envelope(dc_ppp_s2, Lcross,
i="abb2021_pr_dc_mark", j="abb2021_er_dc_mark",
correction='border', nsim=99)
Generating 99 simulations of CSR ...
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.
Done.abb2021pr_abb2021sr_Lcross_dc.csr <- envelope(dc_ppp_s2, Lcross,
i="abb2021_pr_dc_mark", j="abb2021_sr_dc_mark",
correction='border', nsim=99)
Generating 99 simulations of CSR ...
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.
Done.abb2021sr_abb2021er_Lcross_dc.csr <- envelope(dc_ppp_s2, Lcross,
i="abb2021_sr_dc_mark", j="abb2021_er_dc_mark",
correction='border', nsim=99)
Generating 99 simulations of CSR ...
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.
Done.8.7.5 Plot Cross L-Function
- plot() is used to plot the simulated Cross L-Function
layout(matrix(c(1,2,3,4,5,6), 3, 2, byrow = TRUE), widths=c(7,7), heights=c(7,7))
par(mfrow=c(3,2))
plot(abb2019pr_abb2019er_Lcross_dc.csr, . -r ~ r, xlab="distance(m)")
plot(abb2019pr_abb2019sr_Lcross_dc.csr, . -r ~ r, xlab="distance(m)")
plot(abb2019sr_abb2019er_Lcross_dc.csr, . -r ~ r, xlab="distance(m)")
plot(abb2021pr_abb2021er_Lcross_dc.csr, . -r ~ r, xlab="distance(m)")
plot(abb2021pr_abb2021sr_Lcross_dc.csr, . -r ~ r, xlab="distance(m)")
plot(abb2021sr_abb2021er_Lcross_dc.csr, . -r ~ r, xlab="distance(m)")
8.7.6 Conclusions for Second-order Point Patterns Analysis: Cross L-Function
- Recall how we should read the plot:
- Red line represents the theoretical curve
- Black line represents the observed curve
- Grey shaded area represents the envelope derived based on the simulation.
8.7.6.1 Private rooms and Entire home/apt rooms 2019 vs 2021
For Airbnb 2019 Private rooms & Airbnb 2019 Entire home/apt rooms:
- For distance < 250, the black line(observed) is above the red line(theo) and it is beyond the upper confidence envelope, we reject the null hypothesis at 90% confidence interval. We conclude that Airbnb 2019 Private rooms & Airbnb 2019 Entire home/apt rooms are not spatially independent for this distance.
- For 250 < distance < 350, the black line(observed) is within confidence envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnb 2019 Private rooms & Airbnb 2019 Entire home/apt rooms are spatially independent for this distance.
For Airbnb 2021 Private rooms & Airbnb 2021 Entire home/apt rooms:
- For distance < 20, the black line(observed) is above the red line(theo) and it is within the upper confidence envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnb 2021 Private rooms & Airbnb 2021 Entire home/apt rooms are spatially independent for this distance.
- For 20 < distance < 150, the black line(observed) is above the red line(theo) and it is beyond the upper confidence envelope, we reject the null hypothesis at 90% confidence interval. We conclude that Airbnb 2021 Private rooms & Airbnb 2021 Entire home/apt rooms are not spatially independent for this distance.
- For 150 < distance < 220, the black line(observed) is within the confidence envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnb 2021 Private rooms & Airbnb 2021 Entire home/apt rooms are spatially independent for this distance.
8.7.6.2 Private rooms and Shared rooms 2019 vs 2021
For Airbnb 2019 Private rooms & Airbnb 2019 Shared rooms:
- For distance < 250, the black line(observed) is above the red line(theo) and it is beyond the upper confidence envelope, we reject the null hypothesis at 90% confidence interval. We conclude that Airbnb 2019 Private rooms & Airbnb 2019 Shared rooms are not spatially independent for this distance.
- For 250 < distance < 350, the black line(observed) is within confidence envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnb 2019 Private rooms & Airbnb 2019 Shared rooms are spatially independent for this distance.
For Airbnb 2021 Private rooms & Airbnb 2021 Shared rooms:
- For distance < 180, the black line(observed) is above the red line(theo) and it is beyond the upper confidence envelope, we reject the null hypothesis at 90% confidence interval. We conclude that Airbnb 2021 Private rooms & Airbnb 2021 Shared rooms are not spatially independent for this distance.
- For 180 < distance < 250, the black line(observed) is within confidence envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnb 2021 Private rooms & Airbnb 2021 Shared rooms are spatially independent for this distance.
8.7.6.3 Shared rooms and Entire home/apt rooms 2019 vs 2021
For Airbnb 2019 Shared rooms & Airbnb 2019 Entire home/apt rooms:
- For distance < 5, the black line(observed) is above the red line(theo) and it is within the upper confidence envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnb 2019 Shared rooms & Airbnb 2019 Entire home/apt rooms are spatially independent for this distance.
- For 5 < distance < 210 and > 250, the black line(observed) is above the red line(theo) and it is beyond the upper confidence envelope, we reject the null hypothesis at 90% confidence interval. We conclude that Airbnb 2019 Shared rooms & Airbnb 2019 Entire home/apt rooms are not spatially independent for this distance.
- For 210 < distance < 240, the black line(observed) is above the red line(theo) and it is within the confidence envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnb 2019 Shared rooms & Airbnb 2019 Entire home/apt rooms are spatially independent for this distance.
For Airbnb 2021 Shared rooms & Airbnb 2021 Entire home/apt rooms:
- For distance < 5, the black line(observed) is above the red line(theo) and it is within the upper confidence envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnb 2021 Shared rooms & Airbnb 2021 Entire home/apt rooms are spatially independent for this distance.
- For 5 < distance < 150, the black line(observed) is above the red line(theo) and it is beyond the upper confidence envelope, we reject the null hypothesis at 90% confidence interval. We conclude that Airbnb 2021 Shared rooms & Airbnb 2021 Entire home/apt rooms are not spatially independent for this distance.
- For 150 < distance < 220, the black line(observed) is within the 2021 envelope, we cannot reject the null hypothesis at 90% confidence interval. We conclude that Airbnb 2021 Shared rooms & Airbnb 2019 Entire home/apt rooms are spatially independent for this distance.
Comparing the 2019 and 2021 listings by room types,
- For Private rooms and Entire home/apt rooms in 2019 and 2021, the distance that we can reject the null hypothesis has dropped from <250 to <150.
- For Private rooms and Shared rooms in 2019 and 2021, the distance that we can reject the null hypothesis has dropped from <250 to <180. -For Shared rooms and Entire home/apt rooms in 2019 and 2021, the distance that we can reject the null hypothesis has dropped from 5 < distance < 210 to 5 < distance < 150.
- Seeing how all the distances of the spatial independence of room types in 2021 has dropped from 2019, we can conclude that COVID-19 indeed had a significant impact on Airbnb’s business in Singapore’s Downtown Core area.
9. Conclusion for Section A and B
Overall, the spatial patterns observed may be due to how near these Airbnb listings are to the location factors. As Airbnb are commonly used by tourists, tourists may choose to stay in Airbnb that are more convenient for them to travel around Singapore. Hence due to the demand, there are more listings around the Center and South regions of Singapore which are also where the tourist attractions or hotspots are in Singapore.
Additionally, for Singaporeans who may want to find an Airbnb to stay temporarily, likely because of moving houses or having a staycation. Hence, if they were to find an Airbnb to stay in, it is most likely near Town/CBD or the South regions as it is more convenient for them to go to work or again, where the shopping paradise is or places of interest are at.
As for the COVID 19 Impact on Airbnbs, we saw for ourselves that there are generally lesser Airbnb listings in 2021 as compared to 2019. This might be due to the travel restrictions where tourists are not able to enter Singapore anymore and also the COVID-19 regulations where the owners of the Airbnbs are less willing to rent out their rooms because of the fear of the rising cases. Also, since the distances for Downtown Core has dropped significantly as seen from the different room types of Airbnb listings in 2021, based on this sample, we can say that we can safely conclude that COVID-19 has made it much more difficult for Airbnb business to operate in Singapore.
10. References
- Department of Computer Science, Worcester Polytechnic Institute. (2019, March 29). Scalable Kernel Density Estimation-based Local Outlier Detection over Large Data Streams∗. Retrieved from: https://openproceedings.org/2019/conf/edbt/EDBT19_paper_115.pdf
- Gimond, M. (n.d.). Chapter 11 point pattern analysis | Intro to GIS and spatial analysis. Main repos | mgimond.github.io. Retrieved from: https://mgimond.github.io/Spatial/chp11_0.html#distance-based-analysis
- Solymosi, R., & Medina, J. (2021, April 28). Chapter 6 studying spatial point patterns | Crime mapping in R. GitHub Pages. Retrieved from: https://maczokni.github.io/crime_mapping_textbook/studying-spatial-point-patterns.html
- Tin Seong.K (2021, September 06) Hands-on Exercise 4: Spatial Point Patterns Analysis-spatstat methods. Retrieved from: https://is415.netlify.app/hands-on_ex/hands-on_ex04/hands-on_ex04-sppa
- Tin Seong.K (2021, September 13)Hands-on Exercise 5: Analysing Marked Point Patterns. Retrieved from: https://is415.netlify.app/hands-on_ex/hands-on_ex05/hands-on_ex05-markedsppa