<jupyter_start><jupyter_text>_Lambda School Data Science_ Ordinary Least Squares Regression What is Linear Regression?Linear Regression is a statistical model that seeks to describe the relationship between some y variable and one or more x variables. ![Linear Regression](https://upload.wikimedia.org/wikipedia/commons/thumb/3/3a/Linear_regression.svg/1200px-Linear_regression.svg.png)In the simplest case, linear regression seeks to fit a straight line through a cloud of points. This line is referred to as the "regression line" or "line of best fit." This line tries to summarize the relationship between our X and Y in a way that enables us to use the equation for that line to make predictions. Synonyms for "y variable" - Dependent Variable- Response Variable- Outcome Variable - Predicted Variable- Measured Variable- Explained Variable- Label- Target Synonyms for "x variable"- Independent Variable- Explanatory Variable- Regressor- Covariate- Feature  Simple Linear Regresion (bivariate) Making PredictionsSay that we were trying to create a model that captured the relationship between temperature outside and ice cream sales. In Machine Learning our goal is often different that of other flavors of Linear Regression Analysis, because we're trying to fit a model to this data with the intention of making **predictions** on new data (in the future) that we don't have yet. What are we trying to predict?So if we had measured ice cream sales and the temprature outside on 11 different days, at the end of our modeling **what would be the thing that we would want to predict? - Ice Cream Sales or Temperature?**We would probably want to be measuring temperature with the intention of using that to **forecast** ice cream sales. If we were able to successfully forecast ice cream sales from temperature, this might help us know beforehand how much ice cream to make or how many cones to buy or on which days to open our store, etc. Being able to make predictions accurately has a lot of business implications. This is why making accurate predictions is so valuable (And in large part is why data scientists are paid so well). Y Variable IntuitionWe want the thing that we're trying to predict to serve as our **y** variable. This is why it's sometimes called the "predicted variable." We call it the "dependent" variable because our prediction for how much ice cream we're going to sell "depends" on the temperature outside.  X Variable IntuitionAll other variables that we use to predict our y variable (we're going to start off just using one) we call our **x** variables. These are called our "independent" variables because they don't *depend* on y, they "explain" y. Hence they are also referred to as our "explanatory" variables.<jupyter_code>%matplotlib inline
from ipywidgets import interact
from matplotlib.patches import Rectangle
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import seaborn as sns
from sklearn.linear_model import LinearRegression
from sklearn.metrics import mean_absolute_error, mean_squared_error, r2_score
import statsmodels.api as sm
columns = ['Year','Incumbent Party Candidate','Other Candidate','Incumbent Party Vote Share']

data = [[1952,"Stevenson","Eisenhower",44.6],
        [1956,"Eisenhower","Stevenson",57.76],
        [1960,"Nixon","Kennedy",49.91],
        [1964,"Johnson","Goldwater",61.34],
        [1968,"Humphrey","Nixon",49.60],
        [1972,"Nixon","McGovern",61.79],
        [1976,"Ford","Carter",48.95],
        [1980,"Carter","Reagan",44.70],
        [1984,"Reagan","Mondale",59.17],
        [1988,"Bush, Sr.","Dukakis",53.94],
        [1992,"Bush, Sr.","Clinton",46.55],
        [1996,"Clinton","Dole",54.74],
        [2000,"Gore","Bush, Jr.",50.27],
        [2004,"Bush, Jr.","Kerry",51.24],
        [2008,"McCain","Obama",46.32],
        [2012,"Obama","Romney",52.00]]
        
df = pd.DataFrame(data=data, columns=columns)
df
df.plot(x='Year', y='Incumbent Party Vote Share', kind='scatter');
df['Incumbent Party Vote Share'].describe()
target = 'Incumbent Party Vote Share'
df['Prediction'] = df[target].mean()
df['Error'] = df['Prediction'] - df[target]
df['Error'].sum()
df['Absolute Error'] = df['Error'].abs()
df['Absolute Error'].sum()
df['Absolute Error'].mean()
df
mean_absolute_error(y_true=df[target], y_pred=df['Prediction'])
r2_score(y_true=df[target], y_pred=df['Prediction'])<jupyter_output><empty_output><jupyter_text>R Squared:  $R^2$One final attribute of linear regressions that we're going to talk about today is a measure of goodness of fit known as $R^2$ or R-squared. $R^2$ is a statistical measure of how close the data are fitted to our regression line. A helpful interpretation for the $R^2$ is the percentage of the dependent variable that is explained by the model.In other words, the $R^2$ is the percentage of y that is explained by the x variables included in the model. For this reason the $R^2$ is also known as the "coefficient of determination," because it explains how much of y is explained (or determined) by our x varaibles. We won't go into the calculation of $R^2$ today, just know that a higher $R^2$ percentage is nearly always better and indicates a model that fits the data more closely.   Add data<jupyter_code>columns = ['Year','Average Recent Growth in Personal Incomes']

data = [[1952,2.40],
        [1956,2.89],
        [1960, .85],
        [1964,4.21],
        [1968,3.02],
        [1972,3.62],
        [1976,1.08],
        [1980,-.39],
        [1984,3.86],
        [1988,2.27],
        [1992, .38],
        [1996,1.04],
        [2000,2.36],
        [2004,1.72],
        [2008, .10],
        [2012, .95]]
        
growth = pd.DataFrame(data=data, columns=columns)
df = df.merge(growth)
df
feature = 'Average Recent Growth in Personal Incomes'
df.plot(x=feature, y=target, kind='scatter');<jupyter_output><empty_output><jupyter_text>We can see from the scatterplot that these data points seem to follow a somewhat linear relationship. This means that we could probably summarize their relationship well by fitting a line of best fit to these points. Lets do it. The Equation for a LineAs we know a common equation for a line is:\begin{align}y = mx + b\end{align}Where $m$ is the slope of our line and $b$ is the y-intercept. If we want to plot a line through our cloud of points we figure out what these two values should be. Linear Regression seeks to **estimate** the slope and intercept values that describe a line that best fits the data points.<jupyter_code>m = 4
b = 44
df['Prediction'] = m * df[feature] + b
df['Error'] = df['Prediction'] - df[target]
df['Absolute Error'] = df['Error'].abs()
df['Absolute Error'].sum()
df['Absolute Error'].mean()
r2_score(y_true=df[target], y_pred=df['Prediction'])
ax = df.plot(x=feature, y=target, kind='scatter')
df.plot(x=feature, y='Prediction', kind='line', ax=ax);
def regression(m, b):
    df['Prediction'] = m * df[feature] + b
    df['Error'] = df['Prediction'] - df[target]
    df['Absolute Error'] = df['Error'].abs()
    sum_absolute_error = df['Absolute Error'].sum()
    
    title = f'Sum of absolute errors: {sum_absolute_error}'
    ax = df.plot(x=feature, y=target, kind='scatter', title=title, figsize=(7, 7))
    df.plot(x=feature, y='Prediction', kind='line', ax=ax)
    
regression(m=4, b=48)<jupyter_output><empty_output><jupyter_text>Residual Error The residual error is the distance between points in our dataset and our regression line.<jupyter_code>def regression(m, b):
    df['Prediction'] = m * df[feature] + b
    df['Error'] = df['Prediction'] - df[target]
    df['Absolute Error'] = df['Error'].abs()
    sum_absolute_error = df['Absolute Error'].sum()
    
    title = f'Sum of absolute errors: {sum_absolute_error}'
    ax = df.plot(x=feature, y=target, kind='scatter', title=title, figsize=(7, 7))
    df.plot(x=feature, y='Prediction', kind='line', ax=ax)
    
    for x, y1, y2 in zip(df[feature], df[target], df['Prediction']):
        ax.plot((x, x), (y1, y2), color='grey')

regression(m=3, b=46)
interact(regression, m=(-10,10,0.5), b=(40,60,0.5));
df['Square Error'] = df['Error'] **2
def regression(m, b):
    df['Prediction'] = m * df[feature] + b
    df['Error'] = df['Prediction'] - df[target]
    df['Absolute Error'] = df['Error'].abs()
    df['Square Error'] = df['Error'] **2
    sum_square_error = df['Square Error'].sum()
    
    title = f'Sum of square errors: {sum_square_error}'
    ax = df.plot(x=feature, y=target, kind='scatter', title=title, figsize=(7, 7))
    df.plot(x=feature, y='Prediction', kind='line', ax=ax)
        
    xmin, xmax = ax.get_xlim()
    ymin, ymax = ax.get_ylim()
    scale = (xmax-xmin)/(ymax-ymin)

    for x, y1, y2 in zip(df[feature], df[target], df['Prediction']):
        bottom_left = (x, min(y1, y2))
        height = abs(y1 - y2)
        width = height * scale
        ax.add_patch(Rectangle(xy=bottom_left, width=width, height=height, alpha=0.1))
interact(regression, m=(-10,10,0.5), b=(40,60,0.5));
b = 46
ms = np.arange(-10,10,0.5)
sses = []

for m in ms:
    predictions = m * df[feature] + b
    errors = predictions - df[target]
    square_errors = errors ** 2
    sse = square_errors.sum()
    sses.append(sse)
    
hypotheses = pd.DataFrame({'Slope': ms})
hypotheses['Intercept'] = b
hypotheses['Sum of Square Errors'] = sses

hypotheses.sort_values(by='Sum of Square Errors')
hypotheses.plot(x='Slope', y='Sum of Square Errors', 
                title=f'Intercept={b}');
X = df[[feature]]
y = df[target]

X.shape, y.shape
model = LinearRegression()
model.fit(X, y)
model.coef_, model.intercept_
model.predict([[0]])
model.predict([[1]])
model.predict([[2]])
model.predict([[3]])
model.predict(X)
df['Prediction'] = model.predict(X)
ax = df.plot(x=feature, y=target, kind='scatter', title='sklearn LinearRegression')
df.plot(x=feature, y='Prediction', kind='line', ax=ax);
df['Error'] = df['Prediction'] - y
df['Absolute Error'] = df['Error'].abs()
df['Square Error'] = df['Error'] ** 2
df['Square Error'].mean()
mean_squared_error(y_true=y, y_pred=model.predict(X))
np.sqrt(mean_squared_error(y, model.predict(X)))
model.score(X, y)
r2_score(y, model.predict(X))<jupyter_output><empty_output><jupyter_text>Statsmodels https://www.statsmodels.org/dev/examples/notebooks/generated/ols.html<jupyter_code>X
sm.add_constant(X)
model = sm.OLS(y, sm.add_constant(X))
results = model.fit()
print(results.summary())<jupyter_output>OLS Regression Results                                
======================================================================================
Dep. Variable:     Incumbent Party Vote Share   R-squared:                       0.580
Model:                                    OLS   Adj. R-squared:                  0.550
Method:                         Least Squares   F-statistic:                     19.32
Date:                        Mon, 29 Apr 2019   Prob (F-statistic):           0.000610
Time:                                12:42:56   Log-Likelihood:                -42.839
No. Observations:                          16   AIC:                             89.68
Df Residuals:                              14   BIC:                             91.22
Df Model:                                   1                                         
Covariance Type:                    nonrobust                                         
===========================================[...]<jupyter_text>The Anatomy of Linear Regression- Intercept: The $b$ value in our line equation $y=mx+b$- Slope: The $m$ value in our line equation $y=mx+b$. These two values together define our regression line.![Slope and Intercept](http://www.ryanleeallred.com/wp-content/uploads/2018/08/linear-regression-diagram.png)- $\hat{y}$ : A prediction- Line of Best Fit (Regression Line)- Predicted (fitted) Values: Points on our regression line- Observed Values: Points from our dataset- Error: The distance between predicted and observed values.![Residual Error](http://www.ryanleeallred.com/wp-content/uploads/2018/08/residual-or-error.gif)  More Formal Notation![Simple Linear Regression](http://www.ryanleeallred.com/wp-content/uploads/2018/08/simple-regression-formula.png)We have talked about a line of regression being represented like a regular line $y=mx+b$ but as we get to more complicated versions we're going to need to extend this equation. So lets establish the proper terminology.**X** - Independent Variable, predictor variable, explanatory variable, regressor, covariate**Y** - Response variable, predicted variable, measured vairable, explained variable, outcome variable$\beta_0$ - "Beta Naught" or "Beta Zero", the intercept value. This is how much of y would exist if X were zero. This is sometimes represented by the letter "a" but I hate that. So it's "Beta 0" during my lecture.$\beta_1$ - "Beta One" The primary coefficient of interest. This values is the slope of the line that is estimated by "minimizing the sum of the squared errors/residuals" - We'll get to that. $\epsilon$ - "Epsilon" The "error term", random noise, things outside of our model that affect y.  How Does it do it? Minimizing the Sum of the Squared ErrorThe most common method of estimating our $\beta$ parameters  is what's known as "Ordinary Least Squares" (OLS). (There are different methods of arriving at a line of best fit). OLS estimates the parameters that minimize the squared distance between each point in our dataset and our line of best fit. \begin{align}SSE = \sum(y_i - \hat{y})^2\end{align}  Linear Algebra!The same result that is found by minimizing the sum of the squared errors can be also found through a linear algebra process known as the "Least Squares Solution:"![OLS Regression](http://www.ryanleeallred.com/wp-content/uploads/2018/08/OLS-linear-algebra.png)Before we can work with this equation in its linear algebra form we have to understand how to set up the matrices that are involved in this equation.   The $\beta$ vectorThe $\beta$ vector represents all the parameters that we are trying to estimate, our $y$ vector and $X$ matrix values are full of data from our dataset. The $\beta$ vector holds the variables that we are solving for: $\beta_0$ and $\beta_1$Now that we have all of the necessary parts we can set them up in the following equation:\begin{align}y = X \beta + \epsilon\end{align}Since our $\epsilon$ value represents **random** error we can assume that it will equal zero on average.\begin{align}y = X \beta\end{align}The objective now is to isolate the $\beta$ matrix. We can do this by pre-multiplying both sides by "X transpose" $X^{T}$.\begin{align}X^{T}y =  X^{T}X \beta\end{align}Since anything times its transpose will result in a square matrix, if that matrix is then an invertible matrix, then we should be able to multiply both sides by its inverse to remove it from the right hand side. (We'll talk tomorrow about situations that could lead to $X^{T}X$ not being invertible.)\begin{align}(X^{T}X)^{-1}X^{T}y =  (X^{T}X)^{-1}X^{T}X \beta\end{align}Since any matrix multiplied by its inverse results in the identity matrix, and anything multiplied by the identity matrix is itself, we are left with only $\beta$ on the right hand side:\begin{align}(X^{T}X)^{-1}X^{T}y = \hat{\beta}\end{align}We will now call it "beta hat" $\hat{\beta}$ because it now represents our estimated values for $\beta_0$ and $\beta_1$ Lets calculate our $\beta$ coefficients with numpy!<jupyter_code>X = sm.add_constant(df[feature]).values
print('X')
print(X)

y = df[target].values[:, np.newaxis]
print('y')
print(y)

X_transpose = X.T
print('X Transpose')
print(X_transpose)

X_transpose_X = X_transpose @ X
print('X Transpose X')
print(X_transpose_X)

X_transpose_X_inverse = np.linalg.inv(X_transpose_X)
print('X Transpose X Inverse')
print(X_transpose_X_inverse)

X_transpose_y = X_transpose @ y
print('X Transpose y')
print(X_transpose_y)

beta_hat = X_transpose_X_inverse @ X_transpose_y
print('Beta Hat')
print(beta_hat)<jupyter_output>X
[[ 1.    2.4 ]
 [ 1.    2.89]
 [ 1.    0.85]
 [ 1.    4.21]
 [ 1.    3.02]
 [ 1.    3.62]
 [ 1.    1.08]
 [ 1.   -0.39]
 [ 1.    3.86]
 [ 1.    2.27]
 [ 1.    0.38]
 [ 1.    1.04]
 [ 1.    2.36]
 [ 1.    1.72]
 [ 1.    0.1 ]
 [ 1.    0.95]]
y
[[44.6 ]
 [57.76]
 [49.91]
 [61.34]
 [49.6 ]
 [61.79]
 [48.95]
 [44.7 ]
 [59.17]
 [53.94]
 [46.55]
 [54.74]
 [50.27]
 [51.24]
 [46.32]
 [52.  ]]
X Transpose
[[ 1.    1.    1.    1.    1.    1.    1.    1.    1.    1.    1.    1.
   1.    1.    1.    1.  ]
 [ 2.4   2.89  0.85  4.21  3.02  3.62  1.08 -0.39  3.86  2.27  0.38  1.04
   2.36  1.72  0.1   0.95]]
X Transpose X
[[16.    30.36 ]
 [30.36  86.821]]
X Transpose X Inverse
[[ 0.18575056 -0.06495418]
 [-0.06495418  0.03423145]]
X Transpose y
[[ 832.88  ]
 [1669.7967]]
Beta Hat
[[46.24764802]
 [ 3.06052805]]<jupyter_text>Multiple RegressionSimple or bivariate linear regression involves a single $x$ variable and a single $y$ variable. However, we can have many $x$ variables. A linear regression model that involves multiple x variables is known as **Multiple** Regression - NOT MULTIVARIATE!![Multiple Regression](http://www.ryanleeallred.com/wp-content/uploads/2018/08/multiple-regression-model.png)<jupyter_code>df.sort_values(by='Square Error', ascending=False)
"""
Fatalities denotes the cumulative number of American military
fatalities per millions of US population the in Korea, Vietnam,
Iraq and Afghanistan wars during the presidential terms
preceding the 1952, 1964, 1968, 1976 and 2004, 2008 and
2012 elections.

http://www.douglas-hibbs.com/HibbsArticles/HIBBS-PRESVOTE-SLIDES-MELBOURNE-Part1-2014-02-26.pdf
"""

columns = ['Year','US Military Fatalities per Million']

data = [[1952,190],
        [1956,  0],
        [1960,  0],
        [1964,  1],
        [1968,146],
        [1972,  0],
        [1976,  2],
        [1980,  0],
        [1984,  0],
        [1988,  0],
        [1992,  0],
        [1996,  0],
        [2000,  0],
        [2004,  4],
        [2008, 14],
        [2012,  5]]
        
deaths = pd.DataFrame(data=data, columns=columns)
df = df.merge(deaths)
features = ['Average Recent Growth in Personal Incomes', 
            'US Military Fatalities per Million']

target = 'Incumbent Party Vote Share'

X = df[features]
y = df[target]

model = sm.OLS(y, sm.add_constant(X))
results = model.fit()
print(results.summary())
model = LinearRegression()
model.fit(X, y)
print('Intercept:', model.intercept_)
pd.Series(model.coef_, features)
np.sqrt(mean_squared_error(y, model.predict(X)))<jupyter_output><empty_output><jupyter_text>Train / Test Split<jupyter_code>train = df.query('Year < 2008')
test  = df.query('Year >= 2008')
X_train = train[features]
y_train = train[target]
X_test  = test[features]
y_test  = test[target]
X_train.shape, y_train.shape, X_test.shape, y_test.shape
model.fit(X_train, y_train)
model.predict(X_test)
test<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.13315649867374005, 'fraction_numerical': 0.04281167108753316, 'mean_word_length': 3.6741879494173073, 'pattern_counts': {'":': 0, '<': 32, '<?xml version=': 0, '>': 32, 'https://': 2, 'lorem ipsum': 0, 'www.': 7, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 9, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '8', 'dir': 'jupyter-structured-clean-dedup', 'id': '43790', 'license': 'MIT', 'n_tokens_mistral': 6536, 'n_tokens_neox': 5893, 'n_words': 2041, 'path': 'module1-ols-regression/ols-regression.ipynb', 'repo_name': 'biovir3/DS-Unit-2-Sprint-2-Regression', 'size': '18850'}