Create semivariogram in package R gstat
Suppose I have precipitation data that were taken at four meteorological stations during the period 2004-2016. I have passed the data to the database for retrieval in R. My goal is to take the data for each day from that period, and krige reuse those values.
So now my data looks like this, each row corresponds to one of the points, and the columns are in order: lat, long and rainfall_data.
I followed this tutorial: https://rpubs.com/nabilabd/118172 to help me get started. So here is my code:
day_1 <- dbGetQuery(con, "SELECT lat, long, rainfall_data FROM schema.sample")
coordinates(day_1) <- ~lat+long
day_1.vgm <- variogram(rainfall_data~1, day_1)...
My problem starts with the last piece of code, every time I run it, all I get is a null (empty) result (as seen from RStudio). I can't even go to the next step:
day_1.fit <- fit.variogram(day_1.vgm, model=vgm(1, "Sph", 900, 1))
Because when I do this, it throws an error that reads:
Error in fit.variogram (day1.vgm, model = vgm (1, "Sph", 900, 1)): the object must be of the gstatVariogram or variogramCloud class
I know the dataset is SUPER LACKING with only 4 points, and I know this leads to some really bad results, but to what I got, so I'm sticking with it. But regardless of the size of the dataset, it should work if I don't miss something.
If I'm average in Java, then R is a completely foreign language to me (although impossible to learn) and statistics are far from my skill list (I'm an IT guy, not a statistician).
I am doing something wrong, can anyone give me directions? Please help. Thank.
EDIT: The data looks like this:
lat long rainfall_data
7.16 124.21 0.25
8.6 123.35 1
8.43 124.28 125.6
8.15 125.08 4.3
source to share
I would say that it is unwise to try to fit the variogram by just 4 points. However, if you really want it, you can do something like this: -
The error you are getting is that the object day_1.vgm
is NULL, so you need to look at the documentation variogram
. There are parameters, namely width
and cutoff
that need to be changed. For example try the following
day_1.vgm <- variogram(rainfall_data~1, day_1, width = 0.02, cutoff = 1.5)
If you look at the plot of this variogram, it looks something like this:
You are now trying to place the variogram at these points. This way you can use the command fit.variogram
as you have already used. However, be careful, look at the parameters. Let's start with a simple spherical model.
day_1.fit <- fit.variogram(day_1.vgm, model=vgm("Sph", psill = 8000, range = 1))
You will probably get a warning about the singular approach, which is probably due to the very low score.
The installed variogram will look something like this:
You can change the fit options accordingly.
source to share