Curve Fitting

Part 1.  Example: Quadratic Fit to U.S. Population Data

Answer questions 1 through 6 here.

Enter the U.S. population data and plot it:
  (Note: times are measured in years since 1900)

\(times = {0, 10, 20, 30, 40, 50, 60, 70, 80, 90}; 

    pops = {75.996, 91.972, 105.711, 122.775, 131.669, 150.697, 179.323, 
        203.185, 226.546, 248.710}; 
popPairs = Transpose[{times, pops}]; 

    popPlot = 
      ListPlot[popPairs, PlotStyle -> {PointSize[0.02], RGBColor[0, 0, 1]}, 
        PlotRange -> {0, 300}]; \)

Enter the vectors 1, t, t2, and y and form the matrix X:

\(Clear[one, t, y]; 
one\  = \ {1, 1, 1, 1, 1, 1, 1, 1, 1, 1}; 

    t\  = \ times; 
t2\  = \ t^2; 
y = pops; 

    X\  = \ Transpose[{one, \ t, \ t2}]; 
X // MatrixForm\)

Solve the normal equations and enter the values of the parameters a, b, and c below:

\(a\  = \ \(?
b\ =\ ?\)\)
c  =  ?

Define the least squares quadratic and show its plot along with the U.S. population data.

q[t_] := a  +  b*t  +  c*t^2\(quadPlot = 
      Plot[q[t], {t, \(-10\), 100}, PlotStyle -> RGBColor[1, 0, 0], 
        DisplayFunction -> Identity]; 

    Show[popPlot, quadPlot, DisplayFunction -> $DisplayFunction]; \)

Give the computations and answers for steps 9 and 10 below.

Part 2.  Linear Least Squares

Enter the Signal Strength data and plot it:

\(Clear[t, y]; 

    t = {\(-2.0\), \(-1.5\), \(-1.0\), \(-0.5\), 0.0, 0.5, 1.0, 1.5, 2.0, 
        2.5, 3.0}; 

    y = {\(-6.32\), \(-3.23\), 1.62, 3.13, 1.74, \(-0.75\), \(-1.41\), 1.78, 
        8.88, 9.98, 7.10}; 
sigPairs = Transpose[{t, y}]; 

    sigPlot = 
      ListPlot[sigPairs, PlotStyle -> {PointSize[0.02], RGBColor[0, 0, 1]}, 
        PlotRange -> {\(-8\), 12}]; \)

Enter the vectors 1, s1, c1, s2, and c2 and form the matrix X:

\(Clear[one, X]; 
one\  = \ {1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1}; 

    s1 = Sin[t]; 
c1 = Cos[t]; 
s2 = Sin[2*t]; 
c2 = Cos[2*t]; 

    X\  = \ Transpose[{one, s1, c1, s2, c2}]; 
X // MatrixForm\)

Solve the normal equations and enter the values of the parameters a0, a1, b1, a2, and b2 below:

\(a0\  = \ \(?
a1\ =\ ?\)\)\(b1\  = \ \(?
a2\ =\ ?\)\)b2  =   ?

Define the least squares trig polynomial and show its plot along with the U.S. population data.

p[t_] := a0  +  a1*Sin[t] + b1*Cos[t] + a2*Sin[2*t] + b2*Cos[2*t]\(tpolyPlot = 
      Plot[p[t], {t, \(-2\), 3}, PlotStyle -> RGBColor[1, 0, 0], 
        DisplayFunction -> Identity]; 

    Show[sigPlot, tpolyPlot, DisplayFunction -> $DisplayFunction]; \)

For step 3, do the computations and make the residual plot below.
 Discuss the goodness of fit.

Part 3.  Linearization

Enter the vectors 1, t, and y* and form the matrix X:

\(Clear[one, t, y]; 
one\  = \ {1, 1, 1, 1, 1, 1, 1, 1, 1, 1}; 

    t\  = \ times; 
y = pops; 
ystar = Log[y]; 

    X\  = \ Transpose[{one, \ t}]; 
X // MatrixForm\)

Solve the normal equations, find and enter the values of the parameters a and b below:

\(a\  = \ \(?
b\ =\ ?\)\)

Define the exponential function corresponding to your semilog fit and show its plot along with the U.S. population data.

s[t_] := a*E^(b*t)\(sePlot = 
      Plot[s[t], {t, \(-10\), 100}, PlotStyle -> RGBColor[1, 0, 0], 
        DisplayFunction -> Identity]; 

    Show[popPlot, sePlot, DisplayFunction -> $DisplayFunction]; \)

For step 3, do the computations and make the residual plot for your exponential model below.  Discuss the goodness of fit.

For step 4, enter the coefficients a, b, and c you found in Part 1 for the quadratic model and then make the residual plot.  Answer all questions about the comparison between the fits of the two models.

Clear[a, b, c, t];
\(a\  = \ \(?
b\ =\ ?\)\)
c  =  ?

q[t_] := a  +  b*t  +  c*t^2\(resid = pops - q[times]; 
resPairs = Transpose[{times, resid}]\)ListPlot[resPairs, 
      PlotStyle -> {PointSize[0.02], RGBColor[0, 0, 1]}];

Step 5: Power function fit to small data set 
        Enter the data below and make a scatterplot:

\(Clear[t, y]; 
t = {1.0, 2.0, 3.0, 4.0}; 
y = {2.5, 8.0, 19.9, 50.0}; 
    
smallPairs = Transpose[{t, y}]; 

    smallPlot = 
      ListPlot[smallPairs, PlotStyle -> {PointSize[0.02], RGBColor[0, 0, 1]}, 
        PlotRange -> {0, 60}]; \)

Enter the vectors 1, t*, and y* and form the matrix X:

\(Clear[one]; 
one\  = \ {1, 1, 1, 1}; 
tstar\  = \ Log[t]; 

    ystar = Log[y]; 
X\  = \ Transpose[{one, \ tstar}]; 
X // MatrixForm\)

Solve the normal equations, find and enter the values of the parameters a and b below:

\(a\  = \ \(?
b\ =\ ?\)\)

Define the power function corresponding to your semilog fit and show its plot along with the small data set.

pow[t_] := a*t^b\(powPlot = 
      Plot[pow[t], {t, 0, 5}, PlotStyle -> RGBColor[1, 0, 0], 
        DisplayFunction -> Identity]; 

    Show[smallPlot, powPlot, DisplayFunction -> $DisplayFunction]; \)

For step 7, do the computations and make the residual plot for your power model below.  Discuss the goodness of fit.

Part 4.  Pitfalls of Linearization

Answer question 1 below:

Enter the data and the vectors 1, t, and y* and form the matrix X:

\(Clear[one, t, y, ystar]; 
t = {1.0, 2.0, 3.0, 4.0}; 

    y = {2.5, 8.0, 19.9, 50.0}; 
one\  = \ {1, 1, 1, 1}; 
ystar = Log[y]; 

    X\  = \ Transpose[{one, \ t}]; 
X // MatrixForm\)

Solve the normal equations, enter the values of the parameters a and b below and define your exponential model:

\(a\  = \ \(?
b\ =\ ?\)\)exp[t_] := a*E^(b*t)

The following commands make a residual plot of the two different exponential models, the one you just found in blue and the model
 y = 0.921 E^(0.999t) in red.

alt[t_] := 0.921*E^(0.999*t)\(resid0 = y - alt[t]; 
altPairs = Transpose[{t, resid0}]\)\(altPlot = 
      ListPlot[altPairs, PlotStyle -> {PointSize[0.02], RGBColor[1, 0, 0]}, 
        DisplayFunction -> Identity]; 
resid = y - exp[t]; 

    expPairs = Transpose[{t, resid}]\)\(expPlot = 
      ListPlot[expPairs, PlotStyle -> {PointSize[0.02], RGBColor[0, 0, 1]}, 
        DisplayFunction -> Identity]; 

    Show[altPlot, expPlot, DisplayFunction -> $DisplayFunction]; \)

Do the computations and discuss the comparison of the two exponential models below.

Part 5.  Nonlinear Curve Fitting

Fit Power Function to Small Data Set
      Enter the small data set data:

\(Clear[t, y]; 
t = {1.0, 2.0, 3.0, 4.0}; 
y = {2.5, 8.0, 19.9, 50.0}; 
    \)

Step 1:  Do the computations.  Then enter your initial guess (a0,b0) for a and b below:

\(a0\  = \ \(?
b0\ =\ ?\)\)

Enter the definitions of the power model and its partial derviatives:

fun[t_, a_, b_] := a*t^bfuna[t_, a_, b_] := t^bfunb[t_, a_, b_] := a*t^b*Log[t]

Step 2A: Form the necessary vectors and the matrix X:

a = a0b = b0\(ystar = y - fun[t, a, b]; 
fa = funa[t, a, b]; 
fb = funb[t, a, b]; 

    X\  = \ Transpose[{fa, fb}]; 
X // MatrixForm\)

Step 2B and 3:  Solve the normal equations and update parameters a and b

Step 4 and 5:  Repeat steps 2 and 3 using updated values of a and b each time.

Step 6:  Set a0 and b0 to their original values and execute the following loop.

\( (*\ Initialize : \ *) 
a = a0\ \)b = b0\(iter\  = \ 0; 
da = 1.0*10^6; \ db = 1.0*10^6; 

     (*\ Iterate\ until\ changes\ to\ parameters\ are\ tiny\ *) 

    While[\ \((Sqrt[da^2 + db^2] > 10^\(-4\))\)\  && \ \((iter < 20)\), 

      \ \ \ \ \ \ ystar\  = \ y\  - \ fun[t, a, b]; 
\ \ \ \ \ \ 
      fa\  = \ funa[t, a, b]; 
\ \ \ \ \ \ fb\  = \ funb[t, a, b]; 

      \ \ \ \ \ \ X\  = \ Transpose[{fa, \ fb}]; 

      \ \ \ \ \ {da, db}\  = LinearSolve[Transpose[X].X, Transpose[X].ystar]; 
      
\ \ \ \ \ iter\  = \ iter\  + \ 1; 
\ \ \ \ \ a = a + da; 

      \ \ \ \ \ b = b + db; 
\ \ \ \ \ 
      Print["iter= ", iter, "  a= ", a, "  b= ", b]
\ \ ]\)

Step 6 (cont):  Compute the residuals and S for your fit and discuss below.

Step 7:  Fit Logistic Model to U.S. Populations
   Enter the definitions of the logistic model and its partial derviatives:

\(Clear[t, P0, K, r]; 

    fun[t_, P0_, K_, r_] := 
      K*P0/\((\ P0\  + \ \((K\  - \ P0)\)*Exp[\(-r\)*t]\ )\)\)funP0[t_, P0_, K_, r_] = D[ fun[t, P0, K, r],  P0]funK[t_, P0_, K_, r_] = D[ fun[t, P0, K, r],  K]funr[t_, P0_, K_, r_] = D[ fun[t, P0, K, r],  r]

Enter the U.S. population data:

\(Clear[t, y]; 
t = times\)y = pops

Enter the initial values of the parameters:

\(P01 = 78; 
K0 = 700; 
r0 = 0.0168; \)

Step 8A: Form the necessary vectors and the matrix X:

P0 = P01K = K0r = r0\(ystar = y - fun[t, P0, K, r]; 
fP0 = funP0[t, P0, K, r]; 

    fK = funK[t, P0, K, r]; 
fr = funr[t, P0, K, r]; 

    X\  = \ Transpose[{fP0, fK, fr}]; 
X // MatrixForm\)

Step 8B:  Solve the normal equations and update parameters P0, K, and r.

Step 9:  Set P0, K, and r to their original values and execute the following loop.

\( (*\ Initialize : \ *) 
P0 = P01\ \)K = K0r = r0\(iter\  = \ 0; 
dP0 = 10^6; \ dK = 10^6; \ dr = 10^6; 

     (*\ Iterate\ until\ changes\ to\ parameters\ are\ tiny\ *) 

    While[\ \((Sqrt[dP0^2 + dK^2 + dr^2] > 10^\(-4\))\)\  && \ 
        \((iter < 20)\), 
\ \ \ \ \ \ ystar = y - fun[t, P0, K, r]; 

      \ \ \ \ \ \ fP0 = funP0[t, P0, K, r]; 
\ \ \ \ \ \ 
      fK = funK[t, P0, K, r]; 
\ \ \ \ \ \ fr = funr[t, P0, K, r]; 

      \ \ \ \ \ \ X\  = \ Transpose[{fP0, fK, fr}]; 

      \ \ \ \ \ {dP0, dK, dr}\  = 
        LinearSolve[Transpose[X].X, Transpose[X].ystar]; 
\ \ \ \ \ 
      iter\  = \ iter\  + \ 1; 
\ \ \ \ \ P0 = P0 + dP0; 
\ \ \ \ \ 
      K = K + dK; 
	\ \ \ r = r + dr; 
\ \ \ \ \ 
      Print[iter, "  P0= ", P0, "  K= ", K, "  r= ", r]
\ \ ]
      \)

Step 10:  Ignore the warning messages above.  Fill in the values of the parameters P0, K, and r and then plot the logistic fit along with the originial population data:

\(P0 = \(?
K\ =\ ?\)\)r =  ?

\(logiPlot = 
      Plot[fun[t, P0, K, r], {t, \(-10\), 100}, 
        PlotStyle -> RGBColor[1, 0, 0], DisplayFunction -> Identity]; 

    Show[popPlot, logiPlot, DisplayFunction -> $DisplayFunction]; \)

Step 11:  Make the residual plot and discuss the comparison of the fit with the quadratic and the exponential models below.

Part 6.  Summary

Answer the questions here.