Description Usage Arguments Details Examples
Creates a marked point process model object with class "mpp"
.
1 |
data |
a |
gif |
ground intensity function. See topic |
marks |
a |
params |
numeric vector of all model parameters. |
gmap |
|
mmap |
|
TT |
vector of length 2, being the time interval over which the integral of the ground intensity function is to be evaluated. |
Let lambda_g(t|Ht) denote the ground intensity function and f(y|Ht) denote the joint mark densities, where y in Y. The log-likelihood of a marked point process is given by
log L = sum{log lambda_g(ti|Ht)} + sum{log f(yi|Ht)} - integral{lambda_g(t|Ht) dt},
where the summation is taken over those events contained in the interval (TT[1], TT[2])
, and the integral is also taken over that interval. However, all events in the data frame data
before t, even those before TT[1]
, form the history of the process Ht. This allows an initial period for the process to reach a “steady state” or “equilibrium”.
The parameter spaces of the ground intensity function and mark distribution are not necessarily disjoint, and can have common parameters. Hence, when the model parameters are estimated, these relationships must be known, and are specified by the arguments gmap
and mmap
. The mapping expressions can also contain arithmetic expressions. The ith element in the params
argument is addressed in the expressions as params[i]
. Here is an example of a five parameter model, where the gif
has 4 parameters, and the mark distribution has 2, with mappings specified as:
1 2 3 | gmap = expression(c(params[1:3], exp(params[4]+params[5])))
mmap = expression(c(log(params[2]/3), params[5]))
|
Note the inclusion of the combine (c
) function, because the expression
must create a vector of parameters. Care must be taken specifying these expressions as they are embedded directly into the code of various functions.
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 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 | data(Tangshan)
# increment magnitudes a fraction so none are zero
Tangshan[,"magnitude"] <- Tangshan[,"magnitude"] + 0.01
dmagn_mark <- function(x, data, params){
# Gamma distribution
# exponential density when params[7]=0
# See topic "marks" for further discussion
lambda <- etas_gif(data, x[,"time"], params=params[1:5])
y <- dgamma(x[,"magnitude"], shape=1+sqrt(lambda)*params[7],
rate=params[6], log=TRUE)
return(y)
}
TT <- c(0, 4018)
# params <- c(0.0067, 1.1025, 1.0794, 0.0169, 0.9506, 1.9159, 0.4704)
params <- c(0.007, 1.1, 1.08, 0.02, 0.95, 1.92, 0.47)
x <- mpp(data=Tangshan,
gif=etas_gif,
marks=list(dmagn_mark, NULL),
params=params,
gmap=expression(params[1:5]),
mmap=expression(params[1:7]),
TT=TT)
allmap <- function(y, p){
# one to one mapping, all p positive
y$params <- exp(p)
return(y)
}
# Parameters must be positive. Transformed so that nlm
# can use entire real line (no boundary problems, see
# topic "neglogLik" for further explanation).
# Argument "iterlim" has been restricted to 2 to avoid
# excessive time in package checks, set much larger to
# ensure convergence.
z <- nlm(neglogLik, log(params), object=x, pmap=allmap,
print.level=2, iterlim=2, typsize=abs(params))
x1 <- allmap(x, z$estimate)
# print parameter estimates
print(x1$params)
print(logLik(x))
print(logLik(x1))
plot(x1, log=TRUE)
|
iteration = 0
Step:
[1] 0 0 0 0 0 0 0
Parameter:
[1] -4.96184513 0.09531018 0.07696104 -3.91202301 -0.05129329 0.65232519
[7] -0.75502258
Function Value
[1] 1133.502
Gradient:
[1] 3.596469 72.437753 111.530955 64.783572 -441.385097 -15.487562
[7] 9.160502
iteration = 1
Step:
[1] -1.121877e-08 -5.579857e-03 -8.281627e-03 -1.649672e-06 2.535932e-02
[6] 3.634613e-03 -1.288213e-04
Parameter:
[1] -4.96184514 0.08973032 0.06867941 -3.91202466 -0.02593398 0.65595980
[7] -0.75515141
Function Value
[1] 1128.627
Gradient:
[1] -1.295572 -5.179578 -4.366239 -6.288302 87.943830 1.277122 1.222243
iteration = 2
Parameter:
[1] -4.96184514 0.09014136 0.06909796 -3.91202450 -0.03011746 0.65566359
[7] -0.75516078
Function Value
[1] 1128.416
Gradient:
[1] -0.4592951 5.5840238 11.2814851 3.8121255 13.0624118 -0.9085600 2.1820004
Iteration limit exceeded. Algorithm failed.
[1] 0.00700000 1.09432897 1.07154117 0.01999997 0.97033156 1.92642045 0.46993505
[1] -1133.502
[1] -1128.416
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.