piecewiseSEM (version 1.1.1)

partial.resid: Calculates partial residuals for two variables

Description

Extracts partial residuals for y ~ x | Z, where Z represents all other variables present in a structured equation upon which x is conditional.

Usage

partial.resid(.formula = y ~ x, modelList, data, model.control = NULL, 
return.data.frame = TRUE, plotit = TRUE, plotreg = TRUE, plotCI = TRUE)

Arguments

.formula

a formula specifying the partial residuals to investigate.

modelList

a single model or list of regressions representing the structural equation model.

data

a data.frame used to construct the model.

model.control

a list of model control arguments to be passed to the partial residual models.

return.data.frame

whether a data.frame of the partial residuals should be returned. Default is TRUE.

plotit

whether the partial plot of y ~ x | Z should be returned. Default is TRUE.

plotreg

whether the partial regression of y ~ x | Z should also be plotted. Default is TRUE.

plotCI

whether the confidence intervals of the partial regression of y ~ x | Z should also be plotted. Default is TRUE.

Value

• Returns a data.frame where the first column are the partial residuals of y ~ x | Z, and the second column is the partial residuals of x ~ Z | y. If plotit = TRUE, then the function also returns a plot of y ~ x | Z. If plotreg = TRUE, then the plot includes the linear regression of y ~ x | Z.

References

Shipley, Bill. Cause and correlation in biology: a user's guide to path analysis, structural equations and causal inference. Cambridge University Press, 2002.

Examples

# Load model package
library(nlme)

# Load data from Shipley (2013)
data(shipley2013) 

shipley2013.modlist = list(
  
  lme(x2~x1, random = ~x1 | species, data = shipley2013),
  
  lme(x3~x2, random = ~x2 | species, data = shipley2013),
  
  lme(x4~x2, random = ~x2 | species, data = shipley2013),
  
  lme(x5~x3+x4, random = ~x3+x4 | species, data = shipley2013)
  
)

# Get partial residuals of x3 on x5 conditional on x4
resids.df = partial.resid(x5 ~ x3, shipley2013.modlist, shipley2013, 
  list(lmeControl(opt = "optim")))

# Also returns raw residuals values for plotting in other packages
head(resids.df)

## NOT RUN ##

# # Create example data
# set.seed(1)
# 
# example.data = data.frame(
#   y = rnorm(100, 0, 1),
#   x1 = rnorm(100, 10, 50),
#   random = letters[1:5]
# )
# 
# example.data$x2 = example.data$y + runif(100, 0, 5)
# 
# example.data$x3 = example.data$x2 + runif(100, 0, 5)
# 
# # Run regular linear model using lm()
# lm.model = lm(y ~ x1 + x2 + x3, example.data)
# 
# partial.resid(y ~ x3, lm.model, example.data)
# 
# # Works with interactions too
# lm.model2 = lm(y ~ x1 * x2 + x3, example.data)
# 
# partial.resid(y ~ x1 * x2, lm.model2, example.data, return.data.frame = FALSE)
# 
# # Run generalized least squared regression
# library(nlme)
# 
# gls.model = gls(y ~ x1 + x2 + x3, example.data)
# 
# partial.resid(y ~ x3, gls.model, example.data, return.data.frame = FALSE)
# 
# # Run mixed effects model using lme()
# lme.model = lme(y ~ x1 + x2 + x3, random = ~ 1 | random, example.data)
# 
# partial.resid(y ~ x3, lme.model, example.data, return.data.frame = FALSE)
# 
# # Run mixed effects model using lmer()
# library(lme4) 
# 
# lmer.model = lmer(y ~ x1 + x2 + x3 + (1 | random), example.data)
# 
# partial.resid(y ~ x3, lmer.model, example.data, return.data.frame = FALSE)
# 
# # Remove some values of x3 from data.frame
# example.data[c(12, 23, 45), "x3"] = NA
# 
# # Run regular linear model using lm()
# lm.model = lm(y ~ x1 + x2 + x3, example.data)
# 
# sum(!is.na(partial.resid(y ~ x3, lm.model, example.data)$x.resids)) # Should be 97

## END NOT RUN ##