# NOT RUN {
########## Simulating a multivariate damped random walk process
n <- 100
k <- 5
obs.time <- cumsum(rgamma(n, shape = 3, rate = 1))
tau <- 100 + 20 * (1 : 5) #rnorm(k, 0, 5)
sigma <- 0.01 * (1 : 5)
#tau <- c(1 : 5) #rnorm(k, 0, 5)
#sigma <- 0.05 + 0.007 * (0 : 4) #rnorm(k, 0, 0.002)
mu <- 17 + 0.5 * (1 : 5)
rho.m <- matrix(0, k, k)
for(i in 1 : k) {
for(j in 1 : k) {
rho.m[i, j] = 1.1^(-abs(i - j))
}
}
rho <- rho.m[upper.tri(rho.m)]
measure.error.band <- c(0.010, 0.014, 0.018, 0.022, 0.026)
measure.error <- NULL
for(i in 1 : k) {
measure.error <- cbind(measure.error, rnorm(n, measure.error.band[i], 0.002))
}
x <- drw.sim(time = obs.time, n.datasets = 5, measure.error.SD = measure.error,
mu = mu, sigma = sigma, tau = tau, rho = rho)
plot(obs.time, x[, 1], xlim = c(min(obs.time), max(obs.time)), ylim = c(17, 20),
xlab = "time", ylab = "observation")
points(obs.time, x[, 2], col = 2, pch = 2)
points(obs.time, x[, 3], col = 3, pch = 3)
points(obs.time, x[, 4], col = 4, pch = 4)
points(obs.time, x[, 5], col = 5, pch = 5)
########## Simulating a univariate damped random walk process
x <- drw.sim(time = obs.time, n.datasets = 1, measure.error.SD = measure.error[, 1],
mu = mu[1], sigma = sigma[1], tau = tau[1])
plot(obs.time, x)
# }
Run the code above in your browser using DataLab