BaPreStoPro (version 0.1)

estimate,jumpRegression-method: Estimation for regression model dependent on Poisson process

Description

Bayesian estimation of the parameter of the regression model $y_i = f(t_i, N_{t_i}, \theta) + \epsilon_i$ with $N_t\sim Pois(\Lambda(t, \xi)), \epsilon_i\sim N(0,\gamma^2\widetilde{s}(t))$.

Usage

"estimate"(model.class, t, data, nMCMC, propSd, adapt = TRUE, proposal = c("normal", "lognormal"), it.xi = 10)

Arguments

model.class

class of the regression model based on the NHPP including all required information, see jumpRegression-class

t

vector of time points

data

vector of observation variables

nMCMC

length of Markov chain

propSd

vector of proposal variances for $(\theta, \xi)$

adapt

if TRUE (default), proposal variance is adapted

proposal

proposal density for $\theta$: "normal" (default) or "lognormal" (for positive parameters)

it.xi

number of iterations for MH step for $\xi$ inside the Gibbs sampler

Proposal densities

For $\theta$, there is the possibility to choose "normal" or "lognormal". The proposal density for $\xi$ depends on the starting value of $\xi$. If all components are positive, the proposal density is lognormal, and normal otherwise.

References

Heeke, G., S. Hermann, R. Maurer, K. Ickstadt, and C. H. Mueller (2015). Stochastic Modeling and Statistical Analysis of Fatigue Tests on Prestressed Concrete Beams under Cyclic Loadings. SFB 823 discussion paper 25/15.

Examples

t <- seq(0,1, by = 0.01)
model <- set.to.class("jumpRegression", fun = function(t, N, theta) exp(theta[1]*t) + theta[2]*N,
                   parameter = list(theta = c(2, 2), gamma2 = 0.25, xi = c(3, 0.5)),
                   Lambda = function(t, xi) (t/xi[2])^xi[1])
data <- simulate(model, t = t, plot.series = FALSE)
est <- estimate(model, t, data, 1000)
plot(est)
## Not run: 
# # work in progress
# est_hid <- estimate(model, t, data$Y, 1000)
# plot(est_hid)
# ## End(Not run)