airquality
This material has been reproduced and communicated to you by or on behalf of Kaplan Business School pursuant to Part VB of the Copyright Act 1968 (the Act).
The material in this communication may be subject to copyright under the Act. Any further reproduction or communication of this material by you may be the subject of copyright protection under the Act.
Do not remove this notice.
| Week | Topic |
|---|---|
| 1–4 | Variable types. Descriptive statistics. Histograms, boxplots, scatter plots, correlation. |
| 5 | Group presentations. Linear trend line. |
| 6 | Visualisation for time series. |
| 7 | Visualisation for multivariate data. |
| 8 | Data frames. Missing values: imputation and interpolation. Data cleaning. |
| 9 | Aggregates and joining tables. |
| 10–11 | Data visualisation with ggplot. |
| 12 | Individual presentations. |
airquality and the question of what
cor(airquality) returns before you do anything about the missing values. This lesson
answers it, and then works out what to do next.NA means in R and predict how it propagates through a calculation.airquality ships with R. Daily air-quality measurements from New
York, 1 May to 30 September 1973.
| Column | Meaning | Missing |
|---|---|---|
| Ozone | Mean ozone in parts per billion, 13:00 to 15:00 | 37 of 153 |
| Solar.R | Solar radiation in langleys, 08:00 to 12:00 | 7 of 153 |
| Wind | Average wind speed in miles per hour | 0 |
| Temp | Maximum daily temperature in degrees Fahrenheit | 0 |
| Month | Month, 5 to 9 | 0 |
| Day | Day of month, 1 to 31 | 0 |
data(airquality)
aq <- airquality
dim(aq)
#> [1] 153 6A data frame is the object every dataset in this unit lives in. Five operations cover almost everything you will do to one, and all five exist in base R.
Rows are observations. Columns are variables. Each column holds one type, and that is what separates a data frame from a spreadsheet.
students <- data.frame(
name = c("John Smith", "Jane Doe", "Joe Schmo"),
age = c(34, 28, 51),
address = c("123 Main St.", "456 Maple Ave.",
"789 Broadway"),
stringsAsFactors = FALSE
)df <- read.csv("my_data.csv")
write.csv(df, "new_data.csv", row.names = FALSE)
# read.csv arguments worth knowing
read.csv("f.csv",
stringsAsFactors = FALSE,
na.strings = c("", "NA", "-999"),
strip.white = TRUE)install.packages(c("readr", "dplyr", "tidyr"))
library(readr)
library(dplyr)
df <- read_csv("my_data.csv")
write_csv(df, "new_data.csv")read_csv returns a tibble, prints ten rows, and never converts text to a
factor.read.csv returns a plain data frame and is available with no installation.na or na.strings for custom missing markers.| Function | What it answers |
|---|---|
| dim(aq) | How many rows and columns? |
| str(aq) | What type is each column, and what do the first few values look like? |
| head(aq, 10) | Does the top of the file look like data or like a stray header? |
| summary(aq) | Range, quartiles, and the count of NAs per column. |
| colnames(aq) | Are the names usable, or full of spaces? |
| colSums(is.na(aq)) | Exactly how many values are missing where? |
summary(aq$Ozone)
#> Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
#> 1.00 18.00 31.50 42.13 63.25 168.00 37
The final column appears only when the variable contains missing values. It is the fastest missingness check in R, and it comes with the summary you wanted anyway.
str() and summary() before any analysis takes ten seconds and
catches the numeric column that arrived as text, the date that arrived as a factor, and the
column that is entirely NA.| Operation | Base R |
|---|---|
| select | aq[, c("Ozone", "Temp")] |
| filter | aq[aq$Month == 5, ]subset(aq, Month == 5) |
| arrange | aq[order(-aq$Temp), ] |
| mutate | aq$TempC <- (aq$Temp - 32) * 5/9transform(aq, TempC = ...) |
| summarise | aggregate(Ozone ~ Month, aq, mean) |
aq[aq$Month == 5, ] keeps rows where the test is NA and fills them with
NA. subset() drops them. The two give different row counts on incomplete data, so
choose one deliberately.R 4.1 added a pipe to the language itself. It needs no package.
# without a pipe: nested, read inside out
head(aq[order(-aq$Temp), c("Ozone", "Temp")], 5)
# with the native pipe: read left to right
aq |>
subset(Month == 7) |>
transform(hot = Temp > 85) |>
head(3)
# x |> f(y) is exactly f(x, y)
# The value on the left becomes the FIRST argument.
The magrittr pipe %>% predates the native one and behaves almost
identically for these purposes.
library(dplyr)
aq %>%
filter(Month == 7) %>%
select(Ozone, Temp) %>%
arrange(desc(Temp)) %>%
mutate(hot = Temp > 85) %>%
head(3)
# Differences worth knowing:
# %>% allows the placeholder . anywhere
# |> puts the value first, and nowhere else
# %>% needs magrittr or dplyr loaded
# |> is part of the language
Tidy means one variable per column and one observation per row. Neither layout below is more correct than the other; each is tidy for a different question.
sum(duplicated(df)) # how many repeats
df[!duplicated(df), ] # keep the first of each
unique(df) # the same thing
# duplicated on chosen columns only
df[!duplicated(df[, c("id", "date")]), ]
str(df) # what R thinks each is
df$price <- as.numeric(df$price)
df$grp <- as.factor(df$grp)
df$date <- as.Date(df$date, "%d/%m/%Y")
p <- c("$1,200", "$980", "$1,050")
as.numeric(gsub("[$,]", "", p))
#> [1] 1200 980 1050
# splitting a fixed-width field
bd <- c("03151998", "11021985")
substr(bd, 1, 2) #> "03" "11"
substr(bd, 3, 4) #> "15" "02"
substr(bd, 5, 8) #> "1998" "1985"
# splitting on a separator
strsplit("guest_AU", "_")[[1]]
#> [1] "guest" "AU"
as.numeric() on
"$1,200" returns NA with a warning, and you have just created a
missing value that was never missing.aq[aq$Ozone > 100, ] and subset(aq, Ozone > 100) return different numbers of rows. Why?"$1,200". What does as.numeric(price) return?Before you fill a single value, find out how much is missing, where it sits, and whether the pattern is related to anything you can observe.
"" is a value of length zero. It is a string that R will happily sort, count and paste.NULL has length zero and disappears from a vector. NA occupies a position and keeps the length.NaN is the result of an undefined calculation such as 0/0. is.na(NaN) is TRUE, but is.nan(NA) is FALSE.na.rm = TRUE drops it and silently reduces n.x == NA returns NA, never TRUE. The only test that works is is.na(x).x <- c(41, 36, NA, 18)
mean(x) #> NA
mean(x, na.rm = TRUE) #> 31.66667
length(x) #> 4
sum(!is.na(x)) #> 3
is.na(x) #> F F T F
x == NA #> NA NA NA NA
# Distinguishing the four
is.na(NA); is.na(NaN) #> TRUE TRUE
is.nan(NA); is.nan(NaN) #> FALSE TRUE
length(c(1, NULL, 3)) #> 2
length(c(1, NA, 3)) #> 3
# na.rm changes the denominator. Report it.
mean(aq$Ozone, na.rm = TRUE) #> 42.12931
sum(!is.na(aq$Ozone)) #> 116# A missingness matrix in one call.
# t() and the reversed index put day 1 at the top.
image(1:ncol(aq), 1:nrow(aq),
t(is.na(aq))[, nrow(aq):1],
col = c("#dfeaf2", "#CC0000"),
axes = FALSE, xlab = "", ylab = "Day")
axis(1, 1:ncol(aq), colnames(aq), las = 2)
box()
# The same information as text
colSums(is.na(aq))
#> Ozone Solar.R Wind Temp Month Day
#> 37 7 0 0 0 0
# Read this chart BEFORE you decide anything.
# A block and a scatter of single days call for
# completely different treatments.# cells
sum(is.na(aq)) #> 44
prod(dim(aq)) #> 918
100 * 44 / 918 #> 4.793028
# rows
sum(!complete.cases(aq)) #> 42
100 * mean(!complete.cases(aq)) #> 27.45098
# per column
colSums(is.na(aq))
round(100 * colMeans(is.na(aq)), 1)
#> Ozone Solar.R Wind Temp Month Day
#> 24.2 4.6 0.0 0.0 0.0 0.0
barplot(colSums(is.na(aq)), las = 2,
col = "#CC0000", ylab = "Missing values")
# The gap between 4.8% and 27.5% is the whole
# argument against reaching for complete.cases()
# as a first move.table(Ozone = is.na(aq$Ozone),
Solar.R = is.na(aq$Solar.R))
#> Solar.R
#> Ozone FALSE TRUE
#> FALSE 111 5
#> TRUE 35 2
# Expected count under independence
37 * 7 / 153 #> 1.692810
# Which days lose both?
which(is.na(aq$Ozone) & is.na(aq$Solar.R))
#> [1] 5 27
# For more than two columns, count the distinct
# patterns directly.
pat <- apply(is.na(aq), 1, paste, collapse = "")
sort(table(pat), decreasing = TRUE)| Name | What it means | Example in this dataset | Consequence |
|---|---|---|---|
| MCAR | Missing completely at random. Whether a value is missing has nothing to do with anything, observed or not. | An instrument that fails on days chosen by a coin toss. | Deleting incomplete rows is unbiased. You lose precision only. |
| MAR | Missing at random. Whether a value is missing depends on other variables you did record. | Ozone missing far more often in June than in any other month. | Deletion is biased. Imputation that conditions on the observed variable can recover the structure. |
| MNAR | Missing not at random. Whether a value is missing depends on the value itself. | A monitor that saturates and records nothing on the highest-ozone days. | No method in this lesson fixes it. You need the reason for the absence. |
Naming the mechanism is a claim about the world, not a calculation. Say which one you are assuming and why.
# Does missingness depend on an observed variable?
tapply(is.na(aq$Ozone), aq$Month, sum)
#> 5 6 7 8 9
#> 5 21 5 5 1
table(aq$Month)
#> 5 6 7 8 9
#> 31 30 31 31 30
# Compare the other columns on missing vs
# observed days
obs <- !is.na(aq$Ozone)
round(rbind(
observed = colMeans(aq[ obs, 2:4], na.rm = TRUE),
missing = colMeans(aq[!obs, 2:4], na.rm = TRUE)), 2)
#> Solar.R Wind Temp
#> observed 184.80 9.86 77.87
#> missing 189.51 10.26 77.92
# Solar.R, Wind and Temp barely move.
# Month moves a great deal. That is the signal.
boxplot(Temp ~ is.na(Ozone), data = aq,
names = c("recorded", "missing"),
ylab = "Temperature (F)")This is the question Week 7 finished on. The default is to refuse.
cor(aq)["Ozone", "Temp"]
#> [1] NA
cor(aq, use = "complete.obs")["Ozone", "Temp"]
#> [1] 0.6985
cor(aq, use = "pairwise.complete.obs")["Ozone", "Temp"]
#> [1] 0.6984
sum(complete.cases(aq)) #> 111
sum(!is.na(aq$Ozone) & !is.na(aq$Temp)) #> 116
| use = | Behaviour |
|---|---|
| "everything" | The default. Any NA in a pair returns NA for that pair. |
| "complete.obs" | Drops any row with an NA anywhere, then computes every correlation on the same 111 rows. |
| "pairwise.complete.obs" | Uses whatever rows each pair has. Ozone against Temp gets 116 rows; Wind against Temp gets all 153. |
complete.obs and report the n.Here the two settings agree to three decimal places, which happens when the extra 5 rows resemble the other 111. That agreement is a result on this dataset rather than a rule.
x <- c(4, NA, 6). Which expression returns c(FALSE, TRUE, FALSE)?Filling a gap adds no information. Every method below trades one distortion for another, and the trade is your decision to defend.
What fraction of cells, and what fraction of rows? A column that is 60 per cent missing is not a candidate for imputation; it is a candidate for deletion.
Scattered single values, or sustained blocks? A block cannot be filled from its neighbours, because it has none nearby.
Does missingness depend on an observed variable, on the unobserved value itself, or on neither? This is the mechanism from Section 2.
A mean, a correlation, a chart, a model? Mean imputation preserves the mean and damages every correlation. The right method depends on which quantity you are about to report.
Reporting a statistic on the complete cases, with the n stated, is a valid answer. It is the baseline every other method has to beat.
Add a logical column recording which values you filled. It costs one column and it makes every downstream result auditable.
Drops every row containing a missing value anywhere, then proceeds as though the data were complete.
cc <- aq[complete.cases(aq), ]
nrow(cc) #> 111
nrow(aq) #> 153
# na.omit() is the same thing
nrow(na.omit(aq)) #> 111
# Restricted to the columns you actually need
keep <- complete.cases(aq[, c("Ozone", "Temp")])
sum(keep) #> 116
complete.cases() to the columns your analysis uses, not to the whole frame.
On this dataset that recovers 5 days for an Ozone against Temp analysis, at no cost.| Rows kept | 111 of 153 |
| Rows lost | 42, or 27.45 per cent |
| June days kept | 9 of 30 |
A value filled at \(\bar{x}\) contributes \((\bar{x}-\bar{x})^2 = 0\) to the sum and 1 to \(n\). Filling \(m\) gaps therefore multiplies \(s^2\) by \(\frac{n-1}{n+m-1}\), which is always below 1. Here \(n = 116\) and \(m = 37\), giving a factor of 0.757 on the variance and 0.870 on the standard deviation.
oz <- aq$Ozone
# mean imputation
mi <- oz
mi[is.na(mi)] <- mean(oz, na.rm = TRUE)
# median imputation, better for a skewed column
mdi <- oz
mdi[is.na(mdi)] <- median(oz, na.rm = TRUE)
c(observed = sd(oz, na.rm = TRUE),
mean_imp = sd(mi),
med_imp = sd(mdi))
#> observed mean_imp med_imp
#> 32.98788 28.69214 29.05375
mean(mi) #> 42.12931 (unchanged)
mean(mdi) #> 39.55556 (pulled towards the median)
# Why the sd must fall:
# var = sum((x - xbar)^2) / (n - 1)
# each filled value adds 0 to the numerator
# and 1 to n, so the ratio can only shrink.If missingness depends on a variable you observed, condition on it. Fill each gap with the mean of its own group rather than the mean of everything.
# ave() computes a group statistic and returns it
# aligned to the original rows
grp <- ave(oz, aq$Month,
FUN = function(v) mean(v, na.rm = TRUE))
gi <- oz
gi[is.na(gi)] <- grp[is.na(gi)]
round(tapply(oz, aq$Month, mean, na.rm = TRUE), 2)
#> 5 6 7 8 9
#> 23.62 29.44 59.12 59.96 31.45
c(mean = mean(gi), sd = sd(gi))
#> mean sd
#> 40.85098 29.59461
fit <- lm(Ozone ~ Solar.R + Wind + Temp, data = aq)
summary(fit)$r.squared #> 0.6059
pr <- predict(fit, newdata = aq)
ri <- oz
k <- is.na(ri) & !is.na(pr)
ri[k] <- pr[k]
sum(k) #> 35 filled
sum(is.na(ri)) #> 2 still missing
sum(ri < 0, na.rm = TRUE) #> 2 impossible
# Two problems, two fixes:
# 1. The model needs Solar.R, which is itself
# missing on 2 of those days.
# 2. Nothing constrains the prediction to the
# physically possible range.
ri <- pmax(ri, min(oz, na.rm = TRUE))
# Stochastic regression restores the scatter by
# adding a residual drawn from the model.
set.seed(3100)
sri <- oz
sri[k] <- pr[k] + rnorm(sum(k), 0, sigma(fit))Fill each gap with a value drawn from the observed values of the same column. The filled column then has the right marginal distribution.
set.seed(3100)
hi <- oz
hi[is.na(hi)] <- sample(oz[!is.na(oz)],
size = sum(is.na(oz)),
replace = TRUE)
c(observed = sd(oz, na.rm = TRUE), hotdeck = sd(hi))
#> observed hotdeck
#> 32.98788 32.22754
cor(hi, aq$Temp) #> 0.5683
cor(oz, aq$Temp, use = "complete.obs") #> 0.6984
# Nearest-neighbour hot deck: draw from donors
# that resemble the recipient
donors <- aq$Temp > 80 & !is.na(aq$Ozone)
sample(aq$Ozone[donors], 1)
| Spread | Correlation | |
|---|---|---|
| Observed | 32.99 | 0.698 |
| Mean | 28.69 | 0.609 |
| Hot deck | 32.23 | 0.568 |
Hot deck preserves the spread better than any other method here, and damages the correlation with Temperature more than any other method here. The two facts have the same cause: a randomly drawn value has a realistic magnitude and no relationship to the day it lands on.
Hot deck is random, so it needs set.seed() to be reproducible.
Report the seed.
comp <- data.frame(
method = c("listwise", "mean", "median",
"month", "regression", "hot deck"),
n = c(116, 153, 153, 153, 151, 153),
sd = round(c(sd(oz, na.rm = TRUE), sd(mi), sd(mdi),
sd(gi), sd(ri, na.rm = TRUE), sd(hi)), 2),
r = round(c(
cor(oz, aq$Temp, use = "complete.obs"),
cor(mi, aq$Temp), cor(mdi, aq$Temp),
cor(gi, aq$Temp),
cor(ri, aq$Temp, use = "complete.obs"),
cor(hi, aq$Temp)), 3))
comp
#> method n sd r
#> 1 listwise 116 32.99 0.698
#> 2 mean 153 28.69 0.609
#> 3 median 153 29.05 0.601
#> 4 month 153 29.60 0.646
#> 5 regression 151 30.57 0.717
#> 6 hot deck 153 32.23 0.568
# Build this table for your own data before you
# choose. It takes ten lines and it is the whole
# justification for whichever row you pick.One extra column turns an unverifiable dataset into an auditable one.
aq2 <- aq
aq2$Ozone_imputed <- is.na(aq2$Ozone)
aq2$Ozone <- mi
sum(aq2$Ozone_imputed) #> 37
# Now every downstream result can be checked
# with and without the filled values
cor(aq2$Ozone[!aq2$Ozone_imputed],
aq2$Temp[!aq2$Ozone_imputed]) #> 0.6984
cor(aq2$Ozone, aq2$Temp) #> 0.6087
# And every chart can show them
plot(aq2$Temp, aq2$Ozone,
pch = ifelse(aq2$Ozone_imputed, 1, 19),
col = ifelse(aq2$Ozone_imputed,
"#CC0000", "#0072B2"))
legend("topleft", bty = "n",
pch = c(19, 1), col = c("#0072B2", "#CC0000"),
legend = c("recorded", "imputed"))
| How many | 37 of 153 Ozone values, 24.2 per cent |
| Which method | Monthly mean |
| Why that one | Missingness depends on Month, and the monthly means differ by a factor of 2.5 |
| What it cost | Standard deviation falls from 32.99 to 29.60 |
| What is still missing | None for Ozone. 7 Solar.R values were left as NA. |
| Seed, if random | set.seed(3100) |
When rows have an order, the neighbours of a gap are informative. That is the one situation where filling a value can be defended on evidence rather than on convenience.
plot(aq$Ozone, type = "o", pch = 19, cex = 0.5,
col = "#0072B2", xlab = "Day index",
ylab = "Ozone (ppb)")
# where the gaps are
gap <- which(is.na(aq$Ozone))
abline(v = gap, col = "#CC000022", lwd = 2)
# how long each run of missing days is
r <- rle(is.na(aq$Ozone))
r$lengths[r$values]
#> [1] 1 1 3 6 1 2 2 10 1 1 1 2 2 1 1 1 1
max(r$lengths[r$values]) #> 10
# where the longest run starts
start <- cumsum(c(1, head(r$lengths, -1)))
start[r$values & r$lengths == 10] #> 52
aq[52, c("Month", "Day")] #> 6, 21
# Interpolation is only available because the rows
# have an order. Shuffle them and every method in
# this section becomes meaningless.LOCF holds the most recent recorded value until a new one arrives. NOCB does the same backwards.
locf <- function(v) {
for (i in seq_along(v)[-1])
if (is.na(v[i])) v[i] <- v[i - 1]
v
}
nocb <- function(v) {
for (i in rev(seq_along(v))[-1])
if (is.na(v[i])) v[i] <- v[i + 1]
v
}
x <- c(NA, NA, 30, 34, NA, NA, 41, 38, NA)
locf(x) #> NA NA 30 34 34 34 41 38 38
nocb(x) #> 30 30 30 34 41 41 41 38 NA
nocb(locf(x)) #> 30 30 30 34 34 34 41 38 38
oz <- aq$Ozone
f <- locf(oz)
sum(is.na(f)) #> 0
# Ozone starts and ends with a recorded value,
# so LOCF alone leaves nothing missing here.
aq$Ozone[1] #> 41
aq$Ozone[153] #> 20
# The ten-day run from index 52
round(f[50:63])
# ten identical values in the middle of a series
# whose day-to-day variation is large
The controlled test on slide 4.6 measures it. On a smooth seasonal series LOCF reaches a mean absolute error of 3.89 degrees, against 3.43 for linear interpolation. It is the second best of the four methods tested, and it is worse than linear on exactly the long runs where the difference matters most.
r <- rle(is.na(oz))
ok <- rep(r$lengths <= 2, r$lengths)
f2 <- oz
f2[ok] <- locf(oz)[ok]idx <- seq_along(oz)
k <- !is.na(oz)
# straight line between neighbouring points
lin <- approx(idx[k], oz[k], xout = idx, rule = 2)$y
sum(is.na(lin)) #> 0
range(lin) #> 1 168
range(oz, na.rm = TRUE) #> 1 168
# Linear interpolation can never leave the range
# of the observed data, because every filled value
# is a weighted average of two real ones.
# approx() also does the reverse: resampling a
# series onto a new grid
approx(idx[k], oz[k], xout = seq(1, 153, by = 7))
# method = "constant" turns approx() into LOCF
approx(idx[k], oz[k], xout = idx,
method = "constant", rule = 2)$yspl <- spline(idx[k], oz[k], xout = idx)$y
range(spl) #> -102.4 261.1
range(oz, na.rm = TRUE) #> 1 168
sum(spl < 0) #> 7
# Enforcing the physical bound afterwards
spl2 <- pmax(spl, min(oz, na.rm = TRUE))
# A monotone spline cannot overshoot between
# points, which removes most of the damage
spl3 <- spline(idx[k], oz[k], xout = idx,
method = "hyman")$y
range(spl3)
# splinefun() returns a function you can evaluate
# anywhere, rather than a fixed set of values
f <- splinefun(idx[k], oz[k])
f(52:61)
# Use a spline when the underlying quantity really
# is smooth and the gaps are short. Ozone over ten
# days is neither.rule argument.rule = 1 returns NA and rule = 2 carries the nearest recorded value outward as a flat line.rule = 2 does it by assuming the series stopped changing. On airquality both settings give the same answer, because Ozone begins at 41 and ends at 20, with no missing values at either end. Check that before you rely on it.x <- c(NA, NA, 30, 34, NA, NA, 41, 38, NA)
i <- seq_along(x); k <- !is.na(x)
approx(i[k], x[k], xout = i, rule = 1)$y
#> NA NA 30.0 34.0 36.3 38.7 41.0 38.0 NA
approx(i[k], x[k], xout = i, rule = 2)$y
#> 30.0 30.0 30.0 34.0 36.3 38.7 41.0 38.0 38.0
# On airquality the distinction does not bite:
is.na(aq$Ozone[1]) #> FALSE
is.na(aq$Ozone[153]) #> FALSE
sum(is.na(approx(idx[k], oz[k],
xout = idx, rule = 1)$y)) #> 0
# rule can take two values, one per end:
# rule = c(1, 2) -> NA before, carried after
# Default to rule = 1 and decide about the ends
# on purpose. An NA you can see is safer than a
# fabricated value you cannot.# Temperature is complete, so we can hide values
# we already know and score the reconstruction.
tp <- aq$Temp
gaps <- c(20:24, 55:57, 88:95, 120:121)
y <- tp
y[gaps] <- NA
k <- !is.na(y)
L <- approx(idx[k], y[k], xout = idx, rule = 2)$y
S <- spline(idx[k], y[k], xout = idx)$y
F <- locf(y)
M <- y; M[is.na(M)] <- mean(y, na.rm = TRUE)
sapply(list(mean = M, locf = F, linear = L, spline = S),
function(z) mean(abs(z[gaps] - tp[gaps])))
#> mean locf linear spline
#> 8.807000 3.889000 3.426000 5.526000
# This is the only honest way to choose. It costs
# ten lines and it turns a preference into
# evidence you can put in a report.The method is a decision, and a decision has to be written down. What you report matters as much as what you compute.
| Situation | Reach for | Why | Avoid |
|---|---|---|---|
| Under a few per cent of rows incomplete, no pattern found | Listwise deletion | Invents nothing and costs only precision. | Any imputation, which adds risk for no gain |
| Missingness depends on a variable you recorded | Group-wise mean, or regression on that variable | Conditioning on the observed cause removes the bias that deletion would leave. | Overall mean, which ignores the very structure you found |
| You will report a spread, a standard error or a confidence interval | Hot deck, or stochastic regression | Single-value methods shrink the spread by construction, so every interval afterwards is too narrow. | Mean and median imputation |
| You will report a correlation or fit a model | Complete cases, with the n stated | Every imputation distorts the relationship in a direction you have to argue about. | Regression imputation, which strengthens the relationship it assumed |
| Ordered series, gaps of one or two | Linear interpolation | The neighbours genuinely constrain the value, and the result stays inside the observed range. | Spline, whose extra flexibility buys nothing over two points |
| Ordered series, long runs | Leave the NA in place | Ten days of invented values look identical to ten days of measurement in every chart you will draw. | Any method that fills silently |
| Step-like series such as a price or a stock level | LOCF | The value genuinely held until it changed. | Linear interpolation, which invents a ramp |
| The error | Why it happens | What to do instead |
|---|---|---|
| Reading a sentinel as data | A file marks absence with −999, 0, or an empty cell. R reads it as a number and your mean drops. | Declare them at import: na.strings = c("", "NA", "-999"), then check
range() on every numeric column. |
| Deleting first and looking later | na.omit() is the fastest thing to type, so it gets typed before anyone has
counted what it removes. |
Count the rows it would drop and check whether they differ from the rest, before running it. |
| Imputing, then reporting a confidence interval | Single-value imputation shrinks the spread, and the software has no way to know the values were invented. | Report the interval from complete cases, or use a method that restores the variability you removed. |
| Filling a long run silently | One function call fills a one-day gap and a ten-day gap, and the output looks the same either way. | Check the run lengths with rle() first. Set a maximum gap you are willing to
fill and leave the rest as NA. |
NA means unknown. It propagates, and it cannot be compared with
==.attitude and on a dataset of your own.Solar.R rather than Ozone. Its
7 gaps behave differently from Ozone's 37.Wind instead of
Temp. Wind is noisier, and the ranking of the four methods may not survive.The aggregate() calls in this lesson were the summarise step used in passing.
Next week they become the topic, along with what happens when two tables have to be matched on
a key.
Solar.R, and one sentence naming which method you would
use and what it costs.Press T or Escape to close