4 Eksploracyjna analiza danych przestrzennych

library("raster")
library("sp")
library("gstat")
library("dismo")
library("rgeos")
library("geostatbook")
library("deldir")
data(punkty)
data(granica)

4.1 Mapy

4.1.1 Podstawowe terminy | Kontekst przestrzenny

  • Populacja - cały obszar, dla którego chcemy określić wybrane właściwości
  • Próba - zbiór obserwacji, dla których mamy informacje. Inaczej, próba to podzbiór populacji. Zazwyczaj niemożliwe (lub bardzo kosztowne) jest zdobycie informacji o całej populacji. Z tego powodu bardzo cenne jest odpowiednie wykorzystanie informacji z próby.

4.1.2 Mapy punktowe

Eksploracyjna analiza danych przestrzennych w przypadku danych punktowych ma na celu:

  • Wgląd w typ próbkowania danych
  • Sprawdzenie poprawności współrzędnych
  • Sprawdzenie poprawności danych, w tym między innymi określenie danych odstających lokalnie
  • Identyfikacja głównych cech struktury przestrzennej zjawiska (np. trend)

4.2 Próbkowanie

4.2.1 Podstawowe typy próbowania

Istnieje cały szereg typów próbkowania danych przestrzennych. Funkcja spsample() z pakietu sp pozwala na stworzenie kilku typów próbkowania (argument type), między innymi:

  • Regularny (ang.regular)
  • Losowy (ang.random)
  • Losowy stratyfikowany (ang.stratified)
  • Preferencyjny (ang.clustered)

4.2.2 Typ próbowania | Regularny

W regularnym typie próbkowania, kolejne punkty rozłożone są równomiernie na badanym obszarze.

set.seed(225)
regularny <- spsample(granica, 150, type = "regular")
plot(regularny)

4.2.3 Typ próbowania | Losowy

W losowym typie próbkowania każda lokalizacja ma takie samo prawdopodobieństwo wystąpienia. Dodatkowo, każdy punkt jest losowany niezależnie od pozostałych.

set.seed(301)
losowy <- spsample(granica, 150, type = "random")
plot(losowy)

4.2.4 Typ próbowania | Losowy stratyfikowany

Losowy stratyfikowany typ próbkowania polega na podzieleniu analizowanego obszaru na regularne komórki, a następnie dla każdej komórki losowana jest lokalizacja punktu.

set.seed(125)
strat <- spsample(granica, 150, type = "stratified")
plot(strat)

4.2.5 Typ próbowania | Preferencyjny

W preferencyjnym typie próbkowania istnieją obszary, które z jakieś powodu (np. specyficzne wartości analizowanej cechy) są znacznie częściej opróbkowane niż inne.

set.seed(425)
pref <- spsample(granica, 150, type = "clustered", nclusters = 30)
plot(pref)

4.3 Dane lokalnie odstające

Wstępne sprawdzenie poprawności współrzędnych można wykonać poprzez wizualizację danych przestrzennych za pomocą funkcji plot().

plot(granica)
plot(punkty, add = TRUE)

Dane lokalnie odstające oznaczają nietypowe przestrzennie wartości danej cechy. Inaczej mówiąc, może to być niska wartość otoczona wysokimi wartościami lub też wysoka wartość otoczona niskimi wartościami. Dane lokalnie ostające mogą oznaczać zarówno błąd w danych albo wpływ innego czynnika na analizowaną cechę. Przyjrzenie się wartościom analizowanej cechy można wykonać z użyciem funkcji spplot(). Na poniższym przykładzie wyświetlona jest zmienna temp oznaczająca średnią temperaturę dobową.

spplot(punkty, "temp", sp.layout = granica)

Dodatkowo, używając argumentu identify = TRUE można interaktywnie określić numer punktu (numer wiersza w tabeli atrybutów). Działanie tej funkcji polega na wybraniu punktów za pomocą lewego przycisk myszy, a po wybraniu serii punktów klawisza Esc.

spplot(punkty, "temp", identify = TRUE)

4.4 Rozgrupowanie danych

4.4.1 Rozgrupowanie danych

Istnieje szereg metod rozgrupowywania danych, między innymi jest to rozgrupowywanie komórkowe oraz rozgrupowywanie poligonowe. Celem tych metod jest nadanie wag obserwacjom w celu zapewnienia reprezentatywności przestrzennej danych. Poniżej znajdują się przykłady użycia dwóch wersji rozgrupowania komórkowego (rozgrupowywanie komórkowe I oraz rozgrupowywanie komórkowe II) i jedna wersja rozgrupowania poligonowego.

4.4.2 Rozgrupowanie danych

Na potrzeby przykładów związanych z rozgrupowaniem danych w pakiecie geostatbook znajduje się zbiór danych punkty_pref. W tym zbiorze gęściej opróbkowane są niskie wartości temperatury, co powoduje, że średnia dla całego obszaru jest znacznie niższa niż w rzeczywistości.

data(punkty_pref)
spplot(punkty_pref, "temp")

summary(punkty_pref$temp)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   8.233  12.073  14.260  14.281  15.538  25.852

4.4.3 Rozgrupowanie komórkowe I | (ang. cell declustering)

Pierwszy rodzaj rozgrupowania komórkowego polega na:

  1. Stworzeniu siatki dla badanego obszaru
  2. Policzeniu liczby obserwacji w każdym oczku siatki
  3. Nadanie wagi dla każdego punktu, zgodnie ze wzorem:

\[w'_j=\frac{\frac{1}{n_i}}{\text{liczba komorek z danymi}} \cdot n\] , gdzie \(n_i\) to liczba obserwacji w komórce, a \(n\) to łączna liczba obserwacji

4.4.4 Rozgrupowanie komórkowe I | (ang. cell declustering)

data(punkty_pref)
punkty_pref$id <- 1:nrow(punkty_pref)
spplot(punkty_pref, "id", colorkey = TRUE)

data(granica)
siatka_n <- raster(extent(gBuffer(granica, width = 500)))
res(siatka_n) <- c(1000, 1000)
siatka_n[] <- 0
proj4string(siatka_n) <- CRS(proj4string(punkty_pref))
siatka_n <- as(siatka_n, "SpatialPolygonsDataFrame")
siatka_n <- siatka_n[!is.na(siatka_n@data$layer), ]
plot(siatka_n)
plot(punkty_pref, add = TRUE)

punkty_pref$liczebnosc <- rep(0, length(punkty_pref))
siatka_nr <- aggregate(punkty_pref["liczebnosc"], by = siatka_n, FUN = length)
spplot(siatka_nr, "liczebnosc")

liczba <- over(punkty_pref, siatka_nr)
punkty_pref$waga <- ((1 / liczba$liczebnosc) / sum(!is.na(siatka_nr$liczebnosc))) * length(punkty_pref)

spplot(punkty_pref, "waga")

srednia_arytmetyczna <- mean(punkty_pref$temp)
srednia_wazona_c1 <- mean(punkty_pref$temp * punkty_pref$waga, na.rm = TRUE)

4.4.5 Rozgrupowanie komórkowe II | (ang. cell declustering)

Drugi rodzaj rozgrupowania komórkowego polega na:

  1. Stworzeniu siatki dla badanego obszaru
  2. Wykonaniu interpolacji z użyciem funkcji krige() z pakietu gstat. W tym wypadku konieczne jest użycie argumentu nmax = 1, który przypisuje wartość najbliższej obserwacji do każdego oczka siatki.
  3. Waga dla każdego punktu nadawana jest poprzez zliczenie liczby oczek siatki dla konkretnej wartości, a następnie podzielenie tego przez liczbę oczek siatki.
data(punkty_pref)
punkty_pref$id <- 1:nrow(punkty_pref)
spplot(punkty_pref, "id", colorkey = TRUE)

data(granica)
siatka_n <- raster(extent(gBuffer(granica, width = 500)))
res(siatka_n) <- c(100, 100)
siatka_n[] <- 0
proj4string(siatka_n) <- CRS(proj4string(punkty_pref))
siatka_n <- as(siatka_n, "SpatialPointsDataFrame")
siatka_n <- siatka_n[!is.na(siatka_n@data$layer), ]
gridded(siatka_n) <- TRUE
plot(siatka_n, border = TRUE)

out <- krige(id~1, punkty_pref, siatka_n, nmax = 1)
## [inverse distance weighted interpolation]
spplot(out, "var1.pred")

df <- as.data.frame(table(out[[1]]))
df$waga <- df$Freq / sum(df$Freq)
punkty_pref <- sp::merge(punkty_pref, df, by.x = "id", by.y = "Var1")
summary(punkty_pref$waga)
##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
## 0.000336 0.001386 0.002520 0.003802 0.004326 0.043935
spplot(out, "var1.pred", sp.layout = list("sp.points", punkty_pref))

# spplot(punkty_pref, "waga")

srednia_arytmetyczna <- mean(punkty_pref$temp)
srednia_wazona_c2 <- sum(punkty_pref$temp * punkty_pref$waga, na.rm = TRUE)

4.4.6 Rozgrupowanie poligonowe | (ang. polygon declustering)

Rozgrupowanie poligonowe polega na zastosowaniu jednej z metod triangulacji, np. poligonów Voronoi’a:

  1. Dla każdego punktu określany jest poligon.
  2. Wyliczana jest powierzchnia poligonu.
  3. Waga każdego punktu wyliczana jest poprzez podzielenie powierzchni indywidualnych przez powierzchnię całego obszaru, a następnie pomnożenie przez liczbę punktów.

\[w'_j=\frac{area_j}{\sum_{j=1}^{n}area_j} \cdot n\] , gdzie \(area_j\) powierzchnia dla wybranej obserwacji, a \(n\) to łączna liczba obserwacji

4.4.7 Rozgrupowanie poligonowe | (ang. polygon declustering)

data(punkty_pref)
punkty_pref$id <- 1:nrow(punkty_pref)
spplot(punkty_pref, "id", colorkey = TRUE)

v <- voronoi(punkty_pref)
plot(punkty_pref, cex = 0.2, col = "red")
plot(v, add = TRUE)

v$pow <- area(v)
v$waga <- v$pow / sum(v$pow) * length(punkty_pref)

punkty_pref <- sp::merge(punkty_pref, v[c("id", "waga")], by = "id")
spplot(punkty_pref, "waga")

srednia_arytmetyczna <- mean(punkty_pref$temp, na.rm = TRUE)
srednia_wazona_p <- mean(punkty_pref$temp * punkty_pref$waga, na.rm = TRUE)

4.4.8 Porównanie metod rozgrupowania

Średnia wartość temperatury dla badanego obszaru wynosiła 15,59 stopni Celsjusza, jednak w preferencyjnej próbie ta wartość wynosiła 14,28 stopni Celsjusza. Porównując trzy zastosowane metody próbkowania warto zauważyć, że najbliższy wynik uzyskano korzystając z pierwszego rodzaju rozgrupowania komórkowego, która nieznacznie zaniżyła rzeczywistą wartość temperatury. Rozgrupowanie komórkowe II oraz rozgrupowanie poligonowe były w tym przypadku mniej dokładne, wyraźnie zawyżając wartości temperatury.

Średnia arytmetyczna
Populacja 15.5913
Próba 14.2815
Rozgrupowanie komórkowe I 15.3789
Rozgrupowanie komórkowe II 16.0734
Rozgrupowanie poligonowe 16.2511

W przypadku metod rozgrupowania należy jednak pamiętać, że ich wynik zależy od szeregu wprowadzonych parametrów, w szczególności granic badanego obszaru oraz zastosowanej wielkości oczka siatki.