RE: generate AR(1) data

[Mike Miller]

If it wasn't clear -- Ted and I are saying the same thing about how to use
Cholesky factorization but we're using slightly different notation.  All
I'm adding is that there is a chol() function in Octave.  --Mike


On Thu, 17 Jul 2003 address@hidden wrote:

> On 17-Jul-03 Heber Farnsworth wrote:
> > I need to generate a large number of random (Gaussian) vectors that
> > have an AR(1) correlation structure.  Obviously I could do this in a
> > loop but this is slow.  I would rather rotate a vector of IID normal
> > draws so that the result will have the right covariance matrix.  So
> > what I need is the matrix square root of a symmetric toeplitz matrix
> >
> > 1   .5  .25  .125 ...
> > .5   1  .5    .25 ...
> > .25 .5  1   .5 ...
> >
> > etc (here the correlation is .5).  It seems inefficient to try to
> > compute this square root numerically.  Isn't it known in closed form?
>
> You could try the Choleski factorisation of T into T=V'*V where V is
> upper triangular. In this case the top row of V is the same as that
> of T; each of the remaining rows has zeros on the left of the main
> diagonal, the on-diagonal term is sqrt(3/4), and subsequent terms are
> all 1/2 of the one on their left.
>
> In general, where your original T has top row 1, r, r^2, r^3, ... ,
> the top row of V is as in T, and (reading rightwards from the
> diagonal term) each subsequent row is of the form a*(1, r, r^2, r^3, ... )
> where a^2 = 1 - r^2. Terms to the left of the diagonal are zero. If you
> work it out (using a^2 = 1 - r^2) you can easily see that
> for i <= j, [V'*V](i,j) = r^(j-i) and of course V'*V is symmetric.
>
> Once you have V, let U be a ROW vector of independent N(0,1)s.
> The covariance matrix of X = U*V is E(X'*X) = V'*E(U'*U)*V = V'*V = T.
>
> Hope this helps!
> Ted.



-------------------------------------------------------------
Octave is freely available under the terms of the GNU GPL.

Octave's home on the web:  http://www.octave.org
How to fund new projects:  http://www.octave.org/funding.html
Subscription information:  http://www.octave.org/archive.html
-------------------------------------------------------------