[R] State space AR models in R: some examples

Pablo Almaraz palmaraz at gmail.com
Thu Apr 27 16:55:19 CEST 2006


Hi all,

Does anyone have an example of an autoregressive (AR) time-series model 
specified as a state space model in R? That is, I want to go beyond the 
locally linear (constant) model, and fit the following Gaussian AR state 
process model:

Xt = a + (1+b)*Xt-1 + epsilon

,where the model for the observation process is

Yt = Xt + tau

I have information of the tau's (observation variance) for each 
observation in the time-series, and it would be perfect to include this 
information during the fitting routine. I have actually coded this as a 
WinBUGS code (pasted below), but I'm not quite sure it works as it 
should. I would be extremely thanked if anyone could submit an example 
of an R code fitting the above problem. Gibb's sampler for solving the 
problem would be great, I'm not sure whether the Kalman filter would 
work well with only 30 data points (?). Additional details, corrections 
and/or help would probably save my life at least for a while.

Thank you all
Cheers

Pablo

WinBUGS state-space R code:
######
model;
    {

# Parameters and priors
        alpha ~ dnorm(0,0.000001)   # Intrinsic rate of increase
        b ~ dnorm(0,0.000001)
        beta[1] <- b-1   # First-order density-dependence
        sigma ~ dunif(0, 1000)   # State process SD
        isigma2 <- pow(sigma, -2)   # State process 1/var
        # isigma2 ~ dgamma(0.01,0.000001)

# Initial state value
        n.exp[1] ~ dnorm(n[1],tau[1])

# State process model
        for(j in 1:(N-1)){
            n.exp.mu[j+1] <- alpha + b*n.exp[j]  # First-order Gompertz 
model
            n.exp[j+1] ~ dnorm(n.exp.mu[j+1], isigma2)
        }

# Observation process model
        for(j in 1:(N-1)){
            n[j+1] ~ dnorm(n.exp[j+1],tau[j+1])
        }
    }

# Loge-transformed and standardized time-series data

list(N=28,  
n=c(-0.24645, 0.015312, 0.442262, -0.05879, -0.17308, -0.03778, 
0.120961, -0.04383, 0.002507, 0.073278, -0.11684, 0.003657, -0.07375, 
0.05006, -0.04489, -0.00826, -0.06713, 0.682228, 0.032058, -0.33254, 
-0.50432, 0.176914, 0.249793, 0.01672, -0.30581, -0.19617, 0.158579, 
0.185296),
tau=c(2.38351, 2.351379, 49.12811, 10.01703, 11.68982, 3.846619, 
1.999254, 1.6685, 3.011932, 5.661051, 168.2524, 1.581, 25.74985, 
50.29332, 3.03117, 7.65013, 3.376606, 17.34871, 4.215985, 2.455294, 
7.685724, 1.918054, 5.588953, 8.503541, 0.5666, 0.923611, 4.986243, 
10.36613))

# Inits for MC 1
list(alpha = 0.5, b = -1, sigma = 0.5)

# Inits for MC 2
list(alpha = 1, b = 0.01, sigma = 1)

# Inits for MC 3
list(alpha = 0.01, b = 1, sigma = 0.01)

### End (not run)

-- 

Pablo Almaraz García
Estación Biológica de Doñana (CSIC)
Pabellón del Perú, Avda. Mª Luísa s/n
E-41013, Sevilla
SPAIN

E-mail: almaraz[AT]ebd[DOT]csic[DOT]es
webpage: http://www.almaraz.org




More information about the R-help mailing list