4

我正在处理 R 中的一些数据，这些数据由由三个空间维度和一个时间维度组成的四维数组组成：x、y、z、t。对于我的一些分析，我想获取一组空间坐标 x、y、z 的时间维度中的所有数据。到目前为止，我已经使用 which 函数来获取感兴趣的空间位置的索引。但是当我去获取与空间位置相对应的时间维度中的所有相关数据时，我找不到一个优雅的 R 解决方案，而不得不使用 repmat，一个移植的 MATLAB 函数。

a4d <- array(rnorm(10000), rep(10,4)) #x, y, z, t

#arbitrary set of 3d spatial indices x, y, z (here, using high values at first timepoint)
indices <- which(a4d[,,,1] > 2, arr.ind=TRUE)
str(indices)

# int [1:20, 1:3] 10 2 6 5 8 2 6 8 2 10 ...
# - attr(*, "dimnames")=List of 2
# ..$ : NULL
# ..$ : chr [1:3] "dim1" "dim2" "dim3"

#Now, I would like to use these indices to get data x, y, z for all t

#Intuitive, but invalid, syntax (also not clear what the structure of the data would be)
#a4d[indices,]

#Ugly, but working, syntax
library(pracma)

#number of timepoints
nt <- dim(a4d)[4]

#create a 4d lookup matrix
lookup <- cbind(repmat(indices, nt, 1), rep(1:nt, each=nrow(indices)))

#obtain values at each timepoint for indices x, y, z
result <- cbind(lookup, a4d[lookup])

该解决方案对于所述目的可以正常工作，但在概念上看起来很丑陋。理想情况下，我希望最后有一个二维矩阵：索引 x 时间。因此，在这种情况下，查找中有 20 个 x、y、z 坐标和 10 个时间点，一个 20 x 10 矩阵将是理想的，其中行表示每行索引（不需要保留 x、y、z , 值必须）并且每一列都是一个时间点。

在 R 中有没有好的方法来做到这一点？我玩过 do.call("[", list ... etc. 并使用 external 和 prod，但这些并没有像我希望的那样工作。

感谢您的任何建议！迈克尔

4

3 回答 3

7

我认为您正在寻找：

apply(a4d, 4, `[`, indices)

并检查我们的结果是否匹配：

result1 <- matrix(result[,5], ncol = 10)
result2 <- apply(a4d, 4, `[`, indices)
identical(result1, result2)
# [1] TRUE
于 2012-07-17T02:21:36.190 回答
1

我可能错过了一些东西，但你不只是想要a4d[indices[,1],indices[,2],indices[,3],]吗？

于 2012-07-17T01:55:58.337 回答
1

按每个维度单独访问并不像@tilo-wiklund 或我所期望的那样工作。结果不是跨 10 个时间步长的 23 行，而是跨 10 个时间步长的 23x23x23 立方体。

r.idvdim <- a4d[indices[,1],indices[,2],indices[,3],]
r.apply  <- apply(a4d, 4, `[`, indices)
r.cbind  <- matrix(a4d[lookup],ncol=nt)

dim(r.idvdim)     # [1] 23 23 23 10
dim(r.apply)      # [1] 23 10
dim(r.cbind)      # [1] 23 10
于 2012-09-17T13:26:55.383 回答