exploratory-data-analysis.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "95f0a171",
   "metadata": {},
   "source": [
    "(exploratory-data-analysis)=\n",
    "# Exploratory Data Analysis\n",
    "\n",
    "## Introduction\n",
    "\n",
    "This chapter will show you how to use visualisation and transformation to explore your data in a systematic way, a task that data scientists call exploratory data analysis, or EDA for short. EDA is an iterative cycle; you:\n",
    "\n",
    "1.  Generate questions about your data.\n",
    "\n",
    "2.  Search for answers by visualising, transforming, and modelling your data.\n",
    "\n",
    "3.  Use what you learn to refine your questions and/or generate new questions.\n",
    "\n",
    "EDA is not a formal process with a strict set of rules and, during the initial phases of EDA, you should feel free to investigate every idea that occurs to you. Some of these ideas will pan out, and some will be dead ends. As your exploration continues, you will home in on a few particularly productive areas that you'll eventually write up and communicate to others. As you explore your data, you should remember that there are some pitfalls: you should always think about how the data were collected, what might be missing, whether there are quality problems, and be really strict about the differences between correlation and causation (this is a huge topic in itself!).\n",
    "\n",
    "### Prerequisites\n",
    "\n",
    "For doing EDA, we'll use the **pandas**, **skimpy**, and **pandas-profiling** packages. We'll also need **lets-plot** for data visualisation. All of these can be installed via `pip install <packagename>`.\n",
    "\n",
    "As ever, we begin by loading these packages that we'll use:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a3377aa6",
   "metadata": {},
   "outputs": [],
   "source": [
    "from skimpy import skim\n",
    "import pandas as pd\n",
    "from pandas.api.types import CategoricalDtype\n",
    "from lets_plot import *\n",
    "from lets_plot.mapping import as_discrete\n",
    "\n",
    "LetsPlot.setup_html()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "51a55374",
   "metadata": {
    "tags": [
     "remove-cell"
    ]
   },
   "outputs": [],
   "source": [
    "import matplotlib_inline.backend_inline\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# Plot settings\n",
    "plt.style.use(\"https://github.com/aeturrell/python4DS/raw/main/plot_style.txt\")\n",
    "matplotlib_inline.backend_inline.set_matplotlib_formats(\"svg\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e4ddb863",
   "metadata": {},
   "source": [
    "## Questions\n",
    "\n",
    "> \"There are no routine statistical questions, only questionable statistical routines.\" --- Sir David Cox\n",
    "\n",
    "> \"Far better an approximate answer to the right question, which is often vague, than an exact answer to the wrong question, which can always be made precise.\" --- John Tukey\n",
    "\n",
    "Your goal during EDA is to develop an understanding of your data. The easiest way to do this is to use questions as tools to guide your investigation. When you ask a question, the question focuses your attention on a specific part of your dataset and helps you decide which graphs, models, or transformations to make.\n",
    "\n",
    "EDA is fundamentally a creative process. And like most creative processes, the key to asking *quality* questions is to generate a large *quantity* of questions. It is difficult to ask revealing questions at the start of your analysis because you do not know what insights are contained in your dataset. On the other hand, each new question that you ask will expose you to a new aspect of your data and increase your chance of making a discovery. You can quickly drill down into the most interesting parts of your data---and develop a set of thought-provoking questions---if you follow up each question with a new question based on what you find.\n",
    "\n",
    "There is no rule about which questions you should ask to guide your research.\n",
    "However, two types of questions will always be useful for making discoveries within your data.\n",
    "You can loosely word these questions as:\n",
    "\n",
    "1.  What type of variation occurs within my variables?\n",
    "\n",
    "2.  What type of covariation occurs between my variables?\n",
    "\n",
    "The rest of this chapter will look at these two questions.\n",
    "We'll explain what variation and covariation are, and we'll show you several ways to answer each question."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6adeff41",
   "metadata": {},
   "source": [
    "## Variation\n",
    "\n",
    "**Variation** is the tendency of the values of a variable to change from measurement to measurement. You can see variation easily in real life; if you measure any continuous variable twice, you will get two different results. This is true even if you measure quantities that are constant, like the speed of light, because of deficiencies in equipment. Each of your measurements includes a small amount of error that varies from measurement to measurement. Variables can also vary if you measure across different subjects (e.g. the eye colours of different people) or different times (e.g. the energy levels of an electron at different moments). The best way to understand that pattern is to visualize the distribution of the variable's values.\n",
    "\n",
    "We'll start our exploration by visualizing the distribution of weights (`\"carat\"`) of \\~54,000 diamonds from the `diamonds` dataset."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "069caa7c",
   "metadata": {},
   "outputs": [],
   "source": [
    "diamonds = pd.read_csv(\n",
    "    \"https://github.com/mwaskom/seaborn-data/raw/master/diamonds.csv\"\n",
    ")\n",
    "diamonds[\"cut\"] = diamonds[\"cut\"].astype(\n",
    "    CategoricalDtype(\n",
    "        categories=[\"Fair\", \"Good\", \"Very Good\", \"Premium\", \"Ideal\"], ordered=True\n",
    "    )\n",
    ")\n",
    "diamonds[\"color\"] = diamonds[\"color\"].astype(\n",
    "    CategoricalDtype(categories=[\"D\", \"E\", \"F\", \"G\", \"H\", \"I\", \"J\"], ordered=True)\n",
    ")\n",
    "diamonds.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "01f5979c",
   "metadata": {},
   "source": [
    "Since `\"carat\"` is a numerical variable, we can use a histogram:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "97900f58",
   "metadata": {},
   "outputs": [],
   "source": [
    "(ggplot(diamonds, aes(x=\"carat\")) + geom_histogram(binwidth=0.5))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2307ba7c",
   "metadata": {},
   "source": [
    "Now that you can visualize variation, what should you look for in your plots?\n",
    "And what type of follow-up questions should you ask?\n",
    "We've put together a list below of the most useful types of information that you will find in your graphs, along with some follow-up questions for each type of information.\n",
    "The key to asking good follow-up questions will be to rely on your curiosity (What do you want to learn more about?) as well as your skepticism (How could this be misleading?).\n",
    "\n",
    "### Typical values\n",
    "\n",
    "In both bar charts and histograms, tall bars show the common values of a variable, and shorter bars show less-common values.\n",
    "Places that do not have bars reveal values that were not seen in your data.\n",
    "To turn this information into useful questions, look for anything unexpected:\n",
    "\n",
    "-   Which values are the most common?\n",
    "    Why?\n",
    "\n",
    "-   Which values are rare?\n",
    "    Why?\n",
    "    Does that match your expectations?\n",
    "\n",
    "-   Can you see any unusual patterns?\n",
    "    What might explain them?\n",
    "\n",
    "Let's take a look at the distribution of `\"carat\"` for smaller diamonds.\n",
    "\n",
    "Note that when we create `smaller_diamonds`, we're doing it by creating a copy. Otherwise, any changes we made to `smaller_diamonds` would also affect `diamonds` (the two point to the same underlying data in your computer's memory). Sometimes you may want to have a cut still be connected to the original dataset, and sometimes you don't; in this case, we'd like them to be distinct so we use `copy()`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "20d75550",
   "metadata": {},
   "outputs": [],
   "source": [
    "smaller_diamonds = diamonds.query(\"carat < 3\").copy()\n",
    "\n",
    "(ggplot(smaller_diamonds, aes(x=\"carat\")) + geom_histogram(binwidth=0.01))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ba20a0e2",
   "metadata": {},
   "source": [
    "This histogram suggests several interesting questions:\n",
    "\n",
    "-   Why are there more diamonds at whole carats and common fractions of carats?\n",
    "\n",
    "-   Why are there more diamonds slightly to the right of each peak than there are slightly to the left of each peak?\n",
    "\n",
    "Visualisations can also reveal clusters, which suggest that subgroups exist in your data.\n",
    "To understand the subgroups, ask:\n",
    "\n",
    "-   How are the observations within each subgroup similar to each other?\n",
    "\n",
    "-   How are the observations in separate clusters different from each other?\n",
    "\n",
    "-   How can you explain or describe the clusters?\n",
    "\n",
    "-   Why might the appearance of clusters be misleading?\n",
    "\n",
    "Some of these questions can be answered with the data while some will require domain expertise about the data.\n",
    "Many of them will prompt you to explore a relationship *between* variables, for example, to see if the values of one variable can explain the behavior of another variable.\n",
    "We'll get to that shortly."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0626d35a",
   "metadata": {},
   "source": [
    "### Unusual values\n",
    "\n",
    "Outliers are observations that are unusual; data points that don't seem to fit the pattern.\n",
    "Sometimes outliers are data entry errors, sometimes they are simply values at the extremes that happened to be observed in this data collection, and other times they suggest important new discoveries.\n",
    "When you have a lot of data, outliers are sometimes difficult to see in a histogram.\n",
    "For example, take the distribution of the `\"y\"` variable from the diamonds dataset.\n",
    "The only evidence of outliers is the unusually wide limits on the x-axis.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d9d7e995",
   "metadata": {},
   "outputs": [],
   "source": [
    "(ggplot(diamonds, aes(x=\"y\")) + geom_histogram(binwidth=0.5))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "05bcf733",
   "metadata": {},
   "source": [
    "There are so many observations in the common bins that the rare bins are very short, making it very difficult to see them (although maybe if you stare intently at 0 you'll spot something).\n",
    "To make it easy to see the unusual values, we need to zoom to small values of the y-axis with `coord_cartesian()`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ea8f8bf3",
   "metadata": {},
   "outputs": [],
   "source": [
    "(\n",
    "    ggplot(diamonds, aes(x=\"y\"))\n",
    "    + geom_histogram(binwidth=0.5)\n",
    "    + coord_cartesian(ylim=[0, 50])\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ba2b2c79",
   "metadata": {},
   "source": [
    "`coord_cartesian()` also has an `xlim()` argument for when you need to zoom into the x-axis.\n",
    "**Lets-Plot** also has `xlim()` and `ylim()` functions that work slightly differently: they throw away the data outside the limits.\n",
    "\n",
    "This allows us to see that there are three unusual values: 0, \\~30, and \\~60. We pluck them out with **pandas**:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e81ffb55",
   "metadata": {},
   "outputs": [],
   "source": [
    "unusual = diamonds.query(\"y < 3 or y > 20\").loc[:, [\"x\", \"y\", \"z\", \"price\"]]\n",
    "unusual"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c9321f36",
   "metadata": {},
   "source": [
    "The `\"y\"` variable measures one of the three dimensions of these diamonds, in mm.\n",
    "We know that diamonds can't have a width of 0mm, so these values must be incorrect.\n",
    "By doing EDA, we have discovered missing data that was coded as 0, which we never would have found by simply searching for `NA`s.\n",
    "Going forward we might choose to re-code these values as `NA`s in order to prevent misleading calculations.\n",
    "We might also suspect that measurements of 32mm and 59mm are implausible: those diamonds are over an inch long, but don't cost hundreds of thousands of dollars!\n",
    "\n",
    "It's good practice to repeat your analysis with and without the outliers.\n",
    "If they have minimal effect on the results, and you can't figure out why they're there, it's reasonable to omit them, and move on.\n",
    "However, if they have a substantial effect on your results, you shouldn't drop them without justification.\n",
    "You'll need to figure out what caused them (e.g., a data entry error) and disclose that you removed them in your write-up.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "142d21d7",
   "metadata": {},
   "source": [
    "### Exercises\n",
    "\n",
    "1.  Explore the distribution of each of the `x`, `y`, and `z` variables in `diamonds`.\n",
    "    What do you learn?\n",
    "    Think about a diamond and how you might decide which dimension is the length, width, and depth.\n",
    "\n",
    "2.  Explore the distribution of `\"price\"`.\n",
    "    Do you discover anything unusual or surprising?\n",
    "    (Hint: Carefully think about the `binwidth=` keyword argument setting and make sure you try a wide range of values.)\n",
    "\n",
    "3.  How many diamonds are 0.99 carat?\n",
    "    How many are 1 carat?\n",
    "    What do you think is the cause of the difference?\n",
    "\n",
    "4.  Compare and contrast `coord_cartesian()` vs. `xlim()` or `ylim()` when zooming in on a histogram.\n",
    "    What happens if you leave binwidth unset?\n",
    "    What happens if you try and zoom so only half a bar shows?\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c51e8d6f",
   "metadata": {},
   "source": [
    "## Unusual Values\n",
    "\n",
    "If you've encountered unusual values in your dataset, and simply want to move on to the rest of your analysis, you have two options.\n",
    "\n",
    "1.  Drop the entire row with the strange values:\n",
    "\n",
    "    ```python\n",
    "    condition = ((diamonds[\"y\"] < 3) | (diamonds[\"y\"] > 20))\n",
    "    diamonds2 = diamonds.loc[~condition, :]\n",
    "    ```\n",
    "\n",
    "    We don't recommend this option because one invalid value doesn't imply that all the other values for that observation are also invalid.\n",
    "    Additionally, if you have low quality data, by the time that you've applied this approach to every variable you might find that you don't have any data left!\n",
    "\n",
    "2.  Instead, we recommend replacing the unusual values with missing values.\n",
    "    One way to do this, which makes a distinction between data frames that have had the unusual values replaced and the original data, is to make a copy and then set the problematic values to `pd.NA`, **pandas**'s special NA value.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ecf345a7",
   "metadata": {},
   "outputs": [],
   "source": [
    "diamonds2 = diamonds.copy()\n",
    "condition = (diamonds2[\"y\"] < 3) | (diamonds2[\"y\"] > 20)\n",
    "diamonds2.loc[condition, \"y\"] = pd.NA"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d26d922e",
   "metadata": {},
   "source": [
    "It's not obvious where you should plot missing values, so **lets-plot** doesn't include them in the plot:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "15a43255",
   "metadata": {},
   "outputs": [],
   "source": [
    "(ggplot(diamonds2, aes(x=\"x\", y=\"y\")) + geom_point())"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e8ae854b",
   "metadata": {},
   "source": [
    "Other times you want to understand what makes observations with missing values different to observations with recorded values.\n",
    "For example, in the nycflights13 data, missing values in the `\"dep_time\"` variable indicate that the flight was cancelled.\n",
    "So you might want to compare the scheduled departure times for cancelled and non-cancelled times.\n",
    "You can do this by making a new variable, using `is.na()` to check if `\"dep_time\"` is missing."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0a4ea922",
   "metadata": {},
   "outputs": [],
   "source": [
    "url = \"https://raw.githubusercontent.com/byuidatascience/data4python4ds/master/data-raw/flights/flights.csv\"\n",
    "flights = pd.read_csv(url)\n",
    "flights.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6849f4d9",
   "metadata": {},
   "outputs": [],
   "source": [
    "flights2 = flights.assign(\n",
    "    cancelled=lambda x: pd.isna(x[\"dep_time\"]),\n",
    "    sched_hour=lambda x: x[\"sched_dep_time\"] // 100,\n",
    "    sched_min=lambda x: x[\"sched_dep_time\"] % 100,\n",
    "    sched_dep_time=lambda x: x[\"sched_hour\"] + x[\"sched_min\"] / 60,\n",
    ")\n",
    "\n",
    "(\n",
    "    ggplot(flights2, aes(x=\"sched_dep_time\"))\n",
    "    + geom_freqpoly(aes(color=\"cancelled\"), binwidth=1 / 4)\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b97e453a",
   "metadata": {},
   "source": [
    "However this plot isn't great because there are many more non-cancelled flights than cancelled flights.\n",
    "In the next section we'll explore some techniques for improving this comparison.\n",
    "\n",
    "### Exercises\n",
    "\n",
    "1.  What happens to missing values in a histogram?\n",
    "    What happens to missing values in a bar chart?\n",
    "    Why is there a difference in how missing values are handled in histograms and bar charts?\n",
    "\n",
    "2.  Recreate the frequency plot of `scheduled_dep_time` coloured by whether the flight was cancelled or not.\n",
    "    Also facet by the `cancelled` variable.\n",
    "    Experiment with different values of the `scales` variable in the faceting function to mitigate the effect of more non-cancelled flights than cancelled flights.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f63d0b9f",
   "metadata": {},
   "source": [
    "## Covariation\n",
    "\n",
    "If variation describes the behavior *within* a variable, covariation describes the behavior *between* variables.\n",
    "**Covariation** is the tendency for the values of two or more variables to vary together in a related way.\n",
    "The best way to spot covariation is to visualise the relationship between two or more variables, but note that covariation doesn't imply a causal relationship between variables.\n",
    "\n",
    "### A categorical and a numerical variable\n",
    "\n",
    "For example, let's explore how the price of a diamond varies with its quality (measured by `\"cut\"`) using `geom_freqpoly()`:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e1719d8f",
   "metadata": {},
   "outputs": [],
   "source": [
    "(\n",
    "    ggplot(diamonds, aes(x=\"price\"))\n",
    "    + geom_freqpoly(aes(color=\"cut\"), binwidth=500, linewidth=0.75)\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "20387b96",
   "metadata": {},
   "source": [
    "The default appearance of `geom_freqpoly()` is not that useful here because the height, determined by the overall count, differs so much across cuts, making it hard to see the differences in the shapes of their distributions.\n",
    "\n",
    "To make the comparison easier we need to swap what is displayed on the y-axis.\n",
    "Instead of displaying count, we'll display the **density**, which is the count standardised so that the area under each frequency polygon is one."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9388e24b",
   "metadata": {},
   "outputs": [],
   "source": [
    "(\n",
    "    ggplot(diamonds, aes(x=\"price\"))\n",
    "    + geom_density(aes(color=\"cut\", fill=\"cut\"), size=1, alpha=0.2)\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "157f63dd",
   "metadata": {},
   "source": [
    "There’s something rather surprising about this plot - it appears that fair diamonds (the lowest quality) have the highest average price! But maybe that’s because density plots are a little hard to interpret - there’s a lot going on in this plot.\n",
    "\n",
    "A visually simpler plot for exploring this relationship is using side-by-side boxplots."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a3f333a6",
   "metadata": {},
   "outputs": [],
   "source": [
    "(ggplot(diamonds, aes(x=\"cut\", y=\"price\")) + geom_boxplot())"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b9d9ef00",
   "metadata": {},
   "source": [
    "We see much less information about the distribution, but the boxplots are much more compact so we can more easily compare them (and fit more on one plot). It supports the counter-intuitive finding that better quality diamonds are typically cheaper! In the exercises, you’ll be challenged to figure out why.\n",
    "\n",
    "`\"cut\"` is an ordered categorical variable: fair is worse than good, which is worse than very good and so on. Many categorical variables don’t have such an intrinsic order, so you might want to reorder them to make a more informative display. One way to do that is according to the median value, though other options are available.\n",
    "\n",
    "As an example, with the mpg dataset, we might want to look at how highway mileage varies across classes:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6949db81",
   "metadata": {},
   "outputs": [],
   "source": [
    "mpg = pd.read_csv(\n",
    "    \"https://vincentarelbundock.github.io/Rdatasets/csv/ggplot2/mpg.csv\", index_col=0\n",
    ")\n",
    "mpg[\"class\"] = mpg[\"class\"].astype(\"category\")\n",
    "\n",
    "(ggplot(mpg, aes(x=\"class\", y=\"hwy\")) + geom_boxplot())"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "871aaf1c",
   "metadata": {},
   "source": [
    "To make the trend easier to see, we can reorder class based on the median value of `\"hwy\"`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a5b1ed09",
   "metadata": {},
   "outputs": [],
   "source": [
    "(ggplot(mpg) + geom_boxplot(aes(as_discrete(\"class\", order_by=\"..middle..\"), \"hwy\")))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "dde59236",
   "metadata": {},
   "source": [
    "If you have long variable names, geom_boxplot() will work better if you flip it 90°. You can do that by adding `coord_flip()`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "920a4268",
   "metadata": {},
   "outputs": [],
   "source": [
    "(\n",
    "    ggplot(mpg)\n",
    "    + geom_boxplot(aes(as_discrete(\"class\", order_by=\"..middle..\"), \"hwy\"))\n",
    "    + coord_flip()\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "299635df",
   "metadata": {},
   "source": [
    "#### Exercises\n",
    "\n",
    "1.  Use what you've learned to improve the visualisation of the departure times of cancelled vs. non-cancelled flights.\n",
    "\n",
    "2.  Based on EDA, what variable in the diamonds dataset appears to be most important for predicting the price of a diamond?\n",
    "    How is that variable correlated with cut?\n",
    "    Why does the combination of those two relationships lead to lower quality diamonds being more expensive?\n",
    "\n",
    "3.  Create a visualisation of diamond prices vs. a categorical variable from the `diamonds` dataset using `geom_violin()`, then a faceted `geom_histogram()`, then a coloured `geom_freqpoly()`, and then a coloured `geom_density()`.\n",
    "    Compare and contrast the four plots.\n",
    "    What are the pros and cons of each method of visualising the distribution of a numerical variable based on the levels of a categorical variable?\n",
    "\n",
    "4.  If you have a small dataset, it's sometimes useful to use `geom_jitter()` to avoid overplotting to more easily see the relationship between a continuous and categorical variable.\n",
    "    The ggbeeswarm package provides a number of methods similar to `geom_jitter()`.\n",
    "    List them and briefly describe what each one does."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6e6a0b46",
   "metadata": {},
   "source": [
    "### Two categorical variables\n",
    "\n",
    "To visualise the covariation between categorical variables, you'll need to count the number of observations for each combination of levels of these categorical variables. You can do this with a `pd.crosstab()` that we then melt to put it in \"tidy\" format."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "68d330d2",
   "metadata": {},
   "outputs": [],
   "source": [
    "ct_cut_color = pd.melt(\n",
    "    pd.crosstab(diamonds[\"cut\"], diamonds[\"color\"]).reset_index(),\n",
    "    id_vars=[\"cut\"],\n",
    "    value_vars=diamonds[\"color\"].unique(),\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c2ebc0c6",
   "metadata": {},
   "source": [
    "Followed by visualising it with `geom_tile()`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e858cd22",
   "metadata": {},
   "outputs": [],
   "source": [
    "(ggplot(ct_cut_color, aes(x=\"color\", y=\"cut\")) + geom_tile(aes(fill=\"value\")))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8ccd5c1b",
   "metadata": {},
   "source": [
    "### Exercises\n",
    "\n",
    "1. How could you rescale the count dataset above to more clearly show the distribution of cut within color, or colour within cut?\n",
    "\n",
    "2. What different data insights do you get with a segmented bar chart if colour is mapped to the `x` aesthetic and cut is mapped to the fill aesthetic? Calculate the counts that fall into each of the segments.\n",
    "\n",
    "3. Use `geom_tile()` together with **pandas** to explore how average flight departure delays vary by destination and month of year. What makes the plot difficult to read? How could you improve it?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b2f19afc",
   "metadata": {},
   "source": [
    "### Two numerical variables\n",
    "\n",
    "You've already seen one great way to visualise the covariation between two numerical variables: draw a scatterplot with `geom_point()`.\n",
    "You can see covariation as a pattern in the points.\n",
    "For example, you can see a positive association between the carat size and price of a diamond: diamonds with more carats have a higher price.\n",
    "The relationship is exponential."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2afe2535",
   "metadata": {},
   "outputs": [],
   "source": [
    "(ggplot(smaller_diamonds, aes(x=\"carat\", y=\"price\")) + geom_point())"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9d8034f2",
   "metadata": {},
   "source": [
    "(In this section we'll use the `smaller_diamonds` dataset to stay focused on the bulk of the diamonds that are smaller than 3 carats)\n",
    "\n",
    "Scatterplots become less useful as the size of your dataset grows, because points begin to overplot, and pile up into areas of uniform black, making it hard to judge differences in the density of the data across the 2-dimensional space as well as making it hard to spot the trend.\n",
    "You've already seen one way to fix the problem: using the `alpha` aesthetic to add transparency.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b55707a9",
   "metadata": {},
   "outputs": [],
   "source": [
    "(ggplot(smaller_diamonds, aes(x=\"carat\", y=\"price\")) + geom_point(alpha=1 / 20))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "351c22bd",
   "metadata": {},
   "source": [
    "But using transparency can be challenging for very large datasets. In that case, we recommend a *binscatter*, or binned scatterplot. A binned scatterplot divides the conditioning variable, `\"carat\"` in our example, into equally sized bins or quantiles, and then plots the conditional mean of the dependent variable, `\"price\"` in our example, within each bin. Bin scatters often come with confidence intervals too. A good bin scatter package in Python is [**binsreg**](https://nppackages.github.io/binsreg/). However, bin scatters are an advanced topic, and we won't cover them here."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "14158b48",
   "metadata": {},
   "source": [
    "## **pandas** built-in tools for EDA\n",
    "\n",
    "**pandas** has some great options for built-in EDA; in fact we've already seen one of them, `df.info()` which, as well as reporting datatypes and memory usage, also tells us how many observations in each column are 'truthy' rather than 'falsy', ie how many have non-null values.\n",
    "\n",
    "### Exploratory tables and descriptive statistics\n",
    "\n",
    "A small step beyond `.info()` to get tables is to use `.describe()` which, if you have mixed datatypes that include floats, will report some basic summary statistics:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "13079065",
   "metadata": {},
   "outputs": [],
   "source": [
    "diamonds.describe()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "db57dbb7",
   "metadata": {},
   "source": [
    "Although helpful, that sure is hard to read! We can improve this by using the `round()` method too:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b4144440",
   "metadata": {},
   "outputs": [],
   "source": [
    "sum_table = diamonds.describe().round(1)\n",
    "sum_table"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8063b8f5",
   "metadata": {},
   "source": [
    "Published summary statistics tables often list one variable per row, and if your data frame has many variables, `describe()` can quickly get too wide to read easily. You can transpose it using the `T` property (or the `transpose()` method):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cd2f8772",
   "metadata": {},
   "outputs": [],
   "source": [
    "sum_table = sum_table.T\n",
    "sum_table"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e67104d3",
   "metadata": {},
   "source": [
    "Of course, the stats provided in this pre-built table are not very customised. So what do we do to get the table that we actually want? Well, the answer is to draw on the contents of the previous data chapters, particularly the introduction to data analysis. Groupbys, merges, aggregations: use all of them to produce the EDA table that you want.\n",
    "\n",
    "If you're exploring data, you might also want to be able to read everything clearly and see any deviations from what you'd expect quickly. **pandas** has some built-in functionality that styles data frames to help you. These styles persist when you export the data frame to, say, Excel, too.\n",
    "\n",
    "Here's an example that highlights some ways of styling data frames, making use of several features such as: unstacking into a wider format (`unstack`), changing the units (`lambda` function; note that `1e3` is shorthand for `1000` on computers), fill NaNs with unobtrusive strings (`.fillna('-')`), removing numbers after the decimal place (`.style.format(precision=0)`), and adding a caption (`.style.set_caption`)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5afcacbc",
   "metadata": {},
   "outputs": [],
   "source": [
    "(\n",
    "    diamonds.groupby([\"cut\", \"color\"])[\"price\"]\n",
    "    .mean()\n",
    "    .unstack()\n",
    "    .apply(lambda x: x / 1e3)\n",
    "    .fillna(\"-\")\n",
    "    .style.format(precision=2)\n",
    "    .set_caption(\"Sale price (thousands)\")\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d16e0fc1",
   "metadata": {},
   "source": [
    "Although a neater one than we've seen, this is still a drab table of numbers. The eye is not immediately drawn to it!\n",
    "\n",
    "To remedy that, let's take a look at another styling technique: the use of colour. Let's say we wanted to make a table that showed a cross-tabulation between cut and color; that is the counts of objects appearing in both of these fields according to the categories.\n",
    "\n",
    "To perform a cross-tabulation, we'll use the built-in `pd.crosstab()` but we'll ask that the values that appear in the table (counts) be lit up with a heatmap using `style.background_gradient()` too:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "21e65189",
   "metadata": {},
   "outputs": [],
   "source": [
    "pd.crosstab(diamonds[\"color\"], diamonds[\"cut\"]).style.background_gradient(cmap=\"plasma\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b7c6aa71",
   "metadata": {},
   "source": [
    "By default, `background_gradient()` highlights each number relative to the others in its column; you can highlight by row using `axis=1` or relative to all table values using `axis=0`. And of course `plasma` is just one of [many available colormaps](https://matplotlib.org/stable/tutorials/colors/colormaps.html)!\n",
    "\n",
    "```{admonition} Exercise\n",
    "Do a new cross-tabulation using a different colourmap.\n",
    "```\n",
    "\n",
    "Here are a couple of other styling tips for data frames.\n",
    "\n",
    "First, use bars to show ordering:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "bb0162ba",
   "metadata": {},
   "outputs": [],
   "source": [
    "(\n",
    "    pd.crosstab(diamonds[\"color\"], diamonds[\"cut\"])\n",
    "    .style.format(precision=0)\n",
    "    .bar(color=\"#d65f5f\")\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e36d4795",
   "metadata": {},
   "source": [
    "Use `.hightlight_max()`, and similar commands, to show important entries:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5d19072c",
   "metadata": {},
   "outputs": [],
   "source": [
    "pd.crosstab(diamonds[\"color\"], diamonds[\"cut\"]).style.highlight_max().format(\"{:.0f}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "63075da4",
   "metadata": {},
   "source": [
    "You can find a full set of styling commands [here](https://pandas.pydata.org/pandas-docs/stable/user_guide/style.html#Styling)."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d8b6d225",
   "metadata": {},
   "source": [
    "### Exploratory Plotting with **pandas**\n",
    "\n",
    "**pandas** has some built-in plotting options to help you look at data quickly. These can be accessed via `.plot.*` or `.plot()`, depending on the context. Let's make a quick `.plot()` using a dataset on taxis."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b479d5b1",
   "metadata": {},
   "outputs": [],
   "source": [
    "taxis = pd.read_csv(\"https://github.com/mwaskom/seaborn-data/raw/master/taxis.csv\")\n",
    "# turn the pickup time column into a datetime\n",
    "taxis[\"pickup\"] = pd.to_datetime(taxis[\"pickup\"])\n",
    "# set some other columns types\n",
    "taxis = taxis.astype(\n",
    "    {\n",
    "        \"dropoff\": \"datetime64[ns]\",\n",
    "        \"pickup\": \"datetime64[ns]\",\n",
    "        \"color\": \"category\",\n",
    "        \"payment\": \"category\",\n",
    "        \"pickup_zone\": \"string\",\n",
    "        \"dropoff_zone\": \"string\",\n",
    "        \"pickup_borough\": \"category\",\n",
    "        \"dropoff_borough\": \"category\",\n",
    "    }\n",
    ")\n",
    "taxis.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ee971c9c",
   "metadata": {},
   "outputs": [],
   "source": [
    "taxis.info()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2015b1dc",
   "metadata": {},
   "outputs": [],
   "source": [
    "(\n",
    "    taxis.set_index(\"pickup\")\n",
    "    .groupby(pd.Grouper(freq=\"D\"))[\"total\"]\n",
    "    .mean()\n",
    "    .plot(\n",
    "        title=\"Mean taxi fares\",\n",
    "        xlabel=\"\",\n",
    "        ylabel=\"Fare (USD)\",\n",
    "    )\n",
    ");"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b9c0f990",
   "metadata": {},
   "source": [
    "Again, if you can get the data in the right shape, you can plot it. The same function works with multiple lines\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "51e86185",
   "metadata": {},
   "outputs": [],
   "source": [
    "(\n",
    "    taxis.set_index(\"pickup\")\n",
    "    .groupby(pd.Grouper(freq=\"D\"))[[\"fare\", \"tip\", \"tolls\"]]\n",
    "    .mean()\n",
    "    .plot(\n",
    "        style=[\"-\", \":\", \"-.\"],\n",
    "        title=\"Components of taxi fares\",\n",
    "        xlabel=\"\",\n",
    "        ylabel=\"USD\",\n",
    "    )\n",
    ");"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cc227b81",
   "metadata": {},
   "source": [
    "Now let's see some of the other quick `.plot.*` options."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1c4278ac",
   "metadata": {},
   "source": [
    "A bar chart (use `barh` for horizontal orientation; `rot` sets rotation of labels):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "79ceca92",
   "metadata": {},
   "outputs": [],
   "source": [
    "taxis.value_counts(\"payment\").sort_index().plot.bar(title=\"Counts\", rot=0);"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4e6f1cc2",
   "metadata": {},
   "source": [
    "This next one, uses `.plot.hist()` to create a histogram."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5efc5817",
   "metadata": {},
   "outputs": [],
   "source": [
    "taxis[\"tip\"].plot.hist(bins=30, title=\"Tip\");"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3f15bfa7",
   "metadata": {},
   "source": [
    "Boxplot:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0b735d15",
   "metadata": {},
   "outputs": [],
   "source": [
    "(taxis[[\"fare\", \"tolls\", \"tip\"]].plot.box());"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "67017bef",
   "metadata": {},
   "source": [
    "Scatter plot:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "66adada2",
   "metadata": {},
   "outputs": [],
   "source": [
    "taxis.plot.scatter(x=\"fare\", y=\"tip\", alpha=0.7, ylim=(0, None));"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9766d2c4",
   "metadata": {},
   "source": [
    "## Other tools for EDA\n",
    "\n",
    "Between **pandas** and visualisation packages, you have a lot of what you need for EDA. But there are some tools just dedicated to making EDA easier that it's worth knowing about."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "813ec70a",
   "metadata": {},
   "source": [
    "### **skimpy** for summary statistics\n",
    "\n",
    "The **skimpy** package is a light weight tool that provides summary statistics about variables in data frames in the console (rather than in a big HTML report, which is what the other EDA packages in the rest of this chapter too). Sometimes running `.summary()` on a data frame isn't enough, and **skimpy** fills this gap. It also comes with the `clean_columns()` function for cleaning column names that we saw in an earlier chapter. To install **skimpy**, run `pip install skimpy` in the terminal.\n",
    "\n",
    "Let's see **skimpy** in action."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "32796b5f",
   "metadata": {},
   "outputs": [],
   "source": [
    "skim(taxis)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f2810c9e",
   "metadata": {},
   "source": [
    "## Summary\n",
    "\n",
    "In this chapter, you've learned a variety of tools to help you understand the variation within your data.\n",
    "You've seen techniques that work with a single variable at a time and with a pair of variables.\n",
    "This might seem painfully restrictive if you have tens or hundreds of variables in your data, but they're foundation upon which all other techniques are built.\n",
    "\n",
    "In the next chapter, we'll focus on the tools we can use to communicate our results."
   ]
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "9d7534ecd9fbc7d385378f8400cf4d6cb9c6175408a574f1c99c5269f08771cc"
  },
  "jupytext": {
   "cell_metadata_filter": "-all",
   "encoding": "# -*- coding: utf-8 -*-",
   "formats": "md:myst",
   "main_language": "python"
  },
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.13"
  },
  "toc-showtags": true
 },
 "nbformat": 4,
 "nbformat_minor": 5
}