Gcom | R Documentation |
Given a point process model fitted to a point pattern dataset,
this function computes the compensator
of the nearest neighbour distance distribution function G
based on the fitted model
(as well as the usual nonparametric estimates
of G
based on the data alone).
Comparison between the nonparametric and model-compensated G
functions serves as a diagnostic for the model.
Gcom(object, r = NULL, breaks = NULL, ...,
correction = c("border", "Hanisch"),
conditional = !is.poisson(object),
restrict=FALSE,
model=NULL,
trend = ~1, interaction = Poisson(),
rbord = reach(interaction),
ppmcorrection="border",
truecoef = NULL, hi.res = NULL)
object |
Object to be analysed.
Either a fitted point process model (object of class |
r |
Optional.
Vector of values of the argument |
breaks |
This argument is for internal use only. |
correction |
Edge correction(s) to be employed in calculating the compensator.
Options are |
conditional |
Optional. Logical value indicating whether to compute the estimates for the conditional case. See Details. |
restrict |
Logical value indicating whether to compute
the restriction estimator ( |
model |
Optional. A fitted point process model (object of
class |
trend , interaction , rbord |
Optional. Arguments passed to |
... |
Extra arguments passed to |
ppmcorrection |
The |
truecoef |
Optional. Numeric vector. If present, this will be treated as
if it were the true coefficient vector of the point process model,
in calculating the diagnostic. Incompatible with |
hi.res |
Optional. List of parameters passed to |
This command provides a diagnostic for the goodness-of-fit of
a point process model fitted to a point pattern dataset.
It computes different estimates of the nearest neighbour distance
distribution function G
of the
dataset, which should be approximately equal if the model is a good
fit to the data.
The first argument, object
, is usually a fitted point process model
(object of class "ppm"
), obtained from the
model-fitting function ppm
.
For convenience, object
can also be a point pattern
(object of class "ppp"
).
In that case, a point process
model will be fitted to it,
by calling ppm
using the arguments
trend
(for the first order trend),
interaction
(for the interpoint interaction)
and rbord
(for the erosion distance in the border correction
for the pseudolikelihood). See ppm
for details
of these arguments.
The algorithm first extracts the original point pattern dataset
(to which the model was fitted) and computes the
standard nonparametric estimates of the G
function.
It then also computes the model-compensated
G
function. The different functions are returned
as columns in a data frame (of class "fv"
).
The interpretation of the columns is as follows
(ignoring edge corrections):
bord
:the nonparametric border-correction estimate of G(r)
,
\hat G(r) = \frac{\sum_i I\{ d_i \le r\} I\{ b_i > r \}}{\sum_i I\{
b_i > r\}}
where d_i
is the distance from the i
-th data point
to its nearest neighbour, and b_i
is the distance from the
i
-th data point to the boundary of the window W
.
bcom
:the model compensator of the border-correction estimate
{\bf C}\, \hat G(r) = \frac{\int \lambda(u,x) I\{ b(u) > r\} I\{ d(u,x)
\le r\}}{ 1 + \sum_i I\{ b_i > r\} }
where
\lambda(u,x)
denotes the conditional intensity
of the model at the location u
, and d(u,x)
denotes the
distance from u
to the nearest point in x
, while
b(u)
denotes the distance from u
to the boundary of the
windowW
.
han
:the nonparametric Hanisch estimate of G(r)
\hat G(r) = \frac{D(r)}{D(\infty)}
where
D(r) = \sum_i
\frac{ I\{x_i \in W_{\ominus d_i}\} I\{d_i \le r\}
}{
\mbox{area}(W_{\ominus d_i})
}
in which W_{\ominus r}
denotes the erosion of the window
W
by a distance r
.
hcom
:the corresponding model-compensated function
{\bf C} \, G(r) = \int_W
\frac{
\lambda(u,x) I(u \in W_{\ominus d(u)}) I(d(u) \le r)
}{
\hat D(\infty) \mbox{area}(W_{\ominus d(u)}) + 1
}
where d(u) = d(u, x)
is the (‘empty space’)
distance from location u
to the nearest point of x
.
If the fitted model is a Poisson point process, then the formulae above are exactly what is computed. If the fitted model is not Poisson, the formulae above are modified slightly to handle edge effects.
The modification is determined by the arguments
conditional
and restrict
.
The value of conditional
defaults to FALSE
for Poisson models
and TRUE
for non-Poisson models.
If conditional=FALSE
then the formulae above are not modified.
If conditional=TRUE
, then the algorithm calculates
the restriction estimator if restrict=TRUE
,
and calculates the reweighting estimator if restrict=FALSE
.
See Appendix E of Baddeley, Rubak and \Moller (2011).
See also spatstat.options('eroded.intensity')
.
Thus, by default, the reweighting estimator is computed
for non-Poisson models.
The border-corrected and Hanisch-corrected estimates of G(r)
are
approximately unbiased estimates of the G
-function,
assuming the point process is
stationary. The model-compensated functions are unbiased estimates
of the mean value of the corresponding nonparametric estimate,
assuming the model is true. Thus, if the model is a good fit, the mean value
of the difference between the nonparametric and model-compensated
estimates is approximately zero.
To compute the difference between the nonparametric and model-compensated
functions, use Gres
.
A function value table (object of class "fv"
),
essentially a data frame of function values.
There is a plot method for this class. See fv.object
.
, \ege and Jesper \Moller.
Baddeley, A., Rubak, E. and \Moller, J. (2011) Score, pseudo-score and residual diagnostics for spatial point process models. Statistical Science 26, 613–646.
Related functions:
Gest
,
Gres
.
Alternative functions:
Kcom
,
psstA
,
psstG
,
psst
.
Model fitting: ppm
.
fit0 <- ppm(cells ~1) # uniform Poisson
G0 <- Gcom(fit0)
G0
plot(G0)
# uniform Poisson is clearly not correct
# Hanisch estimates only
plot(Gcom(fit0), cbind(han, hcom) ~ r)
fit1 <- ppm(cells, ~1, Strauss(0.08))
plot(Gcom(fit1), cbind(han, hcom) ~ r)
# Try adjusting interaction distance
fit2 <- update(fit1, Strauss(0.10))
plot(Gcom(fit2), cbind(han, hcom) ~ r)
G3 <- Gcom(cells, interaction=Strauss(0.12))
plot(G3, cbind(han, hcom) ~ r)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.