R: Calculate sill, range and nugget from a raster object

I need to calculate the sill, range and nugget from a raster layer. I have explored gstat, usdm packages where one can create variogram however I couln't find a function which given a raster layer will estimate these parameters.In most of the functions these parameters have to be defined eg. krigging.

I have raster data layers for different heights which looks similar to 在这里输入图像描述

I would like get the sill, nugget and range from the parameters of semivariogram fitted to these data layers to create a plot similar to this: 在这里输入图像描述

The original data layers are available here as a multiband tiff. Here is a figure from this paper which further illustrates the concept.

在这里输入图像描述


Using gstat, here is an example:

library(raster)
library(gstat)
demo(meuse, ask = FALSE, echo = FALSE)
set.seed(131) # make random numbers reproducible
# add some noise with .1 variance
meuse.grid$dist = meuse.grid$dist + rnorm(nrow(meuse.grid), sd=sqrt(.1))
r = raster(meuse.grid["dist"])
v = variogram(dist~1, as(r, "SpatialPixelsDataFrame"))

(f = fit.variogram(v, vgm("Sph")))
#   model      psill    range
# 1   Nug 0.09035948    0.000
# 2   Sph 0.06709838 1216.737

f$psill[2] # sill
# [1] 0.06709838

f$range[2] # range
# [1] 1216.737

f$psill[1] # nugget
# [1] 0.09035948

Plug in your own raster for r , and it should work. Change the Sph to fit another variogram model, try plot(v,f) to verify the plot.


This is just a guess. This is how I estimate semi variance

where n is the number of layers which their mean is less than the total mean. m is the total mean across all the layers. r is the mean of each layer that fell below the total mean.

s <- stack("old_gap_.tif")
m <- cellStats(mean(s), stat="mean", na.rm=T) # 0.5620522
r <- m[m < 0.5620522]
sem <- 1/53 * (0.5620522 - r)^2
plot(sem, r)
链接地址: http://www.djcxy.com/p/90790.html

上一篇: 如何调用ListView.renderRow()的highlightRow?

下一篇: R:从栅格对象计算窗台,范围和块