0
$\begingroup$

I am having some issues with Kriging in R, and I was looking for some idea where I am going wrong.

From what I can tell, I done a decent job removing the trend, and I believe my transformed data is more or less isotropic. So I cannot figure out why my prediction SE is so high (see plot below). There are few places where I think I may have tripped up...

  1. I didn't understand the R-code, and I just put some argument in the wrong place.
  2. I am miss informed, and I am not Kriging inappropriately.
  3. I took the log transformation of the response and that is somehow messed up my mapping of the SE.

Thanks for the help! enter image description here enter image description hereenter image description here

enter image description here

library(geoR)
library(gstat)
library(scatterplot3d)
library(rgl)
library(sp)
library(lattice)
library(akima)
library(grid)

############ Begin Functions ###############
# This function is written by Dr. Jon Graham 
#(http://www.math.umt.edu/graham/stat544/)

rotateaxis <- function(xy, theta)
{
  pimult <- (theta * 2 * pi)/360
  newx <- c(cos(pimult), sin(pimult))
  newy <- c( - sin(pimult), cos(pimult))
  XY <- as.matrix(xy) %*% cbind(newx, newy)
  as.data.frame(XY)
}

# function written by Dr. Jon Graham (http://www.math.umt.edu/graham/stat544/)

rose <- function(dat.var,critgam){
  numcases <- length(dat.var$dist)
  numdirec <- length(unique(dat.var$dir.hor))
  rose.len <- rep(0,numdirec)
  rose.azi <- rep(0,numdirec)
  j <- 1
  for (i in 1:(numcases-1)){
    prod <- (dat.var$gamma[i] - critgam)*(dat.var$gamma[i+1] - critgam)
    if (dat.var$dir.hor[i] != dat.var$dir.hor[i+1]) prod <- 0
    if (j > 1){
      if (rose.azi[j-1] == dat.var$dir.hor[i]) prod <- 0
    }
    if (prod < 0){
      rose.len[j] <- dat.var$dist[i] + (critgam - dat.var$gamma[i])*
        (dat.var$dist[i+1] - dat.var$dist[i])/
        (dat.var$gamma[i+1] - dat.var$gamma[i])
      rose.azi[j] <- as.numeric(unique(dat.var$dir.hor)[j])
      j <- j + 1
    }
  }
  rose.azi <- rose.azi*pi/180
  dscale <- max(rose.len)
  rose.len <- rose.len/dscale
  plot(0,0,type="n",axes=F,xlim=c(-1.1,1.1),ylim=c(-1.1,1.1),
       xlab="E-W Direction",ylab="",cex.lab=1.6)
  for (j in 1:(length(rose.len))){
    segments(0,0,(rose.len[j]*sin(rose.azi[j])),(rose.len[j]*cos(rose.azi[j])))
    text((rose.len[j]*sin(rose.azi[j]) + (.1*sin(rose.azi[j]))),
         (rose.len[j]*cos(rose.azi[j]) + (.1*cos(rose.azi[j]))),
         (format(round(rose.len[j]*dscale,2))),cex=1.5)
    segments(0,0,-(rose.len[j]*sin(rose.azi[j])),
             -(rose.len[j]*cos(rose.azi[j])))
  }
  maxx <- max(rose.len*sin(rose.azi))
  text(-(maxx+.25),0,"N-S Direction",srt=90,crt=90,cex=1.6)
  title("Rose Diagram",srt=0,crt=0,cex.main=2.0)
  print("Rose Lengths")
  rose.len*dscale
}

## written by Dr. Jon Graham, website: http://www.math.umt.edu/graham/stat544/panel.gamma0.r

panelgamma0 <- function(x,y,gamma0=gamma0,span=2/3, ...){
  lofit <- loess.smooth(x,y,span=span)
  panel.xyplot(x,y,...)
  panel.xyplot(lofit$x,lofit$y,type="l")
  nump <- round(length(lofit$x)/2)
  dist0 <- approx(lofit$y[1:nump],lofit$x[1:nump],xout=gamma0)$y
  panel.segments(0,gamma0,dist0,gamma0)
  panel.segments(dist0,0,dist0,gamma0)
  parusr <- par()$usr
  panel.text(max(x),min(y),paste("d0=",format(round(dist0,4)),sep=""),adj=1)
}

############### End Functions ##########

### Read data ###

air <-structure(list(easting = c(256L, 256L, 256L, 275L, 275L, 275L, 
                                 275L, 275L, 275L, 275L, 275L, 295L, 295L, 295L, 295L, 295L, 295L, 
                                 295L, 295L, 295L, 295L, 295L, 295L, 314L, 314L, 314L, 314L, 314L, 
                                 314L, 314L, 314L, 314L, 314L, 314L, 314L, 314L, 314L, 314L, 314L, 
                                 314L, 333L, 333L, 333L, 333L, 333L, 333L, 333L, 333L, 333L, 333L, 
                                 333L, 333L, 333L, 333L, 333L, 333L, 333L, 333L, 333L, 333L, 333L, 
                                 333L, 353L, 353L, 353L, 353L, 353L, 353L, 353L, 353L, 353L, 353L, 
                                 353L, 353L, 353L, 353L, 353L, 353L, 353L, 353L, 353L, 353L, 353L, 
                                 353L, 353L, 353L, 353L, 372L, 372L, 372L, 372L, 372L, 372L, 372L, 
                                 372L, 372L, 372L, 372L, 372L, 372L, 372L, 372L, 372L, 372L, 372L, 
                                 372L, 372L, 372L, 372L, 372L, 372L, 372L, 392L, 392L, 392L, 392L, 
                                 392L, 392L, 392L, 392L, 392L, 392L, 392L, 392L, 392L, 392L, 392L, 
                                 392L, 392L, 392L, 392L, 392L, 392L, 392L, 392L, 392L, 392L, 392L, 
                                 392L, 392L, 392L, 392L, 411L, 411L, 411L, 411L, 411L, 411L, 411L, 
                                 411L, 411L, 411L, 411L, 411L, 411L, 411L, 411L, 411L, 411L, 411L, 
                                 411L, 411L, 411L, 411L, 411L, 411L, 411L, 411L, 411L, 411L, 411L, 
                                 411L, 430L, 430L, 430L, 430L, 430L, 430L, 430L, 430L, 430L, 430L, 
                                 430L, 430L, 430L, 430L, 430L, 430L, 430L, 430L, 430L, 430L, 430L, 
                                 430L, 430L, 430L, 430L, 430L, 430L, 450L, 450L, 450L, 450L, 450L, 
                                 450L, 450L, 450L, 450L, 450L, 450L, 450L, 450L, 450L, 450L, 450L, 
                                 450L, 450L, 450L, 450L, 450L, 450L, 450L, 450L, 450L, 450L, 450L, 
                                 450L, 450L, 450L, 469L, 469L, 469L, 469L, 469L, 469L, 469L, 469L, 
                                 469L, 469L, 469L, 469L, 469L, 469L, 469L, 469L, 469L, 469L, 469L, 
                                 469L, 469L, 469L, 469L, 469L, 469L, 489L, 489L, 489L, 489L, 489L, 
                                 489L, 489L, 489L, 489L, 489L, 489L, 489L, 489L, 489L, 489L, 489L, 
                                 489L, 489L, 489L, 489L, 489L, 489L, 489L, 489L, 508L, 508L, 508L, 
                                 508L, 508L, 508L, 508L, 508L, 508L, 508L, 508L, 508L, 508L, 508L, 
                                 508L, 508L, 508L, 508L, 508L, 508L, 508L, 508L, 508L, 508L, 508L, 
                                 508L, 527L, 527L, 527L, 527L, 527L, 527L, 527L, 527L, 527L, 527L, 
                                 527L, 527L, 527L, 527L, 527L, 527L, 527L, 527L, 527L, 527L, 527L, 
                                 527L, 527L, 527L, 527L, 527L, 527L, 547L, 547L, 547L, 547L, 547L, 
                                 547L, 547L, 547L, 547L, 547L, 547L, 547L, 547L, 547L, 547L, 547L, 
                                 547L, 547L, 547L, 547L, 547L, 547L, 547L, 547L, 547L, 566L, 566L, 
                                 566L, 566L, 566L, 566L, 566L, 566L, 566L, 566L, 566L, 566L, 566L, 
                                 566L, 566L, 566L, 566L, 566L, 566L, 586L, 586L, 586L, 586L, 586L, 
                                 586L, 586L, 586L, 586L, 586L, 586L, 586L, 586L, 586L, 586L, 586L, 
                                 586L, 586L, 586L, 586L, 586L, 586L, 586L, 605L, 605L, 605L, 605L, 
                                 605L, 605L, 605L, 605L, 605L, 605L, 605L, 605L, 605L, 605L, 605L, 
                                 605L, 605L, 605L, 605L, 605L, 605L, 605L, 605L, 605L, 605L, 624L, 
                                 624L, 624L, 624L, 624L, 624L, 624L, 624L, 624L, 624L, 624L, 624L, 
                                 624L, 624L, 624L, 624L, 624L, 624L, 624L, 644L, 644L, 644L, 644L, 
                                 644L, 644L, 644L, 644L, 644L, 644L, 663L, 663L, 663L, 663L, 663L, 
                                 663L, 663L, 683L, 683L, 683L, 683L, 683L, 683L, 683L, 683L, 683L, 
                                 683L, 683L, 702L, 702L, 702L, 702L, 702L, 702L, 702L, 702L, 721L, 
                                 721L, 721L, 721L, 721L, 721L, 721L, 741L, 741L, 741L, 741L, 741L, 
                                 741L, 760L, 760L, 760L, 760L, 760L), northing = c(258L, 276L, 
                                                                                   293L, 258L, 276L, 293L, 311L, 329L, 346L, 539L, 557L, 276L, 293L, 
                                                                                   311L, 329L, 346L, 434L, 452L, 469L, 487L, 504L, 539L, 575L, 135L, 
                                                                                   241L, 258L, 276L, 293L, 311L, 329L, 399L, 416L, 434L, 452L, 469L, 
                                                                                   504L, 522L, 539L, 575L, 592L, 100L, 118L, 135L, 153L, 223L, 241L, 
                                                                                   258L, 276L, 293L, 311L, 381L, 399L, 416L, 434L, 469L, 487L, 504L, 
                                                                                   522L, 539L, 557L, 575L, 592L, 65L, 82L, 100L, 118L, 135L, 153L, 
                                                                                   170L, 223L, 241L, 258L, 276L, 293L, 311L, 381L, 399L, 416L, 434L, 
                                                                                   452L, 469L, 487L, 504L, 522L, 539L, 557L, 575L, 65L, 82L, 100L, 
                                                                                   135L, 153L, 170L, 223L, 241L, 258L, 276L, 293L, 311L, 329L, 364L, 
                                                                                   381L, 399L, 416L, 434L, 452L, 469L, 487L, 504L, 522L, 539L, 698L, 
                                                                                   65L, 82L, 118L, 135L, 153L, 170L, 188L, 205L, 223L, 241L, 258L, 
                                                                                   276L, 293L, 311L, 329L, 346L, 364L, 381L, 399L, 416L, 434L, 452L, 
                                                                                   487L, 504L, 522L, 539L, 575L, 592L, 610L, 645L, 65L, 82L, 100L, 
                                                                                   135L, 153L, 170L, 188L, 205L, 223L, 241L, 258L, 276L, 293L, 311L, 
                                                                                   329L, 346L, 364L, 381L, 399L, 416L, 434L, 487L, 504L, 522L, 539L, 
                                                                                   557L, 575L, 627L, 663L, 680L, 65L, 82L, 118L, 135L, 153L, 170L, 
                                                                                   188L, 205L, 223L, 241L, 258L, 276L, 293L, 311L, 329L, 346L, 364L, 
                                                                                   381L, 399L, 452L, 469L, 487L, 504L, 522L, 539L, 557L, 680L, 65L, 
                                                                                   100L, 118L, 135L, 153L, 170L, 188L, 205L, 223L, 241L, 258L, 276L, 
                                                                                   293L, 311L, 329L, 346L, 364L, 381L, 399L, 416L, 434L, 452L, 469L, 
                                                                                   487L, 504L, 522L, 539L, 610L, 645L, 663L, 82L, 153L, 223L, 241L, 
                                                                                   258L, 276L, 293L, 311L, 329L, 346L, 364L, 381L, 399L, 416L, 452L, 
                                                                                   469L, 487L, 504L, 522L, 539L, 557L, 575L, 592L, 680L, 698L, 82L, 
                                                                                   100L, 258L, 276L, 293L, 311L, 329L, 346L, 364L, 381L, 399L, 416L, 
                                                                                   434L, 452L, 469L, 487L, 504L, 522L, 557L, 575L, 645L, 663L, 680L, 
                                                                                   698L, 135L, 258L, 276L, 293L, 311L, 329L, 346L, 364L, 381L, 399L, 
                                                                                   416L, 452L, 469L, 487L, 504L, 522L, 539L, 557L, 575L, 592L, 610L, 
                                                                                   627L, 645L, 663L, 680L, 698L, 135L, 205L, 223L, 258L, 276L, 311L, 
                                                                                   329L, 346L, 364L, 381L, 399L, 416L, 434L, 452L, 469L, 487L, 522L, 
                                                                                   539L, 557L, 575L, 592L, 610L, 627L, 645L, 663L, 680L, 698L, 188L, 
                                                                                   205L, 258L, 276L, 293L, 311L, 329L, 346L, 364L, 381L, 399L, 416L, 
                                                                                   434L, 469L, 487L, 522L, 539L, 557L, 575L, 610L, 627L, 645L, 663L, 
                                                                                   680L, 698L, 258L, 276L, 311L, 329L, 346L, 364L, 381L, 399L, 416L, 
                                                                                   434L, 452L, 469L, 539L, 557L, 627L, 645L, 663L, 680L, 698L, 170L, 
                                                                                   205L, 223L, 258L, 276L, 311L, 346L, 364L, 381L, 399L, 416L, 434L, 
                                                                                   452L, 469L, 487L, 504L, 522L, 539L, 592L, 610L, 627L, 645L, 698L, 
                                                                                   223L, 241L, 258L, 276L, 293L, 311L, 329L, 346L, 364L, 381L, 399L, 
                                                                                   416L, 434L, 452L, 469L, 487L, 504L, 522L, 539L, 557L, 575L, 592L, 
                                                                                   645L, 680L, 698L, 153L, 223L, 241L, 258L, 293L, 311L, 346L, 364L, 
                                                                                   381L, 399L, 416L, 469L, 487L, 504L, 522L, 539L, 557L, 575L, 698L, 
                                                                                   153L, 311L, 329L, 346L, 364L, 381L, 399L, 469L, 522L, 539L, 153L, 
                                                                                   311L, 364L, 381L, 504L, 610L, 645L, 293L, 311L, 329L, 346L, 364L, 
                                                                                   381L, 399L, 416L, 434L, 504L, 539L, 346L, 364L, 381L, 399L, 452L, 
                                                                                   522L, 557L, 627L, 416L, 434L, 452L, 469L, 522L, 539L, 592L, 188L, 
                                                                                   241L, 364L, 452L, 504L, 557L, 487L, 539L, 575L, 592L, 610L), 
                     ammonia = c(103L, 175L, 95L, 297L, 370L, 369L, 121L, 101L, 
                                 137L, 5L, 1L, 67L, 13L, 299L, 163L, 37L, 13L, 13L, 9L, 10L, 
                                 6L, 2L, 2L, 12L, 28L, 230L, 112L, 160L, 346L, 6L, 31L, 19L, 
                                 22L, 28L, 15L, 1L, 1L, 2L, 6L, 6L, 4L, 60L, 143L, 112L, 248L, 
                                 119L, 102L, 394L, 331L, 108L, 25L, 23L, 17L, 32L, 5L, 6L, 
                                 10L, 4L, 4L, 3L, 10L, 1L, 9L, 82L, 92L, 78L, 84L, 119L, 7L, 
                                 220L, 376L, 190L, 169L, 158L, 42L, 73L, 24L, 48L, 16L, 8L, 
                                 8L, 2L, 24L, 17L, 23L, 19L, 20L, 29L, 66L, 78L, 96L, 91L, 
                                 78L, 373L, 486L, 342L, 271L, 404L, 196L, 291L, 193L, 500L, 
                                 95L, 73L, 39L, 42L, 4L, 28L, 58L, 33L, 20L, 1L, 27L, 62L, 
                                 80L, 97L, 63L, 58L, 42L, 99L, 285L, 174L, 124L, 182L, 371L, 
                                 751L, 745L, 458L, 380L, 185L, 96L, 48L, 34L, 10L, 60L, 48L, 
                                 29L, 24L, 70L, 5L, 5L, 5L, 40L, 63L, 47L, 45L, 61L, 77L, 
                                 195L, 178L, 279L, 275L, 208L, 181L, 264L, 177L, 279L, 228L, 
                                 326L, 175L, 49L, 15L, 18L, 60L, 17L, 88L, 58L, 70L, 25L, 
                                 10L, 5L, 5L, 3L, 20L, 43L, 41L, 49L, 82L, 114L, 192L, 240L, 
                                 259L, 244L, 215L, 171L, 151L, 269L, 383L, 317L, 239L, 105L, 
                                 117L, 50L, 136L, 24L, 50L, 75L, 75L, 5L, 13L, 18L, 40L, 12L, 
                                 5L, 32L, 85L, 19L, 79L, 132L, 138L, 312L, 322L, 348L, 251L, 
                                 380L, 269L, 242L, 250L, 200L, 46L, 224L, 164L, 137L, 150L, 
                                 110L, 75L, 10L, 10L, 10L, 20L, 33L, 27L, 116L, 297L, 290L, 
                                 156L, 334L, 280L, 293L, 283L, 324L, 319L, 276L, 223L, 211L, 
                                 182L, 150L, 100L, 100L, 75L, 50L, 50L, 10L, 10L, 20L, 30L, 
                                 178L, 446L, 247L, 113L, 288L, 261L, 301L, 260L, 274L, 280L, 
                                 248L, 199L, 199L, 199L, 166L, 150L, 100L, 100L, 25L, 25L, 
                                 50L, 25L, 10L, 108L, 303L, 145L, 151L, 188L, 252L, 244L, 
                                 237L, 274L, 272L, 203L, 199L, 199L, 183L, 150L, 150L, 150L, 
                                 100L, 150L, 75L, 75L, 75L, 75L, 75L, 50L, 30L, 125L, 36L, 
                                 94L, 111L, 122L, 185L, 261L, 232L, 205L, 236L, 240L, 247L, 
                                 235L, 209L, 180L, 150L, 150L, 150L, 150L, 150L, 150L, 150L, 
                                 100L, 100L, 75L, 50L, 17L, 17L, 150L, 73L, 30L, 126L, 132L, 
                                 208L, 288L, 261L, 241L, 245L, 245L, 245L, 245L, 150L, 150L, 
                                 150L, 150L, 150L, 150L, 150L, 150L, 100L, 50L, 113L, 81L, 
                                 120L, 142L, 182L, 214L, 224L, 202L, 250L, 242L, 245L, 245L, 
                                 150L, 150L, 150L, 150L, 150L, 75L, 50L, 20L, 14L, 14L, 70L, 
                                 101L, 70L, 181L, 131L, 194L, 185L, 247L, 248L, 245L, 245L, 
                                 245L, 245L, 216L, 150L, 150L, 150L, 150L, 150L, 100L, 83L, 
                                 24L, 75L, 43L, 15L, 100L, 183L, 177L, 252L, 281L, 252L, 257L, 
                                 219L, 196L, 245L, 245L, 245L, 245L, 233L, 150L, 150L, 100L, 
                                 50L, 50L, 100L, 5L, 33L, 16L, 9L, 38L, 118L, 190L, 202L, 
                                 202L, 313L, 85L, 215L, 184L, 184L, 202L, 216L, 150L, 50L, 
                                 50L, 5L, 117L, 96L, 116L, 172L, 184L, 39L, 175L, 200L, 150L, 
                                 10L, 83L, 87L, 86L, 150L, 5L, 5L, 66L, 66L, 19L, 37L, 133L, 
                                 96L, 124L, 86L, 91L, 30L, 100L, 56L, 124L, 129L, 129L, 219L, 
                                 91L, 70L, 5L, 118L, 146L, 183L, 211L, 13L, 116L, 40L, 5L, 
                                 10L, 70L, 138L, 250L, 133L, 157L, 212L, 128L, 128L, 28L), 
                     nitroxides = c(3L, 3L, 1L, 3L, 57L, 20L, 54L, 35L, 4L, 6L, 
                                    2L, 2L, 3L, 4L, 3L, 5L, 11L, 6L, 3L, 1L, 3L, 2L, 1L, 36L, 
                                    1L, 22L, 108L, 395L, 376L, 1L, 37L, 3L, 1L, 9L, 3L, 1L, 2L, 
                                    2L, 25L, 10L, 1L, 30L, 6L, 201L, 17L, 8L, 4L, 96L, 275L, 
                                    80L, 5L, 4L, 38L, 11L, 8L, 3L, 5L, 5L, 5L, 3L, 38L, 15L, 
                                    2L, 147L, 9L, 1L, 9L, 50L, 4L, 18L, 36L, 15L, 96L, 145L, 
                                    7L, 4L, 5L, 41L, 6L, 7L, 8L, 4L, 9L, 13L, 8L, 7L, 2L, 1L, 
                                    1L, 5L, 7L, 1L, 11L, 57L, 26L, 125L, 18L, 169L, 255L, 36L, 
                                    21L, 49L, 40L, 30L, 24L, 32L, 1L, 20L, 13L, 11L, 7L, 26L, 
                                    15L, 15L, 15L, 74L, 3L, 49L, 10L, 49L, 34L, 32L, 207L, 86L, 
                                    242L, 557L, 358L, 136L, 159L, 84L, 24L, 13L, 15L, 50L, 25L, 
                                    27L, 11L, 5L, 5L, 2L, 9L, 24L, 20L, 5L, 30L, 3L, 10L, 2L, 
                                    34L, 21L, 21L, 26L, 28L, 24L, 69L, 328L, 325L, 163L, 52L, 
                                    10L, 9L, 3L, 1L, 45L, 30L, 28L, 2L, 2L, 2L, 24L, 30L, 2L, 
                                    2L, 14L, 5L, 1L, 151L, 4L, 24L, 43L, 14L, 78L, 43L, 63L, 
                                    163L, 378L, 64L, 284L, 184L, 108L, 2L, 2L, 11L, 11L, 15L, 
                                    18L, 2L, 5L, 2L, 4L, 5L, 43L, 85L, 30L, 98L, 7L, 12L, 53L, 
                                    49L, 74L, 106L, 109L, 124L, 139L, 503L, 283L, 2L, 25L, 25L, 
                                    3L, 5L, 20L, 1L, 5L, 40L, 1L, 2L, 2L, 2L, 1L, 10L, 54L, 58L, 
                                    94L, 121L, 39L, 139L, 24L, 180L, 307L, 34L, 161L, 36L, 12L, 
                                    3L, 2L, 1L, 1L, 1L, 2L, 1L, 1L, 15L, 14L, 2L, 10L, 117L, 
                                    175L, 88L, 38L, 58L, 17L, 225L, 541L, 88L, 108L, 8L, 5L, 
                                    15L, 15L, 2L, 1L, 12L, 12L, 15L, 2L, 24L, 14L, 1L, 68L, 114L, 
                                    46L, 40L, 45L, 71L, 60L, 111L, 59L, 15L, 2L, 100L, 15L, 14L, 
                                    2L, 12L, 149L, 24L, 7L, 2L, 2L, 2L, 5L, 19L, 63L, 2L, 27L, 
                                    14L, 85L, 48L, 26L, 55L, 71L, 106L, 112L, 4L, 9L, 1L, 2L, 
                                    5L, 2L, 1L, 36L, 36L, 24L, 43L, 9L, 12L, 2L, 24L, 2L, 63L, 
                                    4L, 3L, 97L, 83L, 3L, 8L, 11L, 17L, 35L, 24L, 5L, 9L, 1L, 
                                    7L, 1L, 2L, 15L, 38L, 5L, 2L, 24L, 5L, 24L, 2L, 58L, 58L, 
                                    43L, 7L, 7L, 32L, 33L, 6L, 43L, 5L, 5L, 3L, 76L, 33L, 7L, 
                                    24L, 24L, 5L, 12L, 12L, 2L, 9L, 8L, 35L, 30L, 7L, 20L, 6L, 
                                    1L, 6L, 4L, 1L, 1L, 1L, 1L, 2L, 22L, 1L, 24L, 2L, 1L, 2L, 
                                    2L, 58L, 16L, 32L, 29L, 5L, 4L, 3L, 17L, 11L, 95L, 93L, 7L, 
                                    38L, 7L, 76L, 76L, 5L, 65L, 5L, 2L, 5L, 5L, 1L, 7L, 12L, 
                                    2L, 23L, 6L, 4L, 2L, 3L, 64L, 16L, 5L, 7L, 2L, 5L, 2L, 24L, 
                                    413L, 38L, 2L, 38L, 2L, 4L, 47L, 7L, 3L, 41L, 4L, 5L, 1L, 
                                    5L, 2L, 7L, 2L, 3L, 5L, 15L, 5L, 1L, 27L, 30L, 2L, 8L, 3L, 
                                    4L, 37L, 11L, 1L, 5L, 2L, 35L, 1L, 3L, 3L, 1L, 2L, 2L, 1L, 
                                    42L, 1L, 3L, 3L, 3L, 2L, 1L, 4L, 2L, 14L, 1L, 1L, 2L, 1L, 
                                    1L, 2L, 1L, 12L)), .Names = c("easting", "northing", "ammonia", 
                                                                  "nitroxides"), class = "data.frame", row.names = c(NA, -496L))

air.coord <- air[c("easting", "northing")]
coordinates(air.coord) <- ~ easting + northing
coordinates(air) <- ~ easting + northing

vario.emperical <- variogram(nitroxides ~ easting + northing, data = air)
plot(vario.emperical)

scatterplot3d(air$easting, air$northing, air$nitroxides,
              main="3D Scatterplot",
              zlab="Nitroxides",
              xlab="Easting",
              ylab="Northing",
              type="h")

symbols(air$easting, air$northing, circles=air$nitroxides,
        bg = adjustcolor("blue", alpha=0.2))

trend<-lm(log(nitroxides) ~ easting + northing, data=air)
summary(trend)
plot(trend$residuals,  ylim = c(-4,4))


air3 <-as.data.frame(cbind(easting = air$easting,
                           northing = air$northing,
                           ammonia = air$ammonia,
                           resid = trend$residuals))
coordinates(air3) <- ~ easting + northing

plot(air3$northing, air3$resid, ylim = c(-4,4))
plot(air3$easting, air3$resid, ylim = c(-4,4))
hist(air3$resid)
scatterplot3d(air3$easting, air3$northing, air3$resid,
              main="3D Scatterplot Air3",
              zlab="Response",
              xlab="Easting",
              ylab="Northing",
              type="h")

### Semi-variograms ###
vario <- variogram(resid ~ easting + northing, data = air3)
vario
plot(vario, cex=1.3, pch=16, col=1, xlab="Distance", ylab="Semivariogram")

az <- seq(0, 165, 15)
vario.direction <- variogram(resid ~ easting + northing,
                             data = air3, alpha = az,
                             cutoff = 250, width = 12) 
vario.direction
plot(vario.direction)

## Initial Rose Diagram ## 
rose(vario.direction, 2.5)

### Rotate Axis ###
coord <- coordinates(air3)
air3.rot <- cbind(rotateaxis(coord, -90), resid = air3$resid)  
coordinates(air3.rot) <- ~ newx + newy
az <- seq(0, 165, 15)
vario.direction2 <- variogram(resid ~ newx + newy,
                              data = air3.rot, alpha = az,
                              cutoff = 225, width = 12)
rose(vario.direction2, 2.5)

xyplot(gamma ~ dist|as.factor(dir.hor), data = vario.direction2,
       layout=c(4,3), panel = panelgamma0, gamma0=2.1,
       pch=16, main="Interpolated Directional Ranges",col=1)

## Stretching Rose Diagram ##
ratio <- (130/97)
coords.new <- coords.aniso(coordinates(air3.rot), aniso.pars=c(0, ratio))
head(coords.new)
air4 <- data.frame(anix=coords.new[,1],aniy=coords.new[,2],
                   resid = air3.rot$resid )
coordinates(air4) <- ~anix + aniy

# Rose Diagram After Transformation # 
az <- seq(0, 165, 15)
air4.vario<- variogram(resid ~ anix + aniy, data = air4,
                       width = 12, cutoff = 250, alpha = az)
rose(air4.vario, 2.5)

air4.vario2<- variogram(resid ~ anix + aniy, data = air4,
                        width = 12, cutoff = 250)

model.vario <- fit.variogram(air4.vario2,
                             vgm(70000,"Sph",40,20000, anis= c(0,1/ratio)),
                             fit.method=2)
plot(air4.vario2, model.vario, 
     cex=1.3,main="Final Isotropic Variogram",
     pch=16, ylab = "Semivariogram")

plot(vario.emperical, model.vario, 
     cex = 1.3, main = "Final Isotropic Variogram",
     pch = 16, col = 1, ylab = "Semivariogram")

coords <- coordinates(air)
head(coords)
plot(coords, type="n", xlab="Easting", ylab="Northing",
     main= "Plot of Sample and Prediction Sites")
poly <- chull(coords)
poly <- c(poly, poly[1])
lines(coords[poly,])

poly.in <- polygrid(seq(256, 760, 12), seq(64, 700, 12), coords[poly,])
summary(poly.in)

coordinates(poly.in) <- ~ x+y
points(poly.in)
points(coords, pch=16)
krige.out <- krige(air$nitroxides ~ 1, air,
                   poly.in, model= model.vario)

northpred <- seq(min(air$northing), max(air$northing))
eastpred <- seq(min(air$easting), max(air$easting))

int.obs <- interp(air$easting,
                  air$northing,
                  air$nitroxides,
                  xo=eastpred, yo=northpred)

int.pred <- interp(coordinates(poly.in)[,1],
                   coordinates(poly.in)[,2],
                   krige.out$var1.pred,
                   xo=eastpred, yo=northpred)

int.se <- interp(coordinates(poly.in)[,1],
                 coordinates(poly.in)[,2],
                 sqrt(krige.out$var1.var),
                 xo=eastpred, yo=northpred)

zmin <- min(int.obs$z[!is.na(int.obs$z)], 
            int.pred$z[!is.na(int.pred$z)])
zmax <- max(int.obs$z[!is.na(int.obs$z)],
            int.pred$z[!is.na(int.pred$z)])

image(int.obs, ylab="Northing", xlab="Easting",
      main = "Observed Concentrations",
      zlim=c(zmin, zmax),
      col = rev(heat.colors(24)))

image(int.pred, ylab="Northing", xlab="Easting",
      main="Kriging Predicted Values",
      zlim=c(zmin, zmax), col=rev(heat.colors(24)))

image(int.se, ylab="Northing", xlab="Easting",
      main= "Kriging Standard Errors",
      col=rev(heat.colors(24)))
$\endgroup$

1 Answer 1

3
$\begingroup$

It looks like what is happening is that 3 prediction locations are identical to observation locations ((508,592), (760,592), (508,364)), so their prediction error is zero.

Also, the call to krige uses the original data air, on which no detrending has occurred, and it assumes a constant mean (~ 1), whereas the variogram model is based on the residuals in air4 in which any linear trend in $x$ and $y$ trend has been removed.

$\endgroup$
1
  • $\begingroup$ thank you for the tip. I won't be able to check on this till Friday, but I will let you know the results when I do. $\endgroup$
    – k6adams
    Apr 28, 2015 at 2:43

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.