# January 17

library(tidyverse)

d <- read.csv("http://whlevine.hosted.uark.edu/psyc5133/tv.csv")

d <- d %>%
	mutate(tv.c = tv - mean(tv),
				 ability.c = ability - mean(ability))

# non-moderated model

m1 <- lm(achieve ~ tv.c + ability.c, d)
summary(m1)

# moderated model

m2 <- lm(achieve ~ tv.c * ability.c, d)
summary(m2)

library(interactions)

interact_plot(model = m2, pred = tv.c, modx = ability.c)

# TV no longer has *a* slope; it has many!

# let's center ability at different values and examine the slope of TV
d <- d %>% 
	mutate(ability.low = ability - (mean(ability) - sd(ability)),
				 ability.high = ability - (mean(ability) + sd(ability)))

m2.low.ability  <- lm(achieve ~ tv.c * ability.low, d)
m2.M.ability    <- lm(achieve ~ tv.c * ability.c, d)
m2.high.ability <- lm(achieve ~ tv.c * ability.high, d)

# slopes
coef(m2.low.ability)["tv.c"]   #  0.74
coef(m2.M.ability)["tv.c"]     # -0.36
coef(m2.high.ability)["tv.c"]  # -1.45

# relationship between TV-watching and achievement gets increasingly negative as
# ability gets higher

# regions of significance?

# it's useful to fit a model with uncentered predictors
m2uncentered <- lm(achieve ~ tv * ability, d)
johnson_neyman(model = m2uncentered, pred = tv, modx = ability)

# the slope of TV is significantly positive (weird) for low-ability cases, but
# significantly negative for high-ability cases