| Feature | Example | Description |
|---|---|---|
| Asymmetry |  |
The distribution is not symmetrical. |
| Outliers |  |
Some observations are that are far from the rest. |
| Multimodality |  |
There are more than one "peak" in the observations. |
| Gaps |  |
Some continuous interval that are contained within the range but no observations exists. |
| Heaping |  |
Some values occur unexpectedly often. |
| Discretized |  |
Only certain values are found, e.g. due to rounding. |
| Implausible |  |
Values outside of plausible or likely range. |
ETC5521: Diving Deeply into Data Exploration
Working with a single variable, making transformations, detecting outliers, using robust statistics
Numerical features of a single quantitative variables
A measure of central tendency, e.g. mean, median and mode
A measure of dispersion (also called variability or spread), e.g. variance, standard deviation and interquartile range
There are other measures, e.g. skewness and kurtosis that measures “tailedness”, but these are not as common as the measures of first two
The mean is also the first moment and variance, skewness and kurtosis are second, third, and fourth central moments
Significance tests or hypothesis tests
- Testing for \(H_0: \mu = \mu_0\) vs. \(H_1: \mu \neq \mu_0\) (often \(\mu_0 = 0\))
- The \(t\)-test is commonly used if the underlying data are believed to be normally distributed
2019 Australian Federal Election (1/8)
Context
- There are 151 seats in the House of Representative for the 2019 Australian federal election
- The major parties in Australia are:
- the Coalition, comprising of the:
- Liberal,
- Liberal National (Qld),
- National, and
- Country Liberal (NT) parties, and
- the Australian Labor party
- the Coalition, comprising of the:
- The Greens party is a small but notable party
Source: PRObono
2019 Australian Federal Election (5/8)
What does this graph tell you?
- Majority of the country does not have first preference for the Greens
- Some constituents are slightly more supportive than the others
What further questions does it raise?
Notes:
- Australia uses full-preference instant-runoff voting in single member seats
- Following the full allocation of preferences, it is possible to derive a two-party-preferred figure, where the votes have been allocated between the two main candidates in the election.
- In Australia, this is usually between the candidates from the Coalition parties and the Australian Labor Party.
── Data Summary ────────────────────────
Values
Name tdf3
Number of rows 151
Number of columns 6
_______________________
Column type frequency:
character 3
numeric 3
________________________
Group variables None
── Variable type: character ────────────────────────────────
skim_variable n_missing complete_rate min max empty
1 DivisionID 0 1 3 3 0
2 DivisionNm 0 1 4 15 0
3 State 0 1 2 3 0
n_unique whitespace
1 151 0
2 151 0
3 8 0
── Variable type: numeric ──────────────────────────────────
skim_variable n_missing complete_rate mean sd
1 votes_GRN 0 1 9821. 5581.
2 votes_total 0 1 99925. 9801.
3 perc_GRN 0 1 9.87 5.63
p0 p25 p50 p75 p100 hist
1 2744 6555 8676 11532. 45876 ▇▂▁▁▁
2 51009 96372. 100936 105588 116216 ▁▁▁▇▅
3 2.89 6.43 8.55 11.4 47.8 ▇▂▁▁▁Formulating questions for EDA vs making observations from a plot
- BEFORE plotting or making summaries think broad (open-ended) questions about the distribution of values
- Questions with simple answers (i.e. yes or no) less helpful in encouraging exploration using graphics
- For example,
- What is the distribution of the first preference vote percentages for the Labor party across Australia?
- Is it evenly spread across electorates or are there clusters of popularity?
- AFTER plotting or making summaries think
was this what you expected, are there any surprises . Detail what you learn, and how you should follow up on these observations.
Is the outlying observation the electoral district that won the seat?
Lineup of Greens first preference percentages
Using the exponential distribution as the null says that we expect most electorates to have small tallies for Greens, and only a few electorates will have large, potentially winnable tallies.
NOTE: We’ve already seen the data so we can’t be impartial judges for choosing the most different plot. We can use the null plots to check whether the small mode of moderately high tallies is unusual if the tallies really are samples from an exponential distribution.
library(nullabor)
set.seed(241)
ggplot(lineup(null_dist("perc_GRN", "exp"), tdf3, n=10),
aes(x=perc_GRN)) +
geom_histogram(color = "white", fill = "#00843D", bins = 30) +
facet_wrap(~.sample, ncol=5, scales="free") +
theme(axis.text = element_blank(),
axis.title = element_blank(),
panel.grid.major = element_blank())Robust measure of central tendency
- Mean is a non-robust measure of location.
- Median is the 50% quantile of the observations
- Trimmed mean is the sample mean after discarding observations at the tails.
- Winsorized mean is the sample mean after replacing observations at the tails with the minimum or maximum of the observations that remain.
Both trimmed and Winsorized mean trimmed 20% of the tails.
| Plot | Mean | Median | Trimmed Mean | Winsorized Mean |
|---|---|---|---|---|
| 1 | 0.109 | 0.114 | 0.120 | 0.103 |
| 2 | 0.054 | -0.045 | -0.016 | -0.029 |
| 3 | 1.177 | 0.729 | 0.820 | 0.888 |
| 4 | 0.533 | 0.541 | 0.543 | 0.542 |
| 5 | 0.468 | 0.329 | 0.355 | 0.390 |
| 6 | 5.626 | 6.656 | 5.918 | 5.688 |
Robust measure of dispersion
- Standard deviation or its square, variance, is a popular choice of measure of dispersion but is not robust to outliers
- Standard deviation for sample \(x_1, ..., x_n\) is
\[\sqrt{\sum_{i=1}^n \frac{(x_i - \bar{x})^2}{n - 1}}\]
- Interquartile range difference between 1st and 3rd quartile, more robust measure of spread
- Median absolute deviance (MAD) is even more robust
\[\text{median}(|x_i - \text{median}(x_i)|)\]
Measure of dispersion
|
|||||
|---|---|---|---|---|---|
| Plot | SD | IQR | MAD | Skewness | Kurtosis |
| 1 | 0.90 | 1.19 | 0.87 | -0.072 | 3.0 |
| 2 | 0.99 | 1.41 | 1.08 | 0.358 | 2.2 |
| 3 | 1.33 | 1.18 | 0.79 | 1.944 | 7.2 |
| 4 | 0.29 | 0.45 | 0.34 | -0.126 | 1.8 |
| 5 | 0.47 | 0.50 | 0.34 | 1.691 | 6.4 |
| 6 | 2.78 | 5.36 | 2.98 | -0.351 | 1.7 |
Inference for robust statistics
We have seen the re-sampling methods simulation and permutation used for generating null plots in a lineup. Re-sampling methods can be used with numeric statistics also.
Simulation from distribution, can be used to to check for outliers.
We can also compute how many simulated values are more than the observed which gives a simulation \(p\)-value: 0.61.
For sample means, conventional tests provide a means for assessing what might be observed if different samples were taken.
Bootstrapping the current sample, can be used for robust statistics. If we have a sample of values:
[1] 2 2 3 6 7 7 8 8to bootstrap sample with replacement:
Here’s an example of bootstrapping to get a confidence interval for a median.
[1] "Median: 6.34"[1] "95% CI: ( 4.99 , 9.16 )"2019 Australian Federal Election (7/8)
Where are these electorates?
The width of the boxplot is proportional to the number of electoral districts in the corresponding state (which is roughly proportional to the population).
Outliers
Outliers are observations that are significantly different from the majority.
- Outliers can occur by chance in almost all distributions, but could be indicative of:
- a measurement error,
- a different population, or
- an issue with the sampling process.
Closer look at the boxplot
- Observations that are outside the range of lower to upper fence (1.5 times the box length) are often referred to as outliers.
- Plotting boxplots for data from a skewed distribution will almost always show these “outliers” but these are not necessarily outliers.
- Some definitions of outliers assume a symmetrical population distribution (e.g. in boxplots or observations a certain standard deviations away from the mean) and these definitions are ill-suited for asymmetrical distributions.
- Declaring observations outliers typically requires additional data context.
2019 Australian Federal Election (8/8)
Now what do you notice from this graph that you didn’t notice before?
Only two electoral districts in NT.
And only 3 and 5 electoral districts in ACT and TAS, respectively!
Boxplots requires 5 points!
We should have summarised the number of electoral districts for each state with numerical statistics as a first step.
Also the outlier (yes, safe to call this an outlier!) and the cluster in the Victoria electorates.
Melbourne Housing Prices (2/6)
Is missingness more likely for expensive houses?
- Check with a lineup
- To impute missings other variables will need to be used.
Note: Houses with more than 8 rooms removed. Why?
── Data Summary ────────────────────────
Values
Name df2
Number of rows 63023
Number of columns 13
_______________________
Column type frequency:
character 8
Date 1
numeric 4
________________________
Group variables None
── Variable type: character ────────────────────────────────
skim_variable n_missing complete_rate min max empty
1 Suburb 0 1 3 18 0
2 Address 0 1 7 27 0
3 Type 0 1 1 1 0
4 Method 0 1 1 2 0
5 SellerG 0 1 1 27 0
6 Postcode 0 1 4 4 0
7 Regionname 0 1 16 26 0
8 CouncilArea 0 1 17 30 0
n_unique whitespace
1 380 0
2 57754 0
3 3 0
4 9 0
5 476 0
6 225 0
7 8 0
8 34 0
── Variable type: Date ─────────────────────────────────────
skim_variable n_missing complete_rate min
1 Date 0 1 2016-01-28
max median n_unique
1 2018-10-13 2017-09-03 112
── Variable type: numeric ──────────────────────────────────
skim_variable n_missing complete_rate mean sd
1 Rooms 0 1 3.11 0.958
2 Price 14590 0.768 997898. 593499.
3 Propertycount 0 1 7618. 4424.
4 Distance 0 1 12.7 7.59
p0 p25 p50 p75 p100 hist
1 1 3 3 4 31 ▇▁▁▁▁
2 85000 620000 830000 1220000 11200000 ▇▁▁▁▁
3 39 4380 6795 10412 21650 ▅▇▅▂▁
4 0 7 11.4 16.7 64.1 ▇▆▁▁▁df2 |>
mutate(miss = ifelse(is.na(Price),
"Missing", "Recorded")) |>
count(Rooms, miss) |>
filter(Rooms < 8) |>
group_by(miss) |>
mutate(perc = n / sum(n) * 100) |>
ggplot(aes(as.factor(Rooms), perc, fill = miss)) +
geom_col(position = "dodge") +
scale_fill_viridis_d(begin=0.3, end=0.7) +
labs(x = "Rooms", y = "Percentage", fill = "Price") +
theme(aspect.ratio = 0.8)
Break the association between Rooms and Missing/Not on Price, because the null hypothesis is that there is no difference in missing status for price based on the size of the house. Why?
library(nullabor)
df2_d <- df2 |>
mutate(miss = ifelse(is.na(Price), "Missing", "Recorded")) |>
select(Rooms, miss) |>
filter(Rooms < 8)
df2_l <- lineup(null_permute("miss"), df2_d, n=10, pos=7)
df2_l_agg <- df2_l |>
group_by(.sample) |>
count(Rooms, miss) |>
ungroup() |>
group_by(miss) |>
mutate(perc = n / sum(n) * 100) |>
mutate(Rooms = as.factor(Rooms))
ggplot(df2_l_agg, aes(x=Rooms, y=perc, fill = miss)) +
geom_col(position = "dodge") +
scale_fill_viridis_d(begin=0.3, end=0.7) +
facet_wrap(~.sample, ncol=5) +
theme(legend.position = "none",
axis.text = element_blank(),
axis.title = element_blank(),
panel.grid.major.x = element_blank())Melbourne Housing Prices (3/6)
What can we say from this plot?
- The housing prices are right-skewed
- There appears to be a lot of outlying housing prices (how can we tell?)
Note: We determined that it is likely that more higher price houses have not disclosed the sale price. The distribution of price will need to be checked again after imputation.
Melbourne Housing Prices (4/6)
- The x-axis has been \(\log_{10}\)-transformed in this plot
- The plot appears more symmetrical now
- What is a useful measure of central tendency here?
General rules for transform quantitative variables
Non-shape changing, scaling:
- standardise to mean 0, sd 1
- standardise to min 0, max 1
- z-score
Shape changing, transformations: Remember the ladder of power transformations. (eg transforming left-skewed to more uniform using \(^2\))
Distribution changing: quantile
Some features cannot be fixed: gaps, multimodality, heaping. You need to find some explaining variable.
Some features can be artificially fixed: discreteness. If regularly discretized, add random uniform noise to spread equally between gaps.
Melbourne Housing Prices (6/6)
- Distribution from side-by-side univariate plots shows that higher number of rooms generally are pricier.
- This strata could be responsible for multimodality in price distribution, even though it is not visible in the histogram.
- Accounting for rooms is important.
Hidalgo stamps thickness
Famous historical example
- A stamp collector, Walton von Winkle, bought several collections of Mexican stamps from 1872-1874 and measured the thickness of all of them.
- The different bandwidth for the density plot suggests different possibilities for number of modes.
Which do you think most accurately reflects what’s in the data?
── Data Summary ────────────────────────
Values
Name Hidalgo1872
Number of rows 485
Number of columns 3
_______________________
Column type frequency:
numeric 3
________________________
Group variables None
── Variable type: numeric ──────────────────────────────────
skim_variable n_missing complete_rate mean sd p0
1 thickness 0 1 0.0860 0.0150 0.06
2 thicknessA 195 0.598 0.0922 0.0162 0.068
3 thicknessB 289 0.404 0.0768 0.00508 0.06
p25 p50 p75 p100 hist
1 0.075 0.08 0.098 0.131 ▅▇▃▂▁
2 0.0772 0.092 0.105 0.131 ▇▃▆▃▂
3 0.072 0.078 0.08 0.097 ▁▃▇▁▁hid <- ggplot(Hidalgo1872, aes(x=thickness, y=25)) +
geom_quasirandom(width=15, alpha=0.5, size=1) +
labs(x = "Thickness (0.001 mm)", y = "Density") +
ylim(c(0, 70)) +
theme(aspect.ratio = 0.6)
hid_p1 <- hid +
geom_density(aes(x=thickness), alpha=0.5,
color = "#E16A86", bw = 0.01, linewidth=2,
inherit.aes = FALSE)
hid_p2 <- hid +
geom_density(aes(x=thickness), alpha=0.5,
color = "#E16A86", bw = 0.0075, linewidth=2,
inherit.aes = FALSE)
hid_p3 <- hid +
geom_density(aes(x=thickness), alpha=0.5,
color = "#E16A86", bw = 0.004, linewidth=2,
inherit.aes = FALSE)
hid_p4 <- hid +
geom_density(aes(x=thickness), alpha=0.5,
color = "#E16A86", bw = 0.001, linewidth=2,
inherit.aes = FALSE)
hid_p1 + hid_p2 + hid_p3 + hid_p4 + plot_layout(ncol=2)Olive oil content
What do you see?
Mixture of discreteness and normal shape of continuous values. Why might this happen?
Check if there is a difference in the strata (here 1 thru 9), implying measurement policy differences.
Movie length
Upon further exploration, you can find the two movies that are well over 16 hours long are “Cure for Insomnia”, “Four Stars”, and “Longest Most Meaningless Movie in the World”
We can restrict our attention to films under 3 hours:
- Notice that there is a peak at particular times. Why do you think so?
── Data Summary ────────────────────────
Values
Name movies
Number of rows 58788
Number of columns 24
_______________________
Column type frequency:
character 2
numeric 22
________________________
Group variables None
── Variable type: character ────────────────────────────────
skim_variable n_missing complete_rate min max empty
1 title 0 1 1 121 0
2 mpaa 0 1 0 5 53864
n_unique whitespace
1 56007 0
2 5 0
── Variable type: numeric ──────────────────────────────────
skim_variable n_missing complete_rate mean
1 year 0 1 1976.
2 length 0 1 82.3
3 budget 53573 0.0887 13412513.
4 rating 0 1 5.93
5 votes 0 1 632.
6 r1 0 1 7.01
7 r2 0 1 4.02
8 r3 0 1 4.72
9 r4 0 1 6.37
10 r5 0 1 9.80
11 r6 0 1 13.0
12 r7 0 1 15.5
13 r8 0 1 13.9
14 r9 0 1 8.95
15 r10 0 1 16.9
16 Action 0 1 0.0797
17 Animation 0 1 0.0628
18 Comedy 0 1 0.294
19 Drama 0 1 0.371
20 Documentary 0 1 0.0591
21 Romance 0 1 0.0807
22 Short 0 1 0.161
sd p0 p25 p50 p75
1 23.7 1893 1958 1983 1997
2 44.3 1 74 90 100
3 23350085. 0 250000 3000000 15000000
4 1.55 1 5 6.1 7
5 3830. 5 11 30 112
6 10.9 0 0 4.5 4.5
7 5.96 0 0 4.5 4.5
8 6.45 0 0 4.5 4.5
9 7.59 0 0 4.5 4.5
10 9.73 0 4.5 4.5 14.5
11 11.0 0 4.5 14.5 14.5
12 11.6 0 4.5 14.5 24.5
13 11.3 0 4.5 14.5 24.5
14 9.44 0 4.5 4.5 14.5
15 15.7 0 4.5 14.5 24.5
16 0.271 0 0 0 0
17 0.243 0 0 0 0
18 0.455 0 0 0 1
19 0.483 0 0 0 1
20 0.236 0 0 0 0
21 0.272 0 0 0 0
22 0.367 0 0 0 0
p100 hist
1 2005 ▁▁▃▃▇
2 5220 ▇▁▁▁▁
3 200000000 ▇▁▁▁▁
4 10 ▁▃▇▆▁
5 157608 ▇▁▁▁▁
6 100 ▇▁▁▁▁
7 84.5 ▇▁▁▁▁
8 84.5 ▇▁▁▁▁
9 100 ▇▁▁▁▁
10 100 ▇▁▁▁▁
11 84.5 ▇▂▁▁▁
12 100 ▇▃▁▁▁
13 100 ▇▃▁▁▁
14 100 ▇▁▁▁▁
15 100 ▇▃▁▁▁
16 1 ▇▁▁▁▁
17 1 ▇▁▁▁▁
18 1 ▇▁▁▁▃
19 1 ▇▁▁▁▅
20 1 ▇▁▁▁▁
21 1 ▇▁▁▁▁
22 1 ▇▁▁▁▂ggplot(movies, aes(length)) +
geom_histogram(color = "white") +
labs(x = "Length of movie (minutes)", y = "Frequency") +
theme(aspect.ratio = 0.6)
ggplot(movies, aes(length)) +
geom_histogram(color = "white") +
labs(x = "Length of movie (minutes)", y = "Frequency") +
scale_x_log10() +
theme(aspect.ratio = 0.6)2019 Australian Federal Election
df1 <- read_csv(here::here("data/HouseFirstPrefsByCandidateByVoteTypeDownload-24310.csv"),
skip = 1,
col_types = cols(
.default = col_character(),
OrdinaryVotes = col_double(),
AbsentVotes = col_double(),
ProvisionalVotes = col_double(),
PrePollVotes = col_double(),
PostalVotes = col_double(),
TotalVotes = col_double(),
Swing = col_double()
)
)Lineup of tuberculosis count between sexes
For this problem there is nothing more to learn from a lineup that what can be learned from a conventional hypothesis test of \(H_0: p=0.5\).
Exact binomial test
data: tb_oz_2012$count
number of successes = 314, number of trials = 997,
p-value <2e-16
alternative hypothesis: true probability of success is not equal to 0.5
95 percent confidence interval:
0.29 0.34
sample estimates:
probability of success
0.31 library(nullabor)
set.seed(252)
tb_oz_2012 <- tb |>
filter(iso3 == "AUS",
year == 2012) |>
select(iso3, year, new_sp_m04:new_ep_m65) |>
pivot_longer(cols=new_sp_m04:new_ep_m65,
names_to = "var",
values_to = "count") |>
separate(var, into=c("new", "type", "sexage")) |>
select(-new) |>
filter(!(sexage %in% c("sexunk014", "sexunk04",
"sexunk15plus", "sexunk514",
"f15plus", "m15plus"))) |>
group_by(sexage) |>
summarise(count = sum(count, na.rm=TRUE)) |>
mutate(sex=str_sub(sexage, 1, 1),
age=str_sub(sexage, 2, str_length(sexage))) |>
group_by(sex) |>
summarise(count=sum(count)) |>
ungroup() |>
mutate(sex01 = ifelse(sex=="m", 0, 1)) |>
select(-sex)
ggplot(lineup(null_dist("count", "binom",
list(size=sum(tb_oz_2012$count),
p=0.5)),
tb_oz_2012, n=10),
aes(x=sex01, y=count)) +
geom_col() +
facet_wrap(~.sample, ncol=5) +
theme(axis.text = element_blank(),
axis.title = element_blank(),
panel.grid.major = element_blank())Imputing missings for univariate distributions
Quantitative variable: Simulate from a fitted distribution.
df2 <- df2 |>
mutate(lPrice = log10(Price),
price_miss = ifelse(is.na(Price), "yes", "no"))
df2_smry <- df2 |>
summarise(m = mean(lPrice, na.rm=TRUE),
s = sd(lPrice, na.rm=TRUE))
set.seed(1003)
df2 <- df2 |>
rowwise() |>
mutate(lPrice = ifelse(price_miss == "yes",
rnorm(1, df2_smry$m, df2_smry$s), lPrice)) |>
mutate(Price = ifelse(price_miss == "yes", 10^lPrice, Price))
Categorical variable: Simulate from multinomial.
# A tibble: 7 × 2
age count
<chr> <dbl>
1 1524 27
2 2534 48
3 3544 15
4 4554 11
5 5564 9
6 65 15
7 u 12 [1] 5 3 1 1 2 1 1 2 1 1 1 1# A tibble: 7 × 2
age count
<chr> <dbl>
1 1524 35
2 2534 50
3 3544 16
4 4554 11
5 5564 10
6 65 15
7 u 12imputeMulti library can automate for multiple variables.