[R-sig-Geo] Fitting variogram model using NLS function in R
Edzer Pebesma
edzer.pebesma at uni-muenster.de
Thu Sep 5 08:16:27 CEST 2013
-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1
fit.variogram uses a weighted fit, nls() in the way you specify it
uses a non-weighted fit. In order to get comparable results, I would
try to use the same weights in nls() as fit.variogram() uses.
On 09/04/2013 09:35 PM, Moshood Agba Bakare wrote:
> Dear Dr. Edzer, My R script is not different from the previous
> mail. The only change I made is to specify the fit.method = 1 which
> is weight least squares method going by the gstat manual. The
> script is as follows
>
> ## Fit empirical semivariogram using gstat empvar<-variogram(yield
> ~ 1,canmod.sp,cutoff= 400,width = 25)
>
> ## Spherical variogram model sph.var<-vgm(psill=130, model="Sph",
> range=65, nugget=180)
> sph.mod<-fit.variogram(empvar,model=sph.var,fit.method=1)
> print(sph.mod)
>
> The result obtain is
>
>
> print(sph.mod) model psill range 1 Nug
> 292.7856258929 0.000000000 2 Sph 81.4905260953 417.759633211
>
>
>
>
> # Define Spherical Variogram functions for NLS
>
> sph.vgram <- function(dist, range, psill, nugget){ dist <-
> dist/range nugget + psill*ifelse(dist < 1, (1.5 * dist - 0.5 *
> dist^3),1) }
>
> ## Fitting spherical with NLS
>
> fit.var <- nls(gamma~sph.vgram(dist,range,psill,nugget),data =
> empvar, start=list(psill=130, range=65, nugget=180),trace=T) The
> result obtained from nls is
>
>
> Formula: gamma ~ sph.vgram(dist, range, psill, nugget)
>
> Parameters: Estimate Std. Error t value Pr(>|t|) psill
> 90.70714234 7.54388401 12.02393 2.0404e-08 *** range
> 342.10250070 42.62964315 8.02499 2.1611e-06 *** nugget
> 278.95421782 6.70590124 41.59832 3.2326e-15 *** --- Signif.
> codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>
> Residual standard error: 9.36720447 on 13 degrees of freedom
>
> Number of iterations to convergence: 18 Achieved convergence
> tolerance: 9.79837651e-06
>
>
>
>
> On Wed, Sep 4, 2013 at 11:39 AM, Edzer Pebesma <
> edzer.pebesma at uni-muenster.de> wrote:
>
>
>
> On 09/04/2013 05:50 PM, Moshood Agba Bakare wrote:
>>>> Dear Dr. Edzer, Yes, I read about the various options of
>>>> fit.method in gstat manual. The default value of fit.method
>>>> for fitting variogram model to empirical variogram is 7. That
>>>> is, the number of pair of points divided by squares of the
>>>> distance, h. I refitted the model by setting fit.method =1
>>>> which is weighted least squares so that it would be the same
>>>> with default estimation method of nls function. The
>>>> comparison of the results differ greatly despite using the
>>>> same inital starting value. I am worried! Is there wrong in
>>>> my R script? Please help.
>
> Instead of sharing emotions, you might as well share information.
>
> You don't provide a reproducible example, and you haven't shown
> your R script of the second attempt. In any case, I haven't seen
> you specify weights for nls(), and hence don't share your worries.
>
>>>>
>>>> Thank you Moshood
>>>>
>>>>
>>>> On Wed, Sep 4, 2013 at 4:40 AM, Edzer Pebesma
>>>> <edzer.pebesma at uni-muenster.de
>>>>> wrote:
>>>>
>>>> did you read the documentation of argument fit.method in
>>>> function gstat::fit.variogram?
>>>>
>>>>
>>>> On 09/04/2013 02:00 AM, Moshood Agba Bakare wrote:
>>>>>>> Dear All, I tried to compare the result obtained by
>>>>>>> fitting the spherical variogram model using
>>>>>>> fit.variogram and nls functions. The large difference
>>>>>>> in the results is a great concern for me knowing that
>>>>>>> the two functions use Weighted Least Squares (WLS)
>>>>>>> approach for estimating parameters. The steps taken are
>>>>>>> as follows:
>>>>>>>
>>>>>>> ## Fit empirical semivariogram using gstat
>>>>>>> empvar<-variogram(yield ~ 1,canmod.sp,cutoff= 400,width
>>>>>>> = 25, Cressie=TRUE)
>>>>>>>
>>>>>>> # Fitting Spherical variogram model to sample
>>>>>>> variogram
>>>>>>>
>>>>>>> sph.var<- vgm(psill =130, model = "Sph", range = 65,
>>>>>>> nugget = 180) sph.mod<-fit.variogram(empvar, model =
>>>>>>> sph.var) print(sph.mod)
>>>>>>>
>>>>>>> The result obtain from fitting the spherical model to
>>>>>>> sample variogram is sph.mod model psill
>>>>>>> range 1 Nug 230.917736411 0.0000000000 2 Sph
>>>>>>> 108.323055319 87.6889385431
>>>>>>>
>>>>>>> The non-linear least squares (NLS) approach use by
>>>>>>> default s Gauss-newton algorithm in an iterative search
>>>>>>> process. I used the initial values obtained from the
>>>>>>> empirical variogram above (psill=130, range=65, and
>>>>>>> nugget=180) as starting values for the iterative
>>>>>>> procedure.
>>>>>>>
>>>>>>> ## Define Spherical Variogram functions for NLS
>>>>>>>
>>>>>>> sph.vgram <- function(dist, range, psill, nugget){ dist
>>>>>>> <- dist/range nugget + psill*ifelse(dist < 1, (1.5 *
>>>>>>> dist - 0.5 * dist^3),1) }
>>>>>>>
>>>>>>> ## Fitting spherical with NLS
>>>>>>>
>>>>>>> fit.var <-
>>>>>>> nls(gamma~sph.vgram(dist,range,psill,nugget),data =
>>>>>>> empvar, start=list(psill=130, range=65,
>>>>>>> nugget=180),trace=T)
>>>>>>>
>>>>>>> The result obtained from the nls fitting is
>>>>>>>
>>>>>>> Nonlinear regression model model: gamma ~
>>>>>>> sph.vgram(dist, range, psill, nugget) data: empvar
>>>>>>> psill range nugget 90.7071423 342.1025007
>>>>>>> 278.9542178 residual sum-of-squares: 1140.67875
>>>>>>>
>>>>>>> Number of iterations to convergence: 18 Achieved
>>>>>>> convergence tolerance: 9.79837651e-06
>>>>>>>
>>>>>>> Could anyone explain why large difference in the two
>>>>>>> result? Is the R script for fitting the NLS right? I am
>>>>>>> worried for having such disparity. Thanks while looking
>>>>>>> forward to reading your suggestion, comments and
>>>>>>> advice. Moshood.
>>>>>>>
>>>>>>> [[alternative HTML version deleted]]
>>>>>>>
>>>>>>> _______________________________________________
>>>>>>> R-sig-Geo mailing list R-sig-Geo at r-project.org
>>>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>>>>
>>>>
>>>>>
>>>>> _______________________________________________ R-sig-Geo
>>>>> mailing list R-sig-Geo at r-project.org
>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>>
>>>>>
>>>>
>
>>
>> _______________________________________________ R-sig-Geo mailing
>> list R-sig-Geo at r-project.org
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>>
>
- --
Edzer Pebesma
Institute for Geoinformatics (ifgi), University of Münster
Heisenbergstraße 2, 48149 Münster, Germany. Phone: +49 251
83 33081, Fax: +49 251 8339763 http://ifgi.uni-muenster.de
-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.11 (GNU/Linux)
Comment: Using GnuPG with Thunderbird - http://www.enigmail.net/
iQEcBAEBAgAGBQJSKCG7AAoJEM1OCHCtOnfxVeAH/2CV6xOVfX/4dB97IMd8FGad
PrL/rPzgETRQuleTFrqWPI968r9TCOpDHF66AuCC+XVahs1I8b3vBc320YcwTGj3
+RLtWOUzztt6PzmGlNGldQza6fqv3IVZl32eLdD0g1iEgjyndXzD414weoyRSJYv
rjvNcVjYkT234AudFyiu+yLq1+0ytC4FL9sGGOMQ8P67Cn0kq7PDFPjW8ODyBzwD
2nnh08ekbHG2zbAw9wbWfVnGH+G7i0o8hIwAxaYYMDIzf6h5SEXWwQ3rwy1W9aOc
lPcc3wpPtbYSdfrMivX9dmRdNXdgGfQDKtY6J9aIAiSRnStLSC2RM6B10sHX/CQ=
=h54I
-----END PGP SIGNATURE-----
More information about the R-sig-Geo
mailing list