<jupyter_start><jupyter_text>Fundamentals of Information Systems Python Programming (for Data Science) Master's Degree in Data Science <NAME> (Courtesy of Gabriele Tolomei FIS 2018-2019)<EMAIL>University of Padua, Italy2019/2020  Lecture 12: The Regression Problem - Example (Part 1)  Instructions-  We consider the dataset available at this [link](http://archive.ics.uci.edu/ml/machine-learning-databases/auto-mpg/auto-mpg.data), provided by the [UC Irvine Machine Learning Repository](http://archive.ics.uci.edu/ml/index.php)-  The dataset comes with a [**README**](http://archive.ics.uci.edu/ml/machine-learning-databases/auto-mpg/auto-mpg.names) file, which contains information about the number of instances, the number and type of attributes, as well as the prediction goal.-  To be able to work even without any network connection, I have stored both the dataset and README files locally on my machine.-  "The data concerns city-cycle fuel consumption in miles per gallon, to be predicted in terms of 3 multivalued discrete and 5 continuous attributes." (Quinlan, 1993)  Additional Notes-  The original file contains a mixture of whitespace and tab-separated fields. In order to transform it into a legitimate tab-separated file, I had to run the following shell commands:```bash> TAB=$(printf '\t')> sed "s/ \{2,\}/$TAB/g"  ${NEW_DATASET_FILE}```-  The commands above use **sed** to replace **2 or more** whitespaces with a tab character on **every line** of the original dataset file.-  Linux and Mac OS X systems have **sed** natively installed. Windows users can install it from [here](http://gnuwin32.sourceforge.net/packages/sed.htm)<jupyter_code>import math
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
# Import stats module from scipy, which contains a large number 
# of probability distributions as well as an exhaustive library of statistical functions.
import scipy.stats as stats

%matplotlib inline<jupyter_output><empty_output><jupyter_text>1. Data Collection<jupyter_code># Path to the local dataset file
DATASET_PATH = "./data/auto-mpg-regression/dataset.tsv"
# Load the dataset with Pandas
data = pd.read_csv(DATASET_PATH, sep="\t")
print("Shape of the dataset: {}".format(data.shape))
data.head()
# NOTE: the first line of the file is considered as the header
# Load the dataset with Pandas, this time taking into account
# the fact that there is no header line
data = pd.read_csv(DATASET_PATH, sep="\t", header=None)
print("Shape of the dataset: {}".format(data.shape))
data.head()
# Column names (labels) are not so meaningful and we should use the attribute names
# provided in the README file
# Row index, instead, can be left as it is (i.e., default IndexRange)
columns = ['mpg', 'cylinders', 'displacement', 'horsepower', 'weight', 
           'acceleration', 'model_year', 'origin', 'car_name']
data = pd.read_csv(DATASET_PATH, sep="\t", header=None, names=columns)
print("Shape of the dataset: {}".format(data.shape))
data.head()<jupyter_output>Shape of the dataset: (398, 9)<jupyter_text>Checking for any Missing Values<jupyter_code># Check if there is any missing value in the whole dataset
print("There are missing values in the dataset: {}".
     format(data.isnull().any().any()))
# Check if 'horsepower' has really 6 missing values
print("N. of missing values for attribute 'horsepower': {}".
     format(data.horsepower.isnull().sum()))
# Weird... we know from the README file there must be 6 missing values
# How come we are not able to spot those?
# Possibly, because of the way in which NAs are 'encoded' in the file...
# Let's see the set of values contained in the 'horsepower' column
print(sorted(data.horsepower.unique(), reverse=True))
# Apparently, there are some question mark characters '?' in this column.
# Let's see how many records have '?' in their 'horsepower' column.
# Extract the sub-DataFrame using boolean indexing 
# on the 'horsepower' column and count the corresponding number of matching rows.
print("How many records have 'horsepower=?'?: {}"
      .format(data[data.horsepower == '?'].shape[0]))
data.loc[data.horsepower == "?", :]
# We have therefore found that '?' is a sentinel value used to identify NAs
# Let's reload the dataset using this value as a marker for NAs.
data = pd.read_csv(DATASET_PATH, sep="\t", header=None, 
                   names=columns,
                  na_values={'horsepower':'?'})

# Alternatively, we could simply replace '?' on the loaded dataset with np.nan
# 1. Using 'loc':
# data.loc[data.horsepower == '?', 'horsepower'] = np.nan
# 2. Using 'replace':
# data.horsepower.replace('?', np.nan, inplace=True)
print("Shape of the dataset: {}".format(data.shape))
data.head()
# Let's repeat the same check as above on missing values.
# Check if there is any missing value in the whole dataset
print("There are missing values in the dataset: {}".
     format(data.isnull().any().any()))
# Check if 'horsepower' has really 6 missing values
print("N. of missing values for attribute 'horsepower': {}".
     format(data.horsepower.isnull().sum()))
# Let's have a look at the output of the 'describe()' function.
data.describe(include='all')
print(data["mpg"])<jupyter_output>0      18.0
1      15.0
2      18.0
3      16.0
4      17.0
       ... 
393    27.0
394    44.0
395    32.0
396    28.0
397    31.0
Name: mpg, Length: 398, dtype: float64<jupyter_text>Change the Column Layout<jupyter_code># Just as a convention, I prefer to place the column to be predicted
# as the last one.
columns = data.columns.tolist()
print("Orignal order of columns:\n{}".format(columns))
# Popping out 'mpg' from the list and insert it back at the end.
columns.insert(len(columns), columns.pop(columns.index('mpg')))
print("New order of columns:\n{}".format(columns))
# Let's refactor the DataFrame using this new column index
data = data.loc[:, columns]
data.head()
# Alternatively to 'loc' we can also use 'reindex()'
# data = data.reindex(columns=columns)
# data.head()<jupyter_output><empty_output><jupyter_text>Apply Some Simple Data Transformations<jupyter_code># Suppose we want to convert displacement unit from cubic inch to litre
# There is a useful conversion table which tells us how to do that.
# 1 cubic inch = 0.016387064 litre
CI_TO_LITRE = 0.016387064
data.displacement *= CI_TO_LITRE
data.head()<jupyter_output><empty_output><jupyter_text>2. Data Exploration  2.1 Analysis of Data Distributions: Continous Values-  Let's start visualizing the distributions of the **5 continuous-valued** features:    - **displacement**    - **horsepower** (contains 6 **NA**s)    - **weight**    - **acceleration**<jupyter_code># Create a lambda function which will be applied to each entry
# of the numpy 2-D array of AxesSubplot objects
# x is a reference to an AxesSubplot object
y_labeler = lambda x: x.set_ylabel('density')
# np.vectorize() allows calling the function on each element
y_labeler = np.vectorize(y_labeler)
fig, axes = plt.subplots(2, 2, figsize=(10,12))
y_labeler(axes)
# Create a Figure containing 2x2 subplots
fig, axes = plt.subplots(2, 2, figsize=(10,12))
# Call the vectorized function for labeling all the y-axes
y_labeler(axes)
# Plot 'displacement' on the top-left subplot
sns.distplot(data.displacement, color='#808080', ax=axes[0,0], 
                 hist_kws=dict(edgecolor="#404040", linewidth=1))
# Plot 'horsepower' (limited only to non-NA values) on the top-right subplot
sns.distplot(data.loc[data.horsepower.notnull(), 'horsepower'], 
                 color='#df2020', ax=axes[0,1], 
                 hist_kws=dict(edgecolor="#404040", linewidth=1))
# Plot 'weight' on the bottom-left subplot
sns.distplot(data.weight, color='#0033cc', ax=axes[1,0], 
                 hist_kws=dict(edgecolor="k", linewidth=1))
# Plot 'acceleration' on the bottom-right subplot
sns.distplot(data.acceleration, color='#009933', ax=axes[1,1],
                 hist_kws=dict(edgecolor="#006622", linewidth=1))
# Adjust space between plots
plt.subplots_adjust(wspace=.3, hspace=.3)<jupyter_output><empty_output><jupyter_text>Additional Note on Binning-  If we don't specify the number of **bins** as argument of **sns.distplot** function (i.e.,**bins=None**) the **Freedman-Diaconis** rule is used to devise the _best_ number of bins.-  This rule starts from defining the **width** each bin should have on the basis of the range of values observed, as follows:$$\texttt{bin}\_\texttt{width} = 2 * \frac{\texttt{IQR}}{n^{1/3}}$$where $\texttt{IQR}$ stands for **interquartile range**, namely the length of the interval delimited by the 1st and the 3rd quartile, and $n$ is the **number of observations**.-  Finally, the number of bins is computed as:$$\texttt{bins} = (max−min)/\texttt{bin}\_\texttt{width}$$where, $max$ ($min$) is the **maximum** (**minimum**) value observed.<jupyter_code># Let's produce the boxplots corresponding to the distribution plots above
# Create a Figure containing 2x2 subplots
fig, axes = plt.subplots(2, 2, figsize=(10,12))

sns.boxplot(data.displacement, color='#808080', ax=axes[0,0])
sns.boxplot(data.loc[data.horsepower.notnull(), 'horsepower'], 
                color='#df2020', ax=axes[0,1])
sns.boxplot(data.weight, color='#0033cc', ax=axes[1,0])
sns.boxplot(data.acceleration, color='#009933', ax=axes[1,1])
plt.subplots_adjust(wspace=.3, hspace=.3)<jupyter_output><empty_output><jupyter_text>Spotting Outliers-  Generally speaking, an outlier is an observation that is numerically distant from the rest of the data. -  Boxplots are useful to actually spot any possible **outlier**, as they show the distribution of values that are located within the $\texttt{IQR}$ (i.e., any data point between the 1st and 3rd quartile).-  A boxplot defines also 2 other values, called **fences** or **whiskers** which are used to define outliers (i.e., any data point that is located outside the fences).-  Usually, fences are determined as follows: if $Q_1$ and $Q_3$ represents the 1st and 3rd quartile, respectively, we define $F_\textrm{left}$ and $F_\textrm{right}$ as the left and right fence point, respectively, so that:$$F_\textrm{left} = Q_1 - 1.5 * \texttt{IQR};~~F_\textrm{right} = Q_3 + 1.5 * \texttt{IQR}$$  Few Observations from the Plots-  **displacement** distribution is not uni-modal at all; in fact it is bi-, or possibly, tri-modal. Apparently, this is not affected by any outlier.-  **horsepower** distribution is also not uni-modal, with a small bump around the value of 150Hp. However, here we can spot some possible outliers (i.e., large values over 200Hp).-  **weight** distribution is uni-modal yet **right-skewed** (**positively skewed**), which means the mean is shifted to the right due to the presence of some large values towards the positive direction (not necessarily outliers).-  **acceleration** essentially fits nicely to a **Normal** (**Gaussian**) **distribution** and contains some outliers.<jupyter_code># Check how many outliers we have
# 1. 'horsepower'
hp_q1, hp_q3 = data.loc[data.horsepower.notnull(), 'horsepower'].quantile([.25, .75])
print("1st Quartile of 'horsepower': {:.2f}".format(hp_q1))
print("3rd Quartile of 'horsepower': {:.2f}".format(hp_q3))
hp_IQR = (hp_q3 - hp_q1)
print("IQR of 'horsepower': {:.2f}".format(hp_IQR))
hp_fence_left = hp_q1 - 1.5 * hp_IQR
hp_fence_right = hp_q3 + 1.5 * hp_IQR
print("Fence range: [{:.2f}, {:.2f}]".format(hp_fence_left, hp_fence_right))
print("N. of instances containing outlier of 'horsepower': {}".
      format(data[data.horsepower > hp_fence_right].shape[0]))
# Check how many outliers we have
# 2. 'acceleration'
acc_q1, acc_q3 = data['acceleration'].quantile([.25, .75])
print("1st Quartile of 'acceleration': {:.2f}".format(acc_q1))
print("3rd Quartile of 'acceleration': {:.2f}".format(acc_q3))
acc_IQR = (acc_q3 - acc_q1)
print("IQR of 'acceleration': {:.2f}".format(hp_IQR))
acc_fence_left = acc_q1 - 1.5 * acc_IQR
acc_fence_right = acc_q3 + 1.5 * acc_IQR
print("Fence range: [{:.2f}, {:.2f}]".format(acc_fence_left, acc_fence_right))
print("N. of instances containing outlier of 'acceleration': {}".
      format(data[(data.acceleration < acc_fence_left) | 
                  (data.acceleration > acc_fence_right)].shape[0]))
print("N. of instances containing outlier of 'horsepower' or 'acceleration': {}".
      format(data[(data.horsepower > hp_fence_right) |
                  (data.acceleration > acc_fence_right) |
                  (data.acceleration < acc_fence_left)].shape[0]))
# Let's see if we can spot where those outliers are located
# w.r.t. other features (e.g., cylinders)
fig, axes = plt.subplots(2, 1, figsize=(10,8))

sns.boxplot(x=data.cylinders, y=data.loc[data.horsepower.notnull(), 'horsepower'], 
                color='#df2020', ax=axes[0])
sns.boxplot(x=data.cylinders, y=data.acceleration, color='#009933', ax=axes[1])
plt.subplots_adjust(wspace=.3, hspace=.3)<jupyter_output><empty_output><jupyter_text>Other Few Observations from the Plots-  **horsepower**: outliers seem to occur only on 6- and 8-cylinder vehicles, stronger in the former.-  **acceleration**: outliers here are more evenly distributed over 3- 4- and 8-cylinder vehicles.<jupyter_code># Let's now plot the pairwise relationship between our continuous-valued features
sns.pairplot(data.loc[data.horsepower.notnull(),
                          ['displacement', 'horsepower', 'weight', 
                           'acceleration']],
                 kind="reg",
                 diag_kind='kde', 
                 diag_kws={'shade': True, 'color': '#ff6600'}, 
                 plot_kws={'color': '#ff6600'})<jupyter_output><empty_output><jupyter_text>How Continuous-valued Features Relate to Each Other-  As **displacement** increases, so do **weight** and **horsepower**; also **acceleration** tends to decrease (pretty intuitive!).-  As **horsepower** increases, **acceleration** decreases, whilst **weight** and **displacement** increase as well.-  As **weight** increases, both **horsepower** and **displacement** increases but **acceleration** tends to decrease. Now, this might seem counterintuitive at first but the thing is that there might be a **latent factor** which affects this relationship (e.g., a heavy vehicle is also likely to have more horsepower and therefore this could be the actual reason why we observe such a phenomenon).<jupyter_code># Let's now plot the pairwise relationship between our continuous-valued features 
# this time also considering our target variable 'mpg'
sns.pairplot(data.loc[data.horsepower.notnull(),
                          ['displacement', 'horsepower', 'weight', 'acceleration', 'mpg']],
                 kind="reg",
                 diag_kind='kde', 
                 diag_kws={'shade': True, 'color': '#33cccc'}, 
                 plot_kws={'color': '#33cccc'})<jupyter_output><empty_output><jupyter_text>2.1 Analysis of Data Distributions: Categorical Values-  Let's visualize the distributions of the **3 categorical** features:    - **cylinders**    - **model_year**    - **origin**<jupyter_code># Let's see the frequency counts of each categorical variable
# 'cylinders'
print(data.cylinders.value_counts())
print()
# 'model_year'
print(data.model_year.value_counts())
print()
# 'origin'
print(data.origin.value_counts())
# Let's produce the boxplots corresponding to the distribution plots above
# Create a Figure containing 2x2 subplots
fig, axes = plt.subplots(1, 3, figsize=(8,4))
y_labeler(axes)

# Plot 'cylinders'
sns.distplot(data.cylinders, color='#006699', ax=axes[0], 
                 kde=False, norm_hist=True,
                 hist_kws=dict(edgecolor="#404040", linewidth=1))
# Plot 'model_year'
sns.distplot(data.model_year, color='#a6cc33', ax=axes[1], 
                 kde=False, norm_hist=True,
                 hist_kws=dict(edgecolor="#85a329", linewidth=1))
# Plot 'origin'
sns.distplot(data.origin, color='#cc3399', ax=axes[2], 
                 kde=False, norm_hist=True,
                 hist_kws=dict(edgecolor="#8f246b", linewidth=1))
# Adjust space between plots
plt.subplots_adjust(wspace=.5, hspace=.4)
# For categorical variables, 'countplot' is the way to go
# Create a Figure containing 1x3 subplots
fig, axes = plt.subplots(1, 3, figsize=(12,4))

# Plot 'cylinders'
sns.countplot(data.cylinders, ax=axes[0])
sns.countplot(data.model_year, ax=axes[1])
sns.countplot(data.origin, ax=axes[2])
plt.subplots_adjust(wspace=.5, hspace=.4)
# stripplot is also another useful plot to relate categorical vs. target variable
# Create a Figure containing 1x3 subplots
fig, axes = plt.subplots(1, 3, figsize=(12,4))

sns.stripplot(x=data.cylinders, y=data.mpg, color='#006699', ax=axes[0])
sns.stripplot(x=data.model_year, y=data.mpg, color='#a6cc33', ax=axes[1])
sns.stripplot(x=data.origin, y=data.mpg, color='#cc3399', ax=axes[2])

plt.subplots_adjust(wspace=.5, hspace=.3)
# swarmplot is also another useful plot to relate categorical vs. target variable
# Create a Figure containing 1x3 subplots
fig, axes = plt.subplots(1, 3, figsize=(12,4))

sns.swarmplot(x=data.cylinders, y=data.mpg, color='#006699', ax=axes[0])
sns.swarmplot(x=data.model_year, y=data.mpg, color='#a6cc33', ax=axes[1])
sns.swarmplot(x=data.origin, y=data.mpg, color='#cc3399', ax=axes[2])

plt.subplots_adjust(wspace=.5, hspace=.3)
# boxplot is also another useful plot to relate categorical vs. target variable
# Create a Figure containing 1x3 subplots
fig, axes = plt.subplots(1, 3, figsize=(16,6))

sns.boxplot(x=data.cylinders, y=data.mpg, color='#006699', ax=axes[0])
sns.boxplot(x=data.model_year, y=data.mpg, color='#a6cc33', ax=axes[1])
sns.boxplot(x=data.origin, y=data.mpg, color='#cc3399', ax=axes[2])

plt.subplots_adjust(wspace=.5, hspace=.3)
# Let's overlay boxplot and swarmplot
# Create a Figure containing 1x3 subplots
fig, axes = plt.subplots(1, 3, figsize=(16,6))

sns.boxplot(x=data.cylinders, y=data.mpg, color='#006699', ax=axes[0])
sns.swarmplot(x=data.cylinders, y=data.mpg, color=".5", ax=axes[0])
sns.boxplot(x=data.model_year, y=data.mpg, color='#a6cc33', ax=axes[1])
sns.swarmplot(x=data.model_year, y=data.mpg, color=".5", ax=axes[1])
sns.boxplot(x=data.origin, y=data.mpg, color='#cc3399', ax=axes[2])
sns.swarmplot(x=data.origin, y=data.mpg, color=".5", ax=axes[2])

plt.subplots_adjust(wspace=.5, hspace=.3)<jupyter_output><empty_output><jupyter_text>3. Data Preprocessing (Munging)  Summary of the Issues-  From our exploratory data analysis above, **two** main issues are observed:    1. The presence of **6 missing values** for the attribute **horsepower**    2. The presence of a total of **18 outliers** on the attributes **horsepower** and **acceleration**-  In addition to those, we should also consider how to properly handle different feature's scale as well as the fact that we are in presence of both continuous and categorical attributes.  3.1 Handling Missing Values (NA)-  There are just **6 out of 398** (i.e., approximately **1.5%**) of records containing a missing value.-  Since they do not represent a significant subset of the whole dataset, we can simply drop those records.-  Otherwise, we could mark (i.e., replace) those missing value using one of the strategies discussed (e.g., replace them with the median as computed from observed values).<jupyter_code># Let's go for the second option, i.e., replacing missing values on 'horsepower'
# using the median as computed from the other observations.
# NOTE: here's a classical example where using the mean rather than the median
# might affect the result, as the mean is more sensitive to outliers.
# NOTE: by default, median() does not include NAs in the computation.
# In other words, we don't need to explicitly tell pandas to work on non-NA values:
# data.horsepower[data.horsepower.notnull()].median()
data.horsepower.fillna(data.horsepower.median(), inplace=True)
data.head()
data.describe(include='all')<jupyter_output><empty_output><jupyter_text>3.2 Handling Outliers-  There are **18 outliers** shared between **horsepower** and **acceleration**.-  Like missing values, outliers can be simply discarded as well (i.e., a process which is also known as **trimming** or **truncation**).-  Another approach is called **winsorizing** and consists of replacing outliers with a specified percentile of the data (e.g., a 90% winsorization would see all data below the 5th percentile set to the 5th percentile, and data above the 95th percentile set to the 95th percentile).<jupyter_code># Python can winsorize data using 'scipy.stats' module.
# Example:
a = pd.Series([92, 19, 101, 58, 1053, 91, 26, 78, 10, 13, -40, 101, 86, 85, 15, 89, 89, 28, -5, 41])
print("Length of a: {}".format(a.shape[0]))
print("Mean of a: {}".format(a.mean()))
print("Median of a: {}".format(a.median()))
print("Sorted a: {}".format(np.sort(a)))
q_005, q_95 = a.quantile([0.05, 0.95])
print("5th percentile of a: {:.2f}".format(q_005))
print("95th percentile of a: {:.2f}".format(q_95))
stats.mstats.winsorize(a, limits=0.05, inplace=True)
print("Sorted a: {}".format(np.sort(a)))
# Let's winsorize 'horsepower' and 'acceleration'
stats.mstats.winsorize(data.horsepower, limits=0.0375, inplace=True)
stats.mstats.winsorize(data.acceleration, limits=0.0375, inplace=True)
data.describe(include='all')
# Let's verify the outliers are actually gone
fig, axes = plt.subplots(2, 2, figsize=(10,12))
y_labeler(axes)
sns.boxplot(data.displacement, color='#808080', ax=axes[0,0])
sns.boxplot(data.horsepower, color='#df2020', ax=axes[0,1])
sns.boxplot(data.weight, color='#0033cc', ax=axes[1,0])
sns.boxplot(data.acceleration, color='#009933', ax=axes[1,1])
plt.subplots_adjust(wspace=.3, hspace=.3)<jupyter_output><empty_output><jupyter_text>3.3 Encoding Categorical Features-  Categorical variables are typically stored as text values which represent various traits. -  Some examples include **color** = {"Red", "Yellow", "Blue"), **size** =  ("Small", "Medium", "Large"), etc.-  Many ML algorithms can support categorical values without further manipulation but there are many others that do not. -  Therefore, the analyst is faced with the challenge of figuring out how to turn these text attributes into **numerical values** for further processing.  Approach 1: Label Encoding-  Label encoding is simply converting each value in a column to a number. -  For example, the **model_year** column contains 13 different values. We could choose to encode it like this:```python70 --> 071 --> 172 --> 2...```  Approach 2: One-Hot Encoding-  Label encoding is straightforward but it has the disadvantage that numeric values can be "misinterpreted" by the learning algorithms. -  For example, the value of 0 is obviously less than the value of 4 but is that what we really aim for? For example, does a vehicle from 1974 have "4x" more weight than one from 1971?-  A common alternative approach is called **one hot encoding**. Here, the basic strategy is to convert each category value into a new column and assigns a 1 or 0 (True/False) value to the column. -  This has the benefit of not weighting a value improperly but does have the downside of adding more columns to the data set.<jupyter_code># In pandas we can achieve easily one-hot encoding using the 'get_dummies()' function
categorical_features = ['cylinders', 'model_year', 'origin']
data_with_dummies = pd.get_dummies(data, columns = categorical_features)
data_with_dummies.head()
# Just as a convention, I prefer to place the column to be predicted
# as the last one.
columns = data_with_dummies.columns.tolist()
# Popping out 'mpg' from the list and insert it back at the end.
columns.insert(len(columns), columns.pop(columns.index('mpg')))
# Let's refactor the DataFrame using this new column index
data_with_dummies = data_with_dummies.loc[:, columns]
data_with_dummies.head()
# The categorical variable 'car_name' contains a lot of different values.
# Using one-hot encoding might lead to a very sparse dataset, as we need
# to map a single column to 305 columns!
data_with_dummies.car_name.value_counts()
"""
Three solutions can be designed to tackle with this issue:
1) Just drop the column 'car_name' (i.e., our model won't rely on that feature for prediction)
2) Use one-hot encoding scheme and deal with sparsity data (i.e., possibly leading to overfitting)
3) Trade-off: try to build another column which somehow reduces (i.e., cluster) similar values together
and then apply one-hot encoding.
Let's see how to perform 3)
"""
# Suppose we want to create another column called 'automaker_name', which simply contains
# the name of the automaker, disregarding the model.
# For example, automaker_name('<NAME>') = automaker_name('<NAME>') = 'ford'
data_with_dummies['automaker_name'] = data_with_dummies['car_name'].map(lambda x:
                                                                       x.split(' ')[0])
data_with_dummies.head()
data_with_dummies['car_name'].map(lambda x: x.split(' ')[0]).value_counts()
"""
Let's see how many distinct values we have now for this nvalue_countsrical variable
"""
data_with_dummies['automaker_name'].value_counts()
def sanitize_automaker_name(car_name):
    s = car_name.split(' ')[0]
    if s == 'vw' or s == 'vokswagen':
        return car_name.replace(s,'volkswagen')
    if s == 'chevroelt' or s == 'chevy':
        return car_name.replace(s,'chevrolet')
    if s == 'maxda':
        return car_name.replace(s,'mazda')
    if s == 'mercedes':
        return car_name.replace(s,'mercedes-benz')
    if s == 'toyouta':
        return car_name.replace(s,'toyota')
    return car_name
"""
Use the 'sanitize_automaker_name' function to update 'car_name' values.
"""
data_with_dummies['car_name'] = data_with_dummies['car_name'].map(lambda x:
                                                                 sanitize_automaker_name(x))
data_with_dummies.head()
"""
Re-apply the function on the sanitized version of 'car_name'.
"""
data_with_dummies['automaker_name'] = data_with_dummies['car_name'].map(lambda x: 
                                                                        x.split(' ')[0])
data_with_dummies.head()
"""
Let's see how many distinct values we have now for this new categorical variable.
"""
data_with_dummies['automaker_name'].value_counts()
# Create the set of the top-10 automakers
top_10_automakers = set(data_with_dummies['automaker_name'].value_counts().index[:10])
# Label with 'other' any automaker_name which is not in the list above
data_with_dummies['automaker_name'] = np.where(data_with_dummies.
                                               automaker_name.isin(top_10_automakers), 
                                              data_with_dummies.automaker_name,
                                              'other')
"""
Let's verify we did it right!
"""
data.loc[data_with_dummies.automaker_name == 'other', 'car_name'].map(lambda x: 
                                                                      x.split(' ')[0] not in top_10_automakers).all()
"""
We now categorize (i.e., discretize) 'automaker_name' using the 10 + 1 discrete values above.
"""
categorical_features = ['automaker_name']
data_with_dummies = pd.get_dummies(data_with_dummies, columns = categorical_features)
data_with_dummies.head()
# Just as a convention, I prefer to place the column to be predicted
# as the last one.
columns = data_with_dummies.columns.tolist()
# Popping out 'mpg' from the list and insert it back at the end.
columns.insert(len(columns), columns.pop(columns.index('mpg')))
# Popping out 'automaker_name_other' from the list and insert it after 'automaker_name_volkswagen'.
columns.insert(columns.index('automaker_name_volkswagen'), columns.pop(columns.index('automaker_name_other')))
# Popping out 'car_name' from the list and insert it right before 'mpg'
columns.insert(-1, columns.pop(columns.index('car_name')))
# Let's refactor the DataFrame using this new column index
data_with_dummies = data_with_dummies.loc[:, columns]
data_with_dummies.head()<jupyter_output><empty_output><jupyter_text>3.4 Standardize Feature Scale-  Some learning models are sensitive to different feature scales appearing on the training dataset.-  To overcome this issue, one typically standardize (i.e., normalize) the values of each continuous feature.-  Two main strategies are usually enacted:    -  **min-max** normalization    -  **z-score** standardization<jupyter_code># Let's deep copy our DataFrame again
data_norm_0_1 = data_with_dummies.copy()
"""
The easiest way to normalize a (sub)set of features is as follows.
"""
# 1. Decide which list of features to standardize
features_to_standardize = ['displacement', 'horsepower', 'weight', 'acceleration']

# 2. Select those features (i.e., DataFrame columns) and apply, for example, min-max normalization
data_norm_0_1[features_to_standardize] = ((data_norm_0_1[features_to_standardize] - 
                                          data_norm_0_1[features_to_standardize].min()) 
                                          / (data_norm_0_1[features_to_standardize].max() - 
                                             data_norm_0_1[features_to_standardize].min()))
# 3. Verify the result
data_norm_0_1.head()
# Let's make another deep copy of our DataFrame
data_std = data_with_dummies.copy()
"""
The following three functions are used to standardize/normalize features.
"""
# 1. z_score computes the standard z-score of a feature value x
def z_score(x, mu_X, sigma_X):
    return (x - mu_X)/sigma_X

# 2. min_max computes the normalized value of a feature value x in the range [-1, 1]
def min_max(x, X_min, X_max):
    return (2*x - X_max - X_min)/(X_max - X_min)

# 3. min_max_0_1 computes the normalized value of a feature value x in the range [0, 1]
def min_max_0_1(x, X_min, X_max):
    return (x - X_min)/(X_max - X_min)
"""
We use 'map' to call the z_score function above element-wise
on each Series object: displacement, horsepower, weight, and acceleration.
"""
# 1.a. Compute the mean and std deviation of 'displacement'
mu_displacement = data_std.displacement.mean()
sigma_displacement = data_std.displacement.std()
# 1.b. Call the z_score function on the 'displacement' Series
data_std.displacement = data_std.displacement.map(lambda x: 
                                                  z_score(x, mu_displacement, 
                                                          sigma_displacement))

# 2.a. Compute the mean and std deviation of 'horsepower'
mu_horsepower = data_std.horsepower.mean()
sigma_horsepower = data_std.horsepower.std()
# 2.b. Call the z_score function on the 'horsepower' Series
data_std.horsepower = data_std.horsepower.map(lambda x: 
                                              z_score(x, mu_horsepower, 
                                                      sigma_horsepower))

# 3.a. Compute the mean and std deviation of 'weight'
mu_weight = data_std.weight.mean()
sigma_weight = data_std.weight.std()
# 3.b. Call the z_score function on the 'weight' Series
data_std.weight = data_std.weight.map(lambda x: 
                                      z_score(x, mu_weight, 
                                              sigma_weight))

# 4.a. Compute the mean and std deviation of 'acceleration'
mu_acceleration = data_std.acceleration.mean()
sigma_acceleration = data_std.acceleration.std()
# 4.b. Call the z_score function on the 'acceleration' Series
data_std.acceleration = data_std.acceleration.map(lambda x: 
                                                  z_score(x, mu_acceleration, 
                                                          sigma_acceleration))

data_std.head()
# Let's create yet another deep copy of our DataFrame
data_std_z = data_with_dummies.copy()
"""
This is an even more general solution for standardizing multiple features
in one shot, using customized "normalizing" functions.
"""
def standardized_continuous_features(dataset, feature_names, func=z_score):
    
    if func != z_score and func != min_max and func != min_max_0_1:
        func = z_score
            
    for feature in feature_names:
        print("Standardized feature \"{}\" using [{}] function".format(feature, func.__name__))
        if func == min_max or func == min_max_0_1:
            feature_min = dataset[feature].min()
            feature_max = dataset[feature].max()
            print("Min. = {}".format(feature_min))
            print("Max. = {}".format(feature_max))
            dataset[feature] = dataset[feature].map(lambda x: func(x, feature_min, feature_max))
        else:
            feature_mean = dataset[feature].mean()
            feature_std = dataset[feature].std()
            print("Mean = {}".format(feature_mean))
            print("Std. Deviation = {}".format(feature_std))
            dataset[feature] = dataset[feature].map(lambda x: func(x, feature_mean, feature_std))
    
    return dataset
"""
We call the function defined above on our deep copy, using the list of features
that needs to be standardized.
"""
data_std_z = standardized_continuous_features(data_std_z, features_to_standardize)
data_std_z.head()
"""
Let's verify the two different approaches above lead to the same result."
"""
print((np.abs(data_std[features_to_standardize] - data_std_z[features_to_standardize]) < 0.0001).all())
"""
Standardization and min-max scaling can be also performed using scikit-learn.
"""
# The following is the scikit-learn package which provides
# various preprocessing capabilities
from sklearn import preprocessing

# Standardizing features using z-score
std_scale = preprocessing.StandardScaler().fit(data_with_dummies[features_to_standardize])
data_std = std_scale.transform(data_with_dummies[features_to_standardize])

# Normalizing features using min-max
minmax_scale = preprocessing.MinMaxScaler().fit(data_with_dummies[features_to_standardize])
data_minmax = minmax_scale.transform(data_with_dummies[features_to_standardize])

# NOTE: 'data_std' and 'data_minmax' are numpy's ndarray (i.e., not pandas' DataFrame) objects<jupyter_output><empty_output><jupyter_text>3.5 Feature Selection-  This is a topic that would require an in-depth analysis.-  A very simple approach to select highly discriminant features is given by measuring how each feature correlated with the target, and pick the ones with the highest correlation score.<jupyter_code>"""
Let's compute feature correlation between each feature and our target 'mpg'.
Pandas 'corr()' function when applied to a DataFrame returns the whole correlation matrix.
This is a diagonal matrix having all 1's on its diagonal, whilst the entry (i,j) will contain
the correlation coefficient between feature i and feature j (i != j).
"""
corr_matrix = data_with_dummies.corr()
strong_corr = corr_matrix[(corr_matrix['mpg'].abs() > .5) &
                         (corr_matrix['mpg'] != 1)].loc[:, 'mpg']

print("Strongest correlated features:\n{}".format(strong_corr.sort_values()))
print()
print("Strongest correlated features (absolute values):\n{}".
      format(strong_corr.abs().sort_values(ascending=False)))<jupyter_output>Strongest correlated features:
weight         -0.831741
displacement   -0.804203
horsepower     -0.788402
cylinders_8    -0.647308
origin_1       -0.568192
cylinders_4     0.758259
Name: mpg, dtype: float64

Strongest correlated features (absolute values):
weight          0.831741
displacement    0.804203
horsepower      0.788402
cylinders_4     0.758259
cylinders_8     0.647308
origin_1        0.568192
Name: mpg, dtype: float64
{'fraction_non_alphanumeric': 0.12012775940529692, 'fraction_numerical': 0.0237711636846717, 'mean_word_length': 4.169491525423729, 'pattern_counts': {'":': 0, '<': 69, '<?xml version=': 0, '>': 76, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 40, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '16', 'dir': 'jupyter-structured-clean-dedup', 'id': '666891', 'license': 'MIT', 'n_tokens_mistral': 11833, 'n_tokens_neox': 10785, 'n_words': 3555, 'path': 'Small Projects/Auto Regression/Lecture_17_The_Regression_Problem_Example_(Part_1).ipynb', 'repo_name': 'GleydsonHudson/FIS-Exercises-Unipd', 'size': '35432'}