| Re: pth order auto-regressive process with a specified mean and variance [message #79204] |
Mon, 13 February 2012 09:11  |
David Fanning
Messages: 11724
Registered: August 2001
|
Senior Member |
|
|
Yngvar Larsen writes:
> P.S. This is statistical signal processing, so your math books might
> not be the proper place to look for gibberish-to-english
> translation :)
Ah, yeah, I was chasing one of those lithe modern dance
majors and skipped most of that class! ;-)
Cheers,
David
--
David Fanning, Ph.D.
Fanning Software Consulting, Inc.
Coyote's Guide to IDL Programming: http://www.idlcoyote.com/
Sepore ma de ni thui. ("Perhaps thou speakest truth.")
|
|
|
|
|
|
|
|
| Re: pth order auto-regressive process with a specified mean and variance [message #79210 is a reply to message #79209] |
Mon, 13 February 2012 07:04  |
Yngvar Larsen
Messages: 134
Registered: January 2010
|
Senior Member |
|
|
On Feb 13, 9:19 am, Tom Van Niel <mookiethe...@gmail.com> wrote:
> Hi Guys,
>
> Does anybody have IDL code that simulates a pth order auto-regressive
> process with a specified mean and variance? If so, please let me
> know.
>
Well. Assuming that the parameters of your AR(p) process is
chosen such that the process is wide-sense stationary, this should be
easy:
xmean = 0.9
xvariance = 1.2
;; This example vector of AR(4) coefficients
;; results in a WSS process.
;; Homework: make sure this is the case yourself.
;; (Hint: Roots of characteristic polynomial
;; within unit circle)
;; My model:
;; X_t = \sum_{i=1}^p a_iX_{t-i} + n_t
;; with n_t iid normal.
;; Array A below contains (after reverse)
;; A = [a_{p-1}, ..., a_2, a_1]
A = reverse([2.7607, -3.8106, 2.6535, -0.9238])
p = n_elements(A)
transient = 1000 ; Transient throwaway points
npoints = 10000 + transient
drive_proc = sqrt(xvariance)*randomn(seed, npoints)
ar_proc = fltarr(npoints)
ar_proc[0] = drive_proc[0]
for ii=1, p-1 do $
ar_proc[ii] = drive_proc[ii] + $
total(A[p-ii:*]*ar_proc[0:ii-1])
for ii=p, npoints-1 do $
ar_proc[ii] = drive_proc[ii] + $
total(A*ar_proc[ii-p:ii-1])
;; Remove transient points, where
;; the process isn't WSS yet.
;; Add mean value.
ar_proc = ar_proc[transient:*] + xmean
--
Yngvar
|
|
|
|
| Re: pth order auto-regressive process with a specified mean and variance [message #79300 is a reply to message #79210] |
Mon, 13 February 2012 17:52 |
Tom Van Niel
Messages: 4
Registered: February 2012
|
Junior Member |
|
|
Hi Yngvar,
That is just what I need. At the moment, I'm only using AR(1), so as
long as I keep the parameter < 1, it should remain WSS. If I use
AR(2) or greater, then I'll see if I can figure out the homework you
have assigned :)
David, no need to look up your textbooks when Wikipedia is around,
although it sounds like it brought back some good memories.
Thanks
Yngvar Larsen wrote:
> On Feb 13, 9:19 am, Tom Van Niel <mookiethe...@gmail.com> wrote:
>> Hi Guys,
>>
>> Does anybody have IDL code that simulates a pth order auto-regressive
>> process with a specified mean and variance? If so, please let me
>> know.
>>
>
> Well. Assuming that the parameters of your AR(p) process is
> chosen such that the process is wide-sense stationary, this should be
> easy:
>
> xmean = 0.9
> xvariance = 1.2
> ;; This example vector of AR(4) coefficients
> ;; results in a WSS process.
> ;; Homework: make sure this is the case yourself.
> ;; (Hint: Roots of characteristic polynomial
> ;; within unit circle)
> ;; My model:
> ;; X_t = \sum_{i=1}^p a_iX_{t-i} + n_t
> ;; with n_t iid normal.
> ;; Array A below contains (after reverse)
> ;; A = [a_{p-1}, ..., a_2, a_1]
> A = reverse([2.7607, -3.8106, 2.6535, -0.9238])
> p = n_elements(A)
> transient = 1000 ; Transient throwaway points
> npoints = 10000 + transient
>
> drive_proc = sqrt(xvariance)*randomn(seed, npoints)
> ar_proc = fltarr(npoints)
> ar_proc[0] = drive_proc[0]
> for ii=1, p-1 do $
> ar_proc[ii] = drive_proc[ii] + $
> total(A[p-ii:*]*ar_proc[0:ii-1])
>
> for ii=p, npoints-1 do $
> ar_proc[ii] = drive_proc[ii] + $
> total(A*ar_proc[ii-p:ii-1])
>
> ;; Remove transient points, where
> ;; the process isn't WSS yet.
> ;; Add mean value.
> ar_proc = ar_proc[transient:*] + xmean
>
> --
> Yngvar
|
|
|
|
comp.lang.idl-pvwave archive
Members
Search
Help
Login
Home
