Adds random noise to selected variables to evaluate collinearity. Also suitable for other models than linear regression and for categorical independent variables

1 2 3 4 5 6 7 8 |

`mod` |
A model object, not necessarily type lm |

`pvars` |
Contains an array of variables to be perturbed. Random values are added to the variable, after which the model is re-estimated. |

`prange` |
Contains an array of values determining the magnitude of perturbations. There should be as many |

`ptrans` |
Contains an array of transformations to be applied. Each element should contain valid |

`pfac` |
Contains a list of categorical variables to be perturbed and parameters controlling the reclassification process. The first component must be a factor name. The name for the first component is ignored. Other list components should correspond with options for |

`uniform` |
If uniform=TRUE, random values from a uniform distribution |

`niter` |
Indicates the number of times to re-estimate the model. Default is 100. |

`object` |
a |

`x` |
a |

`dec.places` |
number of decimal places to use when printing |

`full` |
if TRUE, some extra information is printed |

`...` |
arguments to be passed on to or from other methods. Print options for class |

Perturb is a tool for assessing the impact of small random changes (perturbations) to variables on parameter estimates. It is an alternative for collinearity diagnostics such as `vif`

in the car package, `vif`

in the rms package or `colldiag`

in this package. Perturb is particularly useful for evaluating collinearity if interactions are present or nonlinear transformations of variables, e.g. a squared term. Perturb can show how the perturbations affect the estimates of a variable and terms derived from it whereas other collinearity diagnostics simply treat interactions and transformations as regular independent variables. Perturb is not limited to linear regression but can be used for all regression-like models. Perturb can also deal with categorical variables by randomly misclassifying them to assess the impact on parameter estimates.

Perturb works by adding a small random “perturbation” value to selected independent variables, then re-estimating the model. This process is repeated *niter* times, after which a summary of the means, standard deviation, minimum and maximum of the parameter estimates is displayed. If collinearity is a serious problem in the data, then the estimates will be unstable and vary strongly.

Perturb can be used with categorical variables. Categorical variables are reclassified according to a table of reclassification probabilities. Random numbers determine to which category each case is reclassified at each iteration. The reclassification probabilities are specified to make reclassification to the same category highly likely. For ordered variables, short distance reclassification can be made more likely than long distance. See the section on *“Categorical variables”* and `reclassify`

for further details.

An object of class “perturb”. The main result is contained in `coef.table`

, which contains the parameter estimates for each iteration of the perturbation analysis. `summary`

prints the mean, standard deviation, minimum and maximum of `coef.table`

over the iterations. If the option `full`

is specified, `reclassify`

prints additional information on how the reclassification probabilities were derived. `Summary`

also prints information on the transformed model formula.

`coef.table` |
Estimated coefficients for each iteration of the perturbation analysis |

`formula` |
The model formula used |

`pvars` |
The continuous variables perturbed in the analysis |

`prange` |
Magnitude of the perturbations |

`ptrans2` |
The transformations using the temporary variables of the perturbation analysis |

`reclassify.tables` |
objects produced by |

`formula2` |
The model formula using temporary variables |

`distribution` |
“normal” or “uniform”, the distribution from which to draw random numbers |

In a perturbation analysis, categorical variables are reclassified with a high probability of remaining in the same category. This could be accomplished in several ways. `reclassify`

lets you specify the target percentage reclassification, then adjusts the reclassification frequencies so that the expected frequencies of the reclassified variable are the same as those of the original. In addition, reclassify imposes a meaningful pattern of association between the original and the reclassified variable. See `reclassify`

for details.

Categorical variables are specified in `perturb`

as a list in the `pfac`

option. The first (unnamed) component is the factor to be reclassified. The names of following components should be valid reclassify options follows by appropriate values. For example, to reclassify the factor “type” with a 95% probability of reclassify to the same category, use:

`p2<-perturb(m2,pvars=c("income","education"),prange=c(1,1),`

`pfac=list("type",pcnt=95))`

This command will iteratively re-estimate model `m2`

100 times (default). Random values from a normal distribution with a standard deviation of 1 will be added to the variables income and education. `Reclassify`

creates a table of initial reclassification probabilities for `type`

with a 95% probability of reclassification to the same category. This initial table is adjusted and the final reclassification probabilities printed in the output are subsequently used to reclassify `type`

at each iteration.

Use a list of lists to specify a model with more than one reclassification factor, for example:

`pfac=list(list("fegp6",pcnt=95),list("eyrc",pcnt=m1),list("expc",pcnt=m2))`

`q<-perturb(r1,pfac=pfac)`

Perturb can be used with estimation procedures other than `lm`

. On the other hand, collinearity is a result of extreme (multiple) correlation among independent variables. Another option is `vif`

in the rms package `vif`

in the car package, or `colldiag`

, which use only the independent variables of a regression model. This will usually be a faster solution since maximum likelihood procedures require iterative estimation for each iteration of the perturbation analysis. It is possible though that certain procedures are more sensitive to collinearity than `lm`

, in which case perturb could be a better solution.

John Hendrickx John_Hendrickx@yahoo.com

D. Belsley, E. Kuh, and R. Welsch (1980). Regression Diagnostics. Wiley.

Belsley, D.A. (1991). Conditioning diagnostics, collinearity and weak data in regression. New York: John Wiley & Sons.

Hendrickx, John, Ben Pelzer. (2004). Collinearity involving ordered and unordered categorical variables. Paper presented at the RC33 conference in Amsterdam, August 17-20 2004. Available at http://home.wanadoo.nl/john.hendrickx/statres/perturb/perturb.pdf

`reclassify`

, `colldiag`

, `[car]`

`vif`

, `[rms]`

`vif`

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 | ```
library(car)
data(Duncan)
attach(Duncan)
m1<-lm(prestige~income+education)
summary(m1)
anova(m1)
vif(m1)
p1<-perturb(m1,pvars=c("income","education"),prange=c(1,1))
summary(p1)
m2<-lm(prestige~income+education+type)
summary(m2)
anova(m2)
vif(m2)
p2<-perturb(m2,pvars=c("income","education"),prange=c(1,1),pfac=list("type",pcnt=95))
summary(p2)
## Not run:
r1<-lm(ses~fegp6+educyr+eyr+exp2)
summary(r1)
q<-perturb(r1,c("eyr","exp"),c(2.5,2.5),ptrans="exp2<-exp^2")
summary(q)
fegp6<-as.factor(fegp6)
# eyr and exp also as factors
eyrc<-cut(eyr,c(min(eyr),40,50,60,70,80,max(eyr)),include.lowest=T,right=F)
table(eyrc)
expc<-cut(exp,c(0,10,20,max(exp)),include.lowest=T,right=F)
table(expc)
# rough initial reclassification probabilities,
# program will ensure they sum to 100 row-wise
m1<-matrix(0,nlevels(eyrc),nlevels(eyrc))
m1[row(m1)==col(m1)]<-80
m1[abs(row(m1)-col(m1))==1]<-8
m1[abs(row(m1)-col(m1))==2]<-2
m1
m2<-matrix(0,nlevels(expc),nlevels(expc))
m2[row(m2)==col(m2)]<-80
m2[abs(row(m2)-col(m2))==1]<-10
m2[abs(row(m2)-col(m2))==2]<-2
m2
r2<-lm(ses~fegp6+eyrc+expc)
summary(r2)
pfac=list(list("fegp6",pcnt=95),list("eyrc",pcnt=m1),list("expc",pcnt=m2))
q2<-perturb(r2,pfac=pfac,niter=1)
summary(q2)
## End(Not run)
``` |

Questions? Problems? Suggestions? Tweet to @rdrrHQ or email at ian@mutexlabs.com.

All documentation is copyright its authors; we didn't write any of that.