Hi Guys, I've been tryin' to figure something out and darnit, I just can't seem to find a satisfying answer. Maybe one of you guys will help me. You recall the idea of a moving averager built with a tapped delay-line FIR filter like the 8-point moving averager shown below in Figure 1. (All the "b" coefficients are equal to 1/8th.) Each new y(n) output is the average of the current x(n) input sample and the seven previous input samples. x(n)--------b0--->(+)----> y(n) | ^ v | [z^-1] | | | +----b1--->(+) | ^ v | [z^-1] | FIG. 1 | | +----b2--->(+) | ^ ... ... v | [z^-1] | | | +----b7--->(+) You'll also recall that a recursive 8-point moving averager can be built, with a "comb" filter (8 delays and a subtract) followed by a multiply and then a 1st-order integrator as shown below in Figure 2. + x(n)---------->(+)-->(*)-->(+)----+--> y(n) | -^ ^ ^ | v | | | v [z^-8] | 1/8 | [z^-1] | | | | +------>(+) +-----+ FIG. 2 The neat part of the Figure 2 network (a cascaded integrator-comb, CIC, filter) is that it reduces the computational workload of the averager in Figure 1. Other neat parts of Figure 2 are that its structure (except for a longer delay line) remains the same, and its comp. workload remains unchanged, should you want to compute a 100-point moving average! OK, you've heard all of this before. What occurred to me was that it might be useful to figure out how to build a network that computes an "N-point moving variance". That is, a network where each new y(n) output is the variance of the current x(n) input sample and the N-1 previous input samples. But darn, I can't come up with such an "N-point moving variance" network whose structure remains (essentially) the same regardless of the value of N. It's those darned subtractions, needed in the variance computation, that make this process "messy". (I suppose for a signal whose mean is zero, we could eliminate those pesky subtractions, so long as we make sure that N is large number. In that case, it seems to me, we could use the above Figure 2 network with the input being x(n)^2 and the coefficient set to N-1 for an "unbiased" variance estimate.) Have any of you guys built such a "moving variance" estimation network whose output is updated at the input sample rate? Is such a network lurking out there somewhere in the literature of DSP. (I spent time searching the Internet, but found nothing useful.) John Sampson (in April 2001) showed us a way to compute an N-point variance, but that scheme provided a new variance result only after N input samples rather than after each new input sample. Maybe the whole idea of computing a new "moving variance" value for each new input sample is silly. I'm not sure. Thanks & See Ya, [-Rick-]
Question about computing a "moving variance"
Started by ●February 27, 2005
Reply by ●February 27, 20052005-02-27
"Rick Lyons" <r.lyons@_BOGUS_ieee.org> wrote in message news:4221b997.481221828@news.sf.sbcglobal.net...> > Hi Guys, > > I've been tryin' to figure something out and > darnit, I just can't seem to find a satisfying > answer. Maybe one of you guys will help me.detailed stuff snipped. Hello Rick, Yes you can build a moving variance calculator. Basically you will need two of your delay lines where one is fed your data and the other is fed the squares of the data. And you will keep a sum of data and a sum of the squares of the data. What you are going to use is known of the "raw score formula" for variance. See the next link for the formula. It is not too hard to derive. http://www.mnstate.edu/wasson/ed602calcvarraws.htm I hope this answers your question. Clay
Reply by ●February 27, 20052005-02-27
Clay S. Turner wrote:> "Rick Lyons" <r.lyons@_BOGUS_ieee.org> wrote in message > news:4221b997.481221828@news.sf.sbcglobal.net... > >>Hi Guys, >> >> I've been tryin' to figure something out and >>darnit, I just can't seem to find a satisfying >>answer. Maybe one of you guys will help me. > > > detailed stuff snipped. > > Hello Rick, > > Yes you can build a moving variance calculator. Basically you will need two > of your delay lines where one is fed your data and the other is fed the > squares of the data. And you will keep a sum of data and a sum of the > squares of the data. What you are going to use is known of the "raw score > formula" for variance. See the next link for the formula. It is not too hard > to derive. > > http://www.mnstate.edu/wasson/ed602calcvarraws.htm > > I hope this answers your question. > > Clay >I would have suggested this, but Clay got there first. Just pay attention to your data path requirements -- you need to keep every bit of that squared term or your variance result will be incorrect. So for 12-bit unsigned input data you need 24 bits just for the singular squared term, plus whatever you need for the sum (signed is better, of course). -- Tim Wescott Wescott Design Services http://www.wescottdesign.com
Reply by ●February 27, 20052005-02-27
Rick Lyons wrote:> Hi Guys, > > I've been tryin' to figure something out and > darnit, I just can't seem to find a satisfying > answer. Maybe one of you guys will help me. > > You recall the idea of a moving averager built > with a tapped delay-line FIR filter like the > 8-point moving averager shown below in Figure 1. > (All the "b" coefficients are equal to 1/8th.) > Each new y(n) output is the average of the > current x(n) input sample and the seven > previous input samples. > > x(n)--------b0--->(+)----> y(n) > | ^ > v | > [z^-1] | > | | > +----b1--->(+) > | ^ > v | > [z^-1] | FIG. 1 > | | > +----b2--->(+) > | ^ > ... ... > v | > [z^-1] | > | | > +----b7--->(+) > > You'll also recall that a recursive 8-point > moving averager can be built, with a "comb" > filter (8 delays and a subtract) followed > by a multiply and then a 1st-order > integrator as shown below in Figure 2. > + > x(n)---------->(+)-->(*)-->(+)----+--> y(n) > | -^ ^ ^ | > v | | | v > [z^-8] | 1/8 | [z^-1] > | | | | > +------>(+) +-----+ > > FIG. 2 > > The neat part of the Figure 2 network (a cascaded > integrator-comb, CIC, filter) is that > it reduces the computational workload of the averager > in Figure 1. Other neat parts of Figure 2 are that > its structure (except for a longer delay line) > remains the same, and its comp. workload remains > unchanged, should you want to compute a 100-point > moving average! > > OK, you've heard all of this before. What occurred > to me was that it might be useful to figure out how > to build a network that computes an > "N-point moving variance". That is, a network where > each new y(n) output is the variance of the > current x(n) input sample and the N-1 previous input > samples. > > But darn, I can't come up with such an "N-point moving > variance" network whose structure remains > (essentially) the same regardless of the value of N. > It's those darned subtractions, needed in the > variance computation, that make this process "messy". > > (I suppose for a signal whose mean is zero, we could > eliminate those pesky subtractions, so long as we make > sure that N is large number. In that case, it seems to > me, we could use the above Figure 2 network with the > input being x(n)^2 and the coefficient set to > N-1 for an "unbiased" variance estimate.) > > Have any of you guys built such a "moving variance" > estimation network whose output is updated > at the input sample rate? > Is such a network lurking out there somewhere > in the literature of DSP. (I spent time searching > the Internet, but found nothing useful.) > > John Sampson (in April 2001) showed us a way > to compute an N-point variance, but that scheme > provided a new variance result only after N > input samples rather than after each new input sample. > > Maybe the whole idea of computing a new "moving variance" > value for each new input sample is silly. I'm not sure. > > Thanks & See Ya, > [-Rick-]Have you seen a paper by R.J. Hanson called "Stably Updating Mean and Standard Deviation of Data" in Communications of ACM # 18? It shows you how to do a one-pass running variance estimate. I have used it in situations where the data is blocked up and the algorithm is called repeatedly. By saving the state between calls you can get the variance of the cumulative data set each time, not just the current block. But it does not give you a variance on every single sample. Google Scholar might turn up a free copy of the paper; I only have hardcopy. John
Reply by ●February 27, 20052005-02-27
Tim Wescott wrote:> Clay S. Turner wrote: > >> "Rick Lyons" <r.lyons@_BOGUS_ieee.org> wrote in message >> news:4221b997.481221828@news.sf.sbcglobal.net... >> >>> Hi Guys, >>> >>> I've been tryin' to figure something out and >>> darnit, I just can't seem to find a satisfying >>> answer. Maybe one of you guys will help me. >> >> >> >> detailed stuff snipped. >> >> Hello Rick, >> >> Yes you can build a moving variance calculator. Basically you will >> need two of your delay lines where one is fed your data and the other >> is fed the squares of the data. And you will keep a sum of data and a >> sum of the squares of the data. What you are going to use is known of >> the "raw score formula" for variance. See the next link for the >> formula. It is not too hard to derive. >> >> http://www.mnstate.edu/wasson/ed602calcvarraws.htm >> >> I hope this answers your question. >> >> Clay >> > I would have suggested this, but Clay got there first. > > Just pay attention to your data path requirements -- you need to keep > every bit of that squared term or your variance result will be > incorrect. So for 12-bit unsigned input data you need 24 bits just for > the singular squared term, plus whatever you need for the sum (signed is > better, of course). >The "definition" in terms of deviations is usually described as a "two pass" algorithm so one sees the other "definition" in terms of difference of (non-central) moments as a "one pass" algorithm. The one pass algorithm is well know for its numerical instability compared to the two pass algorithm. There is a third way which is one pass and has good numerical stability. It is often called "the method of tentative means" which keeps running values for the count, the mean and a sum of squares about the mean. The basis of the algorithm is that the sum of squares is updated for BOTH the new data and the change in the tentative mean caused by the new data. (With enough arithmetic precision the "instability" can be ignored so this is more an issue of better results with less precision requirements.) Ask google for more details if you can not work it out yourself. The more important point is to know that the method exists.
Reply by ●February 27, 20052005-02-27
On Sun, 27 Feb 2005 11:23:16 -0500, "Clay S. Turner" <Physics@Bellsouth.net> wrote:> >"Rick Lyons" <r.lyons@_BOGUS_ieee.org> wrote in message >news:4221b997.481221828@news.sf.sbcglobal.net... >> >> Hi Guys, >> >> I've been tryin' to figure something out and >> darnit, I just can't seem to find a satisfying >> answer. Maybe one of you guys will help me. > >detailed stuff snipped. > >Hello Rick, > >Yes you can build a moving variance calculator. Basically you will need two >of your delay lines where one is fed your data and the other is fed the >squares of the data. And you will keep a sum of data and a sum of the >squares of the data. What you are going to use is known of the "raw score >formula" for variance. See the next link for the formula. It is not too hard >to derive. > >http://www.mnstate.edu/wasson/ed602calcvarraws.htm > >I hope this answers your question. > >ClayHi Clay, INTERESTING! Thanks for the pointer to the web page. I'll review that page carefully. [-Rick-]
Reply by ●February 27, 20052005-02-27
On Sun, 27 Feb 2005 09:18:21 -0800, Tim Wescott <tim@wescottnospamdesign.com> wrote: (snipped)>I would have suggested this, but Clay got there first. > >Just pay attention to your data path requirements -- you need to keep >every bit of that squared term or your variance result will be >incorrect. So for 12-bit unsigned input data you need 24 bits just for >the singular squared term, plus whatever you need for the sum (signed is >better, of course). > >Tim WescottHi Tim, good point. Thanks. [-Rick-]
Reply by ●February 27, 20052005-02-27
On 27 Feb 2005 09:46:04 -0800, "john" <johns@xetron.com> wrote: (snipped)> > >Have you seen a paper by R.J. Hanson called "Stably Updating Mean and >Standard Deviation of Data" in Communications of ACM # 18? It shows you >how to do a one-pass running variance estimate. I have used it in >situations where the data is blocked up and the algorithm is called >repeatedly. By saving the state between calls you can get the variance >of the cumulative data set each time, not just the current block. But >it does not give you a variance on every single sample. > >Google Scholar might turn up a free copy of the paper; I only have >hardcopy. > >JohnHi John, As a matter of fact, I have seen it. It was the original source of John Sampson's old post (that I mentioned). And yes, a copy of the paper is available on the web. I did a Google search and found it. Thanks John, [-Rick-]
Reply by ●February 27, 20052005-02-27
On Sun, 27 Feb 2005 18:30:54 GMT, Gordon Sande <g.sande@worldnet.att.net> wrote:>(snipped)>> >> Just pay attention to your data path requirements -- you need to keep >> every bit of that squared term or your variance result will be >> incorrect. So for 12-bit unsigned input data you need 24 bits just for >> the singular squared term, plus whatever you need for the sum (signed is >> better, of course). >> > >The "definition" in terms of deviations is usually described as a >"two pass" algorithm so one sees the other "definition" in terms of >difference of (non-central) moments as a "one pass" algorithm. The >one pass algorithm is well know for its numerical instability compared >to the two pass algorithm. There is a third way which is one pass and >has good numerical stability. It is often called "the method of >tentative means" which keeps running values for the count, the mean >and a sum of squares about the mean. The basis of the algorithm is that >the sum of squares is updated for BOTH the new data and the change in >the tentative mean caused by the new data. (With enough arithmetic >precision the "instability" can be ignored so this is more an issue >of better results with less precision requirements.) > >Ask google for more details if you can not work it out yourself. The >more important point is to know that the method exists.Hi Gordon, thanks for the new terminology (new to me) that I can use in my Googling. See Ya, [-Rick-]
Reply by ●February 27, 20052005-02-27
in article 4221b997.481221828@news.sf.sbcglobal.net, Rick Lyons at r.lyons@_BOGUS_ieee.org wrote on 02/27/2005 07:37:> In that case, it seems to > me, we could use the above Figure 2 network with the > input being x(n)^2 and the coefficient set to > N-1 for an "unbiased" variance estimate.i've done some moving "power" implementations before which was about the same as this. (i *did* remove DC with a DC block and kept a running sum of the last N values of (x[n])^2. i think what Clay said is the pretty straight-forward way to do what the statistics texts say (keep a running sum of x[n] and x[n]^2 using a CIC and compute mean and variance from that) and his reference to http://www.mnstate.edu/wasson/ed602calcvarraws.htm reminds me of a question i intended to ask the group a few times from long ago. using some of the language of the webpage with link above, there is the "population variance" and the "sample variance". and there is a slightly different formula for the variance of the population vs. the sample. i understand what the mathematics is about regarding this difference (to *some* point, but there is still a philosophical problem i have). let x[n] be some random variable with a known mean and known variance: mean: E{ x[n] } = integral{ x * p(x) dx} = m variance: E{ (x[n]-m)^2 } = integral{ (x-m)^2 * p(x) dx} = sigma^2 where the integral is over all x and p(x) is the p.d.f. of x[n]. now we take a "sample" of these x[n] (that's that statistician's use of the term, we DSPers would say "N samples of these x[n]") and, from that sample, compute a measure or estimate of the mean and variance in the straight way: (you need a mono-spaced font to view this right.) _ N-1 est. mean: x = (1/N) * SUM{ x[n] } n=0 N-1 _ est. variance: S^2 = (1/N) * SUM{ (x[n] - x)^2 } n=0 N-1 _ = (1/N) * SUM{ (x[n])^2 } - ( x )^2 n=0 this computation is valid for the so-called "population mean" and "population variance" if N represents the number of *all* of the samples of x[n]. in fact, it is not an estimate, but is actually computing the actual mean and variance of the entire set of x[n]. but if N is smaller than the entire population of x[n], then this is a "sample mean" and "sample variance" and there is some risk in getting the value wrong because of the choice of which x[n] you use and which x[n] are left out of the estimation. _ in that case, x and S^2 are both random variables and it turns out that the _ mean of x is the same as the mean of x[n] which is m (cool!). but it turns out that the mean of S^2 (the estimated variance) is a little smaller than sigma^2. _ E{ x } = E{ x[n] } = m but E{ S^2 } = (N-1)/N * E{ (x[n]-m)^2 } = (N-1)/N * sigma^2 so, it's a "biased estimator" and we fix that bias by multiplying E{ S^2 } by N/(N-1). now here's what bothers me: it's the discontinuity of the application of the bias correction. first of all the "biased" estimator works, even for N=1 (the mean is just the single sample x[0] and the variance is zero), but the "unbiased" estimator does not (division by 0). i guess that doesn't bother me too much since, if you had a random variable of non-zero sigma, you have to turn that estimated variance of zero to the non-zero sigma^2 and dividing by zero seems like a plausible way to do it. :-) but what really bothers me is the that if N is even only 1 less than the entire population, then you unbias the variance estimate by multiplying by N/(N-1), but if N is the entire population you don't multiply by N/(N-1). let "P" be the size of the entire population, if N = P-1 you multiply by (P-1)/(P-2) but when N = P, you don't multiply by (P-0)/(P-1). in general you always multiply S^2 by (P - (P-N))/(P - (P-N) - 1) except the case when P-N = 0. then you multiply by 1. can any of you statistics gurus explain why this discontinuity exists? why doesn't the unbiased estimate of variance work for the entire population? -- r b-j rbj@audioimagination.com "Imagination is more important than knowledge."