Nonparametric Estimate of Intensity as Function of a Covariate
Computes a nonparametric estimate of the intensity of a point process, as a function of a (continuous) spatial covariate.
rhohat(object, covariate, ...) ## S3 method for class 'ppp' rhohat(object, covariate, ..., baseline=NULL, weights=NULL, method=c("ratio", "reweight", "transform"), horvitz=FALSE, smoother=c("kernel", "local", "decreasing", "increasing"), subset=NULL, dimyx=NULL, eps=NULL, n = 512, bw = "nrd0", adjust=1, from = NULL, to = NULL, bwref=bw, covname, confidence=0.95, positiveCI) ## S3 method for class 'quad' rhohat(object, covariate, ..., baseline=NULL, weights=NULL, method=c("ratio", "reweight", "transform"), horvitz=FALSE, smoother=c("kernel", "local", "decreasing", "increasing"), subset=NULL, dimyx=NULL, eps=NULL, n = 512, bw = "nrd0", adjust=1, from = NULL, to = NULL, bwref=bw, covname, confidence=0.95, positiveCI) ## S3 method for class 'ppm' rhohat(object, covariate, ..., weights=NULL, method=c("ratio", "reweight", "transform"), horvitz=FALSE, smoother=c("kernel", "local", "decreasing", "increasing"), subset=NULL, dimyx=NULL, eps=NULL, n = 512, bw = "nrd0", adjust=1, from = NULL, to = NULL, bwref=bw, covname, confidence=0.95, positiveCI)
object |
A point pattern (object of class |
covariate |
Either a |
weights |
Optional weights attached to the data points.
Either a numeric vector of weights for each data point,
or a pixel image (object of class |
baseline |
Optional baseline for intensity function.
A |
method |
Character string determining the smoothing method. See Details. |
horvitz |
Logical value indicating whether to use Horvitz-Thompson weights. See Details. |
smoother |
Character string determining the smoothing algorithm. See Details. |
subset |
Optional. A spatial window (object of class |
dimyx,eps |
Arguments controlling the pixel resolution at which the covariate will be evaluated. See Details. |
bw |
Smoothing bandwidth or bandwidth rule
(passed to |
adjust |
Smoothing bandwidth adjustment factor
(passed to |
n, from, to |
Arguments passed to |
bwref |
Optional. An alternative value of |
... |
Additional arguments passed to |
covname |
Optional. Character string to use as the name of the covariate. |
confidence |
Confidence level for confidence intervals. A number between 0 and 1. |
positiveCI |
Logical value.
If |
This command estimates the relationship between point process intensity and a given spatial covariate. Such a relationship is sometimes called a resource selection function (if the points are organisms and the covariate is a descriptor of habitat) or a prospectivity index (if the points are mineral deposits and the covariate is a geological variable). This command uses nonparametric methods which do not assume a particular form for the relationship.
If object
is a point pattern, and baseline
is missing or
null, this command assumes that object
is a realisation of a
point process with intensity function
lambda(u) of the form
lambda(u) = rho(Z(u))
where Z is the spatial
covariate function given by covariate
, and
rho(z) is the resource selection function
or prospectivity index.
A nonparametric estimator of the function rho(z) is computed.
If object
is a point pattern, and baseline
is given,
then the intensity function is assumed to be
lambda(u) = rho(Z(u)) * B(u)
where B(u) is the baseline intensity at location u. A nonparametric estimator of the relative intensity rho(z) is computed.
If object
is a fitted point process model, suppose X
is
the original data point pattern to which the model was fitted. Then
this command assumes X
is a realisation of a Poisson point
process with intensity function of the form
lambda(u) = rho(Z(u)) * kappa(u)
where kappa(u) is the intensity of the fitted model
object
. A nonparametric estimator of
the relative intensity rho(z) is computed.
The nonparametric estimation procedure is controlled by the
arguments smoother
, method
and horvitz
.
The argument smoother
selects the type of estimation technique.
If smoother="kernel"
(the default) or smoother="local"
,
the nonparametric estimator is a smoothing estimator
of rho(z), effectively a kind of density estimator
(Baddeley et al, 2012).
The estimated function rho(z) will be
a smooth function of z.
Confidence bands are also computed, assuming a Poisson point process.
See the section on Smooth estimates.
If smoother="increasing"
or smoother="decreasing"
,
we use the nonparametric maximum likelihood estimator
of rho(z) described by Sager (1982).
This assumes that
rho(z) is either an increasing function of z,
or a decreasing function of z.
The estimated function will be a step function,
increasing or decreasing as a function of z.
See the section on Monotone estimates.
See Baddeley (2018) for a comparison of these estimation techniques.
If the argument weights
is present, then the contribution
from each data point X[i]
to the estimate of rho is
multiplied by weights[i]
.
If the argument subset
is present, then the calculations are
performed using only the data inside this spatial region.
This technique assumes that covariate
has continuous values.
It is not applicable to covariates with categorical (factor) values
or discrete values such as small integers.
For a categorical covariate, use
intensity.quadratcount
applied to the result of
quadratcount(X, tess=covariate)
.
The argument covariate
should be a pixel image, or a function,
or one of the strings "x"
or "y"
signifying the
cartesian coordinates. It will be evaluated on a fine grid of locations,
with spatial resolution controlled by the arguments
dimyx,eps
which are passed to as.mask
.
A function value table (object of class "fv"
)
containing the estimated values of rho
(and confidence limits) for a sequence of values of Z.
Also belongs to the class "rhohat"
which has special methods for print
, plot
and predict
.
Smooth estimators of rho(z) were proposed by Baddeley and Turner (2005) and Baddeley et al (2012). Similar estimators were proposed by Guan (2008) and in the literature on relative distributions (Handcock and Morris, 1999).
The estimated function rho(z) will be a smooth function of z.
The smooth estimation procedure involves computing several density estimates
and combining them. The algorithm used to compute density estimates is
determined by smoother
:
If smoother="kernel"
,
the smoothing procedure is based on
fixed-bandwidth kernel density estimation,
performed by density.default
.
If smoother="local"
, the smoothing procedure
is based on local likelihood density estimation, performed by
locfit
.
The argument method
determines how the density estimates will be
combined to obtain an estimate of rho(z):
If method="ratio"
, then rho(z) is
estimated by the ratio of two density estimates,
The numerator is a (rescaled) density estimate obtained by
smoothing the values Z(y[i]) of the covariate
Z observed at the data points y[i]. The denominator
is a density estimate of the reference distribution of Z.
See Baddeley et al (2012), equation (8). This is similar but not
identical to an estimator proposed by Guan (2008).
If method="reweight"
, then rho(z) is
estimated by applying density estimation to the
values Z(y[i]) of the covariate
Z observed at the data points y[i],
with weights inversely proportional to the reference density of
Z.
See Baddeley et al (2012), equation (9).
If method="transform"
,
the smoothing method is variable-bandwidth kernel
smoothing, implemented by applying the Probability Integral Transform
to the covariate values, yielding values in the range 0 to 1,
then applying edge-corrected density estimation on the interval
[0,1], and back-transforming.
See Baddeley et al (2012), equation (10).
If horvitz=TRUE
, then the calculations described above
are modified by using Horvitz-Thompson weighting.
The contribution to the numerator from
each data point is weighted by the reciprocal of the
baseline value or fitted intensity value at that data point;
and a corresponding adjustment is made to the denominator.
Pointwise confidence intervals for the true value of ρ(z)
are also calculated for each z,
and will be plotted as grey shading.
The confidence intervals are derived using the central limit theorem,
based on variance calculations which assume a Poisson point process.
If positiveCI=FALSE
, the lower limit of the confidence
interval may sometimes be negative, because the confidence intervals
are based on a normal approximation to the estimate of ρ(z).
If positiveCI=TRUE
, the confidence limits are always
positive, because the confidence interval is based on a normal
approximation to the estimate of log(ρ(z)).
For consistency with earlier versions, the default is
positiveCI=FALSE
for smoother="kernel"
and positiveCI=TRUE
for smoother="local"
.
The nonparametric maximum likelihood estimator of a monotone function rho(z) was described by Sager (1982). This method assumes that rho(z) is either an increasing function of z, or a decreasing function of z. The estimated function will be a step function, increasing or decreasing as a function of z.
This estimator is chosen by specifying
smoother="increasing"
or smoother="decreasing"
.
The argument method
is ignored this case.
To compute the estimate of rho(z), the algorithm first computes several primitive step-function estimates, and then takes the maximum of these primitive functions.
If smoother="decreasing"
, each primitive step function
takes the form rho(z) = lambda when z ≤ t,
and rho(z) = 0 when z > t, where
and lambda is a primitive estimate of intensity
based on the data for Z <= t. The jump location t
will be the value of the covariate Z at one of the
data points. The primitive estimate lambda
is the average intensity (number of points divided by area)
for the region of space where the covariate value is less than
or equal to t.
If horvitz=TRUE
, then the calculations described above
are modified by using Horvitz-Thompson weighting.
The contribution to the numerator from
each data point is weighted by the reciprocal of the
baseline value or fitted intensity value at that data point;
and a corresponding adjustment is made to the denominator.
Confidence intervals are not available for the monotone estimators.
Smoothing algorithm by Adrian Baddeley Adrian.Baddeley@curtin.edu.au, Ya-Mei Chang, Yong Song, and Rolf Turner r.turner@auckland.ac.nz.
Nonparametric maximum likelihood algorithm by Adrian Baddeley Adrian.Baddeley@curtin.edu.au.
Baddeley, A., Chang, Y.-M., Song, Y. and Turner, R. (2012) Nonparametric estimation of the dependence of a point process on spatial covariates. Statistics and Its Interface 5 (2), 221–236.
Baddeley, A. and Turner, R. (2005) Modelling spatial point patterns in R. In: A. Baddeley, P. Gregori, J. Mateu, R. Stoica, and D. Stoyan, editors, Case Studies in Spatial Point Pattern Modelling, Lecture Notes in Statistics number 185. Pages 23–74. Springer-Verlag, New York, 2006. ISBN: 0-387-28311-0.
Baddeley, A. (2018) A statistical commentary on mineral prospectivity analysis. Chapter 2, pages 25–65 in Handbook of Mathematical Geosciences: Fifty Years of IAMG, edited by B.S. Daya Sagar, Q. Cheng and F.P. Agterberg. Springer, Berlin.
Guan, Y. (2008) On consistent nonparametric intensity estimation for inhomogeneous spatial point processes. Journal of the American Statistical Association 103, 1238–1247.
Handcock, M.S. and Morris, M. (1999) Relative Distribution Methods in the Social Sciences. Springer, New York.
Sager, T.W. (1982) Nonparametric maximum likelihood estimation of spatial patterns. Annals of Statistics 10, 1125–1136.
See ppm
for a parametric method for the same problem.
X <- rpoispp(function(x,y){exp(3+3*x)}) rho <- rhohat(X, "x") rho <- rhohat(X, function(x,y){x}) plot(rho) curve(exp(3+3*x), lty=3, col=2, add=TRUE) rhoB <- rhohat(X, "x", method="reweight") rhoC <- rhohat(X, "x", method="transform") rhoM <- rhohat(X, "x", smoother="increasing") plot(rhoM, add=TRUE, col=5) fit <- ppm(X, ~x) rr <- rhohat(fit, "y")
Please choose more modern alternatives, such as Google Chrome or Mozilla Firefox.