fitmat3: Fit a Matern-III marked point process model via MCMC sampling

Description Usage Arguments Value Author(s) References See Also Examples

View source: R/fitmat3.R

Description

Markov chain Monte Carlo sampler from Matern-III marked spatial point process

Usage

1
2
3
4
5
6
7
fitmat3(mat3, fname = c("centers.txt", "fingers.txt"), win = owin(c(0,
  1), c(0, 1)), R_clusters = 0.005, R_centers = 0.02, L = 10000,
  N = 10000, J = 100, seed = NULL, resultsName = NULL,
  logPriors = preparePriors(), initialValues = pickInitialValues(),
  GeyerThompson = c("large", "resample", "importance", "simple",
  "MollerAlt"), tildeTheta = c(length(mat3)/area(win), 1),
  candidateVar = c(0.1, 0.2, 0.1, 0.1, 0.1), verbose = TRUE, ...)

Arguments

mat3

A mat3 object, see rmat3. If reading from a data file, set to NA instead.

fname

Will only be read if mat3 is set to NA. File names with the centers and fingers dataset, to be read with read.table. The default action assumes the files "centers.txt" and "fingers.txt", each formatted with a header and three columns corresponding to a center ID (integer), and the spatial location of the center/finger end. The fingers file should have at least one finger for each center ID in the centers file.

win

Square window in which the process is to be sampled. See owin. Defaults to the unit square.

R_clusters

Nuisance parameter for the fingers' inhibition. No finger endings will be closer than R_clusters in Euclidean distance during the process sampling step (hard core inhibiition).

R_centers

Nuisance parameter for the centers' inhibition. At the birth-and-death process, process will take into account how many fingers (attached to other active centers) are currently at R_centers's distance from the candidate center being born.

L

Number of samples for the MCMC sampler. Defaults to 10000.

N

Number of samples for integrating the area of the shadow. Defaults to 10000.

J

Number of updates for the Geyer-Thompson step in estimating the sampling distribution for beta and phi (if GeyerThompson is set to "large", 100 times J). Defaults to 100 for "simple" and "importance" options in GeyerThompson, and 10000 for "large".

seed

Fixes the RNG seed for replication. Defaults to NULL, which does not fix the seed.

resultsName

Writes the results to a file; defaults to NULL, which does not save the results.

logPriors

A list of 5 log-priors for the coefficients, each being a function of the parameter and initial value only. See preparePriors, which produces the default values of logPriors.

initialValues

Some holistic initial values for empirical priors on the parameters. See pickInitialValues, which produces the default values of initialValues.

GeyerThompson

Provides different methods for updating the Geyer and Thompson (1999) step in the computation of the likelihood. The "large" option produces J*100 samples from which the Geyer-Thompson approximation ot the likelihood is built; "importance" uses J samples, but they are updated at every 100 iterations. "simple" works like "importance", but sets phi = 1. "Resample" uses a resample method on existing finger. "MollerAlt" uses an alternative approach of auxiliary variabiles to avoid computing the normalizing constant, described in Moller et al (2006).

tildeTheta

Only needed for the "MollerAlt" option in "GeyerThompson". Corresponds to a vector of length two, with estimates to parameters beta and phi. Defaults to average intensity and 1.

candidateVar

Variability of MCMC candidate selection step, if tuning is necessary.

verbose

Prints acceptance rates for candidates. Defaults to TRUE.

...

further arguments passed to read.table.

Value

fitmat3 returns a list containing at least the following components

parameters

An L by 5 matrix with the samples from beta, phi, gamma, sigma, kappa.

Author(s)

Guilherme Ludwig and Nancy Garcia

References

Garcia, N., Guttorp, P. and Ludwig, G. (2020) "Interacting cluster point process model for epidermal nerve fibers", to appear.

Geyer, C. J., and Thompson, E. A. (1992) "Constrained Monte Carlo maximum likelihood for dependent data." Journal of the Royal Statistical Society, Series B, Vol. 54, pp. 657-699.

Moller, J., Pettitt, A. N., Reeves, R. and Berthelsen, K. K. (2006) "An efficient Markov chain Monte Carlo method for distributions with intractable normalising constants." Biometrika, Vol. 93, No. 2, pp. 451-458.

See Also

rmat3

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
set.seed(1234)
x <- rmat3(70, 2, 5, 0.05, 3)
plot(x)
# Changing default sampling sizes to make it run fast
model <- fitmat3(x, L = 100, N = 1000, J = 2, seed = 1234)
burnin <- 20
layout(matrix(c(1:5,0), ncol=2))
plot(model$parameters[-c(1:burnin),1], type="l", ylab="beta")
abline(h = median(model$parameters[-c(1:burnin),1]))
abline(h = 2.5*25, col = "Red")
plot(model$parameters[-c(1:burnin),2], type="l", ylab="phi")
abline(h = median(model$parameters[-c(1:burnin),2]))
abline(h = 2, col = "Red")
plot(model$parameters[-c(1:burnin),3], type="l", ylab="gamma")
abline(h = median(model$parameters[-c(1:burnin),3]))
abline(h = 20, col = "Red")
plot(model$parameters[-c(1:burnin),4], type="l", ylab="sigma")
abline(h = median(model$parameters[-c(1:burnin),4]))
abline(h = 0.05, col = "Red")
plot(model$parameters[-c(1:burnin),5], type="l", ylab="kappa")
abline(h = median(model$parameters[-c(1:burnin),5]))
abline(h = 3, col = "Red")
# Saving results to file:
## Not run: model <- fitmat3(NA, fname = c("centers.txt", "fingers.txt"),
                          seed = 1234, resultsName = "results.txt")
## End(Not run)

guiludwig/mat3c documentation built on Dec. 2, 2019, 1:32 a.m.