Packages used in this chapter
The packages used in this chapter include:
• lmPerm
• vegan
The following commands will install these packages if they are not already installed:
if(!require(lmPerm)){install.packages("lmPerm")}
if(!require(vegan)){install.packages("vegan")}
Permutation tests for linear models with the lmPerm package
The lmPerm package can be used to test coefficients from a linear model with a permutation test. By default, non-sequential sums of squares (that is, type 2 or type 3 sums of squares) are returned by the anova() call. Model p-value and r-squared can be reported with the summary() call. At the time of writing, I don’t know if there are methods for post-hoc tests. Residuals and predicted values can be asked for from the model object.
This example uses the Lung Capacity data set, but we’ll truncate the data set to the first 200 observations, just to get more interesting results.
LungCapData = read.csv(header=TRUE, stringsAsFactors = TRUE,
file="https://rcompanion.org/documents/LungCapData.csv")
LungCapData = LungCapData[1:200,]
library(lmPerm)
model = lmp(LungCap ~ Age + Smoke + Gender + Smoke:Gender, data=LungCapData)
anova(model)
Analysis of Variance Table
Response: LungCap
Df R Sum Sq R Mean Sq Iter Pr(Prob)
Age 1 968.59 968.59 5000 <2e-16 ***
Smoke 1 3.11 3.11 357 0.2213
Gender 1 17.51 17.51 5000 <2e-16 ***
Smoke:Gender 1 0.10 0.10 51 0.7059
Residuals 195 364.53 1.87
summary(model)
Residual standard error: 1.367 on 195 degrees of freedom
Multiple R-Squared: 0.7377, Adjusted R-squared: 0.7323
F-statistic: 137.1 on 4 and 195 DF, p-value: < 2.2e-16
Permutation tests for linear models with the vegan package
The vegan package is tailor-made for applications in ecology. The adonis2() function can test model parameters with permutation tests for models with either a univariate or multivariate left-hand side (dependent variable). The vegdist() function is used to create a distance matrix for the dependent variable(s). This is then used as the left-hand side of the model formula.
This example uses the Lung Capacity data set, but again we’ll truncate the data set to the first 200 observations to compare results to the lmPerm results.
Here, model2 uses the by="margin" option to ask for non-sequential sums of squares (that is, type 2 or type 3 sums of squares) for the model terms. Note that the p-value for Smoke:Gender is returned, but not for Smoke or Gender individually.
LungCapData = read.csv(header=TRUE, stringsAsFactors = TRUE,
file="https://rcompanion.org/documents/LungCapData.csv")
LungCapData = LungCapData[1:200,]
library(vegan)
Yvar = LungCapData$LungCap
LHS = vegdist(Yvar)
model2 = adonis2(LHS ~ Age + Smoke + Gender + Smoke:Gender, data=LungCapData,
permutations=1e4,
by="margin")
model2
Permutation test for adonis under reduced model
Marginal effects of terms
Df SumOfSqs R2 F Pr(>F)
Age 1 3.5762 0.51068 216.2090 9.999e-05 ***
Smoke:Gender 1 0.0030 0.00043 0.1817 0.7763
Residual 195 3.2254 0.46058
Total 199 7.0029 1.00000 00
The overall p-value and r-squared value for the overall model can be obtained by fitting a model with the by=NULL option.
model3 = adonis2(LHS ~ Age + Smoke + Gender + Smoke:Gender, data=LungCapData,
permutations=1e4,
by=NULL)
model3
Permutation test for adonis under reduced model
Df SumOfSqs R2 F Pr(>F)
Model 4 3.7775 0.53942 57.094 9.999e-05 ***
Residual 195 3.2254 0.46058
Total 199 7.0029 1.00000
Tests for sequentially-added terms will be reported with the by="term" option. Note that the p-value for all of Smoke, Gender, and Smoke:Gender is returned.
model4 = adonis2(LHS ~ Age + Smoke + Gender + Smoke:Gender, data=LungCapData,
permutations=1e4,
by="term")
model4
Permutation test for adonis under reduced model
Terms added sequentially (first to last)
Df SumOfSqs R2 F Pr(>F)
Age 1 3.5988 0.51390 217.5745 9.999e-05 ***
Smoke 1 0.0207 0.00296 1.2531 0.2629
Gender 1 0.1549 0.02212 9.3664 0.0004 ***
Smoke:Gender 1 0.0030 0.00043 0.1817 0.7797
Residual 195 3.2254 0.46058
Total 199 7.0029 1.00000