Homework in Preparation for the Midterm

Donna Odhiambo

1. Use Newton’s method to find a root of the non-linear function: f(x)=3-2(e^(x/2)) starting from initial guess (Subscript[x, 0])=-1.548

f[x_]:=3-2(E^(x/2)); x0=-1.548; r=Log[9/4];

i. true value of root: r=ln((9/4))
ii. Draw the first tangent line to the graph of f, through the point ((Subscript[x, 0]), f((Subscript[x, 0])))

Plot[{f[x], f'[x0](x-x0)+f[x0]},{x,-5,5},PlotLegends->"Expressions", Epilog->{PointSize[Medium],Point[{x0,f[x0]}],Point[{r,0}]}]

iii. Find the root (Subscript[x, 1]) of the tangent line function, (Subscript[T, 1])((Subscript[x, 1]))=0 and use this as an improved guess. 
iv. Iterate this process four more times to produce a table of estimates (Subscript[x, 0]),(Subscript[x, 1]),(Subscript[x, 2]),(Subscript[x, 3]),(Subscript[x, 4]),(Subscript[x, 5]).

guess=x0;
list={};
n=0;
While[n<5,
newguess=Solve[ f'[guess](x-guess)+f[guess]==0,x];
guess=x/.newguess[[1]];
AppendTo[list,guess];
n++]
Grid[{{"(Subscript[x, 0]):",x0},{"(Subscript[x, 1]):",list[[1]]},{"(Subscript[x, 2]):",list[[2]]},{"(Subscript[x, 3]):",list[[3]]},{"(Subscript[x, 4]):",list[[4]]},{"(Subscript[x, 5]):",list[[5]]},{"r:",N[r]}},Frame->All,Alignment->"."]

v. Discuss:
The estimates from (Subscript[x, 1]) on decrease as they approach the real root r. The last estimate, (Subscript[x, 5]), and root r have a difference of 7.66803×(10^(-6)) which is close to zero, however if another iteration is done and estimate (Subscript[x, 6]) is found, the difference between that and r would be 1.46997×(10^(-16)) which is closer to zero, and enough so to state that the root has been found with satisfactory precision.

2. Use non-linear regression to fit a Hollings Type III model, Scaled logistic model, and von Bertalanffy model to the cadaver cumulative data.

Clear["Global`*"]
bodies={2,14,31,14,20,11,7,4,7,1,1,2} ;
hours=Table[i,{i,2,13}]; 
cumulative=Accumulate[bodies];
ybar=Mean[cumulative];
data=Table[{hours[[i]],cumulative[[i]]},{i,1,12}];
p1=ListPlot[data,PlotRange->{{0,15},{0,120}},AxesLabel->{"bodies","hours"}]

i. Hollings Type III model: f(x)=((ax^2)/((b^2)+(x^2)))

Clear[a,b]
Manipulate[
Module[
{},
p2=Plot[((a (x^2))/((b^2)+(x^2))),{x,0,15}];
Show[p1,p2]
]
,{{a,133},100,150,1}
,{{b,4.75},0,50,.01}
]

NLM1=NonlinearModelFit[data, ((a (x^2))/((b^2)+(x^2))),{{a,133},{b,4.75}},x]

Show[p1, Plot[NLM1[x],{x,0,15}]]
r1=ListPlot[NLM1["FitResiduals"],PlotLabel->"Hollings Type III Model",AxesLabel->{"Fitted Value", "Residual"},PlotRange->{Automatic,{-15,10}}];

NLM1["ParameterConfidenceIntervalTable"]
NLM1["ParameterTable"]

yhat=Table[NLM1[i], {i,hours}];
residuals=cumulative-yhat;
rsquared1=1-residuals.residuals/((cumulative-ybar).(cumulative-ybar));
n=Length[hours];
p=2; (* number of parameters estimated *)
rsquaredadj1=1-(n-1)/(n-p) (1-rsquared1);

ii. Scaled Logistic Model: f(x)=s (((e^(a+bx))/(1+(e^(a+bx)))))

Clear[a,b,s]
Manipulate[
Module[
{},
p3=Plot[s(((Exp[a+b x])/(1+Exp[a+b x]))),{x,0,15}];
Show[p1,p3]
]
,{{a,-4.00},-10,0,0.01}
,{{b,0.74},0,1,.01}
,{{s,114},110,115,1}
]

NLM2=NonlinearModelFit[data,s(((Exp[a+b x])/(1+Exp[a+b x]))),{{a,-4.00},{b,0.74},{s,144}},x]

Show[p1, Plot[NLM2[x],{x,0,15}]]
r2=ListPlot[NLM2["FitResiduals"],PlotLabel->"Scaled Logistic Model",AxesLabel->{"Fitted Value", "Residual"},PlotRange->{Automatic,{-15,10}}];

NLM2["ParameterConfidenceIntervalTable"]
NLM2["ParameterTable"]

yhat=Table[NLM2[i], {i,hours}];
residuals=cumulative-yhat;
rsquared2=1-residuals.residuals/((cumulative-ybar).(cumulative-ybar));
n=Length[hours];
p=3; (* number of parameters estimated *)
rsquaredadj2=1-(n-1)/(n-p) (1-rsquared2);

iii. von Bertalanffy model: f(x)=((a(1-(e^(-k(a-d)x))))^((1/(1-d))))

Clear[a,d,k]
Manipulate[
Module[
{},
p4=Plot[a(1-Exp[-k(a-d)x])^(1/(1-d)),{x,0,15}];
Show[p1,p4]
]
,{{a,114},110,120,1}
,{{d,0.81},0,0.99,.01}
,{{k,0.004},0,0.009}
]

NLM3=NonlinearModelFit[data,a(1-Exp[-k(a-d)x])^(1/(1-d)),{{a,114},{d,0.81},{k,0.004}},x]

Show[p1, Plot[NLM3[x],{x,0,15}]]
r3=ListPlot[NLM3["FitResiduals"],PlotLabel->"von Bertalanffy Model",AxesLabel->{"Fitted Value", "Residual"},PlotRange->{Automatic,{-15,10}}];

NLM3["ParameterConfidenceIntervalTable"]
NLM3["ParameterTable"]

yhat=Table[NLM3[i], {i,hours}];
residuals=cumulative-yhat;
rsquared3=1-residuals.residuals/((cumulative-ybar).(cumulative-ybar));
n=Length[hours];
p=3; (* number of parameters estimated *)
rsquaredadj3=1-(n-1)/(n-p) (1-rsquared3);

iv. Compare all models to Dr. Long’s model
-Long’s Model: f(x)=((ae)^(-(b/(x^c))))

NLM=NonlinearModelFit[data,a Exp[-(b/x^c)],{{a,120},{b,26.2772},{c,2.39415}},x]

Show[p1, Plot[NLM[x],{x,0,15}]]
r4=ListPlot[NLM["FitResiduals"],PlotLabel->"Long's Model",AxesLabel->{"Fitted Value", "Residual"},PlotRange->{Automatic,{-15,10}}];

NLM["ParameterConfidenceIntervalTable"]
NLM["ParameterTable"]

yhat=Table[NLM[i], {i,hours}];
residuals=cumulative-yhat;
rsquared=1-residuals.residuals/((cumulative-ybar).(cumulative-ybar));
n=Length[hours];
p=3; (* number of parameters estimated *)
rsquaredadj=1-(n-1)/(n-p) (1-rsquared);

COMPARISONS

Below are residual plots of all models generated as well as the adjusted (R^2) values for the fits. The model that gives the best fit to the data is Dr. Long’s model, followed by the von Bertalanffy model, the scaled logistic model, then lastly the Hollings type III model. The residuals for both the Hollings type III model and the scaled logistic model do not look entirely random and at values of 0.963 and 0.983, these models have the lowest adjusted (R^2). The von Bertalanffy model fares much better at 0.995 for its adjusted (R^2) and a residual plot that shows less patterned behaviour, however Dr. Long’s model had an adjusted (R^2) of 0.997 along with a residual plot that looks completely random. While some models of this cadaver data clearly fare better than others, no model calculated was irrelevant, as shown by each model’s parameter’s confidence interval excluding zero, the p-values being less than 0.05, and the t-statistic values being extreme in nature.

GraphicsGrid[{{r4,r1},{r2,r3}},Frame->True]

Grid[{{"Model","(R^2)","Adjusted (R^2)"},{"Long", rsquared, rsquaredadj},{"Hollings Type III", rsquared1,rsquaredadj1},{"Scaled Logistic",rsquared2,rsquaredadj2},{"von Bertalanffy", rsquared3,rsquaredadj3}},Frame->All]