<jupyter_start><jupyter_text>Notebook inspired by the lecture _Data Analysis with Python_ by Cognitive Class / IBM<jupyter_code>import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
from sklearn.linear_model import LinearRegression
%matplotlib inline 
path="data/clean/auto.csv"
df = pd.read_csv(path)<jupyter_output><empty_output><jupyter_text>Simple linear regressionFor this example, we want to look at how highway-mpg can help us predict car price. Using simple linear regression, we will create a linear function with "highway-L/100km" as the predictor variable and the "price" as the response variable.<jupyter_code>lm = LinearRegression()

X = df[['highway-L/100km']]
Y = df['price']

lm.fit(X,Y)

Yhat=lm.predict(X)
Yhat[0:5]

print(f'Yhat = {lm.intercept_} + {lm.coef_} X')<jupyter_output>Yhat = -14617.843054664594 + [3458.68511314] X<jupyter_text>Let's do the same exercize with the variable 'engine-size'.<jupyter_code>lm = LinearRegression()

X = df[['engine-size']]
Y = df['price']

lm.fit(X,Y)

Yhat=lm.predict(X)
Yhat[0:5]

print(f'Yhat = {lm.intercept_} + {lm.coef_} X')<jupyter_output>Yhat = -7963.338906281049 + [166.86001569] X<jupyter_text>Multiple Linear RegressionFrom the previous analysis we know that other good predictors of price could be:- Horsepower- Curb-weight- Engine-size- Highway-mpg<jupyter_code>Z = df[['horsepower', 'curb-weight', 'engine-size', 'highway-L/100km']]
lm.fit(Z, df['price'])
print('intercept:', lm.intercept_)
print('coefs:', lm.coef_)<jupyter_output>intercept: -14385.634549360107
coefs: [ 36.6377371    3.50038215  85.37370862 500.51979785]<jupyter_text>As we saw above, we should get a final linear function with the structure:$$Yhat = a + b_1 X_1 + b_2 X_2 + b_3 X_3 + b_4 X_4$$**Price** = -15678.742628061467 + 52.65851272 x **horsepower** + 4.69878948 x **curb-weight** + 81.95906216 x **engine-size** + 33.58258185 x **highway-L/100km**  Model Evaluation using VisualizationA good way to visualize the variance of the data is to use a residual plot. What is a residual?The difference between the observed value (y) and the predicted value (Yhat) is called the residual (e). When we look at a regression plot, the residual is the distance from the data point to the fitted regression line. So what is a residual plot?A residual plot is a graph that shows the residuals on the vertical y-axis and the independent variable on the horizontal x-axis. What do we pay attention to when looking at a residual plot?We look at the spread of the residuals: if the points in a residual plot are randomly spread out around the x-axis, then a linear model is appropriate for the data. Why is that? Randomly spread out residuals means that the variance is constant, and thus the linear model is a good fit for this data.<jupyter_code>width = 12
height = 10
plt.figure(figsize=(width, height))
sns.residplot(df['highway-L/100km'], df['price'])
plt.show()<jupyter_output><empty_output><jupyter_text>What is this plot telling us?We can see from this residual plot that the residuals are not randomly spread around the x-axis, which leads us to believe that maybe a non-linear model is more appropriate for this data. Multiple Linear RegressionHow do we visualize a model for Multiple Linear Regression? This gets a bit more complicated because you can't visualize it with regression or residual plot.One way to look at the fit of the model is by looking at the distribution plot: we can look at the distribution of the fitted values that result from the model and compare it to the distribution of the actual values.First lets make a prediction<jupyter_code>Y_hat = lm.predict(Z)
plt.figure(figsize=(width, height))


ax1 = sns.distplot(df['price'], hist=False, color="r", label="Actual Value")
sns.distplot(Yhat, hist=False, color="b", label="Fitted Values" , ax=ax1)


plt.title('Actual vs Fitted Values for Price')
plt.xlabel('Price (in dollars)')
plt.ylabel('Proportion of Cars')

plt.show()
plt.close()<jupyter_output><empty_output><jupyter_text>We can see that the fitted values are reasonably close to the actual values, since the two distributions overlap a bit. However, there is definitely some room for improvement.  Polynomial Regression and PipelinesPolynomial regression is a particular case of the general linear regression model or multiple linear regression models.We get non-linear relationships by squaring or setting higher-order terms of the predictor variables.There are different orders of polynomial regression:Quadratic - 2nd order$$Yhat = a + b_1 X^2 +b_2 X^2 $$Cubic - 3rd order$$Yhat = a + b_1 X^2 +b_2 X^2 +b_3 X^3\\$$Higher order:$$Y = a + b_1 X^2 +b_2 X^2 +b_3 X^3 ....\\$$We saw earlier that a linear model did not provide the best fit while using highway-mpg as the predictor variable. Let's see if we can try fitting a polynomial model to the data instead.<jupyter_code>def PlotPolly(model, independent_variable, dependent_variable, Name):
    x_new = np.linspace(2, 16, 100)
    y_new = model(x_new)

    plt.plot(independent_variable, dependent_variable, '.', x_new, y_new, '-')
    plt.title('Polynomial Fit with Matplotlib for Price ~ Length')
    ax = plt.gca()
    ax.set_facecolor((0.898, 0.898, 0.898))
    fig = plt.gcf()
    plt.xlabel(Name)
    plt.ylabel('Price of Cars')

    plt.show()
    plt.close()

x = df['highway-L/100km']
y = df['price']

# Here we use a polynomial of the 3rd order (cubic) 
f = np.polyfit(x, y, 3)
p = np.poly1d(f)
print(p)
# Plot the function
PlotPolly(p, x, y, 'highway-L/100km')<jupyter_output><empty_output><jupyter_text>We can perform a polynomial transform on multiple features. First, we import the module:<jupyter_code>from sklearn.preprocessing import PolynomialFeatures
pr=PolynomialFeatures(degree=2)
Z_pr=pr.fit_transform(Z)

print(Z.shape)
print(Z_pr.shape)<jupyter_output>(201, 4)
(201, 15)<jupyter_text>The original data is of 201 samples and 4 features, after the transformation, there 201 samples and 15 features.  PipelineData Pipelines simplify the steps of processing the data. We use the module Pipeline to create a pipeline. We also use StandardScaler as a step in our pipeline.<jupyter_code>from sklearn.pipeline import Pipeline
from sklearn.preprocessing import StandardScaler<jupyter_output><empty_output><jupyter_text>We create the pipeline, by creating a list of tuples including the name of the model or estimator and its corresponding constructor.<jupyter_code>Input=[('scale',StandardScaler()), ('polynomial', PolynomialFeatures(include_bias=False)), ('model',LinearRegression())]
pipe=Pipeline(Input)
pipe.fit(Z,y)
ypipe=pipe.predict(Z)

plt.figure(figsize=(width, height))

ax1 = sns.distplot(df['price'], hist=False, color="r", label="Actual Value")
sns.distplot(ypipe, hist=False, color="b", label="Fitted Values" , ax=ax1)

plt.title('Actual vs Fitted Values for Price')
plt.xlabel('Price (in dollars)')
plt.ylabel('Proportion of Cars')

plt.show()
plt.close()<jupyter_output><empty_output><jupyter_text>We can see that we quite improve the result in comparison with the non polynomial model.  Measures for In-Sample EvaluationTwo very important measures that are often used in Statistics to determine the accuracy of a model are:- R-squared- Mean Squared Error (MSE) R-squaredR squared, also known as the coefficient of determination, is a measure to indicate how close the data is to the fitted regression line.The value of the R-squared is the percentage of variation of the response variable (y) that is explained by a linear model. Mean Squared Error (MSE)The Mean Squared Error measures the average of the squares of errors, that is, the difference between actual value (y) and the estimated value (ŷ). Model 1: Simple Linear Regression<jupyter_code>from sklearn.metrics import mean_squared_error
#Compute the R^2 score for the simple linear regression using highway_L/100km
lm = LinearRegression()

X = df[['highway-L/100km']]
Y = df['price']

lm.fit(X, Y)
print('The R-square is: ', lm.score(X, Y))<jupyter_output>The R-square is:  0.6417894513258818<jupyter_text>We can say that ~ 64.18% of the variation of the price is explained by this simple linear model "highway_L/100km_fit".Let's calculate the MSE<jupyter_code>Yhat=lm.predict(X)
mse = mean_squared_error(df['price'], Yhat)
print('The mean square error of price and predicted value is: ', mse)<jupyter_output>The mean square error of price and predicted value is:  22510543.777085222<jupyter_text>Model 2: Multiple Linear Regression<jupyter_code>Z = df[['horsepower', 'curb-weight', 'engine-size', 'highway-L/100km']]
lm.fit(Z, Y)
print('The R-square is: ', lm.score(Z, Y))

Yhat=lm.predict(Z)
mse = mean_squared_error(Y, Yhat)
print('The mean square error of price and predicted value is: ', mse)<jupyter_output>The R-square is:  0.811811561534475
The mean square error of price and predicted value is:  11826072.956532082<jupyter_text>We can say that ~ 80.896 % of the variation of price is explained by this multiple linear regression "multi_fit". Model 3: Polynomial Fit<jupyter_code>from sklearn.metrics import r2_score

x = df['highway-L/100km']
y = df['price']

# Here we use a polynomial of the 3rd order (cubic) 
f = np.polyfit(x, y, 3)
p = np.poly1d(f)

r_squared = r2_score(y, p(x))
print('The R-square value is: ', r_squared)
print('The mean square error of price and predicted value is: ', mean_squared_error(y, p(x)))<jupyter_output>The R-square value is:  0.6753581660812171
The mean square error of price and predicted value is:  20401030.180019405
{'fraction_non_alphanumeric': 0.10054633326329061, 'fraction_numerical': 0.042130699726833365, 'mean_word_length': 5.41010101010101, 'pattern_counts': {'":': 0, '<': 49, '<?xml version=': 0, '>': 49, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, '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': '14', 'dir': 'jupyter-structured-clean-dedup', 'id': '121255', 'license': 'MIT', 'n_tokens_mistral': 3386, 'n_tokens_neox': 3040, 'n_words': 1122, 'path': 'IBM_DA0101EN/Model Development.ipynb', 'repo_name': 'QuentinElGuay/notebooks', 'size': '9518'}