spate.mcmc: MCMC algorithm for fitting the model.

In spate: Spatio-Temporal Modeling of Large Data Using a Spectral SPDE Approach

MCMC algorithm for fitting the model.

Description

MCMC algorithm for fitting the model.

Usage

spate.mcmc(y,coord=NULL,lengthx=NULL,lengthy=NULL,Sind=NULL,n=NULL,
          IncidenceMat=FALSE,x=NULL,SV=c(rho0=0.2,sigma2=0.1,
          zeta=0.25,rho1=0.2,gamma=1,alpha=0.3,muX=0,muY=0,tau2=0.005),
          betaSV=rep(0,dim(x)[1]),RWCov=NULL,parh=NULL,tPred=NULL,
          sPred=NULL,P.rho0=Prho0,P.sigma2=Psigma2,P.zeta=Pzeta,P.rho1=Prho1,
          P.gamma=Pgamma,P.alpha=Palpha,P.mux=Pmux,P.muy=Pmuy,P.tau2=Ptau2,
          lambdaSV=1,sdlambda=0.01,P.lambda=Plambda,DataModel="Normal",
          DimRed=FALSE,NFour=NULL,indEst=1:9,Nmc=10000,BurnIn =1000,
          path=NULL,file=NULL,SaveToFile=FALSE,PlotToFile=FALSE,
          FixEffMetrop=TRUE,saveProcess=FALSE,Nsave=200,seed=NULL,
          Padding=FALSE,adaptive=TRUE,NCovEst=500,BurnInCovEst=500,
          MultCov=0.5,printRWCov=FALSE,MultStdDevLambda=0.75,
          Separable=FALSE,Drift=!Separable,Diffusion=!Separable,
          logInd=c(1,2,3,4,5,9),nu=1,plotTrace=TRUE,
          plotHist=FALSE,plotPairs=FALSE,trueVal=NULL,
          plotObsLocations=FALSE,trace=TRUE,monitorProcess=FALSE,
          tProcess=NULL,sProcess=NULL)

Arguments

y	Observed data in an T x N matrix with columns and rows corresponding to time and space (observations on a grid stacked into a vector), respectively. By default, at each time point, the observations are assumed to lie on a square grid with each axis scaled so that it has unit length.
coord	If specified, this needs to be a matrix of dimension N x 2 with coordinates of the N observation points. Observations in 'y' can either be on a square grid or not. If not, the coordinates of each observation point need to be specified in 'coord'. According to these coordinates, each observation location is then mapped to a grid cell. If 'coord' is not specified, the observations in 'y' are assumed to lie on a square grid with each axis scaled so that it has unit length.
lengthx	Use together with 'coord' to specify the length of the x-axis. This is usefull if the observations lie in a rectangular area instead of a square. The length needs to be at least as large as the largest x-distance in 'coord.
lengthy	Use together with 'coord' to specify the length of the y-axis. This is usefull if the observations lie in a rectangular area instead of a square. The length needs to be at least as large as the largest y-distance in 'coord.
Sind	Vector of indices of grid cells where observations are made, in case, the observation are not made at every grid cell. Alternatively, the coordinates of the observation locations can be specfied in 'coord'.
n	Number of point per axis of the square into which the points are mapped. In total, the process is modeled on a grid of size n*n.
IncidenceMat	Logical; if 'TRUE' an incidence matrix relating the latent process to observation locations is used. This is only recommended to use when the observations are relatively low-dimensional and when the latent process is modeled in a reduced dimensional space as well.
x	Covariates in an array of dimensions p x T X N, where p denotes the number of covariates, T the number of time points, and N the number of spatial points.
SV	Starting values for parameters. Parameters for the SPDE in the following order: rho_0, sigma^2, zeta, rho_1, gamma, alpha, mu_x, mu_y, tau^2. rho_0 and sigma^2 are the range and marginal variance of the Matern covariance funtion for the innovation term epsilon. zeta is the damping parameter. rho_1, gamma, and alpha parametrize the diffusion matrix with rho_1 being a range parameter, gamma and alpha determining the amount and the direction, respectively, of anisotropy. mu_x and mu_y are the two components of the drift vector. tau^2 denotes the nugget effect or measurment error.
betaSV	Starting values for regression coefficients.
RWCov	Covariance matrix of the proposal distribution used in the random walk Metropolis-Hastings step for the hyper-parameters.
parh	Only used in prediction mode. If 'parh' is not 'NULL', this indicates that 'spate.mcmc' is used for making predictions at locations (tPred,sPred) instead of applying the traditional MCMC algorithm. In case 'parh' is not 'NULL', it is a Npar x Nsim matrix containing Nsim samples from the posterior of the Npar parameters. This argument is used by the wrapper function 'spate.predict'.
tPred	Time points where predictions are made.This needs to be a vector if predictions are made at multiple times. For instance, if T is the number of time points in the data 'y', then tPred=c(T+1, T+2) means that predictions are made at time 'T+1' and 'T+2'. This argument is used by the wrapper function 'spate.predict'.
sPred	Vector of indices of grid cells (positions of locations in the stacked spatial vector) where predictions are made. This argument is used by the wrapper function 'spate.predict'.
P.rho0	Function specifying the prior for rho0.
P.sigma2	Function specifying the prior for sigma2.
P.zeta	Function specifying the prior for zeta.
P.rho1	Function specifying the prior for rho1.
P.gamma	Function specifying the prior for gamma.
P.alpha	Function specifying the prior for alpha.
P.mux	Function specifying the prior for mux.
P.muy	Function specifying the prior for muy.
P.tau2	Function specifying the prior for tau2.
lambdaSV	Starting value for transformation parameter lambda in the Tobit model.
sdlambda	Standard deviation of the proposal distribution used in the random walk Metropolis-Hastings step for lambda.
P.lambda	Function specifying the prior for lambda.
DataModel	Specifies the data model. "Normal" or "SkewTobit" are available options.
DimRed	Logical; if 'TRUE' dimension reduction is applied. This means that not the full number (n*n) of Fourier functions is used but rather only a reduced dimensional basis of dimension 'NFour'.
NFour	If 'DimRed' is 'TRUE', this specifies the number of Fourier functions.
indEst	A vector of numbers specifying which for which parameters the posterior should be computed and which should be held fix (at their starting value). If the corresponding to the index of rho_0, sigma^2, zeta, rho_1, gamma, alpha, mu_x, mu_y, tau^2 is present in the vector, the parameter will be estimated otherwise not. Default is indEst=1:9 which means that one samples from the posterior for all parameters.
Nmc	Number of MCMC samples.
BurnIn	Length of the burn-in period.
path	Path, in case plots and / or the spateMCMC object should be save in a file.
file	File name, in case plots and / or the spateMCMC object should be save in a file.
SaveToFile	Indicates whether the spateMCMC object should be save in a file.
PlotToFile	Indicates whether the MCMC output analysis plots should be save in a file.
FixEffMetrop	The fixed effects, i.e., the regression coefficients, can either be sampled in a Gibbs step or updated together with the hyperparameters in the Metropolis-Hastings step. The latter is the default and recommended option since correlations between fixed effects and the random process can result in slow mixing.
saveProcess	Logical; if 'TRUE' samples from the posterior of the latent spatio-temporal process xi are saved.
Nsave	Number of samples from the posterior of the latent spatio-temporal process xi that should be save.
seed	Seed for random generator.
Padding	Indicates whether padding is applied or not. If the range parameters are large relative to the domain, this is recommended since otherwise spurious periodicity can occur.
adaptive	Indicates whether an adaptive Metropolis-Hastings algorithm is used or not. If yes, the proposal covariance matrix 'RWCov' is adaptively estimated during the algorithm and tuning does not need to be done by hand.
NCovEst	Minimal number of samples to be used for estimating the proposal matrix.
BurnInCovEst	Burn-in period for estimating the proposal matrix.
MultCov	Numeric used as multiplier for the adaptively estimated proposal cocariance matrix 'RWCov' of the hyper-parameters. I.e., the estimated covariance matrix is multiplied by 'MultCov'.
printRWCov	Logical, if 'TRUE' the estimated proposal cocariance matrix is printed each time.
MultStdDevLambda	Numeric used as multiplier for the adaptively estimated proposal standard deviation of the Tobit transformation parameter lambda. I.e., the estimated standard deviation is multiplied by 'MultStdDevLambda'.
Separable	Indicates whether a separable model, i.e., no transport / drift and no diffusion, should be estimated.
Drift	Indicates whether a drift term should be included.
Diffusion	Indicates whether a diffusion term should be included.
logInd	Indicates which parameters are sampled on the log-scale. Default is logInd=c(1, 2, 3, 4, 5, 9) corresponding to rho_0, sigma2, zeta, rho_1, gamma, and tau^2.
nu	Smoothness parameter of the Matern covariance function for the innovations. By default this equals 1 corresponding to the Whittle covariance function.
plotTrace	Indicates whether trace plots are made.
plotHist	Indicates whether histograms of the posterior distributions are made.
plotPairs	Indicates whether scatter plots of the hyper-parameters and the regression coefficients are made.
trueVal	In simulations, true values can be supplied for comparison with the MCMC output.
plotObsLocations	Logical; if 'TRUE' the observations locations are ploted together with the grid cells.
trace	Logical; if 'TRUE' tracing information on the progress of the MCMC algorithm is produced.
monitorProcess	Logical; if 'TRUE' in addition to the trace plots of the hyper-parameters, the mixing properties of the latent process xi=Phi*alpha is monitored. This is done by plotting the current sample of the process. More specifically, the time series at locations 'sProcess' and the spatial fieldd at time points 'tProcess'.
tProcess	To be secified if 'monitorProcess=TRUE'. Time points at which spatial fields of the sampled process should be plotted.
sProcess	To be secified if 'monitorProcess=TRUE'. Locations at which time series of the sampled process should be plotted.

Value

The function returns a 'spateMCMC' object with, amongst others, the following entries

Post	Matrix containing samples from the posterior of the hyper-parameters and the regression coefficient
xiPost	Array with samples from the posterior of the spatio-temporal process
RWCov	(Estimated) proposal covariance matrix

Author(s)

Fabio Sigrist

Examples

##Specify hyper-parameters
par <- c(rho0=0.1,sigma2=0.2,zeta=0.5,rho1=0.1,gamma=2,alpha=pi/4,muX=0.2,muY=-0.2,tau2=0.01)
##Simulate data
spateSim <- spate.sim(par=par,n=20,T=20,seed=4)
w <- spateSim$w

##Below is an example to illustrate the use of the MCMC algorithm.
##In practice, more samples are needed for a sufficiently large effective sample size.

##The following takes a couple of minutes.
##Load the precomputed object some lines below to save time.
##spateMCMC <- spate.mcmc(y=w,x=NULL,SV=c(rho0=0.2,sigma2=0.1,
##              zeta=0.25,rho1=0.2,gamma=1,alpha=0.3,muX=0,muY=0,tau2=0.005),
##              RWCov=diag(c(0.005,0.005,0.05,0.005,0.005,0.001,0.0002,0.0002,0.0002)),
##              Nmc=10000,BurnIn=2000,seed=4,Padding=FALSE,plotTrace=TRUE,NCovEst=500,
##              BurnInCovEst=500,trueVal=par,saveProcess=TRUE)
##spateMCMC
##plot(spateMCMC.fit,true=par,postProcess=TRUE)

##Instead of waiting, you can also use this precomputed object
data("spateMCMC")
spateMCMC
plot(spateMCMC,true=par,medianHist=FALSE)

spate documentation built on Oct. 3, 2023, 5:09 p.m.