plot3d - r surface3d




R: Plotten einer 3D-Oberfläche aus x, y, z (4)

Stellen Sie sich vor, ich habe eine 3-Säulen-Matrix
x, y, z wobei z eine Funktion von x und y ist.

Ich weiß, wie man ein "Streudiagramm" dieser Punkte mit plot3d(x,y,z)

Aber wenn ich stattdessen eine Oberfläche haben möchte, muss ich andere Befehle verwenden, wie z. B. "surface3d". Das Problem ist, dass es nicht die gleichen Eingaben akzeptiert wie plot3d, für das es eine Matrix zu benötigen scheint

(nº elements of z) = (n of elements of x) * (n of elements of x)

Wie kann ich diese Matrix bekommen? Ich habe es mit dem Befehl interp versucht, so wie ich es mache, wenn ich Konturplots verwenden muss.

Wie kann ich eine Fläche direkt aus x, y, z plotten, ohne diese Matrix zu berechnen? Wenn ich zu viele Punkte hätte, wäre diese Matrix zu groß.

Prost


Sie können die Funktion outer() , um sie zu generieren.

Sehen Sie sich die Demo für die Funktion persp() , eine persp() zum Zeichnen von Perspektivzeichnungen für Oberflächen.

Hier ist ihr erstes Beispiel:

x <- seq(-10, 10, length.out = 50)  
y <- x  
rotsinc <- function(x,y) {
    sinc <- function(x) { y <- sin(x)/x ; y[is.na(y)] <- 1; y }  
    10 * sinc( sqrt(x^2+y^2) )  
}

z <- outer(x, y, rotsinc)  
persp(x, y, z)

Das gleiche gilt für die surface3d() :

require(rgl)  
surface3d(x, y, z)

Sie könnten mit Lattice suchen. In diesem Beispiel habe ich ein Raster definiert, über das ich z ~ x, y plotten möchte. Es sieht ungefähr so ​​aus. Beachten Sie, dass der Großteil des Codes nur eine 3D-Form erstellt, die ich mithilfe der Drahtmodellfunktion plotten kann.

Die Variablen "b" und "s" könnten x oder y sein.

require(lattice)

# begin generating my 3D shape
b <- seq(from=0, to=20,by=0.5)
s <- seq(from=0, to=20,by=0.5)
payoff <- expand.grid(b=b,s=s)
payoff$payoff <- payoff$b - payoff$s
payoff$payoff[payoff$payoff < -1] <- -1
# end generating my 3D shape


wireframe(payoff ~ s * b, payoff, shade = TRUE, aspect = c(1, 1),
    light.source = c(10,10,10), main = "Study 1",
    scales = list(z.ticks=5,arrows=FALSE, col="black", font=10, tck=0.5),
    screen = list(z = 40, x = -75, y = 0))

Wenn Ihre x- und y-Koordinaten nicht auf einem Raster liegen, müssen Sie Ihre x-, y-, z-Fläche auf eine Ebene interpolieren. Sie können dies mit Kriging mit einem der Geostatistik-Pakete (geor, gstat, andere) oder einfachere Techniken wie inverse Distanzgewichtung tun.

Ich nehme an, die Interp-Funktion, die Sie erwähnen, stammt aus dem Akima-Paket. Beachten Sie, dass die Ausgabematrix unabhängig von der Größe Ihrer Eingabepunkte ist. Sie könnten 10000 Punkte in Ihrer Eingabe haben und diese auf ein 10x10 Raster interpolieren, wenn Sie möchten. Standardmäßig führt akima :: interp auf ein 40x40 Raster:

> require(akima) ; require(rgl)
> x=runif(1000)
> y=runif(1000)
> z=rnorm(1000)
> s=interp(x,y,z)
> dim(s$z)
[1] 40 40
> surface3d(s$x,s$y,s$z)

Das wird stachelig und Müll wegen seiner zufälligen Daten aussehen. Hoffentlich sind deine Daten nicht!


rgl ist großartig, aber ein wenig experimentieren, um die Achsen richtig zu machen.

Wenn Sie viele Punkte haben, nehmen Sie eine Stichprobe von ihnen und zeichnen Sie dann die resultierende Fläche. Sie können mehrere Oberflächen basierend auf Stichproben aus den gleichen Daten hinzufügen, um festzustellen, ob der Stichprobenprozess Ihre Daten erheblich beeinträchtigt.

Also, hier ist eine ziemlich schreckliche Funktion, aber es tut, was ich denke, dass Sie es wollen (aber ohne das Sampling). Bei einer Matrix (x, y, z), bei der z die Höhe ist, werden sowohl die Punkte als auch eine Fläche aufgetragen. Einschränkungen bestehen darin, dass für jedes (x, y) Paar nur ein z vorhanden sein kann. So verursachen Flugzeuge, die sich selbst umkehren, Probleme.

Der plot_points = T wird die einzelnen Punkte, aus denen die Oberfläche gemacht wird, plot_points = T - dies ist nützlich, um zu überprüfen, dass die Oberfläche und die Punkte tatsächlich zusammentreffen. Das plot_contour = T wird ein 2D-Konturdiagramm unter der 3D-Visualisierung darstellen. Setzen Sie die Farbe auf " rainbow , um schöne Farben zu erhalten, bei allen anderen wird sie grau dargestellt. Sie können die Funktion jedoch ändern, um eine benutzerdefinierte Palette zu erhalten. Das macht den Trick für mich sowieso, aber ich bin mir sicher, dass es aufgeräumt und optimiert werden kann. Die verbose = T druckt eine Menge der Ausgabe aus, die ich verwende, um die Funktion zu debuggen, wenn und wenn es bricht.

plot_rgl_model_a <- function(fdata, plot_contour = T, plot_points = T, 
                             verbose = F, colour = "rainbow", smoother = F){
  ## takes a model in long form, in the format
  ## 1st column x
  ## 2nd is y,
  ## 3rd is z (height)
  ## and draws an rgl model

  ## includes a contour plot below and plots the points in blue
  ## if these are set to TRUE

  # note that x has to be ascending, followed by y
  if (verbose) print(head(fdata))

  fdata <- fdata[order(fdata[, 1], fdata[, 2]), ]
  if (verbose) print(head(fdata))
  ##
  require(reshape2)
  require(rgl)
  orig_names <- colnames(fdata)
  colnames(fdata) <- c("x", "y", "z")
  fdata <- as.data.frame(fdata)

  ## work out the min and max of x,y,z
  xlimits <- c(min(fdata$x, na.rm = T), max(fdata$x, na.rm = T))
  ylimits <- c(min(fdata$y, na.rm = T), max(fdata$y, na.rm = T))
  zlimits <- c(min(fdata$z, na.rm = T), max(fdata$z, na.rm = T))
  l <- list (x = xlimits, y = ylimits, z = zlimits)
  xyz <- do.call(expand.grid, l)
  if (verbose) print(xyz)
  x_boundaries <- xyz$x
  if (verbose) print(class(xyz$x))
  y_boundaries <- xyz$y
  if (verbose) print(class(xyz$y))
  z_boundaries <- xyz$z
  if (verbose) print(class(xyz$z))
  if (verbose) print(paste(x_boundaries, y_boundaries, z_boundaries, sep = ";"))

  # now turn fdata into a wide format for use with the rgl.surface
  fdata[, 2] <- as.character(fdata[, 2])
  fdata[, 3] <- as.character(fdata[, 3])
  #if (verbose) print(class(fdata[, 2]))
  wide_form <- dcast(fdata, y ~ x, value_var = "z")
  if (verbose) print(head(wide_form))
  wide_form_values <- as.matrix(wide_form[, 2:ncol(wide_form)])  
  if (verbose) print(wide_form_values)
  x_values <- as.numeric(colnames(wide_form[2:ncol(wide_form)]))
  y_values <- as.numeric(wide_form[, 1])
  if (verbose) print(x_values)
  if (verbose) print(y_values)
  wide_form_values <- wide_form_values[order(y_values), order(x_values)]
  wide_form_values <- as.numeric(wide_form_values)
  x_values <- x_values[order(x_values)]
  y_values <- y_values[order(y_values)]
  if (verbose) print(x_values)
  if (verbose) print(y_values)

  if (verbose) print(dim(wide_form_values))
  if (verbose) print(length(x_values))
  if (verbose) print(length(y_values))

  zlim <- range(wide_form_values)
  if (verbose) print(zlim)
  zlen <- zlim[2] - zlim[1] + 1
  if (verbose) print(zlen)

  if (colour == "rainbow"){
    colourut <- rainbow(zlen, alpha = 0)
    if (verbose) print(colourut)
    col <- colourut[ wide_form_values - zlim[1] + 1]
    # if (verbose) print(col)
  } else {
    col <- "grey"
    if (verbose) print(table(col2))
  }


  open3d()
  plot3d(x_boundaries, y_boundaries, z_boundaries, 
         box = T, col = "black",  xlab = orig_names[1], 
         ylab = orig_names[2], zlab = orig_names[3])

  rgl.surface(z = x_values,  ## these are all different because
              x = y_values,  ## of the confusing way that 
              y = wide_form_values,  ## rgl.surface works! - y is the height!
              coords = c(2,3,1),
              color = col,
              alpha = 1.0,
              lit = F,
              smooth = smoother)

  if (plot_points){
    # plot points in red just to be on the safe side!
    points3d(fdata, col = "blue")
  }

  if (plot_contour){
    # plot the plane underneath
    flat_matrix <- wide_form_values
    if (verbose) print(flat_matrix)
    y_intercept <- (zlim[2] - zlim[1]) * (-2/3) # put the flat matrix 1/2 the distance below the lower height 
    flat_matrix[which(flat_matrix != y_intercept)] <- y_intercept
    if (verbose) print(flat_matrix)

    rgl.surface(z = x_values,  ## these are all different because
                x = y_values,  ## of the confusing way that 
                y = flat_matrix,  ## rgl.surface works! - y is the height!
                coords = c(2,3,1),
                color = col,
                alpha = 1.0,
                smooth = smoother)
  }
}

Das add_rgl_model führt den gleichen Job ohne die Optionen aus, überlagert jedoch eine Oberfläche auf dem vorhandenen 3D-Plot.

add_rgl_model <- function(fdata){

  ## takes a model in long form, in the format
  ## 1st column x
  ## 2nd is y,
  ## 3rd is z (height)
  ## and draws an rgl model

  ##
  # note that x has to be ascending, followed by y
  print(head(fdata))

  fdata <- fdata[order(fdata[, 1], fdata[, 2]), ]

  print(head(fdata))
  ##
  require(reshape2)
  require(rgl)
  orig_names <- colnames(fdata)

  #print(head(fdata))
  colnames(fdata) <- c("x", "y", "z")
  fdata <- as.data.frame(fdata)

  ## work out the min and max of x,y,z
  xlimits <- c(min(fdata$x, na.rm = T), max(fdata$x, na.rm = T))
  ylimits <- c(min(fdata$y, na.rm = T), max(fdata$y, na.rm = T))
  zlimits <- c(min(fdata$z, na.rm = T), max(fdata$z, na.rm = T))
  l <- list (x = xlimits, y = ylimits, z = zlimits)
  xyz <- do.call(expand.grid, l)
  #print(xyz)
  x_boundaries <- xyz$x
  #print(class(xyz$x))
  y_boundaries <- xyz$y
  #print(class(xyz$y))
  z_boundaries <- xyz$z
  #print(class(xyz$z))

  # now turn fdata into a wide format for use with the rgl.surface
  fdata[, 2] <- as.character(fdata[, 2])
  fdata[, 3] <- as.character(fdata[, 3])
  #print(class(fdata[, 2]))
  wide_form <- dcast(fdata, y ~ x, value_var = "z")
  print(head(wide_form))
  wide_form_values <- as.matrix(wide_form[, 2:ncol(wide_form)])  
  x_values <- as.numeric(colnames(wide_form[2:ncol(wide_form)]))
  y_values <- as.numeric(wide_form[, 1])
  print(x_values)
  print(y_values)
  wide_form_values <- wide_form_values[order(y_values), order(x_values)]
  x_values <- x_values[order(x_values)]
  y_values <- y_values[order(y_values)]
  print(x_values)
  print(y_values)

  print(dim(wide_form_values))
  print(length(x_values))
  print(length(y_values))

  rgl.surface(z = x_values,  ## these are all different because
              x = y_values,  ## of the confusing way that 
              y = wide_form_values,  ## rgl.surface works!
              coords = c(2,3,1),
              alpha = .8)
  # plot points in red just to be on the safe side!
  points3d(fdata, col = "red")
}

Also mein Ansatz wäre, zu versuchen, es mit all Ihren Daten zu tun (ich plotten leicht Oberflächen erzeugt von ~ 15k Punkten). Wenn das nicht funktioniert, nimm mehrere kleinere Proben und zeichne sie alle gleichzeitig mit diesen Funktionen auf.







rgl