Univeral kriging
Hi there, good routines already exist in several R libraries for doing this. They will save you trouble in developing bug-free and efficient code. Check out gstat and maybe fields. If you still want to implement kriging manually, these scripts, which I developed for a class exercise, might be a useful starting place. They are not there permanently, probably through the end of 2009: http://www.msu.edu/~ashton/classes/866/notes/lect18/manual_krig.R http://www.msu.edu/~ashton/classes/866/notes/lect18/covmodels.R yours, Ashton
On Monday 26 October 2009 07:56:36 TANKISO THEJANE wrote:
R geo-helpers
i am new to r and would like you to help me as i am trying to go through
basic steps. I need to compute a predicted value at point (3,3) and its
variance given la set of locations and value. i want to make the code
generic for any size such that i would be able to estimate universal krige
at each point of grid (0,..,6) (0,..6).
the code
g = function(h) { return(0.1+3*(1-exp(-h/4)))}## variogram of exponential
form x = rbind(c(1,1),c(2,5),c(4,1),c(5,4))##4 observed locations
x = rbind(x,c(3,3)) # Row five is the point at which we want to estimate.
d = sqrt(x[,1]^2+x[,2]^2)##distance of the locations from the origin
s = c(1,9,5,12)## values at the corresponding 4 locations
## pairsv = pair semivariogram
pairsv = function(i,j) {
p1 = x[i,]
p2 = x[j,]
d = sqrt(sum((p1-p2)^2))
return(g(d))
}
G = matrix(0,4,4)
for(i in 1:4) for(j in 1:4) G[i,j] = pairsv(i,j)
X = matrix(0,4,3)
for(i in 1:4) for(j in 1:3) X[i,j] = d[i]^(j-1)
##Asymmetric matrix
Gu = rbind(cbind(G,X), cbind(t(X),matrix(0,3,3)))
x
gu = (1:7)*0
gu
for(i in 1:4) gu[i] = pairsv(i,5)
for(i in 1:3) gu[i+4] = d[5]^(i-1)
lu = inv(Gu) %*% gu
[[alternative HTML version deleted]]
Ashton Shortridge Associate Professor ashton at msu.edu Dept of Geography http://www.msu.edu/~ashton 235 Geography Building ph (517) 432-3561 Michigan State University fx (517) 432-1671