Main Content

compare

Class: GeneralizedLinearMixedModel

Compare generalized linear mixed-effects models

Description

example

results = compare(glme,altglme) returns the results of a likelihood ratio test that compares the generalized linear mixed-effects models glme and altglme. To conduct a valid likelihood ratio test, both models must use the same response vector in the fit, and glme must be nested in altglme. Always input the smaller model first, and the larger model second.

compare tests the following null and alternate hypotheses:

  • H0: Observed response vector is generated by glme.

  • H1: Observed response vector is generated by model altglme.

results = compare(glme,altglme,Name,Value) returns the results of a likelihood ratio test using additional options specified by one or more Name,Value pair arguments. For example, you can check if the first input model, glme, is nested in the second input model, altglme.

Input Arguments

expand all

Generalized linear mixed-effects model, specified as a GeneralizedLinearMixedModel object. For properties and methods of this object, see GeneralizedLinearMixedModel.

You can create a GeneralizedLinearMixedModel object by fitting a generalized linear mixed-effects model to your sample data using fitglme. To conduct a valid likelihood ratio test on two models that have response distributions other than normal, you must fit both models using the 'ApproximateLaplace' or 'Laplace' fit method. Models with response distributions other than normal that are fitted using 'MPL' or 'REMPL' cannot be compared using a likelihood ratio test.

Alternative generalized linear mixed-effects model, specified as a GeneralizedLinearMixedModel object. altglme be must fit to the same response vector as glme, but with different model specifications. glme must be nested in altglme, such that you can obtain glme from altglme by setting some of the model parameters of altglme to fixed values such as 0.

You can create a GeneralizedLinearMixedModel object by fitting a generalized linear mixed-effects model to your sample data using fitglme. To conduct a valid likelihood ratio test on two models that have response distributions other than normal, you must fit both models using the 'ApproximateLaplace' or 'Laplace' fit method. Models with response distributions other than normal that are fitted using 'MPL' or 'REMPL' cannot be compared using a likelihood ratio test.

Name-Value Arguments

Specify optional pairs of arguments as Name1=Value1,...,NameN=ValueN, where Name is the argument name and Value is the corresponding value. Name-value arguments must appear after other arguments, but the order of the pairs does not matter.

Before R2021a, use commas to separate each name and value, and enclose Name in quotes.

Indicator to check nesting between two models, specified as the comma-separated pair consisting of 'CheckNesting' and either true or false. If 'CheckNesting' is true, then compare checks if the smaller model glme is nested in the larger model altglme. If the nesting requirements are not satisfied, then compare returns an error. If 'CheckNesting' is false, then compare does not perform this check.

Example: 'CheckNesting',true

Output Arguments

expand all

Results of the likelihood ratio test, returned as a table with two rows. The first row is for glme, and the second row is for altglme. The columns of results contain the following.

Column NameDescription
ModelName of the model
DFDegrees of freedom
AICAkaike information criterion for the model
BICBayesian information criterion for the model
LogLikMaximized log likelihood for the model
LRStatLikelihood ratio test statistic for comparing altglme and glme
deltaDFDF for altglme minus DF for glme
pValuep-value for the likelihood ratio test

Examples

expand all

Load the sample data.

load mfr

This simulated data is from a manufacturing company that operates 50 factories across the world, with each factory running a batch process to create a finished product. The company wants to decrease the number of defects in each batch, so it developed a new manufacturing process. To test the effectiveness of the new process, the company selected 20 of its factories at random to participate in an experiment: Ten factories implemented the new process, while the other ten continued to run the old process. In each of the 20 factories, the company ran five batches (for a total of 100 batches) and recorded the following data:

  • Flag to indicate whether the batch used the new process (newprocess)

  • Processing time for each batch, in hours (time)

  • Temperature of the batch, in degrees Celsius (temp)

  • Categorical variable indicating the supplier of the chemical used in the batch (supplier)

  • Number of defects in the batch (defects)

The data also includes time_dev and temp_dev, which represent the absolute deviation of time and temperature, respectively, from the process standard of 3 hours at 20 degrees Celsius.

Fit a fixed-effects-only model using newprocess, time_dev, temp_dev, and supplier as fixed-effects predictors. Specify the response distribution as Poisson, the link function as log, and the fit method as Laplace. Specify the dummy variable encoding as 'effects', so the dummy variable coefficients sum to 0.

FEglme = fitglme(mfr,'defects ~ 1 + newprocess + time_dev + temp_dev + supplier','Distribution','Poisson','Link','log','FitMethod','Laplace','DummyVarCoding','effects');

Fit a second model that uses the same fixed-effects predictors, response distribution, link function, and fit method. This time, include a random-effects intercept grouped by factory, to account for quality differences that might exist due to factory-specific variations.

The number of defects can be modeled using a Poisson distribution

defectsijPoisson(μij)

This corresponds to the generalized linear mixed-effects model

log(μij)=β0+β1newprocessij+β2time_devij+β3temp_devij+β4supplier_Cij+β5supplier_Bij+bi,

where

  • defectsij is the number of defects observed in the batch produced by factory i during batch j.

  • μij is the mean number of defects corresponding to factory i (where i=1,2,...,20) during batch j (where j=1,2,...,5).

  • newprocessij, time_devij, and temp_devij are the measurements for each variable that correspond to factory i during batch j. For example, newprocessij indicates whether the batch produced by factory i during batch j used the new process.

  • supplier_Cij and supplier_Bij are dummy variables that use effects (sum-to-zero) coding to indicate whether company C or B, respectively, supplied the process chemicals for the batch produced by factory i during batch j.

  • biN(0,σb2) is a random-effects intercept for each factory i that accounts for factory-specific variation in quality.

glme = fitglme(mfr,'defects ~ 1 + newprocess + time_dev + temp_dev + supplier + (1|factory)','Distribution','Poisson','Link','log','FitMethod','Laplace','DummyVarCoding','effects');

Compare the two models using a theoretical likelihood ratio test. Specify 'CheckNesting' as true, so compare returns a warning if the nesting requirements are not satisfied.

results = compare(FEglme,glme,'CheckNesting',true)
results = 
    THEORETICAL LIKELIHOOD RATIO TEST

    Model     DF    AIC       BIC       LogLik     LRStat    deltaDF    pValue    
    FEglme    6     431.02    446.65    -209.51                                   
    glme      7     416.35    434.58    -201.17    16.672    1          4.4435e-05

Since compare did not return an error, the nesting requirements are satisfied. The small p-value indicates that compare rejects the null hypothesis that the observed response vector is generated by the model FEglme, and instead accepts the alternate model glme. The smaller AIC and BIC values for glme also support the conclusion that glme provides a better fitting model for the response.

More About

expand all