Re: moving average function
- To: mathgroup at smc.vnet.net
- Subject: [mg126328] Re: moving average function
- From: "Oleksandr Rasputinov" <oleksandr_rasputinov at ymail.com>
- Date: Tue, 1 May 2012 05:23:49 -0400 (EDT)
- Delivered-to: l-mathgroup@mail-archive0.wolfram.com
- References: <jnlj94$n3k$1@smc.vnet.net>
On Mon, 30 Apr 2012 09:42:12 +0100, Robert McHugh <rtmphone09 at gmail.com>
wrote:
> Below is a moving average function that has the following features:
> 1. returns a list with the same length as the original length
> 2. provides a reasonable estimate for averages on the "sides" of the
> window.
>
> Have failed to figure out how to do this with ListConvolve and
> ListCorrelate, so I submit this with the hope that others can
> recommend how it might be improved. Also searched this website for
> alternatives but didn't find any that met the above criteria.
>
> I was motivated to do this in order to keep my code free of handling
> special cases related to the edges of the widow size. Note that in
> one particular case, I have data measured every minute and would like
> to compare the results of using averaging the data over window sizes
> of 61, 121, and 181.
>
> Recommendations for how to improve the function or alternatives are
> appreciated.
> Bob.
>
>
> movingAverageBalanced[list_List, nAvg_Integer?OddQ ] :=
> Module[{nHang, middle, left, right, all},
> nHang = (nAvg - 1)/2;
>
> middle = MovingAverage[list, nAvg];
>
> left = Total[ Take[list, 2 # - 1]] /(2 # - 1) & /@ Range[nHang];
> right =
> Reverse[Total[ Take[list, -( 2 # - 1)]] /(2 # - 1) & /@
> Range[nHang]];
> all = Join[left, middle, right] ;
> Return[all];
> ]
>
> Example
> listTest = {a, b, c, d, e, f, g, h, i, j, k};
> r = movingAverageBalanced[listTest, 5];
> r // TableForm
>
> {
> {a},
> {1/3 (a + b + c)},
> {1/5 (a + b + c + d + e)},
> {1/5 (b + c + d + e + f)},
> {1/5 (c + d + e + f + g)},
> {1/5 (d + e + f + g + h)},
> {1/5 (e + f + g + h + i)},
> {1/5 (f + g + h + i + j)},
> {1/5 (g + h + i + j + k)},
> {1/3 (i + j + k)},
> {k}
> }
>
This is almost a special case of the Savitzky-Golay filter, code for which
I posted here recently (see
http://forums.wolfram.com/mathgroup/archive/2012/Feb/msg00036.html). The
main difference is in the treatment of the endpoints; while you chose to
use a shortened window in dealing with these, I used a skewed window, i.e.
taking the endpoint values as given by the best-fitting n-th order
polynomial. The latter will obviously give worse results for the simple
moving average where n = 0 (i.e., the endpoints would be approximated by
their mean value) but can provide better results for moving polynomial
windows with larger n. For example, the effect of a 5-point quadratic
smoothing kernel is given by:
SavitzkyGolayFilter[
{a, b, c, d, e, f, g, h, i, j, k},
2, 2, 0
] // Factor
i.e.,
{1/35 (31 a + 9 b - 3 c - 5 d + 3 e),
1/35 (9 a + 13 b + 12 c + 6 d - 5 e),
1/35 (-3 a + 12 b + 17 c + 12 d - 3 e),
1/35 (-3 b + 12 c + 17 d + 12 e - 3 f),
1/35 (-3 c + 12 d + 17 e + 12 f - 3 g),
1/35 (-3 d + 12 e + 17 f + 12 g - 3 h),
1/35 (-3 e + 12 f + 17 g + 12 h - 3 i),
1/35 (-3 f + 12 g + 17 h + 12 i - 3 j),
1/35 (-3 g + 12 h + 17 i + 12 j - 3 k),
1/35 (-5 g + 6 h + 12 i + 13 j + 9 k),
1/35 (3 g - 5 h - 3 i + 9 j + 31 k)}
One may also take derivatives in this way; again using a 5-point kernel, a
piecewise linear approximation to the first derivative is:
SavitzkyGolayFilter[
{a, b, c, d, e, f, g, h, i, j, k},
2, 2, 1
] // Factor
or,
{1/70 (-34 a + 3 b + 20 c + 17 d - 6 e),
1/35 (-12 a - b + 5 c + 6 d + 2 e),
1/10 (-2 a - b + d + 2 e),
1/10 (-2 b - c + e + 2 f),
1/10 (-2 c - d + f + 2 g),
1/10 (-2 d - e + g + 2 h),
1/10 (-2 e - f + h + 2 i),
1/10 (-2 f - g + i + 2 j),
1/10 (-2 g - h + j + 2 k),
1/35 (-2 g - 6 h - 5 i + j + 12 k),
1/70 (6 g - 17 h - 20 i - 3 j + 34 k)}
and so on.