pompp_model: Build a model to be used in the 'pompp' fitting function

In pompp: Presence-Only for Marked Point Process

Build a model to be used in the pompp fitting function

Description

Constructor for pompp_model-class objects, built to facilitate the use of the fitting function. The output of this function has the necessary signature for the fit_pompp function to start the model fit.

Usage

pompp_model(
  po,
  intensitySelection,
  observabilitySelection,
  marksSelection,
  coordinates,
  intensityLink = "logit",
  observabilityLink = "logit",
  initial_values = 1,
  joint_prior = prior(beta = NormalPrior(rep(0, length(intensitySelection) + 1), 10 *
    diag(length(intensitySelection) + 1)), delta = NormalPrior(rep(0,
    length(observabilitySelection) + 1), 10 * diag(length(observabilitySelection) + 1)),
    lambdaStar = GammaPrior(1e-10, 1e-10), marksMean = NormalPrior(1, 100),
    marksPrecision = GammaPrior(0.001, 0.001)),
  verbose = TRUE
)

Arguments

po	A matrix whose rows represent the presence-only data and the columns the covariates observed at each position.
intensitySelection	Either a numeric or character vector and represents the selection of covariates used for the intensity set. If numeric it is the positions of the columns and if character, the names of the columns.
observabilitySelection	Either a numeric or character vector and represents the selection of covariates used for the observability set. If numeric it is the positions of the columns and if character, the names of the columns.
marksSelection	Either a numeric or character vector and represents the selection of column used for the marks. If numeric it is the position of the column and if character, the name of the column.
coordinates	Either a numeric or character vector and represents the columns to be used for the coordinates. If numeric it is the positions of the columns and if character, the names of the columns. They must be in longitude, latitude order.
intensityLink	A string to inform what link function the model has with respect to the intensity covariates. Current version accepts 'logit'.
observabilityLink	A string to inform what link function the model has with respect to the observabilitycovariates. Current version accepts 'logit'.
initial_values	Either a single integer, a single pompp_initial-class or a list containing pompp_initial-class objects. The length of the list will inform the model how many independent chains will be run. If an integer, that many initial values will be randomly generated.
joint_prior	A pompp_prior object.
verbose	Set to FALSE to suppress all messages to console.

Value

A pompp_model object with the requested slots. It is ready to be used in the fit_pompp function.

See Also

initial, prior and fit_pompp.

Examples

# Let us simulate some data to showcase the creation of the model.
beta <- c(-1, 2)
delta <- c(3, 4)
lambdaStar <- 1000
gamma <- 2
mu <- 5

total_points <- rpois(1, lambdaStar)
random_points <- cbind(runif(total_points), runif(total_points))

# Find covariate values to explain the species occurrence.
# We give them a Gaussian spatial structure.
Z <- MASS::mvrnorm(1, rep(0, total_points), 3 * exp(-as.matrix(dist(random_points)) / 0.2))

# Thin the points by comparing the retaining probabilities with uniforms
# in the log scale to find the occurrences
occurrences <- log(runif(total_points)) <= -log1p(exp(-beta[1] - beta[2] * Z))
n_occurrences <- sum(occurrences)
occurrences_points <- random_points[occurrences,]
occurrences_Z <- Z[occurrences]

# Find covariate values to explain the observation bias.
# Additionally create a regular grid to plot the covariate later.
W <- MASS::mvrnorm(1, rep(0, n_occurrences), 2 * exp(-as.matrix(dist(occurrences_points)) / 0.3))
S <- MASS::mvrnorm(1, rep(0, n_occurrences), 2 * exp(-as.matrix(dist(occurrences_points)) / 0.3))

# Find the presence-only observations.
marks <- exp(mu + S + rnorm(n_occurrences, 0, 0.3))
po_sightings <- log(runif(n_occurrences)) <= -log1p(exp(-delta[1] - delta[2] * W - gamma * S))
n_po <- sum(po_sightings)
po_points <- occurrences_points[po_sightings, ]
po_Z <- occurrences_Z[po_sightings]
po_W <- W[po_sightings]
po_marks <- marks[po_sightings]

# Now we create the model
model <- pompp_model(po = cbind(po_Z, po_W, po_points, po_marks),
  intensitySelection = 1, observabilitySelection = 2,
  marksSelection = 5, coordinates = 3:4,
  intensityLink = "logit", observabilityLink = "logit",
  initial_values = 2, joint_prior = prior(
    NormalPrior(rep(0, 2), 10 * diag(2)),
    NormalPrior(rep(0, 3), 10 * diag(3)),
    GammaPrior(1e-4, 1e-4),
    NormalPrior(0, 100), GammaPrior(0.001, 0.001)))
# Check how it is.
model

pompp documentation built on Dec. 28, 2022, 2:24 a.m.