ppdist | R Documentation |
Based on an arbitrary matrix of "distances" between the points of two point patterns xi and eta, this function computes versions of the transport-transform distance between xi and eta.
ppdist( dmat, penalty = 1, type = c("tt", "rtt", "TT", "RTT"), ret_matching = FALSE, p = 1, precision = NULL )
dmat |
a matrix specifying in its (i,j)-th entry the distance from the i-th point of xi to the j-th point of eta. |
penalty |
a positive number. The penalty for adding/deleting points. |
type |
either |
ret_matching |
logical. Shall the optimal point matching be returned? |
p |
a number >0. The matching is chosen such that the |
precision |
a small positive integer value. The precisions of the computations, which
are currently performed in integers. After correcting for the penalty, |
The transport-transform (TT) distance gives the minimal total cost for “morphing”
the pattern xi into the pattern eta by way of shifting points (at costs
specified in dmat
) and adding or deleting points (each at cost penalty
).
The total cost is determined as
(sum_{j=1}^n c_j^p)^{1/p},
where c_j denotes the cost for the jth individual operation and n is the cardinality of the larger point pattern.
The relative transport-transform (RTT) metric is exactly the same, but the sum in the total cost is divided by the larger cardinality:
(1/n * sum_{j=1}^n c_j^p)^{1/p}.
The TT- and RTT-metrics form an umbrella concept that includes the OSPA and Spike Time metrics frequently used in the literature. See Müller, Schuhmacher and Mateu (2020) for details.
The corresponding distance between the point patterns if ret_matching
is FALSE
.
Otherwise a list with components dist
containing
this distance and two vectors target1, target2
of integers, where
target
i specifies the indices of the points in the other pattern
that the points of the i-th pattern are matched to and NA
every
time a point is deleted.
There may be a minus in front of an index, where
-j
indicates that the corresponding pairing with point j
would be over a distance of more than 2^{1/p} * penalty. This is
equivalent to saying that the corresponding point of the first pattern
is deleted and the j-th point of the second pattern is added.
Note that having more than one minus implies that the matching is non-unique.
Dominic Schuhmacher schuhmacher@math.uni-goettingen.de
Raoul Müller, Dominic Schuhmacher and Jorge Mateu (2020).
Metrics and Barycenters for Point Pattern Data.
Statistics and Computing 30, 953-972.
doi: 10.1007/s11222-020-09932-y
# small example # ------------- set.seed(181230) xi <- spatstat.random::rpoispp(20) eta <- spatstat.random::rpoispp(20) dmat <- spatstat.geom::crossdist(xi,eta) res <- ppdist(dmat, penalty=1, type="rtt", ret_matching=TRUE, p=1) plotmatch(xi, eta, dmat, res, penalty=1, p=1) res$dist # for comparison: ospa-distance computation from spatstat: res_ospa <- spatstat.geom::pppdist(xi,eta,"spa") res_ospa$distance # exactly the same as above because nothing gets cut off # same example, but with a much smaller penalty for adding/deleting points # --------------------------------------------------------------- res <- ppdist(dmat, penalty=0.1, type="rtt", ret_matching=TRUE, p=1) plotmatch(xi, eta, dmat, res, penalty=0.1, p=1) # dashed lines indicate points that are deleted and re-added at new position # grey segments on dashed lines indicate the cost of deletion plus re-addition res$dist # for comparison: ospa-distance computation from spatstat # (if things do get cut off, we have to ensure that the cutoff distances # are the same, thus cutoff = 2^(1/p) * penalty): res_ospa <- spatstat.geom::pppdist(xi,eta,"spa",cutoff=0.2) res_ospa$distance # NOT the same as above res_ospa$distance - abs(xi$n-eta$n) * 0.1 / max(xi$n,eta$n) # the same as above # a larger example # --------------------------------------------------------------- set.seed(190203) xi <- spatstat.random::rpoispp(2000) eta <- spatstat.random::rpoispp(2000) dmat <- spatstat.geom::crossdist(xi,eta) res <- ppdist(dmat, penalty = 0.1, type = "rtt", ret_matching = TRUE, p = 1) res$dist # takes about 2-3 seconds
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.