Exploratory Data Analysis

Look before you conclude

Margriet de Wit asks a reasonable question: "Do customers on a dynamic contract use more energy than the others?" The query that answers it takes a minute. Whether the answer can be trusted depends on what is in the data, and nobody at Kestrel has looked at the new analytics copy yet.

Exploratory data analysis (EDA) is that look: examining a dataset with summary numbers and simple charts to learn its structure, its quality and its surprises, before testing any specific idea. The approach was championed by the statistician John Tukey, whose 1977 book gave it its name, and whose box plot and five-number summary are still its basic tools. EDA is open-ended on purpose. It is where you find the errors that would otherwise end up in a board report, and the questions worth asking next.

A useful order of questions:

1. What is one row? One meter on one day, one invoice, one customer. Is that true for every row?
2. Is it complete? All periods, all meters, no duplicates, no rows without a parent?
3. Are the values plausible? Ranges, units, negative numbers, impossible dates, spelling variants.
4. What does each variable look like? Center, spread, shape, extremes.
5. How do variables relate? And what else could explain that relationship?

The examples use the Kestrel sample data from Analytical SQL: daily meter readings for September 2026, with contracts and customers.

Profiling: the first ten minutes

Data profiling is the systematic first pass over a table: for every column, count the rows, the distinct values and the missing values, and find the minimum and maximum. It is quick, it is the same for every dataset, and it finds most data-quality problems before any analysis starts.

SELECT count(*)                     AS rows_total,
       count(DISTINCT meter_id)     AS meters,
       count(DISTINCT reading_date) AS days,
       min(reading_date)            AS first_day,
       max(reading_date)            AS last_day,
       min(kwh_import)              AS min_kwh,
       max(kwh_import)              AS max_kwh,
       count(*) - count(kwh_export) AS export_nulls
FROM daily_readings;

┌────────────┬────────┬───────┬────────────┬────────────┬──────────────┬──────────────┬──────────────┐
│ rows_total │ meters │ days  │ first_day  │  last_day  │   min_kwh    │   max_kwh    │ export_nulls │
│   int64    │ int64  │ int64 │    date    │    date    │ decimal(8,3) │ decimal(8,3) │    int64     │
├────────────┼────────┼───────┼────────────┼────────────┼──────────────┼──────────────┼──────────────┤
│        360 │     12 │    30 │ 2026-09-01 │ 2026-09-30 │       -4.210 │       91.300 │          210 │
└────────────┴────────┴───────┴────────────┴────────────┴──────────────┴──────────────┴──────────────┘

Twelve meters times thirty days is 360, and there are 360 rows. It is tempting to stop here. But the minimum is negative, which a meter that only measures consumption cannot produce, and the maximum of 91.300 kWh is about ten times what most households use in a day. The 210 missing export values are expected: seven meters have no solar panels, and the export is NULL ("not measured") rather than 0 for them. Profiling the grain, one row per meter per day, shows that the total of 360 is a coincidence:

SELECT meter_id, count(*) AS readings
FROM daily_readings
GROUP BY meter_id
HAVING count(*) <> 30
ORDER BY meter_id;

┌──────────┬──────────┐
│ meter_id │ readings │
│ varchar  │  int64   │
├──────────┼──────────┤
│ M-30005  │       31 │
│ M-30010  │       29 │
└──────────┴──────────┘

M-30005 has a duplicated day and M-30010 a missing one; together they cancel out in the total. Every profiling check that matters compares against an expectation ("30 per meter"), not against a single total.

Summary statistics and the five-number summary

The mean and median and the standard deviation are explained in Math for IT. For exploration, Tukey's five-number summary is the more robust description: the minimum, the first quartile (Q1, the 25th percentile), the median, the third quartile (Q3, the 75th percentile) and the maximum. It shows the center, the spread of the middle half and the extremes in five numbers, and a single bad value can only move the minimum or the maximum. For meter M-30001:

SELECT min(kwh_import)                  AS min,
       quantile_cont(kwh_import, 0.25)  AS q1,
       median(kwh_import)               AS median,
       quantile_cont(kwh_import, 0.75)  AS q3,
       max(kwh_import)                  AS max
FROM daily_readings
WHERE meter_id = 'M-30001';

┌──────────────┬──────────────┬──────────────┬──────────────┬──────────────┐
│     min      │      q1      │    median    │      q3      │     max      │
│ decimal(8,3) │ decimal(8,3) │ decimal(8,3) │ decimal(8,3) │ decimal(8,3) │
├──────────────┼──────────────┼──────────────┼──────────────┼──────────────┤
│        6.650 │        7.520 │        8.626 │        9.232 │       91.300 │
└──────────────┴──────────────┴──────────────┴──────────────┴──────────────┘

Half of this household's days fall between 7.520 and 9.232 kWh. The distance between the quartiles is the interquartile range (IQR), here 1.712 kWh: the spread of the middle 50% of the values, unaffected by the extremes. A maximum of 91.300 next to a Q3 of 9.232 is a signal by itself.

Compare what the same bad reading does to the mean and the standard deviation. The Python statistics module computes both, for the month as recorded and with the reading corrected to 9.130:

import statistics as st

# Meter M-30001, 1-30 September 2026 (kWh per day)
kwh = [6.650, 8.490, 7.730, 6.970, 10.117, 9.240, 7.290, 91.300, 8.370, 7.610,
       6.850, 9.958, 9.081, 7.170, 9.010, 8.250, 7.490, 9.330, 9.799, 8.922,
       7.050, 8.890, 8.130, 7.370, 9.210, 9.640, 8.763, 9.530, 8.770, 8.010]
fixed = [9.130 if value == 91.300 else value for value in kwh]

for label, data in (("as recorded", kwh), ("corrected", fixed)):
    print(f"{label:<12} mean {st.mean(data):6.3f}  median {st.median(data):6.3f}"
          f"  stdev {st.stdev(data):6.3f}")

print("quartiles, exclusive:", [round(q, 3) for q in st.quantiles(kwh, n=4)])
print("quartiles, inclusive:", [round(q, 3) for q in st.quantiles(kwh, n=4, method="inclusive")])

as recorded  mean 11.166  median  8.627  stdev 15.168
corrected    mean  8.427  median  8.627  stdev  1.014
quartiles, exclusive: [7.46, 8.627, 9.262]
quartiles, inclusive: [7.52, 8.627, 9.232]

One value in thirty raises the mean by a third and the standard deviation fifteenfold; the median does not move. The last two lines show a smaller trap: quartiles are computed with different methods. Python's default ("exclusive") gives a Q1 of 7.460, the "inclusive" method gives 7.520, the same as DuckDB's quantile_cont and pandas' default. Neither is wrong (see methods that interpolate), but two reports that use different tools can disagree in the second decimal. State the method when the exact value matters. (The medians differ in the third decimal for a similar reason: the exact median is 8.6265; DuckDB keeps the column's DECIMAL(8,3) type and shows 8.626, Python rounds the float to 8.627.)

Looking at distributions

A summary compresses a column into a few numbers; a histogram shows its whole shape. Counting the household readings in bins of 2 kWh takes one query:

SELECT floor(r.kwh_import / 2) * 2 AS bin_start, count(*) AS days
FROM daily_readings r
JOIN meters m     USING (meter_id)
JOIN contracts c  USING (contract_id)
JOIN customers cu USING (customer_id)
WHERE cu.segment = 'household'
GROUP BY bin_start
ORDER BY bin_start;

┌───────────┬───────┐
│ bin_start │ days  │
│  double   │ int64 │
├───────────┼───────┤
│      -6.0 │     1 │
│       4.0 │    26 │
│       6.0 │    68 │
│       8.0 │    92 │
│      10.0 │    42 │
│      12.0 │    10 │
│      90.0 │     1 │
└───────────┴───────┘

The bulk of the readings form one hump around 8 to 10 kWh, slightly stretched to the right, which is what household consumption usually looks like. Two values sit far outside it. When you read a distribution, describe its center, its spread, its shape (one peak or two, symmetric or skewed) and anything that does not belong. With businesses included the picture would have a second cluster far to the right (the bakery and the printing company use 40 to 72 kWh on a weekday), which is a reason to look at segments separately rather than all meters at once. Which chart suits which question, including box plots that show the five-number summary, is on Choosing the Right Chart.

Finding outliers: flag, then investigate

Tukey's rule turns the IQR into a test: a value is a possible outlier when it lies more than 1.5 × IQR below Q1 or above Q3. These limits are called fences. Because each meter has its own normal level, compute the fences per meter:

WITH fences AS (
  SELECT meter_id,
         quantile_cont(kwh_import, 0.25) AS q1,
         quantile_cont(kwh_import, 0.75) AS q3
  FROM daily_readings
  GROUP BY meter_id
)
SELECT r.meter_id, r.reading_date, r.kwh_import,
       round(f.q1 - 1.5 * (f.q3 - f.q1), 3) AS low_fence,
       round(f.q3 + 1.5 * (f.q3 - f.q1), 3) AS high_fence
FROM daily_readings r
JOIN fences f USING (meter_id)
WHERE r.kwh_import < f.q1 - 1.5 * (f.q3 - f.q1)
   OR r.kwh_import > f.q3 + 1.5 * (f.q3 - f.q1)
ORDER BY r.meter_id, r.reading_date;

┌──────────┬──────────────┬──────────────┬───────────────┬───────────────┐
│ meter_id │ reading_date │  kwh_import  │   low_fence   │  high_fence   │
│ varchar  │     date     │ decimal(8,3) │ decimal(12,3) │ decimal(12,3) │
├──────────┼──────────────┼──────────────┼───────────────┼───────────────┤
│ M-30001  │ 2026-09-08   │       91.300 │         4.952 │        11.800 │
│ M-30009  │ 2026-09-21   │       -4.210 │         3.860 │         7.860 │
│ M-30010  │ 2026-09-27   │       12.310 │         7.495 │        11.603 │
└──────────┴──────────────┴──────────────┴───────────────┴───────────────┘

Three flags, three different stories:

• M-30001, 91.300. A recording error: the decimal point is one place off (probably 9.130). Confirm it against the meter's register reading, correct it at the source and note the correction.
• M-30009, −4.210. A consumption meter cannot run backward; this is an error too, but the right value is unknown. Exclude it from averages, report it to the meter data team, and say in the report that you did.
• M-30010, 12.310 on Sunday 27 September. Probably real: households use more at weekends, and this one is a little above its fence. Removing it would make the data look tidier and less true.

Relationships: plot before you calculate

For two numeric variables the usual summary is the correlation coefficient r and a fitted straight line. In 1973 the statistician Francis Anscombe published four small datasets to show how little those numbers say on their own. Known as Anscombe's quartet, the four sets have practically the same means, variances, correlation and regression line, yet look nothing alike:

import statistics as st

x = [10, 8, 13, 9, 11, 14, 6, 4, 12, 7, 5]
quartet = {
    "I":   (x, [8.04, 6.95, 7.58, 8.81, 8.33, 9.96, 7.24, 4.26, 10.84, 4.82, 5.68]),
    "II":  (x, [9.14, 8.14, 8.74, 8.77, 9.26, 8.10, 6.13, 3.10, 9.13, 7.26, 4.74]),
    "III": (x, [7.46, 6.77, 12.74, 7.11, 7.81, 8.84, 6.08, 5.39, 8.15, 6.42, 5.73]),
    "IV":  ([8, 8, 8, 8, 8, 8, 8, 19, 8, 8, 8],
            [6.58, 5.76, 7.71, 8.84, 8.47, 7.04, 5.25, 12.50, 5.56, 7.91, 6.89]),
}

print("set  mean_x  var_x  mean_y  var_y     r  line")
for name, (xs, ys) in quartet.items():
    fit = st.linear_regression(xs, ys)
    print(f"{name:<4} {st.mean(xs):6.2f} {st.variance(xs):6.2f} {st.mean(ys):7.2f}"
          f" {st.variance(ys):6.2f} {st.correlation(xs, ys):5.2f}"
          f"  y = {fit.intercept:.2f} + {fit.slope:.3f}x")

set  mean_x  var_x  mean_y  var_y     r  line
I      9.00  11.00    7.50   4.13  0.82  y = 3.00 + 0.500x
II     9.00  11.00    7.50   4.13  0.82  y = 3.00 + 0.500x
III    9.00  11.00    7.50   4.12  0.82  y = 3.00 + 0.500x
IV     9.00  11.00    7.50   4.12  0.82  y = 3.00 + 0.500x

Set I is what the numbers suggest. Set II is a curve, so a straight line is the wrong model. In set III one outlier tilts the line away from an otherwise perfect relationship. In set IV there is no relationship at all except one point. The correlations agree to two decimals (0.816 for the first three and 0.817 for the fourth, to three). A modern version, the Datasaurus Dozen (Matejka and Fitzmaurice, 2017), has thirteen datasets, one shaped like a dinosaur, that share their means, standard deviations and correlation to two decimals (truncated, not rounded). The lesson of both: always plot the data.

Correlation, confounding and coincidence

Back to Margriet's question. Average consumption per meter per day, by tariff:

SELECT c.tariff,
       count(DISTINCT r.meter_id)   AS meters,
       round(avg(r.kwh_import), 2)  AS avg_kwh_per_day
FROM daily_readings r
JOIN meters m    USING (meter_id)
JOIN contracts c USING (contract_id)
GROUP BY c.tariff
ORDER BY c.tariff;

┌───────────┬────────┬─────────────────┐
│  tariff   │ meters │ avg_kwh_per_day │
│  varchar  │ int64  │     double      │
├───────────┼────────┼─────────────────┤
│ day_night │      3 │           20.32 │
│ dynamic   │      3 │           22.86 │
│ fixed     │      3 │           12.07 │
│ variable  │      3 │            6.57 │
└───────────┴────────┴─────────────────┘

Dynamic customers seem to use the most. But tariff is not the only difference between these groups: each of the dynamic, day_night and fixed tariffs includes one business, and businesses use several times as much as households. Split by segment as well:

SELECT cu.segment, c.tariff,
       count(DISTINCT r.meter_id)   AS meters,
       round(avg(r.kwh_import), 2)  AS avg_kwh_per_day
FROM daily_readings r
JOIN meters m     USING (meter_id)
JOIN contracts c  USING (contract_id)
JOIN customers cu USING (customer_id)
GROUP BY cu.segment, c.tariff
ORDER BY cu.segment, c.tariff;

┌───────────┬───────────┬────────┬─────────────────┐
│  segment  │  tariff   │ meters │ avg_kwh_per_day │
│  varchar  │  varchar  │ int64  │     double      │
├───────────┼───────────┼────────┼─────────────────┤
│ business  │ day_night │      1 │           39.85 │
│ business  │ dynamic   │      1 │           51.83 │
│ business  │ fixed     │      1 │           15.44 │
│ household │ day_night │      2 │           10.55 │
│ household │ dynamic   │      2 │            8.37 │
│ household │ fixed     │      2 │           10.36 │
│ household │ variable  │      2 │            6.17 │
└───────────┴───────────┴────────┴─────────────────┘

Among households, dynamic customers use less than fixed and day_night customers (8.37 against 10.36 and 10.55 kWh per day). The high overall figure came from one bakery. Segment here is a confounding variable: a third factor that influences both things being compared (which tariff a customer chooses, and how much energy it uses), and so creates or hides an apparent relationship between them. The usual remedy is to compare within groups that share the confounder: households with households.

Notice also that the variable tariff lost a meter between the two queries: 3 meters without the customer join, 2 with it. The third belongs to contract CT-20014, whose customer does not exist, and the inner join dropped it without a warning.

Two honest caveats belong next to this table in any report. Each household group has two meters, far too few to conclude anything (see Statistics for Analysts). And the household fixed figure still contains the 91.300 error, which lifts it by almost 1.4 kWh per day: dirty data is a confounder too.

A spurious correlation is a correlation between two variables that are not causally connected: they move together by coincidence, or because both follow a common driver such as time or size. If Kestrel's number of installed smart meters and the number of calls to the service desk both rise over a year simply because the customer base grows, they correlate strongly; calls per customer are what Daniel Mensah should compare. The same happens with any two series that trend upward over the same months. Before believing a correlation, ask what else changed at the same time, compare rates instead of totals, and, where the decision matters, test the cause directly with an experiment (see A/B tests).

Simpson's paradox

Confounding at its most striking is Simpson's paradox: a trend that holds in most or all subgroups reverses, or disappears, when the groups are combined. It happens when the subgroups differ in size and in their baseline, so that the combined figure mostly reflects the mix of the groups, not the comparison you wanted.

The classic case is graduate admissions at the University of California, Berkeley, in 1973, analyzed by Bickel, Hammel and O'Connell in Science (1975). Over the whole university 44% of male and 35% of female applicants were admitted, which looked like discrimination against women. The six largest departments tell a different story:

# Berkeley 1973, six largest departments: (applied, admitted)
men = {"A": (825, 512), "B": (560, 353), "C": (325, 120),
       "D": (417, 138), "E": (191, 53), "F": (373, 22)}
women = {"A": (108, 89), "B": (25, 17), "C": (593, 202),
         "D": (375, 131), "E": (393, 94), "F": (341, 24)}

def rate(applied, admitted):
    return 100 * admitted / applied

def total(group):
    return [sum(values) for values in zip(*group.values())]

print("dept  men %  women %  women applied")
for dept in men:
    print(f"{dept:4} {rate(*men[dept]):6.1f} {rate(*women[dept]):8.1f} {women[dept][0]:14}")
print(f"all  {rate(*total(men)):6.1f} {rate(*total(women)):8.1f} {total(women)[0]:14}")

dept  men %  women %  women applied
A      62.1     82.4            108
B      63.0     68.0             25
C      36.9     34.1            593
D      33.1     34.9            375
E      27.7     23.9            393
F       5.9      7.0            341
all    44.5     30.4           1835

In four of the six departments women were admitted at a higher rate than men, yet over the six departments together 30.4% of women were admitted against 44.5% of men. Women applied mostly to departments C to F, which admitted few applicants of either sex; men applied mostly to A and B, which admitted most. The combined rate measures where people applied more than how they were treated. (The six-department totals are not the same as the whole-university figures above; do not mix the two.)

The effect was described by Karl Pearson and George Udny Yule around 1900 and by Edward Simpson in 1951; the name "Simpson's paradox" came later. In daily analytics it shows up whenever a mix changes: average handling time rises because the service desk got more complex tickets, not because anyone got slower; average consumption per customer drops because many small households joined. Whenever a total and its parts seem to disagree, look at the mix.

From exploration to findings

Exploration produces a list, and the list is a deliverable: it tells the next person what the data can and cannot answer. For the Kestrel readings, the first hour produced this:

Each row is also a candidate for an automated check, so the problem is caught the next time without anyone looking: see Data Pipelines on data tests, and Data Types and Data Quality on cleaning. The questions that remain, such as whether the difference between tariffs is real or chance, are for Statistics for Analysts.