1.1 Research written in code

Python is a free, open-source programming language for working with data: cleaning it, analysing it, and turning it into tables, graphs, and statistical models. It was created in the late 1980s by Guido van Rossum in the Netherlands, released in 1991, and named after the British comedy troupe Monty Python. Today, Python is the most popular programming language in the world, powering scientific research in universities, national laboratories, hospitals, tech companies, and governments across every continent.

Four properties make Python worth learning for a researcher:

1. It is built for clear reading: Python was designed from the ground up to look like clear, readable English. Its syntax is clean and uncluttered, making it accessible to beginners who have never programmed before.
2. It is free and open source: You, your students, and anyone reviewing your thesis can install and run your code without buying expensive proprietary licences.
3. It makes research reproducible: Your written script records every cleaning choice, transformation, and statistical test. Anyone can re-run your script on the raw data and reproduce your exact numbers, tables, and figures. Chapter 17 explores this idea in depth.
4. An unmatched scientific ecosystem: Thousands of scientists have contributed free, battle-tested libraries: pandas for spreadsheets and data tables, matplotlib and seaborn for publication-quality charts, scipy and statsmodels for hypothesis testing and regression, and scikit-learn for machine learning. Whatever statistical method your thesis needs, a Python library already exists for it.

The price of code is that Python asks you to type rather than click. At first, typing feels unfamiliar. Within a few days, however, it becomes second nature—and it pays you back every time you need to re-run an analysis, update a dataset, or generate new figures for your thesis defense.

1.2 Installing Python and VS Code

Working with Python involves two tools:

1. Python is the engine that executes your commands.
2. An editor (or Integrated Development Environment, IDE) is the program in which you drive it: you write Python code there, run it, and inspect your tables and charts.

You will almost never open Python by itself; you open your editor, and the editor talks to Python for you.

The most popular, modern editor for scientific Python is Visual Studio Code (VS Code), developed by Microsoft. It is free, fast, and available for Windows, macOS, and Linux.

1. Install Python: Download Python (version 3.10 or newer) from python.org. On Windows, ensure you check the box that says “Add Python to PATH” during the installer setup.
2. Install VS Code: Download and install VS Code from code.visualstudio.com.
3. Install the Python Extension: Open VS Code, click the Extensions icon on the left sidebar (four squares), search for Python (by Microsoft), and click Install.

Appendix A walks through this setup step by step for Windows, macOS, and Linux, with troubleshooting for common path issues.

1.3 The editor layout

When you open VS Code with a Python file, the window is organized into functional areas, described in Table 1.1.

You can type code directly into the Python Interactive Terminal (or REPL). When you type 2 + 2 and press Enter, Python evaluates the expression and prints 4.

1.4 Calculations

The simplest use of Python is as a calculator. Suppose you slept 6.5, 7, and 5.5 hours on three consecutive nights last week. Your average sleep is:

(6.5 + 7 + 5.5) / 3

6.333333333333333

Python follows standard mathematical order of operations (PEMDAS): parentheses first, followed by exponents, multiplication and division, and finally addition and subtraction. Without the parentheses, Python would divide only the last number by 3:

6.5 + 7 + 5.5 / 3

15.333333333333334

The standard arithmetic operators in Python are: * + for addition * - for subtraction * * for multiplication * / for division * ** for exponents/powers (for example, 2 ** 3 evaluates to 8; note that in Python, ^ is a bitwise operator, not power!)

2 ** 3

8

1.5 Storing values in variables

Typing numbers repeatedly is slow and error-prone. Instead, you can store a value under a memorable name. In Python, this stored value is held in a variable, created using the assignment operator =:

nights = 3

When you run this line, Python produces no visible output, but it now remembers that nights holds the value 3. You can now use nights anywhere you would use the number:

nights

3

nights * 7

21

To store several values together, you can place them in a list, written inside square brackets [ and ], separated by commas:

sleep = [6.5, 7.0, 5.5]
sleep

[6.5, 7.0, 5.5]

Notice how much clearer and self-documenting your analysis becomes when the numbers have meaningful names:

sum(sleep) / nights

6.333333333333333

Variable names in Python follow clear conventions: * Names start with a letter or underscore, and may contain letters, numbers, and underscores (e.g., sleep_hours, student_01). * Names cannot contain spaces. Use an underscore _ to separate words (this is called snake_case, the standard Python style). * Python is strictly case-sensitive: Sleep, sleep, and SLEEP are three completely different variables. * Choose descriptive names: sleep_hours tells you immediately what the data represents; single letters like x or temp will leave you guessing when you review your code months later.

If you assign a new value to an existing variable name, Python overwrites the old value:

nights = 4
nights

4

1.6 Types of data

Research data is not only numbers. A study records measurements (hours of sleep), text categories (a student’s faculty), and yes-or-no status flags (whether a student attended a workshop). Python keeps track of the data type of every value, because the type determines what mathematical and logical operations are permitted.

Table 1.2 lists the core types you will encounter in research.

The built-in function type() inspects the type of any object:

hours = 6.5
faculty = "Education"
invited = True

print(type(hours))
print(type(faculty))
print(type(invited))

<class 'float'>
<class 'str'>
<class 'bool'>

Text must always be enclosed in quotes (single '...' or double "..."). If you omit the quotes, Python assumes you are referring to a variable name:

faculty = Education

---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
Cell In[12], line 1
----> 1 faculty = Education

NameError: name 'Education' is not defined

This error—NameError: name 'Education' is not defined—is one of the most frequent beginner mistakes. The error message is helpful and specific: Python searched your computer’s memory for a variable named Education and found none.

import numpy as np

sleep_with_gap = [6.5, np.nan, 5.5]
sleep_with_gap

[6.5, nan, 5.5]

1.7 Functions and methods

You have already used several built-in functions: sum(), len(), type(), and print(). A function takes one or more inputs, performs a procedure, and returns a result. You call a function by writing its name followed by parentheses (), enclosing the inputs inside:

print(sum(sleep))
print(max(sleep))
print(len(sleep))

19.0
7.0
3

Inputs passed to a function are called arguments. If a function accepts multiple arguments, they are separated by commas. For example, round() takes the number you want to round, and the number of decimal digits:

round(6.333333, 1)

6.3

You can also pass arguments by their name (called keyword arguments):

round(6.333333, ndigits=2)

6.33

In Python, objects also have methods. A method is a function that belongs directly to an object, called with a dot .. For example, strings have methods to change their casing:

faculty = "Education"
print(faculty.lower())
print(faculty.upper())

education
EDUCATION

help(round)

1.8 Writing scripts

Commands typed in an interactive prompt are temporary. For real research work that needs to be defended, documented, and published, you write your code in a script: a plain text file ending in .py.

1. In your editor, create a new file named analysis.py.
2. Type your code, with one instruction per line.
3. Add comments using the hash symbol #. Python completely ignores everything after # on that line; comments are written for human readers (your future self, your supervisor, and your peer reviewers) to explain why a particular decision was made:

# Sleep recorded across three nights (hours per night)
sleep_hours = [6.5, 7.0, 5.5]

# Calculate the mean sleep duration
average_sleep = sum(sleep_hours) / len(sleep_hours)

# Display result rounded to one decimal place for thesis table
print(f"Average sleep: {round(average_sleep, 1)} hours")

Average sleep: 6.3 hours

4. To run the script, press Ctrl+Enter (or Cmd+Enter on macOS) to run the selected lines in the interactive window, or run the entire script in the terminal with python analysis.py.

1.9 Managing research projects

A thesis involves dozens of files: survey exports, cleaned datasets, analysis scripts, generated plots, and chapter drafts.

The gold standard for scientific computing is to create a dedicated project directory on your computer (for example, graduate-wellbeing-thesis/). Keep all code and data self-contained within this folder:

graduate-wellbeing-thesis/
├── data/              <- Raw survey data (read-only)
├── scripts/           <- Your Python analysis scripts (.py)
├── figures/           <- Exported charts and graphs
└── thesis_chapter.qmd <- Your written manuscript

By organizing your work with relative paths (e.g. data/students.csv), your analysis is fully portable: if you zip the folder and email it to your supervisor, or push it to GitHub, the code runs on their computer without any file path errors.

1.10 Installing and importing packages

Python includes a rich standard library, but its true power comes from third-party packages (also called libraries). Using a package is a two-step process:

1. Install the package once onto your computer using Python’s package manager, pip, in your system terminal:

   BASH Copy

   pip install pandas seaborn data2thesis

2. Import the package in each script where you wish to use it:

   PYTHON Copy

   import pandas as pd
   import data2thesis as dt

A helpful analogy is a workshop: running pip install purchases a specialized power tool and places it in your workshop cabinet. You only need to buy it once. Writing import takes the tool off the shelf and places it on your workbench for today’s session.

1.11 A first look at Elaf’s data

With Python and the companion package installed, we can load the Graduate Wellbeing Study:

import data2thesis as dt

students = dt.load_students()

The students object is a pandas DataFrame: a two-dimensional table of rows and columns, similar to an Excel sheet or SPSS dataset, but optimized for fast computation.

To see the dimensions of the table, check its .shape attribute:

students.shape

(600, 14)

The output (600, 14) tells us immediately that Elaf’s study contains 600 students (rows) and 14 variables (columns).

To inspect the first five rows of the dataset, call the .head() method:

students.head()

To see the names of all 14 variables:

list(students.columns)

['student_id',
 'supervisor_id',
 'age',
 'gender',
 'faculty',
 'programme',
 'study_mode',
 'employment',
 'has_children',
 'lives_away',
 'financial_worry',
 'workshop',
 'workshop_sessions',
 'considering_dropout']

Each column is described in the study codebook (dt.load_students?).

To select a single variable from the DataFrame, write the DataFrame name followed by the column name in brackets and quotes: students['age'].

A natural first question in Elaf’s thesis: What is the average age of graduate students entering the study?

In pandas, we call the .mean() method directly on the column:

students['age'].mean()

np.float64(29.724540901502504)

Notice something important: in pandas, .mean() is designed from the ground up for real-world empirical research: it calculates the average across all known observations while safely ignoring missing values by default (skipna=True).

The average age of the students in the study is 29.7 years.

We can also check whether any student had a missing age using .isna().sum():

students['age'].isna().sum()

np.int64(1)

Only 1 student’s age was missing from the initial intake survey. Chapter 3 shows how to identify why it was missing and how to clean it.

Next, we can generate frequency tables for categorical variables using the .value_counts() method. For example, how many students participated from each academic faculty?

students['faculty'].value_counts()

faculty
Health Sciences     154
Education           148
Social Sciences     116
Humanities           95
Natural Sciences     87
Name: count, dtype: int64

And how many were invited to the wellbeing workshop?

students['workshop'].value_counts()

workshop
Not invited    300
Invited        300
Name: count, dtype: int64

Exactly 300 students—precisely 50% of the sample—were invited to the workshop. That is not an accident: the workshop invitations were assigned using random allocation. Chapter 5 explores why random assignment is the cornerstone of causal inference in empirical research.

import seaborn as sns

penguins = sns.load_dataset("penguins")
print(penguins.shape)
penguins.head(3)

(344, 7)

print(penguins['species'].value_counts())
print(f"Average body mass: {penguins['body_mass_g'].mean():.1f} grams")

species
Adelie       152
Gentoo       124
Chinstrap     68
Name: count, dtype: int64
Average body mass: 4201.8 grams

1.12 Chapter review

1.13 Exercises

Work through these short exercises to test your understanding. The Chapter 1 Playground provides interactive browser versions of these exercises with hints and complete solutions.

1. A supervisor sleeps 8, 7.5, 6, and 7 hours on four consecutive nights. Store these four numbers in a list named supervisor_sleep, calculate the mean sleep duration, and round the result to one decimal place.
2. Check the data types of "600" and 600 using Python’s type() function. Explain why adding 1 to one succeeds while adding 1 to the other raises a TypeError.
3. Explain the difference between Python’s built-in None and pandas/numpy’s NaN. How does pandas .mean() handle missing values by default?
4. Using students = dt.load_students(), find the age of the youngest student and the oldest student in the study using .min() and .max().
5. Using .value_counts(), find how many students are enrolled full-time versus part-time (the variable is study_mode).
6. Reflect on an analysis you have previously performed in Excel, SPSS, or another point-and-click tool. List three specific analytical decisions that a thesis examiner could not verify without an explicit code script.

1.14 Further reading

• McKinney, W. (2022). Python for Data Analysis (3rd ed.). O’Reilly Media. (Written by the creator of pandas).
• VanderPlas, J. (2016). Python Data Science Handbook. O’Reilly Media.
• The official VS Code Python Tutorial.

2.1 Data as measurement

Almost all empirical research data can be conceptualized as a rectangular grid: each row represents a case (or observation), and each column represents a variable. A case is the individual unit under observation—a student, a patient, an experimental plot, or an institution. A variable is a measured property recorded across all cases, such as chronological age, anxiety score, or study hours. Every cell in the grid contains a single value representing the measurement of that variable for that specific case. This structure, known as the data matrix, forms the backbone of statistical modeling.

What constitutes a case depends on the research design, defining the unit of analysis. Elaf’s wellbeing study involves multiple units of analysis across different tables. In the students table, each row represents an individual student (\(N = 600\)). In the semesters table, however, each row represents a single student in a specific academic semester, meaning each student appears up to four times:

semesters.head(5)

The first four rows all correspond to student S0001, documenting their trajectory across semesters 1 through 4. Confusing the unit of analysis is one of the most frequent errors in quantitative research. A demographic question concerning baseline enrollment must count each student exactly once. Conversely, an inquiry into GPA trajectories across time requires repeated-measures observations, coupled with statistical techniques that account for within-subject dependency (?@sec-ch10). In ?@sec-ch05, we explore how research hypotheses dictate the appropriate unit of analysis.

2.2 Variables and their values

In Python’s data science ecosystem, a single variable’s values are held in a Pandas Series (or a NumPy 1D array). A Series is an indexed, one-dimensional array capable of holding numbers, text, or boolean flags:

sleep = pd.Series([6.5, 7.0, 5.5, 8.0, 6.0])
sleep

0    6.5
1    7.0
2    5.5
3    8.0
4    6.0
dtype: float64

A Series enforces a single homogeneous data type (dtype). If you inadvertently mix numbers with text, Pandas converts the entire Series to an object (or string) dtype to accommodate the text:

mixed = pd.Series([6.5, "seven", 5.5])
mixed

0      6.5
1    seven
2      5.5
dtype: object

Notice that the numbers 6.5 and 5.5 have become text strings. They can no longer be averaged or plotted numerically. This behavior is vital to understand: if a single survey respondent types "seven" or "N/A" into a numeric field in a web form, the entire column imports into Python as text. We will encounter this exact issue when inspecting the raw Excel export later in this chapter.

sleep * 60

0    390.0
1    420.0
2    330.0
3    480.0
4    360.0
dtype: float64

sleep < 7

0     True
1    False
2     True
3    False
4     True
dtype: bool

print("Count of short nights:", (sleep < 7).sum())
print("Proportion of short nights:", (sleep < 7).mean())

Count of short nights: 3
Proportion of short nights: 0.6

print("First night:", sleep.iloc[0])
print("First and third nights:\n", sleep.iloc[[0, 2]])

First night: 6.5
First and third nights:
 0    6.5
2    5.5
dtype: float64

sleep[sleep < 7]

0    6.5
2    5.5
4    6.0
dtype: float64

sleep.iloc[1] = 7.5
sleep

0    6.5
1    7.5
2    5.5
3    8.0
4    6.0
dtype: float64

2.3 Categories (Categorical Data)

Research datasets frequently contain qualitative categories: faculty division, experimental arm, gender identity, or Likert response levels. Pandas represents categorical variables using the specialized category dtype (similar to categorical variables in statistics or value labels in SPSS). A categorical Series displays text labels but stores them internally as compact integers mapped to a defined list of categories:

faculty = pd.Series(
    ["Education", "Humanities", "Education", "Health Sciences"],
    dtype="category"
)
print(faculty)
print("Categories:", faculty.cat.categories)

0          Education
1         Humanities
2          Education
3    Health Sciences
dtype: category
Categories (3, str): ['Education', 'Health Sciences', 'Humanities']
Categories: Index(['Education', 'Health Sciences', 'Humanities'], dtype='str')

By default, categories are arranged in alphanumeric order. However, when categories possess an intrinsic hierarchy (such as employment intensity or Likert scales), you can define an ordered category using pd.Categorical:

employment = pd.Series(
    pd.Categorical(
        ["None", "Full-time job", "Part-time job", "None"],
        categories=["None", "Part-time job", "Full-time job"],
        ordered=True
    )
)
employment.value_counts(sort=False)

None             2
Part-time job    1
Full-time job    1
Name: count, dtype: int64

Specifying explicit category orders ensures that subsequent summary tables and plots display categories logically from lowest to highest rather than alphabetically.

2.4 A table of cases: The DataFrame

Most empirical data is organized as a two-dimensional table where each row represents a case and each column represents a variable. In Python, this table is a Pandas DataFrame. Each column is a Series sharing a common index. While each column must have a uniform data type, different columns across the DataFrame can store different types.

You can construct a DataFrame manually from a Python dictionary:

pilot = pd.DataFrame({
    "student": ["S1", "S2", "S3", "S4"],
    "faculty": pd.Categorical(["Education", "Humanities", "Education", "Health Sciences"]),
    "sleep": [6.5, 7.5, 5.5, 8.0],
    "invited": [True, False, True, False]
})
pilot

The most comprehensive way to inspect a DataFrame’s architecture is the .info() method, which summarizes column names, non-null counts, and data types:

pilot.info()

<class 'pandas.DataFrame'>
RangeIndex: 4 entries, 0 to 3
Data columns (total 4 columns):
 #   Column   Non-Null Count  Dtype   
---  ------   --------------  -----   
 0   student  4 non-null      str     
 1   faculty  4 non-null      category
 2   sleep    4 non-null      float64 
 3   invited  4 non-null      bool    
dtypes: bool(1), category(1), float64(1), str(1)
memory usage: 228.0 bytes

pilot["sleep"]

0    6.5
1    7.5
2    5.5
3    8.0
Name: sleep, dtype: float64

# Select the second row across all columns
pilot.iloc[1, :]

# Select all rows for two specific columns
pilot.loc[:, ["student", "sleep"]]

# Filter rows where sleep is under 7 hours
pilot.loc[pilot["sleep"] < 7, :]

pilot["short_sleep"] = pilot["sleep"] < 7
pilot

students.info()

<class 'pandas.DataFrame'>
RangeIndex: 600 entries, 0 to 599
Data columns (total 14 columns):
 #   Column               Non-Null Count  Dtype  
---  ------               --------------  -----  
 0   student_id           600 non-null    str    
 1   supervisor_id        600 non-null    str    
 2   age                  599 non-null    float64
 3   gender               600 non-null    str    
 4   faculty              600 non-null    str    
 5   programme            600 non-null    str    
 6   study_mode           600 non-null    str    
 7   employment           600 non-null    str    
 8   has_children         600 non-null    str    
 9   lives_away           600 non-null    str    
 10  financial_worry      562 non-null    float64
 11  workshop             600 non-null    str    
 12  workshop_sessions    600 non-null    int64  
 13  considering_dropout  600 non-null    str    
dtypes: float64(2), int64(1), str(11)
memory usage: 65.8 KB

students["faculty"] = students["faculty"].astype("category")
students["faculty"].cat.categories

Index(['Education', 'Health Sciences', 'Humanities', 'Natural Sciences',
       'Social Sciences'],
      dtype='str')

hs_phd = students[
    (students["faculty"] == "Health Sciences") & (students["programme"] == "PhD")
]
print("Health Sciences PhD candidates:", len(hs_phd))

Health Sciences PhD candidates: 44

2.5 Matrices and Dictionaries

Beyond Series and DataFrames, two additional data structures frequently appear in quantitative workflows: 2D numerical matrices and dictionaries.

first_semester = semesters[semesters["semester"] == 1]
metric_vars = first_semester[["gpa", "sleep_hours", "study_hours", "wellbeing"]]

# Compute bivariate correlation matrix
correlations = metric_vars.corr()
round(correlations, 2)

correlations.loc["sleep_hours", "wellbeing"]

np.float64(0.45630448626446257)

student = {
    "id": "S0001",
    "programme": "Master's",
    "sleep": [7.2, 8.4, 7.8, 7.2]
}
student["sleep"]

[7.2, 8.4, 7.8, 7.2]

import scipy.stats as stats

result = stats.ttest_1samp(first_semester["sleep_hours"], popmean=7.0, nan_policy="omit")
print("t-statistic:", round(result.statistic, 2))
print("p-value:", result.pvalue)

t-statistic: -12.18
p-value: 1.331700973081752e-30

2.6 Importing research data

Different software packages store tabular data in distinct binary or text formats. Python provides dedicated reader functions for every standard format, each parsing the input directly into a clean Pandas DataFrame. The sample files below are installed with data2thesis and located via dt.example_file():

path = dt.example_file("semesters.csv")
semester_data = pd.read_csv(path)
semester_data.head(3)

raw_path = dt.example_file("wellbeing_raw.xlsx")
raw_data = pd.read_excel(raw_path)
print("Dimensions (rows, columns):", raw_data.shape)

Dimensions (rows, columns): (615, 64)

raw_data.iloc[:5, :6]

import pyreadstat

spss_path = dt.example_file("wellbeing.sav")
spss_df, meta = pyreadstat.read_sav(spss_path)

spss_df.loc[:2, ["student_id", "gender", "faculty", "stress_1"]]

print("Question text for stress_1:", meta.column_names_to_labels["stress_1"])

Question text for stress_1: I feel unable to control important things in my studies.

spss_labeled = pd.read_spss(spss_path)
spss_labeled["faculty"].value_counts()

faculty
Health Sciences     154
Education           148
Social Sciences     116
Humanities           95
Natural Sciences     87
Name: count, dtype: int64

2.7 Automated codebook generation

In an empirical dataset, a cryptic column header like stress_1 or a numeric value like 4 is scientifically uninterpretable on its own. It acquires meaning only when accompanied by documentation specifying the original prompt (“I feel unable to control important things in my studies”), the response options (1 = Strongly Disagree to 5 = Strongly Agree), and missing data conventions. A codebook provides this reference, serving as the bridge between participant questionnaires and statistical analysis.

Because pyreadstat parses embedded SPSS metadata, we can programmatically assemble a comprehensive codebook DataFrame:

codebook = pd.DataFrame({
    "variable": meta.column_names,
    "label": meta.column_labels,
    "dtype": [str(spss_df[col].dtype) for col in meta.column_names]
})
codebook.head(10)

For categorical items, the mapping between underlying integer codes and verbal labels is preserved in meta.variable_value_labels:

print("Value labels for faculty:", meta.variable_value_labels.get("faculty"))

Value labels for faculty: {1.0: 'Education', 2.0: 'Health Sciences', 3.0: 'Humanities', 4.0: 'Natural Sciences', 5.0: 'Social Sciences'}

Exporting this table into an appendix provides an automated, error-free documentation trail for your thesis committee and future secondary analysts.

2.8 Exporting data

Once raw data has been cleaned, validated, or transformed, the resulting DataFrame can be exported into standard file formats:

# Export to CSV (suppressing automatic integer index column)
students.to_csv("students_clean.csv", index=False)

# Export to Excel workbook
students.to_excel("students_clean.xlsx", index=False)

# Export to SPSS system file with pyreadstat
import pyreadstat
pyreadstat.write_sav(students, "students_clean.sav")

Setting index=False is essential when saving CSV and Excel files; otherwise, Pandas inserts an unlabelled leading column of row numbers into the file.

[!IMPORTANT] Preserve Raw Data Immutability: Never overwrite primary data files. Save processed datasets under descriptive filenames (e.g., students_clean.csv), ensuring that your Python script can recreate the clean data from raw inputs at any time.

2.9 Mapping SPSS and Excel workflows to Python

Moving from point-and-click tools to scripted programming involves mapping familiar conceptual operations to Python commands (Table 2.2).

The paramount advantage of Python is transparency. In Excel or SPSS Data View, manual edits and filter clicks are transient and undocumented. In Python, the input dataset remains untouched while your script provides an indelible, auditable record of every transformation.

2.10 Reading error messages and tracebacks

Python tracebacks are diagnostic reports explaining precisely where and why execution halted. When an error occurs, scroll to the final line of the traceback to read the specific exception type and description. Four exceptions account for the vast majority of beginner stumbling blocks:

1. NameError: name 'slep' is not defined A variable name is misspelled, or the code defining it was not executed prior to this line:

print(slep)

---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
Cell In[33], line 1
----> 1 print(slep)

NameError: name 'slep' is not defined

2. AttributeError or KeyError Referencing a non-existent column name inside brackets raises a KeyError (often caused by capitalisation mismatches):

students["Faculty"]

---------------------------------------------------------------------------
KeyError                                  Traceback (most recent call last)
File /opt/hostedtoolcache/Python/3.12.15/x64/lib/python3.12/site-packages/pandas/core/indexes/base.py:3816, in Index.get_loc(self, key)
   3815 try:
-> 3816     return self._engine.get_loc(casted_key)
   3817 except KeyError as err:

File pandas/_libs/index.pyx:172, in pandas._libs.index.IndexEngine.get_loc()
--> 172 'Could not get source, probably due dynamically evaluated source code.'

File pandas/_libs/index.pyx:201, in pandas._libs.index.IndexEngine.get_loc()
--> 201 'Could not get source, probably due dynamically evaluated source code.'

File pandas/_libs/hashtable_class_helper.pxi:7726, in pandas._libs.hashtable.PyObjectHashTable.get_item()
-> 7726 'Could not get source, probably due dynamically evaluated source code.'

File pandas/_libs/hashtable_class_helper.pxi:7734, in pandas._libs.hashtable.PyObjectHashTable.get_item()
-> 7734 'Could not get source, probably due dynamically evaluated source code.'

KeyError: 'Faculty'

The above exception was the direct cause of the following exception:

KeyError                                  Traceback (most recent call last)
Cell In[34], line 1
----> 1 students["Faculty"]

File /opt/hostedtoolcache/Python/3.12.15/x64/lib/python3.12/site-packages/pandas/core/frame.py:4330, in DataFrame.__getitem__(self, key)
   4326 
   4327         if is_single_key:
   4328             if self.columns.nlevels > 1:
   4329                 return self._getitem_multilevel(key)
-> 4330             indexer = self.columns.get_loc(key)
   4331             if is_integer(indexer):
   4332                 indexer = [indexer]  # type: ignore[assignment]
   4333         else:

File /opt/hostedtoolcache/Python/3.12.15/x64/lib/python3.12/site-packages/pandas/core/indexes/base.py:3823, in Index.get_loc(self, key)
   3818     if isinstance(casted_key, slice) or (
   3819         isinstance(casted_key, abc.Iterable)
   3820         and any(isinstance(x, slice) for x in casted_key)
   3821     ):
   3822         raise InvalidIndexError(key) from err
-> 3823     raise KeyError(key) from err
   3824 except TypeError:
   3825     # If we have a listlike key, _check_indexing_error will raise
   3826     #  InvalidIndexError. Otherwise we fall through and re-raise
   3827     #  the TypeError.
   3828     self._check_indexing_error(key)

KeyError: 'Faculty'

3. TypeError: unsupported operand type(s) Attempting mathematical operations on incompatible data types (for example, multiplying text strings):

"6.5" * "2"

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
Cell In[35], line 1
----> 1 "6.5" * "2"

TypeError: can't multiply sequence by non-int of type 'str'

4. FileNotFoundError: No such file or directory Python searched the current working directory but could not locate the target filename. Verify your relative path against the folder opened in VS Code.

2.11 Descriptive summary of cohort measurements

With the data structures in place, Pandas provides .describe() for a rapid multi-variable statistical overview of central tendency, dispersion, and range:

first_semester[
    ["gpa", "sleep_hours", "study_hours", "wellbeing"]
].describe().round(2)

This single command confirms that in semester 1, students slept an average of 6.48 hours (standard deviation = 0.94), well below recommended benchmarks. In ?@sec-ch06, we expand these basic summaries into publication-ready descriptive tables.

import seaborn as sns

penguins = sns.load_dataset("penguins")
penguins.info()
penguins[["bill_length_mm", "bill_depth_mm", "body_mass_g"]].describe().round(1)

<class 'pandas.DataFrame'>
RangeIndex: 344 entries, 0 to 343
Data columns (total 7 columns):
 #   Column             Non-Null Count  Dtype  
---  ------             --------------  -----  
 0   species            344 non-null    str    
 1   island             344 non-null    str    
 2   bill_length_mm     342 non-null    float64
 3   bill_depth_mm      342 non-null    float64
 4   flipper_length_mm  342 non-null    float64
 5   body_mass_g        342 non-null    float64
 6   sex                333 non-null    str    
dtypes: float64(4), str(3)
memory usage: 18.9 KB

2.12 Chapter review

2.13 Exercises

Work through these exercises in your browser using the Chapter 2 Playground or download the starter project for VS Code.

1. Create a Series with participant ages: 24, 31, 28, 45, 26. Use boolean indexing to filter ages over 30, and compute how many participants exceed 30.
2. Construct an ordered Categorical Series from ["Agree", "Disagree", "Neutral", "Agree"] with categories ordered logically from "Disagree" to "Neutral" to "Agree". Verify the ordering using .value_counts(sort=False).
3. In students, use boolean indexing with & to count how many participants study part-time (study_mode == "Part-time") AND have dependent children (children == "Yes").
4. Add a boolean column over_30 to students that evaluates to True for students older than 30. Check how many missing values exist using .isna().sum().
5. Import wellbeing.sav with pyreadstat and extract the full survey question text mapped to the variable support_3.
6. Contrast the unit of analysis in students versus counselling_visits. Formulate one distinct research question appropriate for each dataset.

2.14 Further reading

• Python for Data Analysis (3rd ed.) by Wes McKinney (McKinney 2022): the authoritative guide to Pandas indexing, Series arithmetic, and data wrangling by the creator of Pandas.
• Official documentation for pandas.pydata.org and pyreadstat.

3.1 What clean data means

Data can be considered clean when it faithfully represents the empirical phenomena under investigation. This standard encompasses five distinct qualities:

1. Validity: Every recorded value conforms to the permissible domain of its variable (e.g., chronological age cannot be 250 years, and daily sleep cannot exceed 24 hours).
2. Accuracy: Values capture what respondents actually intended, ensuring system codes (such as 99 for skipped items) are not treated as genuine quantitative answers.
3. Completeness: Data captures all eligible participants, and omitted observations are explicitly marked as missing (NaN) rather than masquerading as numbers.
4. Consistency: Categorical responses follow a single standardized nomenclature (e.g., "F", "female", and "Female " are normalized into "Female").
5. Uniqueness: Each participant appears exactly once per unit of analysis, with zero test records and zero duplicate submissions.

Violating these criteria distorts empirical conclusions. The impact of unhandled missing-data codes is immediately obvious. Suppose five students answer an item on a 1-to-5 scale, and one student leaves the question blank, which the web software logs as 99:

answers = pd.Series([4, 3, 5, 99, 2])
print("Uncleaned mean:", answers.mean())

clean_answers = answers.replace(99, np.nan)
print("Cleaned mean:", clean_answers.mean())

Uncleaned mean: 22.6
Cleaned mean: 3.5

Leaving 99 in the calculation produces a mean of 22.6, a mathematical impossibility on a 5-point scale. When 99 is properly converted to NaN, Pandas excludes the missing entry by default, yielding the true average of 3.5. In large datasets, such distortions are often smaller and harder to detect visually, making automated cleaning pipelines essential. Duplications exert a similar distorting effect: duplicating a participant gives them double weight in regressions and artificially narrows confidence intervals.

Three fundamental rules govern scientific data cleaning:

• Raw files remain immutable: Never edit primary data files directly in Excel. Raw files must remain intact, exactly as received from instruments or survey platforms. All transformations are executed via code.
• Every decision is documented: Your Python script serves as an auditable record of all exclusions, recodings, imputations, and composite score derivations.
• Inspect results at every step: Validate table shapes, unique identifiers, and value distributions after each operation, because many data-cleaning bugs produce no syntax errors.

3.2 Tidy data principles in Python

Tabular data should be organized to make analysis straightforward. The standard architecture in modern data science is tidy data (Wickham 2014): 1. Each variable forms a column. 2. Each observation forms a row. 3. Each type of observational unit forms a table.

Consider Table 3.1, which records two semesters of wellbeing scores for three students side by side, as spreadsheet exports often do.

While compact for human reading, this layout spreads a single construct (wellbeing) across multiple columns, while embedding a second variable (semester) inside column headers. In contrast, Table 3.2 presents the identical measurements in tidy (long) format.

In tidy format, each measurement occupies its own row, and computing grouped summaries—such as mean wellbeing per semester—becomes straightforward. The students table in data2thesis is tidy at the participant level (\(N = 600\)), whereas semesters is tidy at the semester observation level (\(N = 2,326\)). Uncurated survey exports, by contrast, frequently pack longitudinal observations side by side. Much of data preparation consists of reshaping data into tidy tables.

3.3 The Pandas manipulation toolkit

In Python, tabular manipulation centers on Pandas. To keep transformations readable and eliminate intermediate variables, we use method chaining: calling sequential methods on a DataFrame enclosed within parentheses (...).

Let’s illustrate on a small pilot cohort:

pilot = pd.DataFrame({
    "student": ["S1", "S2", "S3", "S4", "S5"],
    "faculty": ["Education", "Humanities", "Education", "Health Sciences", "Humanities"],
    "sleep":   [6.5, 7.5, 5.5, 8.0, 6.0],
    "stress":  [3.2, 2.1, 4.5, 1.8, 3.9]
})
pilot

# Boolean indexing
pilot[pilot["sleep"] < 7]

# Multiple conditions combined with & (AND)
pilot[(pilot["sleep"] < 7) & (pilot["faculty"] == "Education")]

# Matching against a collection using .isin()
pilot[pilot["faculty"].isin(["Education", "Humanities"])]

pilot.query("sleep < 7 and faculty == 'Education'")

# Keep only student ID and sleep hours
pilot[["student", "sleep"]]

# Drop the faculty column
pilot.drop(columns=["faculty"])

# Sort by stress in descending order
pilot.sort_values("stress", ascending=False)

pilot = pilot.assign(
    sleep_minutes = pilot["sleep"] * 60,
    short_sleep   = pilot["sleep"] < 7
)
pilot[["student", "sleep", "sleep_minutes", "short_sleep"]]

conditions = [
    pilot["stress"] >= 4.0,
    pilot["stress"] >= 2.5
]
choices = ["High", "Medium"]

pilot["stress_level"] = np.select(conditions, choices, default="Low")
pilot[["student", "stress", "stress_level"]]

# Overall cohort summary
print("Mean sleep across all pilot cases:", pilot["sleep"].mean())

# Grouped summary by faculty
faculty_summary = (
    pilot
    .groupby("faculty")["sleep"]
    .agg(mean_sleep="mean", count="count")
    .reset_index()
)
faculty_summary

Mean sleep across all pilot cases: 6.7

pilot["faculty"].value_counts().reset_index()

wellbeing_trajectory = (
    semesters
    .merge(students[["student_id", "workshop"]], on="student_id", how="left")
    .groupby(["semester", "workshop"])["wellbeing"]
    .agg(mean_wellbeing="mean", n_students="count")
    .reset_index()
    .sort_values(["semester", "workshop"])
)
wellbeing_trajectory

3.4 Reshaping data between wide and long

Data is often collected in wide format but analyzed in long format. In wide format, repeated measurements appear side by side across columns. In long format, each measurement occupies its own row.

Consider the wide table from Table 3.1:

wide = pd.DataFrame({
    "student":     ["S1", "S2", "S3"],
    "wellbeing_1": [58, 71, 64],
    "wellbeing_2": [63, 70, 69]
})
wide

To convert wide columns into tidy rows, use pd.melt():

long = pd.melt(
    wide,
    id_vars=["student"],
    value_vars=["wellbeing_1", "wellbeing_2"],
    var_name="semester",
    value_name="wellbeing"
)
# Clean the semester column to extract the integer number
long["semester"] = long["semester"].str.replace("wellbeing_", "").astype(int)
long = long.sort_values(["student", "semester"]).reset_index(drop=True)
long

The arguments specify: - id_vars: Columns to preserve as identifiers (student). - value_vars: Columns to unpivot into rows. - var_name: Name for the new column storing original column headers (semester). - value_name: Name for the new column storing the data values (wellbeing).

Conversely, when formatting summary tables for a thesis chapter, wide presentation is often clearer for human readers. Use .pivot() to spread long data across columns:

wide_presentation = long.pivot(index="student", columns="semester", values="wellbeing").reset_index()
wide_presentation

We can apply this to our longitudinal trajectory summary, spreading the workshop groups into adjacent columns:

(wellbeing_trajectory
 .assign(mean_wellbeing=lambda df: df["mean_wellbeing"].round(1))
 .pivot(index="semester", columns="workshop", values="mean_wellbeing")
 .reset_index()
)

3.5 Joining relational tables

In empirical research, data is frequently split across multiple tables sharing a common key: participant demographics in one file, longitudinal outcomes in another, and supervisor attributes in a third. In Pandas, relational tables are combined using pd.merge().

Consider two simple tables:

people = pd.DataFrame({
    "student": ["S1", "S2", "S3"],
    "faculty": ["Education", "Humanities", "Education"]
})

scores = pd.DataFrame({
    "student":   ["S1", "S1", "S2", "S4"],
    "wellbeing": [58, 63, 71, 66]
})

A left join (how="left") preserves all rows from the primary table, matching attributes from the secondary table where the keys align:

scores.merge(people, on="student", how="left")

Student S4 has a score but no corresponding demographic record in people, so their faculty is assigned NaN.

To find cases that have no matching counterpart in another table (an anti-join), use boolean negation ~ with .isin():

# Find students in 'people' who have no records in 'scores'
missing_scores = people[~people["student"].isin(scores["student"])]
missing_scores

Student S3 enrolled in the study but contributed zero scores.

In Elaf’s research cohort, merging connects supervisory rank (from supervisors) to student dropout intentions (from students):

supervisors = dt.supervisors

by_rank = (
    students
    .merge(supervisors[["supervisor_id", "rank"]], on="supervisor_id", how="left")
    .assign(dropout_flag = lambda df: (df["considering_dropout"] == "Yes").astype(int))
    .groupby("rank")
    .agg(
        dropout_rate=("dropout_flag", "mean"),
        n_students=("student_id", "count")
    )
    .reset_index()
)
by_rank

Similarly, set filtering identifies students who dropped out of the longitudinal study by Semester 3:

sem3_students = semesters[semesters["semester"] == 3]["student_id"]
attrition_students = students[~students["student_id"].isin(sem3_students)]

print("Students who left before Semester 3:", len(attrition_students))
print("\nDropout thoughts among those who left:")
print(attrition_students["considering_dropout"].value_counts())

Students who left before Semester 3: 37

Dropout thoughts among those who left:
considering_dropout
Yes    22
No     15
Name: count, dtype: int64

Of the 37 students who departed prior to Semester 3, the vast majority had indicated on their baseline questionnaire that they were considering discontinuing their studies. Participant attrition is rarely random; tracking who dropped out is essential context for interpreting later semesters (?@sec-ch06).

3.6 Cleaning an uncurated survey export

We are now equipped to sanitize the raw wellbeing export (wellbeing_raw.xlsx), resolving the real-world data issues introduced in ?@sec-ch02.

raw = pd.read_excel(dt.example_file("wellbeing_raw.xlsx"))
print("Raw table shape (rows, columns):", raw.shape)

Raw table shape (rows, columns): (615, 64)

print(raw["Q3_Gender"].value_counts(dropna=False).head(8))
print("\nSample test records:")
print(raw[raw["Q1_Student ID"].str.upper().str.startswith("TEST", na=False)][["Q1_Student ID", "Q2_Age", "Q3_Gender"]].head(3))

Q3_Gender
Female     203
Male       152
male        78
female      56
M           46
Female      30
F           28
Male        21
Name: count, dtype: int64

Sample test records:
  Q1_Student ID  Q2_Age Q3_Gender
0          TEST    99.0    Female
1          test     NaN       NaN
2         TEST2    30.0      Male

# Filter out test accounts
no_test = raw[~raw["Q1_Student ID"].str.upper().str.startswith("TEST", na=False)]

# Drop technical response ID and eliminate duplicate submissions
responses = no_test.drop(columns=["Response ID"]).drop_duplicates()
print("Cleaned participant count:", len(responses))

Cleaned participant count: 600

# Construct questionnaire item names: stress_1..6, burnout_1..6, support_1..6, satisfaction_1..4
item_names = (
    [f"stress_{i}" for i in range(1, 7)] +
    [f"burnout_{i}" for i in range(1, 7)] +
    [f"support_{i}" for i in range(1, 7)] +
    [f"satisfaction_{i}" for i in range(1, 5)]
)

rename_dict = {
    "Q1_Student ID": "student_id",
    "Q2_Age": "age",
    "Q3_Gender": "gender",
    "Q4_Faculty": "faculty",
    "Q5_Programme": "programme",
    "Q6_Study mode": "study_mode",
    "Q7_Paid work": "employment",
    "Q8_Children": "has_children",
    "Q9_Moved away from family": "lives_away",
    "Q10_How worried are you about money? (1-5)": "financial_worry",
    "Workshop group": "workshop",
    "Workshop sessions attended": "workshop_sessions",
    "Y1_Considered leaving?": "considering_dropout"
}

# Add item column mappings Q11_1 through Q11_22
for i, name in enumerate(item_names, 1):
    rename_dict[f"Q11_{i}"] = name

responses = responses.rename(columns=rename_dict)

# Standardize gender
responses["gender"] = np.where(
    responses["gender"].astype(str).str.strip().str.lower().str.startswith("f"),
    "Female",
    "Male"
)

# Standardize faculty using regex classification
fac_lower = responses["faculty"].astype(str).str.lower()
responses["faculty"] = np.select(
    [
        fac_lower.str.contains("educ"),
        fac_lower.str.contains("health"),
        fac_lower.str.contains("humanit"),
        fac_lower.str.contains("social"),
        fac_lower.str.contains("natural|^sciences?$")
    ],
    [
        "Education",
        "Health Sciences",
        "Humanities",
        "Social Sciences",
        "Natural Sciences"
    ],
    default="Other"
)

# Standardize degree programme and study mode
responses["programme"] = np.where(
    responses["programme"].astype(str).str.lower().str.contains("ph"),
    "PhD",
    "Master's"
)
responses["study_mode"] = np.where(
    responses["study_mode"].astype(str).str.lower().str.contains("part"),
    "Part-time",
    "Full-time"
)

# Standardize binary Yes/No indicators
for col in ["has_children", "lives_away", "considering_dropout"]:
    responses[col] = np.where(
        responses[col].astype(str).str.strip().str.lower().str.startswith("y"),
        "Yes",
        "No"
    )

print("Standardized gender counts:")
print(responses["gender"].value_counts())
print("\nStandardized faculty counts:")
print(responses["faculty"].value_counts())

Standardized gender counts:
gender
Female    312
Male      288
Name: count, dtype: int64

Standardized faculty counts:
faculty
Health Sciences     154
Education           148
Social Sciences     116
Humanities           95
Natural Sciences     87
Name: count, dtype: int64

# Replace missing codes (99, -9) in survey items
for col in item_names + ["financial_worry"]:
    responses[col] = responses[col].replace([99, -9, "99", "-9"], np.nan)
    responses[col] = pd.to_numeric(responses[col], errors="coerce").astype("Int64")

responses["workshop_sessions"] = pd.to_numeric(responses["workshop_sessions"], errors="coerce").astype("Int64")

sample_text = pd.Series(["7 hrs", "6.5", "6,5"])
cleaned_numbers = sample_text.str.replace(",", ".").str.extract(r"([\d\.]+)", expand=False).astype(float)
print(cleaned_numbers.tolist())

[7.0, 6.5, 6.5]

# Parse numeric age and set impossible values (> 100) to NaN
responses["age"] = pd.to_numeric(
    responses["age"].astype(str).str.extract(r"([\d\.]+)", expand=False),
    errors="coerce"
)
responses.loc[responses["age"] > 100, "age"] = np.nan

print("Participants with age > 100:", (responses["age"] > 100).sum())
print("Missing age count:", responses["age"].isna().sum())

Participants with age > 100: 0
Missing age count: 1

questionnaire_clean = responses[["student_id"] + item_names].copy()

# Create scale scores
scores = pd.DataFrame({"student_id": questionnaire_clean["student_id"]})

# Reverse-code stress_4
reversed_stress_4 = 6 - questionnaire_clean["stress_4"]

# Compute row-wise means across scale items (skipping missing items)
scores["stress_score"] = pd.concat([
    questionnaire_clean[["stress_1", "stress_2", "stress_3", "stress_5", "stress_6"]],
    reversed_stress_4.rename("stress_4")
], axis=1).mean(axis=1)

scores["burnout_score"] = questionnaire_clean[[f"burnout_{i}" for i in range(1, 7)]].mean(axis=1)
scores["support_score"] = questionnaire_clean[[f"support_{i}" for i in range(1, 7)]].mean(axis=1)
scores["satisfaction_score"] = questionnaire_clean[[f"satisfaction_{i}" for i in range(1, 5)]].mean(axis=1)

scores.head()

import re

sem_cols = [c for c in responses.columns if re.search(r"_S[1-4]$", c)]
sem_wide = responses[["student_id"] + sem_cols]

# Reshape wide columns into tidy rows
sem_long = pd.wide_to_long(
    sem_wide,
    stubnames=["GPA", "Sleep", "Study", "Exercise", "Caffeine", "Meetings", "Wellbeing"],
    i="student_id",
    j="semester",
    sep="_S"
).reset_index()

sem_long = sem_long.rename(columns={
    "GPA": "gpa",
    "Sleep": "sleep_hours",
    "Study": "study_hours",
    "Exercise": "exercise_days",
    "Caffeine": "caffeine_mg",
    "Meetings": "supervisor_meetings",
    "Wellbeing": "wellbeing"
})

# Standardize numeric dtypes and decimal commas
sem_long["sleep_hours"] = (
    sem_long["sleep_hours"]
    .astype(str)
    .str.replace(",", ".", regex=False)
    .str.extract(r"([\d\.]+)", expand=False)
    .astype(float)
)

for col in ["gpa", "study_hours", "exercise_days", "caffeine_mg", "supervisor_meetings", "wellbeing"]:
    sem_long[col] = pd.to_numeric(
        sem_long[col].astype(str).str.extract(r"([\d\.]+)", expand=False),
        errors="coerce"
    )

# Filter out impossible measurement values
sem_long.loc[sem_long["sleep_hours"] > 24, "sleep_hours"] = np.nan
sem_long.loc[sem_long["study_hours"] > 168, "study_hours"] = np.nan
sem_long.loc[sem_long["gpa"] > 4.0, "gpa"] = np.nan

# Drop unobserved rows where all semester measurements are missing (attrition)
measure_cols = ["gpa", "sleep_hours", "study_hours", "exercise_days", "caffeine_mg", "supervisor_meetings", "wellbeing"]
semesters_clean = sem_long.dropna(subset=measure_cols, how="all").copy()
semesters_clean = semesters_clean.sort_values(["student_id", "semester"]).reset_index(drop=True)

print("Cleaned longitudinal semester observations:", len(semesters_clean))

Cleaned longitudinal semester observations: 2326

# Verify questionnaire items match canonical data exactly
canonical_q = dt.questionnaire.sort_values("student_id").reset_index(drop=True)
q_test = questionnaire_clean.sort_values("student_id").reset_index(drop=True)

diffs = sum(
    not (q_test[c].fillna(-999).astype(int) == canonical_q[c].fillna(-999).astype(int)).all()
    for c in item_names
)
print("Questionnaire item discrepancies:", diffs)

# Verify semesters records match canonical data exactly
canonical_sem = dt.semesters.sort_values(["student_id", "semester"]).reset_index(drop=True)
max_sem_diff = max(
    (semesters_clean[c] - canonical_sem[c]).abs().max()
    for c in measure_cols
)
print("Maximum absolute difference across semester measurements:", max_sem_diff)

Questionnaire item discrepancies: 0
Maximum absolute difference across semester measurements: 0.0

# Load airquality via data2thesis or build sample dataframe
air_sample = pd.DataFrame({
    "Ozone": [41.0, 36.0, 12.0, 18.0, np.nan, 28.0, 23.0, 19.0, 8.0, np.nan],
    "Month": [5, 5, 5, 5, 5, 5, 5, 5, 5, 5],
    "Temp": [67, 72, 74, 62, 56, 66, 65, 59, 61, 69]
})

air_summary = (
    air_sample
    .groupby("Month")
    .agg(
        mean_ozone=("Ozone", "mean"),
        missing_days=("Ozone", lambda s: s.isna().sum())
    )
    .reset_index()
)
air_summary

3.7 Chapter review

3.8 Exercises

Complete these exercises in your browser using the Chapter 3 Playground or download the starter project for VS Code.

1. Using students, calculate the average chronological age for each faculty, sorted from oldest to youngest (handling missing values with default skipna=True).
2. Add a categorical column sleep_category to semesters that assigns "Short" when sleep_hours < 6.0, "Recommended" when sleep_hours >= 7.0, and "Moderate" otherwise. Count records per category.
3. Using semesters, create a wide summary table displaying mean GPA with one row per student and one column per semester (Semesters 1 to 4) using .pivot().
4. Merge semesters with students and compare average semester GPA between full-time and part-time students.
5. Explain why item stress_4 must be reverse-coded before computing stress_score, and describe what would happen to empirical correlations if reverse-coding were omitted.
6. For each of the five qualities of clean data, describe a specific data-quality challenge from your own academic discipline and the corresponding Python cleaning step required.

3.9 Further reading

• Python for Data Analysis (3rd ed.) by Wes McKinney (McKinney 2022): chapters on “Data Cleaning and Preparation” and “Data Wrangling: Join, Combine, and Reshape”.
• “Tidy Data” (Wickham 2014): the seminal paper formalizing tabular data organization.
• Official Pandas guide to Reshaping and Pivot Tables.

4.1 What a graph is for

Summary statistics compress empirical distributions into a handful of numerical parameters (means, standard deviations, correlations). That compression can conceal radically different data distributions. The statistician Francis Anscombe demonstrated this vulnerability with four synthetic datasets constructed in 1973 (Anscombe 1973), bundled in Seaborn as anscombe:

anscombe = sns.load_dataset("anscombe")

summary = (
    anscombe
    .groupby("dataset")
    .agg(
        mean_x=("x", "mean"),
        sd_x=("x", "std"),
        mean_y=("y", "mean"),
        sd_y=("y", "std")
    )
    .round(2)
)
summary

Across all four datasets, the summary statistics are practically identical: \(\bar{x} = 9.0\), \(s_x = 3.32\), \(\bar{y} = 7.50\), and \(s_y = 2.03\). Bivariate regression yields identical linear fits (\(y = 3.0 + 0.5x\)) with identical correlation coefficients (\(r = 0.82\)). Judged solely by numerical summaries, the four datasets appear indistinguishable.

Yet Figure 4.1 demonstrates that their underlying relationships are entirely dissimilar.

PYTHON Copy

g = sns.lmplot(
    data=anscombe, x="x", y="y", col="dataset", col_wrap=2,
    ci=None, height=2.5, aspect=1.3,
    scatter_kws={"s": 35, "color": "#1f77b4"},
    line_kws={"color": "grey", "linewidth": 1.5}
)
g.set_axis_labels("X", "Y")
plt.tight_layout()
plt.show()

Only Dataset I represents the bivariate linear association assumed by Pearson’s \(r\). Dataset II is a non-linear quadratic curve, Dataset III is a perfect linear relationship corrupted by one outlier, and Dataset IV is an unassociated vertical cluster whose correlation is produced entirely by a single high-leverage point. An analyst relying strictly on summary statistics would draw identical—and incorrect—conclusions across all four. Visual inspection prior to modeling is not optional; it is a prerequisite for valid inference.

4.2 Graphical perception: How vision decodes data

Graphs function by translating numbers into visual encodings: position, length, angle, area, slope, and color. Human perception does not decode these attributes with equal precision. In pioneering experiments, William Cleveland and Robert McGill (Cleveland and McGill 1984) evaluated how accurately observers judge quantitative differences across visual channels, establishing an empirical hierarchy:

1. Position along a common aligned scale (most accurate)
2. Length
3. Slope and Angle
4. Area
5. Color saturation and hue (least accurate)

Effective scientific graphics map primary research comparisons to position or length, reserving color and area for secondary grouping.

The perceptual difference is readily apparent when comparing identical proportions displayed as a pie chart versus a bar chart (Figure 4.2).

PYTHON Copy

shares = pd.DataFrame({
    "group": ["A", "B", "C", "D", "E"],
    "percent": [23, 21, 20, 19, 17]
})

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(7.5, 3.2))

# Pie chart (angle / area)
ax1.pie(shares["percent"], labels=shares["group"], startangle=90, colors=sns.color_palette("mako", 5))
ax1.set_title("Pie chart (Angle / Area)", fontsize=11)

# Bar chart (aligned length)
sns.barplot(data=shares, x="group", y="percent", color="#4c72b0", ax=ax2)
ax2.set(title="Bar chart (Aligned Length)", xlabel="Category", ylabel="Percentage (%)")
ax2.set_ylim(0, 26)

plt.tight_layout()
plt.show()

In the pie chart, the five slices appear virtually identical; determining the rank order requires scrutinizing legends and mentally comparing angles. In the bar chart, the monotonic decline from category A (23%) to E (17%) is evident immediately because comparisons occur along a common baseline. Consequently, scientific reporting favors bar charts and dot plots over pie charts.

4.3 Selecting the right visualization

The choice of graph is governed by the research question and the levels of measurement of the variables under study (Table 4.1).

4.4 Core statistical visualizations with Seaborn

To explore Elaf’s wellbeing cohort, we merge Semester 1 measurements with baseline student background demographics:

first_sem = (
    semesters
    .query("semester == 1")
    .merge(students, on="student_id", how="left")
)

4.4.1 Histograms and density curves

The initial question concerning any continuous variable addresses its empirical distribution: central tendency, spread, skewness, and the presence of extreme outliers. A histogram partitions the continuous range into discrete bins, displaying the count of observations falling within each interval:

PYTHON Copy

fig, ax = plt.subplots(figsize=(6, 3.5))

sns.histplot(
    data=first_sem, x="sleep_hours", binwidth=0.5,
    color="#4c72b0", edgecolor="white", ax=ax
)
ax.set(xlabel="Hours of sleep per night", ylabel="Number of students")
plt.tight_layout()
plt.show()

Figure 4.3 displays a roughly symmetrical, unimodal distribution centered near 6.5 hours. Exactly r round(100 * share_short, 1)% of graduate students sleep less than the clinically recommended 7 hours during their first semester.

4.4.2 Bar charts for categorical counts

A bar chart plots frequency distributions across categorical factor levels:

PYTHON Copy

fig, ax = plt.subplots(figsize=(6, 3.5))

sns.countplot(
    data=students, x="faculty", color="#4c72b0",
    order=students["faculty"].value_counts().index, ax=ax
)
ax.set(xlabel="Faculty", ylabel="Student Count")
plt.xticks(rotation=20, ha="right")
plt.tight_layout()
plt.show()

Because bar lengths represent quantities proportionally, the numerical baseline must always start at zero.

4.4.3 Box plots for group comparisons

A box plot summarizes five key distribution parameters: the median (center line), interquartile range (IQR, box edges covering the 25th to 75th percentiles), whiskers (extending to \(1.5 \times \text{IQR}\)), and individual points for potential outliers. It is the gold standard for comparing continuous metrics across qualitative groups:

PYTHON Copy

fig, ax = plt.subplots(figsize=(5, 3.5))

sns.boxplot(
    data=first_sem, x="study_mode", y="wellbeing",
    palette=["#4c72b0", "#55a868"], ax=ax, width=0.5
)
ax.set(xlabel="Study Mode", ylabel="Wellbeing Index (0 to 100)")
plt.tight_layout()
plt.show()

Part-time students report a slightly lower median wellbeing score than full-time students, but the distributions overlap substantially. In ?@sec-ch08, we use inferential ANOVA to test whether this group discrepancy is statistically significant.

4.4.4 Scatter plots with regression trend lines

The relationship between two continuous variables—such as sleep hours and academic performance—is evaluated via a scatter plot. When plotting 600 participants, individual markers often overlap (overplotting). We mitigate this by setting transparency (alpha=0.4) and overlaying a linear regression trend line with sns.regplot():

PYTHON Copy

fig, ax = plt.subplots(figsize=(6, 3.8))

sns.regplot(
    data=first_sem, x="sleep_hours", y="gpa",
    scatter_kws={"alpha": 0.35, "color": "#4c72b0"},
    line_kws={"color": "#c44e52", "linewidth": 2},
    ax=ax
)
ax.set(xlabel="Hours of sleep per night", ylabel="Grade Point Average (GPA)")
plt.tight_layout()
plt.show()

The regression line slopes upward: students sleeping more hours tend to achieve modestly higher GPAs. However, observations scatter broadly around the line (\(r \approx 0.30\)), indicating that sleep explains only a fraction of total academic variance.

4.4.5 Longitudinal trajectory line plots

Elaf’s research design tracked cohorts across four consecutive semesters. A line plot illustrates how group means evolve over time:

PYTHON Copy

wellbeing_trajectory = (
    semesters
    .merge(students[["student_id", "workshop"]], on="student_id", how="left")
    .groupby(["semester", "workshop"])["wellbeing"]
    .mean()
    .reset_index()
)

fig, ax = plt.subplots(figsize=(6.5, 3.8))

sns.lineplot(
    data=wellbeing_trajectory, x="semester", y="wellbeing",
    hue="workshop", marker="o", markersize=7, linewidth=2,
    palette={"Invited": "#2ca02c", "Not invited": "#7f7f7f"},
    ax=ax
)
ax.set(
    title="Cohort Wellbeing Trajectory Over Two Years",
    xlabel="Academic Semester",
    ylabel="Mean Wellbeing (0 to 100)",
    xticks=[1, 2, 3, 4]
)
ax.legend(title="Workshop Group", frameon=True)
plt.tight_layout()
plt.show()

Figure 4.7 captures the core narrative of Elaf’s intervention study: groups start at identical baselines in Semester 1, invited students pull ahead following the intervention in Semester 2, and the gap stabilizes over subsequent terms.

4.5 Preparing publication-quality figures

Exploratory graphics suffice during initial data exploration, but explanatory figures destined for a thesis or journal article require publication-standard polish.

4.5.1 Annotations and reference lines

Adding explicit reference lines and annotations directs reader attention to clinically meaningful thresholds:

PYTHON Copy

fig, ax = plt.subplots(figsize=(6.5, 3.8))

sns.histplot(
    data=first_sem, x="sleep_hours", binwidth=0.5,
    color="#4c72b0", edgecolor="white", ax=ax
)
ax.axvline(7.0, color="#c44e52", linestyle="--", linewidth=1.8)
ax.text(
    7.1, ax.get_ylim()[1] * 0.9,
    "Recommended: 7.0 hrs",
    color="#c44e52", fontweight="semibold", verticalalignment="top"
)
ax.set(
    title="Baseline Sleep Duration",
    xlabel="Hours of sleep per night",
    ylabel="Number of students"
)
plt.tight_layout()
plt.show()

4.6 Honest graphs and visual ethics

A visualization can be mathematically accurate while remaining visually misleading. Distortion frequently stems from unexamined software defaults:

1. Truncated vertical baselines on bar charts: Because a bar encodes magnitude through length, truncating the baseline exaggerates differences.
2. Extreme aspect ratios: Altering aspect ratios artificially steepens or flattens longitudinal trends.
3. Chartjunk and 3D effects: Pseudo-3D perspectives distort angular perception and distance.
4. Hiding underlying distribution variance: Reporting only summary bar heights masks critical within-group heterogeneity.

4.6.1 Redesigning a misleading figure

Suppose an investigator evaluates whether mean wellbeing differs across faculties. A careless chart with an artificially truncated axis appears in Figure 4.9:

PYTHON Copy

faculty_means = (
    first_sem
    .groupby("faculty")["wellbeing"]
    .mean()
    .reset_index()
)

fig, ax = plt.subplots(figsize=(6.5, 3.5))
sns.barplot(data=faculty_means, x="faculty", y="wellbeing", palette="Blues_d", ax=ax)
ax.set_ylim(58, 62.5)  # Truncated baseline!
ax.set(title="Misleading: Truncated Baseline (58 to 62.5)", xlabel="Faculty", ylabel="Mean Wellbeing")
plt.xticks(rotation=20, ha="right")
plt.tight_layout()
plt.show()

In Figure 4.9, the Education faculty bar appears triple the height of Social Sciences. Yet the actual group means differ by barely 2.5 points on a 100-point scale. Furthermore, the bar chart completely conceals individual variation.

In Figure 4.10, the display is redesigned into an honest, informative jittered strip plot with overlaid group means:

PYTHON Copy

order = faculty_means.sort_values("wellbeing")["faculty"]

fig, ax = plt.subplots(figsize=(7, 3.8))

# Jittered individual student scores
sns.stripplot(
    data=first_sem, y="faculty", x="wellbeing",
    order=order, color="grey", alpha=0.3, jitter=0.2, size=5, ax=ax
)

# Overlay group means as prominent diamonds
sns.pointplot(
    data=first_sem, y="faculty", x="wellbeing",
    order=order, color="#c44e52", markers="D", scale=0.8,
    errorbar=None, join=False, ax=ax
)

ax.set(
    title="Truthful: Student Wellbeing Distribution by Faculty",
    xlabel="Wellbeing Index (0 to 100)",
    ylabel=None,
    xlim=(0, 100)  # Full theoretical scale!
)
plt.tight_layout()
plt.show()

Figure 4.10 communicates the true scientific story: individual student scores vary widely between 25 and 95 across all departments, and average differences between faculties are trivial by comparison.

4.7 Advanced scientific visualizations

g = sns.lmplot(
    data=first_sem, x="sleep_hours", y="gpa", col="programme",
    scatter_kws={"alpha": 0.35, "color": "#4c72b0"},
    line_kws={"color": "#c44e52"},
    height=3.5, aspect=1.1
)
g.set_axis_labels("Hours of Sleep per Night", "GPA")
g.set_titles(col_template="{col_name} Candidates")
plt.tight_layout()
plt.show()