2

I have a raster stack with 11 layers that I am running through a triply nested lapply(). I expect 44 layers to be returned, but I'm only seeing 4!

Reproducible example:

# Load packages
library(raster)
library(pbapply)
library(magrittr)

# Reading in the data
dhm <- readRDS(gzcon(url('http://jacksonvoelkel.com/geog597/labs/DHM.RDS'))) 

# Looping through some terrain functions
terrain_opts <- c('slope'
                  ,'TPI'
                  ,'TRI'
                  ,'roughness'
                  ,'aspect'
                  #,'flowdir'
)

# Create terrain stack
dhm_terrain <- pblapply(terrain_opts, function(t) {
  r_terrain <- pblapply(c(4, 8), function(n) {
    r_neigh <- raster::terrain(x = dhm, opt = t, unit = 'radians', neighbors = n)
    names(r_neigh) <- paste0(names(dhm), '_terrain_', t, '_nb_', n)
    return(r_neigh)
  }) %>% stack()
  return(r_terrain)
}) %>% stack()

# Attach the original dhm to terrain stack
dhm_terrain %<>% stack(dhm)

Here's where I'm running into issues:

# Apply boundaries filter to terrain stack
terrain_bound <- pblapply(nlayers(dhm_terrain), function(i) {
  t_stack <- pblapply(c("inner", "outer"), function(t) {
    d_stack <-  pblapply(c(4, 8), function(d) {
      r_bound <- raster::boundaries(x = dhm_terrain[[i]], directions = d, type = t, asNA=F, classes=T)
      names(r_bound) <- paste0(names(dhm_terrain[[i]]), '_', d, '_', t)
      return(r_bound)
      }) %>% stack() 
    return(d_stack)
    }) %>% stack()
  return(t_stack)    
  }) %>% stack()
  • My dhm_terrain has 10 layers, not 11. Were you referring to somethign else? – Spacedman Mar 14 at 23:04
  • Thank you for the heads up. I missed a line in my reprex. I've edited my example to include the 11th layer. – philiporlando Mar 14 at 23:06
4

This:

terrain_bound <- lapply(nlayers(dhm_terrain), function(i) {

should be:

terrain_bound <- lapply(1:nlayers(dhm_terrain), function(i) {

(or use pblapply - I took it out in order to see what i was inside the function and the printing was easier without the progress bar updating).

  • I feel dumb now, but thank you Spacedman! Also, I really enjoyed your datacamp series on spdep. – philiporlando Mar 14 at 23:12
  • 2
    Heh, thanks, and we've all done this so don't feel dumb. The trick is to not give up and debug smartly - all it took was message("i=",i) inside the function to see that i wasn't going anywhere. – Spacedman Mar 14 at 23:14
  • I will add the message() function to my tool belt. Many thanks! – philiporlando Mar 14 at 23:15

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.