1
$\begingroup$

Many statistical packages provide a "F statistic" and the corresponding p value in their output of a multiple linear regression.

I don't find that in the parameters of the Mathematica function LinearModelFit. All I can find are the t-statistics for the individual coefficients but no general statistic related to the regression as a whole.

I also don't find how to test interactions between predictor variables.

Where can I find the general F statistic and see the interactions between my predictor variables?

$\endgroup$
  • $\begingroup$ Hi ! All in all - good, but this isn't a question - please, edit it to reflect what you are trying to say. $\endgroup$ – Sektor Sep 7 '14 at 19:33
  • $\begingroup$ For interactions you may want to have a look at the ANOVA package. $\endgroup$ – Sjoerd C. de Vries Sep 7 '14 at 21:48
  • $\begingroup$ @SjoerdC.deVries degrees of freedom are calculated differently in ANOVA than in LinearModelFit, which changes the MS and thus F ratios. What can I really infer from ANOVA[] for the interactions in my multiple regression? $\endgroup$ – Sulli Sep 16 '14 at 14:19
  • $\begingroup$ ANOVAs are meant to examine the difference between groups or treatments and their interactions. Not mutch to do with regression. I guess I had been glancing over your question and was triggered by the one sentence about interactions. $\endgroup$ – Sjoerd C. de Vries Sep 16 '14 at 20:07
6
$\begingroup$

Assume that you have assigned the variable "model" as your LinearModelFit result. Then you can get the F statistic and its p-value with:

 model[{"ANOVATableFStatistics", "ANOVATablePValues"}]

For interactions, you can include them when you build your data. For example, assume that you have two independent variables. Build your data list as:

data = Table[{x1[[i]], x2[[i]], x1[[i]] x2[[i]], y[[i]]}, {i, 1, Length[y]}]

Now, run your regression:

model = LinearModelFit[data, {x1, x2, x3}, {x1, x2, x3}]

and examine the results:

model["ParameterTable"]

You can calculate the overall F statistic for the model by averaging the individual F statistics:

Mean[model["ANOVATableFStatistics"]]

To get the p value is a bit of a hassle, but this will do it:

1 - CDF[FRatioDistribution[Total[model["ANOVATableDegreesOfFreedom"][[1;;Length[model["ANOVATableFStatistics"]]]]], model["ANOVATableDegreesOfFreedom"][[-2]]], Mean[model["ANOVATableFStatistics"]]]

and here is a function to give both the F statistic and its p value as a list {f,p}:

fStat[m_FittedModel] := {Mean[m["ANOVATableFStatistics"]],
   1 - CDF[
     FRatioDistribution[
      Total[m["ANOVATableDegreesOfFreedom"][[1 ;; 
          Length[m["ANOVATableFStatistics"]]]]], 
      m["ANOVATableDegreesOfFreedom"][[-2]]], 
     Mean[m["ANOVATableFStatistics"]]]};,

There may be easier ways to get these values, and I'd love to see them, but this is what I was able to figure out. I don't understand why MMA doesn't have these as properties.

$\endgroup$
  • $\begingroup$ model[{"ANOVATableFStatistics"}] gives F statistics for each predictor variable, but what I would like is the overall F statistic, usually computed with the formula F = MS Model / MS Error (and the p value associated) $\endgroup$ – Sulli Sep 9 '14 at 11:56
  • $\begingroup$ thanks a lot! Yeah I'd love to see them too... $\endgroup$ – Sulli Sep 14 '14 at 18:50
  • $\begingroup$ With a F value of 570.523, I get a p value of 0. I guess it's because my F value is really really high, but how should I report this result? Prob > F = 0.000 ? $\endgroup$ – Sulli Sep 16 '14 at 13:55
  • $\begingroup$ I also see that adding the interactions to the model has changed the coefficients I had without interactions. How do I know which coefficients I should keep? $\endgroup$ – Sulli Sep 16 '14 at 14:55
  • $\begingroup$ Yes, with a high F like that I would expect the p-value to be approximately zero. Of course adding interactions will change the other coefficients, as will adding any additional variable. I would generally keep the set of variables that give you the highest F. IIRC, this is how stepwise regressions work. $\endgroup$ – Tim Mayes Sep 16 '14 at 18:35
2
$\begingroup$

Using the data from the first example in LinearModelFit one can construct the overall model F-statistic and the associated P-value.

(* Data and fit *)
data = Map[{#[[1]], #[[2]], #[[3]], 1.2 + (3.7 + RandomReal[{-1, 1}]) #[[1]] - 2 #[[2]] + 
      23.4 #[[3]]} &, RandomReal[10, {100, 3}]];
lm = LinearModelFit[data, {x, y, z}, {x, y, z}];

df = lm["ANOVATableDegreesOfFreedom"];
df = {Total[df[[1 ;; -3]]], df[[-2]]}; (* Model df and error df *)
ss = lm["ANOVATableSumsOfSquares"] ;
ss = {Total[ss[[1 ;; -3]]], ss[[-2]]}; (* Model sum of squares and error sum of squares *)

(* Construct overall model F-ratio *)
overallFratio = (ss[[1]]/df[[1]])/(ss[[2]]/df[[2]])

(* Associated P-value *)
1 - CDF[FRatioDistribution[df[[1]], df[[2]]], overallFratio]

The ANOVA table with all of the individual F-tests is found with

lm["ANOVATable"]

\begin{array}{l|lllll} \text{} & \text{DF} & \text{SS} & \text{MS} & \text{F-Statistic} & \text{P-Value} \\ \hline x & 1 & 10010.1 & 10010.1 & 886.643 & \text{2.789559773927423$\grave{ }$*${}^{\wedge}$-50} \\ y & 1 & 885.605 & 885.605 & 78.4423 & \text{4.2431663209054734$\grave{ }$*${}^{\wedge}$-14} \\ z & 1 & 424541. & 424541. & 37603.7 & \text{2.4832751695208124$\grave{ }$*${}^{\wedge}$-126} \\ \text{Error} & 96 & 1083.83 & 11.2899 & \text{} & \text{} \\ \text{Total} & 99 & 436521. & \text{} & \text{} & \text{} \\ \end{array}

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.