Ahhhh, the lightbulb just went off! That makes perfect sense. Thanks! Ben
On Fri, Jan 17, 2020 at 10:45 AM Bede-Fazekas ?kos <bfalevlist at gmail.com> wrote:
Hello Ben,
I think you are absolutely right about "raster's implementation of `[[`
is different than base R". But I don't agree on your interpretation on
the logical indexing of rasters ('first logical index is used'). It is
not the first one. The logical vector is coerced to integer, and c(TRUE,
TRUE, TRUE) is treated as c(1, 1, 1), while c(TRUE, FALSE, TRUE) is
treated as c(1, 0, 1). This is why it cause a 'not a valid subset' error
(there is no sense of searching for the 0th layer of the RasterStack).
If the first logical index was used, c(TRUE, TRUE, TRUE) and c(TRUE,
FALSE, TRUE) would give exactly the same result, i.e a RasterLayer
created from the first layer of the RasterStack.
Have a nice weekend,
?kos
2020.01.17. 15:40 keltez?ssel, Ben Tupper ?rta:
Hi,
Thanks for this. I think the root of my muddle is the mish-mash model
of a raster that I have in my head. Depending upon what is most
convenient, I sometimes view a raster as a multi-dimensional array and
at other times as a list of single layer matrices. If we step back
from logical indexing and use integers it is easier to identify
raster's slight variation on base R recursive indexing. The example
below uses integer indices on a list and on a raster.
Back to logical indexing, in a way it makes perfect sense as just the
first logical index is used; just as if() does. But what is different
is that it uses that first logical index for each element in the index
vector. That's why logo[[c(TRUE, TRUE, TRUE)]] yields red.1, red.2
and red.3.
Thanks again and cheers,
Ben
library(raster)
x <- list(red = "R", green = "G", blue = "B")
logo <- raster::brick(system.file("external/rlogo.grd", package="raster"))
# `[` index a list, get a list back (with NULL for not found)
x[c(1,3)]
# $red
# [1] "R"
#
# $green
# [1] "G"
# `[` index a raster, get a matrix back (or vector for single layer)
logo[c(1,3)]
# red green blue
# [1,] 255 255 255
# [2,] 255 255 255
# `[[` recursive index of a list fails
x[[c(1,3)]]
# Error in x[[c(1, 3)]] : subscript out of bounds
# `[[` index of a raster yields raster
# so raster's implementation of `[[` is different than base R
logo[[c(1,3)]]
# class : RasterStack
# dimensions : 77, 101, 7777, 2 (nrow, ncol, ncell, nlayers)
# resolution : 1, 1 (x, y)
# extent : 0, 101, 0, 77 (xmin, xmax, ymin, ymax)
# crs : +proj=merc +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
# names : red, blue
# min values : 0, 0
# max values : 255, 255
On Fri, Jan 17, 2020 at 7:38 AM Bede-Fazekas ?kos <bfalevlist at gmail.com> wrote:
Dear Ben, Although I cannot answer your question on why logical subsetting was not implemented in package raster, there is a very easy workaround: logo[[(1:nlayers(logo))[c(TRUE, TRUE, TRUE)]]] logo[[(1:nlayers(logo))[c(TRUE, FALSE, TRUE)]]] Also note that in case of lists '[[' does recursive indexing, and this type of logical indexing you are asking about works only with '['. HTH, ?kos Bede-Fazekas Hungarian Academy of Sciences 2020.01.16. 17:50 keltez?ssel, Ben Tupper ?rta:
Hi, I usually avoid using logical indexes with multilayer rasters (stacks and bricks), but I have never understood why indexing ala `[[` with logicals isn't supported. Below is an example that shows how it doesn't work like the typical indexing with logicals. I'm not asking to have logical indices considered (although it would be handy), but instead I really just want to understand why it is the way it is. I scanned over the introductory vignette (https://rspatial.org/raster) and issues (https://github.com/rspatial/raster/issues) but found nothing there. Thanks and cheers, Ben ### START library(raster) logo <- raster::brick(system.file("external/rlogo.grd", package="raster")) # red-red-red logo[[c(TRUE, TRUE, TRUE)]] # class : RasterStack # dimensions : 77, 101, 7777, 3 (nrow, ncol, ncell, nlayers) # resolution : 1, 1 (x, y) # extent : 0, 101, 0, 77 (xmin, xmax, ymin, ymax) # crs : +proj=merc +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 # names : red.1, red.2, red.3 # min values : 0, 0, 0 # max values : 255, 255, 255 # red-red logo[[c(TRUE, TRUE)]] # class : RasterStack # dimensions : 77, 101, 7777, 2 (nrow, ncol, ncell, nlayers) # resolution : 1, 1 (x, y) # extent : 0, 101, 0, 77 (xmin, xmax, ymin, ymax) # crs : +proj=merc +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 # names : red.1, red.2 # min values : 0, 0 # max values : 255, 255 # red logo[[c(TRUE)]] # class : RasterLayer # band : 1 (of 3 bands) # dimensions : 77, 101, 7777 (nrow, ncol, ncell) # resolution : 1, 1 (x, y) # extent : 0, 101, 0, 77 (xmin, xmax, ymin, ymax) # crs : +proj=merc +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 # source : /usr/lib64/R/library/raster/external/rlogo.grd # names : red # values : 0, 255 (min, max) logo[[c(TRUE, FALSE, TRUE)]] #Error in .local(x, ...) : not a valid subset #sessionInfo() # R version 3.5.1 (2018-07-02) # Platform: x86_64-redhat-linux-gnu (64-bit) # Running under: CentOS Linux 7 (Core) # # Matrix products: default # BLAS/LAPACK: /usr/lib64/R/lib/libRblas.so # # locale: # [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8 # [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8 LC_PAPER=en_US.UTF-8 LC_NAME=C # [9] LC_ADDRESS=C LC_TELEPHONE=C LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C # # attached base packages: # [1] stats graphics grDevices utils datasets methods base # # other attached packages: # [1] raster_3.0-7 sp_1.3-2 # # loaded via a namespace (and not attached): # [1] compiler_3.5.1 rgdal_1.4-8 tools_3.5.1 yaml_2.2.0 Rcpp_1.0.3 codetools_0.2-15 # [7] grid_3.5.1 lattice_0.20-35 ### END
_______________________________________________ R-sig-Geo mailing list R-sig-Geo at r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Ben Tupper Bigelow Laboratory for Ocean Science West Boothbay Harbor, Maine http://www.bigelow.org/ https://eco.bigelow.org