1
$\begingroup$

I am applying the LinearModelFit to a 90 x 5 matrix, in which the 4th variable is a nominal variable (i.e., x4 is a categorical variable with 3 levels: a, b, and c levels) and I found that the design matrix is not correctly specified in Mathematica and thus produce incorrect results due to ill-conditioning. I found that one trick is to specify all nominal variables first and before continuous variables. I appreciate if these findings are correct or not.

Here is an example data (90x5 matrix is reduced to 10x5 matrix for illustration):

SeedRandom[1234];

dataXY = Insert[ RandomReal[{0, 1}, {9, 4}] // Transpose, 
   Join[Table[a, 3], Table[b, 3], Table[c, 3]],
   4 (*4th column is nominal variable*)] // Transpose

lmf1=LinearModelFit[dataXY,{x1,x2,x3,x4},{x1,x2,x3,x4},NominalVariables -> x4]

lmf1[{"ANOVATable","ParameterTable","VarianceInflationFactors"}]

Map[MatrixForm,lmf1[{"DesignMatrix","Response"}]]

Note that the regression estimates for x4(a), x4(b), x4(c) are really big and they are identical. Also, the VarianceInflationFactors are close to infinity, indicating the data are singular. The design matrix contains all 3 levels of x4 so that it causes ill-conditioning and thus incorrect regression results. The design matrix should contain only 2 levels (e.g., b and c) and do not include one level (e.g., a) in order to avoid redundant levels.

After much much struggle, I found that if the nominal variables are specified as the 1st set of input variables before continuous variable, then the design matrix does not contain all three levels and thus computation seems correct.

SeedRandom[1234];

dataXY2 = Insert[RandomReal[{0, 1},{9, 4}]//Transpose, Join[Table[a,3],Table[b,3],Table[c,3]],
   1 ] // Transpose

lmf2 = LinearModelFit[dataXY2, {x1, x2, x3, x4}, {x1, x2, x3, x4}, 
  NominalVariables -> x1]

lmf2[{"ANOVATable", "ParameterTable", "VarianceInflationFactors"}]

Map[MatrixForm, lmf2[{"DesignMatrix", "Response"}]]

The documentation did not indicate that the nominal variable must be specified as the first sets of arguments before continuous variables. This is really critical information if someone wants to do regression-type analyses. I have not checked other commands but it is likely that these methods are also affected such as GeneralizedLinearModelFit, LogitModelFit, ProbitModelFit, NonLinearModelFit and even LeastSquares and LinearSolve. I checked the data with two nominal variables that are the 1st and 2nd variables and found the same results: nominal variables must be specified first and before the continuous variables.

Sangdon

$\endgroup$

1 Answer 1

0
$\begingroup$

I think this could be a reportable bug but the conditions where the wrong analysis is produced might be a bit more complicated.

Consider generating 4 "identical" datasets the way you did but have the position of the nominal variable change:

SeedRandom[1234];

dataXY1 = Insert[RandomReal[{0, 1}, {9, 4}] // Transpose, Join[Table[a, 3], Table[b, 3], Table[c, 3]], 1] // Transpose;
dataXY2 = Transpose[{dataXY1[[All, 2]], dataXY1[[All, 1]], dataXY1[[All, 3]], dataXY1[[All, 4]], dataXY1[[All, 5]]}];
dataXY3 = Transpose[{dataXY1[[All, 2]], dataXY1[[All, 3]], dataXY1[[All, 1]], dataXY1[[All, 4]], dataXY1[[All, 5]]}];
dataXY4 = Transpose[{dataXY1[[All, 2]], dataXY1[[All, 3]], dataXY1[[All, 4]], dataXY1[[All, 1]], dataXY1[[All, 5]]}];

Now with dataXY3 (with the nominal predictor in the 3rd position) look at 4 potential sets of LinearModelFit that would be expected to have identical parameter estimates:

lmf31 = LinearModelFit[dataXY3, {x3, x1, x2, x4}, {x1, x2, x3, x4}, NominalVariables -> x3];
lmf32 = LinearModelFit[dataXY3, {x1, x3, x2, x4}, {x1, x2, x3, x4}, NominalVariables -> x3];
lmf33 = LinearModelFit[dataXY3, {x1, x2, x3, x4}, {x1, x2, x3, x4}, NominalVariables -> x3];
lmf34 = LinearModelFit[dataXY3, {x1, x2, x4, x3}, {x1, x2, x3, x4}, NominalVariables -> x3];
#["ParameterTable"] & /@ {lmf31, lmf32, lmf33, lmf34} // TableForm

ParameterTables when order of occurrence changes

One can see that the first two work fine but the second two do not.

A modified set of conditions might be that the nominal variable needs to occur before two continuous variables are given as basis functions.

The same rule applies for any other position of the nominal variable. (I did not try multiple nominal predictors.)

If one includes the option IncludeConstantBasis -> False (and that is appropriate for this particular case), then the nominal predictor variable can be in any position (other than the very last position as that location is reserved for the response variable).

You should report this to Wolfram, Inc.

$\endgroup$
1
  • $\begingroup$ Thank you very much. I reported this to Wolfram. $\endgroup$
    – SDL
    Commented May 28, 2021 at 12:39

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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