Fit Fourier Models
About Fourier Series Models
A Fourier series describes a periodic function as a sum of sine and cosine functions. You can separate an arbitrary periodic function into simple components by using a Fourier series. These components are easy to integrate, differentiate, and analyze. For this reason, Fourier series are often used to approximate periodic signals.
Fourier series are represented in several forms. Curve Fitting Toolbox™ uses the trigonometric Fourier series form
where a0 models a constant (intercept) term in the data and is associated with the i = 0 cosine term, w is the fundamental frequency of the signal, and n is the number of terms (harmonics). Curve Fitting Toolbox supports Fourier series regression for 1 ≤ n ≤ 8.
For more information about Fourier series, refer to Fourier Analysis and Filtering.
Fit Fourier Model Interactively in Curve Fitter App
This example shows how to use the Curve Fitter app to fit a Fourier model to data.
Load the sound signal sample data.
load gong.mat
The variables y
and Fs
contain sound signal and frequency data for a gong ring, respectively. Create a sound clip by storing the first 1000 elements of y
in a vector named gongClip
.
gongClip = y(1:1000);
To calculate the time corresponding to each element in gongClip
, divide the index of the elements by Fs
.
t = [1:1000]./Fs;
Open the Curve Fitter app from the command line.
curveFitter
Alternatively, on the Apps tab, in the Math, Statistics and Optimization group, click Curve Fitter.
In the Curve Fitter app, select the data variables for the fit. On the Curve Fitter tab, in the Data section, click Select Data. In the Select Fitting Data dialog box, select t
as the X data value and gongClip
as the Y data value.
The app plots the data points as you select variables. By default, the app fits a polynomial to the data. To fit a Fourier model, click Fourier
in the Fit Type section of the Curve Fitter tab.
The app fits a single-term Fourier model.
The fitted one-term Fourier model is a periodic function with a simple oscillatory behavior. The Results panel displays the general equation for the model, fitted coefficient estimates with 95% confidence intervals, fundamental frequency, and goodness-of-fit statistics.
The fitted one-term Fourier model has a root mean square error (RMSE) of 0.1996. To compare the one-term Fourier model with a Fourier model that has four terms, select 4
for Number of terms in the Fit Options panel. The app fits a Fourier model with four terms to the data.
The fitted four-term Fourier model has more complex oscillatory behavior than the one-term Fourier model. An RMSE of 0.1685 for the four-term model indicates that four terms predict the sound data more accurately than one. However, the plot shows that some of the data points in gongClip
are outside of the range of the four-term model.
Export the fitted four-term Fourier model to the workspace by clicking Export in the Export section and then selecting Export to Workspace
. In the dialog box, uncheck the second and third options. Store the fit in the variable name in the box next to the first option.
You can listen to the sound data in gongClip
by using the function sound
.
sound(gongClip,Fs)
pause(2) % Allow gongClip to play before executing next line
To get the sound data for the Fourier model approximation of gongClip
, use feval
to evaluate gongFourierModel
at the times in t
. Play the approximated sound data.
gongClipApprox = feval(gongFourierModel,t); sound(gongClipApprox,Fs)
The two clips have the same approximate average tone. However, the approximated sound data does not have as many fluctuations in tone as the sound data in gongClip
.
Fit Fourier Model at the Command Line
This example shows how to fit a Fourier model to data using the fit
function.
Fit Two-Term Fourier Model
Load the El Niño-Southern Oscillation (ENSO) data.
load enso;
The variable pressure
contains data for the averaged atmospheric pressure difference between Easter Island, Chile and Darwin, Australia. The variable month
contains data for the month in which each pressure difference occurred.
Plot pressure
against month
.
plot(month,pressure)
The pressure data oscillates between 0 and 18, which indicates that it can be described by a Fourier series.
Fit a two-term Fourier model by using the Fourier library model. Specify the model type as fourier
followed by the number of terms. Save the goodness-of-fit statistics for later comparison.
[f2,gof2] = fit(month,pressure,"fourier2")
f2 = General model Fourier2: f2(x) = a0 + a1*cos(x*w) + b1*sin(x*w) + a2*cos(2*x*w) + b2*sin(2*x*w) Coefficients (with 95% confidence bounds): a0 = 10.63 (10.23, 11.03) a1 = 2.923 (2.27, 3.576) b1 = 1.059 (0.01593, 2.101) a2 = -0.5052 (-1.086, 0.07532) b2 = 0.2187 (-0.4202, 0.8576) w = 0.5258 (0.5222, 0.5294)
gof2 = struct with fields:
sse: 1.1230e+03
rsquare: 0.4279
dfe: 162
adjrsquare: 0.4103
rmse: 2.6329
f2
is a cfit
object containing the general formula, coefficient estimates with 95% confidence bounds, and fundamental frequency for the fit w
. The confidence bounds on a2
and b2
cross zero, so not enough evidence exists to conclude that they differ from zero or that the fitted model differs from a one-term Fourier model. The root mean square error (RMSE) of 2.6329 is useful for comparing the accuracy of f2
to the accuracy of other fits.
To calculate the period from the fundamental frequency, use the formula T = 2*pi/w
.
w = f2.w
w = 0.5258
T = 2*pi/w
T = 11.9497
The period of the fitted two-term Fourier model is approximately 12 months, or one year.
Plot f2
with a scatter plot of the data.
plot(f2,month,pressure)
The shape of f2
is similar to the shape of a one-term Fourier model, and the oscillation peaks approximately once every 12 months.
Fit Seven-Term Fourier Model
Fit a seven-term Fourier model to the data. Save the goodness-of-fit statistics.
[f7,gof7] = fit(month,pressure,"fourier7")
f7 = General model Fourier7: f7(x) = a0 + a1*cos(x*w) + b1*sin(x*w) + a2*cos(2*x*w) + b2*sin(2*x*w) + a3*cos(3*x*w) + b3*sin(3*x*w) + a4*cos(4*x*w) + b4*sin(4*x*w) + a5*cos(5*x*w) + b5*sin(5*x*w) + a6*cos(6*x*w) + b6*sin(6*x*w) + a7*cos(7*x*w) + b7*sin(7*x*w) Coefficients (with 95% confidence bounds): a0 = 10.63 (10.28, 10.97) a1 = 0.5669 (0.08285, 1.051) b1 = 0.1969 (-0.29, 0.6838) a2 = -1.203 (-1.687, -0.7189) b2 = -0.8085 (-1.307, -0.31) a3 = 0.9323 (0.4325, 1.432) b3 = 0.7599 (0.2622, 1.258) a4 = -0.6653 (-1.149, -0.1817) b4 = -0.2038 (-0.6995, 0.292) a5 = -0.02913 (-0.5129, 0.4547) b5 = -0.3701 (-0.8566, 0.1164) a6 = -0.04841 (-0.5437, 0.4469) b6 = -0.1367 (-0.6286, 0.3552) a7 = 2.812 (2.19, 3.433) b7 = 1.333 (0.4017, 2.264) w = 0.07527 (0.07478, 0.07576)
gof7 = struct with fields:
sse: 768.3656
rsquare: 0.6086
dfe: 152
adjrsquare: 0.5700
rmse: 2.2483
f7
contains several coefficients with confidence bounds that cross zero, so not enough evidence exists to conclude that their corresponding terms increase the accuracy of the fitted Fourier model. The RMSE of 2.2483 is smaller than the RMSE error of f2
, which confirms that the seven-term Fourier model predicts the pressure more accurately than the two-term Fourier model.
To calculate the period from the fundamental frequency, use the formula T = 2*pi/w
to calculate the period.
w = f7.w
w = 0.0753
T = (2*pi)/w
T = 83.4745
The period of the fitted seven-term Fourier model is approximately 83 months, or roughly seven years. The amplitude of the fitted coefficients determines which terms contribute most to the predicted value of the pressure difference.
The period of a sinusoid of the form sin(Ax)
or cos(Ax)
is given by the formula T = 2*pi/|A|
. a7
and b7
are the largest coefficients.
T = 2*pi/(w*7)
T = 11.9249
The period of the terms corresponding to a7
and b7
is approximately 12 months, indicating that the annual cycle is the strongest.
Use the same formula to calculate the periods of the following terms:
The terms
a1
andb1
have a period of 7 years each.The terms
a2
andb2
have a period of 3.5 (7/2) years each. Thea2
andb2
coefficients have larger magnitude than a1 and b1, so the 3.5-year cycle contributes more to the predicted value of the pressure difference than the 7-year cycle.The terms
a3
andb3
are strong, indicating 2.3-year (7/3) cycle.
Smaller terms such as a6
, b6
, a5
, and b5
are less important for the fit.
Plot f7
with a scatter plot of the data.
plot(f7,month,pressure)
The seven-term Fourier model oscillates in a more complex pattern and captures a wider range of values in the pressure difference than the one-term Fourier model. The cycle repeats approximately every 84 months, or 7 years. Typically, the El Niño warming happens at irregular intervals of two to seven years, and lasts nine months to two years. The average period length is five years. The model results reflect some of these periods.
Set Start Points
The fit
function uses the data
input argument to calculate optimized start points for the coefficient and fundamental frequency calculations. Fourier series models are particularly sensitive to start points, and the optimized values might be accurate for only a few terms in the associated equations. You can override the optimized start points by specifying the StartPoint
name-value argument.
The extreme values in the scatter plot of the data suggest that a four-year cycle might be present. To confirm this suggestion, set the start point of the fundamental frequency to the value corresponding to a period of eight years, or 96 months. An eight-year period for the fitted Fourier model increases the period of the terms a2
and b2
from 3.5
to 4
.
w_8 = (2*pi)/96
w_8 = 0.0654
Find the index of the fundamental frequency in the cell vector of f7
coefficient names by using the coeffnames
function.
coeffnames(f7)
ans = 16×1 cell
{'a0'}
{'a1'}
{'b1'}
{'a2'}
{'b2'}
{'a3'}
{'b3'}
{'a4'}
{'b4'}
{'a5'}
{'b5'}
{'a6'}
{'b6'}
{'a7'}
{'b7'}
{'w' }
The fundamental frequency is in the last entry of the vector of coefficient names. Create a vector of coefficient values from the coefficients of f7
, and replace the value for the fundamental frequency with the value corresponding to an eight-year period.
coeffs = coeffvalues(f7); coeffs(:,end) = w_8
coeffs = 1×16
10.6262 0.5669 0.1969 -1.2031 -0.8085 0.9323 0.7599 -0.6653 -0.2038 -0.0291 -0.3701 -0.0484 -0.1367 2.8120 1.3330 0.0654
Fit a seven-term Fourier model to the pressure difference data using the coefficients with the new value for the fundamental frequency as the start point. Save the goodness-of-fit statistics.
[f7_8,gof7_8] = fit(month,pressure,"fourier7",StartPoint=coeffs)
f7_8 = General model Fourier7: f7_8(x) = a0 + a1*cos(x*w) + b1*sin(x*w) + a2*cos(2*x*w) + b2*sin(2*x*w) + a3*cos(3*x*w) + b3*sin(3*x*w) + a4*cos(4*x*w) + b4*sin(4*x*w) + a5*cos(5*x*w) + b5*sin(5*x*w) + a6*cos(6*x*w) + b6*sin(6*x*w) + a7*cos(7*x*w) + b7*sin(7*x*w) Coefficients (with 95% confidence bounds): a0 = 10.58 (10.05, 11.1) a1 = 0.3286 (-0.4339, 1.091) b1 = -0.05917 (-0.7884, 0.6701) a2 = -0.8667 (-1.738, 0.004258) b2 = 1.094 (0.2819, 1.906) a3 = -0.4524 (-1.232, 0.3272) b3 = -0.3117 (-1.099, 0.4753) a4 = 0.181 (-0.7949, 1.157) b4 = 0.5806 (-0.1796, 1.341) a5 = 0.03263 (-0.7174, 0.7827) b5 = -0.2299 (-0.9767, 0.5169) a6 = 0.3726 (-0.39, 1.135) b6 = -0.2745 (-1.165, 0.6161) a7 = 0.4309 (-0.491, 1.353) b7 = -0.3547 (-1.316, 0.6062) w = 0.06795 (0.06519, 0.0707)
gof7_8 = struct with fields:
sse: 1.6851e+03
rsquare: 0.1416
dfe: 152
adjrsquare: 0.0568
rmse: 3.3296
The coefficients of f7_8
are slightly shifted from the f7
coefficients. The higher RMSE for f7_8
indicates that f7
is a better fit for the data. Plot both fits to visually compare the models.
plot(f7_8,month,pressure) hold on plot(f7, 'b') hold off legend("Data","f7_8","f7")
The plot shows that f7
captures the variation in the pressure difference data more accurately than f7_8
.
Display Fourier Fit Iterations
An alternative to specifying additional options using name-value arguments is to pass a fitoptions
object to the fit
function. To view the available options for a Fourier model fit, pass the model name as an input argument to the fitoptions
function.
fitoptions("fourier7")
ans = Normalize: 'off' Exclude: [] Weights: [] Method: 'NonlinearLeastSquares' Robust: 'Off' StartPoint: [1×0 double] Lower: [1×0 double] Upper: [1×0 double] Algorithm: 'Trust-Region' DiffMinChange: 1.0000e-08 DiffMaxChange: 0.1000 Display: 'Notify' MaxFunEvals: 600 MaxIter: 400 TolFun: 1.0000e-06 TolX: 1.0000e-06
Create a fitoptions
object, and specify to display the output after each iteration.
optionsf7 = fitoptions("fourier7",Display="iter")
options = Normalize: 'off' Exclude: [] Weights: [] Method: 'NonlinearLeastSquares' Robust: 'Off' StartPoint: [1×0 double] Lower: [1×0 double] Upper: [1×0 double] Algorithm: 'Trust-Region' DiffMinChange: 1.0000e-08 DiffMaxChange: 0.1000 Display: 'Iter' MaxFunEvals: 600 MaxIter: 400 TolFun: 1.0000e-06 TolX: 1.0000e-06
optionsf7
is a fitoptions
object that contains the options for a seven-term Fourier model fit.
To view the iteration steps involved in creating f7
, fit another seven-term Fourier model using the options in optionsf7
.
f7_iter = fit(month,pressure,"fourier7",optionsf7)
Norm of First-order Iteration Func-count f(x) step optimality CG-iterations 0 2 768.41 1.93e+03 1 4 768.366 2.2176e-05 69.1 0 2 6 768.366 7.94962e-07 2.48 0 Success, but fitting stopped because change in residuals less than tolerance (TolFun).
f7_iter = General model Fourier7: f7_iter(x) = a0 + a1*cos(x*w) + b1*sin(x*w) + a2*cos(2*x*w) + b2*sin(2*x*w) + a3*cos(3*x*w) + b3*sin(3*x*w) + a4*cos(4*x*w) + b4*sin(4*x*w) + a5*cos(5*x*w) + b5*sin(5*x*w) + a6*cos(6*x*w) + b6*sin(6*x*w) + a7*cos(7*x*w) + b7*sin(7*x*w) Coefficients (with 95% confidence bounds): a0 = 10.63 (10.28, 10.97) a1 = 0.5669 (0.08285, 1.051) b1 = 0.1969 (-0.29, 0.6838) a2 = -1.203 (-1.687, -0.7189) b2 = -0.8085 (-1.307, -0.31) a3 = 0.9323 (0.4325, 1.432) b3 = 0.7599 (0.2622, 1.258) a4 = -0.6653 (-1.149, -0.1817) b4 = -0.2038 (-0.6995, 0.292) a5 = -0.02913 (-0.5129, 0.4547) b5 = -0.3701 (-0.8566, 0.1164) a6 = -0.04841 (-0.5437, 0.4469) b6 = -0.1367 (-0.6286, 0.3552) a7 = 2.812 (2.19, 3.433) b7 = 1.333 (0.4017, 2.264) w = 0.07527 (0.07478, 0.07576)
To investigate the Fourier model fit further, you can experiment with specifying different options available for the NonlinearLeastSquares
fitting algorithm. See fitoptions
for more information.
See Also
Apps
Functions
fit
|fittype
|fitoptions