Description Usage Arguments Value Author(s) References See Also Examples
Markov chain Monte Carlo sampler from Matern-III marked spatial point process
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, ...)
|
mat3 |
A mat3 object, see |
fname |
Will only be read if mat3 is set to NA. File names with the centers
and fingers dataset, to be read with |
win |
Square window in which the process is to be sampled. See |
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 |
initialValues |
Some holistic initial values for empirical priors on the parameters.
See |
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 |
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. |
Guilherme Ludwig and Nancy Garcia
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.
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)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.