revision of regression package
- To: mathgroup at yoda.physics.unc.edu
- Subject: revision of regression package
- From: Gossett <gossett at bethel.edu>
- Date: Mon, 3 Aug 1992 10:20:40 -0500
This revision of the package regress.m fixes problems with vertical and
horizontal data sets.
BeginPackage["Regression`"]
Regress::usage = "Regress[Data,Predictions,X_Plot_Range,Y_Plot_Range] does a linear regression on the data points. Data is either a list with the y values, or a list of {x,y} pairs. Predictions is an optional list of x values at which the least squares function is to be evaluated. X_Plot_Range and Y_Plot_Range are lists with an estimated lower and upper value for the axis limits. For example, Regress[{2,4,5,6,8,9},{3,7},{0,10},{0,10}] or
Regress[{{-1,3},{3,4},{4,4}},{-2,5},{0,6.25}]"
(* This package may be freely used for non-commercial purposes. It has been
designed for educational purposes, but may be useful in other contexts.
Copyright June 1992
Dr. Eric Gossett
Bethel College
3900 Bethel Drive
St. Paul, MN 55112
gossett at bethel.edu
Version 1.1 (Fixes problems with vertical and horizontal data sets).
*)
Begin["`Private`"]
Unprotect[Regress]
Regress[Indata_List/;(Length[Dimensions[Indata]]<3),Predictions_List:{},
Domain_List,Range_List] :=
Module[{Data,L,Pmax,Pmin,NumData,Xvals,Yvals,CalcTable,
Xsum,Ysum,XYsum,X2sum,preds},
(* Finalize the data array: get length, add x coordinates if needed *)
NumData = Length[Indata];
If[Length[Dimensions[Indata]]==1,
Data = Table[{i,Indata[[i]]},{i,1,NumData}],
Data = Indata];
(* Check for vertical lines *)
Xvals = Transpose[Data][[1]];
Yvals = Transpose[Data][[2]];
If[Union[Flatten[Outer[(Chop[#1-#2])&,Xvals,Xvals]]] == {0},
Print["The data forms a vertical line"]; Return[]];
(* Compute the least squares line *)
L = Chop[Fit[Data,{1,x},x]];
(* Echo the original data *)
Print[" "];
Print["The Data set is:"];
Print[" "];
Print[TableForm[Data,TableAlignments->{Right,Right}]];
Print[" "];
(* Print the equation of the least squares line *)
If[SameQ[Head[L],Plus],
Print["The equation of the line is y = ",L[[2]][[1]]," x + ",
L[[1]] ],
If[FreeQ[L,x],Print["The equation of the line is y = ",L],
Print["The equation of the line is y = ",L[[1]]," x"]
]
];
Print[" "];
(* Plot the data and the least squares line *)
Pmax = Max[Map[(#[[2]])&,Data],
{L/.x->Domain[[1]],L/.x->Domain[[2]]},Range];
Pmin = Min[Map[(#[[2]])&,Data],
{L/.x->Domain[[1]],L/.x->Domain[[2]]},Range];
Show[ListPlot[Data,PlotRange->{Domain,{Pmin,Pmax}},
PlotStyle->{PointSize[.0125]},
DisplayFunction->Identity],
Plot[L,{x,Domain[[1]],Domain[[2]]},PlotRange->{Domain,{Pmin,Pmax}},
DisplayFunction->Identity],
DisplayFunction->$DisplayFunction];
(* Print any requested estimates using the least squares line *)
If[Length[Predictions] >0,
Print["The predicted values are:"];
Print[" "];
preds=Transpose[{Predictions,Table[L/.x->Predictions[[i]],
{i,1,Length[Predictions]}]}];
preds = Prepend[preds,{"X","Y"}];
Print[TableForm[preds,TableAlignments->{Right,Right}]];
];
(* Print the hand calculation table *)
Print[" "];
Print["The hand calculations are:"];
Print[" "];
CalcTable = Table[{Xvals[[n]],Yvals[[n]],Xvals[[n]] Yvals[[n]],
Xvals[[n]]^2},{n,1,NumData}];
Xsum = Apply[Plus,Xvals];
Ysum = Apply[Plus,Yvals];
XYsum = Apply[Plus,Part[Transpose[CalcTable],3]];
X2sum = Apply[Plus,Part[Transpose[CalcTable],4]];
CalcTable = Prepend[CalcTable,{"X","Y","XY","X^2"}];
CalcTable = Append[CalcTable,{"-----","-----","-----","-----"}];
CalcTable = Append[CalcTable,{Xsum,Ysum,XYsum,X2sum}];
CalcTable = Append[CalcTable,{Xsum^2," "," "," "}];
Print[TableForm[CalcTable,TableAlignments->{Right,Right}]];
] (* end Regress *)
End[]
Protect[Regress]
EndPackage[]