naginterfaces.library.blgm.lm_design_matrix¶
- naginterfaces.library.blgm.lm_design_matrix(hform, hddesc, dat, hxdesc)[source]¶
lm_design_matrix
generates a design matrix from a data matrix and model description.Note: this function uses optional algorithmic parameters, see also:
optset()
,optget()
.For full information please refer to the NAG Library document for g22yc
https://www.nag.com/numeric/nl/nagdoc_29.2/flhtml/g22/g22ycf.html
- Parameters
- hformHandle
A G22 handle to the internal data structure containing a description of the model as returned in by
lm_formula()
.- hddescHandle
A G22 handle to the internal data structure containing a description of the data matrix, as returned in by
lm_describe_data()
.- datfloat, array-like, shape
The data matrix, . By default , the th value for the th variable, for , for , should be supplied in .
If the option ‘Storage Order’, described in
lm_describe_data()
, is set to ‘VAROBS’, should be supplied in .- hxdescHandle, modified in place
On entry: must be set to a null Handle, alternatively an existing G22 handle may be supplied in which case this function will destroy the supplied G22 handle as if
handle_free()
had been called.On exit: holds a G22 handle to the internal data structure containing a description of the design matrix, . You must not change the G22 handle other than through the functions in submodule
blgm
.
- Returns
- xfloat, ndarray, shape
The design matrix, . By default , the th value for the th column, for , for , is returned in .
If the option ‘Storage Order’, described in
lm_formula()
, is set to ‘VAROBS’, is returned in .
- Other Parameters
- ‘Formula’str
This option returns a verbose formula string describing the model, , used to create the design matrix. This formula will only contain variable names, the operators ‘’ and ‘’ and any contrast identifiers present.
- ‘Min Number of Columns’int
This option returns the minimum number of columns required to hold the design matrix, . In most cases . The one exception is when = 71, that is the size of was too small but the data matrix given in can be used as the design matrix. In this case, and holds the number of columns that would be required if only the relevant parts of were copied into a new array.
- ‘Number of Columns’int
This option returns , the number of columns in the design matrix.
- ‘Number of Observations’int
This option returns , the number of observations in the design matrix.
- ‘Storage Order’str
This option returns how the design matrix, , is stored in .
If , , the value for the th variable of the th observation of the design matrix is stored in .
If , , the value for the th variable of the th observation of the design matrix is stored in .
It should be noted that ‘Storage Order’ is not writeable. If you wish to change the storage order of the design matrix you need to change ‘Storage Order’ in as described in Other Parameters for lm_formula prior to calling
lm_design_matrix
.
- Raises
- NagValueError
- (errno )
has not been initialized or is corrupt.
- (errno )
is not a G22 handle as generated by
lm_formula()
.- (errno )
A variable name used when creating is not present in .
Variable name: .
- (errno )
has not been initialized or is corrupt.
- (errno )
is not a G22 handle as generated by
lm_describe_data()
.- (errno )
On entry, column of the data matrix, , is not consistent with information supplied in , .
- (errno )
On entry, and .
Constraint: .
- (errno )
On entry, and .
Constraint: .
- (errno )
On entry, and .
Constraint: .
- (errno )
On entry, and .
Constraint: .
- (errno )
On entry, is not a null Handle or a recognised G22 handle.
- Warns
- NagAlgorithmicWarning
- (errno )
The model contains categorical variables, but no intercept or main effects terms have been requested.
Please check the design matrix returned matches the model you require.
- (errno )
Column of the data matrix, , required rounding more than expected when being treated as a categorical variable, .
- Notes
lm_design_matrix
generates a design matrix from a data matrix and a model description. Design matrices encapsulate the observed values of the independent variables and the required model in a form that can be used by many of the model fitting functions available in the NAG Library, for example those in submodulecorreg
.Notation
Let denote a data matrix with observations on independent variables, denoted by , for . If is a categorical variable, let denote the number of levels associated with it. If is a binary, ordinal or continuous variable, let .
Let denote the th value of .
Let denote a model made up of one or more terms, denoted by . Each term consists of either a main effect or an interaction and hence can be described using one or more variable names and the interaction operator ‘’. The operator ‘’ is used to denote the addition of a term to the model. Therefore, denotes a model with three terms, the first two terms being the main effects for variables and and the last term the interaction between them. For simplicity we reorder the terms of the model by the number of variables in them, so main effects come first, then two-way interactions, then three-way interactions etc. By default it is assumed that the model contains a mean effect (or intercept term), if the mean effect is excluded, this will be denoted by ‘’, so is a model with one term and a mean effect and is the same model with the mean effect dropped.
lm_design_matrix
generates an design matrix, , from and .Dummy Variables
When constructing a design matrix, we cannot work directly with categorical variables. Categorical variables must first be recoded into dummy variables. A categorical variable requires dummy variables. Let denote an matrix of dummy variables for defined as
where is the th column of and is the th element of .
For a binary, ordinal or continuous variable, .
Full Design Matrix
Given a model, , and the matrices of dummy variables constructing the full design matrix is trivial. Each term is processed in order and
If term is a main effect, that is for some , is copied into .
If term is a two-way interaction, that is , for some , then
Loop over .
Loop over .
Add a column to corresponding to the element-wise product of and .
Higher interaction terms are handled in a similar manner as the two-way interactions by adding columns constructed from multiplying all combinations of the columns of the corresponding s that correspond to the variables involved. In all cases, the variables towards the right hand side of a term are iterated over the quickest.
Contrasts
Using the full design matrix in an analysis can result in an overparameterized model. This is due to often not being of full rank as the sum of all the dummy variables for a particular variable is a vector of ones. This source of overparameterization can be alleviated by using a design matrix where (some) dummy variables are replaced by contrasts. For a categorical variable the contrasts are a set of functionally independent linear combinations of the dummy variables.
Whilst the choice of contrasts used in term will affect the individual model coefficients (parameters), it has no effect on the overall contribution of .
For a given variable , the contrasts can be represented by an matrix, . The rows of correspond to a particular value of and the columns correspond to the values to use in the design matrix.
Six types of contrast are available in
lm_design_matrix
; two types of treatment contrasts, two types of sum contrasts, Helmert contrasts and polynomial contrasts. Unless specified otherwise, the contrasts used bylm_design_matrix
are treatment contrasts relative to the first level. See the description of the option ‘Contrast’ inlm_formula()
for ways of changing the contrasts used.Treatment Contrasts
Treatment contrasts are taken relative to either the first or last level of the variable. For example, if ,
would be the contrast matrix for using treatment contrasts relative to the first level. The contrast matrix obtained when using treatment contrasts relative to the last level is similar, but the row of zeros appears at the bottom and all other rows are shifted up one.
Strictly speaking, the term contrast implies that each row in the contrast matrix sums to zero. That is not the case for treatment contrasts, however they are included as this coding is commonly used in practice.
Sum Contrasts
Sum contrasts are similar to treatment contrasts and again can be taken relative to the first or last level of the variable. Unlike treatment contrasts, sum contrasts effectively constrain the coefficients related to the variable to sum to zero. For example, if ,
would be the contrast matrix for using treatment contrasts relative to the last level. The contrast matrix obtained when using treatment contrasts relative to the first level is similar, but the row of s appears at the top and all other rows are shifted down one.
Helmert Contrasts
With Helmert contrasts level of the variable is compared with the average effect of all previous levels. For example, if ,
would be the contrast matrix for using Helmert contrasts.
Polynomial Contrasts
With polynomial contrasts the entries in the columns of correspond in linear, quadratic, cubic, quartic, etc. terms to a hypothetical underlying numeric variable that takes equally spaced values at each level. For example, if ,
would be the contrast matrix for using polynomial contrasts.
When Contrasts Can Be Used
Depending on the specifics of the model, , it may not be possible to always replace the dummy variables with contrasts for all variables in all terms and retain the same model. A simple example of this is a data matrix, , with four observations and two variables which have two and three levels respectively. This data matrix might look something like:
For the sake of argument, assume that our model contains the main effect for each variable, but does not contain a mean effect (or intercept term). So using the notation established earlier, . The full design matrix, , for this data matrix and model would be
However, is not of full rank (and hence is overparameterized) because the sum of the first two columns is a vector of ones as is the sum of the last three columns.
In order to alleviate this we might try constructing where the dummy variables have been replaced by contrasts. Assuming treatment contrasts, relative to the first level, we would have
However, using makes an implicit assumption that the expected value of the dependent variable (the quantity being modelled) is zero when and . This assumption was not made when we used and hence the two design matrices are not equivalent. One solution would be to use dummy variables for and contrasts for , which would result in a design matrix, of
Using would give an equivalent model to using .
The algorithm used by
lm_design_matrix
to decide which variables, in which terms, can be coded as contrasts and which need to be coded as dummy variables is described below.Suppose is any variable that appears in term , let denote the term obtained by dropping from . For example, if , . In this context, the empty term is taken to be the mean effect (or intercept term). We say that appears in if there exists a term , , that contains all of the variables appearing in . In most cases , but this is not required. Note, as stated earlier, the terms in are ordered by the number of variables in them.
A variable, in term is coded by contrasts if appears in and by dummy variables otherwise. It is, therefore, possible for variable to be coded by contrasts in some terms and dummy variables in others within the same .
The above rule assumes the presence of a mean effect. If no such effect is present in the model, the main effect of the first categorical variable is coded by dummy variables to compensate. If no main effects appear in the model, the warning = 14 is returned.
A longer description and informal proof that the resulting is a suitable design matrix for the model of interest can be found in module two of Chambers and Hastie (1992).
Mean Effect
The mean effect (or intercept term) is included in a design matrix by adding a column of ones as the first column of . However, many model fitting functions in the NAG Library handle the mean effect as a special case and do not require it to be explicitly added to the design matrix. Therefore, by default,
lm_design_matrix
does not explicitly add the mean effect to the design matrix. This behaviour can be changed via the option ‘Explicit Mean’ inlm_formula()
.
- References
Chambers, J M and Hastie, T J, 1992, Statistical Models in S, Wadsworth and Brooks/Cole Computer Science Series