Mathematica Stack Exchange is a question and answer site for users of Mathematica. Join them; it only takes a minute:

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer
  3. The best answers are voted up and rise to the top

I've done my best to search the help and google White's test ( in Mathematica. Can anyone point me to a function or Mathematica reference for this or do you have to code it up yourself if you want it?

Thanks, Peter

share|improve this question
while i could give an answer there are examples to do the White test and it's alternative, the Breusch-Pagan. Please have a look at LinearModelFit under Applications. There you'll find examples for White's heteroskedasticity-consistent covariance estimate and Breusch-Pagan test statistics. – Stefan Mar 1 '14 at 21:30
@stefan, thanks I didn't see them buried in the examples. – pjc42 Mar 2 '14 at 5:02
up vote 9 down vote accepted

I was playing around a little with what I've learned from rm-rf.

Here's a small package which encapsulates the White and the Breusch-Pagan test. The math is correct, but I wanted to have a kind of statistical distribution behaviour, where you can query certain properties using a string.

I'd appreciate to have overloaded the LinearModelFit, but decided to do it differently, since I was concerned on creating more trouble.


FitModel::usage = "";


FitModel[data_List, coords_List] :=
    With[{lm = LinearModelFit[data, coords, coords]},
        dispatch[{"BreuschPagan" -> BreuschPagan[data, lm, coords],
                  "WhitesHeteroskedasticity" -> WhitesHeteroskedasticity[lm, coords]}]
dispatch[list_][field_] := field /. list
dispatch[list_]["Properties"] := list /. Rule[field_, _] :> field
dispatch /: ReplaceAll[fields_, dispatch[list_]] := fields /. list
Format[dispatch[list_], StandardForm] := 
    HoldForm[dispatch]["<" <> ToString@Length@list <> ">"]

BreuschPagan[data_List, lm_, coords_List] := Module[{lm2},
    lm2 = LinearModelFit[Block[{newdata = data},
    newdata[[All, -1]] = lm["FitResiduals"]^2; newdata], coords, coords];
    With[{sqResids = 
         lm["FitResiduals"]^2}, (Variance[sqResids] (Length[data] - 1) - 
           Total[lm2["FitResiduals"]^2]) / 2 / (Total[sqResids]/Length[data])^2]

WhitesHeteroskedasticity[lm_, coords_List] := 
    Module[{desmat, resids},
        {desmat, resids} = lm[{"DesignMatrix", "FitResiduals"}];
        With[{inv = Inverse[Transpose[desmat].desmat], 
            xresid = resids*desmat},


To test it:

data = Flatten[ Table[{x, y, z, 1.2 x - 3.4 y + 10 z + RandomReal[10]}, {x, 
 RandomReal[10, 3]}, {y, RandomReal[10, 3]}, {z, 
 RandomReal[10, 3]}], 2];

f = FitModel[data, {x, y, z}];
    => yields with actual data: 0.264212

f["WhitesHeteroskedasticity"] // MatrixForm
=> yields:
(2.95722    -0.262317   -0.294977   -0.0841742
 -0.262317  0.108253    0.022681    -0.0249064
 -0.294977  0.022681    0.0542495   -0.0110319
 -0.0841742 -0.0249064  -0.0110319  0.0409517)

I really hope this helps.

share|improve this answer
brilliant, thanks – pjc42 Mar 2 '14 at 5:03

Your Answer


By posting your answer, you agree to the privacy policy and terms of service.

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