[R-sig-Geo] Advice on GLS on residual variogram?

[ONKELINX, Thierry]

Dear Lyndon,

Have a look at the normalised residuals when creating a variogram. Those
take the correlation structure into account. And have a look at the
parameters of the correlation structure. Are they from a similar
magnitude as you would expect from the OLS variogram? If not, rerun the
GLS with starting values for range and sill based on the OLS variogram.
I had some strange results in the past with range parameters being
smaller than the smallest distance between two points. Supplying
reasonable starting values yielded better results.

HTH,

Thierry

------------------------------------------------------------------------
----
ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek
team Biometrie & Kwaliteitszorg
Gaverstraat 4
9500 Geraardsbergen
Belgium

Research Institute for Nature and Forest
team Biometrics & Quality Assurance
Gaverstraat 4
9500 Geraardsbergen
Belgium

tel. + 32 54/436 185
Thierry.Onkelinx at inbo.be
www.inbo.be

To call in the statistician after the experiment is done may be no more
than asking him to perform a post-mortem examination: he may be able to
say what the experiment died of.
~ Sir Ronald Aylmer Fisher

The plural of anecdote is not data.
~ Roger Brinner

The combination of some data and an aching desire for an answer does not
ensure that a reasonable answer can be extracted from a given body of
data.
~ John Tukey
  

> -----Oorspronkelijk bericht-----
> Van: r-sig-geo-bounces at stat.math.ethz.ch 
> [mailto:r-sig-geo-bounces at stat.math.ethz.ch] Namens Lyndon Estes
> Verzonden: maandag 14 juni 2010 18:15
> Aan: r-sig-geo
> Onderwerp: [R-sig-Geo] Advice on GLS on residual variogram?
> 
> Hello,
> 
> I am attempting to implement regression kriging, following 
> methods recommended by Hengl et al. (2007; and others), for % 
> soil organic carbon over South Africa. I was hoping to ask 
> for advice as to whether the results I am getting from GLS make sense.
> 
> For background, I have 3300+ soil profiles providing A 
> horizon OC %, and I have derived 8 spatial predictors 
> including slope, solar radiation, a topographic moisture
> index, etc. These have been transformed using principal
> components analysis (in ArcGIS).
> 
> My question concerns variograms resulting from the first 
> component of the methodology, which is to find the 
> coefficients and residuals of a GLS model as follows:
> 
> 1. Use OLS to fit my model:
> 
> oc<-read.dbf('~/oc/ocdat.dbf')
> oc_2<-as.data.frame(oc[,c(2,5:6,21:27,207:225)])
> 
> gls.all<-gls(log(CTOP)~PCB1+PCB2+PCB3+PCB4+PCB5+PCB6+PCB8,data=oc_2)
> # Note: GLS is OLS if correlation structure is not specified
> 
> 2. Find the GLS coefficients (and residuals) using an
> appropriate autocorrelation structure. In this case, fitting 
> a variogram to the OLS residuals suggested a spherical 
> autocorrelation structure:
> 
> gls.all.update<-update(gls.all,correlation=corSpher(form=~X+Y,
> nugget=TRUE))
> 
> Step 2 took a long time to complete, given my dataset--an 
> overnight run (not sure how many hours though) using 64-bit 
> R2.10.1 on a MacBook Pro with a 2.4 GHZ processor and 4 GB of 
> RAM. I am not sure that the results make sense, however, as 
> the GLS shows greater autocorrelation in the residuals then 
> the original OLS residuals. The following produced the
> illustration of the plotted residual variograms posted here 
> (http://sites.google.com/site/ldemisc/variogram):
> 
> # Fit variograms to residuals of OLS and GLS
> oc.var<-variogram(residuals(gls.all)~1,oc_2)
> oc.var.update<-variogram(residuals(gls.all.update)~1,oc_2)
> 
> # Create variograms plots
> oc.var.update.pl<-plot(oc.var.update,main="GLS residuals") 
> oc.var.pl<-plot(oc.var,main="OLS residuals")
> 
> # Display variograms side-by-side
> oc.var.pl$x.limits<-oc.var.update.pl$x.limits
> oc.var.pl$y.limits<-oc.var.update.pl$y.limits
> print(oc.var.pl,split=c(1,1,2,1),more=TRUE)
> print(oc.var.update.pl,split=c(2,1,2,1),more=FALSE)
> 
> Does it seem sensible that GLS residuals show a stronger 
> degree of spatial autocorrelation (with no sign of a sill) 
> than OLS? For comparison with another spatially 
> autocorrelated dataset, I used the meuse dataset (following 
> an example from Hengl 2009) with the models:
> 
> data(meuse)
> coordinates(meuse)=~x+y
> meu.ols<-gls(log(zinc)~dist.m+ffreq,meuse)
> meu.gls<-update(meu.ols,correlation=corExp(form=~x+y))
> 
> Plotting the variograms:
> 
> zinc.vgmOLS<-variogram(residuals(meu.ols)~1,meuse)
> ols.vgm.pl<-plot(zinc.vgmOLS,main="OLS plot")
> 
> zinc.vgm.GLS<-variogram(residuals(meu.gls)~1,meuse)
> gls.vgm.pl<-plot(zinc.vgm.GLS,main="GLS plot")
> 
> # Display variograms side-by-side
> ols.vgm.pl$x.limits<-gls.vgm.pl$x.limits
> ols.vgm.pl$y.limits<-gls.vgm.pl$y.limits
> print(ols.vgm.pl,split=c(1,1,2,1),more=TRUE)
> print(gls.vgm.pl,split=c(2,1,2,1),more=FALSE)
> 
> In contrast, this shows nearly identical spatial 
> autocorrelation patterns for OLS and GLS.
> 
> I would greatly appreciate any advice regarding the 
> (seemingly) unusual variogram results, and/or clearing up of 
> any misunderstandings I might have.
> 
> Thanks in advance.
> 
> Regards, Lyndon
> 
> 
> Hengl, T., G.B.M. Heuvelink, and D.G. Rossiter. 2007. About
> regression-kriging: From equations to case studies. Computers 
> & Geosciences 33, no. 10 (October): 1301-1315.
> 
> Hengl, T. 2009. A Practical Guide to Geostatistical Mapping.
> Luxembourg: Office for Official Publications of the European 
> Communities.
> 
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> 

Druk dit bericht a.u.b. niet onnodig af.
Please do not print this message unnecessarily.

Dit bericht en eventuele bijlagen geven enkel de visie van de schrijver weer 
en binden het INBO onder geen enkel beding, zolang dit bericht niet bevestigd is
door een geldig ondertekend document. The views expressed in  this message 
and any annex are purely those of the writer and may not be regarded as stating 
an official position of INBO, as long as the message is not confirmed by a duly 
signed document.