Fits point process models by maximising the profile likelihood, profile pseudolikelihood, profile composite likelihood or AIC.
profilepl(s, f, ..., aic=FALSE, rbord=NULL, verbose = TRUE, fast=TRUE)An object of class "profilepl". There are methods
for plot,
as.ppm,
fitin
and
parameters for objects of this class.
The components of the object include
Best-fitting model
The data frame s
Row index of the best-fitting parameters in s
To extract the best fitting model you can also use as.ppm.
Data frame containing values of the irregular parameters over which the criterion will be computed.
Function (such as Strauss)
that generates an interpoint interaction object, given
values of the irregular parameters.
Data passed to ppm to fit the model.
Logical value indicating whether to find the parameter values
which minimise the AIC (aic=TRUE) or maximise the
profile likelihood (aic=FALSE, the default).
Radius for border correction (same for all models). If omitted, this will be computed from the interactions.
Logical value indicating whether to print progress reports.
Logical value indicating whether to use a faster, less accurate model-fitting technique when computing the profile pseudolikelihood. See Section on Speed and Accuracy.
Computation of the profile pseudolikelihood can be time-consuming.
We recommend starting with a small experiment in which
s contains only a few rows of values. This will indicate
roughly the optimal values of the parameters.
Then a full calculation using more finely
spaced values can identify the exact optimal values.
It is normal that the procedure appears to slow down at the end. During the computation of the profile pseudolikelihood, the model-fitting procedure is accelerated by omitting some calculations that are not needed for computing the pseudolikelihood. When the optimal parameter values have been identified, they are used to fit the final model in its entirety. Fitting the final model can take longer than computing the profile pseudolikelihood.
If fast=TRUE (the default), then additional shortcuts are taken
in order to accelerate the computation of the profile
log pseudolikelihood. These shortcuts mean that the values of the profile
log pseudolikelihood in the result ($prof)
may not be equal to the values that would be obtained if the model was
fitted normally. Currently this happens only for the area interaction
AreaInter. It may be wise to do a small experiment with
fast=TRUE and then a definitive calculation with fast=FALSE.
Adrian Baddeley Adrian.Baddeley@curtin.edu.au, Rolf Turner rolfturner@posteo.net and Ege Rubak rubak@math.aau.dk
The model-fitting function ppm fits point process
models to point pattern data. However,
only the ‘regular’ parameters of the model can be fitted by
ppm. The model may also depend on ‘irregular’
parameters that must be fixed in any call to ppm.
This function profilepl is a wrapper which finds the values of the
irregular parameters that give the best fit.
If aic=FALSE (the default),
the best fit is the model which maximises the
likelihood (if the models are Poisson processes) or maximises
the pseudolikelihood or logistic likelihood.
If aic=TRUE then the best fit is the model which
minimises the Akaike Information Criterion AIC.ppm.
The argument s must be a data frame whose columns contain
values of the irregular parameters over which the maximisation is
to be performed.
An irregular parameter may affect either the interpoint interaction or the spatial trend.
in a call to ppm, the argument interaction
determines the interaction between points. It is usually given
as a call to a function such as Strauss. The
arguments of this call are irregular parameters.
For example, the interaction radius parameter \(r\) of the Strauss
process, determined by the argument r
to the function Strauss, is an irregular parameter.
the form of the spatial trend is specified by a model formula
involving the names of spatial covariates. Spatial covariates may
take many forms; in particular they may be functions in the
R language (written by the user) of the form function(x,y)
or function(x,y, a, b, c, ...); then the extra named arguments
a,b,c are irregular parameters. The model formula includes
only the name of the function, without parentheses and
without its arguments. When calling ppm we would
pass the irregular parameters through the argument covfunargs.
in ppm
a covariate is permitted to be a single numerical value
(and is interpreted as a covariate which takes the same value
at all locations).
So for example ppm(X ~ pmin(x, thresh), covariates=list(thresh=3))
and ppm(X ~ I(x < thresh), covariates=list(thresh=3))
are permissible; the threshold value thresh can be
regarded as an irregular parameter.
The argument f determines the interaction
for each model to be fitted. It would typically be one of the functions
Poisson,
AreaInter,
BadGey,
DiggleGatesStibbard,
DiggleGratton,
Fiksel,
Geyer,
Hardcore,
LennardJones,
OrdThresh,
Softcore,
Strauss or
StraussHard.
Alternatively it could be a function written by the user.
Columns of s which match the names of arguments of f
will be interpreted as interaction parameters. Columns of s
which match the names of variables in the model formula
will be interpreted as spatially constant covariates. Other columns will be
interpreted as trend parameters.
The data frame s must provide values for each argument of
f, except for the optional arguments, which are those arguments of
f that have the default value NA.
To find the best fit,
each row of s will be taken in turn. The parameters
in this row are used as follows:
Parameters which match the name of arguments of f
will be interpreted as interaction parameters and passed to f,
resulting in an interaction object.
Any parameter in s which matches the name of a
variable in the model formula will be passed to ppm
as a covariate consisting of a single number.
Any remaining parameters will be interpreted as trend parameters
(irregular parameters of covariate functions)
and passed to ppm through the argument
covfunargs.
Note that the additional arguments ... may include an
argument called covfunargs which provides the values of some
of the trend parameters. In that case, the trend parameters
given in the current row of s will be combined with those given
in the argument covfunargs with conflicts resolved by
giving priority to the entries in s. This makes it possible
to fix some of the irregular parameters and vary the others.
Then ppm will be applied to the data ...
using this interaction object, along with covariates
and covfunargs.
This results in a fitted point process model.
The value of the log pseudolikelihood or AIC from this model is stored.
After all rows of s have been processed in this way, the
row giving the maximum value of log pseudolikelihood will be found.
The object returned by profilepl contains the profile
pseudolikelihood (or profile AIC) function,
the best fitting model, and other data.
It can be plotted (yielding a
plot of the log pseudolikelihood or AIC values against the irregular
parameters) or printed (yielding information about the best fitting
values of the irregular parameters).
In general, f may be any function that will return
an interaction object (object of class "interact")
that can be used in a call to ppm. Each argument of
f must be a single value.
Baddeley, A. and Turner, R. (2000) Practical maximum pseudolikelihood for spatial point patterns. Australian and New Zealand Journal of Statistics 42, 283--322.
plot.profilepl
human <- interactive()
# one irregular parameter
if(human) {
rr <- data.frame(r=seq(0.05,0.15, by=0.01))
} else {
rr <- data.frame(r=c(0.05,0.1,0.15))
}
ps <- profilepl(rr, Strauss, cells)
ps
plot(ps)
# two irregular parameters
if(human) {
rs <- expand.grid(r=seq(0.05,0.15, by=0.01),sat=1:3)
} else {
rs <- expand.grid(r=c(0.07,0.12),sat=1:2)
}
pg <- profilepl(rs, Geyer, cells)
pg
as.ppm(pg)
## more information
summary(pg)
# \donttest{
# multitype pattern with a common interaction radius
RR <- data.frame(R=seq(0.03,0.05,by=0.01))
MS <- function(R) { MultiStrauss(radii=diag(c(R,R))) }
pm <- profilepl(RR, MS, amacrine ~marks)
# }
## Easiest way to estimate a threshold
tdf <- data.frame(thresh=seq(0, 153, by=if(human) 1 else 3))
pt <- profilepl(tdf, Poisson, nztrees ~ I(x < thresh))
Run the code above in your browser using DataLab