36

3 列の行列
x、y、z があるとします。ここで、z は x と y の関数です。

これらの点の「散布図」をプロットする方法を知っています plot3d(x,y,z)

しかし、代わりにサーフェスが必要な場合は、surface3d などの他のコマンドを使用する必要があります。問題は、plot3d と同じ入力を受け入れないことです。

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

このマトリックスを取得するにはどうすればよいですか? 等高線図を使用する必要がある場合と同様に、コマンド interp を使用してみました。

この行列を計算せずに、x、y、z から直接サーフェスをプロットするにはどうすればよいですか? ポイントが多すぎると、このマトリックスが大きくなりすぎます。

乾杯

4

5 に答える 5

31

x 座標と y 座標がグリッド上にない場合は、x、y、z サーフェスを 1 つに補間する必要があります。これは、任意の地理統計パッケージ (geoR、gstat など) を使用したクリギング、または逆距離重み付けなどのより単純な手法で行うことができます。

あなたが言及した「interp」機能はakimaパッケージのものだと思います。出力行列は、入力ポイントのサイズとは無関係であることに注意してください。入力に ​​10000 ポイントがあり、必要に応じてそれを 10x10 グリッドに補間できます。デフォルトでは、akima::interp は 40x40 グリッド上でそれを行います:

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)

ランダムなデータであるため、それはとがってゴミに見えます。うまくいけば、あなたのデータはそうではありません!

于 2010-10-29T12:54:41.730 に答える
9

関数outer()を使用して生成できます。

関数 のデモをご覧くださいpersp()。これは、サーフェスの透視図を描画する基本グラフィックス関数です。

最初の例は次のとおりです。

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)

同じことが次の場合にも当てはまりますsurface3d()

require(rgl)  
surface3d(x, y, z)
于 2010-10-29T12:53:07.833 に答える
7

ラティスを使用して見ることができます。この例では、z~x,y をプロットしたいグリッドを定義しました。こんな感じです。ほとんどのコードは、ワイヤーフレーム関数を使用してプロットする 3D 形状を構築しているだけであることに注意してください。

変数 "b" と "s" は x または y です。

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))
于 2013-06-17T06:11:47.220 に答える
5

rgl素晴らしいですが、軸を正しくするには少し実験が必要です。

点がたくさんある場合は、それらからランダムなサンプルを取得して、結果のサーフェスをプロットしてみませんか。すべて同じデータからのサンプルに基づいていくつかのサーフェスを追加して、サンプリングのプロセスがデータにひどく影響しているかどうかを確認できます。

したがって、これはかなり恐ろしい関数ですが、あなたがやりたいと思っていることを実行します(ただし、サンプリングはありません)。z が高さである行列 (x, y, z) を指定すると、点と面の両方がプロットされます。制限は、(x,y) ペアごとに 1 つの z しか存在できないことです。そのため、自分自身の上にループバックするプレーンは問題を引き起こします。

は、サーフェスが作成される個々のポイントをプロットします。plot_points = Tこれは、サーフェスとポイントが実際に一致することを確認するのに役立ちます。はplot_contour = T、3D ビジュアライゼーションの下に 2D 等高線図をプロットします。それ以外の場合はグレーに設定rainbowされますが、関数を変更してカスタム パレットを指定できます。とにかくこれでうまくいきますが、整理して最適化できると確信しています。は、関数が壊れたときに関数をデバッグするために使用する多くのverbose = T出力を出力します。

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)
  }
}

オプションなしで同じジョブを実行しますが、サーフェスを既存のadd_rgl_model3dplot にオーバーレイします。

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")
}

したがって、私のアプローチは、すべてのデータでそれを実行することです (約 15k ポイントから生成されたサーフェスを簡単にプロットします)。それがうまくいかない場合は、いくつかの小さなサンプルを取り、これらの関数を使用して一度にプロットします。

于 2012-10-08T10:03:17.307 に答える
4

もう遅いかもしれませんが、Spacedman に従って、 duplicate="strip" などのオプションを試しましたか?

x=runif(1000)
y=runif(1000)
z=rnorm(1000)
s=interp(x,y,z,duplicate="strip")
surface3d(s$x,s$y,s$z,color="blue")
points3d(s)
于 2016-07-29T10:55:20.193 に答える