First step: let’s load in the data and required packages. We’ll use the L2 English Lexical Decision dataset.

library(tidyverse)

lex <- read.csv("/Users/Thomas/Library/CloudStorage/OneDrive-YorkUniversity/LING 3300/Datasets/L2_English_Lexical_Decision_Data.csv", header = TRUE)

Correlation

Correlations measure the strength of association between two variables. Several types of correlation exist; in this course, we have mainly discussed the most common type, which is Pearson’s correlation. Pearson’s correlation measures the strength to which two continuous variables are linearly related to each other. Correlations range from -1 to +1, with 0 indicating no relationship between the two variables, +1 indicating a perfect positive relationship, and -1, a perfect negative relationship.

As an example, let’s consider two of the measures in this dataset. Remember that the workerID column is the anonymized code for each participant, and the domLang column is for the participant’s dominant language. In the RT column we have the number of milliseconds it took for that participant to see a word (in the word column) and respond on a keyboard whether or not that word is a real word in English. Other information that the researchers obtained for each participant is how confident they are in their English (confident), how frequently they speak English (hoursDay), how many long they’d studied English (yearsStudy), whether they live in an English speaking country (engCountry), and how many months they’d lived in an English speaking county (monthsLive).

Perhaps we wonder whether it’s the case that the longer an L2 speaker of English had lived in an English speaking country, the faster they’d react to the words. That seems like a reasonable hypothesis.

Let’s start by looking at, say, all the data we have from German speakers. We’ll filter the data so it only includes speakers whose dominant language is German. We’ll also only consider the data for which the reaction time is lower than 2000ms, since reaction times over two seconds are probably erroneous.

german <- lex %>%
  filter(domLang == "German") %>%
  filter(RT < 2000)

It’s a good idea to visualize our data. Let’s see if we can make a scatterplot that will show us the relationship between two continuous variables: monthsLive and RT.

ggplot(german) + 
  geom_point(aes(x = monthsLive, y = RT))

Oh bother. That’s not what we want. This is treating monthsLive as if it’s a categorical variable, not a continuous variable! And it’s putting the numbers in alphabetical order (like… 120 comes after 12 because they both start with ‘12’), which is… not how numbers work. Why could that be?

typeof(german$monthsLive)
## [1] "character"

Well that explains it. R is treating this column as a list of discrete character strings, not a bunch of numbers! We need to tell R that it should be numeric, even if that means turning whatever missing numbers in this column into “NA” values.

german$monthsLive <- as.numeric(german$monthsLive)

Let’s try again to get ourselves a scatterplot.

ggplot(german) + 
  geom_point(aes(x = monthsLive, y = RT))
## Warning: Removed 100 rows containing missing values or values outside the scale range
## (`geom_point()`).

Note that we get a warning just because it’s telling us it can’t plot the data where monthsLive is NA. That’s to be expected. But the weird thing now is that we have these big vertical lines of dots rather than single points. This is because we’re plotting all of the raw observations, and each participant responded to a lot of different words. It would make more sense for us to find the mean for each participant’s reaction time (the “by-participant means”) and then plot the relationship again. This is particularly important because a major assumption of correlation testing is that your data points are independent – that is, you only have one data point per participant.

german_means <- german %>%
  group_by(workerID, monthsLive) %>% 
  summarize(mean_RT = mean(RT)) %>%
  drop_na() # Let's drop all the rows where there's an NA so we have a clean dataset
## `summarise()` has regrouped the output.
## ℹ Summaries were computed grouped by workerID and monthsLive.
## ℹ Output is grouped by workerID.
## ℹ Use `summarise(.groups = "drop_last")` to silence this message.
## ℹ Use `summarise(.by = c(workerID, monthsLive))` for per-operation grouping
##   (`?dplyr::dplyr_by`) instead.
ggplot(german_means) + 
  geom_point(aes(x = monthsLive, y = mean_RT))

Alright! Now we can see the relationship between how many months a German speaker has lived in an English-speaking country, and how quickly on average they react to the stimulus in this experiment. By eyeballing it, it sort of looks like there’s a negative relationship as we’d hypothesized. Let’s actually test what the correlation is.

cor(german_means$monthsLive, german_means$mean_RT)
## [1] -0.1940825

This correlation is negligible in magnitude, and negative in direction. We can use the geom_smooth() function to draw a best-fit line for this data. It will automatically come with a 95% confidence interval shaded around the line.

ggplot(german_means) + 
  geom_point(aes(x = monthsLive, y = mean_RT)) +
  geom_smooth(aes(x = monthsLive, y = mean_RT),
              method = "lm",
              color = "red")
## `geom_smooth()` using formula = 'y ~ x'

Regression: A brief intro

Now let’s run a linear regression so that we can find the slope of this line as well as its y-intercept. We’ll put our predicted value (from the y-axis) on the left of the ~ and we’ll put the predictor (from the x-axis) on the right.

fit <- lm(mean_RT ~ monthsLive, german_means) # First we fit a model
summary(fit) # Then we run a summary() of that model
## 
## Call:
## lm(formula = mean_RT ~ monthsLive, data = german_means)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -213.30  -96.21  -46.99   74.21  423.54 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 786.9631    26.5470  29.644   <2e-16 ***
## monthsLive   -0.2258     0.1553  -1.454    0.152    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 146.8 on 54 degrees of freedom
## Multiple R-squared:  0.03767,    Adjusted R-squared:  0.01985 
## F-statistic: 2.114 on 1 and 54 DF,  p-value: 0.1518

From this output, the key info is in the Coefficients table. You can see that the y-intercept is 786.96, and that the slope of the line is -0.23. Our line intersects the y-axis (where x=0) at 786.96. This means that when a German has lived in an English speaking country for 0 months, we can predict (based on our sample) that on average they should have about a 786.96 ms reaction time to rating English words as real or not. And then for every month that they live in an English speaking country, their reaction time drops 0.23 ms – that is, they get 0.23 ms faster at the task. In terms of the line, this means that for one unit (month) moved right along the x-axis, you move down 0.23 units (ms) along the y-axis.

Don’t confuse the slope with the correlation of -0.19 we found above! Remember that a correlation quantifies how spread out the points are from the best-fit line, and needs to be between -1 and 1, with weaker correlations closer to 0. The slope can be of any magnitude (that is, way bigger or smaller than 1 or -1) and describes how steep the line is. What the two have in common is that they’ll both be negative if the line is going down toward the right and positive if going up toward the right. But it’s totally possible to have a very low slope with high correlation, or a very high slope with a low correlation. The magnitude of the slope will depend on the scale of the x and y axis units (e.g. if you measure your reaction time in seconds rather than ms, or your residence time in years rather than months, you’ll have much smaller numbers, and thus the slope coefficient will change despite the relationship between the two variables remaining the same), while the magnitude of the correlation will depend on how perfectly the data fit along a line (so this won’t be affected by changing the scales used to measure the variables).

You may also be interested in the R-squared value, which will always range from 0 to 1. Here we’re given two types of R-squared: Multiple R-squared and Adjusted R-squared. We usually care more about Adjusted R-squared, because it uses some fancy math that takes into account the number of predictors and sample size of our data, while Multiple R-squared is just the fraction of total variation explained by the independent variables without any other fancy math. In this case, the Adjusted R-squared is 0.0198, which is ridiculously low: it means that in this sample, the months lived in an English speaking country accounts for approximately 2% of the variation in mean reaction time. That’s not very much!

If you’re thinking that the R-squared seems to measure something very similar to correlation, you’d be correct! Check this out:

cor(german_means$monthsLive, german_means$mean_RT) # The Pearson's r
## [1] -0.1940825
(cor(german_means$monthsLive, german_means$mean_RT))^2 # The Pearson's r, squared
## [1] 0.037668

Look back at the Multiple R-squared shown in your linear regression model output. We can see that if you square the Pearson’s r, you get the Multiple R-squared! This means that the magnitude of the correlation and the magnitude of the R-squared are related to each other – but because the R-squared is squared, it will always be a positive number and won’t tell us anything about the direction of the relationship between the variables, while the plain old correlation will be positive if the slope it’s describing is positive, and negative if the slope is negative.


Practice with correlations and regressions

What we did above was get the overall mean of each German speaker’s RTs, and then we found the correlation with monthsLive. But perhaps it’s important to look separately at the conditions where the words were fake versus real – like maybe people were overall faster at responding to real words, or maybe living in an English speaking country for a long time only helps at speeding up RTs for real words and not fake ones.

  1. Read in the lexical decision dataset. Get a dataframe with just the data from the German-dominant speakers, and filter out all of the RTs that are above 2000ms.

  2. Now recode the condition column (or make a new column) such that if condition is ‘1’, it’s now marked as “real”, and if it has ‘2’, it’s marked as “fake”. This will make it more transparent which condition is which. You may want to use an ifelse() function here, though there are other ways to do this too.

  3. Make a dataframe that finds the by-participant, by-condition RT means. As above, you’ll want to also include monthsLive in your group_by() step so that we also have that information in the new dataframe. It would also be a good idea to drop the NA values, as above.

  4. Use ggplot to show the relationship between monthsLive (on the x-axis) and by-participant mean RT (on the y-axis). Facet over condition such that we end up with two side-by-side scatterplots, with each including the best-fit line.

  5. Consider the relationship between monthsLive and by-participant mean reaction times to fake words. Find the correlation as well as the y-intercept and slope of the best-fit line as shown in your plot.

  6. Consider the relationship between monthsLive and by-participant mean reaction times to real words. Find the correlation as well as the y-intercept and slope of the best-fit line as shown in your plot.