```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE, warning = FALSE, message = FALSE)
set.seed(84322)
```

## Introduction

Chapter 15 discussed mediation analysis, a way of understanding the process of an effect from one variable to another. Mediation is built on linear regression as we've discussed but combines paths to create a single conceptual model.

1. Let's start by loading the `tidyverse` package, the `furniture` package, the `lavaan` and the `MarginalMediation` package.
```{r}
library(tidyverse)
library(furniture)
library(lavaan)
library(MarginalMediation)
```
2. Import the `FacialBurns` data set found in `lavaan` into R.
```{r}
data(FacialBurns)
str(FacialBurns)
```

## Mediation

3. We hypothesize that anxiety (`HADS`) mediates the effect of `TBSA` on `Selfesteem`. Let's test this out first with `lavaan`. We will create the model object with two regression models and then we define the indirect effect and the total effect from the coefficients of the two regressions.
```{r}
med_model <- "
   ## a path
   HADS ~ a*TBSA
   ## b and c' paths
   Selfesteem ~ b*HADS + c*TBSA

   ## Define Indirect and Total Effects
   ind := a * b
   tot := a * b + c
"
sem(med_model, data = FacialBurns, se = "bootstrap", bootstrap = 500) %>% 
  summary()
```
4. The output gives us the a, b, and c paths (as we defined them) along with the defined parameters `ind` and `tot`. Is the indirect path (from TBSA -> HADS -> Self-esteem) significant? What does it mean here?
5. Let's get the standardized indirect and total effects.
```{r}
sem(med_model, data = FacialBurns) %>% 
  summary(standardized = TRUE)
```
6. Interpet the standardized indirect effect. How much of the total effect is explained by the process through anxiety?
7. We can get similar results using the `MarginalMediation` package that uses linear regression, although this provides us with confidence intervals instead of p-values. Notably, the standardized results will be slightly different given that `lavaan` standardizes all the outcomes and predictors whereas `MarginalMediation` only standardizes the outcome (since there can also be categorical predictors in many mediation models).
```{r}
bcpath <- glm(Selfesteem ~ HADS + TBSA, data = FacialBurns)
apath <- glm(HADS ~ TBSA, data = FacialBurns)

mma(bcpath, apath,
    ind_effects = c("TBSA-HADS"))
```
8. Because each part of the mediation is a linear model, we can assess the assumptions as we've done before.
```{r}
par(mfrow = c(2,2))
plot(bcpath)
educ7610::diagnostics(bcpath) %>% head(20)
```
9. Do you see any issues or extreme values from the figures? What about the diagnostic statistics (for space, we are only showing the first 20 observations)?

## Mediation for Categorical Outcomes (Marginal Mediation Analysis)

10. A new approach to when the the mediator and/or outcome is categorical is called Marginal Mediation Analysis. To use this method, we will use the function `mma()` from `MarginalMediation`. To use it, we fit two separate models using `glm()`. But first, we are going to dichotomize the anxiety variable (just to demonstrate the method--don't do this in practice).
```{r}
FacialBurns <- FacialBurns %>%
  mutate(HADS_binary = cut_number(HADS, 2))

bcpath <- glm(Selfesteem ~ HADS_binary + TBSA, data = FacialBurns)
apath <- glm(HADS_binary ~ TBSA, data = FacialBurns, family = binomial(link = "logit"))

mma(bcpath, apath,
    ind_effects = c("TBSA-HADS_binary"))
```
11. Both the unstandardized and Standardized effects are reported. Interpret the standardized indirect effect here.


## Moderated-Mediation (Conditional Process Analysis)

12. There are times that a combination of moderation and mediation are of interest. This can be a challenge but it is possible. Let's use `lavaan` to see if `sex` moderates the effect from anxiety to self-esteem. 
```{r}
FacialBurns <- FacialBurns %>%
  mutate(HADS_Sex = HADS * Sex)
med_model <- "
   ## a path
   HADS ~ a*TBSA
   ## b and c' paths
   Selfesteem ~ b*HADS + c*TBSA + b2 * HADS_Sex

   ## Define Indirect effects for each sex
   ind_sex1 := a * b
   ind_sex2 := a * b + a * b2
"
sem(med_model, data = FacialBurns, se = "bootstrap", bootstrap = 500) %>% 
  summary()
```
13. Did sex moderate the b path (see the p-value associated with the `b2` path)?
14. For more information on moderated mediation see [nickmichalak.com](http://nickmichalak.com/blog_entries/2018/nrg02/nrg02.html)

## Conclusion

Mediation analysis is an important part of our field. It helps us understand how one variable affects another. Because each part of a mediation is essentially a linear regression (or a GLM), the things we have learned throughout the course apply to each piece of the entire conceptual model.

