Agregowanie punktów do siatki za pomocą R.

14

Mam pytanie dotyczące agregacji przestrzennej w R. Próbuję agregować punktowy zestaw danych do siatki. Nie jestem jednak pewien, jak to zrobić, ponieważ mam niewielkie doświadczenie w tego typu sprawach. Miałem nadzieję, że ktoś z was może mieć jakieś przydatne wskazówki / możliwe rozwiązanie.

Moim punktem obserwacyjnym jest zestaw danych zawierający dane georeferencyjne dotyczące konfliktów w Afryce (patrz www.acleddata.com). Punkty są georeferencyjne ze współrzędnymi szerokości i długości geograficznej i zawierają dane dotyczące rodzaju i czasu zdarzenia. To, co chcę zrobić, to zsumować te punkty do siatki 1x1 stopnia.

Zatem komórka siatki powinna zawierać informacje o punktach danych, jeśli zdarzenie miało miejsce w tej komórce siatki. Ostatecznym produktem tego powinna być ramka danych lub coś, co mogę wyeksportować do pliku csv, ponieważ dane są przeznaczone do wykorzystania w zestawie danych panelu do analizy statystycznej.

Do tej pory załadowałem i wydrukowałem dane oraz plik kształtu za pomocą poniższego kodu. Uważam, że powinienem użyć funkcji over z pakietu sp do agregacji, ale nie wiem jak. Mam nadzieję, że jeden z was może pomóc.

Kod, którego do tej pory użyłem, można znaleźć tutaj z odpowiednim wynikiem wizualnym tam .

Sugestie dotyczące zrobienia tego w QGIS są również mile widziane.

rok
źródło
Jest to szybka, prosta operacja, wymagająca jedynie małej arytmetyki. Ale w jakim formacie chcesz otrzymać wynik? „CSV” sugeruje jedynie, że powinna to być tabela relacyjna, ale stanowi to problem: po agregacji każda komórka potencjalnie będzie odpowiadać różnej liczbie punktów. Zazwyczaj wybierasz jedną z dwóch opcji: albo wyprowadzasz jeden rekord na punkt (w tym identyfikator jego zawierającej komórki), albo wyprowadzasz jeden rekord na komórkę i dołączasz niektóre statystyczne podsumowania zawartych w nim punktów. Czego potrzebujesz?
whuber
1
Przepraszam, nie podałem tego. Potrzebuję jednego rekordu na komórkę . Używam pliku csv do tworzenia danych panelu w formacie roku komórkowego .
końca roku

Odpowiedzi:

13

Pobrane dane zawierają pewne szczere błędy lokalizacyjne, więc pierwszą rzeczą do zrobienia jest ograniczenie współrzędnych do rozsądnych wartości:

data.df <- read.csv("f:/temp/All_Africa_1997-2011.csv", header=TRUE, sep=",",row.names=NULL)
data.df <- subset(data.df, subset=(LONGITUDE >= -180 & LATITUDE >= -90))

Obliczenie współrzędnych i identyfikatorów komórki siatki jest jedynie kwestią obcięcia liczb dziesiętnych od wartości szerokości i długości geograficznej. (Mówiąc bardziej ogólnie, w przypadku dowolnych rastrów najpierw wyśrodkuj je i skaluj, aby dopasować do rozmiaru komórki, skróć dziesiętne, a następnie przeskaluj i ponownie wyśrodkuj z powrotem do ich pierwotnej pozycji, jak pokazano w poniższym kodzie ji.) Możemy połączyć te współrzędne w unikalne identyfikatory, dołączając je do wejściowej ramki danych i zapisz rozszerzoną ramkę danych jako plik CSV. Będzie jeden rekord na punkt:

ji <- function(xy, origin=c(0,0), cellsize=c(1,1)) {
  t(apply(xy, 1, function(z) cellsize/2+origin+cellsize*(floor((z - origin)/cellsize))))
}
JI <- ji(cbind(data.df$LONGITUDE, data.df$LATITUDE))
data.df$X <- JI[, 1]
data.df$Y <- JI[, 2]
data.df$Cell <- paste(data.df$X, data.df$Y)

Zamiast tego możesz chcieć wyników podsumowujących zdarzenia w każdej komórce siatki. Aby to zilustrować, obliczmy liczby na komórkę i wyprowadzamy je, jeden rekord na komórkę:

counts <- by(data.df, data.df$Cell, function(d) c(d$X[1], d$Y[1], nrow(d)))
counts.m <- matrix(unlist(counts), nrow=3)
rownames(counts.m) <- c("X", "Y", "Count")
write.csv(as.data.frame(t(counts.m)), "f:/temp/grid.csv")

W przypadku innych podsumowań zmień functionargument w obliczeniach counts. (Ewentualnie użyj arkusza kalkulacyjnego lub oprogramowania bazy danych, aby podsumować pierwszy plik wyjściowy według identyfikatora komórki).

W celu sprawdzenia, niech map liczniki za pomocą siatki centrów zlokalizować symbole na mapach. (Punkty położone na Morzu Śródziemnym, w Europie i na Oceanie Atlantyckim mają podejrzane lokalizacje: podejrzewam, że wiele z nich wynika z pomieszania szerokości i długości geograficznej w procesie wprowadzania danych.)

count.max <- max(counts.m["Count",])
colors = sapply(counts.m["Count",], function(n) hsv(sqrt(n/count.max), .7, .7, .5))
plot(counts.m["X",] + 1/2, counts.m["Y",] + 1/2, cex=sqrt(counts.m["Count",]/100),
     pch = 19, col=colors,
     xlab="Longitude of cell center", ylab="Latitude of cell center",
     main="Event counts within one-degree grid cells")

Mapa Afryki

Ten przepływ pracy jest teraz

  • Dokładnie udokumentowane (za pomocą samego Rkodu),

  • Powtarzalne (przez ponowne uruchomienie tego kodu),

  • Rozszerzalny (poprzez modyfikację kodu w oczywisty sposób), oraz

  • Dość szybko (cała operacja zajmuje mniej niż 10 sekund na przetworzenie tych 53052 obserwacji).

Whuber
źródło
Kod jest doskonale powtarzalny. Mam jeszcze jedno pytanie. Zamiast podsumowania, jak dołączyć informacje z pliku danych wejściowych do komórki w utworzonej siatce?
końca roku
1
Nie jest to możliwe w przypadku tabeli wyjściowej , ponieważ pełna informacja dla komórek ma zmienną długość. Właściwy sposób na zapisanie, czyli pierwszą formę danych wyjściowych, którą wykazałem: jeden rekord na punkt z atrybutem identyfikatora komórki. Jeden z tych dwóch formatów - tabele na punkt i na komórkę - będzie oczekiwany przez każdy program statystyczny, którego używasz.
whuber
1
Ach ok. Rozumiem, co masz na myśli. Wystarczy utworzyć siatkę dla wszystkich komórek i połączyć ją. Dzięki za pomoc.
Horseoftheyear
3

Cóż, to czego potrzebujesz to podstawowe tak zwane „łączenie przestrzenne”, które dopasowuje do siebie dwa pliki kształtów i przydziela sumę (liczbę) do wynikowej tabeli atrybutów. Jeśli szukasz „Spatial Join in R”, znajdziesz tu wiele przykładów nawet tutaj na GIS.Stackexchange. Szybko przejrzałem Google i znalazłem na przykład ten kod opublikowany na liście mailingowej.

Jeśli chcesz uzyskać sprzężenie atrybutu przestrzennego w QGIS, wykonaj następujące czynności:

  • Zapisz kształty jako pliki .shp (polecenie writeOGR z pakietu rgdal)
  • Załaduj je do QGIS. Odtwórz siatkę wektorów za pomocą wtyczki MMQGIS (Utwórz -> Utwórz warstwę siatki) z odpowiednim skalowaniem.
  • Użyj narzędzia „Dołącz atrybuty” z menu Wektor -> Zarządzanie danymi. Wybierz atrybut swojej warstwy punktowej (może to być prosta kolumna reprezentująca wartości PRAWDA (1) lub FAŁSZ (0) dla różnych zdarzeń konfliktu).
  • Wybierz siatkę i Sumuj wszystkie wystąpienia i wykonaj. Następnie obciąłbym również twoją siatkę kształtem kontynentu afrykańskiego.

Jeśli Łączenie jakoś się nie powiedzie (nie działa dla mnie za każdym razem), trzymaj się SEXTANTE i poszukaj przybornika SAGA, który ma również bardzo dobre funkcje łączenia.

Kulik
źródło
Chociaż jest to rozwiązanie, jest to szczególnie złożone i nieefektywne, biorąc pod uwagę, że sumowanie punktów do siatki jest tylko kwestią kilku prostych operacji arytmetycznych, w których Rwyróżnia się. Korzystanie z plików kształtów, rgdalQGIS i Sextante przypomina trochę zalecenie, aby ktoś wynajął nowoczesny zautomatyzowany zakład przemysłowy w celu połączenia dwóch desek :-).
whuber
Spróbuję tego podejścia w ten weekend. W niedalekiej przyszłości mógłbym chcieć łączyć ze sobą różne pliki kształtów, aby było to przydatne. Dzięki za wkład i sugestie.
Horseoftheyear
@ whuber: To prawda, ale jeśli chcesz rozpowszechniać i stylizować swoje wyniki, oczywistym wyborem jest plik kształtu. Niemniej jednak ładny przykład R.
Curlew
W końcu spróbowałem. Problem w tym podejściu polega na tym, że sumuje wszystkie obserwacje wielokąta. Chociaż idealnie chcę zachować informacje o różnych wydarzeniach w czasie. Ale możliwe, że zrobiłem coś złego.
Horseoftheyear