Skip to content

Commit

Permalink
Revert "Add tidySEM to lesson 1"
Browse files Browse the repository at this point in the history
This reverts commit b3f7a76.
  • Loading branch information
mattansb committed Feb 28, 2021
1 parent b3f7a76 commit a7d64d4
Show file tree
Hide file tree
Showing 6 changed files with 291 additions and 213 deletions.
129 changes: 0 additions & 129 deletions 01 path analysis/Appendix - advanced plotting.R

This file was deleted.

74 changes: 74 additions & 0 deletions 01 path analysis/Exercise 1 - solution.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,74 @@
# Exercise ----------------------------------------------------------------

income_psych <- read.csv("income_psych.csv")
library(lavaan)

## 1. Look at the plot in "plot2.png". Fit this model with `lavaan`.

model <- "
# We need to build a regression equation for each endogenous variable:
anxiety ~ shyness
mood_neg ~ anxiety
income ~ anxiety + mood_neg + shyness
"

fit <- sem(model, data = income_psych)
summary(fit, standardize = TRUE)


library(semPlot)
semPaths(fit)
# See `advanced plotting` for making "plot2.png"






## 2. Compute all the paths in this model from *anxiety* to *income*, and
## the total of these paths.
model <- "
anxiety ~ a * shyness
mood_neg ~ b * anxiety
income ~ c * anxiety + d * mood_neg + e * shyness
# Define paths:
direct := c
ind_shyness := a * e
ind_mood := b * d
total := direct + ind_shyness + ind_mood
"

fit <- sem(model, data = income_psych)
summary(fit, standardize = TRUE)
# Seems like the total path is not significant.
# Looks like a case of supression... Some paths are +, some are -...
# They LOOK like they cancel each other out!

# However...
# The total is significant in the standerdized solution!! *
standardizedsolution(fit, output = "text")

## - Why is the std total not equal exactly to the real correlation?
cor(income_psych$anxiety, income_psych$income)
# Theres a diff of (-0.001). Because the model isn't *just-identified*.
# We are missing the arrow from `shyness` to `mood_neg`.


## - is it very different?
# The diff is of (-0.001). Is this a lot? Is it significant?
# We will see next time how to answer these questions!


# * -----------------------------------------------------------------------

# Why does this happen?
# In SEM, standard errors (SE) are NOT invariant - that is, data
# transformation do not affect SEs the same way they do estimates (the
# parameters).
# However, Chis-square atatistics ARE invariant - so this is not an issue
# in model comparison (which we will learn next week).
#
# You can read more (mathy), here:
# https://doi.org/10.1007/BF01065882
# (Thanks to Michael Neale for his explanation)
184 changes: 184 additions & 0 deletions 01 path analysis/advanced plotting.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,184 @@
library(lavaan)
library(tidySEM)
library(dplyr)

income_psych <- read.csv("income_psych.csv")

mediation_model <- '
income ~ anxiety + mood_neg
mood_neg ~ anxiety
'
fit <- sem(mediation_model, data = income_psych)

## See below for how "plot1.png" and "plot2.png" were made.

## Based on:
## https://cjvanlissa.github.io/tidySEM/articles/Plotting_graphs.html


# Changing the layout -----------------------------------------------------

lay <- get_layout(NA, "mood_neg", NA,
NA, NA, NA,
"anxiety", NA, "income",
rows = 3)




graph_data <- prepare_graph(fit, layout = lay, label = "est_std")

plot(graph_data)




# Adding var labels -------------------------------------------------------

# We can add labels to our variables (our "nodes") by chaging the "label" column
# in the "nodes" part of the plot:

(N <- nodes(graph_data))

nodes(graph_data) <- N %>%
mutate(label = c("Anxiety", "Income", "Negative Mood"))

plot(graph_data)



# Adding / changing arrow labels ------------------------------------------


(tab <- table_results(fit, columns = c("label", "est_std", "confint_std", "pval"), digits = 2))

(E <- edges(graph_data))

edges(graph_data) <- E %>%
mutate(est_std = tab$est_std,
label = paste0(round(est_std, 2), " ",tab$confint))

plot(graph_data)





# Even more ---------------------------------------------------------------

E <- edges(graph_data)

edges(graph_data) <- E %>%
mutate(
# color lines by sign:
colour = case_when(
arrow == "both" ~ "black",
est_std < 0 ~ "red",
est_std > 0 ~ "green"
),
# mark which is sig
linetype = c("dashed","dashed","solid", "dashed", "solid", "solid")
)


N <- nodes(graph_data)

nodes(graph_data) <- N %>%
mutate(size = c(1,3,1))

plot(graph_data)





# Making "plot1" ----------------------------------------------------------

mediation_model <- '
income ~ anxiety + mood_neg
mood_neg ~ anxiety
'

fit <- sem(mediation_model, data = income_psych)



lay <- get_layout(NA, "mood_neg", NA,
"anxiety", NA, "income",
rows = 2)

graph_data <- prepare_graph(fit, layout = lay)

E <- edges(graph_data)
E$label <- ""
E$connect_from[3] <- "top"
E$connect_to[3] <- "left"
E$connect_from[2] <- "right"
E <- E[1:3,]
edges(graph_data) <- E



(N <- nodes(graph_data))
N$label <- c("Anxiety", "Income", "Negative Mood")
nodes(graph_data) <- N





(p1 <- plot(graph_data))

# because the result it a ggplot, we will use ggsave to save it
ggplot2::ggsave(p1, filename = "plot1.png", width = 6, height = 3)





# Making "plot2" ----------------------------------------------------------


mediation_model <- '
income ~ anxiety + mood_neg + shyness
mood_neg ~ anxiety
anxiety ~ shyness
'

fit <- sem(mediation_model, data = income_psych)






lay <- get_layout(NA, NA, "mood_neg", NA,
NA, "anxiety", NA, "income",
"shyness", NA, NA, NA,
rows = 3)

graph_data <- prepare_graph(fit, layout = lay)

E <- edges(graph_data)
E$label <- ""
E$connect_from[5] <- "top"
E$connect_to[5] <- "left"
E$connect_from[4] <- "top"
E$connect_to[4] <- "left"
E$connect_to[3] <- "bottom"
E <- E[1:5,]
edges(graph_data) <- E

plot(graph_data)

(N <- nodes(graph_data))
N$label <- c("Anxiety", "Income", "Negative Mood", "Shyness")
N$node_xmin[3] <- N$node_xmin[3]-0.5
N$node_xmax[3] <- N$node_xmax[3]+0.5
nodes(graph_data) <- N


(p2 <- plot(graph_data))

# because the result it a ggplot, we will use ggsave to save it
ggplot2::ggsave(p2, filename = "plot2.png", width = 6, height = 3)
Loading

0 comments on commit a7d64d4

Please sign in to comment.