- endemic
right-hand side formula for the exponential (Cox-like
multiplicative) endemic component. May contain offsets (to be marked
by the special function offset
). If omitted or ~0
there will be no endemic component in the model. A type-specific
endemic intercept can be requested by including the term
(1|type)
in the formula.
- epidemic
formula representing the epidemic model for the event-specific
covariates (marks) determining infectivity. Offsets are not
implemented here. If omitted or ~0
there will be no epidemic
component in the model.
- siaf
spatial interaction function. Possible specifications are:
NULL
or missing, corresponding to
siaf.constant()
, i.e. spatially homogeneous
infectivity independent of the distance from the host
a list as returned by siaf
or, more commonly,
generated by a predefined interaction function such as
siaf.gaussian
as in Meyer et al. (2012) or
siaf.powerlaw
as in Meyer and Held (2014).
The latter requires unique event locations, possibly after random
tie-breaking (untie
) or imputation of
interval-censored locations.
siaf.exponential
is a simpler alternative.
a numeric vector corresponding to the knots of a step
function, i.e. the same as siaf.step(knots)
If you run into “false convergence” with a non-constant
siaf
specification, the numerical accuracy of the cubature
methods is most likely too low (see the control.siaf
argument).
- tiaf
temporal interaction function. Possible specifications are:
NULL
or missing, corresponding to
tiaf.constant()
, i.e. time-constant infectivity
a list as returned by tiaf
or by a
predefined interaction function such as
tiaf.exponential
a numeric vector corresponding to the knots of a step
function, i.e. the same as tiaf.step(knots)
- qmatrix
square indicator matrix (0/1 or FALSE
/TRUE
) for
possible transmission between the event types. The matrix will be
internally converted to logical
. Defaults to the \(Q\) matrix
specified in data
.
- data
an object of class "epidataCS"
.
- subset
an optional vector evaluating to logical indicating a subset of
data$events
to keep. Missing values are
taken as FALSE
. The expression is evaluated in the context of the
data$events@data
data.frame
, i.e. columns of this
data.frame
may be referenced directly by name.
- t0, T
events having occurred during (-Inf;t0] are regarded as part of the
prehistory \(H_0\) of the process. Only events that occurred in the
interval (t0; T] are considered in the likelihood.
The time point t0
(T
) must
be an element of data$stgrid$start
(data$stgrid$stop
).
The default time range covers the whole spatio-temporal grid
of endemic covariates.
- na.action
how to deal with missing values in data$events
? Do not use
na.pass
. Missing values in the spatio-temporal grid
data$stgrid
are not accepted.
- start
a named vector of initial values for (a subset of) the parameters.
The names must conform to the conventions of twinstim
to be
assigned to the correct model terms. For instance,
"h.(Intercept)"
= endemic intercept,
"h.I(start/365)"
= coefficient of a linear time
trend in the endemic component, "h.factorB"
=
coefficient of the level B of the factor variable factor
in
the endemic predictor, "e.(Intercept)"
= epidemic intercept,
"e.VAR"
= coefficient of the epidemic term VAR
,
"e.siaf.2"
= second siaf
parameter,
"e.tiaf.1"
= first tiaf
parameter.
Elements which don't match any of the model parameters are ignored.
Alternatively, start
may also be a named list with elements
"endemic"
or "h"
, "epidemic"
or "e"
,
"siaf"
or "e.siaf"
, and "tiaf"
or "e.tiaf"
,
each of which containing a named numeric vector with the term labels
as names (i.e. without the prefix "h."
, "e."
, etc).
Thus, start=list(endemic=c("(Intercept)"=-10))
is equivalent
to start=c("h.(Intercept)"=-10)
.
- partial
logical indicating if a partial likelihood similar to the approach
by Diggle et al. (2010) should be used (default is FALSE
).
Note that the partial likelihood implementation is not well tested.
- epilink
a character string determining the link function to be used for the
epidemic
linear predictor of event marks. By default, the
log-link is used. The experimental alternative
epilink = "identity"
(for use by epitest
) does
not guarantee the force of infection to be positive. If this leads
to a negative total intensity (endemic + epidemic), the point
process is not well defined (the log-likelihood will be
NaN
).
- control.siaf
a list with elements "F"
and "Deriv"
, which are lists
of extra arguments passed to the functions siaf$F
and
siaf$Deriv
, respectively.
These arguments control the accuracy of the cubature routines from
package polyCub involved in non-constant siaf
specifications, e.g., the bandwidth of the midpoint rule
polyCub.midpoint
, the number of Gaussian
quadrature points for polyCub.SV
, or the relative
tolerance of integrate
in polyCub.iso
.
For instance, siaf.gaussian(F.adaptive = TRUE)
uses the
midpoint-cubature polyCub.midpoint
with an
adaptive bandwidth of eps=adapt*sd
to numerically integrate the
kernel \(f(\bold{s})\), and the default adapt
value (0.1)
can be overwritten by setting control.siaf$F$adapt
.
However, the default version siaf.gaussian()
as well as siaf.powerlaw()
and friends use
polyCub.iso
and thus accept control arguments for the
standard integrate
routine (such as rel.tol
)
via control.siaf$F
and control.siaf$Deriv
.
This argument list is ignored in the case
siaf=siaf.constant()
(which is the default if siaf
is
unspecified).
- optim.args
an argument list passed to optim
, or NULL
, in
which case no optimization will be performed
but the necessary functions will be returned in a list (similar to
what is returned if model = TRUE
).
Initial values for the parameters may be given as list element
par
in the order (endemic, epidemic, siaf, tiaf)
.
If no initial values are provided, crude estimates will be used for
the endemic intercept and the Gaussian kernel, -9 for the epidemic
intercept, and zeroes for the remaining parameters.
Any initial values given in the start
argument
take precedence over those in par
.
Note that optim
receives the negative log-likelihood for
minimization (thus, if used, optim.args$control$fnscale
should be
positive). The hessian
argument defaults to TRUE
, and
in the control
list, trace
ing is enabled with
REPORT=1
by default. By setting
optim.args$control$trace = 0
, all output from the
optimization routine is suppressed.
For the partial
likelihood, the analytic score function and
the Fisher information are not implemented and the default is to use
robust method="Nelder-Mead"
optimization.
There may be an extra component fixed
in the
optim.args
list, which determines which parameters should stick
to their initial values. This can be specified by a
logical vector of the same length as the par
component, by an
integer vector indexing par
or by a character vector following
the twinstim
naming conventions. Furthermore, if
isTRUE(fixed)
, then all parameters are fixed at their initial
values and no optimization is performed.
Importantly, the method
argument in the optim.args
list may also be "nlminb"
, in
which case the nlminb
optimizer is used. This is also the
default for full likelihood inference.
In this case, not only the score function but also the
expected Fisher information can be used during optimization (as
estimated by what Martinussen and Scheike (2006, p. 64) call the
“optional variation process”, or see Rathbun (1996, equation
(4.7))). In our experience this gives better convergence than
optim
's methods.
For method="nlminb"
, the following parameters of the
optim.args$control
list may be named like for
optim
and are renamed appropriately:
maxit
(-> iter.max
), REPORT
(-> trace
,
default: 1), abstol
(-> abs.tol
), and
reltol
(-> rel.tol
, default: 1e-6
).
For nlminb
, a logical hessian
argument (default:
TRUE
) indicates if the negative expected Fisher
information matrix should be used as the Hessian during optimization
(otherwise a numerical approximation is used).
Similarly, method="nlm"
should also work but is not
recommended here.
- finetune
logical indicating if a second maximisation should be performed with
robust Nelder-Mead optim
using the resulting
parameters from the first maximisation as starting point. This
argument is only considered if partial = FALSE
and the
default is to not conduct a second maximization (in most cases this
does not improve upon the MLE).
- model
logical indicating if the model environment should be kept with the
result, which is required for
intensityplot
s and
R0(..., trimmed = FALSE)
.
Specifically, if model=TRUE
, the return value will have the
evaluation environment set as its environment
,
and the returned functions
element will contain the
log-likelihood function (or partial log-likelihood function, if
partial = TRUE
), and optionally the score and the expected
Fisher information functions (not for the partial likelihood, and
only if siaf
and tiaf
provide the necessary
derivatives).
Note that fitted objects with a model environment might consume
quite a lot of memory since they contain the data
.
- cumCIF
logical (default: FALSE
) indicating whether to
calculate the fitted cumulative ground intensity at event times.
This is the residual process, see residuals.twinstim
.
- cumCIF.pb
logical indicating if a progress bar should be shown
during the calculation of cumCIF
. Defaults to do so in an
interactive R session, and will be FALSE
if cores != 1
.
- cores
number of processes to use in parallel operation. By default
twinstim
runs in single-CPU mode. Currently, only the
multicore-type of parallel computing via forking is
supported, which is not available on Windows, see
mclapply
in package parallel.
Note that for a memoised siaf.step
kernel,
cores=1
is fixed internally since parallelization would slow
down model fitting significantly.
- verbose
logical indicating if information should be printed during
execution. Defaults to TRUE
.