# R/etas_gif.R In PtProcess: Time Dependent Point Process Modelling

#### Documented in etas_gif

```etas_gif <-
function(data, evalpts, params, TT = NA, tplus = FALSE)
{
#----------------------------------------------------------------------
# calculate the conditional intensity at given times as a function
# of past event history up to each time.
#
# data           = converted catalogue data (time, magnitude matrix)
# eval.times     = time instances at which conditional intensity is to be
#                  calculated
# params         = vector of parameters, with names attributes
#                  described below
#  mu            = background rate
#  A             = proportionality constant
#  alpha         = sensitivity coefficient
#  CC            = time constant
#  P             = decay coefficient
#  TT		 = time to evaluate, if not NA, only integral term evaluated
#----------------------------------------------------------------------
mu <- params[1]
A <- params[2]
alpha <- params[3]
CC <- params[4]
P <- params[5]
#------------------------------------------------------------
if (is.null(data)) data <- cbind(time = Inf, magnitude = 0)
magnitudes <- data[, "magnitude"]
times <- data[, "time"]
if(any(is.na(TT))) {
#----------------------------------------------------
#  Calculate the cumulative sum, S(t), say. Can only
#  use apply if calculating conditional intensity for
#  more than one time otherwise simply calculate directly
#  if TT is not missing, then just evaluate the integral term
#---------------------------------------------------
if(is.vector(evalpts)) eval.times <- evalpts
else eval.times <- evalpts[, "time"]
if(length(eval.times) > 1) {
St <- apply(t(eval.times), 2, function(z, magnitudes,
times, mu, A, P, CC, alpha, tplus)
{
if (!tplus) use <- times < z
else use <- times <= z
if(sum(use) > 0)
sum(exp(alpha * magnitudes[use]) * (1 + (z -
times[use])/CC)^( - P))
else 0
}
, magnitudes = magnitudes, times = times, mu = mu,
A = A, P = P, CC = CC, alpha = alpha,
tplus = tplus)
}
else {
if (!tplus) use <- times < eval.times
else use <- times <= eval.times
if(sum(use) > 0)
St <- sum(exp(alpha * magnitudes[use]) * (1 + (
eval.times - times[use])/CC)^( - P))
else St <- 0
}
ci <- mu + A * St
}
else {
#  if TT is not missing, then just evaluate the integral term
S <- c(NA, NA)
for(Ii in 1:2) {
Z0 <- (times < TT[Ii])
if(any(Z0)) {
if(P != 1)
S[Ii] <- sum(exp(alpha * data[, "magnitude"][
Z0]) * (1 - (1 + (TT[Ii] - times[Z0])/CC)^( -
P + 1)))/(P - 1)
else S[Ii] <- sum(exp(alpha * data[, "magnitude"
][Z0]) * log(1 + (TT[Ii] - times[Z0])/CC))
}
else S[Ii] <- 0
}
ci <- mu * (TT[2] - TT[1]) + CC * A * (S[2] - S[1])
}
names(ci) <- NULL
return(ci)
}
attr(etas_gif, "rate") <- "decreasing"
```

## Try the PtProcess package in your browser

Any scripts or data that you put into this service are public.

PtProcess documentation built on Nov. 17, 2017, 7:12 a.m.