Saving results in automap R package for time series data

Question

I'm planning to interpolate monthly rainfall data from 1980 and 2010. I'm using following code:

require(gstat)
require(sp)
require(automap)
require(rgdal)
require(raster)

for(i in 1:ncol(data@data)) {
  kriging_new <- autoKrige(data[,i],data,grid)
  plot(kriging_new)
 }

I want to save predicted values (kriging estimates), kriging prediction variance, kriging standard deviation for each i and save it as a img file. Can anyone help me to set up this within the loop.

Answer 1

You can use either rgdal or raster (with an additional step of coersion) to export the prediction, prediction variance or prediction standard deviation as rasters. The sums of squares is stored as a vector where the fit and experimental variogram models are data.frame objects. You will have to attribute a column and output them into a single data.frame or write each one independently using something like write.csv.

Here is a worked example that writes results to the current working directory.

library(sp)
library(automap)
library(rgdal)
library(raster)

setwd("D:/TMP")

data(meuse)
coordinates(meuse) =~ x+y
data(meuse.grid)
gridded(meuse.grid) =~ x+y

vars <- c("cadmium", "copper", "lead", "zinc")
sum.squares <- vector()
var.model <- data.frame() 
  for(i in vars) {
    kriging_new <- autoKrige(meuse@data[,i]~1, meuse, meuse.grid)
    sum.squares <- append(sum.squares, kriging_new$sserr)
    kriging_new$var_model <- data.frame(y=i,kriging_new$var_model)
    var.model <- rbind(var.model, kriging_new$var_model)
    writeGDAL(kriging_new$krige_output["var1.pred"], 
      paste(paste(i, "pred",  sep="_"), "img", sep="."))   
    writeGDAL(kriging_new$krige_output["var1.var"], 
      paste(paste(i, "var", sep="_"), "img", sep="."))   
   } 

var.model

names(sum.squares) <- vars
  print(sum.squares)

r <- raster("cadmium_pred.img")
plot(r)