On Sat, Oct 9, 2010 at 10:36 PM, <[email protected]> wrote: > On Sat, Oct 9, 2010 at 10:01 PM, Charles R Harris > <[email protected]> wrote: >> >> >> On Sat, Oct 9, 2010 at 7:47 PM, <[email protected]> wrote: >>> >>> I'm trying to see whether I can do this without reading the full manual. >>> >>> Is it intended that fromroots normalizes the highest order term >>> instead of the lowest? >>> >>> >>> >>> import numpy.polynomial as poly >>> >>> >>> p = poly.Polynomial([1, -1.88494037, 0.0178126 ]) >>> >>> p >>> Polynomial([ 1. , -1.88494037, 0.0178126 ], [-1., 1.]) >>> >>> pr = p.roots() >>> >>> pr >>> array([ 0.53320748, 105.28741219]) >>> >>> poly.Polynomial.fromroots(pr) >>> Polynomial([ 56.14003571, -105.82061967, 1. ], [-1., 1.]) >>> >>> >>> >>> renormalizing >>> >>> >>> p2 = poly.Polynomial.fromroots(pr) >>> >>> p2/p2.coef[0] >>> Polynomial([ 1. , -1.88494037, 0.0178126 ], [-1., 1.]) >>> >>> >>> this is, I think what I want to do, invert roots that are >>> inside/outside the unit circle (whatever that means >>> >>> >>> pr[np.abs(pr)<1] = 1./pr[np.abs(pr)<1] >>> >>> p3 = poly.Polynomial.fromroots(pr) >>> >>> p3/p3.coef[0] >>> Polynomial([ 1. , -0.54270529, 0.0050643 ], [-1., 1.]) >>> >> >> Wrong function ;) You defined the polynomial by its coefficients. What you >> want to do is > > My coefficients are from a lag-polynomial in time series analysis > (ARMA), and they really are the (estimated) coefficients. It is > essentially the same as the model for scipy.signal.lfilter. > I just need to check the roots to see whether the process is > stationary and invertible. If one of the two lag-polynomials (moving > average) has roots on the wrong side of the unit circle, then I can > invert them. >
I have just been doing np.roots(np.r_[1,ma_params]) and the MA coefficients are invertible if np.abs(np.roots(np.r_[1,ma_params])) < 1 That is the solution to (L**q + L**(q-1)*macoef[0] + ... + macoef[q-1]) = 0 are inside the unit circle Gretl, for example, uses the inverse of this so that the MA coefficients are invertible if (1 + macoef[0]*L + macoef[1]*L**2 + ... + macoef[q-1]*L**q) are outside the unit circle, which uses the complex inverse 1/np.roots(np.r_[1,ma_params]) or np.roots(np.r_[1,ma_params][::-1]) with the roots in reverse order of the original coefficients. The AR coefficients are the same, except the parameters should be subtracted. Ie., np.roots(np.r_[1,-ar_params]) should be inside the unit circle for stationarity. Skipper > I'm coding from memory of how this is supposed to work, so maybe I'm > back to RTFM and RTFTB (TB=text book). > > (I think what I really would need is a z-transform, but I don't have > much of an idea how to do this on a computer) > > Thanks, the main thing I need to do is check the convention or > definition for the normalization. And as btw, I like that the coef are > in increasing order > e.g. seasonal differencing multiplied with 1 lag autoregressive > poly.Polynomial([1.,0,0,-1])*poly.Polynomial([1,0.8]) > > (I saw your next message: Last time I played with function > approximation, I didn't figure out what the basis does, but it worked > well without touching it) > > Josef > >> >> In [1]: import numpy.polynomial as poly >> >> In [2]: p = poly.Polynomial.fromroots([1, -1.88494037, 0.0178126 ]) >> >> In [3]: p >> Out[3]: Polynomial([ 0.03357569, -1.90070346, 0.86712777, 1. ], >> [-1., 1.]) >> >> In [4]: p.roots() >> Out[4]: array([-1.88494037, 0.0178126 , 1. ]) >> >> Chuck >> >> >> _______________________________________________ >> NumPy-Discussion mailing list >> [email protected] >> http://mail.scipy.org/mailman/listinfo/numpy-discussion >> >> > _______________________________________________ > NumPy-Discussion mailing list > [email protected] > http://mail.scipy.org/mailman/listinfo/numpy-discussion > _______________________________________________ NumPy-Discussion mailing list [email protected] http://mail.scipy.org/mailman/listinfo/numpy-discussion
