[R] AR(1) with an error term arima.sim parameter question

[Rolf Turner]

Please see below.

On 10/12/14 20:21, Michael Selevan wrote:
> Hello,
>
> I am attempting to plot an AR(1) model with a standard deviation and I am a
> little confused as how to do that. I have been looking through the
> interwebs and some documentation and I see that there is potentially a few
> different ways to do this.
>
> First, simply using the documentation I came up with the command
>
> arima.sim(n=10, list(ar=0.8), innov=rnorm(10, sd=0.2))
>
> which would give me the standard deviation I want. Or I believe that to be
> the case. However, after some more searching and googling, I saw an example
> where someone used this as a means of adding the AR error term
>
> error.model=function(n){rnorm(n, sd=0.2)}
>
> y = arima.sim(n=10, list(ar=0.8), innov=rnorm(10, sd=0.2), rand.gen=
> error.model)
> Now, I am a little confused by this. Would having the error term in the
> innov parameter as well as the rand.gen be redundant? What would be the
> expected differences between the two? Should only 1 be used?
>
> Just looking for some clarification. Been searching and havent found too
> many examples that explicitly state how to add the error term to an AR(1)
> model.

It's a little bit subtle, but in a way that's not too important.

There is, in addition to "innov" a starting innovations vector 
"start.innov" that is needed.  If either innov or start.innov is not 
supplied their values get supplied by rand.gen().  So in your second 
call to arima.sim() ***start.innov*** is being supplied by rand.gen()
(but ***innov*** will be taken to be equal to the argument supplied.

In your first call, where rand.gen() is not specified (and start.innov
is not specified), the supplied value of innov will be used and 
start.innov will be produced by the *default* value of rand.gen()
which is rnorm(), you'll get rnorm(n.start,0,1).

Thus in your first call, the starting innovations will be done with a 
different standard deviation than the other innovations.  Which is 
probably not what you want.

Hence the second call is correct --- but it *is* kind of redundant and 
confusing to supply "innov" as well as rand.gen().  The code would be
clearer if "innov" were dispensed with and it was just left to 
rand.gen() to do the work.

The following is not important, but it might be mystifying:  If you 
leave out "innov" you will get a different result --- even if you set a 
seed for the random number generators a priori.  E.g.:

# 1.
set.seed(42)
innov <- rnorm(10,0,0.2)
error.model=function(n){rnorm(n, sd=0.2)}
y1 <- arima.sim(n=10, list(ar=0.8), innov=innov,
                 rand.gen=error.model)

# 2.
set.seed(42)
error.model=function(n){rnorm(n, sd=0.2)}
y2 <- arima.sim(n=10, list(ar=0.8),rand.gen=error.model)

The vectors y1 and y2 are (surprisingly until you think carefully) 
different.

This is because for y1, innov.start is generated *after* innov is
generated, whereas for y2 innov.start is generated *before* innov is
generated.  The first entry of innov for y1 will be the same as the
first entry of innov.start for y2.  So the sequence of innovations is
different.

Bottom line:  I would recommend *not* using the "innov" argument and
just specifying rand.gen() to get the standard deviations that you want.

HTH

cheers,

Rolf Turner

-- 
Rolf Turner
Technical Editor ANZJS