Model fitting and hypothesis testing

Dear all,

I have started using MATLAB statistics toolbox and have a few problems in model fitting. I have read the manual for glm, however I'm still confused with my own data. I have 3 groups of patients,with two nuisance covariates and one dependent variable. I would like to compare the measured dependent variable among three groups, controlling for two nuisance covariates. Statistics toolbox provides several functions for model fitting which is a little confusing for me. I'm confused whether I should start with glmfit or GeneralizedLinearModel and how I should proceed, whether with coeftest, glmval . I would be grateful if someone point me where I should begin and proceed.

-Arman

Answer 1

[Tom Lane]

The GeneralizedLinearModel interface is intended to be a step forward over glmfit, though the latter is still there. The new version, for example, will take a grouping variable as a predictor without your having to create a dummy variable for it. This version produces an object that provides methods for plotting and for testing the coefficients, among other things.

There's still a place for using glmfit, though. For instance if you wanted to do bootstrapping and didn't need the features of GeneralizedLinearModel, you'd find that glmfit is faster because it avoids any overhead in computing various parts of the fit object.

[Arman]

Thanks Tom. Would you point me in the right direction for hypothesis testing in GeneralizedLinearModel as well? An example is mentioned in my above post. I'm a little puzzled and was not able to find an example similar to my study in the documentations, so if you only show me what you would do in the mentioned example it would be great.

[Tom Lane]

Try this:

load carsmall
d = dataset(MPG,Weight,Cylinders);
d.Cylinders = ordinal(d.Cylinders);
m = GeneralizedLinearModel.fit(d,'MPG~Weight*Cylinders')
devianceTest(m)
coefTest(m,[1 0 0 0 0 0])
coefTest(m,[0 0 0 0 1 1])

The first commands fit a generalized linear model to some data including a categorical predictor. It turns out this is a normal linear model but we're using the GeneralizedLinearModel class. The devianceTest function tests for significance of the model as a whole. The first coefTest function tests for the significance of the intercept, and you should see that the p-value matches the one shown in the coefficients table display.

The second coefTest function tests for the significance of the interaction term. This is similar to your case, where the term involves a three-level grouping variable and it is represented by two coefficients in the model.

I believe this is somewhat similar to your example. Here Cylinders is a categorical predictor similar to your three-group predictor. Weight could play the role of a nuisance variable.

[Arman]

Thanks a lot for your answer!

Answer 2

[John]

Hi Tom,

I'm having a bit of a problem understanding what the coefTest function is doing (and can't seem to find where the root code is to look at it). I ran a fitlm model on some data and got the following stats matrix:

                          tStat       pValue    
    (Intercept)            -7.0023    1.1537e-10
    Early                   2.4744      0.014615
    Correct               -0.56112       0.57567
    Bins                    8.0292    4.7528e-13
    Early:Correct          -3.0349     0.0028993
    Early:Bins             -2.8373     0.0052691
    Correct:Bins           0.64341       0.52107
    Early:Correct:Bins      3.4799    0.00068014

I had someone replicate this same linear regression in STATA and got the same values. Now I want to test the significance between Early x Bins and Correct x Bins. When you do this in STATA (test Bins_Correct = Bins_Early) you get F=12.12 p=.0007. When I run coefTest(stats1,[0 0 0 0 0 1 1 0]), which seems like it should be equivalent, I get p=.2075. When I run coefTest(stats1,[0 0 0 0 0 1 1 1]) I get .0072.

What input do I need to get the equivalent hypothesis test value from STATA? Or at a minimum where do I find the coefTest m-file so I can see what it's doing?

Thanks!

[Tom Lane]

I don't have your data. If you want to test whether two things are equal, this is the same as testing that their difference is 0, so try setting up the vector to specify a difference. Maybe [0 0 0 0 0 1 -1 0].

You can look at the implementation of the method for an object by typing something like

edit LinearModel/coefTest

[Alexis]

Hi Tom, can you explain how to set up the contrast matrix for 2- or 3-way interactions? This isn't covered in the documentation examples at all. Thanks!

Answer 3

[John]

Your suggestion nailed it. Thanks Tom!!