MTH 420-520 LAB 2

In this lab you will explore different aspects of least squares.

  • 420 turns in A, (C or D), and (E or G).
  • The omegas are encouraged to solve more, and/or attempt B, F.


Groups alpha and mu should follow demonstrations.
DEMONSTRATION 14:05-14:10 Objects and rotations
  1. Set up an object (curve) as a matrix of points. Each column is a point. Recall that the cardioid curve has the parametric equation r=1-cos(theta).
    theta=linspace(0,2*pi);
    cardiod = [cos(theta).*(1-cos(theta));sin(theta).*(1-cos(theta))];
    object = cardioid;
    rot =- [1 -1; 1 1];
    plot(object(1,:),object(2,:)); %% plots the x coordinates againts y coordinates of all points.
    object = a*object; plot(object(1,:),object(2,:)); %% plots the rotated object
    Extra: try playing also with limacon r=1+cos(theta), or other cool shapes of your choice.

PROJECTS: linear least squares
  1. Given the following data t=[0.5, 1.5, 2, 5]; b = [3.1 6.8 10.3 25]; setup the matrix A and right hand side b, and solve the least squares problem for x=[a;c], fitting the data to the linear model b=at+c
    (i) Set up the problem by hand and get the solution using MATLAB. (As part of your solution describe the process of finding the critical points of the function phi(a,c)). Once you have the coefficients a,c, plot the curve and mark the points.
    Here is how you would plot if a=5,c=-1.7; plot(t,b,'*',t,5*t-1.7,'-');legend('data','linear fit');
  2. MATLAB function polyfit can be used to achieve the same as what you did in A. Leanr how and compare to A.

PROJECTS: Exploring nonlinear least squares
  1. Derive by hand calculations the equations that will give you the coefficients a,d,c of a quadratic function fit to the data from A
    b=at2 + dt +c
  2. Derive by hand the equations that would fit the data from A to an exponential model
    b=ceat
    How would you solve this problem ? Describe the challenges. Can you use linear LSQ ? What if you knew a and wanted c ? What if you knew c and wanted a ? How can you combine the two ?

    As an example, b could be the concentration in time t of microbes that are growing with rate a with initial condition c. Extra: provide details why this exponential model might be a good model. Can you think of an appropriate differential equation corresponds to microbe growth ?

    Extra: learn how to use lsqnonlin to solve C and D.

Yaw, pitch, and roll [OL Ex. 12.2/p159]
  1. In this example we explore both: nonlinear data fitting and linear (actually, orthogonal) transformations.

    First, we experiment in 2D. Consider object = [1 0 0;0 0 2 ]
    Plot to see what it represents plot(object(1,:),object(2,:));axis ([ -2 2 -2 2]);axis equal
    Set up the rotation matrix with some chosen angle phi
    rot = [ cos(phi) sin(phi); -sin(phi) cos(phi)];
    In each step transform your object and plot to see how it rotates.
    object = rot*object;

    Now assume that you are given the initial (above) and the final (possibly slightly perturbed) position of the object
    objectnew = [0.9557 0.0965 0.9648; -0.3381 0.0158 1.8977]
    but that you do not know the angle phi of rotation. You also know that there was a perturbation during the rotation so that there is no exact solution (how would you find the exact solution?). Set up the equations to solve this problem for the best value of phi with least squares. (What if there was no perturbation ?)
  2. Solve the problem with lsqnonlin.

    Make sure you set up bounds for the angle phi. Take advantage of options that you can set up for lsqnonlin. Explore sensitivity and accuracy of the routine. Extra: Learn more about it: what is the underlying method ?

  3. Now we move on to 3D.

    Recall the setup of Euler angles of rotation. Yaw: rotation in xy plane. Roll: rotation in yz plane. Pitch: in xz plane.
    Pick some angles phi, psi, theta. Set up the matrices (3 x 3) corresponding to each of these rotations, and calculate the matrix Q which composes them together. Play with showing the object in 3D (use plot3)

    Describe how you would solve the problem analogous to E in 3D. You can start with the object object=[ 1 0 0; 0 0 2; -1 -2 -3];
    Solve the problem in 3D (for all three angles) when As newobject you can use the "true" object rotated with phi=pi/3;psi=pi/3;theta=pi/3, which you then slightly perturb.