Bug in surf.ls (spatial package)?
On Tue, 10 Apr 2007, ONKELINX, Thierry wrote:
-----Oorspronkelijk bericht----- Van: Prof Brian Ripley [mailto:ripley at stats.ox.ac.uk] Verzonden: dinsdag 10 april 2007 16:20 Aan: ONKELINX, Thierry CC: r-sig-geo at stat.math.ethz.ch Onderwerp: Re: Bug in surf.ls (spatial package)? On Tue, 10 Apr 2007, ONKELINX, Thierry wrote:
Dear listers, When fitting a trend surface I noticed that surf.ls()
generates other
values for the model parameters than an ordinary lm().
Since both fit,
in my opinion, the same model to the same data. So I would suppect that both would generate the same parameter values. Or am I missing something?
You are missing a great deal, including the posting guide. Who said the ...$beta were the parameters of the least-squares surface in the raw coordinates? Certainly not the help page, and this is explained in the book this software supports which the posting guide asked you to consult.
First of all, I'd like to apologize for not having every book or paper mentioned in the packages I use. We have a lot of books and papers, but not all... Secondly, ?surf.ls says that "beta: the coefficients". ?lm gives "coefficients: a named vector of coefficients". Neither says that they are the parameters for the raw data. But those for lm() are. So I, wrongly, presumed those for surf.ls() where for the raw data too.
Please do consult the book reference. In both surf.ls() in spatial and equivalent functions in gstat, the coordinates are first rescaled to [-1,1] (p. 420 in the fourth edition of MASS), because at higher orders, powers of at least the latitude term become very large. For smallish study areas, a computer cannot distinguish between I(y^3) and I(y^4) when y is in meters north or south of the Equator. If you try this, you'll see that lm() drops the coefficients for terms that become indistinguishable, but surf.ls() handles this numerical problem by rescaling.
Some sample code.
library("gstat")
data(jura)
library("spatial")
surface.surf <- surf.ls(np = 1, x = prediction.dat$Xloc, y =
prediction.dat$Yloc, z = prediction.dat$Ni) surface.lm <-
lm(Ni ~ Xloc
+ Yloc, data = prediction.dat) surface.surf$beta [1] 18.828803 4.823595 -2.375896 surface.lm$coef (Intercept) Xloc Yloc 15.5140337 2.2466673 -0.9299007 Thanks, Thierry
----------------------------------------------------------------------
-- ---- ir. Thierry Onkelinx Instituut voor natuur- en bosonderzoek / Reseach Institute
for Nature
and Forest Cel biometrie, methodologie en kwaliteitszorg / Section biometrics, methodology and quality assurance Gaverstraat 4 9500 Geraardsbergen Belgium tel. + 32 54/436 185
Thierry.Onkelinx at inbo.be
www.inbo.be Do not put your faith in what statistics say until you have
carefully
considered what they do not say. ~William W. Watt A statistical analysis, properly conducted, is a delicate dissection of uncertainties, a surgery of suppositions. ~M.J.Moroney
-- Brian D. Ripley, ripley at stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UK Fax: +44 1865 272595
_______________________________________________ R-sig-Geo mailing list R-sig-Geo at stat.math.ethz.ch https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Roger Bivand Economic Geography Section, Department of Economics, Norwegian School of Economics and Business Administration, Helleveien 30, N-5045 Bergen, Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43 e-mail: Roger.Bivand at nhh.no