Description
We consider fitting experimental data to a known function defined by several parameters.
The particular scenario is the ubiquitous decaying sinusoidal function indicative of damped
harmonic motion:
W(t) = A0 + Ae−t/τ sin [(ω + αt)t + ϕ] (1)
The unknown parameters are A0, A, τ, ω, α, ϕ. The first term A0 is some fixed offset, A is the
zero-time sinusoidal amplitude, τ is a characteristic amplitude decay time, and the sinusoidal
term includes a phase offset ϕ and a linearly-chirped frequency ω + αt. If α = 0 we have a
simple damped harmonic oscillator equation.
As a specific example, consider the scenario in the figure which describes how to measure
the damped oscillation of an otherwise stationary fluid sphere. The spherical object is in the
path of a collimated (uniform intensity) laser beam which is focused on a photodiode. The
object scatters light out of the beam path (creates a shadow) so the intensity seen by the
photodiode is a function of the cross sectional area of the sphere. The photodiode outputs
a current proportional to the light intensity which can be determined by a voltage reading
across a large resistor. The spherical object is initially driven into sinusoidal volume or shape
oscillations which decay in amplitude due to natural viscous forces. The parameters τ , ω
and α are important system quantities that tell us about fluid properties (viscosity, density,
pressure, etc.) and interface properties (surface tension, impurities, etc.).
object
collimated
laser
beam lens
photodiode
I(t)
R V (t)
Typical data is as shown below, where the voltage is sampled at discrete times and contains
some measurement noise.
P1-1
The goal is to determine the parameters which provide the best fit to the discrete data
set {(tk, Vk)}
n
k=1. That is, for any particular data point (tk, Vk) we wish parameter choices
for which
Vk ≈ W(tk) = A0 + Ae−tk/τ sin [(ω + αtk)tk + ϕ] .
As an optimization problem, we want our parameter choices to be good for all data points.
We can take a least squares approach and minimize the following objective:
f(A0, A, τ, ω, α, ϕ) = Xn
k=1
(W(tk) − Vk)
2
(2)
=
Xn
k=1
A0 + Ae−tk/τ sin [(ω + αtk)tk + ϕ] − Vk
2
. (3)
The gradient of the objective function is
∇f =
∂f
∂A0
∂f
∂A
∂f
∂τ
∂f
∂ω
∂f
∂α
∂f
∂ϕ
= 2Xn
k=1
(W(tk) − Vk)
∂V (tk)
∂A0
∂V (tk)
∂A
∂V (tk)
∂τ
∂V (tk)
∂ω
∂V (tk)
∂α
∂V (tk)
∂ϕ
= 2Xn
k=1
(W(tk) − Vk)
1
Sk
Atkτ
−2Sk
AtkCk
At2
kCk
ACk
. (4)
where Sk = exp(−tk/τ ) sin [(ω + αtk)tk + ϕ] and Ck = exp(−tk/τ ) cos [(ω + αtk)tk + ϕ].
An alternative way of expressing this gradient is
∇f = 2∇W(t) (W(t) − V ), (5)
where W(t) ∈ R
n
is the vector of amplitudes predicted by the model, V ∈ R
n
is the
corresponding vector of data amplitudes, and ∇W(t) ∈ R
6×n
is the matrix whose columns
are the gradient vectors of the model evaluated at each time.
Complete the following tasks.
1. Obtain the data set from the course website.
2. Compose a function that returns the objective (f) and gradient (∇f) for a given set
of input parameters (A0, A, τ, ω, α, ϕ) based on the given data ({tk, Vk}).
3. Devise and employ a method for determining good initial values for the unknown
parameters (A0, A, τ, ω, α, ϕ) based on the values given in the data set. In other words,
try and automate selection of initial values that should work for any similar data set.
4. Compose general code for solving unconstrained problems in R
n by the method of
gradient descent using appropriate line search concepts. Solve the data fitting problem
using your code.
P1-2
5. Compose general code for solving unconstrained problems in R
n by conjugate gradient
using an update method of your choice. Solve the data fitting problem using your code.
6. Compose general code for solving unconstrained problems in R
n by the quasi-Newton
method using a BFGS update. Solve the data fitting problem using your code.
7. Compare your solutions and explore relative efficiency and any other interesting features that you find.
P1-3