Plotten und Einfärben von Daten auf einem unregelmäßigen Raster

Question

Ich habe Daten in der Form (x, y, z), wobei x und y nicht auf einem regelmäßigen Raster liegen. Ich möchte eine 2D-Färbekarte dieser Daten anzeigen, wobei die Intensität (zum Beispiel Graustufen) auf die Variable z abgebildet wird. Eine offensichtliche Lösung besteht darin, (siehe unten) auf einem regelmäßigen Raster zu interpolieren,

d <- data.frame(x=runif(1e3, 0, 30), y=runif(1e3, 0, 30))
d$z = (d$x - 15)^2 + (d$y - 15)^2

library(akima)
d2 <- with(d, interp(x, y, z, xo=seq(0, 30, length = 30),
                     yo=seq(0, 30, length = 50), duplicate="mean"))

pal1 <- grey(seq(0,1,leng=500))
with(d2, image(sort(x), sort(y), z, useRaster=TRUE, col = pal1))
points(d$x, d$y, col="white", bg=grey(d$z/max(d$z)), pch=21, cex=1,lwd=0.1)

Allerdings geht hierbei die Information des initialen Gitters (Positionen der Punkte mit tatsächlichen Daten) verloren, das an bestimmten Stellen sehr fein oder grob sein könnte. Mein bevorzugter Ansatz wäre eine Delaunay-Verkachelung mit Dreiecken, die die tatsächliche Position und Dichte der ursprünglichen Datenpunkte genau darstellt.

Die Lösung sollte idealerweise

• die Erstellung der Verkachelung außerhalb der Plot-Funktion berechnen, sodass die resultierenden Polygone mit ggplot2, lattice oder Basisgrafiken geplottet werden können

• schnell sein. In meinem Beispiel in der Realität (~1e5 Punkte) kann die Berechnung der Verkachelung über deldir wirklich langsam sein.

Mit "Verkachelung" meine ich entweder Delaunay-Dreiecke oder Voronoi-Diagramme, obwohl meine Vorliebe für ersteres gilt. Dies bringt jedoch die zusätzliche Komplexität mit sich, die Farbe jedes Dreiecks basierend auf den Originaldatenpunkten zu interpolieren.

Answer 1

Hier ist eine Lösung auf Basis von dirichlet aus dem Paket maptools,

d <- data.frame(x=runif(1e3, 0, 30), y=runif(1e3, 0, 30))
d$z = (d$x - 15)^2 + (d$y - 15)^2

library(spatstat) 
library(maptools)

W <- ripras(df, shape="rectangle") 
W <- owin(c(0, 30), c(0, 30)) 
X <- as.ppp(d, W=W) 
Y <- dirichlet(X) 
Z <- as(Y, "SpatialPolygons") 
plot(Z, col=grey(d$z/max(d$z)))

Ich bin mir immer noch unsicher darüber, wie man die Polygone aus dieser SpatialPolygons-Klasse extrahiert.

Außerdem würde ich gerne erfahren, ob es einen einfachen Weg gibt, die "korrekten" Farben für die zugehörige Delaunay-Triangulation zu erzeugen.

Answer 2

Hier ist eine Gitterlösung mit deldir

d <- data.frame(x=runif(1e3, 0, 30), y=runif(1e3, 0, 30))
d$z = (d$x - 15)^2 + (d$y - 15)^2

pal1 <- grey(seq(0,1,leng=500))
library(latticeExtra)

 levelplot(z~x*y, data=d,
           panel = function(...) panel.voronoi(..., points=FALSE),
           interpolate=TRUE,
           col.regions = colorRampPalette(pal1)(1e3), cut=1e3)