0

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)
  }
} 
...
majom
  • 7,863
  • 7
  • 55
  • 88
Lince202
  • 143
  • 10
  • Can you include the code that doesn't work as it may only be something minor. – steveb Aug 07 '16 at 14:53
  • @steveb i edit the post with the code – Lince202 Aug 07 '16 at 15:02
  • Taking into account how wou wrote the code, you should try to make `df_interp` a list and include it in the loop to ensure that the information of each iteration is stored (in an element of this list specified by an iterator) and not overwritten. – majom Aug 07 '16 at 21:21
  • @majom i understand your advice... but I need a dataframe in output not a list... maybe it would be better using apply family function to do that... what do you think about? – Lince202 Aug 07 '16 at 22:03
  • Your approach was Ok. Just saving the results in list had to be done, See my answer below, – majom Aug 07 '16 at 22:26

1 Answers1

0

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.

majom
  • 7,863
  • 7
  • 55
  • 88
  • your code is ok!! only `df_interp_idw <- bind_rows(lst_interp_idw, .id = "interact") %>% separate(interact, c("id_rep", "perc", "f"), sep = "\\_")` doesn't work... i want in output a big dataframe like this: `Source: local data frame [...x 8] id_rep perc f x y E_pred idp nmax` – Lince202 Aug 08 '16 at 15:14
  • It depends how you want to prepare the data. "a big dataframe" is not sufficient as you have multiple alternatives how the same values have been imputed and you for sure want to be able to distinguish those. However, once you know how to prepare the data, it should be easy for you to adapt your command to the list of lists. – majom Aug 08 '16 at 15:17
  • yes... indeed i add two columns to your code: `... E_idw= idw(E ~ 1, X, grid, nmax = i, idp = j) %>% as.data.frame() df_i <- select(E_idw, x,y,E_pred=var1.pred) df_i$Idp <- rep(i, nrow(df_i)) df_i$Nmax <- rep(j, nrow(df_i)) df_interp[[k]] <- df_i ...` so i can have memory of each value of idp and nmax... now i want to bind_rows all these lists and have a single dataframe with 8 columns (id_rep perc f x y E_pred idp nmax)... i hope it's clear what i mean... – Lince202 Aug 08 '16 at 15:31
  • Did you try to use `as.data.table()`instead of às.data.frame`and then apply the command `rbindlist()`: http://stackoverflow.com/questions/15673550/why-is-rbindlist-better-than-rbind – majom Aug 08 '16 at 15:59