parallel mediation with more than 2 mediators [closed]

i am an absolute noob when it comes to R programming but, alas, I have a college course I have to attend (and pass) and i am really, coming from social sciences, struggling. we're doing this assignment where we have to do a parallel mediation model/analysis with more than 2 mediators, dataset of choice (i have a dataset but in class we've only done it with one mediator)... any nice soul that could help out with an "for dummies" r script? i have this script from the lectures but struggling to replicate it with more mediators... I get the basic condition influences result, condition influences mediator than by doing so the result but not with two more mediators. much obliged to anyone willing to help. i suck at this, i know

## assign a name to that dataset as "sample.data2" sample.data2 <- read_spss("pmi.sav") sample.data2 <- data.table(sample.data2)

## regression predicting reported reactions model1 <- lm(reaction ~ cond, data = sample.data2) summary(model1)

model.M <- lm(pmi ~ cond, data = sample.data2) model.Y <- lm(reaction ~ cond + pmi, data = sample.data2) summary(model.M) summary(model.Y)

## sobel test ## model.M = regression model predicting mediator, with "apath" being the name of the IV ## model.Y = regression model predicting DV, with "bpath" being the name of the mediator sobel.test <- function(model.M, apath, model.Y, bpath) { a <- coef(model.M)[apath] sa <- summary(model.M)$$coef[apath, 2] b <- coef(model.Y)[bpath] sb <- summary(model.Y)$$coef[bpath, 2] indirect <- a*b

z <- indirect/sqrt(b^2*sa^2 + a^2*sb^2 + sa^2*sb^2)
pval <- 2*pnorm(-abs(z));

return(c(indirect = indirect, z = z, pval = pval))


}

sobel.test(model.M, "cond", model.Y, "pmi")

## 10,000 resample example ## function to get empirical indirect effect indirect.effect <- function(model.M, apath, model.Y, bpath) { a <- coef(model.M)[apath] b <- coef(model.Y)[bpath]

indirect <- a*b
indirect


}

indirect.effect(model.M, "cond", model.Y, "pmi")

rep.indirect <- sapply(1:10000, function(i) { resample <- sample.data2[sample(1:123, replace = T), ] model.M.resample <- lm(formula(model.M), data = resample) model.Y.resample <- lm(formula(model.Y), data = resample) indirect.effect(model.M.resample, "cond", model.Y.resample, "pmi") }, simplify = "array")

quantile(rep.indirect, c(0.025, 0.975))

## plot for 10,000 replication of indirect effects hist(rep.indirect, breaks = 100, col = "gray", main = "", xlab = "Indirect effects") abline(v = quantile(rep.indirect, 0.025), col = "red", lty = 2) ## red dotted (lty = 2) line, lower 95th percentile abline(v = quantile(rep.indirect, 0.975), col = "red", lty = 2) ## red dotted (lty = 2) line, upper 95th percentile

## ------------------------ ## ## indirect effect est in R ## ## ------------------------ ##

## function to get empirical indirect effect indirect.effect <- function(model.M, apath, model.Y, bpath) { a <- coef(model.M)[apath] b <- coef(model.Y)[bpath]

indirect <- a*b
indirect


}

indirect.effect(model.M, "cond", model.Y, "pmi")

## function to make inferecences for indirect effect boot.indirect <- function(data, i, model.M, apath, model.Y, bpath) {

## resample data n times
resample <- data[i,]

## refit regression model on resampled data
model.M.resample <- lm(formula(model.M), data = resample)
model.Y.resample <- lm(formula(model.Y), data = resample)

a <- coef(model.M.resample)[apath] ## a path
b <- coef(model.Y.resample)[bpath] ## b path
c <- coef(model.Y.resample)[apath] ## c path

## effect decomposition
## indirect effect
indirect.boot <- a*b
direct.boot <- c
return(c(direct.boot = direct.boot, indirect.boot = indirect.boot))


}

require(boot) boot_indirect_effect <- boot(data = sample.data2, statistic = boot.indirect, model.M = "model.M", apath = "cond", model.Y = "model.Y", bpath = "pmi", R = 5000, parallel = "multicore", ncpus = 8)

## inference for direct and indirect effect boot.ci(boot_indirect_effect, index = 1, type = "perc") boot.ci(boot_indirect_effect, index = 2, type = "perc")

## more complex model ## Model of M1 model.m1 <- lm(pmi ~ cond + gender + age, data = sample.data2) summary(model.m1)

## Model of M2 model.m2 <- lm(import ~ cond + gender + age, data = sample.data2) summary(model.m2)

## Model of Y model.Y <- lm(reaction ~ cond + pmi + import + gender + age, data = sample.data2) summary(model.Y)

## custom function for inference boot.indirect2 <- function(data, i) {

## resample data based on index i
resample.data <- data[i, ]

## re-estimate the regression models
model.m1.resample <- lm(formula(model.m1), data = resample.data)
model.m2.resample <- lm(formula(model.m2), data = resample.data)
model.Y.resample <- lm(formula(model.Y), data = resample.data)

## extract coefficients and compute indirect effects
a1 <- coef(model.m1.resample)["cond"]
a2 <- coef(model.m2.resample)["cond"]
b1 <- coef(model.Y.resample)["pmi"]
b2 <- coef(model.Y.resample)["import"]
c <- coef(model.Y.resample)["cond"]

indirect_pmi <- a1*b1
indirect_import <- a2*b2
direct <- c

## output vector, name - values
out <- c(indirect_pmi = indirect_pmi,
indirect_import = indirect_import,
total_indirect = indirect_pmi + indirect_import,
direct = direct)
return(out)


}

boot_indirect_effect2 <- boot(data = sample.data2, statistic = boot.indirect2, R = 10000, parallel = "multicore", ncpus = 8)

## inference for direct and indirect effect boot.ci(boot_indirect_effect2, index = 1, type = "perc") boot.ci(boot_indirect_effect2, index = 2, type = "perc") boot.ci(boot_indirect_effect2, index = 3, type = "perc")

## comparing two specific indirect effects boot_indirect_effect2$$t0[1] - boot_indirect_effect2$$t0[2] test <- boot_indirect_effect2$$t[,1] - boot_indirect_effect2$$t[,2] quantile(test, c(0.025, 0.975))

## Multicategorical variables ## read in "web.sav" data and ## assign a name to that dataset as "sample.data3" sample.data3 <- read_spss("web.sav")

multicat.model <- lm(attitude ~ d1 + d2, data = sample.data3) summary(multicat.model)

multicat.model.M <- lm(inter ~ d1 + d2, data = sample.data3) multicat.model.Y <- lm(attitude ~ d1 + d2 + inter, data = sample.data3) summary(multicat.model.M) summary(multicat.model.Y)

boot.indirect.multicat <- function(data, i) {

## resample data
resample.data <- data[i, ]

## re-estimate model of M and Y based on resampled data
multicat.M.rspl <- lm(formula(multicat.model.M), data = resample.data)
multicat.Y.rspl <- lm(formula(multicat.model.Y), data = resample.data)

## recover coefficients
a1 <- coef(multicat.M.rspl)["d1"]
a2 <- coef(multicat.M.rspl)["d2"]
b <- coef(multicat.Y.rspl)["inter"]
cp1 <- coef(multicat.Y.rspl)["d1"]
cp2 <- coef(multicat.Y.rspl)["d2"]

## effect decomposition
relative.ind.d1 <- a1*b
relative.ind.d2 <- a2*b
relative.direct.d1 <- cp1
relative.direct.d2 <- cp2

out <- c(relative.ind.d1 = relative.ind.d1,
relative.ind.d2 = relative.ind.d2,
relative.direct.d1 = relative.direct.d1,
relative.direct.d2 = relative.direct.d2)

return(out)


}

boot_indirect_effect3 <- boot(data = sample.data3, statistic = boot.indirect.multicat, R = 10000, parallel = "multicore", ncpus = 8)

## inference for direct and indirect effect boot_indirect_effect3$t0 boot.ci(boot_indirect_effect3, index = 1, type = "perc") boot.ci(boot_indirect_effect3, index = 2, type = "perc") boot.ci(boot_indirect_effect3, index = 3, type = "perc") boot.ci(boot_indirect_effect3, index = 4, type = "perc") ## omnibus test of total effect ## total effect model without no group dummies multicat.model.null <- lm(attitude ~ 1, data = sample.data3) summary(multicat.model.null) ## total effect model with group dummies (we already have it) summary(multicat.model) ## changes in R-sq (as a whole) summary(multicat.model)$$r.squared - summary(multicat.model.null)$$r.squared ## anova (whether R-sq change is significant?) anova(multicat.model, multicat.model.null)  ## omnibus test of direct effect ## direct effect model without no group dummies multicat.model.null2 <- lm(attitude ~ inter, data = sample.data3) summary(multicat.model.null2) ## direct effect model with group dummies (we already have it) summary(multicat.model.Y) ## changes in R-sq (as a whole) summary(multicat.model.Y)$$r.squared - summary(multicat.model.null2)$$r.squared ## anova (whether R-sq change is significant?) anova(multicat.model.null2, multicat.model.Y)  ## omnibus test of indirect effect boot.omnibus.multicat <- function(data, i) { ## resample data resample.data <- data[i, ] ## re-estimate model of M and Y based on resampled data multicat.M.null.rspl <- lm(inter ~ 1, data = resample.data) multicat.M.rspl <- lm(formula(multicat.model.M), data = resample.data) multicat.Y.rspl <- lm(formula(multicat.model.Y), data = resample.data) ## recover quantity of interest rsq.change <- summary(multicat.M.rspl)$$r.squared - summary(multicat.M.null.rspl)$$r.squared b <- coef(multicat.Y.rspl)["inter"] omnibus.indirect.multicat <- rsq.change*b return(omnibus.indirect.multicat)  } boot_omnibus.test <- boot(data = sample.data3, statistic = boot.omnibus.multicat, R = 10000, parallel = "multicore", ncpus = 8) ## inference for direct and indirect effect boot_omnibus.test$t0 boot.ci(boot_omnibus.test, index = 1, type = "perc")

closed as off-topic by Michael Chernick, Peter Flom♦Apr 24 at 11:11

This question appears to be off-topic. The users who voted to close gave this specific reason:

• "This question appears to be off-topic because EITHER it is not about statistics, machine learning, data analysis, data mining, or data visualization, OR it focuses on programming, debugging, or performing routine operations within a statistical computing platform. If the latter, you could try the support links we maintain." – Michael Chernick, Peter Flom
If this question can be reworded to fit the rules in the help center, please edit the question.