[R] Simulate an AR(1) process via distributions? (without specifying a model specification)

[Prof Brian Ripley]

We provide arima.sim() to simulate from stationary ARIMA processes, and I 
had guessed you were trying to use that.

But what your code did was to simulate a non-stationary Gaussian AR(1), 
fit a marginal distribution to it as if it were iid (it is neither 'i') 
and then simulate iid samples from the fitted marginal and fit an ARIMA 
model to that (and arima() assumes a stationary Gaussian ARIMA process).

I could not reproduce the first step: the AR(1) simulation exploded to 
numbers around 10^175 and arima() got numerical overflow (the sum of 
squares overflowed to Inf).

You wrote x = 0.05 + 0.64*x(t-1) + e but simulated
x = 0.05 + 1.50*x(t-1) + e.  Which did you intend?  If the former, the 
marginal of a stationary Gaussian AR(1) is Gaussian, so you don't need a 
GLD to model it.

I don't understand what you are trying to do, but you are violating the 
assumptions of all of the tools you are using.  This is not an R issue, 
so please seek local statistical expertise to help you.


On Wed, 28 Nov 2007, Pedro.Rodriguez at sungard.com wrote:

> Thanks Prof. Ripley.
>
> My apologies for not including the code.
>
> Below I illustrate my point using the GLD package.
>
> Thank you very much for your time.
>
> Kind Regards,
>
> Pedro N. Rodriguez
>
>
> # Code begins
>
> # Simulate an ar(1) process
> # x = 0.05 + 0.64*x(t-1) + e
>
> # Create the vector x
> 	x	<- vector(length=1000)
>
> #simulate the own risk
> 	e	<- rnorm(1000)
>
> #Set the coefficient
> 	beta 	<- 1.50
>
> # set an initial value
> 	x[1] 	<- 5
>
> #Fill the vector x
> 	for(i in 2:length(x))
> 	{
> 		x[i]	<- 0.05 + beta*x[i-1] + e[i]
> 	}
>
> #Check the AR(1)
> 	simulated_data_ar <- arima(x,order=c(1,0,0))
> 	simulated_data_ar
>
> #Using the G Lambda Distribution to fit the distribution.
> 	library(gld)
> 	resul1		<- starship(x,optim.method="Nelder-Mead")
> 	lambdas1	<- resul1$lambda
>
> #Plot the Distribution
> 	plotgld(lambdas1[1],lambdas1[2],lambdas1[3],lambdas1[4])
>
> #Random Deviates from GLD
> 	x_sim		<-
> rgl(1000,lambdas1[1],lambdas1[2],lambdas1[3],lambdas1[4])
>
> #Fit an AR(1)
> 	gld_simulated	<- arima(x_sim,order=c(1,0,0))
> 	gld_simulated
>
> #Code ends
>
>
> -----Original Message-----
> From: Prof Brian Ripley [mailto:ripley at stats.ox.ac.uk]
> Sent: Wednesday, November 28, 2007 11:37 AM
> To: Rodriguez, Pedro
> Cc: r-help at stat.math.ethz.ch
> Subject: Re: [R] Simulate an AR(1) process via distributions? (without
> specifying a model specification)
>
> On Wed, 28 Nov 2007, Pedro.Rodriguez at sungard.com wrote:
>
>> Is it possible to simulate an AR(1) process via a distribution?
>
> Any distribution *of errors*, yes.  Of the process values, not in
> general.
>
>> I have simulated an AR(1) process the usual way (that is, using a
> model
>> specification and using the random deviates in the error), and used
> the
>> generated time series to estimate 3- and 4-parameter distributions
> (for
>> instance, GLD). However, the random deviates generated from these
>> distributions do not follow the specified AR process.
>
> How do you know that?  Please give us the reproducible example we asked
> for (in the posting guide, at the bottom of every message), and we
> should
> be able to explain it to you.
>
>

-- 
Brian D. Ripley,                  ripley at stats.ox.ac.uk
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford,             Tel:  +44 1865 272861 (self)
1 South Parks Road,                     +44 1865 272866 (PA)
Oxford OX1 3TG, UK                Fax:  +44 1865 272595