Hi, I simulate an ARMA (3,3) in Matlab/Simulink: a_t-0.9391a_(t-1)+2.8763a_(t-2)-2.9372a_(t-3)=0.0376e_t+0.1127e_(t-1) +0.1127e_(t-2)+0.0376e_(t-3) First, it is transformed: a_t=0.9391a_(t-1)-2.8763a_(t-2)+2.9372a_(t-3)+0.0376e_t+0.1127e_(t-1) +0.1127e_(t-2)+0.0376e_(t-3) According to the source paper, e_t~N(0, 0.01^2), I find that the output a_t has extremely huge magnitude signal after 1000 or so noise samples. The delay units initial values are zero, but I think it is not a big problem to this IIR yet for 0 value now. Could you tell me what is wrong? Thanks

# What is wrong with my ARMA implementation (IIR)?

Started by ●June 24, 2014

Reply by ●June 24, 20142014-06-24

On Tue, 24 Jun 2014 16:19:32 -0700, fl wrote:> Hi, > I simulate an ARMA (3,3) in Matlab/Simulink: > > a_t-0.9391a_(t-1)+2.8763a_(t-2)-2.9372a_(t-3)=0.0376e_t+0.1127e_(t-1) > +0.1127e_(t-2)+0.0376e_(t-3) > > First, it is transformed: > a_t=0.9391a_(t-1)-2.8763a_(t-2)+2.9372a_(t-3)+0.0376e_t+0.1127e_(t-1) > +0.1127e_(t-2)+0.0376e_(t-3) > > According to the source paper, e_t~N(0, 0.01^2), I find that the output > a_t has extremely huge magnitude signal after 1000 or so noise samples. > > The delay units initial values are zero, but I think it is not a big > problem to this IIR yet for 0 value now. > > Could you tell me what is wrong? > > ThanksAside from the fact that the specified filter is grossly unstable, what do you think the problem might be? A check that's so quick that I'll do it for free is to look at the coefficient of the state with the highest lag (-2.9372, in your case). This coefficient is the product of all the poles. If its absolute value is greater than 1, then your filter absolutely positively must be unstable (feel free to work that out as a formal proof -- it's not hard). The more its absolute value exceeds 1, the worse case you're in. -- Tim Wescott Wescott Design Services http://www.wescottdesign.com

Reply by ●June 24, 20142014-06-24

On Tue, 24 Jun 2014 18:40:03 -0500, Tim Wescott wrote:> On Tue, 24 Jun 2014 16:19:32 -0700, fl wrote: > >> Hi, >> I simulate an ARMA (3,3) in Matlab/Simulink: >> >> a_t-0.9391a_(t-1)+2.8763a_(t-2)-2.9372a_(t-3)=0.0376e_t+0.1127e_(t-1) >> +0.1127e_(t-2)+0.0376e_(t-3) >> >> First, it is transformed: >> a_t=0.9391a_(t-1)-2.8763a_(t-2)+2.9372a_(t-3)+0.0376e_t+0.1127e_(t-1) >> +0.1127e_(t-2)+0.0376e_(t-3) >> >> According to the source paper, e_t~N(0, 0.01^2), I find that the output >> a_t has extremely huge magnitude signal after 1000 or so noise samples. >> >> The delay units initial values are zero, but I think it is not a big >> problem to this IIR yet for 0 value now. >> >> Could you tell me what is wrong? >> >> Thanks > > Aside from the fact that the specified filter is grossly unstable, what > do you think the problem might be? > > A check that's so quick that I'll do it for free is to look at the > coefficient of the state with the highest lag (-2.9372, in your case). > This coefficient is the product of all the poles. If its absolute value > is greater than 1, then your filter absolutely positively must be > unstable (feel free to work that out as a formal proof -- it's not > hard). The more its absolute value exceeds 1, the worse case you're in.When you figure out what the problem is with your numbers, we'll explain to you why you should always break your IIR filters down into 1 and 2 pole sections in cascade. -- Tim Wescott Wescott Design Services http://www.wescottdesign.com

Reply by ●June 25, 20142014-06-25

On Tuesday, June 24, 2014 7:19:32 PM UTC-4, fl wrote:> Hi, > > I simulate an ARMA (3,3) in Matlab/Simulink: > > > > a_t-0.9391a_(t-1)+2.8763a_(t-2)-2.9372a_(t-3)=0.0376e_t+0.1127e_(t-1) > > +0.1127e_(t-2)+0.0376e_(t-3) > > > > First, it is transformed: > > a_t=0.9391a_(t-1)-2.8763a_(t-2)+2.9372a_(t-3)+0.0376e_t+0.1127e_(t-1) > > +0.1127e_(t-2)+0.0376e_(t-3) > > > > According to the source paper, e_t~N(0, 0.01^2), I find that the output a_t > > has extremely huge magnitude signal after 1000 or so noise samples. > > > > The delay units initial values are zero, but I think it is not a big problem to this IIR yet for 0 value now. > > > > Could you tell me what is wrong? > > > > ThanksThanks Tim. The ARMA filter is from a paper: J. R. Statist. Soc. B (2000) 62, Part 3, pp. 493-508. Mixture Kalman filters I check the coefficients which are the same with the paper. It said that it emulates a fast flat fading channel. Have you noticed that the small var of N(0, 0.01^2)? I have no further idea to confirm the coef now.

Reply by ●June 25, 20142014-06-25

On Tuesday, June 24, 2014 7:19:32 PM UTC-4, fl wrote:> Hi, > > I simulate an ARMA (3,3) in Matlab/Simulink: > > > > a_t-0.9391a_(t-1)+2.8763a_(t-2)-2.9372a_(t-3)=0.0376e_t+0.1127e_(t-1) > > +0.1127e_(t-2)+0.0376e_(t-3) > > > > First, it is transformed: > > a_t=0.9391a_(t-1)-2.8763a_(t-2)+2.9372a_(t-3)+0.0376e_t+0.1127e_(t-1) > > +0.1127e_(t-2)+0.0376e_(t-3) > > > > According to the source paper, e_t~N(0, 0.01^2), I find that the output a_t > > has extremely huge magnitude signal after 1000 or so noise samples. > > > > The delay units initial values are zero, but I think it is not a big problem to this IIR yet for 0 value now. > > > > Could you tell me what is wrong? > > > > ThanksThanks Tim. The ARMA filter is from a paper: J. R. Statist. Soc. B (2000) 62, Part 3, pp. 493-508. Mixture Kalman filters I check the coefficients which are the same with the paper. It said that it emulates a fast flat fading channel. Have you noticed that the small var of N(0, 0.01^2)? I have no further idea to confirm the coef now.

Reply by ●June 25, 20142014-06-25

On Tuesday, June 24, 2014 7:19:32 PM UTC-4, fl wrote:> Hi, > > I simulate an ARMA (3,3) in Matlab/Simulink: > > > > a_t-0.9391a_(t-1)+2.8763a_(t-2)-2.9372a_(t-3)=0.0376e_t+0.1127e_(t-1) > > +0.1127e_(t-2)+0.0376e_(t-3) > > > > First, it is transformed: > > a_t=0.9391a_(t-1)-2.8763a_(t-2)+2.9372a_(t-3)+0.0376e_t+0.1127e_(t-1) > > +0.1127e_(t-2)+0.0376e_(t-3) > > > > According to the source paper, e_t~N(0, 0.01^2), I find that the output a_t > > has extremely huge magnitude signal after 1000 or so noise samples. > > > > The delay units initial values are zero, but I think it is not a big problem to this IIR yet for 0 value now. > > > > Could you tell me what is wrong? > > > > ThanksThanks Tim. The ARMA filter is from a paper: J. R. Statist. Soc. B (2000) 62, Part 3, pp. 493-508. Mixture Kalman filters I check the coefficients which are the same with the paper. It said that it emulates a fast flat fading channel. Have you noticed that the small var of N(0, 0.01^2)? I have no further idea to confirm the coef now.

Reply by ●June 25, 20142014-06-25

On Tuesday, June 24, 2014 7:19:32 PM UTC-4, fl wrote:

Reply by ●June 25, 20142014-06-25

On Wed, 25 Jun 2014 10:21:18 -0700, fl wrote:> On Tuesday, June 24, 2014 7:19:32 PM UTC-4, fl wrote: >> Hi, >> >> I simulate an ARMA (3,3) in Matlab/Simulink: >> >> >> >> a_t-0.9391a_(t-1)+2.8763a_(t-2)-2.9372a_(t-3)=0.0376e_t+0.1127e_(t-1) >> >> +0.1127e_(t-2)+0.0376e_(t-3) >> >> >> >> First, it is transformed: >> >> a_t=0.9391a_(t-1)-2.8763a_(t-2)+2.9372a_(t-3)+0.0376e_t+0.1127e_(t-1) >> >> +0.1127e_(t-2)+0.0376e_(t-3) >> >> >> >> According to the source paper, e_t~N(0, 0.01^2), I find that the output >> a_t >> >> has extremely huge magnitude signal after 1000 or so noise samples. >> >> >> >> The delay units initial values are zero, but I think it is not a big >> problem to this IIR yet for 0 value now. >> >> >> >> Could you tell me what is wrong? >> >> >> >> Thanks > > On the description of the ARMA model. It said that "In the > communications literature, this is called a (low pass) Butterworth > filter of order 3 with cut-off frequency 0.01. It is normalized to have > a stationary variance 1." > > I create an IIR, order 3, cut-off frequency 0.01 in Matlab. It works > well. I have their poles/zeros. How to get difference equations from > polse/zeros? > From IIR order, cut-off frequency, the difference equation is unique?Any good text on basic DSP practices will tell you how to do that. This paper touches the subject, but not in the usual approach: http:// wescottdesign.com/articles/zTransform/z-transforms.html. Search on "Z Transforms and Software" to see how to go from a transfer function to a difference equation and thence code. In general the procedure is to go from poles and zeros to a collection of transfer functions of no more than 2nd order (and 1st-order if you can pair up a real pole and a real zero). Going from a transfer function to a difference equation is NOT a 1:1 mapping. There are an infinite number of ways to do it, half a dozen popular forms, and probably a dozen or so that you'll stumble across in the literature if you're well-read in DSP and related fields. Each approach has its advantages and disadvantages. Wikipedia, unfortunately, lets us down on this one, but this page may give you enough search terms to get a start: http://en.wikipedia.org/wiki/ Digital_biquad_filter. Or, if you have access to a technical library, find the section where the DSP books live, pick one up and look up "direct form" in the index, then look at the author's discussion. Repeat until you find a book you like, check it out, and read it. -- Tim Wescott Wescott Design Services http://www.wescottdesign.com