Pearson Correlation in R: Complete Step-by-Step Tutorial

Pearson's correlation is probably the first inferential technique you learn in a psychology degree and the one you'll use most throughout your academic career. It measures the strength and direction of the linear relationship between two continuous variables, returning a coefficient r between -1 and +1. Sounds simple, but mistakes happen constantly: skipping assumption checks, confusing correlation with causation, reporting results halfway. In this tutorial you'll learn to do it properly in R, from loading data to having your write-up ready for the manuscript.

R is particularly good for this because cor.test() gives you r, the p-value, and the confidence interval in a single call. No bootstrap tricks needed. Let's get to it.

If you end up comparing more than two variables at once, don't build the correlation table by hand: the APA 7 correlation table generator turns your matrix of r and p values into the table with means, SD and significance asterisks your journal expects.

When to use Pearson's correlation

Pearson is appropriate when you want to assess whether a linear relationship exists between two continuous (or at least interval-level) variables. For example: you want to know whether self-esteem and life satisfaction move together in a sample of 150 university students. That's Pearson. If one of your variables is ordinal (a single 5-point Likert item treated as a rank) or the relationship is monotonic but not linear, the alternative is Spearman, and I'll show you how to run it at the end of this tutorial.

Throughout, we'll use a concrete example: the relationship between self-esteem (Rosenberg Self-Esteem Scale, 10-40) and life satisfaction (Diener's SWLS, 5-35) in 150 students.

Assumptions of Pearson's correlation

Before computing anything, you need to verify that your data meet the assumptions. Skipping this is one of the most common mistakes in statistics and can invalidate your conclusions. Here's what to check.

Continuous variables

Both variables must be continuous (or at least interval-level). Total scores from validated psychometric scales like the Rosenberg or the SWLS are conventionally treated as continuous. If one of your variables is genuinely ordinal (a single Likert item, for instance), use Spearman.

Linear relationship

Pearson only captures linear relationships. If the association between your variables is U-shaped or curvilinear, r can come out close to zero even when there's a strong association. The most reliable check is a scatter plot. If the point cloud follows a reasonably straight trend, you're good.

Bivariate normality

Strictly speaking, Pearson assumes bivariate normality (that the joint distribution of both variables is normal). In practice, with samples of 30 or more, the significance test is fairly robust to moderate violations thanks to the central limit theorem. Evaluate each variable's normality separately using histograms or Q-Q plots. If either shows extreme skewness or severe kurtosis, consider a transformation or switch to Spearman. The article on how to verify statistical assumptions covers this in depth.

No extreme outliers

Outliers have a disproportionate effect on Pearson's r because the coefficient is based on means and standard deviations, which are sensitive to extreme values. A single atypical case can dramatically inflate or deflate r. Inspect the scatter plot for points that clearly break away from the cloud. If you find any, figure out whether they're coding errors, belong to a different population, or are genuine values, and report results both with and without those cases.

Pearson correlation in R: the full workflow

Time for code. We'll walk through the entire workflow: loading data, exploring, checking assumptions visually, computing the correlation, and visualizing the matrix.

Step 1: Load and explore the data

First, import your data and take a quick look. If you're working with a CSV (the most common scenario), here's what you need:

# Load data
data <- read.csv("self_esteem_data.csv")

# First rows: make sure columns look right
head(data)

# Basic descriptives
summary(data[, c("self_esteem", "life_satisfaction")])

# Standard deviations (summary doesn't include them)
sapply(data[, c("self_esteem", "life_satisfaction")], sd, na.rm = TRUE)

Check for out-of-range values. Self-esteem should be between 10 and 40 (Rosenberg) and life satisfaction between 5 and 35 (SWLS). If you see a 99 or a -1, you have a coding problem.

Step 2: Scatter plot with ggplot2

The scatter plot is your main tool for checking linearity and spotting outliers. Don't skip this step.

library(ggplot2)

ggplot(data, aes(x = self_esteem, y = life_satisfaction)) +
  geom_point(alpha = 0.6, color = "#2f6b57", size = 2) +
  geom_smooth(method = "lm", se = TRUE, color = "#101a14",
              linewidth = 1) +
  labs(
    x = "Self-Esteem (Rosenberg)",
    y = "Life Satisfaction (SWLS)",
    title = "Relationship between self-esteem and life satisfaction"
  ) +
  theme_minimal(base_size = 13)

What you're looking for: a point cloud that follows a reasonably straight trend (linearity) with no isolated points far from the rest (outliers). The regression line and the grey confidence band give you a quick visual reference. If the cloud looks curved, Pearson isn't your test.

Step 3: Check normality

Histograms and Q-Q plots let you assess each variable's normality separately. With ggplot2 and base R you can generate both easily:

# Histograms
library(gridExtra)

p1 <- ggplot(data, aes(x = self_esteem)) +
  geom_histogram(bins = 20, fill = "#2f6b57", color = "white",
                 alpha = 0.8) +
  labs(title = "Self-Esteem", x = NULL) +
  theme_minimal()

p2 <- ggplot(data, aes(x = life_satisfaction)) +
  geom_histogram(bins = 20, fill = "#2f6b57", color = "white",
                 alpha = 0.8) +
  labs(title = "Life Satisfaction", x = NULL) +
  theme_minimal()

grid.arrange(p1, p2, ncol = 2)

# Q-Q plots (base R, quick)
par(mfrow = c(1, 2))
qqnorm(data$self_esteem, main = "Q-Q Self-Esteem")
qqline(data$self_esteem, col = "#2f6b57", lwd = 2)
qqnorm(data$life_satisfaction, main = "Q-Q Life Satisfaction")
qqline(data$life_satisfaction, col = "#2f6b57", lwd = 2)
par(mfrow = c(1, 1))

If the Q-Q points line up reasonably well with the diagonal, normality is acceptable. You don't need perfect normality, just no extreme skewness or heavy tails.

Step 4: Compute the correlation with cor.test()

Here's the core of it. cor.test() is your function:

# Pearson correlation with 95% CI
result <- cor.test(data$self_esteem, data$life_satisfaction,
                   method = "pearson")
result

The output gives you everything in one block:

# Pearson's product-moment correlation
#
# data:  data$self_esteem and data$life_satisfaction
# t = 6.48, df = 148, p-value = 1.2e-09
# alternative hypothesis: true correlation is not equal to 0
# 95 percent confidence interval:
#  0.3345  0.5872
# sample estimates:
#       cor
# 0.4700

Let's break this down:

• r = .47: moderate-to-large positive correlation.
• t = 6.48, df = 148: the t statistic and degrees of freedom (N - 2).
• p-value = 1.2e-09: that's p < .001. The correlation is significant.
• 95% CI [.33, .59]: the range of plausible values for the population correlation.

If you need to extract individual values programmatically for reporting:

# Extract components
result$estimate   # r
result$p.value    # p-value
result$conf.int   # 95% CI
result$statistic  # t
result$parameter  # degrees of freedom

# r-squared
result$estimate^2  # proportion of shared variance

Step 5: Directional hypothesis (one-tailed)

If your hypothesis is directional (you predicted a positive correlation and justified this a priori), you can request a one-tailed test:

# One-tailed: hypothesis of positive correlation
cor.test(data$self_esteem, data$life_satisfaction,
         method = "pearson",
         alternative = "greater")

The p-value will be half the two-tailed value. Only use this if you genuinely had a directional hypothesis before looking at the data. Deciding on the direction after seeing results is HARKing, and that's not okay.

Step 6: Correlation matrix

When you have more than two variables and want to explore all bivariate correlations at once, use cor() for the matrix and corrplot to visualize it:

# Select variables
variables <- data[, c("self_esteem", "life_satisfaction",
                       "anxiety", "social_support")]

# Correlation matrix
cor_matrix <- cor(variables, use = "pairwise.complete.obs")
round(cor_matrix, 2)

The cor() function gives you the matrix but no p-values. To get a table with both r and p for each pair, the Hmisc package is very handy:

library(Hmisc)

# Matrix with r and p-values
rcorr(as.matrix(variables))

And for the visual version, corrplot is the standard:

library(corrplot)

corrplot(cor_matrix, method = "color", type = "upper",
         addCoef.col = "black", tl.col = "black",
         tl.srt = 45, diag = FALSE,
         col = colorRampPalette(c("#e74c3c", "white",
                                   "#2f6b57"))(200))

The heat map lets you spot patterns quickly: clusters of strongly intercorrelating variables may suggest latent factors. Much more readable than a table of numbers once you have five or more variables.

Bonus: Spearman as an alternative

If your data don't meet Pearson's assumptions (the relationship isn't linear, there are severe outliers, or a variable is ordinal), Spearman is your alternative. The R syntax is identical, just change the method argument:

# Spearman correlation
cor.test(data$self_esteem, data$life_satisfaction,
         method = "spearman")

Spearman works with ranks instead of raw values, making it resistant to outliers and free from the linearity requirement. If Pearson and Spearman give you similar results, that's a good sign: your correlation is robust.

Interpreting the results

You have the R output in front of you. Now you need to interpret it properly. Four key elements. (And if you would rather check it against your own case, copy what your console returns and it breaks it down element by element.)

The r coefficient

The r value tells you the strength and direction of the linear relationship. An r = .47 means a positive relationship: higher self-esteem goes with higher life satisfaction. For magnitude, Cohen's (1988) conventions are the standard in psychology:

• r = .10: small effect
• r = .30: medium effect
• r = .50: large effect

Our r = .47 falls in the medium-to-large range. But don't apply these mechanically: in some areas of psychology, correlations of .20 are substantively important, while in others .50 might be modest. Your research context matters.

The p-value

The p-value tells you the probability of getting an r this extreme (or more) if the population correlation were zero. A p < .001 means it's extremely unlikely the observed correlation is due to chance. But watch out: the p-value depends heavily on sample size. With N = 500, correlations as low as .09 are already significant, and that doesn't mean they're practically relevant.

The confidence interval

The 95% CI (in our example, [.33, .59]) tells you the range of plausible values for the population correlation. It's far more informative than the p-value alone because it shows the precision of your estimate. A narrow CI means precise estimation; a wide one says you need more data. The article on confidence intervals goes deeper into this.

The coefficient of determination r-squared

Squaring r gives you the proportion of shared variance between the two variables. If r = .47, then r-squared = .22, meaning self-esteem accounts for roughly 22% of the variability in life satisfaction (and vice versa). This is crucial for judging practical relevance. An r = .20 looks modest, but r-squared = .04 tells you only 4% of the variance is shared. That puts things in perspective.

Correlation matrix: how to read it

When you analyze more than two variables, the output is a correlation matrix: a symmetric table where each cell holds the r between that pair of variables. The diagonal is always 1 (each variable correlates perfectly with itself) and the lower half mirrors the upper half.

When reading a matrix, focus on three things: (1) the magnitude of the coefficients, to identify which pairs are most strongly associated; (2) the direction, positive or negative; and (3) overall patterns, like clusters of variables that strongly intercorrelate (which may suggest latent factors). If your matrix has many variables, the corrplot heat map we built earlier is far more readable than a table of numbers.

APA 7 reporting format

APA 7 requires reporting r, degrees of freedom, the p-value, and ideally the confidence interval. Here's the structure:

Formatting details that reviewers always catch: the r coefficient is italicized and reported without a leading zero (.47, not 0.47). Degrees of freedom go in parentheses after r and equal N - 2. The p-value is also italicized. For a comprehensive guide, check the APA 7 reporting tutorial.

A handy R trick: you can build the APA text directly in code so you don't mistype values:

# Generate APA text automatically
r_val <- round(result$estimate, 2)
df_val <- result$parameter
p_val <- result$p.value
ci_low <- round(result$conf.int[1], 2)
ci_high <- round(result$conf.int[2], 2)

cat(sprintf("r(%d) = %.2f, p < .001, 95%% CI [%.2f, %.2f]",
            df_val, r_val, ci_low, ci_high))

Common mistakes with Pearson's correlation

Confusing correlation with causation

Self-esteem and life satisfaction correlating at .47 doesn't mean self-esteem causes life satisfaction. The relationship could run in reverse, be bidirectional, or be driven by a third variable (like social support). Establishing causation requires experimental designs or, at minimum, longitudinal studies with proper controls. The article on mediation and moderation explores this further.

Failing to check linearity

Running Pearson without first looking at a scatter plot is like driving with your eyes closed. If the relationship is curvilinear, r will severely underestimate the true association. Always inspect the plot first.

Using leading zeros when reporting r and p

In APA format, statistics that cannot exceed 1 in absolute value (r, p, standardized betas) are reported without a leading zero. It's .47, not 0.47. It's p < .001, not p < 0.001. Minor detail, but reviewers notice it instantly.

Ignoring effect size

Reporting only significance without contextualizing it with r-squared or Cohen's conventions is incomplete. An r = .12 may be significant with N = 300, but r-squared = .014 means barely 1.4% of the variance is shared. Significance tells you the effect probably isn't zero; effect size tells you whether it matters.

Dichotomizing continuous variables

Converting a continuous variable into two groups (high vs low self-esteem via median split) to run a t-test instead of a correlation drastically reduces statistical power and wastes information. If both variables are continuous, use the correlation.

If you're running correlations for your thesis and need help checking assumptions, interpreting the output, or writing the results section in APA format, my statistical consulting service can help. We go through the analysis together and I deliver the write-up ready for your manuscript.