Title: | Presence-Only for Marked Point Process |
---|---|
Description: | Inspired by Moreira and Gamerman (2022) <doi:10.1214/21-AOAS1569>, this methodology expands the idea by including Marks in the point process. Using efficient 'C++' code, the estimation is possible and made faster with 'OpenMP' <https://www.openmp.org/> enabled computers. This package was developed under the project PTDC/MAT-STA/28243/2017, supported by Portuguese funds through the Portuguese Foundation for Science and Technology (FCT). |
Authors: | Guido Alberti Moreira [cre, aut] |
Maintainer: | Guido Alberti Moreira <[email protected]> |
License: | GPL (>= 3) |
Version: | 0.1.3 |
Built: | 2024-11-11 07:30:27 UTC |
Source: | CRAN |
Generic class for the beta and delta parameters.
## S4 method for signature 'BetaDeltaPrior' show(object) ## S4 method for signature 'BetaDeltaPrior' print(x, ...) ## S3 method for class 'BetaDeltaPrior' print(x, ...)
## S4 method for signature 'BetaDeltaPrior' show(object) ## S4 method for signature 'BetaDeltaPrior' print(x, ...) ## S3 method for class 'BetaDeltaPrior' print(x, ...)
object |
The BetaDeltaPrior object. |
x |
The BetaDeltaPrior object. |
... |
Ignored. |
show
and print
: The invisible object.
family
The family of distributions of the prior.
Objects of this class is the output of the "covariates_importance" object
from the pompp_fit-class
. It can be plotted which uses
the graphics
package. The print
method
gives a point-wise estimation, the same seen in the bacplot
method.
Both plot
and boxplot
methods use the posterior distribution
of the importance.
## S3 method for class 'covariates_importance' print(x, component = "intensity", ...) ## S3 method for class 'covariates_importance' plot( x, component = "intensity", y = "importance", quantiles = c(0.025, 0.5, 0.975), ... ) ## S3 method for class 'covariates_importance' barplot(height, component = "intensity", y, ...) ## S3 method for class 'covariates_importance' boxplot(x, component = "intensity", ...)
## S3 method for class 'covariates_importance' print(x, component = "intensity", ...) ## S3 method for class 'covariates_importance' plot( x, component = "intensity", y = "importance", quantiles = c(0.025, 0.5, 0.975), ... ) ## S3 method for class 'covariates_importance' barplot(height, component = "intensity", y, ...) ## S3 method for class 'covariates_importance' boxplot(x, component = "intensity", ...)
x |
The |
component |
Either |
... |
Other parameters passed to |
y |
Either |
quantiles |
A 2- or 3-simensional vector with the desired quantiles
specified. If 3-dimensiona, the middle point is drawn as a dot when the
|
height |
The |
Objects of this class have two matrices where the Monte Carlo samples on the rows and parameters on the columns. One matrix is for the intensity importance and the other for the observability importance.
The invisible object.
Nothing is returned. Plot is called and drawn on the configured device.
A barplot. See barplot
for details. If component is selected
as "both"
, only the second barplot is returned.
A boxplot. See boxplot
for details. If component is selected
as "both"
, only the second boxplot is returned.
The model uses a data augmentation scheme to avoid performing approximations on the likelihood function.
fit_pompp( object, background, mcmc_setup = list(iter = 5000), verbose = TRUE, ... ) ## S4 method for signature 'pompp_model,matrix' fit_pompp( object, background, neighborhoodSize = 20, mcmc_setup, verbose = TRUE, area = 1, cores = parallel::detectCores(), ... ) checkFormatBackground(object, background)
fit_pompp( object, background, mcmc_setup = list(iter = 5000), verbose = TRUE, ... ) ## S4 method for signature 'pompp_model,matrix' fit_pompp( object, background, neighborhoodSize = 20, mcmc_setup, verbose = TRUE, area = 1, cores = parallel::detectCores(), ... ) checkFormatBackground(object, background)
object |
Either a |
background |
A matrix where the rows are the grid cells for the studied
region and the columns are the covariates. |
mcmc_setup |
A list containing |
verbose |
Set to |
... |
Parameters passed on to specific methods.
|
neighborhoodSize |
Number of neighbors to use in the NNGP method. |
area |
A positive number with the studied region's area. |
cores |
Number of cores to pass to OpenMP. |
The background is kept outside of the
An object of class "pompp_fit"
.
pompp_model
and pompp_fit-class
.
beta <- c(-1, 2) # Intercept = -1. Only one covariate delta <- c(3, 4) # Intercept = 3. Only one covariate lambdaStar <- 1000 gamma <- 2 mu <- 5 total_points <- rpois(1, lambdaStar) random_points <- cbind(runif(total_points), runif(total_points)) grid_size <- 50 # Find covariate values to explain the species occurrence. # We give them a Gaussian spatial structure. reg_grid <- as.matrix(expand.grid(seq(0, 1, len = grid_size), seq(0, 1, len = grid_size))) Z <- MASS::mvrnorm(1, rep(0, total_points + grid_size * grid_size), 3 * exp(-as.matrix(dist(rbind(random_points, reg_grid))) / 0.2)) Z1 <- Z[1:total_points]; Z2 <- Z[-(1:total_points)] # 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] * Z1)) n_occurrences <- sum(occurrences) occurrences_points <- random_points[occurrences,] occurrences_Z <- Z1[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 + grid_size * grid_size), 2 * exp(-as.matrix(dist(rbind(occurrences_points, reg_grid))) / 0.3)) W1 <- W[1:n_occurrences]; W2 <- W[-(1:n_occurrences)] 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] * W1 - gamma * S)) n_po <- sum(po_sightings) po_points <- occurrences_points[po_sightings, ] po_Z <- occurrences_Z[po_sightings] po_W <- W1[po_sightings] po_marks <- marks[po_sightings] jointPrior <- prior( NormalPrior(rep(0, 2), 10 * diag(2)), # Beta NormalPrior(rep(0, 3), 10 * diag(3)), # Delta GammaPrior(0.00001, 0.00001), # LambdaStar NormalPrior(0, 100), GammaPrior(0.001, 0.001) # Marks ) 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 = jointPrior) bkg <- cbind(Z2, W2, reg_grid) # Create background # Be prepared to wait a long time for this fit <- fit_pompp(model, bkg, neighborhoodSize = 20, area = 1, mcmc_setup = list(burnin = 1000, iter = 2000), cores = 1) summary(fit) # Rhat upper CI values are above 1.1. More iterations are needed...
beta <- c(-1, 2) # Intercept = -1. Only one covariate delta <- c(3, 4) # Intercept = 3. Only one covariate lambdaStar <- 1000 gamma <- 2 mu <- 5 total_points <- rpois(1, lambdaStar) random_points <- cbind(runif(total_points), runif(total_points)) grid_size <- 50 # Find covariate values to explain the species occurrence. # We give them a Gaussian spatial structure. reg_grid <- as.matrix(expand.grid(seq(0, 1, len = grid_size), seq(0, 1, len = grid_size))) Z <- MASS::mvrnorm(1, rep(0, total_points + grid_size * grid_size), 3 * exp(-as.matrix(dist(rbind(random_points, reg_grid))) / 0.2)) Z1 <- Z[1:total_points]; Z2 <- Z[-(1:total_points)] # 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] * Z1)) n_occurrences <- sum(occurrences) occurrences_points <- random_points[occurrences,] occurrences_Z <- Z1[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 + grid_size * grid_size), 2 * exp(-as.matrix(dist(rbind(occurrences_points, reg_grid))) / 0.3)) W1 <- W[1:n_occurrences]; W2 <- W[-(1:n_occurrences)] 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] * W1 - gamma * S)) n_po <- sum(po_sightings) po_points <- occurrences_points[po_sightings, ] po_Z <- occurrences_Z[po_sightings] po_W <- W1[po_sightings] po_marks <- marks[po_sightings] jointPrior <- prior( NormalPrior(rep(0, 2), 10 * diag(2)), # Beta NormalPrior(rep(0, 3), 10 * diag(3)), # Delta GammaPrior(0.00001, 0.00001), # LambdaStar NormalPrior(0, 100), GammaPrior(0.001, 0.001) # Marks ) 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 = jointPrior) bkg <- cbind(Z2, W2, reg_grid) # Create background # Be prepared to wait a long time for this fit <- fit_pompp(model, bkg, neighborhoodSize = 20, area = 1, mcmc_setup = list(burnin = 1000, iter = 2000), cores = 1) summary(fit) # Rhat upper CI values are above 1.1. More iterations are needed...
Constructor for GammaPrior-class
objects
GammaPrior(shape, rate)
GammaPrior(shape, rate)
shape |
A positive number. |
rate |
A positive number. |
A GammaPrior
object with adequate slots.
This is used to represent the prior for lambdaStar individually. It still needs to be joined with the prior for Beta and Delta to be used in a model.
## S4 method for signature 'GammaPrior' names(x) ## S4 method for signature 'GammaPrior' x$name ## S4 replacement method for signature 'GammaPrior' x$name <- value ## S4 method for signature 'GammaPrior' show(object) ## S4 method for signature 'GammaPrior' print(x, ...) ## S3 method for class 'GammaPrior' print(x, ...)
## S4 method for signature 'GammaPrior' names(x) ## S4 method for signature 'GammaPrior' x$name ## S4 replacement method for signature 'GammaPrior' x$name <- value ## S4 method for signature 'GammaPrior' show(object) ## S4 method for signature 'GammaPrior' print(x, ...) ## S3 method for class 'GammaPrior' print(x, ...)
x |
The GammaPrior object. |
name |
The requested slot. |
value |
New value. |
object |
The GammaPrior object. |
... |
Ignored. |
names
: A character vector with the prior parameters.
`$`
The requested slot's value.
`$<-`
: The new object with the updated slot.
show
and print
: The invisible object.
shape
The shape parameter of the Gamma distribution.
rate
The rate parameter of the Gamma distribution.
GammaPrior(0.0001, 0.0001)
GammaPrior(0.0001, 0.0001)
Helper function to create a valid set of initial values to be used with the fit_pompp function.
initial( beta = numeric(), delta = numeric(), lambdaStar = numeric(), marksMean = numeric(), marksPrecision = numeric(), random = FALSE )
initial( beta = numeric(), delta = numeric(), lambdaStar = numeric(), marksMean = numeric(), marksPrecision = numeric(), random = FALSE )
beta |
Either a vector or a single integer. The vector is used if the initial values are provided and the integer is used as the vector size to be randomly generated. |
delta |
Either a vector or a single integer. The vector is used if the initial values are provided and the integer is used as the vector size to be randomly generated. |
lambdaStar |
A positive number. |
marksMean |
Any real number. If random, defines the mean of the random value. |
marksPrecision |
A positive number. If random, defines the mean of the random value. |
random |
A logical value. If |
A pompp_initial
object. It can be used in the
fit_pompp
function by itself, but must be in a list if multiple
initial values are supplied. Initial values can be combined by adding them
(with the use of +
).
# Let us create initial values for a model with, say, 3 intensity covariates # and 4 observability covariates. We add an initial values for both these # cases due to the intercepts. # This first one is in1 <- initial(rep(0, 4), c(0, 2, -1, -2, 3), 100, 0, 1) # Then we initalize some randomly. in2 <- initial(4, 5, 100, 0, 1, random = TRUE) # We can even multiply the random one to generate more. Let us join them all # to include in a model. initial_values <- in1 + in2 * 3 # 4 chains are initialized.
# Let us create initial values for a model with, say, 3 intensity covariates # and 4 observability covariates. We add an initial values for both these # cases due to the intercepts. # This first one is in1 <- initial(rep(0, 4), c(0, 2, -1, -2, 3), 100, 0, 1) # Then we initalize some randomly. in2 <- initial(4, 5, 100, 0, 1, random = TRUE) # We can even multiply the random one to generate more. Let us join them all # to include in a model. initial_values <- in1 + in2 * 3 # 4 chains are initialized.
Generic class for the LambdaStar parameters.
## S4 method for signature 'LambdaStarPrior' show(object)
## S4 method for signature 'LambdaStarPrior' show(object)
object |
The LambdaStarPrior object. |
show
and print
: The invisible object.
family
The family of distributions of the prior.
Constructor for NormalPrior-class
objects
NormalPrior(mu, Sigma)
NormalPrior(mu, Sigma)
mu |
The mean vector for the Normal distribution. |
Sigma |
The covariance matrix for the Normal distribution. |
Matrix Sigma must be square and positive definite. Its dimensions must match mu's length.
A NormalPrior
object with adequate slots.
NormalPrior(rep(0, 10), diag(10) * 10)
NormalPrior(rep(0, 10), diag(10) * 10)
This is used to represent the prior for Beta and Delta individually. They still need to be joined to be used in a model.
## S4 method for signature 'NormalPrior' names(x) ## S4 method for signature 'NormalPrior' x$name ## S4 replacement method for signature 'NormalPrior' x$name <- value ## S4 method for signature 'NormalPrior' show(object) ## S4 method for signature 'NormalPrior' print(x, ...) ## S3 method for class 'NormalPrior' print(x, ...)
## S4 method for signature 'NormalPrior' names(x) ## S4 method for signature 'NormalPrior' x$name ## S4 replacement method for signature 'NormalPrior' x$name <- value ## S4 method for signature 'NormalPrior' show(object) ## S4 method for signature 'NormalPrior' print(x, ...) ## S3 method for class 'NormalPrior' print(x, ...)
x |
The NormalPrior object. |
name |
The requested slot. |
value |
New value. |
object |
The NormalPrior object. |
... |
Ignored. |
names
: A character vector with the prior parameters.
`$`
: The requested slot's value.
`$<-`
: The new object with the updated slot.
show
and print
: The invisible object.
mu
The mean vector for the prior.
Sigma
The covariance matrix for the prior.
Objects of this class are the main objects of this package. They contain much information about the fitted model.
## S4 method for signature 'pompp_fit' show(object) ## S4 method for signature 'pompp_fit' print(x, ...) ## S3 method for class 'pompp_fit' print(x, ...) ## S4 method for signature 'pompp_fit' summary(object, ...) ## S3 method for class 'pompp_fit' summary(object, ...) ## S4 method for signature 'pompp_fit' names(x) ## S3 method for class 'pompp_fit' names(x) ## S4 method for signature 'pompp_fit,ANY,ANY' x[[i]] ## S4 method for signature 'pompp_fit' x$name ## S4 method for signature 'pompp_fit' as.array(x, ...) ## S3 method for class 'pompp_fit' as.array(x, ...) ## S4 method for signature 'pompp_fit' as.matrix(x, ...) ## S3 method for class 'pompp_fit' as.matrix(x, ...) ## S4 method for signature 'pompp_fit' as.data.frame(x, row.names = NULL, optional = FALSE, ...) ## S3 method for class 'pompp_fit' as.data.frame(x, row.names = NULL, optional = FALSE, ...) ## S4 method for signature 'pompp_fit,pompp_fit' e1 + e2 ## S4 method for signature 'pompp_fit' c(x, ...)
## S4 method for signature 'pompp_fit' show(object) ## S4 method for signature 'pompp_fit' print(x, ...) ## S3 method for class 'pompp_fit' print(x, ...) ## S4 method for signature 'pompp_fit' summary(object, ...) ## S3 method for class 'pompp_fit' summary(object, ...) ## S4 method for signature 'pompp_fit' names(x) ## S3 method for class 'pompp_fit' names(x) ## S4 method for signature 'pompp_fit,ANY,ANY' x[[i]] ## S4 method for signature 'pompp_fit' x$name ## S4 method for signature 'pompp_fit' as.array(x, ...) ## S3 method for class 'pompp_fit' as.array(x, ...) ## S4 method for signature 'pompp_fit' as.matrix(x, ...) ## S3 method for class 'pompp_fit' as.matrix(x, ...) ## S4 method for signature 'pompp_fit' as.data.frame(x, row.names = NULL, optional = FALSE, ...) ## S3 method for class 'pompp_fit' as.data.frame(x, row.names = NULL, optional = FALSE, ...) ## S4 method for signature 'pompp_fit,pompp_fit' e1 + e2 ## S4 method for signature 'pompp_fit' c(x, ...)
object |
A pompp_fit object. |
x |
A pompp_fit object. |
... |
Ignored in this version. |
i |
The requested slot. |
name |
The requested slot. |
row.names |
NULL or a character vector giving the row names for the data frame. Missing values are not allowed. |
optional |
logical. If TRUE, setting row names and converting column
names to syntactic names is optional. See help('as.data.frame') for more.
Leaving as |
e1 |
A pompp_fit object. |
e2 |
A pompp_fit object with the same slots (except for initial
values) as |
show
and print
: The invisible object.
summary
: A matrix with the summary statistics of the
fit. It is also printed in the print
method. The summary can be
treated as a matrix, such as retrieving rows/columns and creating tables
with the xtable
package.
names
: A character vector with the available options
for the `$`
and `[[`
methods.
`$`
and `[[`
: The requested slot.
Available options are not necessarily the class slots, and can be checked
with the names
method.
as.array
: An array
with dimensions I x C x P,
where I stands for number of iterations, C for number of chains and P for
total number of parameters. P is actually larger than the number of
parameters in the model, as the the generated sizes of the latent processes
and the log-posterior are also included. This is organized so that is ready
for the bayesplot
package functions.
as.matrix
: The dimension of the output is
I * C x (P + 2), where I stands for number of iterations, C for number of
chains and P for total number of parameters. P is actually larger than the
number of parameters in the model, as the generated sizes of the latent
processes and the log-posterior are also included.
Two extra columns are included to indicate to which chain and to which iteration that draw belongs.
as.data.frame
: The dimension of the output is
I*C x P + 2, where I stands for number of iterations, C for number of chains
and P for total number of parameters. P is actually larger than the number
of parameters in the model, as the generated sizes of the latent processes
and the log-posterior are also included.
Two extra columns are included to indicate to which chain and to which
iteration that draw belongs. This is to facilitate the use of plotting
results via the ggplot2
package if desired.
If row.names
is left at NULL
then row names are created as
CcIi where c is the chain and i is the iteration of that row.
+
: A new pompp_fit
object where the chains
are combined into a new multi-chain object. This can be used if chains are
run in separate occasions or computers to combine them into a single object
for analysis.
c
: A new pompp_fit
object where the chains
are combined into a new multi-chain object. The +
method is
used for that, with the same arguments restrictions and results.
fit
The actual fit from the model. It is an object of class
mcmc.list
, as generated from the coda
package.
original
The model used to generate the chains, an object with class
pompp_model
.
backgroundSummary
A small summary of the original background covariates. This is to ensure that continuing the chains will use the identical background matrix. Only the summary is kept for storage efficiency.
area
A positive number indicating the area measure of the region being studied.
parnames
The names of the parameters. If the model used selects the covariates with column names, they are replicated here. If they are the column indexes, names are generated for identification.
mcmc_setup
The original mcmc setup used.
Class for the initial values for the MCMC for the pompp package
## S4 method for signature 'pompp_initial' names(x) ## S4 method for signature 'pompp_initial' x$name ## S4 method for signature 'pompp_initial,ANY' e1 + e2 ## S4 method for signature 'list,pompp_initial' e1 + e2 ## S4 method for signature 'pompp_initial,list' e1 + e2 ## S4 method for signature 'pompp_initial,numeric' e1 * e2 ## S4 method for signature 'numeric,pompp_initial' e1 * e2 ## S4 method for signature 'pompp_initial' show(object) ## S4 method for signature 'pompp_initial' print(x, ...) ## S3 method for class 'pompp_initial' print(x, ...)
## S4 method for signature 'pompp_initial' names(x) ## S4 method for signature 'pompp_initial' x$name ## S4 method for signature 'pompp_initial,ANY' e1 + e2 ## S4 method for signature 'list,pompp_initial' e1 + e2 ## S4 method for signature 'pompp_initial,list' e1 + e2 ## S4 method for signature 'pompp_initial,numeric' e1 * e2 ## S4 method for signature 'numeric,pompp_initial' e1 * e2 ## S4 method for signature 'pompp_initial' show(object) ## S4 method for signature 'pompp_initial' print(x, ...) ## S3 method for class 'pompp_initial' print(x, ...)
x |
The pompp_initial object. |
name |
The requested slot. |
e1 |
A pompp_initial object. |
e2 |
Another pompp_initial object or a list with pompp_initial objects for + and a positive integer for *. e1 and e2 can be switched (+ and * are commutative). |
object |
A pompp_initial object. |
... |
Currently unused. |
names
: A character vector with the initialized
parameter names.
`$`
: The requested initial value (in case of
LambdaStar) or values (in case of Beta or Delta).
+
: A list with the objects. Useful to start the
fit_pompp
function, as it requires a list of initial values.
*
: A list with e2
random initial values.
show
and print
: The invisible object.
beta
Initial values for beta.
delta
Initial values for delta.
lambdaStar
Initial values for lambdaStar.
tag
Indicates the source of the initial values.
pompp
fitting functionConstructor 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.
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 )
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 )
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
|
joint_prior |
A |
verbose |
Set to |
A pompp_model
object with the requested slots. It is ready
to be used in the fit_pompp
function.
# 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
# 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
The model includes the presence-only data, all selected variables, the link
functions for and
, the initial values and the prior
distribution.
## S4 method for signature 'pompp_model' names(x) ## S4 method for signature 'pompp_model' x$name ## S4 replacement method for signature 'pompp_model' x$name <- value ## S4 method for signature 'pompp_model' show(object) ## S4 method for signature 'pompp_model' print(x, ...) ## S3 method for class 'pompp_model' print(x, ...)
## S4 method for signature 'pompp_model' names(x) ## S4 method for signature 'pompp_model' x$name ## S4 replacement method for signature 'pompp_model' x$name <- value ## S4 method for signature 'pompp_model' show(object) ## S4 method for signature 'pompp_model' print(x, ...) ## S3 method for class 'pompp_model' print(x, ...)
x |
The pompp_model object. |
name |
The requested slot. |
value |
New value. |
object |
The pompp_model object. |
... |
Currently unused. |
names
: A character vector with possible options
for the `$`
and `$<-`
methods.
`$`
: The requested slot's value.
`$<-`
: The new object with the updated slot.
show
and print
: The invisible object.
po
The matrix containing the covariates values for the data.
intensityLink
A string informing about the chosen link for the
intensity covariates. Current acceptable choice is only "logit"
.
intensitySelection
A vector containing the indexes of the selected
intensity columns in the po
matrix.
observabilityLink
A string informing about the chosen link for the
observability covariates. Current acceptable choice is only "logit"
.
observabilitySelection
A vector containing the indexes of the selected
observability columns in the po
matrix.
marksSelection
A single value containing the index of the selected
marks column in the po
matrix.
coordinates
A vector of two values containing the column positions
of the longitude and latitude in the po
matrix.
init
A list with objects of class pompp_initial
indicating
the initial values for each chain. The length of this list tells the program
how many chains are requested to be run.
prior
An object of class pompp_prior
which indicates the
joint prior distribution for the model parameters.
iSelectedColumns
If the intensity covariates selection was made with the name of the columns, they are stored in this slot.
oSelectedColumns
If the observability covariates selection was made with the name of the columns, they are stored in this slot.
mSelectedColumns
If the marks selection was made with the name of the column, it is stored in this slot.
pompp_initial-class
and
pompp_prior-class
and pompp_model
Objects of this class are the joining of independent priors for Beta, Delta
and LambdaStar. They can be used in the fit_pompp
function.
## S4 method for signature 'pompp_prior' names(x) ## S4 method for signature 'pompp_prior' x$name ## S4 method for signature 'pompp_prior' show(object) ## S4 method for signature 'pompp_prior' print(x, ...) ## S3 method for class 'pompp_prior' print(x, ...) ## S4 method for signature 'pompp_prior' x$name ## S4 replacement method for signature 'pompp_prior' x$name <- value
## S4 method for signature 'pompp_prior' names(x) ## S4 method for signature 'pompp_prior' x$name ## S4 method for signature 'pompp_prior' show(object) ## S4 method for signature 'pompp_prior' print(x, ...) ## S3 method for class 'pompp_prior' print(x, ...) ## S4 method for signature 'pompp_prior' x$name ## S4 replacement method for signature 'pompp_prior' x$name <- value
x |
The pompp_prior object. |
name |
The requested slot. |
object |
The pompp_prior object. |
... |
Ignored. |
value |
New value. |
names
: A character vector with the model parameters
names.
`$`
: The requested slot's value.
`$<-`
: The new object with the updated slot.
beta
An object of a class which inherits the BetaDeltaPrior
S4
class with the appropriate Beta prior.
delta
An object of a class which inherits the BetaDeltaPrior
S4
class with the appropriate Delta prior.
lambdaStar
An object of a class which inherits the
LambdaStarPrior
S4 class with the appropriate LambdaStar prior.
marksMean
An object of S4 class NormalPrior
with the chosen
prior for the marks mean
marksPrecision
An object of S4 class GammaPrior
with the chosen
prior for the marks precision
Constructor for pompp_prior
objects, which is used in the
pompp_fit
function. The generated prior is so that Beta, Delta
and LambdaStar are indepdendent a priori.
prior(beta, delta, lambdaStar, marksMean, marksPrecision)
prior(beta, delta, lambdaStar, marksMean, marksPrecision)
beta |
An S4 object whose class inherits from |
delta |
An S4 object whose class inherits from |
lambdaStar |
An S4 object whose class inherits from |
marksMean |
An S4 object of class |
marksPrecision |
An S4 object of class |
A pompp_prior
object with the adequate slots. It is ready to
be included in a model via the pompp_model
function.
fit_pompp
, NormalPrior
,
GammaPrior
and pompp_model
.
# Let us say there are 3 intensity covariates and 4 observability covariates. # One more element is included in both sets due to the intercepts. new_prior <- prior( NormalPrior(rep(0, 4), 10 * diag(4)), NormalPrior(rep(0, 5), 10 * diag(5)), GammaPrior(0.0001, 0.0001), NormalPrior(0, 100), GammaPrior(0.001, 0.001) )
# Let us say there are 3 intensity covariates and 4 observability covariates. # One more element is included in both sets due to the intercepts. new_prior <- prior( NormalPrior(rep(0, 4), 10 * diag(4)), NormalPrior(rep(0, 5), 10 * diag(5)), GammaPrior(0.0001, 0.0001), NormalPrior(0, 100), GammaPrior(0.001, 0.001) )