markcorr: Mark Correlation Function
In spatstat.explore: Exploratory Data Analysis for the 'spatstat' Family

markcorr

R Documentation

Mark Correlation Function

Description

Estimate the marked correlation function of a marked point pattern.

Usage

markcorr(X, f = function(m1, m2) { m1 * m2}, r=NULL,
         correction=c("isotropic", "Ripley", "translate"),
         method="density", ..., weights=NULL,
         f1=NULL, normalise=TRUE, fargs=NULL, internal=NULL)

Arguments

`X`	The observed point pattern. An object of class `"ppp"` or something acceptable to `as.ppp`.
`f`	Optional. Test function `f` used in the definition of the mark correlation function. An R function with at least two arguments. There is a sensible default.
`r`	Optional. Numeric vector. The values of the argument `r` at which the mark correlation function `k_f(r)` should be evaluated. There is a sensible default.
`correction`	A character vector containing any selection of the options `"isotropic"`, `"Ripley"`, `"translate"`, `"translation"`, `"none"` or `"best"`. It specifies the edge correction(s) to be applied. Alternatively `correction="all"` selects all options.
`method`	A character vector indicating the user's choice of density estimation technique to be used. Options are `"density"`, `"loess"`, `"sm"` and `"smrep"`.
`...`	Arguments passed to the density estimation routine (`density`, `loess` or `sm.density`) selected by `method`.
`weights`	Optional. Numeric weights for each data point in `X`. A numeric vector, a pixel image, or a `function(x,y)`. Alternatively, an `expression` to be evaluated to yield the weights; the expression may involve the variables `x,y,marks` representing the coordinates and marks of`X`.
`f1`	An alternative to `f`. If this argument is given, then `f` is assumed to take the form `f(u,v)=f_1(u)f_1(v)`.
`normalise`	If `normalise=FALSE`, compute only the numerator of the expression for the mark correlation.
`fargs`	Optional. A list of extra arguments to be passed to the function `f` or `f1`.
`internal`	Do not use this argument.

Details

By default, this command calculates an estimate of Stoyan's mark correlation k_{mm}(r) for the point pattern.

Alternatively if the argument f or f1 is given, then it calculates Stoyan's generalised mark correlation k_f(r) with test function f.

Theoretical definitions are as follows (see Stoyan and Stoyan (1994, p. 262)):

For a point process X with numeric marks, Stoyan's mark correlation function k_{mm}(r), is

k_{mm}(r) = \frac{E_{0u}[M(0) M(u)]}{E[M,M']}

where E_{0u} denotes the conditional expectation given that there are points of the process at the locations 0 and u separated by a distance r, and where M(0),M(u) denote the marks attached to these two points. On the denominator, M,M' are random marks drawn independently from the marginal distribution of marks, and E is the usual expectation.
For a multitype point process X, the mark correlation is

k_{mm}(r) = \frac{P_{0u}[M(0) M(u)]}{P[M = M']}

where P and P_{0u} denote the probability and conditional probability.
The generalised mark correlation function k_f(r) of a marked point process X, with test function f, is

k_f(r) = \frac{E_{0u}[f(M(0),M(u))]}{E[f(M,M')]}

The test function f is any function f(m_1,m_2) with two arguments which are possible marks of the pattern, and which returns a nonnegative real value. Common choices of f are: for continuous nonnegative real-valued marks,

f(m_1,m_2) = m_1 m_2

for discrete marks (multitype point patterns),

f(m_1,m_2) = 1(m_1 = m_2)

and for marks taking values in [0,2\pi),

f(m_1,m_2) = \sin(m_1 - m_2)

Note that k_f(r) is not a “correlation” in the usual statistical sense. It can take any nonnegative real value. The value 1 suggests “lack of correlation”: if the marks attached to the points of X are independent and identically distributed, then k_f(r) \equiv 1. The interpretation of values larger or smaller than 1 depends on the choice of function f.

The argument X must be a point pattern (object of class "ppp") or any data that are acceptable to as.ppp. It must be a marked point pattern.

The argument f determines the function to be applied to pairs of marks. It has a sensible default, which depends on the kind of marks in X. If the marks are numeric values, then f <- function(m1, m2) { m1 * m2} computes the product of two marks. If the marks are a factor (i.e. if X is a multitype point pattern) then f <- function(m1, m2) { m1 == m2} yields the value 1 when the two marks are equal, and 0 when they are unequal. These are the conventional definitions for numerical marks and multitype points respectively.

The argument f may be specified by the user. It must be an R function, accepting two arguments m1 and m2 which are vectors of equal length containing mark values (of the same type as the marks of X). (It may also take additional arguments, passed through fargs). It must return a vector of numeric values of the same length as m1 and m2. The values must be non-negative, and NA values are not permitted.

Alternatively the user may specify the argument f1 instead of f. This indicates that the test function f should take the form f(u,v)=f_1(u)f_1(v) where f_1(u) is given by the argument f1. The argument f1 should be an R function with at least one argument. (It may also take additional arguments, passed through fargs).

The argument r is the vector of values for the distance r at which k_f(r) is estimated.

This algorithm assumes that X can be treated as a realisation of a stationary (spatially homogeneous) random spatial point process in the plane, observed through a bounded window. The window (which is specified in X as Window(X)) may have arbitrary shape.

Biases due to edge effects are treated in the same manner as in Kest. The edge corrections implemented here are

isotropic/Ripley: Ripley's isotropic correction (see Ripley, 1988; Ohser, 1983). This is implemented only for rectangular and polygonal windows (not for binary masks).
translate: Translation correction (Ohser, 1983). Implemented for all window geometries, but slow for complex windows.

Note that the estimator assumes the process is stationary (spatially homogeneous).

The numerator and denominator of the mark correlation function (in the expression above) are estimated using density estimation techniques. The user can choose between

"density": which uses the standard kernel density estimation routine density, and works only for evenly-spaced r values;
"loess": which uses the function loess in the package modreg;
"sm": which uses the function sm.density in the package sm and is extremely slow;
"smrep": which uses the function sm.density in the package sm and is relatively fast, but may require manual control of the smoothing parameter hmult.

If normalise=FALSE then the algorithm will compute only the numerator

c_f(r) = E_{0u} f(M(0),M(u))

of the expression for the mark correlation function. In this case, negative values of f are permitted.

Value

A function value table (object of class "fv") or a list of function value tables, one for each column of marks.

An object of class "fv" (see fv.object) is essentially a data frame containing numeric columns

`r`	the values of the argument `r` at which the mark correlation function `k_f(r)` has been estimated
`theo`	the theoretical value of `k_f(r)` when the marks attached to different points are independent, namely 1

together with a column or columns named "iso" and/or "trans", according to the selected edge corrections. These columns contain estimates of the mark correlation function k_f(r) obtained by the edge corrections named.

Author(s)

\spatstatAuthors

References

Stoyan, D. and Stoyan, H. (1994) Fractals, random shapes and point fields: methods of geometrical statistics. John Wiley and Sons.

Examples

    # CONTINUOUS-VALUED MARKS:
    # (1) Spruces
    # marks represent tree diameter
    # mark correlation function
    ms <- markcorr(spruces)
    plot(ms)

    # (2) simulated data with independent marks
    
      X <- rpoispp(100)
      X <- X %mark% runif(npoints(X))
      Xc <- markcorr(X)
      plot(Xc)
    
    
    # MULTITYPE DATA:
    # Hughes' amacrine data
    # Cells marked as 'on'/'off'
    X <- if(interactive()) amacrine else amacrine[c(FALSE, TRUE)]
    # (3) Kernel density estimate with Epanecnikov kernel
    # (as proposed by Stoyan & Stoyan)
    M <- markcorr(X, function(m1,m2) {m1==m2},
                  correction="translate", method="density",
                  kernel="epanechnikov")
    # Note: kernel="epanechnikov" comes from help(density)

    # (4) Same again with explicit control over bandwidth
    
      M <- markcorr(X, 
                  correction="translate", method="density",
                  kernel="epanechnikov", bw=0.02)
      # see help(density) for correct interpretation of 'bw'
   
   

   # weighted mark correlation
   X <- if(interactive()) betacells else betacells[c(TRUE,FALSE)]
   Y <- subset(X, select=type)
   a <- marks(X)$area
   v <- markcorr(Y, weights=a)

spatstat.explore documentation built on Oct. 22, 2024, 9:07 a.m.