<jupyter_start><jupyter_text>Exploratory Data Analysis in Action - Data preparationIn this notebook we apply some useful techniques for cleaning and organizing our data set. **Import statements**<jupyter_code>import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
%matplotlib inline<jupyter_output><empty_output><jupyter_text>**Global settings**<jupyter_code>pd.options.display.max_rows = 999
pd.options.display.max_columns = 100
plt.rcParams["figure.figsize"] = [15,6]<jupyter_output><empty_output><jupyter_text>Cleaning and organizing data ![](./_img/Time_data_science.png)Source: [Gil Press (2016)](https://www.forbes.com/sites/gilpress/2016/03/23/data-preparation-most-time-consuming-least-enjoyable-data-science-task-survey-says/55852a146f63) **Load data set**<jupyter_code>PATH = "./data/"
df_raw = pd.read_csv(PATH + "operations.csv", low_memory=False)<jupyter_output><empty_output><jupyter_text>Data cleaning[Data cleansing or data cleaning](https://en.wikipedia.org/wiki/Data_cleansing) is the process of **detecting and correcting (or removing) corrupt or inaccurate records** from a record set, table, or database and refers to **identifying incomplete, incorrect, inaccurate or irrelevant parts of the data and then replacing, modifying, or deleting the dirty or coarse data**.<jupyter_code>df_clean = df_raw.copy()
df_clean.sample(10)<jupyter_output><empty_output><jupyter_text>Open issues: Deal with incomplete, incorrect, inaccurate or irrelevant parts of the data-  Identify    -  Correct-  Spatial subsetting the data set  Dealing with incomplete data (`NaN`)Missing values in data sets are a well-known problem as nearly everywhere, where data is measured and recorded, issues with missing values occur. Various reasons lead to missing values: values may not be measured, values may be measured but get lost or values may be measured but are considered unusable. Missing values can lead to problems, because often further data processing and analysis steps rely on complete data sets. Therefore missing values need to be replaced with reasonable values. In statistics this process is called **imputation**.When faced with the problem of missing values it is important to understand the mechanism that causes missing data. Such an understanding is useful, as it may be employed as background knowledge for selecting an appropriate imputation strategy.  **Check for `NaN` **Note that in many cases missing values are assigned special characters, such as `-999`, `NA`, `k.A.` etc.; hence, you as a data analyst are responsible for taking appropriate action.<jupyter_code>df_clean.shape[0]
df_clean.isnull().sum()
df_clean.notnull().sum()<jupyter_output><empty_output><jupyter_text>**Strategies to deal with missing data in Python**In general there are many options to consider when imputing missing values, for example:* A constant value that has meaning within the domain, such as 0, distinct from all other values.* A value from another randomly selected record.* A mean, median or mode value for the column.* A value estimated by another predictive model.There are some libraries implementing more or less advanced missing value imputation strategies such as * [`statsmodels`](http://www.statsmodels.org/dev/imputation.html) ([Multiple Imputation with Chained Equations (MICE)](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3074241/))* [`fancyimpute`](https://github.com/iskandr/fancyimpute) (matrix completion and imputation algorithms)* [`scikit-learn`](http://scikit-learn.org/stable/modules/generated/sklearn.preprocessing.Imputer.html) (mean, median, most frequent)* [`pandas`](https://pandas.pydata.org/pandas-docs/stable/missing_data.html) ([`fillna`](https://pandas.pydata.org/pandas-docs/stable/generated/pandas.DataFrame.fillna.html), [`interpolate`](https://pandas.pydata.org/pandas-docs/stable/generated/pandas.Series.interpolate.html) methods) **Our working strategy to deal with missing data**_Owing to the fact that the amount of missing values in our data set is considerable high and that probably our domain knowledge with respect to World War II history is limited, we simply remove all features of the data set where more the 50% of the data is missing._  > **Challenge**: Write a function (or a script) which takes in our data set and a threshold value of 0.5 and which returns the column names of features that include more data than the given threshold. We then use the returned column names for subsetting our data set.<jupyter_code>def features_above_threshold(df, threshold=0.5):
    ## your code here ...
    pass
# %load ./src/_solutions/features_above_threshold.py<jupyter_output><empty_output><jupyter_text>**Apply `features_above_threshold` on our data set**<jupyter_code>cols2keep = features_above_threshold(df_clean, threshold=0.5)
cols2keep<jupyter_output><empty_output><jupyter_text>**Reassign subsetted data set to the variable `df_clean`**<jupyter_code>df_clean = df_clean[cols2keep]
df_clean.shape
df_clean.sample(3)
df_clean.isnull().sum()<jupyter_output><empty_output><jupyter_text>Dealing with incorrect data (looking for outliers) * **Check that mission dates are between 1939 and 1945**<jupyter_code>df_clean["Mission Date"].dtype
df_clean.loc[0:15, "Mission Date"]<jupyter_output><empty_output><jupyter_text>> **Challenge**: Reassign the `Mission Date` variable, a `string` object with o a time-aware `datetime` object.    _Hint: pandas comes with a_ `.to_datetime` _method_<jupyter_code>## your code here ...
# %load ./src/_solutions/infere_mission_data.py
df_clean["Mission Date"].head()<jupyter_output><empty_output><jupyter_text>Once we have a `datetime` object we can use it to extract the year, month and day assign those to a unique columns in out data set (`year`, `month` and `day`)<jupyter_code>df_clean["year"] = df_clean["Mission Date"].dt.year
df_clean["month"] = df_clean["Mission Date"].dt.month
df_clean["day"] = df_clean["Mission Date"].dt.day
df_clean[["Mission Date", "year", "month", "day"]].sample(10)<jupyter_output><empty_output><jupyter_text>Finally we ask for the minimum and the maximum value of the `year` column.<jupyter_code>print(df_clean["year"].min())
print(df_clean["year"].max())<jupyter_output><empty_output><jupyter_text>* __Check if the values in the altitude column (`Altitude (Hundreds of Feet)`) is within reasonable bounds__  > __Challenge__: Add a new column to the data set (`Altitude (meters)`) by transforming th values in the `Altitude (Hundreds of Feet)` from feet to meters.   _Hint: You may use the `apply` method._$$1\;\text{foot} =  0.3048\; \text{meter}$$<jupyter_code>## your code here ...
df_clean["Altitude (meters)"] = None 
# %load ./src/_solutions/altitude_in_meters.py
df_clean["Altitude (meters)"].describe()
df_clean["Altitude (meters)"].max()<jupyter_output><empty_output><jupyter_text>Considering the output from above and answers to the question "_At what altitude could WW2 planes fly?_" on [Quora](https://www.quora.com/At-what-altitude-could-WW2-planes-fly-What-was-the-record-altitude-during-the-war-and-what-planes-could-fly-at-a-high-altitude-Could-they-fly-at-the-same-altitude-as-these-days-or-not), indicating that flight heights above 40,000 feet (approx. 12,000 m) occurred rarely.> __Challenge__: Limit the reported altitudes to 15,000 m and replace  outliers with `np.nan`.<jupyter_code>max_height = 15000
print((df_clean["Altitude (meters)"] > max_height).sum())
## your code here ...
# %load ./src/_solutions/limit_altitude.py
# check replacement
print((df_clean["Altitude (meters)"] > max_height).sum())

df_clean["Altitude (meters)"].describe()<jupyter_output><empty_output><jupyter_text>* __Check if the values in the target location columns (`Target Latitude` and `Target Longitude`) is within reasonable bounds__.Note that the location is given in geographic coordinates, hence, the values should be within -90 to 90 degree latitude and -180 to 180 degree longitude.<jupyter_code>df_clean[['Target Latitude', 'Target Longitude']].describe()
df_clean.plot.scatter( x='Target Longitude', y='Target Latitude');<jupyter_output><empty_output><jupyter_text>Considering the table and the plot above there are some deviations from the expected values. > __Challenge__: Replace geographical coordinates outside the natural limits [-90, 90] degrees latitude and [-180, 180] degrees longitude with `np.nan`.<jupyter_code>## your code here ...
# %load ./src/_solutions/limit_coordinates.py
# check result
df_clean.plot.scatter( x='Target Longitude', y='Target Latitude');<jupyter_output><empty_output><jupyter_text>Dealing with irrelevant dataOur data set, as well as other data sets we encounter, sometimes come with features which are irrelevant for the purpose our research question. Hence, to reduce the data set size for the purpose of better readability as well as memory issues, among other, we may drop data columns. However, once again domain knowledge helps to decide which columns are of interest or not.<jupyter_code>df_clean.columns
df_clean.head()<jupyter_output><empty_output><jupyter_text>Exploring the data set above we decide to drop the columns `Theater of Operations`, `Air Force`, `Target ID`, `Target Priority`, `Altitude (Hundreds of Feet)`, `Source ID`.<jupyter_code>cols2exlude = ['Theater of Operations', 'Air Force', 'Target ID', 
               'Target Priority', 'Altitude (Hundreds of Feet)', 'Source ID']<jupyter_output><empty_output><jupyter_text>One way to achieve this task with Python is to use the powerful [list comprehension](https://docs.python.org/3/tutorial/datastructures.htmllist-comprehensions).      df_clean = df_clean[[c for c in df_clean.columns if c not in cols2exlude]]  Another way is to leverage the pandas functionality and use the `drop` method. Note that for the purpose of memory and computation efficiency in many cases pandas returns a view of the object, rather than a copy. Hence, if to make a permanent change we have to assign/reassign the object to a variable:    df_clean = df_clean.drop(cols2exlude, axis=1)or use the `inplace=True` argument:    df_clean.drop(cols2exlude, axis=1, inplace=True)<jupyter_code>df_clean.drop(cols2exlude, axis=1, inplace=True)
# check operation
df_clean.columns<jupyter_output><empty_output><jupyter_text>Spatial subsetting the data set_Note: In the subsequent cells we load Python library for spatial data analysis, such as `shapely`, `fiona`,`geopandas`, `cartopy` and `folium`. Make sure that you have installed the [GDAL bindings](http://www.gdal.org/index.html) on your computer._ > __Challenge__: Install the modules `geopandas`, `cartopy` and `folium` into your environment.<jupyter_code>df_subset = df_clean.copy()<jupyter_output><empty_output><jupyter_text>Some of you may have already noticed that the data set contains data of aerial bombing during World War II for all araound the world. In this tutorial however, we want to focus on Europe. Hence, we are going to subset the data set accordingly.Owing to the manifold of spatial data representation, working with spatial data becomes sometimes more involved. In this section we make use third party libraries, such as [GeoPandas](http://geopandas.org/index.html) and [shapely](http://toblerity.org/shapely/), which abstract away many algorithmic or computational issues related to spatial data processing and plotting by integrating the workhorses of geospatial computing, such as [GEOS](http://trac.osgeo.org/geos/), [GDAL](http://www.gdal.org/), [OGR](http://gdal.org/1.11/ogr/) and [proj.4](http://proj4.org/), among others. **Transform the variables `Target Latitude` and `Target Longitude` to spatial coordinates**<jupyter_code>df_subset[['Target Latitude', 'Target Longitude']].head() 
from shapely.geometry import Point
geometry = [Point(xy) for xy in zip(df_subset['Target Longitude'], df_subset['Target Latitude'])]
geometry[0:5]<jupyter_output><empty_output><jupyter_text>**Use the GeoPandas to make a pandas `DataFrame` spatially aware.**[GeoPandas](http://geopandas.org/index.html) extends the datatypes used by pandas to allow spatial operations on geometric types. Geometric operations are performed by [shapely](http://toblerity.org/shapely/). GeoPandas further depends on [fiona](http://toblerity.org/fiona/README.html) for file access and descartes and [matplotlib](https://matplotlib.org/) for plotting.It combines the capabilities of pandas and shapely, providing geospatial operations in pandas and a high-level interface to multiple geometries to shapely.<jupyter_code>import geopandas as gpd
gdf = gpd.GeoDataFrame(df_subset, geometry=geometry)
gdf.head()<jupyter_output><empty_output><jupyter_text>**Make sure that for every entry we have a valid spatial coordinates**<jupyter_code>print(gdf.shape)
# subset only vaild spatial coordinates
gdf = gdf.loc[gdf[['Target Longitude', 'Target Latitude']].notnull().all(axis = 1)]
print(gdf.shape)
gdf[['Target Longitude', 'Target Latitude']].isnull().sum()<jupyter_output><empty_output><jupyter_text>**Assign a spatial coordinate reference system (`crs`) to our GeoPandas object**In general the CRS may be defined in several ways, for example the CRS may be defined as [Well-known text (WKT)](https://en.wikipedia.org/wiki/Well-known_text) format, or [JSON](https://en.wikipedia.org/wiki/JSON) format, or [GML](https://en.wikipedia.org/wiki/Geography_Markup_Language) format, or in the [Proj4](https://en.wikipedia.org/wiki/PROJ.4) format, among many others.The Proj4 format is a generic, string-based description of a CRS. It defines projection types and parameter values for particular projections. For instance the Proj4 format string for the [European Terrestrial Reference System 1989 (ETRS89)](https://en.wikipedia.org/wiki/European_Terrestrial_Reference_System_1989) is:    +proj=longlat +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +no\_defsWith respect to the enormous amount of existing CRS the [International Association of Oil & Gas Producers (IOGP)](https://en.wikipedia.org/wiki/International_Association_of_Oil_%26_Gas_Producers), formerly known as **_European Petroleum Survey Group (EPSG)_**, built a collection of definitions for global, regional, national and local coordinate reference systems and coordinate transformations, the [EPSG Geodetic Parameter Dataset](http://www.epsg.org/). Within this collection each particular coordinate reference systems gets an unique integer identifier, commonly denoted as EPSG. For instance, the EPSG identifier for the the latest revision of the [World Geodetic System (WGS84)](https://en.wikipedia.org/wiki/World_Geodetic_System) is simply [4326](http://spatialreference.org/ref/epsg/4326/).A nice look up page for different coordinate reference systems is found [here](https://epsg.io/) and a fancy visualization of many prominent map projections is found [here](https://bl.ocks.org/mbostock/raw/3711652/).<jupyter_code>gdf.crs = {'init' :'epsg:4326'}<jupyter_output><empty_output><jupyter_text>Context matters: Load _Natural Earth countries_ dataset, bundled with GeoPandas[Natural Earth](http://www.naturalearthdata.com/) is a public domain map dataset available at 1:10m, 1:50m, and 1:110 million scales. Featuring tightly integrated vector and raster data, with Natural Earth you can make a variety of visually pleasing, well-crafted maps with cartography or GIS software. A subset comes bundled with GeoPandas and is accessible from the `gpd.datasets` module. We’ll use it as a helpful global base layer map.<jupyter_code>world = gpd.read_file(gpd.datasets.get_path('naturalearth_lowres'))
world.head(2)
world.crs
world.plot(facecolor='lightgray');<jupyter_output><empty_output><jupyter_text>**Combine world map and the Aerial Bombing data set**<jupyter_code>base = world.plot(facecolor='lightgray')
gdf.plot(ax=base, marker='o', color='red', markersize=10, alpha=0.01);<jupyter_output><empty_output><jupyter_text>As indicated above we want to focus on Europe. For that we first want to know how many times the targets in our data set we in Germany, France or Italy.<jupyter_code>percent = (int(np.round(gdf["Target Country"].
                        isin(["GERMANY", "FRANCE", "ITALY"]).sum()/gdf.shape[0] * 100)))

print("Approximately {}% of the targets in our data set are located in Germany, France or Italy.".format(percent))<jupyter_output><empty_output><jupyter_text>Then, we want to restrict our analysis to data points, which refer to an area within Europe. Although it is not straightforward to define Europe as an entity, in terms of geography, politics or sphere of cultural identity, we define Europe as an area between the coordinates  $$\text{33.0 to 73.5°N and 27.0°W to 45.0°E.}$$In order to represent that area spatially; we construct a `Polygon` object to represent the [bounding box](https://en.wikipedia.org/wiki/Minimum_bounding_box) of Europe.<jupyter_code>from shapely.geometry import Polygon
# generate geopandas object
poly_europe = gpd.GeoSeries([Polygon([(-27,33), (45,33), (45,73.5), (-27,73.5)])])
bb_europe = gpd.GeoDataFrame({'geometry': poly_europe})
# assign crs
bb_europe.crs = {'init':'epsg:4326'}
bb_europe<jupyter_output><empty_output><jupyter_text>**Plot world map and bounding box of Europe**<jupyter_code>base=world.plot(facecolor='lightgray')
bb_europe.plot(ax=base, alpha=0.5);<jupyter_output><empty_output><jupyter_text>**Subset (intersect) the GeoPandas `DataFrame` with the bounding box of Europe**<jupyter_code>gdf_europe = gpd.sjoin(gdf, bb_europe, how="inner", op='intersects').drop("index_right", axis=1)
print(gdf_europe.shape)
gdf_europe.sample(5)
gdf_europe.plot();<jupyter_output><empty_output><jupyter_text>**In order to keep the spatial context we extract the area of Europe from the world map**<jupyter_code>europe = gpd.overlay(world, bb_europe, how='intersection')
europe.crs = {'init': 'epsg:4326'}
base = europe.plot(facecolor='lightgray')
gdf_europe.plot(ax=base, marker='o', color='red', markersize=5, alpha=0.01);<jupyter_output><empty_output><jupyter_text>Ready!So now we have everything we need, a cleaned and subsetted data set `gdf_europe` and a spatial representation of Europe, `europe`, both in form of GeoPandas `DataFrame` objects.For further usage we write the GeoPandas `GeoDataFrame` objects to disk. `GeoDataFrames` can be exported to many different standard formats using the `GeoDataFrame.to_file()` method. Howerver, for the purpose of this tutorial we serialize the data by applying the [`pickle` module](https://docs.python.org/3/library/pickle.html).<jupyter_code>import pickle
pickle.dump(gdf_europe, open("./data/gdf_europe.p", "wb"))
pickle.dump(europe, open("./data/europe.p", "wb"))<jupyter_output><empty_output><jupyter_text>(Optional) Advanced spatial plotting using cartopy<jupyter_code>import os
import sys

# add the 'src' directory as one where we can import modules
src_dir = os.path.abspath(os.path.join(os.getcwd(), 'src'))
sys.path.append(src_dir)
print(src_dir)
import helper_funcs as hf
hf.cuteplot(gpd_df=gdf_europe, crs="Orthographic", title="Map projection: Orthographic")<jupyter_output><empty_output><jupyter_text>Note that the following loop may take some time, as each plot consists of 100,000+ data points.<jupyter_code>for proj in ["Mollweide", "Robinson", "PlateCarree", "UTM32N"]:
    hf.cuteplot(gpd_df=gdf_europe, crs=proj, title="Map projection: " + proj)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.10931195033626487, 'fraction_numerical': 0.013347128815312985, 'mean_word_length': 6.300226586102719, 'pattern_counts': {'":': 0, '<': 137, '<?xml version=': 0, '>': 145, 'https://': 21, 'lorem ipsum': 0, 'www.': 8, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '34', 'dir': 'jupyter-structured-clean-dedup', 'id': '79631', 'license': 'MIT', 'n_tokens_mistral': 6141, 'n_tokens_neox': 5802, 'n_words': 1925, 'path': 'Workshop_Data_Science/Exploratory Data Analysis WW2 - Data preparation.ipynb', 'repo_name': 'chrisbilliard/python', 'size': '19330'}