SavitzkyGolayFilter - Maple Help
For the best experience, we recommend viewing online help using Google Chrome or Microsoft Edge.

Online Help

All Products    Maple    MapleSim


SignalProcessing

  

SavitzkyGolayFilter

  

apply Savitzky-Golay filter to a signal

 

Calling Sequence

Parameters

Options

Description

Examples

References

Compatibility

Calling Sequence

SavitzkyGolayFilter( signal, degree, frameradius, options )

SavitzkyGolayFilter( signal, degree, weights, options )

Parameters

signal

-

1-D rtable or list of real-valued data

degree

-

nonnegative integer for the degree of the polynomials used for filtering

frameradius

-

positive integer for the radius of the frame used for filtering

weights

-

2-D rtable or list of positive real-valued weights for the least-squares fitting

Options

• 

extend: Keyword true or false. Specifies if the signal is to be extended on the left and right by frameradius in order to improve the quality of the filtered signal. The default is false.

• 

extrapolation: Keyword periodic, polynomial, or spline. Specifies how the signal is to be extended when extend=true. The default is polynomial.

• 

filteredplotoptions: List of additional plot options used to create the plot of the filtered signal. The default is [].

• 

plotoptions: List of additional plot options used to create the combined plot of the unfiltered and filtered signals. The default is [].

• 

recurrencevariables: List of three unique names for the index, old value, and new value used when constructing the recurrence relation for the filter. The default is [n,x,y].

• 

splinepoints: Positive integer other than 1 or the keyword infinity. Specifies the maximum number of data points used to perform spline extrapolation. The default is 10.

• 

timerange: Range of the form t1..t2, where t1 and t2 are numeric values. Specifies the time range for the plot and Vector of times. The default is 1..numelems(signal).

• 

unfilteredplotoptions: List of additional plot options used to create the plot of the unfiltered signal. The default is [].

• 

output: The type of output. Let m be the signal size, d be the degree, r be the frame radius, and l=2r+1 be the frame length. The supported options are:

– 

coefficients: The float[8] Vector, of size l, storing the coefficients of the filter.

– 

convolution: The float[8] Vector, of size m2r when extend=false and m when extend=true, storing the convolution of the filter.

– 

derivatives: The float[8] Matrix, of size m by d, storing the derivatives of the signal up to degree d.

– 

differentiator: The float[8] Matrix, of size l by d+1, storing the differentiator of the filter.

– 

filteredsignal: The float[8] Vector, of size m, storing the filtered signal.

– 

hankel: The float[8] Matrix, of size d+1 by d+1, storing the Hankel Matrix of the filter.

– 

nrr: The numeric value corresponding to the noise-reduction ratio of the filtered signal.

– 

plot: The plot of the unfiltered and filtered signals.

– 

projection: The float[8] Matrix, of size l by l, storing the projection of the filter.

– 

recurrence: The linear recurrence relation corresponding to the filter.

– 

times: The float[8] Vector with the time values determined by timerange.

– 

unfilteredsignal: The float[8] Vector, of size m, storing the unfiltered signal.

– 

vandermonde: The float[8] Matrix, of size l by d+1, storing the Vandermonde Matrix of the filter.

– 

weights: The float[8] Vector, of size l, storing the weights.

– 

record: Returns a record with the previous options.

– 

list of any of the above options: Returns an expression sequence with the corresponding outputs, in the same order. The default is filteredsignal.

Description

• 

The SavitzkyGolayFilter command applies the Savitzky-Golay Filter to a 1-D signal with real-valued data, which is useful for smoothing noisy data and estimating derivatives. The method is also known as the Least-Squares Smoothing, Polynomial Smoothing Filter, and Locally Weighted Scatterplot Smoothing (LOWESS).

• 

The filter works by applying polynomial fitting to a succession of frames/windows containing a fixed number of points. More specifically, for each interior data point of the signal, say xk which corresponds to time tk, the frame consists of the points xkr,..,xk+r, where r is the frame radius. The total number of points in the frame is the frame length, l=2r+1, and the index k satisfies r+1k and km+r, where m is the size of the signal. The filtered value yk of xk is determined by the value of the best-fit polynomial of degree d through the points in the frame at time tk.

• 

To determine the filtered signal, Maple does the following:

1. 

Suppose the frame consists of r nodes ni spaced uniformly between −1 and 1.

2. 

Let p be the polynomial of degree d having coefficients c with respect to the basis defined by the columns of a Vandermonde Matrix V, with dimensions l by d+1. The elements of V are such that the value of p at node ni is element i of y=V·c. The Vandermonde Matrix used by Maple has elements defined by Vi,j=ir1rj1.

3. 

It can be shown that the optimal choice of c is given by the the middle row of the Projection Matrix P=F·VTypesetting:-_Hold%T, where the Differentiator Matrix is given by F=V·Typesetting:-_Hold%^H,1 and the Hankel Matrix is given by H=VTypesetting:-_Hold%T·W·V. More specifically, if x is the signal over the frame, the objective function to minimize is given by φ=xV·cTypesetting:-_Hold%T·W·xV·c.

4. 

For each point of the unfiltered signal, say xk with r+1k and km+r, the filtered value yk is found using a frame centered at xk and computing the value of the best-fit polynomial at the center node using the above method. All the filtered values together between indices r+1 and m+r, namely the steady-state signal, are computed as the valid convolution of the unfiltered signal and the coefficients Vector.

5. 

To find the filtered signal on the left and right of the steady-state signal, namely the input-on transient and input-off transient, respectively, let u be the first r elements of the unfiltered signal x in reverse order, v be the last r elements of x in reverse order, Q be the first r rows of P in reverse order, and R be the last r rows of P in reverse order. The input-on transient is given by R·u and the input-off transient is given by Q·v.

• 

When extend=false, the original signal is extended on the left and right by r so that the overall signal has size m+2r. This is done so that the resulting steady-state of the filtered signal is of size m, which will be used as the filtered version of the original signal.

• 

The noise-reduction ratio is computed as the sum of the squares of the coefficients. Note that the sum of the coefficients is equal to 1.

• 

The Derivatives Matrix is found by convolving the unfiltered signal with the appropriately weighted columns of the Differentiator Matrix. When extend=false the original unfiltered signal is used with the convolution shape being same. On the other hand, when extend=true, the extended unfiltered signal is used with the convolution shape being valid.

• 

The frame length l and radius r satisfy l=2r+1. Moreover, the size of the signal m must satisfy lm, and the degree of the filtering polynomials d must satisfy d<2r.

• 

When a weights container is passed, its size corresponds to l and must be an odd integer larger than 1. Internally, the weights will be normalized with respect to root mean square.

• 

The signal and weights rtables cannot have an indexing function and must use Fortran ordering and rectangular storage.

• 

Maple attempts to coerce signal and weights to an rtable of float[8] datatype, and an error is thrown if this is not possible. For this reason, it is most efficient for the passed containers to be rtables of this datatype.

• 

The SavitzkyGolayFilter command is not thread safe.

Examples

withSignalProcessing&colon;

Example 1

XGenerateSignal&ExponentialE;0.5tsin2t&comma;t&equals;0..2Pi&comma;200&comma;&apos;noisedeviation&apos;&equals;0.05

d3

d3

(1)

r25

r25

(2)

SavitzkyGolayFilterX&comma;d&comma;r

SavitzkyGolayFilterX&comma;d&comma;r&comma;&apos;output&apos;&equals;&apos;plot&apos;

Example 2

• 

Consider the following signal:

fsint&plus;1cos3t5

fsint+cos3t5

(3)

t10

t10

(4)

t22Pi

t22π

(5)

m200

m200

(6)

T&comma;XGenerateSignalf&comma;t&equals;t1..t2&comma;m&comma;&apos;noisedeviation&apos;&equals;0.2&comma;&apos;output&apos;&equals;&apos;times&apos;&comma;&apos;signal&apos;&colon;

• 

Use of non-constant weights can improve the smoothing:

d2

d2

(7)

r10

r10

(8)

Wseq1..r&plus;1&comma;seqr..1&comma;1&colon;

SavitzkyGolayFilterX&comma;d&comma;r&comma;&apos;timerange&apos;&equals;t1..t2&comma;&apos;plotoptions&apos;&equals;&apos;title&apos;&equals;Unweighted Savitzky-Golay Filtered Signal&comma;&apos;output&apos;&equals;&apos;plot&apos;

SavitzkyGolayFilterX&comma;d&comma;W&comma;&apos;timerange&apos;&equals;t1..t2&comma;&apos;plotoptions&apos;&equals;&apos;title&apos;&equals;Weighted Savitzky-Golay Filtered Signal&comma;&apos;output&apos;&equals;&apos;plot&apos;

Example 3

• 

The Savitzky-Golay Filter can be used to determine the derivatives of a signal. For example, consider the following signal:

gt&rarr;2cost&plus;sin3t

gt2cost+sin3t

(9)

t12Pi

t12π

(10)

t22Pi

t22π

(11)

m500

m500

(12)

T&comma;XGenerateSignalg&comma;t1..t2&comma;m&comma;&apos;includefinishtime&apos;&equals;&apos;false&apos;&comma;&apos;output&apos;&equals;&apos;times&apos;&comma;&apos;signal&apos;&colon;

• 

Here, we will compare the un-extended and extended cases with the actual derivatives:

d3

d3

(13)

r20

r20

(14)

Y1SavitzkyGolayFilterX&comma;d&comma;r&comma;&apos;timerange&apos;&equals;t1..t2&comma;&apos;extend&apos;&equals;&apos;false&apos;&comma;&apos;output&apos;&equals;&apos;derivatives&apos;&colon;

Y2SavitzkyGolayFilterX&comma;d&comma;r&comma;&apos;timerange&apos;&equals;t1..t2&comma;&apos;extend&apos;&equals;&apos;true&apos;&comma;&apos;extrapolation&apos;&equals;&apos;periodic&apos;&comma;&apos;output&apos;&equals;&apos;derivatives&apos;&colon;

dgVectord&comma;k&rarr;Dkg&colon;

ZMatrixm&comma;d&comma;i&comma;j&rarr;dgjTi&comma;&apos;datatype&apos;&equals;&apos;float8&apos;&colon;

• 

As we can see, both versions are accurate but the extended one is more accurate at the ends of the time interval:

dataplotT&comma;Y1&comma;&apos;style&apos;&equals;&apos;line&apos;&comma;&apos;legend&apos;&equals;seqsprintf%s derivative&comma;converti&comma;&apos;english&apos;&comma;&apos;ordinal&apos;&comma;i&equals;1..d&comma;&apos;title&apos;&equals;Savitzky-Golay Derivatives Without Extension&comma;&apos;view&apos;&equals;&apos;DEFAULT&apos;&comma;minZ..maxZ

dataplotT&comma;Y2&comma;&apos;style&apos;&equals;&apos;line&apos;&comma;&apos;legend&apos;&equals;seqsprintf%s derivative&comma;converti&comma;&apos;english&apos;&comma;&apos;ordinal&apos;&comma;i&equals;1..d&comma;&apos;title&apos;&equals;Savitzky-Golay Derivatives With Extension&comma;&apos;view&apos;&equals;&apos;DEFAULT&apos;&comma;minZ..maxZ

• 

Using relative root mean square error, we can quantify the accuracy:

seqRelativeRootMeanSquareErrorY1..&comma;j&comma;Z..&comma;j&comma;j&equals;1..d

0.343537779779587105&comma;0.393736670852185722&comma;0.940391637995256491

(15)

seqRelativeRootMeanSquareErrorY2..&comma;j&comma;Z..&comma;j&comma;j&equals;1..d

0.0104018750855234855&comma;0.158908529185620623&comma;0.130136176920078145

(16)
• 

As expected, though, the unextended version is very accurate away from the endpoints:

seqRelativeRootMeanSquareErrorY1r&plus;1..r&plus;1&comma;j&comma;Zr&plus;1..r&plus;1&comma;j&comma;j&equals;1..d

0.0102865451746987236&comma;0.159079567157564300&comma;0.130110621919648534

(17)

Example 4

• 

The linear recurrence relation corresponding to the filter can be returned:

req1SavitzkyGolayFilterseq7&plus;1k&comma;k&equals;1..25&comma;2&comma;2&comma;&apos;recurrencevariables&apos;&equals;&apos;i&apos;&comma;&apos;u&apos;&comma;&apos;v&apos;&comma;&apos;output&apos;&equals;&apos;recurrence&apos;

req1vi=0.0857142857142857ui2+0.342857142857143ui1+0.485714285714286ui+0.342857142857143ui+10.0857142857142857ui+2

(18)
• 

In terms of fractions:

req2subsindetsreq1&comma;&apos;float&apos;&comma;identify

req2vi=3ui235+12ui135+17ui35+12ui+1353ui+235

(19)

Example 5

• 

And finally, an example included merely because it looks cool:

SavitzkyGolayFilterArray1..100&comma;k&rarr;1k&comma;2&comma;15&comma;&apos;output&apos;&equals;&apos;plot&apos;

References

  

Orfanidis, Sophocles J. Introduction to Signal Processing. Pearson Education, Inc. (2009).

Compatibility

• 

The SignalProcessing[SavitzkyGolayFilter] command was introduced in Maple 2023.

• 

For more information on Maple 2023 changes, see Updates in Maple 2023.

See Also

CurveFitting[Lowess]

SignalProcessing

SignalProcessing[Convolution]

SignalProcessing[DifferentiateData]

SignalProcessing[Welch]

Statistics[Lowess]