**Usage:** This tutorial accompanies [Discovering Statistics Using R and RStudio](https://www.discovr.rocks/) [@field_discovering_2023] by [Andy Field](https://en.wikipedia.org/wiki/Andy_Field_(academic)). It contains material from the book so there are some copyright considerations but I offer them under a [Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License](http://creativecommons.org/licenses/by-nc-nd/4.0/). Tl;dr: you can use this tutorial for teaching and non-profit activities but please don't meddle with it or claim it as your own work.
r cat_space(fill = blu) Welcome to the discovr space pirate academy
Hi, welcome to discovr space pirate academy. Well done on embarking on this brave mission to planet r rproj()s, which is a bit like Mars, but a less red and more hostile environment. That's right, more hostile than a planet without water. Fear not though, the fact you are here means that you can master r rproj(), and before you know it you'll be as brilliant as our pirate leader Mae Jemstone (she's the badass with the gun). I am the space cat-det, and I will pop up to offer you tips along your journey.
On your way you will face many challenges, but follow Mae's system to keep yourself on track:
r bmu(height = 1.5) This icon flags materials for teleporters. That's what we like to call the new cat-dets, you know, the ones who have just teleported into the academy. This material is the core knowledge that everyone arriving at space academy must learn and practice. For accessibility, these sections will also be labelled with [(1)]{.alt}.
r user_visor(height = 1.5) Once you have been at space pirate academy for a while, you get your own funky visor. It has various modes. My favourite is the one that allows you to see everything as a large plate of tuna. More important, sections marked for cat-dets with visors goes beyond the core material but is still important and should be studied by all cat-dets. However, try not to be disheartened if you find it difficult. For accessibility, these sections will also be labelled with [(2)]{.alt}.
r user_astronaut(height = 1.5) Those almost as brilliant as Mae (because no-one is quite as brilliant as her) get their own space suits so that they can go on space pirate adventures. They get to shout RRRRRR really loudly too. Actually, everyone here gets to should RRRRRR really loudly. Try it now. Go on. It feels good. Anyway, this material is the most advanced and you can consider it optional unless you are a postgraduate cat-det. For accessibility, these sections will also be labelled with [(3)]{.alt}.
It's not just me that's here to help though, you will meet other characters along the way:
r alien(height = 1.5) aliens love dropping down onto the planet and probing humanoids. Unfortunately you'll find them probing you quite a lot with little coding challenges. Helps is at hand though.
r robot(height = 1.5)bend-R is our coding robot. She will help you to try out bits of r rproj() by writing the code for you before you encounter each coding challenge.
r bug(height = 1.5) we also have our friendly alien bugs that will, erm, help you to avoid bugs in your code by highlighting common mistakes that even Mae Jemstone sometimes makes (but don't tell her I said that or my tuna supply will end).
Also, use hints and solutions to guide you through the exercises (Figure 1).
Figure 1: In a code exercise click the hints button to guide you through the exercise.
The tutorials are self-contained (you practice code in code boxes). However, so you get practice at working in r rstudio() I strongly recommend that you create an Quarto document within an r rstudio() project and practice everything you do in the tutorial in the Quarto document, make notes on things that confused you or that you want to remember, and save it. Within this document you will need to load the relevant packages and data.
Packages
This tutorial uses the following packages:
afex [@singmann_afex_2020]
here [@R-here]
Hmisc [@R-Hmisc] is loaded by ggplot2
knitr [@R-knitr]
It also uses these tidyverse packages [@R-tidyverse; @tidyverse2019]: dplyr [@R-dplyr], forcats [@R-forcats], ggplot2 [@wickhamGgplot2ElegantGraphics2016], and readr [@R-readr].
Coding style
There are (broadly) two styles of coding:
Explicit: Using this style you declare the package when using a function: package::function(). For example, if I want to use the mutate() function from the package dplyr, I will type dplyr::mutate(). If you adopt an explicit style, you don't need to load packages at the start of your Quarto document (although see below for some exceptions).
Concise: Using this style you load all of the packages at the start of your Quarto document using library(package_name), and then refer to functions without their package. For example, if I want to use the mutate() function from the package dplyr, I will use library(dplyr) in my first code chunk and type the function as mutate() when I use it subsequently.
Coding style is a personal choice. The Google r rproj() style guide and tidyverse style guide recommend an explicit style, and I use it in teaching materials for two reasons (1) it helps you to remember which functions come from which packages, and (2) it prevents clashes resulting from using functions from different packages that have the same name. However, even with this style it makes sense to load tidyverse because the dplyr and ggplot2 packages contain functions that are often used within other functions and in these cases explicit code is difficult to read. Also, no-one wants to write ggplot2:: before every function from ggplot2.
You can use either style in this tutorial because all packages are pre-loaded. If working outside of the tutorial, load the tidyverse package (and any others if you're using a concise style) at the beginning of your Quarto document:
library(tidyverse)
Data
To work outside of this tutorial you need to download the following data files:
Set up an r rstudio() project in the way that I recommend in this tutorial, and save the data files to the folder within your project called [data]{.alt}. Place this code in the first code chunk in your Quarto document:
This code reads in the data and converts the variables strategy, looks and personality to be factors (categorical variable). It also uses fct_relevel() to set the order of the levels of each factor.
For the profile picture data use the code below. The variables rel_status and profile_pic are converted to factors and the reference category for each is set using fct_relevel():
The main examples in this tutorial are from [@field_discovering_2023]. A big discussion in magazines seems to be the factors that get you a relationship, for example, the relative importance of looks, personality, and dating strategies (whether you should 'treat them mean to keep them keen' and all that stuff). Scientists have looked at these issues too. For example, the top three most highly rated attributes of a partner in teenagers are reliability, honesty and kindness [@ha_effects_2010]. Beyond that, in the same study boys tended to rate attractiveness slightly higher than girls, and girls rate a sense of humour more highly than boys (although both are ranked in the top 10 by both sexes). With regard to dating strategies, @dai_when_2014 suggest that if someone is committed to pursuing a relationship with a person who plays hard to get, they will find that person more desirable but less likeable.
Imagine a scientist designed a study to look at the interplay between looks, personality and dating strategies on evaluations of a date. She set up a speed-dating night with nine tables at which there sat a 'date'. All the dates were stooges selected to vary in their attractiveness (high, average, and low), you might remember from a previous lecture that this could perhaps be gauged with facial symmetry (symmetrical faces tend to be perceived as more attractive). The dates were trained to vary their charisma (high, average, low), and also the strategy they employ during the conversation (normal or playing hard to get). The stooge dates were trained before the study to act charismatically to varying degrees, and also how to act in a way that made them seem unobtainable (hard to get) or not. As such, across the nine dates/stooges there were three attractive dates, one of whom acted charismatically, one who acted normally (average) and another who acted like they write statistics books, and likewise for the three average-looking dates and the three low attractiveness dates.
Therefore, each participant attending the speed-dating night was exposed to all combinations of attractiveness and charisma (these are repeated measures). (There were separate sets of nine male stooges and nine females so that those attending could meet 'dates' of whichever sex interested them.) Upon arrival participants were randomly assigned a green or red sticker. For the participants with the red sticker the stooges played hard to get (unobtainable) and for those with a green sticker they acted normally. Over the course of a few nights 20 people attended, spent 5 minutes with each of the nine 'dates' and then rated how much they'd like to have a proper date with each one as a percentage (100% = 'I'd pay a large sum of money for their phone number', 0% = 'I'd pay a large sum of money for a plane ticket to get me as far away from them as possible').
r alien() Alien coding challenge
View the data in [date_tib]{.alt}.
date_tib
Note that there are four variables: the id, which is a character variable (note the under the name), strategy (normal or hard to get), looks (high, average, low) and personality (high, average or low charisma), which are all factor variables (note the under the name) and date (the rating of how much the participant would like a date with the stooge), which is numeric and has the data type 'double' (note the under the name).
The data are in tidy format, which means that each row represents an instance of the outcome variable and the columns code information about each instance, for example, which participant the instance is related to and the characteristics of each of the 9 'dates' they met. Consequently each participant occupies 9 rows of the tibble (because each participant contributes 9 instances of the outcome variable, date).
The variables strategy, looks and personality are factors (categorical variable), so having read the data file and converted these variables to factors it's a good idea to check the order of the levels for each one.
r alien() Alien coding challenge
Using what you've learnt in previous tutorials check the order of the levels of the variables strategy, looks and personality.
# use this function:levels()
# Remember that to access a variable you use:
name_of_tibble$name_of_variable
You'll find that the factor levels are in a sensible order:
For strategy, the baseline category is 'normal', which makes sense.
For looks and personality, the factor levels are ordered low, average, high. This makes sense because low is the baseline and the levels go in ascending order.
These orders are because I set the data up within this tutorial. Working outside of the tutorial you might need to manually set the order of levels for each factor using forcats::fct_relevel() as described in [Data]{.alt} part of this tutorial.
r bmu() Exploring the data [(1)]{.alt}
We can create means, confidence intervals and plots in the usual way using ggplot2 and dplyr but for this tutorial we will use the methods built into the afex and emmeans packages.
r user_visor() Contrasts [(2)]{.alt}
Before we fit the model we need to think about what contrasts would be sensible. The variable, looks, has three conditions: high, average and low attractiveness. It makes sense to compare the high and low conditions to the average because the average person represents the norm (although it wouldn't be wrong to, for example, compare high and average attractiveness to low attractiveness).
Similarly, personality has a category, average charisma, that represents the norm. So, we could use this as the reference/control against which to compare our two extremes (when the date displays high and low charisma).
To implement contrasts we can use the emmeans::contrast() function, which we have used before to set our own bespoke contrasts. This function also has some built in methods, that will make our life easier. You might remember from [discovr_15]{.alt} that the function takes this general form:
In which we replaced [my_emm_object]{.alt} with an object that contains the estimated marginal means for the effect of interest, and we replaced [my_contrasts]{.alt} with the name of the object containing our bespoke contrasts. To use the built in contrasts, we use the [method]{.alt} argument instead of [my_contrasts]{.alt}:
Be default the method is set to [eff]{.alt}. Here are some of the different settings we can use:
[eff]{.alt}: each level/category of a variable is compared to the average of all levels/categories.
[pairwise]{.alt}: each level/category is compared to all other levels/categories in turn.
[poly]{.alt}: polynomial contrasts (i.e. tests for an ordered trend in category means).
[trt.vs.ctrl]{.alt}: compares each level/category to a declared reference category (by default the first category). Use [ref = x]{.alt} to make [x]{.alt} the reference category. For example. [ref = 2]{.alt} will compare each level/category to the second level.
[consec]{.alt}: compares each level/category (except the first) to the previous category.
Given that we want to set contrasts for both looks and personality that compare each category to average, and we know from the previous exercise that 'average' is the second factor level for both variables, we could set contrasts for both variables using [trt.vs.ctrl]{.alt} and [ref = 2]{.alt}. Remember this as we work through the tutorial.
r user_visor() Fitting a model to a mixed design [(2)]{.alt}
We can fit a model that predict the rating of the date from strategy, looks and personality using the afex package, which we met in discovr_13 and discovr_15.
We can use the aov_4() function in much the same way as in discovr_15. Remember that the aov_4() function has the following format:
afex::aov_4(outcome ~ predictors + (rm_predictors|id_variable), data = my_tib)
`r cat_space()` **Tip**
Remember that we can specify all main effects and their interactions using `*`. For example, `strategy*looks*personality` will introduce the main effect of **strategy**, the main effect of **looks**, the main effect of **personality**, all of the two-way interaction terms ($\text{strategy} \times \text{looks}$, $\text{strategy} \times \text{personality}$, and $\text{looks} \times \text{personality}$), and the three-way interaction ($\text{strategy} \times \text{looks} \times \text{personality}$).
In the current design, we can specify the predictors and their interaction using [strategy*looks*personality]{.alt}, therefore, we replace [predictors]{.alt} in the code above with [strategy*looks*personality]{.alt}. We also need to specify which predictors are repeated measures by specifying which variables are nested within the participant. Remembering that the variable that uniquely identifies the participants is [id]{.alt}, we'd replace [(rm_predictors\|id_variable)]{.alt} with [(looks*personality\|id)]{.alt} to indicate that the variables looks, personality and their interaction are (repeated measures) nested within the variable id.
r robot() Code example
Remembering that the outcome variable is date, and the tibble containing the data is called [date_tib]{.alt}, we can put it all together to fit the model with this code:
# fit the model (replace the xs):
date_afx <- afex::aov_4(xxxxxx ~ xxxx*xxxxx*xxxxx + (xxxxx*xxxxx|xxxxx), data = xxxxx)
# fit the model:
date_afx <- afex::aov_4(date ~ strategy*looks*personality + (looks*personality|id), data = date_tib)
date_afx #this shows us the model
Working down from the top of the table, we find a non-significant effect of strategy but significant effects (the [p.value]{.alt} is less than 0.05) of looks, the strategy × looks interaction, personality, the strategy × personality interaction, the looks × personality interaction and the strategy × looks × personality interaction. Everything, basically. You wouldn't normally be interested in main effects when there are significant interactions. In fact, the only effect we're really interested in is the significant strategy × looks × personality interaction, but for completeness of your understanding we'll interpret each effect in turn, starting with the main effect of strategy.
r user_visor() The main effect of strategy [(2)]{.alt}
quiz(caption ="The main effect of strategy (level 2)",
question("How would you interpret the effect of **strategy**?",
answer("There was a non-significant effect of strategy because the *p*-value is 0.949, which is greater than 0.05. This effect tells us that if we ignore all other variables, ratings were equivalent regardless of whether the date adopted a hard to get persona or not.", correct =T),
answer("There was a significant effect of strategy because the *p*-value is 0.949, which is greater than 0.05. This effect tells us that if we ignore all other variables, ratings were depended on whether the date adopted a hard to get persona or not.", message ="A *p* value *less* than 0.05 is typically interpretted as a significant effect."),
correct ="Correct - well done!",
incorrect ="Sorry, this answer is incorrect.",
random_answer_order =TRUE,
allow_retry =T
)
)
r user_visor() The main effect of looks [(2)]{.alt}
The main output from the afex model is reproduced below:
date_afx
quiz(caption ="The main effect of looks (level 2)",
question("How would you interpret the effect of **looks** in the output above?",
answer("There was a significant main effect of looks, *F*(1.92, 34.62) = 423.73, *p* < 0.001, which means that if we ignore all other variables, ratings of attractive, average and unattractive dates differed.", correct =T),
answer("There was a non-significant main effect of looks, *F*(1.92, 34.62) = 423.73, *p* < 0.001, which means that if we ignore all other variables, ratings of attractive, average and unattractive dates were equivalent.", message ="A *p* value *less* than 0.05 is typically interpretted as a *significant* effect."),
correct ="Correct - well done!",
incorrect ="Sorry, this answer is incorrect.",
random_answer_order =TRUE,
allow_retry =T
)
)
To pick the effect apart we can first get the estimated marginal means using emmeans(), which we used in both which we used in [discovr_13]{.alt} and which we used in [discovr_15]{.alt}.
r robot() Code example
Remember that the emmeans() function takes this general form when we have repeated measures:
my_new_emm <- emmeans::emmeans(my_afx_model, ~predictor(s), model ="multivariate")
In which [my_new_emm]{.alt} is whatever we choose to call the saved means, [my_afx_model]{.alt} is the name of the afex model we have created, in our case [date_afx]{.alt}, and [predictor(s)]{.alt} is a list of the predictors for which we want means.
r alien() Alien coding challenge
Obtain the estimated marginal means for the effect of looks and store them as [looks_emm]{.alt}
# replace the Xs
looks_emm <- emmeans::emmeans(xxxxxx, ~xxxxxx, model ="xxxxxx")
looks_emm <- emmeans::emmeans(date_afx, ~looks, model ="multivariate")
looks_emm # shows us the means
From this table you can see that as attractiveness of the date increases from low to high, the desirability to have a date with them increases too. Contrasts will help us to understand exactly what's going on.
r robot() Code example
We mentioned earlier that a reasonable contrast would be to use the average attractiveness as a reference category and compare the low and high attractiveness dates to this 'norm'. We also saw that we could achieve this using emmeans::contrast() and defining the method as [trt.vs.ctrl]{.alt} and the reference category as the second category (average), [ref = 1]{.alt}. Remember that we feed into the function the estimated marginal means that we just calculated [looks_emm]{.alt}. Putting this all together, we'd execute:
These contrasts tell us (in combination with the means themselves that) highly attractive dates were rated significantly higher than average dates, r report_em(looks_con), and averagely attractive dates were rated significantly higher than those low on attractiveness, r report_em(looks_con, row = 2).
r user_visor() The main effect of personality [(2)]{.alt}
The main output from the afex model is reproduced below:
date_afx
quiz(caption ="The main effect of personality (level 2)",
question("How would you interpret the effect of **personality** in the output above?",
answer("There was a significant main effect of charisma, *F*(1.87, 33.62) = 328.25, *p* < 0.001, which tells us that if we ignore all other variables, ratings for highly charismatic, a bit charismatic and dull dates differed.", correct =T),
answer("There was a non-significant main effect of charisma, *F*(1.87, 33.62) = 328.25, *p* < 0.001, which tells us that if we ignore all other variables, ratings for highly charismatic, a bit charismatic and dull dates were equivalent.", message ="A *p* value *less* than 0.05 is typically interpretted as a *significant* effect."),
correct ="Correct - well done!",
incorrect ="Sorry, this answer is incorrect.",
random_answer_order =TRUE,
allow_retry =T
)
)
Again, we could pick the effect apart with the estimated marginal means and contrasts. Remember from earlier, that we want to set the same contrasts: that is compare the average category of personality against the low and high categories. Therefore, the code will be the same as we just used for the variable looks except that:
We should give our estimated marginal means object a different name.
We replace the variable [looks]{.alt} with [personality]{.alt} in the code
r alien() Alien coding challenge
Obtain the estimated marginal means for the effect of looks and store them as [pers_emm]{.alt}. Adapt the code from the previous section and use the hints if you get stuck.
# replace the Xs
pers_emm <- emmeans::emmeans(xxxxxx, ~xxxxxx, model ="xxxxxx")
pers_emm <- emmeans::emmeans(date_afx, ~personality, model ="multivariate")
pers_emm # shows us the means# now get the contrasts:
emmeans::contrast(xxxxxxx, method ="xxxxxxx", ref =2, adjust ="holm")
pers_emm <- emmeans::emmeans(date_afx, ~personality, model ="multivariate")
pers_emm # shows us the means# now get the contrasts
emmeans::contrast(pers_emm, method ="trt.vs.ctrl", ref =2, adjust ="holm")
The contrasts and means tell us that dates with low charisma were rated significantly lower than those acting with average charisma, r report_em(pers_con), and dates acting with high charisma were rated significantly higher than dates with average charisma, r report_em(pers_con, row = 2).
r user_visor() The interaction between strategy and looks [(2)]{.alt}
The main output from the afex model is reproduced below:
date_afx
The output above tells us that dating strategy significantly interacted with the attractiveness of the date, r report_afx(date_tbl, row = 3): the profile of ratings across dates of different attractiveness was different depending on whether or not they played hard to get.
To pick the effect apart we can again inspect the estimated marginal means using emmeans(). We use the same approach as before, but remember that in this code:
my_new_emm <- emmeans::emmeans(my_afx_model, ~predictor(s), model ="multivariate")
We now replace [\~predictor(s)]{.alt} with the predictors involved in the interaction of interest by listing them withinc(). For example, for the strategy*looks interaction we'd specify [c("strategy", "looks")]{.alt} (note we drop the [\~]{.alt} at the beginning).
r alien() Alien coding challenge
Obtain the estimated marginal means for the strategy*looks interaction and store them as [strat_looks_emm]{.alt}
# replace the Xs
strat_looks_emm <- emmeans::emmeans(xxxxxx, c(xxxxxx), model ="xxxxxx")
strat_looks_emm <- emmeans::emmeans(date_afx, c("strategy", "looks"), model ="multivariate")
strat_looks_emm # shows us the means
quiz(caption ="The interaction of strategy and looks (level 2)",
question("Using the means, which of these statements about the **strategy × looks** interaction are true? (Select ALL that apply)",
answer("Playing hard to get led to higher ratings (than acting normally) for attractive dates, made no difference for average looking dates and led to lower ratings in unattractive dates", correct =T),
answer("The level of attractiveness had less impact when dates acted normally than when they played hard to get", correct =T),
answer("The level of attractiveness had more impact when dates acted normally than when they played hard to get."),
answer("The effect of playing hard to get (compared to acting normally) was similar regardless of how attractive the date was."),
answer("The effect of playing hard to get (compared to acting normally) was similar for attractive and unattractive dates but different for average looking dates."),
correct ="Correct - well done! The means show that for average-looking dates it doesn't make a difference whether they played hard to get. For attractive dates, ratings were higher when the date played hard to get than when they didn't and for unattractive dates the opposite was true – ratings were lower when dates played hard to get. In short, playing hard to get only has an effect at the extremes of attractiveness.",
incorrect ="You have at least one incorrect answer. The means show that for average-looking dates it doesn't make a difference whether they played hard to get. For attractive dates, ratings were higher when the date played hard to get than when they didn't and for unattractive dates the opposite was true – ratings were lower when dates played hard to get. In short, playing hard to get only has an effect at the extremes of attractiveness.",
random_answer_order =TRUE,
allow_retry =T
)
)
We can get contrasts for this interaction term. To use the emmeans::contrast() function on an interaction term we have to use the [interaction]{.alt} argument to specify the contrast for each variable in the interaction. In general, then:
In which we replaced [my_emm_object]{.alt} with an object that contains the estimated marginal means for the effect of interest, [interaction = c(variable_a = "contrast", variable_b = "contrast" ...)]{.alt} lists each variable and the contrast you want to apply to it and [ref = x]{.alt} is used to set the reference category for the contrasts (if relevant).
In previous contrasts we have used [trt.vs.ctrl]{.alt} as the contrast for the variables looks and personality, so we should continue to do so. We may as well use [trt.vs.ctrl]{.alt} for strategy too, because it has only two categories, so [trt.vs.ctrl]{.alt} will compare these categories. Therefore, we can use [interaction = c(strategy = "trt.vs.ctrl", looks = "trt.vs.ctrl")]{.alt} to specify the contrast for the strategy*looks interaction. We must also remember to feed the estimated marginal means into the function, which we craeted above and stored as [strat_looks_emm]{.alt}.
r robot() Code example
Putting this together, the code for the strategy*looks contrasts will be:
The first contrast for the interaction term shows the effect of unattractive dates compared to average dates, comparing playing hard to get to normal. This contrast is highly significant, r report_em(strat_looks_con), suggesting that the increased interest in average-looking dates compared to unattractive dates is significantly greater when dates played hard to get than when they acted normally.
The second contrast, which compares playing hard to get to normal in high attractive dates compared to average-looking dates is also highly significant, r report_em(strat_looks_con, row = 2). This contrast tells us that the decreased interest in average-looking dates compared to attractive dates is significantly more when dates played hard to get than when they acted normally.
r user_visor() The interaction between strategy and personality [(2)]{.alt}
The main output from the afex model is reproduced below:
date_afx
The output reproduced above reminds us that the dating strategy significantly interacted with how charismatic the date was,
r report_afx(date_tbl, row = 5). This effect means that the profile of ratings across dates of different levels of charisma was influenced by the dating strategy employed.
We can adapt the code from the previous section to get the estimated marginal means and contrasts for this interaction. The main changes will be:
We should call the object storing the estimated marginal means something different to the previous section.
We need to replace [looks]{.alt} with [personality]{.alt} in the code.
r alien() Alien coding challenge
Create an object called [strat_pers_emm]{.alt} containing the estimated marginal means for the strategy*personality interaction. Also obtain the contrasts as in the previous section.
# get the estimated marginal means:
strat_pers_emm <- emmeans::emmeans(xxxx, c("xxxxx", "xxxx"), model ="xxxxxx")
# get the estimated marginal means:
strat_pers_emm <- emmeans::emmeans(date_afx, c("strategy", "personality"), model ="multivariate")
# display the means
xxxxx
# get the estimated marginal means:
strat_pers_emm <- emmeans::emmeans(date_afx, c("strategy", "personality"), model ="multivariate")
# display the means
strat_pers_emm
# get contrasts
emmeans::contrast(xxxxxxx, interaction =c(xxxxxx ="xxxxxx", xxxx ="xxxxxxx"), ref = x, adjust ="holm" )
strat_pers_emm <- emmeans::emmeans(date_afx, c("strategy", "personality"), model ="multivariate")
strat_pers_emm # display the means
emmeans::contrast(strat_pers_emm, interaction =c(strategy ="trt.vs.ctrl", personality ="trt.vs.ctrl"), ref =2, adjust ="holm")
quiz(caption ="The strategy × personality interaction (level 2)",
question("Using the means, which of these statements about the **strategy × personality** interaction are true? (Select ALL that apply)",
answer("Playing hard to get led to lower ratings (than acting normally) for charismatic dates, made no difference for dates with average charisma and led to higher ratings in dates who behaved like dullards", correct =T),
answer("The level of charisma had more impact when dates acted normally than when they played hard to get", correct =T),
answer("The level of charisma had less impact when dates acted normally than when they played hard to get."),
answer("The effect of playing hard to get (compared to acting normally) was similar regardless of how charismatic the date was."),
answer("The effect of playing hard to get (compared to acting normally) was similar for dates who acted highly charismatic and those who acted like dullards but different for dates who acted with average charisma."),
correct ="Correct - well done! The means show that for dates with normal amounts of charisma the dating strategy they adopted had little impact. However, highly charismatic dates were rated higher when they acted normally compared to when they played hard to get.",
incorrect ="You have at least one incorrect answer. The means show that for dates with normal amounts of charisma the dating strategy they adopted had little impact. However, highly charismatic dates were rated higher when they acted normally compared to when they played hard to get. Conversely, dull dates were rated higher when they played hard to get compared to when they acted normally. Although interest in dating decreases as charisma decreases, this decreases is less pronounced when the dates play hard to get.",
random_answer_order =TRUE,
allow_retry =T
)
)
The decreased interest in low charisma dates compared to averagely charismatic dates is significantly less when dates play hard to get than when they act normally, r report_em(strat_pers_con).
The increased interest in highly charismatic dates compared to averagely charismatic dates found when dates acted normally is significantly more than when the dates played hard to get, r report_em(strat_pers_con, row = 2).
r user_visor() The interaction between looks and personality [(2)]{.alt}
The main output from the afex model is reproduced below:
date_afx
The reproduced output reminds us that there was a significant looks × personality interaction, r report_afx(date_tbl, row = 6). This effect tells us that the profile of ratings across dates of different levels of charisma was different for attractive, average and unattractive dates.
We can (again) adapt the code from the previous two sections to get the estimated marginal means and contrasts for this interaction. The main changes will be:
We should call the object storing the estimated marginal means something different to the previous section.
We need to specify [looks]{.alt} and [personality]{.alt} as the variables in the interaction.
r alien() Alien coding challenge
Create an object called [looks_pers_emm]{.alt} containing the estimated marginal means for the looks*personality interaction. Also obtain the contrasts as in the previous section.
# get the estimated marginal means:
looks_pers_emm <- emmeans::emmeans(xxxxxx, c("xxxxx", "xxxxx"), model ="xxxxxx")
# get the estimated marginal means:
looks_pers_emm <- emmeans::emmeans(date_afx, c("looks", "personality"), model ="multivariate")
# display the means
xxxxx
# get the estimated marginal means:
looks_pers_emm <- emmeans::emmeans(date_afx, c("looks", "personality"), model ="multivariate")
# display the means
looks_pers_emm
# get contrasts
emmeans::contrast(xxxxxxx, interaction =c(xxxxxx ="xxxxxx", xxxx ="xxxxxxx"), ref = x, adjust ="holm" )
looks_pers_emm <- emmeans::emmeans(date_afx, c("looks", "personality"), model ="multivariate")
looks_pers_emm # display the means
emmeans::contrast(looks_pers_emm, interaction =c(looks ="trt.vs.ctrl", personality ="trt.vs.ctrl"), ref =2, adjust ="holm")
To answer the following questions you might it useful to plot the means (because there are a lot of them). Do this with the code below:
afex::afex_plot(date_afx,
x ="looks",
trace ="personality",
mapping ="colour"
) +theme_minimal()
quiz(caption ="The looks × personality interaction (level 2)",
question("Using the means, which of these statements about the **looks × personality** interaction are true? (Select ALL that apply)",
answer("For attractive dates, having high or average charisma resulted in similar (high) ratings, but being a dullard led to lower ratings", correct =T),
answer("For average-looking dates ratings increased as charisma went from low to high", correct =T),
answer("For unattractive dates being a dullard or having average charisma results in similar (low) ratings, but being highly charismatic led to higher ratings", correct =T),
answer("The effect of charisma was similar for attractive and unattractive dates, but different for average-looking dates."),
answer("Attractive dates were rated higher regardless of charisma."),
answer("The effect of charisma is present when comparing attractive and average-looking dates but not when comparing unattractive and average-looking dates."),
answer("For dates acting with high charisma, attractiveness had no impact on ratings"),
answer("For dates acting with low charisma, attractiveness had no impact on ratings"),
correct ="Correct - well done! The graph shows the mean ratings of dates of different levels of attractiveness when the date also had high levels of charisma (blue line), average charisma (green line) and low charisma (red line). Look first at the difference between attractive and average-looking dates. The interest in highly charismatic dates doesn't change (the blue line is more or less flat between these two points), but for dates with average charisma or low charisma interest levels decline (the red and green lines slope down). If you have lots of charisma you can get away with being average-looking and people will still want to date you. Now look at the difference between average-looking and unattractive dates. A different pattern is observed: for dates with low charisma there is little difference between unattractive and average-looking people (the red line is flat) but for those with average or high charisma, there is a decline in interest if they are low attractive (the green and blue lines slope down).",
incorrect ="You have at least one incorrect answer. The graph shows the mean ratings of dates of different levels of attractiveness when the date also had high levels of charisma (blue line), average charisma (green line) and low charisma (red line). Look first at the difference between attractive and average-looking dates. The interest in highly charismatic dates doesn't change (the blue line is more or less flat between these two points), but for dates with average charisma or low charisma interest levels decline (the red and green lines slope down). If you have lots of charisma you can get away with being average-looking and people will still want to date you. Now look at the difference between average-looking and unattractive dates. A different pattern is observed: for dates with low charisma there is little difference between unattractive and average-looking people (the red line is flat) but for those with average or high charisma, there is a decline in interest if they are low attractive (the green and blue lines slope down).",
random_answer_order =TRUE,
allow_retry =T
)
)
The contrasts help to pick apart this interaction.
The first contrast addresses the question 'is the difference between low charisma and average charisma the same for unattractive people and average-looking people? For average-looking dates, ratings were higher when they had average charisma than when they had low charisma, but for unattractive dates the ratings were roughly the same regardless of the level of charisma. This contrast is highly significant, r report_em(looks_pers_con).
The second contrast asks the question 'is the difference between low charisma and average charisma the same for attractive people and average-looking people? We can again focus on the relevant part of the interaction graph. Interest in attractive dates was higher when they had average charisma (green) than when they had low charisma (red); the same is also true for average-looking dates. The two lines are fairly parallel, which is reflected in the non-significant contrast, r report_em(looks_pers_con, row = 2). It seems that as the attractiveness of dates is reduced there is a decline in interest both when charisma is average and when it is low.
The third contrast asks 'is the difference between high charisma and average charisma the same for low attractive dates and average-looking dates?' Interest in dating decreases from average-looking dates to unattractive ones in dates with both high and average charisma dates; however, this fall is significantly greater in the low-charisma dates (the red line is slightly steeper than the blue), r report_em(looks_pers_con, row = 3). As dates' attractiveness is reduced there is a significantly greater decline in interest when dates have average charisma compared to when they have a lot.
The final contrast asks 'is the difference between high charisma and average charisma the same for attractive and average-looking dates?' The best way to understand this contrast is to focus on the relevant bit of the interaction graph. Interest (as indicated by high ratings) in attractive dates was the same regardless of whether they had high or average charisma; however, for average-looking dates, there was more interest when that person had high charisma rather than average charisma. This contrast is highly significant, r report_em(looks_pers_con, row = 4), and tells us that as attractiveness is reduced there is a significantly greater decline in interest when charisma is average compared to when it is high.
r user_visor() The interaction between strategy, looks and personality [(2)]{.alt}
The three-way interaction will blow your mind! The main output from the afex model is reproduced below:
date_afx
The significant strategy × looks × personality interaction, r report_afx(date_tbl, row = 7), tells us whether the looks × personality interaction described above is the same (statistically) when dates played hard to get compared to when they didn't.
To pick the effect apart we can, take a very deep breath and, inspect the estimated marginal means using emmeans(). We use the same approach as before, but adapt the code from the previous section to include the third predictor:
three_way_emm <- emmeans::emmeans(date_afx, c("strategy", "looks", "personality"), model ="multivariate")
three_way_emm # shows us the means
r alien() Alien coding challenge
Obtain the estimated marginal means for the strategy*looks*personality interaction and store them as [three_way_emm]{.alt}
# replace the Xs
three_way_emm <- emmeans::emmeans(xxxxxx, c(xxxxxx, xxxxx, xxxxx), model ="xxxxxx")
three_way_emm <- emmeans::emmeans(date_afx, c("strategy", "looks", "personality"), model ="multivariate")
three_way_emm # shows us the means
We now have 18 means to wrap our heads around, which is a lot, so lets plot them, using afex_plot(). The code below plots looks on the x-axis, personality as different coloured lines and strategy is split across two panels.
quiz(caption ="The three-way interaction (level 2)",
question("Look at the interaction graph and compare the graph for the dates why acted normally to the one for the dates who played hard to get. The following statements highligh the *differences* between these two graphs, which of them are true? (Select ALL that apply)",
answer("For attractive dates, charisma had no effect on ratings when dates played hard to get, but when they acted normally dullards were rated much lower", correct =T),
answer("For average-looking dates, ratings increased as charisma went from low to high regardless of the dating strategy employed", correct =T),
answer("For unattractive dates charisma had a very small effect on ratings when dates played hard to get, but when they acted normally highly charismatic dates were rated much higher", correct =T),
answer("There are no differences in the patterns of means between the graphs depicting dates who acted normally and those that played hard to get."),
answer("The patterns of means between the graphs depicting dates who acted normally and those that played hard to get differ only for attractive dates."),
answer("For dates acting with high charisma, attractiveness had a similar impact for dates who acted normally and those that played hard to get"),
answer("For dates acting with average charisma, attractiveness had a greater impact for dates who acted normally than for those that played hard to get"),
correct ="Correct - well done! The graph for dates who played hard to get shows that when dates are attractive, a high interest was expressed regardless of charisma levels (the blue, green and red lines meet). At the opposite end of the attractiveness scale, when a date is unattractive, regardless of charisma, very little interest is expressed (ratings are all low). If the date plays hard to get, the only time charisma makes a difference is if the date is average-looking, in which case high charisma (blue) boosts interest, having low charisma (red) reduces interest, and having average charisma leaves things somewhere in between. The take-home message is that playing hard to get only works if you're averagely attractive: if you're highly charismatic it will boost your appeal but its a disastrous strategy if you're dull. The picture when the date doesn't play hard to get (acts normally) is different. If someone has high levels of charisma then what they look like won't affect interest in them (the blue line is relatively flat). At the other extreme, if the date has low charisma, then they will have very little interest expressed in them regardless of how attractive they are (the red line is relatively flat). The only time attractiveness makes a difference is when someone has an average amount of charisma (the green line), in which case being attractive boosts interest, and being unattractive reduces it. If you don't play hard to get and you are averagely attractive then you can influence others interest you with your charisma.",
incorrect ="You have at least one incorrect answer. The graph for dates who played hard to get shows that when dates are attractive, a high interest was expressed regardless of charisma levels (the blue, green and red lines meet). At the opposite end of the attractiveness scale, when a date is unattractive, regardless of charisma, very little interest is expressed (ratings are all low). If the date plays hard to get, the only time charisma makes a difference is if the date is average-looking, in which case high charisma (blue) boosts interest, being a dullard (red) reduces interest, and having some charisma leaves things somewhere in between. The take-home message is that playing hard to get only works if you're averagely attractive: if you're highly charismatic it will boost your appeal but its a disastrous strategy if you're dull. The picture when the date doesn't play hard to get (acts normally) is different. If someone has high levels of charisma then what they look like won't affect interest in them (the blue line is relatively flat). At the other extreme, if the date has low charisma, then they will have very little interest expressed in them regardless of how attractive they are (the red line is relatively flat). The only time attractiveness makes a difference is when someone has an average amount of charisma (the green line), in which case being attractive boosts interest, and being unattractive reduces it. If you don't play hard to get and you are averagely attractive then you can influence others interest you with your charisma.",
random_answer_order =TRUE,
allow_retry =T
)
)
We can again get contrasts for this interaction term. For the two-way interactions we used code of this form (I'm using the example of the looks*personality interaction):
For the three-way interaction we need to add strategy to the list of variables and assign it a contrast type. Like we did in previous interactions, we may as well use [trt.vs.ctrl]{.alt} as a contrast for strategy has only two categories, so [trt.vs.ctrl]{.alt} will compare these categories. We also need to replace [looks_pers_emm]{.alt} with the estimated means object for the three-way interaction that we just created [three_way_emm]{.alt}.
The first contrast was not significant, r report_em(three_way_con). This contrast looks at the effect of Strategy when comparing low attractive dates to average-looking dates, when displaying low charisma relative average charisma. Interest in unattractive dates was the same regardless of whether they had average charisma or low (the green and red dots are in the same place). Interest in average-looking dates was greatest when they had average charisma rather than when they had low charisma (the green dot is higher than the red). This pattern of results is very similar when dates played hard to get and when they didn't.
The second contrast explores effect of high attractive dates relative to average-looking ones, when displaying low charisma compared to average charisma, when dates played hard to get relative to when they didn't. The contrast is significant, r report_em(three_way_con, row = 2), which reflects the fact that the pattern of means is different when dates played hard to get compared to when they didn't. First, if we look at average-looking dates, more interest was expressed when the date has average charisma than when they have low charisma, and this is true whether or not dates played hard to get (the distance between the green and red lines is about the same in the two dating strategy groups). So, the difference created by playing hard to get doesn't appear to be here. Now look at attractive dates. When dates played hard to get the interest in the date is high regardless of their charisma (the lines meet). However, when dates acted normally interest in dating an attractive person is much lower if they are a dullard (the red dot is much lower than the green). Another way to look at it is that for dates with average charisma, the reduction in interest as attractiveness goes down is about the same regardless of whether they played hard to get (the green lines have the same slope). However, for dates displaying low charisma, the decrease in interest if these dates are average-looking rather than attractive is much more dramatic if they play hard to get (the red line is steeper in the hard to get group).
The third contrast was also significant, r report_em(three_way_con, row = 3). This contrast compares low attractive dates to average-looking ones, when displaying high charisma relative to average charisma, when dates played hard to get relative to when they didn't. First, let's look at when dates played hard to get. As attractiveness goes down, so does interest when the date has high charisma and when they have average charisma (the slopes of the blue and green lines are similar). So, regardless of charisma, there is a similar reduction in interest as attractiveness declines. Now let's look at when the dates acted normally. The picture is quite different: when charisma is high, there is no decline in interest as attractiveness falls (the blue line is flat); however, when charisma is average, interest is lower in an unattractive date than in an average-looking date (the green line slopes down). Another way to look at it is that for dates with average charisma, the reduction in interest as attractiveness goes down is about the same regardless of whether dates play hard to get (the green lines have similar slopes). However, for dates who have high charisma, the decrease in interest if these dates are unattractive rather than average-looking is much more dramatic when dates played hard to get than when they didn't (the blue line is steeper for dates that played hard to get).
The final contrast shows the effect of high attractive dates relative to average-looking ones, when they display high charisma compared average charisma, when dates played hard to get relative to when they didn't, r report_em(three_way_con, row = 4). This contrast seems to show that interest in dating (as indicated by high ratings) attractive dates was the same regardless of whether they had high or average charisma (the green and blue dots are in a similar place). However, for average-looking dates, there was more interest when that person had high charisma rather than average charisma (the green dot is lower than the blue dot). The non-significance of this contrast indicates that this pattern of results is very similar when dates played hard to get and when they didn't.
`r pencil()` **Report**`r rproj()`
Greenhouse-Geisser corrected degrees of freedom are reported throughout. The main effect of strategy was not significant, `r report_afx(date_tbl, row = 1)`, but the main effects of looks, `r report_afx(date_tbl, row = 2)`, and personality, `r report_afx(date_tbl, row = 4)`, were. These effects were superseded by the following significant interactions: strategy $\times$ looks, `r report_afx(date_tbl, row = 3)`, strategy $\times$ personality, `r report_afx(date_tbl, row = 5)`, and personality $\times$ looks, `r report_afx(date_tbl, row = 6)`. These interactions were also superseded by the significant strategy $\times$ personality $\times$ looks interaction, `r report_afx(date_tbl, row = 7)`. Contrasts were used to break down this interaction.
The first contrast looked at the effect of strategy when comparing low attractive dates to average-looking dates, when displaying low charisma relative average charisma. The contrast was not significant, `r report_em(three_way_con)`, suggesting that interest in unattractive dates was the same regardless of whether they had average charisma or low (the green and red dots are in the same place) and that this pattern was similar for the two dating strategies.
The second contrast explored the effect of high attractive dates relative to average-looking ones, when displaying low charisma compared to average charisma, when dates played hard to get relative to when they did not. The contrast was significant, `r report_em(three_way_con, row = 2)`. For average-looking dates, more interest was expressed when the date has average charisma than when they had low charisma, and this was true whether or not dates played hard to get. For attractive dates, when they played hard to get the interest in them was high regardless of their charisma. However, when attractive dates acted normally interest in them was much lower if they had low charisma.
The third contrast was also significant, `r report_em(three_way_con, row = 3)`. This contrast compared low attractive dates to average-looking ones, when displaying high charisma relative to average charisma, when dates played hard to get relative to when they did not. When dates played hard to get, as attractiveness decreased so did interest in the date regardless of charisma. However, when the dates acted normally, at high levels of charisma there was no decline in interest as attractiveness falls; however, when charisma is average, interest was lower for unattractive dates than for average-looking dates.
The final contrast showed the effect of high attractive dates relative to average-looking ones, when they display high charisma compared average charisma, when dates played hard to get relative to when they did not. This contrast was not significant, `r report_em(three_way_con, row = 4)`. This contrast suggests that interest in dating attractive dates was unaffected by whether they had high or average charisma. However, for average-looking dates, there was more interest when that person had high charisma rather than average charisma. The non-significance of this contrast indicates that this pattern of results was similar for the two dating strategies.
**A message from Mae Jemstone:**
No-one reads these do they? I'll write this properly when it's not 1am (again ...).
I need sleep sooooo badly.
Resources {data-progressive="FALSE"}
Statistics
The tutorials typically follow examples described in detail in @field_discovering_2023. That book covers the theoretical side of the statistical models, and has more depth on conducting and interpreting the models in these tutorials.
If any of the statistical content doesn't make sense, you could try my more introductory book An adventure in statistics [@fieldAdventureStatisticsReality2016].
There are free lectures and screencasts on my YouTube channel.
R for data science by @wickhamDataScience2017 is an open-access book by the creator of the tidyverse (Hadley Wickham). It covers the tidyverse and data management.
ModernDive is an open-access textbook on r rproj() and r rstudio().
I'm extremely grateful to Allison Horst for her very informative blog post on styling learnr tutorials with CSS and also for sending me a CSS template file and allowing me to adapt it. Without Allison, these tutorials would look a lot worse (but she can't be blamed for my colour scheme).
**Usage:** This tutorial accompanies [Discovering Statistics Using R and RStudio](https://www.discovr.rocks/) [@field_discovering_2023] by [Andy Field](https://en.wikipedia.org/wiki/Andy_Field_(academic)). It contains material from the book so there are some copyright considerations but I offer them under a [Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License](http://creativecommons.org/licenses/by-nc-nd/4.0/). Tl;dr: you can use this tutorial for teaching and non-profit activities but please don't meddle with it or claim it as your own work.
