Transforms and Filters for Stochastic Processes (779449), страница 5
Текст из файла (страница 5)
.. ~ ( N - p + l )z(N - 1) ...z(N-p)X =7X =.(p)....(l)(5.182)132Chapter 5. Tkansforms and Filters for Stochastic ProcessesandWe haveand?Lf)(l) = X H x .(5.185)Covariance Method. The cowariance method takes into account the prediction errors in steady state only and yields an unbiased estimate of theautocorrelation matrix. In this case X and X are defined asX =[;i;]x ( N - 1) .. .
x ( N - p ):::X@)(5.186)and(5.187)The equation to be solved is(C")R,,a=h-+L;")(I),(5.188)=XHz.(5.190)wherer+")(l)xxNote that kg")is not a Toeplitz matrix, so that solving (5.188) is much morecomplex than solving (5.181) via the Levinson-Durbin recursion. However, thecovariance method has the advantage of being unbiased; we haveE{RE"'} = R x x .(5.191)1335.7.
Estimation of Autocorrelation Sequences and Power Spectral DensitiesEstimation of Autocorrelation Sequencesand Power Spectral Densities5.75.7.1Estimation of Autocorrelation SequencesIn the following, we will discuss methods for estimating the autocorrelationsequence of random processes from given sample values z(n),n = 0 , . . . , N - 1.We start the discussion with the estimatel?!,(m) = -NcN-In-l"x * ( n ) z(n+ m),(5.192)n=Owhichis the same astheestimate ?iic)(m),usedin theautocorrelationmethod explained in the last section. As can easily be verified, the estimate?:,(m) is biased with meanE {?!,(m)} =eNr,,(m).(5.193)However, sincelim E {?!,(m)} = r,,(m),N+CCthe estimate is asymptotically unbiased.
The triangular windowoccurs in (5.193) is known as the Bartlett window.The variance of the estimate can be approximated as [77]cl o ov 4 e x (m11lr,,(n)12+ rZx(n - m ) r z z ( n+ m ) .(5.194)vthat(5.195)n=-ccThus, as N+ m, the variance tends to zero:(5.196)Such an estimate is said to be consistent. However, although consistency isgiven, we cannot expect good estimates for large m as long as N is finite,because the bias increases as Iml + N .Unbiased Estimate. An unbiased estimate of the autocorrelation sequenceis given byN-In-l(5.197)134Chapter 5. Tkansforms and FiltersStochasticfor ProcessesThe variance of the estimate can be approximated as [77]var[+:, (m)]McccN(N - lmD2+l~,,(n)1~ rZ,(n - m) r,,(n + m ) .
(5.199)As N + CO,this giveslim var[?&(m)] + 0,N+cc(5.200)which means that ?,",(m) is a consistent estimate. However, problems arisefor large m as long as N is finite, because the variance increases for Iml + N.5.7.2Non-Parametric Estimation of PowerSpectralDensitiesIn many real-world problems one is interested in knowledge about the powerspectral density of the data to be processed. Typically, only a finite set ofobservations).(Xwith n = 0 , 1 , .
. . ,N-l is available. Since the power spectraldensity is the Fourier transform of the autocorrelation sequence, and sincewe have methods for the estimation of the autocorrelation sequence, it is alogical consequence to look at the Fourier transforms of these estimates. Westart with ?:,(m). The Fourier transform of +!,(m) will be denoted asWe know that ?!,(m) is a biased estimateof the true autocorrelationsequencer z z ( m ) ,so that we canconclude that the spectrumP,,(eJ") is a biasedestimate of the true power spectral density S,,(eJW).In order to be explicit,let us recall thatwith wB(m) being the Bartlett window; i.e.(5.203)5.7.
Estimation of Autocorrelation Sequences and Power Spectral Densities135In the spectral domain, we have(5.204)where W B ( e j w ) is the Fourier transform of w ~ ( mgiven) by(5.205)Thus, E {P,,(ej")} is a smoothed version of the true power spectral densityS,,(ej"), where smoothing is carried out with the Fourier transform of theBartlett window.A second way of computing P,,(eJ") is to compute the Fourier transformof ).(Xfirst and to derive P,,(ej") from X(ej"). By inserting (5.192) into(5.201) and rearranging the expression obtained, we get(5.206)In the form (5.206) P,,(eJ") is known as the periodogram.Another wayof deriving an estimate of the power spectral density is toconsider the Fourier transform of the estimate ?:,(m). We use the notationQ,,(eJW) for this type of estimate:cN-lQ z z ( e j w )=m=-(N-l)?:,(m)e-jwm.(5.207)136Chapter 5.
Tkansforms and FiltersStochasticfor ProcessesThe expected value iscN-lE{Q,,(ej")}=,!{+&(m)}e-jwmm=-(N-l)N-l(5.208)m=-(N-l)Mm=-mwhere w R ( m ) is the rectangular windoww R ( m )=andWR(ej")Iml 5 N{ 1, forotherwise,-1(5.209)0,is its Fourier transform:(5.210)This means that although +&(m)is an unbiased estimate of T,, ( m ) ,thequantity Q z z ( e j w ) is a biased estimate of S,,(ej"). The reason for this is thefact that only a finite number of taps of the autocorrelation sequenceis used inthe computation of Q s s ( e j w ) . The mean E { Q s s ( e j w ) } is a smoothed versionof S,,(eJW), where smoothingis carried out with the Fourier transform of therectangular window.As N + co both estimates ?:,(m) and ?&(m)become unbiased.
Thesame holds for P,,(ejw) and Q,,(ej"), so that both estimates of the powerspectral density are asymptotically unbiased. The behavior of the variance ofthe estimatesis different. While the estimatesof the autocorrelation sequencesare consistent, those of the power spectral density are not. For example, fora Gaussian process ).(Xwith power spectral density S Z z ( e J W ) the,varianceof the periodogram becomeswhich yieldslimvarN+CC[ P,, (ej") 1 = S,",(ej").Thus, the periodogram does not give a consistent estimate ofproof of (5.211) is straightforward and is omitted here.(5.212)S,,(ej").The1375.7. Estimation of Autocorrelation Sequences and Power Spectral DensitiesUse of the DFT or FFT forComputing the Periodogram.
Sincethe periodogram is computed from the Fourier transform of the finite datasequence, it can be efficiently evaluated at a discrete set of frequencies byusing the FFT.Given a length-N sequence ~ ( n we) , may consider a length-NDFT, resulting inP z z ( e J W k=)112N-lx ( n ) e-jZr'lNI)(5.213)with w k = 27rk/N. In many applications, the obtained number of samplesof P z z ( e j w )may be insufficient in order to draw a clear picture of the periodogram. Moreover, the DFT length may be inconvenient for computation,because no powerful FFT algorithm is at hand for the given length.
Theseproblemscanbe solvedby extending the sequence ~ ( nwith)zeros to anarbitrary length N' 2 N . This procedure is known as zero padding. We obtain(5.214)with w k = 27rk/N'. The evaluation of (5.214) is typically carried out via theFFT.Bartlett Method. Various methods have been proposed for achieving consistent estimates of the power spectral density. The Bartlett method does thisby decomposing the sequence ).(Xinto disjoint segments of smaller lengthand taking the ensemble average of the spectrum estimates derived from thesmaller segments.
With~ ( ~ )=( nz(n)+ iM), i = 0,1,. . . , K - 1,n = 0 , 1 , . . . , M - 1, (5.215)we get the K periodograms: : : 1C-(ejw) = MX ( i ) (n)e-jwn,i = 0,1,. . . , K- 1.(5.216)The Bartlett estimate then is(5.217)138Chapter 5. Tkansforms and FiltersStochasticfor Processeswith W B( e~J w )being the Fourier transform of the length-M Bartlettwindow.Assuming a Gaussian process z(n), the variance becomesvar[P,,(eBjw11 = -1v a r [ ~ , , ( e j ~ ) ]=Ksin (W M )K(5.219)Thus, as N , M , K + 00, the variance tends to zero and the estimate isconsistent. For finite N , the decomposition of ).(Xinto K sets results ina reduced variance, butthe bias increases accordingly andthespectrumresolution decreases.Blackman-Tukey Method.
Blackman and Tukey proposed windowing theestimated autocorrelation sequence prior the Fourier transform [8]. The argument is that windowing allows us to reduce the influence of the unreliableestimates of the autocorrelation sequence for large m. Denoting the windowand its Fourier transform as w(m)and W ( e j w ) respectively,,the estimate canbe written asN-lP,,BT (ej wC-w(m)?;,(m) e--jwm.(5.220)m=-(N-l)In the frequency domain, this means that(5.221)The windowW(.)should be chosen such thatW(ejw) >oVW(5.222)in order to ensure that PLT(ejw)is positive for all frequencies.The expected value of PLT(ejw)is most easily expressed in the formE {P,,BT ( ej w,>N-l=Cw ( m ) w B ( m )r,,(m) e-jwm.(5.223)m=-(N-l)Provided that w ( m ) is wide with respect to r z z ( m )and narrow with respectto W B ( ~ )the, expected value can be approximated asE { p , T ( e j w ) }= w(0) S,,(ejw).(5.224)Thus, in order to achieve an asymptotically unbiased estimate, the windowshould satisfyw ( 0 ) = 1.(5.225)5.7.
Estimation of Autocorrelation Sequences and Power Spectral Densities139For a symmetric window w ( m ) = W(-m) the variance can be estimated as [8]This approximation is based on the assumption that W(ej") iswide withrespect to W ~ ( e j " and) narrow with respect to the variations of S,,(ej").Welch Method. In the Welch method [l621 the data is divided into overlapping blocks+z(~)(,) = ~ ( ni D ) , i = 0 , 1 , . . . , K - 1,n = 0 , 1 , . . .
, M - 1 (5.227)with D 5 M . For D = M we approach the decomposition in the Bartlettmethod. For D < M we have more segments than in the Bartlett method.Each block is windowedprior to computation of the periodogram, resultingin K spectral estimatesThe factor a is chosen asc1 M-1a=w2(m) = - M m=OM 21r1"S&,(ejw) dw,(5.229)-Twhich means that the analysis is carried out with aenergy. Taking the average yields the final estimatewindow of normalized(5.230)The expected value becomeswith-M -Ml -lIn the spectral domain, this can be rewritten as140Chapter 5.
Tkansforms and FiltersStochasticfor Processes(5.234)With increasing N and M , SEW(ej("- ")) becomes narrow with respect toS z Z ( d v )and the expected value tends toThis shows that the Welch method is asymptotically unbiased.For a Gaussian process, the variance of the estimate is(5.236)If no overlap is considered (D = M ) , the expression reduces toWvar[P,,(ej w) I =-var[v$(ejw)11KM-1~ : ~ ( e j ~ ) .K(5.237)For k + 00 the variance approaches zero, which shows that the Welch methodis consistent.Various windows with different properties are known for the purpose ofspectral estimation.
In the following, a brief overview is given.Hanning Window.w(n) =0.5-0.5~0s0,n = 0 , 1 , ..., N - l(5.238)otherwise.Hamming Window.w(n) =0.54 - 0.46cos0,n = 0,1,. . . ,N - 1(5.239)otherwise.Blackman Window.w(n) =+0.08 cos,n = (),l,...,N-1otherwise.(5.240)5.7. Estimation of Autocorrelation Sequences and Power Spectral Densities14114Time IndexFigure 5.5. Window functions.Figure 5.5 shows the windows, and Figure 5.6shows their magnitudefrequency responses. The spectrum of the Bartlett window is positive for allfrequencies, which also means that the bias due to the Bartlett window isstrictly positive. The spectra of the HanningandHamming window haverelatively large negative side lobes, so that the estimated power spectraldensity may have a negative bias in the vicinity of large peaks in S,, (ej'").The Blackman window is a compromise between the Bartlett and the Hanning/Hamming approaches.5.7.3Parametric Methods in Spectral EstimationParametric methods in spectral estimation have been the subject of intensiveresearch, and many different methods have been proposed.