rinterpolationgstat

IDW parameters in R


I want to perform IDW interpolation using R using the idw command from the gstat package. I have this data:

#settings
library(gstat)
library(dplyr)
library(sp)
library(tidyr)

id_rep <- rep(c(1,2), 20)
f <- rep(c(930,930.2), each=20)
perc <- rep(c(90, 80), each=10)
x <- sample(1:50, 40)
y <- sample(50:100, 40)
E <- runif(40)
df <- data.frame(id_rep, perc, x,y, f, E)
df_split <- split(df, list(df$id_rep, df$perc, df$f), drop = TRUE, sep="_")

#grid
x.range <- range(df$x)
y.range <- range(df$y)

grid <- expand.grid(x = seq(x.range[1], x.range[2], by=1), 
                       y = seq(y.range[1], y.range[2], by=1))
coordinates(grid) <- ~x + y

#interpolation
lst_interp_idw <- lapply(df_split, function(X) {

   coordinates(X) <- ~x + y
    E_idw <- idw(E~ 1, X, grid, idp=1, nmax=3) %>% as.data.frame()

   df_interp <- select(E_idw, x,y,E_pred=var1.pred)
   df_interp
})

  df_interp_idw <- bind_rows(lst_interp_idw, .id = "interact") %>%
  separate(interact, c("id_rep", "perc", "f"), sep = "\\_")

Now I want to perform each run with different idp and nmax parameters within certain values​ (idp from 1 to 3 by 0.5, and nmax 3 to 6 by 1) and get out a data frame with columns for each combination of idp and nmax values. I try with two for loops but it doesn't work.

EDIT the code that doesn't work is:

idp = seq(from = 1, to = 3, by = 0.5)
nmax = seq(from = 3, to = 6, by = 1)

...
for(i in idp) {
  for(j in nmax)
{ E_idw= idw(E ~ 1, X, grid, nmax = i, idp = j)
  }
} 
...

Solution

  • Here is a way how to store the result of every iteration in a list.

    #settings
    #install.packages("gstat")
    library(gstat)
    library(dplyr)
    library(sp)
    library(tidyr)
    
    id_rep <- rep(c(1,2), 20)
    f <- rep(c(930,930.2), each=20)
    perc <- rep(c(90, 80), each=10)
    x <- sample(1:50, 40)
    y <- sample(50:100, 40)
    E <- runif(40)
    df <- data.frame(id_rep, perc, x,y, f, E)
    df_split <- split(df, list(df$id_rep, df$perc, df$f), drop = TRUE, sep="_")
    
    #grid
    x.range <- range(df$x)
    y.range <- range(df$y)
    
    grid <- expand.grid(x = seq(x.range[1], x.range[2], by=1), 
                        y = seq(y.range[1], y.range[2], by=1))
    coordinates(grid) <- ~x + y
    
    # ==============================================
    # NEW function
    # ==============================================
    
    idp = seq(from = 1, to = 3, by = 0.5)
    nmax = seq(from = 3, to = 6, by = 1)
    
    #interpolation
    lst_interp_idw <- lapply(df_split, function(X) {
    
      coordinates(X) <- ~x + y
    
      df_interp <- vector(length(idp)*length(nmax), mode = "list" )
    
      k <- 0
    
      for(i in idp) {
    
        for(j in nmax) {
    
          print(paste(i, j))
    
          # Iterator
          k <- k + 1
    
          E_idw= idw(E ~ 1, X, grid, nmax = i, idp = j) %>% as.data.frame()
    
          df_interp[[k]] <- select(E_idw, x,y,E_pred=var1.pred)
    
        }
      }
    
      return(df_interp)
    })
    
    # ==============================================
    

    Some plausibility checks (lapply is applied to 8 list elements and 20 variations are calculated):

    length(lst_interp_idw) # 8
    length(lst_interp_idw[[1]]) #20
    length(lst_interp_idw[[1]]) #20
    

    It should be easy for you to adapt the last line of your code

    df_interp_idw <- bind_rows(lst_interp_idw, .id = "interact") %>%
      separate(interact, c("id_rep", "perc", "f"), sep = "\\_")
    

    to format the output in the desired format. This highly depends on how you want to present the different interpolation alternatives.