Statystyki centrograficzne stanowią przestrzenną modyfikację statystyk opisowych w klasycznej statystyce. Są podstawową formą opisu rozkładu przestrzennego danych punktowych. Statystyki centrograficzne oblicza się za pomocą funkcji z pakietu sfdep.
Żródło: https://mgimond.github.io/Spatial/chp11_0.html#centrography
Średnia centralna wyznaczona jest przez punkt o współrzędnych \(\overline{x_c}\), \(\overline{y_c}\) stanowiących średnią arytmetyczną długości i szerokości geograficznej obiektów rozproszonych w przestrzeni.
\[\overline{x_c} = \frac{\sum_{i=1}^{n} x_i}{n}\]
\[\overline{y_c} = \frac{\sum_{i=1}^{n} y_i}{n}\]
Statystyki centrograficzne (średnia centralna, odległość standardowa) mogą być obliczone na podstawie ramki danych zawierającej współrzędne x oraz y. Funkcja st_coordinates() z pakietu sf pozwala na wydobycie współrzędnych z obiektu punktowego klasy sf oraz stworzenie ramki danych zawierającej dwie kolumny: x oraz y.
library(sf)Linking to GEOS 3.10.2, GDAL 3.4.1, PROJ 8.2.1; sf_use_s2() is TRUEp2019 = read_sf("data/przestepstwa_2019.gpkg")
#granica miasta Poznań
pzn = read_sf("data/pzn_borders.gpkg")
p2019_coords = st_coordinates(p2019)
head(p2019_coords) X Y
[1,] 361379.5 503546.3
[2,] 361691.3 503694.6
[3,] 360347.3 506216.0
[4,] 353099.3 504204.1
[5,] 361715.2 503305.9
[6,] 360260.9 510415.0Średnią centralną obliczamy poprzez wyliczenie średniej wartości dla kolumny X oraz Y.
# mean center
mc = apply(p2019_coords, 2, mean)
mc X Y
358716.4 506363.5 Wizualizacja średniej centralnej wymaga utworzenia obiektu przestrzennego w oparciu o współrzędne punktu. W tym celu wykorzystanie zostanie funkcja st_point() oraz st_sfc() z pakietu sf.
#Przekształcenie obiektu mc na obiekt przestrzenny klasy sfc_POINT
mc_point = st_sfc(st_point(mc), crs = st_crs(p2019))plot(mc_point)Odległość standardowa określa średnią odległość puntków (\(x_{i}\), \(y_{i}\)) od punktu centralnego (centroidu, (\(\overline{x_c}\), \(\overline{y_c}\))) i obliczana jest wg wzoru:
\[d = \sqrt{\frac{\sum_{i=1}^{n} d_{ic}^2}{n}}\]
\[d_{ic}^2 = (x_i - \overline{x_c})^2 + (y_i - \overline{y_c})^2 \]
W powyższym wzorze
# standard distance
sd = sqrt(sum((p2019_coords[,1] - mc[1])^2 + (p2019_coords[,2] - mc[2])^2) / nrow(p2019_coords))
sd[1] 3886.935Odległość standardowa graficznie przedstawiana jest jako koło o środku w punkcie wyznaczonym przez średnią centralną oraz promieniu równym odległości standardowej. Do stworzenia obiektu przestrzennego przedstawiającego zasięg odległości standardowej wykorzystana zostanie funkcja st_buffer() z pakietu sf.
#wyznaczenie strefy buforowej wokół średniej centralnej o promieniu równym odległości standardowej
sd_buffer = st_buffer(mc_point, sd)library(tmap)
tm_shape(pzn) +
tm_borders() +
tm_shape(p2019) +
tm_dots() +
tm_shape(mc_point) +
tm_symbols(fill = "red") +
tm_shape(sd_buffer) +
tm_borders(col = "red", lwd = 2)
Elipsa odchylenia standardowego pozwala na określenia kierunków rozrzutu obiektów (obserwacji) w przestrzeni. W R może być obliczona za pomocą funkcji std_dev_ellipse() z pakietu sfdep.
Statystyki centrograficzne w R oblicza sie z wykorzystaniem funkcji z pakietu sfdep. Pakiet sfdep pozwala na obliczanie różnych statystyk centrograficznych w oparciu o dane punktowe wczytane za pomocą pakietu sf:
center_mean() - średnia centralnastd_distance() - odległość standardowastd_dev_ellipse() - elipsa odchylenia standardowegoW przykładzie 1 wykorzystano dane dotyczące przestępczości w Poznaniu (plik przestepstwa_2019.gpkg) wraz z granicą miasta Poznania (plik pzn_borders.gpkg). Do podsumowania rozkładu przestępczości wykorzystano statystyki centrograficzne. Plik przestepstwa_2019.gpkg zawiera tylko kolumnę geom bez dodatkowych atrybutów.
library(sf)
# dane punktowe
p2019 = read_sf("data/przestepstwa_2019.gpkg")
#granica miasta Poznań
pzn = read_sf("data/pzn_borders.gpkg")head(p2019)Simple feature collection with 6 features and 0 fields
Geometry type: POINT
Dimension: XY
Bounding box: xmin: 353099.3 ymin: 503305.9 xmax: 361715.2 ymax: 510415
Projected CRS: ETRF2000-PL / CS92
# A tibble: 6 × 1
geom
<POINT [m]>
1 (361379.5 503546.3)
2 (361691.3 503694.6)
3 (360347.3 506216)
4 (353099.3 504204.1)
5 (361715.2 503305.9)
6 (360260.9 510415)library(tmap)
tm_shape(pzn) +
tm_borders() +
tm_shape(p2019) +
tm_dots()
Statystyki centrograficzne w R można wyznaczyć używając funkcji z pakietu sfdep:
center_mean() - średnia centralnalibrary(pracma)
library(sfdep)
#srednia centralna
mc_point2 = center_mean(p2019)std_distance() - odległość standardowalibrary(pracma)
library(sfdep)
# odległość standardowa
sd2 = std_distance(p2019) #zwraca wartość z odległością
sd2_buffer = st_buffer(mc_point2, sd2) #zwraca obiekt przestrzennystd_dev_ellipse() - elipsa odchylenia standardowegolibrary(pracma)
library(sfdep)
# funkcja std_dev_ellipse zwraca parametry elipsy, st_ellipse tworzy na ich podstawie obiekt liniowy
p2019_e = std_dev_ellipse(p2019)
p2019_e = st_ellipse(p2019_e, sx = p2019_e$sx, sy = p2019_e$sy)tm_shape(pzn) +
tm_borders() +
tm_shape(p2019) +
tm_dots() +
tm_shape(mc_point) +
tm_dots(fill = "red", size = 2) +
tm_shape(sd2_buffer) +
tm_borders(col = "red", lwd = 3) +
tm_shape(p2019_e) +
tm_lines(col = "darkgreen", lwd = 3)
W drugim przykładzie wykorzystano dane pozyskane z BDOT10k z warstwy przedstawiającej budynki (OT_BUBD_A) oraz granice administracyjne (OT_ADJA_A) do przeanalizowania rozkładu przestrzennego szkół podstawowych w Poznaniu.
Granica miasta Pozania została pozyskana z warstwy OT_ADJA_A.
granica = read_sf("data/bdot10k/PL.PZGiK.308.BDOT10k.3064__OT_ADJA_A.gpkg")
granica = granica[granica$RODZAJ == "powiat", "geometria"]Warstwa OT_BUBD_A jest warstwą poligonową zawierającą różne typy budynków.
budynki <- read_sf("data/bdot10k/PL.PZGiK.308.BDOT10k.3064__OT_BUBD_A.gpkg")head(budynki)Simple feature collection with 6 features and 19 fields
Geometry type: POLYGON
Dimension: XY
Bounding box: xmin: 359155.6 ymin: 500467.5 xmax: 362596.6 ymax: 502295.9
Projected CRS: ETRF2000-PL / CS92
# A tibble: 6 × 20
TERYT LOKALNYID PRZESTRZENNAZW WERSJA POCZATEKWERSJIOBIEKTU
<chr> <chr> <chr> <dttm> <dttm>
1 3064 3c25a636-5303-… PL.PZGiK.308.… 2023-12-31 12:00:00 2023-12-31 12:00:00
2 3064 12e66037-0d58-… PL.PZGiK.308.… 2023-12-31 12:00:00 2023-12-31 12:00:00
3 3064 021dfd88-65df-… PL.PZGiK.308.… 2023-12-31 12:00:00 2023-12-31 12:00:00
4 3064 5b5acb18-3c3a-… PL.PZGiK.308.… 2023-12-31 12:00:00 2023-12-31 12:00:00
5 3064 a4b459f2-551c-… PL.PZGiK.308.… 2023-12-31 12:00:00 2023-12-31 12:00:00
6 3064 3fddfa4f-6f68-… PL.PZGiK.308.… 2023-12-31 12:00:00 2023-12-31 12:00:00
# ℹ 15 more variables: OZNACZENIEZMIANY <chr>,
# ZRODLODANYCHGEOMETRYCZNYCH <chr>, KATEGORIAISTNIENIA <chr>, UWAGI <chr>,
# INFORMACJADODATKOWA <chr>, KODKARTO10K <chr>, SKROTKARTOGRAFICZNY <chr>,
# KODKST <chr>, FUNKCJAOGOLNABUDYNKU <chr>, PRZEWAZAJACAFUNKCJABUDYNKU <chr>,
# LICZBAKONDYGNACJI <dbl>, NAZWA <chr>, FSBUD <chr>, IDEGIB <chr>,
# geometria <POLYGON [m]>Informacje o typach budynków zapisane są w dwóch kolumnach: FUNKCJAOGOLNABUDYNKU oraz PRZEWAZAJACAFUNKCJABUDYNKU.
unique(budynki$FUNKCJAOGOLNABUDYNKU) [1] "budynki szpitali i inne budynki opieki zdrowotnej"
[2] "budynki transportu i łączności"
[3] "budynki mieszkalne"
[4] "budynki produkcyjne, usługowe i gospodarcze dla rolnictwa"
[5] "zbiorniki, silosy i budynki magazynowe"
[6] "budynki przemysłowe"
[7] "budynki oświaty, nauki i kultury oraz budynki sportowe"
[8] "budynki handlowo-usługowe"
[9] "pozostałe budynki niemieszkalne"
[10] "budynki biurowe" funcja_szczegolowa_budynku = unique(budynki$PRZEWAZAJACAFUNKCJABUDYNKU)Szkoły podstawowe to kategoria określona jako przeważająca funkcja budynku (kolumna PRZEWAZAJACAFUNKCJABUDYNKU).
szkoly_podstawowe = budynki[budynki$PRZEWAZAJACAFUNKCJABUDYNKU%in%c("szkoła podstawowa"),]Warstwa BUBD_A jest warstwą poligonową. Aby zamienić budynki (poligony) na punkty należy użyć funkcji st_centroid() z pakietu sf.
szkoly_pkt = st_centroid(szkoly_podstawowe[, "geometria"])st_write(szkoly_pkt, dsn = "data/out_poznan_szkoly.gpkg", layer = "szkoly", delete_dsn = TRUE)Deleting source `data/out_poznan_szkoly.gpkg' using driver `GPKG'
Writing layer `szkoly' to data source
`data/out_poznan_szkoly.gpkg' using driver `GPKG'
Writing 84 features with 0 fields and geometry type Point.st_write(granica, dsn = "data/out_poznan_szkoly.gpkg", layer = "granica")Writing layer `granica' to data source
`data/out_poznan_szkoly.gpkg' using driver `GPKG'
Writing 1 features with 0 fields and geometry type Polygon.Do wizualizacji lokalizacji szkół podstawowych w mieście Poznań wykorzystano pakiet tmap().
library(tmap)
tm_shape(granica) +
tm_borders(col = "black") +
tm_shape(szkoly_pkt) +
tm_dots(col = "black", size = 0.7)
Wykorzystując funkcje z pakietu sfdep zostały obliczone statystyki centrograficzne: średnia centralna (funkcja center_mean()), odległość standardowa (funkcja std_distance()), elipsa odchylenia standardowego (std_dev_ellipse()).
library(pracma)
library(sfdep)
#srednia centralna
mc_point2 = center_mean(szkoly_pkt)
# odległość standardowa
sd2 = std_distance(szkoly_pkt) #zwraca wartość z odległością
sd2_buffer = st_buffer(mc_point2, sd2) #zwraca obiekt przestrzenny
# funkcja std_dev_ellipse zwraca parametry elipsy, st_ellipse tworzy na ich podstawie obiekt liniowy
szkoly_pkt_e = std_dev_ellipse(szkoly_pkt)
szkoly_pkt_e = st_ellipse(szkoly_pkt_e, sx = szkoly_pkt_e$sx, sy = szkoly_pkt_e$sy)tm_shape(granica) +
tm_borders(col = "black") +
tm_shape(szkoly_pkt) +
tm_dots(fill = "black", size = 0.7) +
tm_shape(mc_point2) +
tm_dots(fill = "red", size = 2) +
tm_shape(sd2_buffer) +
tm_borders(col = "red", lwd = 3) +
tm_shape(szkoly_pkt_e) +
tm_lines(col = "darkgreen", lwd = 3)