Hello all
I have done quite a bit of searching including in the archives of this list and I cannot find an answer to this question. The script below calculates a sample spacetime variogram from water levels measured in wells over about 40 years. The variogram plots look fine but I dont understand the values in the timelag field of the output, it comes out as decimals, with a max value of 0.52. If I have formatted the input data correctly shouldn't the time lags be counted in months as the time values of the data points are "month year", just as the lags for the spatial coordinates are in meters? Do the decimal values in the timefield of the variogram ST output mean the time lag as a fraction of the entire time series of data? I have also formatted the input data as the full date of measurement, including the day, the time lags then result in days (I think) but with a maximum of about 190. But the total time series of data spans about 40 years. I would like to use the temporal data spanning a couple of tens of years to inform the sample variogram. I have tried various values for the tlags parameter in variogram ST but that has only furthered my confusion. I am probably missing something obvious about how variogramST works and/or how to format the temporal data correctly for input to the function. Any advice is appreciated. Geoff Rawling dropbox link to two csv files for script https://www.dropbox.com/s/w2vcelj4sr5d779/doi_wells.zip?dl=0 # script to calculate experimental spacetime variogram library(sp) library(spacetime) library(gstat) library(rgdal) library(ggplot2) library(xts) # read in data of interest, doi #read in data comma separated text files of well locs and wls doi = read.table("doi.csv", header = TRUE, sep = ",", stringsAsFactors = FALSE) wells = read.table("wells.csv", header = TRUE, sep = ",", stringsAsFactors = FALSE) # basic heat map of wl meas pts thru time d <- ggplot(doi, aes(mnth_ct, Well_SID )) + scale_fill_gradient(low="black", high="red") d + geom_bin2d(binwidth = c(3, 1)) # simplify to create index matrix for input STSDF, each unique well site gets a number, each unique time of meas gets a number # n x 2 matrix, each row is a wl meas, 1st column is well site number, second column is meas time number index <- as.matrix(doi[,c("Well_SID","T_SID")]) # sort the well locs by increasing x, then y wells <- wells[order(wells$LCC_East, wells$LCC_North),] # plot to check data locations, symbolize w well id ggplot(doi, aes(LCC_East,LCC_North)) + geom_point() ggplot(doi, aes(LCC_East,LCC_North)) + geom_text(aes(label=Well_SID)) rownames(wells) <- wells[,1] wells <- wells[,c(2,3)] wells = SpatialPoints(wells, proj4string=CRS("+proj=lcc +lat_0=0 +lon_0=-106 +lat_1=33 +lat_2=45 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs ")) # time of water level measurements # DateL is from lubridate time <- as.Date(doi[,"DateL"]) # simplify time to month and year timeYM <- as.yearmon(time) str(timeYM) summary(timeYM) # PCA trend surface residuals at spatial points and time of measurement tsr <- doi[,c("PCA.tsR","Point_ID")] # make spacetime object STSDF # STSDF(sp, time, data, index, endTime = delta(time)) SO <- STSDF(wells, time, tsr, index) timeSO <- SO@time str(timeSO) str(SO) # calc exp st variogram exp_st_vgm <- variogramST(PCA.tsR~1, SO, width = 2000, cutoff = 20000, assumeRegular=FALSE, na.omit = TRUE) plot(exp_st_vgm)# plot variograms plot(exp_st_vgm, wireframe = T) # plot variograms plot(exp_st_vgm, map=F) # plot variograms -- ********************************************* Geoffrey Rawling Senior Field Geologist NM Bureau of Geology and Mineral Resources [hidden email] 1015 Tijeras Avenue NW, Suite 200 Albuquerque, NM, 87102 505-366-2535 phone 505-999-8307 cell ********************************************* [[alternative HTML version deleted]] _______________________________________________ R-sig-Geo mailing list [hidden email] https://stat.ethz.ch/mailman/listinfo/r-sig-geo |
Free forum by Nabble | Edit this page |