Timothy_Handley at nps.gov wrote:> A couple weeks ago I posted a message on this topic to r-help, the response
> was that this seemed like odd behavior, and that I ought to post it to one
> of the developer lists. I posted to r-sig-mixed-models, but didn't get
any
> response. So, with good intentions, I decided to try posting once more, but
> to this more general list.
>
> The goal is (1) FYI, to make you aware of this issue, in case it is an
> error in gls (2) For my information, in case I have made an error, in the
> hope that one of you folks might be able to correct me. Thanks in advance
> for your time.
>
> The issue is in 2 parts.
>
> (A) I've used gls to fit a model with two fixed effects and a corExp
> object. By my count, this fitting process estimates 5 parameters:
> (Intercept), l10area, newx, range, and nugget. With 118 total df, there
> should be 118 - 5 = 113 residual df. However, the output from summary.gls
> reports 115 residual degrees of freedom. Is this an error in summary or
> gls, or is there an error in my count?
Those for the correlation structure are not counted for these residuals
as you can find in
nlme:::print.summary.gls
that contains the line
cat("Degrees of freedom:", dd[["N"]], "total;",
dd[["N"]] -
dd[["p"]], "residual\n")
> (B) Summary.gls reports logLik=-273.6. Using my count of 5 estimated
> parameters, the AIC should be -2*(-273.6) + 2*5 = 557.2. However,
> summary.gls reports an AIC of 559.2. If one works backwards from the
> reported AIC of 559.2, it seems that gls believes it has estimated 6
> parameters in the fitting process. Again, is this an error in gls, or an
> error on my part?
1 additional was used for the estimation of the variance. Accordingly
nlme:::logLik.gls
contains the line
attr(val, "df") <- p +
length(coef(object[["modelStruct"]])) + 1
The AICs should be comparable at least within R and if others think 5
rather than 6 should be used its still fine since the difference will be
there in all models to be compared.
Best wishes,
Uwe Ligges
> Copied from R terminal:
>
>> summary(sppl.i.ex)
> Generalized least squares fit by maximum likelihood
> Model: all.all.rch ~ l10area + newx
> Data: gtemp
> AIC BIC logLik
> 559.167 575.7911 -273.5835
>
> Correlation Structure: Exponential spatial correlation
> Formula: ~x + y | area
> Parameter estimate(s):
> range nugget
> 15.4448835 0.3741476
>
> Coefficients:
> Value Std.Error t-value p-value
> (Intercept) 7.621306 0.7648135 9.964921 0.0000
> l10area 6.332931 0.5589199 11.330659 0.0000
> newx 0.066535 0.0204417 3.254857 0.0015
>
> Correlation:
> (Intr) l10are
> l10area -0.605
> newx 0.358 -0.024
>
> Standardized residuals:
> Min Q1 Med Q3 Max
> -3.0035983 -0.5990432 -0.2226852 0.5113270 2.4444263
>
> Residual standard error: 2.820337
> Degrees of freedom: 118 total; 115 residual
>
> Tim Handley
> Fire Effects Monitor
> Santa Monica Mountains National Recreation Area
> 401 W. Hillcrest Dr.
> Thousand Oaks, CA 91360
> 805-370-2347
>
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel