Statistical Model Analysis

When fitting models to data, it is often useful to analyze how well the model fits the data and how well the fitting meets the assumptions of the model. For a number of common statistical models, this is accomplished in Mathematica by way of fitting functions that construct FittedModel objects.

Object for fitted model information.

FittedModel objects can be evaluated at a point or queried for results and diagnostic information. Diagnostics vary somewhat across model types. Available model fitting functions fit linear, generalized linear, and nonlinear models.

Functions that generate FittedModel objects.

The major difference between model fitting functions such as LinearModelFit and functions such as Fit and FindFit is the ability to easily obtain diagnostic information from the FittedModel objects. The results are accessible without refitting the model.

Fitting options relevant to property computations can be passed to FittedModel objects to override defaults.

Typical data for these model-fitting functions takes the same form as data in other fitting functions such as Fit and FindFit.

Data specifications.

Linear Models

Linear models with assumed independent normally distributed errors are among the most common models for data. Models of this type can be fitted using the LinearModelFit function.

Linear model fitting.

Linear models have the form , where is the fitted or predicted value, the are parameters to be fitted, and the are functions of the predictor variables . The models are linear in the parameters . The can be any functions of the predictor variables. Quite often the are simply the predictor variables .

Options for model specification and for model analysis are available.

Options for LinearModelFit.

The Weights option specifies weight values for weighted linear regression. The NominalVariables option specifies which predictor variables should be treated as nominal or categorical. With NominalVariables->All, the model is an analysis of variance (ANOVA) model. With NominalVariables->{x₁, ..., x_i-1, x_i+1, ..., x_n} the model is an analysis of covariance (ANCOVA) model with all but the predictor treated as nominal. Nominal variables are represented by a collection of binary variables indicating equality and inequality to the observed nominal categorical values for the variable.

ConfidenceLevel, VarianceEstimatorFunction, and WorkingPrecision are relevant to the computation of results after the initial fitting. These options can be set within LinearModelFit to specify the default settings for results obtained from the FittedModel object. These options can also be set within an already constructed FittedModel object to override the option values originally given to LinearModelFit.

IncludeConstantBasis, LinearOffsetFunction, NominalVariables, and Weights are relevant only to the fitting. Setting these options within an already constructed FittedModel object will have no further impact on the result.

A major feature of the model-fitting framework is the ability to obtain results after the fitting. The full list of available results can be obtained using .

The properties include basic information about the data, fitted model, and numerous results and diagnostics.

Properties related to data and the fitted function.

The property gives the fitted parameter values . is the fitted function and gives the fitted function as a pure function. gives the list of functions , with being the constant 1 when a constant term is present in the model. The is the design or model matrix for the data. gives the list of the response or values from the original data.

Types of residuals.

Residuals give a measure of the pointwise difference between the fitted values and the original responses. gives the differences between the observed and fitted values . and are scaled forms of the residuals. The standardized residual is , where is the estimated error variance, is the diagonal element of the hat matrix, and is the weight for the data point. The studentized residual uses the same formula with replaced by , the variance estimate omitting the data point.

Properties related to the sum of squared errors.

gives a formatted analysis of variance table for the model. gives the numeric entries in the table and the remaining properties give the elements of columns in the table so individual parts of the table can easily be used in further computations.

Properties and diagnostics for parameter estimates.

gives the covariance between fitted parameters. The matrix is , where is the variance estimate, is the design matrix, and is the diagonal matrix of weights. is the associated correlation matrix for the parameter estimates. is equivalent to the square root of the diagonal elements of the covariance matrix.

and contain information about the individual parameter estimates, tests of parameter significance, and confidence intervals.

In[14]:=

The Estimate column of these tables is equivalent to . The -statistics are the estimates divided by the standard errors. Each -value is the two-sided -value for the -statistic and can be used to assess whether the parameter estimate is statistically significantly different from 0. Each confidence interval gives the upper and lower bounds for the parameter confidence interval at the level prescribed by the ConfidenceLevel option. The various and properties can be used to get the columns or the unformatted array of values from the table.

is used to measure the multicollinearity between basis functions. The inflation factor is equal to , where is the coefficient of variation from fitting the basis function to a linear function of the other basis functions. With IncludeConstantBasis->True, the first inflation factor is for the constant term.

gives the eigenvalues, condition indices, and variance partitions for the nonconstant basis functions. The Index column gives the square root of the ratios of the eigenvalues to the largest eigenvalue. The column for each basis function gives the proportion of variation in that basis function explained by the associated eigenvector. gives the values in the variance partitioning for all basis functions in the table.

Properties related to influence measures.

Pointwise measures of influence are often employed to assess whether individual data points have a large impact on the fitting. The hat matrix and catcher matrix play important roles in such diagnostics. The hat matrix is the matrix such that , where is the observed response vector and is the predicted response vector. gives the diagonal elements of the hat matrix. is the matrix such that , where is the fitted parameter vector.

gives the DFFITS values that provide a measure of influence of each data point on the fitted or predicted values. The DFFITS value is given by , where is the hat diagonal and is the studentized residual.

gives the DFBETAS values that provide measures of influence of each data point on the parameters in the model. For a model with parameters, the element of is a list of length with the value giving the measure of the influence of data point on the parameter in the model. The vector can be written as , where is the , element of the catcher matrix.

gives the Cook distance measures of leverage. The Cook distance is given by , where is the standardized residual.

The element of is given by and the value is equal to , where is the single deletion variance.

The Durbin-Watson -statistic is used for testing the existence of a first-order autoregressive process. The -statistic is equivalent to , where is the residual.

Properties of predicted values.

Tabular results for confidence intervals are given by and . These include the observed and predicted responses, standard error estimates, and confidence intervals for each point. Mean prediction confidence intervals are often referred to simply as confidence intervals and single prediction confidence intervals are often referred to as prediction intervals.

Mean prediction intervals give the confidence interval for the mean of the response at fixed values of the predictors and are given by , where is the quantile of the Student distribution with degrees of freedom, is the vector of basis functions evaluated at fixed predictors, and is the estimated covariance matrix for the parameters. Single prediction intervals provide the confidence interval for predicting at fixed values of the predictors, and are given by , where is the estimated error variance.

and give formulas for mean and single prediction confidence intervals as functions of the predictor variables.

Goodness-of-fit measures.

Goodness-of-fit measures are used to assess how well a model fits or to compare models. The coefficient of determination is the ratio of the model sum of squares to the total sum of squares. penalizes for the number of parameters in the model and is given by .

and are likelihood-based goodness-of-fit measures. Both are equal to times the log-likelihood for the model plus , where is the number of parameters to be estimated including the estimated variance. For is , and for is .

Generalized Linear Models

The linear model can be seen as a model with each response value being an observation from a normal distribution with mean value . The generalized linear model extends to models of the form , with each assumed to be an observation from a distribution of known exponential family form with mean , and being an invertible function over the support of the exponential family. Models of this sort can be obtained via GeneralizedLinearModelFit.

Generalized linear model fitting.

The invertible function is called the link function and the linear combination is referred to as the linear predictor. Common special cases include the linear regression model with the identity link function and Gaussian or normal exponential family distribution, logit and probit models for probabilities, Poisson models for count data, and gamma and inverse Gaussian models.

The error variance is a function of the prediction and is defined by the distribution up to a constant , which is referred to as the dispersion parameter. The error variance for a fitted value can be written as , where is an estimate of the dispersion parameter obtained from the observed and predicted response values, and is the variance function associated with the exponential family evaluated at the value .

Logit and probit models are common binomial models for probabilities. The link function for the logit model is and the link for the probit model is the inverse CDF for a standard normal distribution . Models of this type can be fitted via GeneralizedLinearModelFit with ExponentialFamily->"Binomial" and the appropriate LinkFunction or via LogitModelFit and ProbitModelFit.

Logit and probit model fitting.

Parameter estimates are obtained via iteratively reweighted least squares with weights obtained from the variance function of the assumed distribution. Options for GeneralizedLinearModelFit include options for iteration fitting such as PrecisionGoal, options for model specification such as LinkFunction, and options for further analysis such as ConfidenceLevel.

Options for GeneralizedLinearModelFit.

The options for LogitModelFit and ProbitModelFit are the same as for GeneralizedLinearModelFit except that ExponentialFamily and LinkFunction are defined by the logit or probit model and so are not options to LogitModelFit and ProbitModelFit.

ExponentialFamily can be , , , , , or . Binomial models are valid for responses from 0 to 1. Poisson models are valid for non-negative integer responses. Gaussian or normal models are valid for real responses. Gamma and inverse Gaussian models are valid for positive responses. Quasi-likelihood models define the distributional structure in terms of a variance function such that the log of the quasi-likelihood function for the data point is given by . The variance function for a model can be optionally set via ExponentialFamily->{"QuasiLikelihood", "VarianceFunction"->fun}, where fun is a pure function to be applied to fitted values.

DispersionEstimatorFunction defines a function for estimating the dispersion parameter . The estimate is analogous to in linear and nonlinear regression models.

ExponentialFamily, IncludeConstantBasis, LinearOffsetFunction, LinkFunction, NominalVariables, and Weights all define some aspect of the model structure and optimization criterion and can only be set within GeneralizedLinearModelFit. All other options can be set either within GeneralizedLinearModelFit or passed to the FittedModel object when obtaining results and diagnostics. Options set in evaluations of FittedModel objects take precedence over settings given to GeneralizedLinearModelFit at the time of the fitting.

Properties related to data and the fitted function.

gives the parameter estimates for the basis functions. gives the fitted function , and gives the linear combination . gives the list of functions , with being the constant 1 when a constant term is present in the model. is the design or model matrix for the basis functions.

Properties related to dispersion and model deviances.

Deviances and deviance tables generalize the model decomposition given by analysis of variance in linear models. The deviance for a single data point is , where is the log-likelihood function for the fitted model. gives a list of the deviance values for all data points. The sum of all deviances gives the model deviance. The model deviance can be decomposed as sums of squares, which are in an ANOVA table for linear models. The full model is the model whose predicted values are the same as the data.

In[31]:=

As with sums of squares, deviances are additive. The Deviance column of the table gives the increase in the model deviance when the given basis function is added. The Residual Deviance column gives the difference between the model deviance and the deviance for the submodel containing all previous terms in the table. For large samples, the increase in deviance is approximately distributed with degrees of freedom equal to that for the basis function in the table.

is the deviance for the null model, the constant model equal to the mean of all observed responses for models including a constant, or if a constant term is not included.

As with , a number of properties are included to extract the columns or unformatted array of entries from .

Types of residuals.

is the list of residuals, differences between the observed and predicted responses. Given the distributional assumptions, the magnitude of the residuals is expected to change as a function of the predicted response value. Various types of scaled residuals are employed in the analysis of generalized linear models.

If and are the deviance and residual for the data point, the deviance residual is given by . The Pearson residual is defined as , where is the variance function for the exponential family distribution. Standardized deviance residuals and standardized Pearson residuals include division by , where is the diagonal of the hat matrix. values combine deviance and Pearson residuals. The likelihood residual is given by .

provide a transformation of the residuals toward normality, so a plot of these residuals should be expected to look roughly like white noise. The Anscombe residual can be written as .

gives the residuals from the last step of the iterative fitting. The working residual can be obtained as evaluated at .

Properties and diagnostics for parameter estimates.

gives the covariance between fitted parameters and is very similar to the definition for linear models. With CovarianceEstimatorFunction->"ExpectedInformation", the expected information matrix obtained from the iterative fitting is used. The matrix is , where is the design matrix and is the diagonal matrix of weights from the final stage of the fitting. The weights include both weights specified via the Weights option and the weights associated with the distribution's variance function. With CovarianceEstimatorFunction->"ObservedInformation", the matrix is given by , where is the observed Fisher information matrix, which is the Hessian of the log-likelihood function with respect to parameters of the model.

is the associated correlation matrix for the parameter estimates. is equivalent to the square root of the diagonal elements of the covariance matrix. and contain information about the individual parameter estimates, tests of parameter significance, and confidence intervals. The test statistics for generalized linear models asymptotically follow normal distributions.

Properties related to influence measures.

and extend the leverage measures from linear regression to generalized linear models. The hat matrix from which the diagonal elements are extracted is defined using the final weights of the iterative fitting.

The Cook distance measures of leverage are defined as in linear regression with standardized residuals replaced by standardized Pearson residuals. The Cook distance is given by , where is the standardized Pearson residual.

Properties of predicted values.

Goodness-of-fit measures.

is the log-likelihood for the fitted model. and are penalized log-likelihood measures , where is the log-likelihood for the fitted model, is the number of parameters estimated including the dispersion parameter, and is for and for for a model of data points. is given by , where is the log-likelihood for the null model.

A number of the goodness-of-fit measures generalize from linear regression as either a measure of explained variation or as a likelihood-based measure. is given by . is a scaled version of Cox and Snell's measure . involves the ratio of log-likelihoods , and adjusts by penalizing for the number of parameters . uses the sum of squares interpretation of and is given as , where is the residual and is the mean of the responses .

is equal to , where the are Pearson residuals.

Nonlinear Models

A nonlinear least-squares model is an extension of the linear model where the model need not be a linear combination of basis function. The errors are still assumed to be independent and normally distributed. Models of this type can be fitted using the NonlinearModelFit function.

Nonlinear model fitting.

Nonlinear models have the form , where is the fitted or predicted value, the are parameters to be fitted, and the are predictor variables. As with any nonlinear optimization problem, a good choice of starting values for the parameters may be necessary. Starting values can be given using the same parameter specifications as for FindFit.

Options for model fitting and for model analysis are available.

Options for NonlinearModelFit.

General numeric options such as AccuracyGoal, Method, and WorkingPrecision are the same as for FindFit.

The Weights option specifies weight values for weighted nonlinear regression. The optimal fit is for a weighted sum of squared errors.

All other options can be relevant to computation of results after the initial fitting. They can be set within NonlinearModelFit for use in the fitting and to specify the default settings for results obtained from the FittedModel object. These options can also be set within an already constructed FittedModel object to override the option values originally given to NonlinearModelFit.

Properties related to data and the fitted function.

Basic properties of the data and fitted function for nonlinear models behave like the same properties for linear and generalized linear models with the exception that returns a rule as is done for the result of FindFit.

Many diagnostics for nonlinear models extend or generalize concepts from linear regression. These extensions often rely on linear approximations or large sample approximations.

Types of residuals.

As in linear regression, gives the differences between the observed and fitted values , and and are scaled forms of these differences.

The standardized residual is , where is the estimated error variance, is the diagonal element of the hat matrix, and is the weight for the data point, and the studentized residual is obtained by replacing with the single deletion variance . For nonlinear models a first-order approximation is used for the design matrix, which is needed to compute the hat matrix.

Properties related to the sum of squared errors.

provides a decomposition of the variation in the data attributable to the fitted function and to the errors or residuals.

The uncorrected total sums of squares gives the sum of squared responses, while the corrected total gives the sum of squared differences between the responses and their mean value.

Properties and diagnostics for parameter estimates.

gives the approximate covariance between fitted parameters. The matrix is , where is the variance estimate, is the design matrix for the linear approximation to the model, and is the diagonal matrix of weights. is the associated correlation matrix for the parameter estimates. is equivalent to the square root of the diagonal elements of the covariance matrix.

and contain information about the individual parameter estimates, tests of parameter significance, and confidence intervals obtained using the error estimates.

Curvature diagnostics.

The first-order approximation used for many diagnostics is equivalent to the model being linear in the parameters. If the parameter space near the parameter estimates is sufficiently flat, the linear approximations and any results that rely on first-order approximations can be deemed reasonable. Curvature diagnostics are used to assess whether the approximate linearity is reasonable. is a table of curvature diagnostics.

and are scaled measures of the normal and tangential curvatures of the parameter spaces at the best-fit parameter values. is a scaled measure of the radius of curvature of the parameter space at the best-fit parameter values. If the normal and tangential curvatures are small relative to the value of the , the linear approximation is considered reasonable. Some rules of thumb suggest comparing the values directly, while others suggest comparing with half the .

Properties related to influence measures.

The hat matrix is the matrix such that , where is the observed response vector and is the predicted response vector. gives the diagonal elements of the hat matrix. As with other properties, uses the design matrix for the linear approximation to the model.

The element of is equivalent to , where is the number of data points, is the number of parameters, is the hat diagonal, is the variance estimate for the full dataset, and is the residual.

Properties of predicted values.

Tabular results for confidence intervals are given by and . These results are analogous to those for linear models obtained via LinearModelFit, again with first-order approximations used for the design matrix.

and give functions of the predictor variables.

Goodness-of-fit measures.

, , , and are all direct extensions of the measures as defined for linear models. The coefficient of determination is , where is the residual sum of squares and is the uncorrected total sum of squares. The coefficient of determination does not have the same interpretation as the percentage of explained variation in nonlinear models as it does in linear models because the sum of squares for the model and for the residuals do not necessarily sum to the total sum of squares. penalizes for the number of parameters in the model and is given by .

and are equal to times the log-likelihood for the model plus , where is the number of parameters to be estimated including the estimated variance. For is , and for is .