Interpolation von Geodaten auf der Oberfläche einer Kugel

Question

Ich habe einen Datensatz mit Lat/Lon-Koordinaten und einem entsprechenden 0/1-Wert für jede Geolokation (4 bis 200+ Datenpunkte). Jetzt möchte ich die Lücken interpolieren und Farben auf der Oberfläche des Globus basierend auf den Interpolationsergebnissen hinzufügen. Das Hauptproblem, das ich habe, ist "rund um den Globus" zu interpolieren, denn derzeit mache ich in einem Flugzeug, was offensichtlich nicht funktioniert.

Meine Daten

set.seed(41) 
n <- 5 
s <- rbind(data.frame(lon = rnorm(n, 0, 180), 
         lat = rnorm(n, 90, 180), 
         value = 0), 
      data.frame(lon = rnorm(n, 180, 180), 
         lat = rnorm(n, 90, 180), 
         value = 1)) 
s$lon <- s$lon %% 360 -180 
s$lat <- s$lat %% 180 -90 
s_old <- s 

Visualize Datenpunkte

library(sp) 
library(rgdal) 
library(scales) 
library(raster) 
library(dplyr) 

par(mfrow=c(2,1), mar=c(0,0,0,0)) 
grd <- expand.grid(lon = seq(-180,180, by = 20), 
        lat = seq(-90, 90, by=10)) 
coordinates(grd) <- ~lon + lat 
gridded(grd) <- TRUE 
plot(grd, add=F, col=grey(.8)) 

coordinates(s) = ~lon + lat 
points(s, col=s$value + 2, pch=16, cex=.6) 

bivariate in der Ebene interpolieren

Gegenwärtig wird die bivariate Spline-Interpolation direkt auf den Lat/Lon-Koordinaten unter Verwendung des akima Pakets ausgeführt. Dies funktioniert, berücksichtigt jedoch nicht, dass die Koordinaten lat/lon auf einer Kugel liegen.

nx <- 361 
ny <- 181 
xo <- seq(-180, 179, len=nx) 
yo <- seq(-90, 89, len=ny) 
xy <- as.data.frame(coordinates(s)) 
int <- akima:::interp(x = xy$lon, y = xy$lat, z = s$value, 
         extrap = T, 
         xo = xo, yo = yo, 
         nx = nx, ny=100, 
         linear = F) 
z <- int$z 
# correct for out of range interpolations values 
z[z < 0] <- 0 
z[z > 1] <- 1 

grd <- expand.grid(lon = seq(-180,180, by = 20), 
        lat = seq(-90, 90, by=10)) 
coordinates(grd) <- ~lon + lat 
gridded(grd) <- TRUE 
plot(grd, add=F, col=grey(.8)) 

## create raster image 
r <- raster(nrows=ny, ncols=nx, crs='+proj=longlat', 
      xmn=-180, xmx=180, ymn=-90, ymx=90) 
values(r) <- as.vector(z) 

# tweaking of color breaks 
colors <- alpha(colorRampPalette(c("red", "yellow", "green"))(21), .4) 
br <- seq(0.3, 0.7, len=20) 
image(xo, yo, z, add = T, col = colors, breaks=c(-.1, br, 1.1)) 
points(s, col=s$value + 2, pch=16, cex=.6) 

Offensichtlich ist dies für eine Kugel nicht funktioniert, wie die linke Seite auf die rechte Seite nicht übereinstimmt. Auf einer Kugel sollte die Interpolation nahtlos sein.

Welche Ansätze kann ich auf einer Kugel in R interpolieren?

Answer 1

Sie können Abstände zwischen Punkten und Gitter selbst berechnen und dann Ihre eigene Interpolation verwenden. Zum Beispiel ist unten eine inverse Distanzinterpolation für Ihr Datenbeispiel.

generieren Daten

library(sp) 
library(rgdal) 

# Data 
set.seed(41) 
n <- 5 
s <- rbind(data.frame(lon = rnorm(n, 0, 180), 
         lat = rnorm(n, 90, 180), 
         value = 0), 
      data.frame(lon = rnorm(n, 180, 180), 
         lat = rnorm(n, 90, 180), 
         value = 1)) 
s$lon <- s$lon %% 360 -180 
s$lat <- s$lat %% 180 -90 
s_old <- s 

erstellen Raster für Raster-Interpolation

## create raster image 
r <- raster(nrows=ny, ncols=nx, crs='+proj=longlat', 
      xmn=-180, xmx=180, ymn=-90, ymx=90) 

berechnen Abstände zwischen den Punkten und Raster

Funktion spDists in der Bibliothek sp die Großkreisentfernung verwenden, wenn Koordinaten nicht projiziert . Dies bedeutet, dass die Entfernung zwischen zwei Punkten am kürzesten ist.

# Distance between points and raster 
s.r.dists <- spDists(x = coordinates(s), y = coordinates(r), longlat = TRUE) 

Interpolation auf einer Kugel unter Verwendung von inverser Interpolation Abstand

Hier schlagen I mit Strom 2 (idp=2) unter Verwendung der klassischen inversen Abstand Interpolation interpoliert werden. Sie können die Berechnung ändern, wenn Sie eine andere Leistung oder lineare Interpolation wünschen oder wenn Sie mit einer begrenzten Anzahl von Nachbarn interpolieren möchten.

# Inverse distance interpolation using distances 
# pred = 1/dist^idp 
idp <- 2 
inv.w <- (1/(s.r.dists^idp)) 
z <- (t(inv.w) %*% matrix(s$value))/apply(inv.w, 2, sum) 

r.pred <- r 
values(r.pred) <- z 

plotten dann die Ergebnisse

# tweaking of color breaks 
colors <- alpha(colorRampPalette(c("red", "yellow", "green"))(21), .4) 
br <- seq(0.3, 0.7, len=20) 
plot(r.pred, col = colors, breaks=c(-.1, br, 1.1), legend=F) 
points(s, col=s$value + 2, pch=16, cex=.6)