Transformations

Centering, standardising, correlation and logarithms

Dr Massimiliano Canzi

Transforming variables

Did you know: the blue whale is the largest animal known to have lived. The biggest accurately measured one, a female, was 29.5 m long and weighed 180 tonnes.

Sources: Encyclopaedia Britannica; WWF

Today

  • Why we transform variables
  • Linear transformations: units, centering, standardising
  • Correlation: Pearson’s r, and what it hides
  • Logarithms: for predictors and for outcomes
  • A note on the word “normalising”

Last week

tidy(m)[, 1:3]
# A tibble: 2 × 3
  term        estimate std.error
  <chr>          <dbl>     <dbl>
1 (Intercept) 137.      2.00    
2 body_mass_g   0.0153  0.000467

Each extra gram of body mass goes with 0.0153 mm more flipper, and a penguin of 0 g would have flippers 136.7 mm long.

The model is fine, but its numbers are awkward: a slope too small to picture, and an intercept that describes an impossible penguin. Today we fix both without touching the fit.

Why transform?

Goal Example Transformation
Coefficients that are easier to read slope per kg, intercept for an average penguin units, centering
Predictors on a common scale body mass (g) against bill length (mm) standardising
A shape the linear model can handle skewed reaction times, word frequencies logarithm

The data stay the same. A transformation changes the scale on which we describe them.

Linear vs non-linear

A linear transformation adds, subtracts, multiplies or divides by a constant. A non-linear one, such as the logarithm, changes the shape of the distribution.

Simulated reaction times on four scales. The three linear ones are the same histogram with different numbers on the axis. Only the log changes the shape.

Linear transformations

Did you know: a nursing blue whale calf gains up to 90 kg a day: that is 3.75 kg an hour, or 0.09 tonnes a day. Same calf, different units.

Source: Encyclopaedia Britannica

Changing units

pengs <- pengs |> mutate(body_mass_kg = body_mass_g / 1000)

coef(lm(flipper_length_mm ~ body_mass_g,  data = pengs))
 (Intercept)  body_mass_g 
136.72955927   0.01527592 
coef(lm(flipper_length_mm ~ body_mass_kg, data = pengs))
 (Intercept) body_mass_kg 
   136.72956     15.27592 

Divide the predictor by 1000, and the slope is multiplied by 1000: 15.3 mm per kilogram.

The line and its R² (0.76) stay the same; the numbers are just easier to read.

Where is the intercept?

The intercept is where the line crosses x = 0: far outside the data.

Centering

Subtract the mean, so that 0 means “a penguin of average body mass”:

pengs <- pengs |>
  mutate(body_mass_c = body_mass_g - mean(body_mass_g))

coef(lm(flipper_length_mm ~ body_mass_c, data = pengs))
 (Intercept)  body_mass_c 
200.91520468   0.01527592 
mean(pengs$flipper_length_mm)
[1] 200.9152

The slope does not change. The intercept is now the mean flipper length: the predicted value for an average penguin.

Centering

The data and the line are unchanged, and only the x axis has moved. The intercept (white) sits in the middle of the data.

Centre on any meaningful value

The mean is not the only choice. Subtract any value that makes 0 meaningful:

pengs <- pengs |> mutate(body_mass_4kg = body_mass_g - 4000)
coef(lm(flipper_length_mm ~ body_mass_4kg, data = pengs))
  (Intercept) body_mass_4kg 
 197.83322170    0.01527592 

The intercept is now the predicted flipper length of a 4 kg penguin.

In linguistic data:

Centred predictor The intercept describes
age - 18 an 18-year-old speaker
year - 2000 the year 2000
trial - 1 the first trial

Standardising

A z-score: subtract the mean, then divide by the standard deviation.

z = (x − mean) / SD
pengs <- pengs |>
  mutate(body_mass_z = (body_mass_g - mean(body_mass_g)) / sd(body_mass_g))

coef(lm(flipper_length_mm ~ body_mass_z, data = pengs))
(Intercept) body_mass_z 
  200.91520    12.25059 

The slope is now the change in flipper length for one standard deviation more body mass.

What a z-score says

A z-score says how many standard deviations a value lies from the mean: z = 0 is average, z = -2 is two SDs below it. From lecture 3, about 95% of values in a normal distribution lie between z = -2 and z = 2.

Comparing predictors

Which predicts flipper length more strongly: body mass (g) or bill length (mm)?

c(coef(lm(flipper_length_mm ~ body_mass_g,    data = pengs))[2],
  coef(lm(flipper_length_mm ~ bill_length_mm, data = pengs))[2])
   body_mass_g bill_length_mm 
    0.01527592     1.69006185 

Grams and millimetres cannot be compared. Standardise both. scale() computes z-scores, but returns a matrix, so wrap it in as.numeric():

pengs <- pengs |> mutate(bill_length_z = as.numeric(scale(bill_length_mm)))
c(coef(lm(flipper_length_mm ~ body_mass_z,   data = pengs))[2],
  coef(lm(flipper_length_mm ~ bill_length_z, data = pengs))[2])
  body_mass_z bill_length_z 
    12.250590      9.227034 

Per standard deviation, body mass goes with a larger change in flipper length.

Speaker normalisation

A classic linguistic use of z-scores: vowel formants.

Speakers with longer vocal tracts have lower formants overall. Raw F1 and F2 mix up which vowel and which speaker.

set.seed(8)
vowel_means <- tribble(
  ~vowel, ~F1, ~F2,
  "i",    300, 2300,
  "a",    750, 1300,
  "u",    320,  800
)
formants <- crossing(speaker = c("A", "B"), vowel_means, token = 1:20) |>
  mutate(k  = if_else(speaker == "A", 1, 1.2),   # speaker B: shorter vocal tract
         F1 = F1 * k + rnorm(n(), 0, 40),
         F2 = F2 * k + rnorm(n(), 0, 120))

Raw formants

The same vowels from two speakers land in different places.

Lobanov normalisation

z-score each formant within each speaker (Lobanov, 1971). With .by from lecture 2:

formants_z <- formants |>
  mutate(across(c(F1, F2), \(x) as.numeric(scale(x))),
         .by = speaker)

.by = speaker: each speaker’s F1 and F2 are z-scored using that speaker’s own mean and SD.

After normalisation

After normalisation, the vowels overlap across speakers.

Discussion

Lobanov normalisation removes the differences between speakers, so that vowels can be compared across them.

  • When could those differences be exactly what you want to study?
  • What would you lose by normalising in that case?

Linear transformations: summary

predictor 0 means intercept slope slope per R2
body_mass_g 0 g 136.730 0.015 1 g 0.76
body_mass_kg 0 kg 136.730 15.276 1 kg 0.76
body_mass_c average mass 200.915 0.015 1 g 0.76
body_mass_4kg 4 kg 197.833 0.015 1 g 0.76
body_mass_z average mass 200.915 12.251 1 SD 0.76

Linear transformations change what the numbers mean, while the fit stays the same.

Correlation

Did you know: the heart of a blue whale that washed ashore in Newfoundland, Canada, in 2014 weighed about 180 kg. It is now on display at the Royal Ontario Museum in Toronto.

Source: Royal Ontario Museum

Pearson’s r

The correlation coefficient r measures how closely two variables follow a straight line.

  • From -1 (perfect negative) through 0 (no linear relation) to +1 (perfect positive)
  • It has no units: it is the same whether you measure in g, kg or tonnes

cor() in R

cor(penguins$flipper_length_mm, penguins$body_mass_g)
[1] NA

Missing values again. Use only complete pairs:

cor(penguins$flipper_length_mm, penguins$body_mass_g, use = "complete.obs")
[1] 0.8712018

Per group, with summarise() and .by:

pengs |> summarise(r = cor(flipper_length_mm, body_mass_g), .by = species)
# A tibble: 3 × 2
  species       r
  <fct>     <dbl>
1 Adelie    0.468
2 Gentoo    0.703
3 Chinstrap 0.642

Correlation and regression

Standardise both variables, and the regression slope is the correlation:

pengs <- pengs |> mutate(flipper_z = as.numeric(scale(flipper_length_mm)))

coef(lm(flipper_z ~ body_mass_z, data = pengs))
  (Intercept)   body_mass_z 
-1.011386e-15  8.712018e-01 
cor(pengs$flipper_length_mm, pengs$body_mass_g)
[1] 0.8712018

And in simple regression, R² = r²:

cor(pengs$flipper_length_mm, pengs$body_mass_g)^2
[1] 0.7589925
glance(m)$r.squared
[1] 0.7589925

Same r, different data

Anscombe’s quartet (1973), built into R as anscombe:

Four data sets with the same means, the same regression line and the same correlation, r = 0.82. Only one of them looks like a straight-line relationship. Always plot first.

Correlation by group

Bill length and bill depth, from last week’s warning:

bills <- drop_na(penguins, bill_length_mm, bill_depth_mm)
cor(bills$bill_length_mm, bills$bill_depth_mm)
[1] -0.2350529
bills |> summarise(r = cor(bill_length_mm, bill_depth_mm), .by = species)
# A tibble: 3 × 2
  species       r
  <fct>     <dbl>
1 Adelie    0.391
2 Gentoo    0.643
3 Chinstrap 0.654

Negative overall, positive within every species: Simpson’s paradox again. A correlation describes only the data you give it, and says nothing about the process behind them.

Correlation is not causation

In children, shoe size and vocabulary size are correlated: children with bigger feet tend to know more words.

Bigger feet do not teach words. Both grow with a third variable: age.

A correlation, or a regression slope, tells us that two variables go together. Whether one causes the other is a question for the design of the study.

Discussion

Suppose children who are read to at home every day know more words than children who are not.

  • Does reading aloud cause vocabulary growth? What else could explain the correlation?
  • What kind of study would let you decide?

Logarithms

Did you know: blue whale calls can reach 188 decibels. Decibels are a logarithmic scale: every 10 dB more means ten times the sound intensity.

Source: WWF

Magnitudes

Some linguistic variables span many orders of magnitude. Word frequency is the classic example:

set.seed(5)
words <- tibble(rank = 1:2000,
                freq = round(1e6 / rank)) |>        # Zipf-like frequencies
  mutate(rt = exp(6.75 - 0.06 * log10(freq) + rnorm(n(), 0, 0.18)))

range(words$freq)
[1] 5e+02 1e+06

The most frequent word is 2000 times more frequent than the least frequent one. rt: a simulated lexical decision time for each word, in ms.

What is a logarithm?

A logarithm answers one question: how many times do I multiply the base by itself to get this number? log10(1000) is 3, because 10 × 10 × 10 = 1000.

So the log turns multiplying into adding: each step of 1 on the log scale is one more multiplication by the base.

log10(1000)        # 10^3 = 1000
[1] 3
log2(8)            # 2^3  = 8
[1] 3
log(exp(2))        # natural log, base e = 2.718...
[1] 2
exp(log(50))       # exp() undoes log()
[1] 50

Which base?

Function Base One unit means
log10() 10 a tenfold increase
log2() 2 a doubling
log() e ≈ 2.718 the default in R and in most models

The base changes only the numbers on the axis, and the shape stays the same: log10(x) is just log(x) / log(10). Choose the base that is easiest to explain.

Zeros and negatives

log(0)
[1] -Inf
log(-5)
[1] NaN

The log is only defined for positive numbers.

Frequency counts often contain zeros (a word that never occurs in a corpus). A common fix is to add 1 first:

log1p(c(0, 9, 99))     # same as log(c(0, 9, 99) + 1)
[1] 0.000000 2.302585 4.605170

Logging a predictor

On the raw scale, most words are squeezed against zero and a straight line fits badly. On the log scale, the relationship is close to linear.

Logging a predictor

m_freq <- lm(rt ~ log10(freq), data = words)
coef(m_freq)
(Intercept) log10(freq) 
  876.69485   -48.27585 

Each tenfold increase in frequency (10 → 100 → 1000) goes with responses about 48 ms faster.

Frequency effects in language are usually multiplicative: going from 10 to 100 matters about as much as going from 100,000 to 1,000,000.

Logging the outcome

Reaction times are positive and right-skewed (lecture 3). So are the residuals of a model fitted to raw RTs:

After the log, the residuals are much more symmetric: closer to what the linear model assumes (lecture 4).

Reading a log outcome

m_log <- lm(log(rt) ~ log10(freq), data = words)
coef(m_log)
(Intercept) log10(freq) 
 6.78462766 -0.06845576 
exp(coef(m_log)[2])
log10(freq) 
  0.9338348 

On the log scale, effects multiply. Each tenfold increase in frequency multiplies the RT by 0.934: responses about 6.6% faster.

Back to milliseconds with exp():

exp(predict(m_log, newdata = tibble(freq = c(10, 1000, 100000))))
       1        2        3 
825.6508 720.0066 627.8799 

When to log

  • The variable is positive and right-skewed: RTs, durations, frequencies, counts
  • Effects are plausibly multiplicative: “x% faster” rather than “x ms faster”
  • The residuals of the raw model are skewed or funnel-shaped

The log is not a cure-all. Report which variables you transformed and why, and interpret the coefficients on the transformed scale.

The word “normalising”

A slippery word

“Normalising” is used for several different operations:

What people mean What it does
z-scoring subtracts the mean and divides by the SD
min-max scaling squeezes values into the range 0 to 1
a normalising transformation makes a distribution look more normal, e.g. a log
speaker normalisation removes vocal tract differences in phonetics, e.g. Lobanov
x <- c(4, 8, 15, 16, 23, 42)
(x - min(x)) / (max(x) - min(x))     # min-max scaling
[1] 0.0000000 0.1052632 0.2894737 0.3157895 0.5000000 1.0000000

In your write-up, avoid the word on its own. Say exactly what you did: “F1 and F2 were z-scored within speaker (Lobanov, 1971)”. See also Winter (2019), section 5.7.

Summary

  • Units and centering change what the intercept and slope mean; the fit stays the same
  • Standardising puts variables on a common scale: slopes per SD
  • Correlation is the slope between two z-scored variables: unitless, from −1 to 1
  • The same r can come from very different data, so always plot
  • Logarithms turn multiplicative effects into additive ones, and tame skew
  • Name your transformation precisely, instead of just writing “normalised”

Slides and code

QR code linking to these slides

These slides

Code, data and all lectures