# Agenda

- Introduction to regression analysis

    - Dependent variable, independent variables
    - Prediction, forecasting, inference
    - Types of regression analysis
    

- Introduction to linear regression

    - Discussion of Continuous, ordinal and nominal data
    - Simple linear regression
    - Formalization of the problem
    - Objective Function
    - Using Python for linear regression
    - Linear regression with more than one predictors
    

# Regression Analysis


Regression analysis is a __set of statistical processes__ for estimating the relationship between a __dependent variable__ and one or more __independent variables__.

- Goal: Understand the distribution of the response variable (Y) conditioned on predictors (X).

- Can be used for prediction or forecasting or inference.


__Prediction__

- Make educated guesses about a population (whole) from data (sample)

__Forecast__

- Make systematic and educated guesses about future from knowledge gain from past

__Inference__

Is there any casual relationship between a dependent variable and independent variables in a fix dataset.

Source: Statistical Learning From Regression Perspective 2nd edition p3

## Types of Regression Analysis

- Linear Regression

- Binomial Regression

- Poisson Regression

- Robust Regression

- Time Series 

- And many other variations

# Linear Regression

## Measurement Scales

__Nominal Scales (Categorical Data)__

- They are labels more than scales

- Examples: Male, Female, Republican, Democrat, Independent

- Note that they can show up as numerical values too:{"Math":1, "Economics": 2, "Philosophy": 3}.

__Ordinal Scales__

- Orders - Rank

- Military: commander< captain < admiral

- Another example: Pain Level: 0 < 1 < 2 < 3 < 4 < ... < 10

__Interval Scales__

- Differences makes sense but ratios don't!

- Example: 10&deg;F, 20&deg;F, 40&deg;F, 80&deg;F

__Ratio Scales__

- There is an absolute zero point.

- 5 seconds 10 seconds

- 5 meters 10 meters

- 5 lb, 10 lb

## Simple Linear Regression

__Understanding the task__

For this part of the lecture, we will follow the textbook [Introduction to Statistical Learning with R, 7th edition.](https://faculty.marshall.usc.edu/gareth-james/ISL/ISLR%20Seventh%20Printing.pdf) (ISLR)


```python
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline

```


```python
## load data is also supporting code
import load_data
```


```python
# read data to a dataframe
data = load_data.data()[['TV', 'Sales']]

print(data.shape)
data.head()
```

    (200, 2)





<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>TV</th>
      <th>Sales</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>1</th>
      <td>230.1</td>
      <td>22.1</td>
    </tr>
    <tr>
      <th>2</th>
      <td>44.5</td>
      <td>10.4</td>
    </tr>
    <tr>
      <th>3</th>
      <td>17.2</td>
      <td>9.3</td>
    </tr>
    <tr>
      <th>4</th>
      <td>151.5</td>
      <td>18.5</td>
    </tr>
    <tr>
      <th>5</th>
      <td>180.8</td>
      <td>12.9</td>
    </tr>
  </tbody>
</table>
</div>



__About this dataset:__ The dataset contains statistics about the sales of a product in 200 different markets, together with advertising budgets in each of these markets for different media channels: TV, radio and newspaper.

The sales are in thousands of units and the budget is in thousands of dollars.

[For more informatio about this dataset](https://cran.r-project.org/web/packages/ISLR/ISLR.pdf)

__Some interesting questions we can try to answer__

- Is there any relationship between advertising budget and sales?

- If there is any relation how strong is the relationship?

- Which media contribute to sales?

For more read ISLR p60.


```python
## plot data
data.plot(kind = 'scatter', x = 'TV', y= 'Sales')
plt.show()
```


```python
## utils is a py file that I added a homemade functions
from utils import draw_line

tv = data.TV.values
sales = data.Sales.values

```


```python
## draw line will draw a line with a given slope and intercept
draw_line(tv, sales, slope= 0.07, intercept = 6.8, figsize = (10,6) )
```



__Formalization of the Process__

- In simple linear regression we assume that if we would have population of both X and y variables then we would see the following relation:

$$ Y = \beta_{0} + \beta_{1}X + \epsilon$$
 
- $ \beta_{0}, \beta_{1}$ are parameters of the model and called the __intercept__ and __coefficient__ of the linear model respectively. 

- $\epsilon$ is the irreducible error term. Depend on the problem at hand we might assume that these errors are coming from measurement mistakes, personal beliefs, recording errors, etc.

__Note__

The same equation above can appear in different notations. We choose to stick to the notation in ISLR.

Some alternatives: 

$$ Y = \theta_{0} + \theta_{1}X + \epsilon$$
 
 
$$ Y = w_{0} + w_{1}X + \epsilon$$


## Assumptions

* Assume $\mbox{Sales} \approx \mbox{TV}$  linearly.

* There is an initial value (baseline) of the Sales. 

* Data might not be fit on a line but the errors are random.

* Observations(rows) are independent of each other

__Task?__

Let's discuss

__Task?__


- Finding $\beta_{0}$, $\beta_{1}$ that gives us the best line!

Is research problem specific enough?

__Task?__


- Finding $\beta_{0}$, $\beta_{1}$ that gives the line which minimizes the prediction error!

Is research problem specific enough?

__Task?__


- Finding $\beta_{0}$, $\beta_{1}$ that gives the line which minimizes the sum of individual errors!

Is research problem specific enough?

__Task?__


- Finding $\beta_{0}$, $\beta_{1}$ that gives the line which minimizes the sum of squares of individual errors!

This approach is called __Least Squares Approach__

__The same statement just fancier__

__Formula (Model)__

$$ Y = \beta_{0} + \beta_{1}X $$


__Objective Function (Cost Function)__

$$ J(\beta_{0}, \beta_{1}) = \frac{1}{2N}\sum_{i = 1}^{N} \left( e_{i} \right)^{2} $$

__Goal:  minimize $J(\beta_{0}, \beta_{1})$ with respect to $\beta_{0}, \beta_{1}$__

$$ e_i = y_{i} - (\beta_{0} + \beta_{1}x_{i})$$



__Full Picture__




```python
import matplotlib.image as mpimg 

f, axes = plt.subplots(nrows= 1, ncols = 3)
for i in range(3):
    axes[i].axis('off')
    if i==0:
        img = mpimg.imread(f'images/simple{i}.jpg')
    else:
        img = mpimg.imread(f'images/simple{i}.png')
    axes[i].imshow(img)
plt.gcf().set_size_inches(20, 20) 
```

__Solution__

$$ \boldsymbol{\beta} = (X^{T}X)^{-1}X^{T}\boldsymbol{Y}$$


__Finding $\beta$'s with Python__


```python
## Equation 3.4 from ISLR

def least_squares(X,y):
    x_bar = X.mean()
    y_bar = y.mean()
    b1_num = (X - x_bar).dot(y - y_bar)
    b1_denom = ((X - x_bar)**2).sum()
    b1 = b1_num /b1_denom
    b0 = y_bar - b1*x_bar
    J = 0.5 *((y - (b0 + b1*X))**2).mean()
    return  b0,b1, J

least_squares(data.TV, data.Sales)
```




    (7.032593549127698, 0.047536640433019736, 5.256326457828378)




```python
## Four Step workflow with sklearn
## 1. import the tool
from sklearn.linear_model import LinearRegression

## 2. instantiate the tool you would like to use.
lr = LinearRegression()
X = data[['TV']]
y = data.Sales
## 3. fit the model
lr.fit(X, y)

## 4. check-out the results or make predictions
print(lr.intercept_, lr.coef_)
```

    7.032593549127693 [0.04753664]



```python
## R_squared Score

lr.score(X, y)
```




    0.611875050850071



$R^{2}:$ The proportion of the variability in $Y$ that can be explained using X with the linear model



## More than one Variables

- When we have more than one variable our linear model will look like as:

$$ Y = \beta_{0} + \beta_{1} X_{1} + \beta_{2}X_{2} + \cdots \beta_{p}X_{p} + \epsilon$$

- Cost function (objective function) become:

\begin{align}
RSS &= \sum\limits_{i=1}^{N} (y_{i} - \hat{y}_{i})^{2} \\
&= \sum\limits_{i=1}^{N} (y_{i} - \beta_{0} - \beta_{1} x_{i1}  - \beta_{2}x_{i2} - \cdots  - \beta_{p}x_{ip})^{2}
\end{align}

__The same statement just fancier__

__Formula (Model)__

$$ Y = \beta_{0} + \beta_{1} X_{1} + \beta_{2}X_{2} + \cdots \beta_{p}X_{p} + \epsilon$$

__Objective Function (Cost Function)__

$$ J(\beta_{0}, \beta_{1}, \cdots, \beta_{p}) = \frac{1}{2N}\sum_{i = 1}^{N} \left( e_{i} \right)^{2} $$

__Goal:  minimize $J(\beta_{0}, \beta_{1}, \cdots, \beta_{p})$ with respect to $\beta_{0}, \beta_{1}, \cdots, \beta_{p}$__


```python
## This time we will use statsmodels api

import statsmodels.api as sm
data = load_data.data()
X = data[['TV', 'Radio', 'Newspaper']]
y = data.Sales
X = sm.add_constant(X)
mod = sm.OLS(y, X, hasconst= True)
res = mod.fit()
res.summary()
```




<table class="simpletable">
<caption>OLS Regression Results</caption>
<tr>
  <th>Dep. Variable:</th>          <td>Sales</td>      <th>  R-squared:         </th> <td>   0.897</td>
</tr>
<tr>
  <th>Model:</th>                   <td>OLS</td>       <th>  Adj. R-squared:    </th> <td>   0.896</td>
</tr>
<tr>
  <th>Method:</th>             <td>Least Squares</td>  <th>  F-statistic:       </th> <td>   570.3</td>
</tr>
<tr>
  <th>Date:</th>             <td>Thu, 01 Oct 2020</td> <th>  Prob (F-statistic):</th> <td>1.58e-96</td>
</tr>
<tr>
  <th>Time:</th>                 <td>18:41:41</td>     <th>  Log-Likelihood:    </th> <td> -386.18</td>
</tr>
<tr>
  <th>No. Observations:</th>      <td>   200</td>      <th>  AIC:               </th> <td>   780.4</td>
</tr>
<tr>
  <th>Df Residuals:</th>          <td>   196</td>      <th>  BIC:               </th> <td>   793.6</td>
</tr>
<tr>
  <th>Df Model:</th>              <td>     3</td>      <th>                     </th>     <td> </td>   
</tr>
<tr>
  <th>Covariance Type:</th>      <td>nonrobust</td>    <th>                     </th>     <td> </td>   
</tr>
</table>
<table class="simpletable">
<tr>
      <td></td>         <th>coef</th>     <th>std err</th>      <th>t</th>      <th>P>|t|</th>  <th>[0.025</th>    <th>0.975]</th>  
</tr>
<tr>
  <th>const</th>     <td>    2.9389</td> <td>    0.312</td> <td>    9.422</td> <td> 0.000</td> <td>    2.324</td> <td>    3.554</td>
</tr>
<tr>
  <th>TV</th>        <td>    0.0458</td> <td>    0.001</td> <td>   32.809</td> <td> 0.000</td> <td>    0.043</td> <td>    0.049</td>
</tr>
<tr>
  <th>Radio</th>     <td>    0.1885</td> <td>    0.009</td> <td>   21.893</td> <td> 0.000</td> <td>    0.172</td> <td>    0.206</td>
</tr>
<tr>
  <th>Newspaper</th> <td>   -0.0010</td> <td>    0.006</td> <td>   -0.177</td> <td> 0.860</td> <td>   -0.013</td> <td>    0.011</td>
</tr>
</table>
<table class="simpletable">
<tr>
  <th>Omnibus:</th>       <td>60.414</td> <th>  Durbin-Watson:     </th> <td>   2.084</td>
</tr>
<tr>
  <th>Prob(Omnibus):</th> <td> 0.000</td> <th>  Jarque-Bera (JB):  </th> <td> 151.241</td>
</tr>
<tr>
  <th>Skew:</th>          <td>-1.327</td> <th>  Prob(JB):          </th> <td>1.44e-33</td>
</tr>
<tr>
  <th>Kurtosis:</th>      <td> 6.332</td> <th>  Cond. No.          </th> <td>    454.</td>
</tr>
</table><br/><br/>Warnings:<br/>[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.



__Q__ Is there a relationship between predictors and response?

__Answer__ Check the F-stats

__Q__ Which factors are important?

This is called variable selection and not a very straightforward process.

- Ideally we would check all possible prediction combinations and then compare them with AIC, BIC, $C_P$

__Q__: Predictions


```python
print(res.predict(X.loc[13].values), y.loc[13])

X.loc[13]
```

    [10.57712073] 9.2





    const         1.0
    TV           23.8
    Radio        35.1
    Newspaper    65.9
    Name: 13, dtype: float64




```python
residuals = res.resid

y_pred = res.predict(X)
plt.scatter(y_pred, residuals)
plt.axhline(y=0, color = 'red', label = '0')

plt.xlabel('predicted values')

plt.ylabel('residuals')

plt.tight_layout()
```

__Q__ Model Fit

- We can check $R^{2}$ or adjusted $R^{2}$

######  Play Ground


```python
import numpy as np
ncols = 14
nrows = 15
X = np.zeros(shape = (nrows,ncols))
for i in range(ncols):
    X[:,i] = np.random.normal(size = nrows)
    
y = np.random.normal(size = nrows)
```


```python
import statsmodels.api as sm

mod = sm.OLS(y, X, hasconst= True)
res = mod.fit()
print(res.summary())
```

                                OLS Regression Results                            
    ==============================================================================
    Dep. Variable:                      y   R-squared:                       1.000
    Model:                            OLS   Adj. R-squared:                  0.997
    Method:                 Least Squares   F-statistic:                     420.4
    Date:                Thu, 01 Oct 2020   Prob (F-statistic):             0.0382
    Time:                        18:59:08   Log-Likelihood:                 44.707
    No. Observations:                  15   AIC:                            -61.41
    Df Residuals:                       1   BIC:                            -51.50
    Df Model:                          13                                         
    Covariance Type:            nonrobust                                         
    ==============================================================================
                     coef    std err          t      P>|t|      [0.025      0.975]
    ------------------------------------------------------------------------------
    x1             0.4817      0.018     27.328      0.023       0.258       0.706
    x2            -0.9031      0.112     -8.093      0.078      -2.321       0.515
    x3            -0.4771      0.055     -8.682      0.073      -1.175       0.221
    x4             0.1321      0.042      3.158      0.195      -0.399       0.663
    x5            -0.8227      0.049    -16.873      0.038      -1.442      -0.203
    x6            -0.6451      0.057    -11.380      0.056      -1.366       0.075
    x7            -0.1962      0.130     -1.515      0.371      -1.842       1.450
    x8             0.6703      0.024     27.698      0.023       0.363       0.978
    x9             0.5629      0.068      8.334      0.076      -0.295       1.421
    x10            0.1525      0.115      1.330      0.411      -1.305       1.610
    x11            1.0894      0.140      7.800      0.081      -0.685       2.864
    x12            0.2359      0.028      8.485      0.075      -0.117       0.589
    x13            0.8089      0.067     12.014      0.053      -0.047       1.664
    x14           -0.1957      0.049     -4.003      0.156      -0.817       0.425
    ==============================================================================
    Omnibus:                        1.174   Durbin-Watson:                   2.182
    Prob(Omnibus):                  0.556   Jarque-Bera (JB):                0.798
    Skew:                          -0.533   Prob(JB):                        0.671
    Kurtosis:                       2.623   Cond. No.                         37.5
    ==============================================================================
    
    Warnings:
    [1] Standard Errors assume that the covariance matrix of the errors is correctly specified.


__Important Note__

- As we mentioned in the beginning of the lecture linear regression is a basic approach but not necessarily simple one. 

- Especially when the assumptions do not hold things results might not be reliable. 

- For more details on this please read ISLR 3.3.3 very carefully.

## Further Reading 

- [What is the parameter copy_X in sklearn](https://stackoverflow.com/questions/52541717/what-is-the-parameter-copy-x-in-sklearn-linear-regression)

- https://newonlinecourses.science.psu.edu/stat501/node/252/

- https://towardsdatascience.com/simple-and-multiple-linear-regression-in-python-c928425168f9

- https://www.statsmodels.org/stable/examples/notebooks/generated/ols.html

- https://lagunita.stanford.edu/assets/courseware/v1/e4f2f4823252ead8c1a68722e3415c28/c4x/HumanitiesSciences/StatLearning/asset/linear_regression-handout.pdf

- [Statsmodels API](https://www.statsmodels.org/stable/discretemod.html)

"Some of the figures in this presentation are taken from "An Introduction to Statistical Learning, with applications in R"  (Springer, 2013) with permission from the authors: G. James, D. Witten,  T. Hastie and R. Tibshirani " 

d705a854b63d93edd80e1f98696002b8797513bb
ipynb
Jupyter Notebook
week-6/Regression.ipynb
jaredfincke/UMBC_Data601
0507d287391f8e857f2f4aa52e0131f4f893879a
MIT
2020-09-03T20:35:50.000Z
2020-11-05T00:57:25.000Z
week-6/Regression.ipynb
jaredfincke/UMBC_Data601
0507d287391f8e857f2f4aa52e0131f4f893879a
MIT
2020-10-15T20:36:36.000Z
2020-10-22T02:31:45.000Z
week-6/Regression.ipynb
jaredfincke/UMBC_Data601
0507d287391f8e857f2f4aa52e0131f4f893879a
MIT
2020-09-03T20:34:39.000Z
2021-01-28T16:05:56.000Z
Qwen/Qwen-72B
1. YES
2. YES