rst | R Documentation |

This function simulates spatial and temporal random effects with mean zero. The method is described in Algorithm 3.1 of Rue & Held 2015.

rst( n = 1, type = c("s", "t", "st")[1], type.s = "ICAR", type.t = c("RW1", "RW2")[2], Amat = NULL, n.t = NULL, scale.model = TRUE )

`n` |
sample size |

`type` |
type of random effects: temporal (t), spatial (s), or spatial-temporal (st) |

`type.s` |
type of spatial random effect, currently only ICAR is available |

`type.t` |
type of temporal random effect, currently only RW1 and RW2 are available |

`Amat` |
adjacency matrix for the spatial regions |

`n.t` |
number of time points for the temporal random effect |

`scale.model` |
logical indicator of whether to scale the random effects to have unit generalized variance. See Sørbye 2013 for more details |

a matrix (for spatial or temporal) or a three-dimensional array (for spatial-temporal) of the random effects.

Zehang Richard Li

Rue, H., & Held, L. (2005). *Gaussian Markov random fields: theory and applications*. CRC press.

Sørbye, S. H. (2013). *Tutorial: Scaling IGMRF-models in R-INLA*. Department of Mathematics and Statistics, University of Tromsø.

## Not run: data(DemoMap) ## Spatial random effects out <- rst(n=10000, type = "s", Amat = DemoMap$Amat) # To verify the mean under the conditional specification mean(out[,1] - apply(out[,c(2,3,4)], 1, mean)) mean(out[,2] - apply(out[,c(1,3)], 1, mean)) mean(out[,3] - apply(out[,c(1,2,4)], 1, mean)) mean(out[,4] - apply(out[,c(1,3)], 1, mean)) ## Temporal random effects (RW1) out <- rst(n=1, type = "t", type.t = "RW1", n.t = 200, scale.model = FALSE) par(mfrow = c(1,2)) plot(1:dim(out)[2], out, col = 1, type = "l", xlab = "Time", ylab = "Random effects") # verify the first order difference is normally distributed first_diff <- diff(as.numeric(out[1,])) qqnorm(first_diff ) abline(c(0,1)) ## Temporal random effects (RW2) out <- rst(n=1, type = "t", type.t = "RW2", n.t = 200, scale.model = FALSE) par(mfrow = c(1,2)) plot(1:dim(out)[2], out, col = 1, type = "l", xlab = "Time", ylab = "Random effects") # verify the second order difference is normally distributed first_diff <- diff(as.numeric(out[1,])) second_diff <- diff(first_diff) qqnorm(second_diff) abline(c(0,1)) ## Spatial-temporal random effects out <- rst(n=1, type = "st", type.t = "RW2", Amat = DemoMap$Amat, n.t = 50) dimnames(out) par(mfrow = c(1,1)) plot(1:dim(out)[3], out[1,1,], col = 1, type = "l", ylim = range(out), xlab = "Time", ylab = "Random effects") for(i in 2:4) lines(1:dim(out)[3], out[1,i,], col = i) legend("bottomright", colnames(DemoMap$Amat), col = c(1:4), lty = rep(1,4)) ## End(Not run)

Embedding an R snippet on your website

Add the following code to your website.

For more information on customizing the embed code, read Embedding Snippets.