DSC365: Introduction to Data Science

class: center, middle, inverse, title-slide

.title[
# DSC365: Introduction to Data Science
]
.author[
### Linear Models
]
.date[
### March 3, 2026
]

---

### Announcements

**Lab 4**

- Due Thursday March 5th at 11:59 pm in Blueline

**Quiz 2**

- Last 30 minutes of class Thursday March 5th 
- Topics: Data Communication and Statistical Foundations/Linear Models

**Lab 5**: work day in class on Thursday March 5th

- Due Tuesday March 24th at 11:59 pm in Blueline

**Mini Project 2**

- Due Thursday March 19th at 11:59 in Blueline

---
### Language of Models

Sometimes we want to know whether any variables are related, thus we need to fit a statistical model.

- We will focus on linear models, but many other types of models too!

.pull-left[
<img src="linear-models_files/figure-html/unnamed-chunk-2-1.png" alt="An example of a scatterplot that has a linear relationship between the two variables." />
].pull-right[
<img src="linear-models_files/figure-html/unnamed-chunk-3-1.png" alt="An example of a scatterplot that has a non-linear relationship between the two variables. The relationship looks like a parabola." />

]

---
### Fit a model

Sometimes we want to know whether any variables are related, thus we need to fit a statistical model.

__Statistical models__:

There are two things we can do with fitting a model:

---
### Fit a model

Today's lecture (and mini-project 2) focuses on interpretation, while mini-project 3 will focus on prediction.

---
### Vocabulary

- **Response variable**:

- **Explanatory variables**: 
  
---
### Benefits and Drawbacks of Models

- **Benefits**:

- **Drawbacks**:

---
### Looking for Potential Correlation

__Example__: Review NYC flight data. What variables, if any, are related to flight arrival delays? Consider a random sample of 1000 flights from NYC to Chicago in 2013.

**Focus**: if the time of the flights (`hour`) is related to arrival delay (`arr_delay`)?
+ A popular model we can use is a Linear Model.

---
### Linear model

__Population Linear model__: The relation between the observation `$Y$` and independent variables `$X_1,..., X_p$` is formulated as 
`$$Y = \beta_0 + \beta_1X_1 + ... + \beta_pX_p + \epsilon$$`

__Sample Linear model__: `$$Y = \hat{\beta_0} + \hat{\beta_1}X_1 + ... + \hat{\beta_p}X_p$$`

---
### Linear Model: One Numerical Predictor

``` r
model = lm(arr_delay ~ hour, data = Chicago1000)
summary(model)
```

```
## 
## Call:
## lm(formula = arr_delay ~ hour, data = Chicago1000)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -69.14 -25.96 -10.22   8.21 389.28 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -19.5500     4.2317  -4.620 4.34e-06 ***
## hour          2.0384     0.3094   6.587 7.23e-11 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 46.13 on 998 degrees of freedom
## Multiple R-squared:  0.04167,	Adjusted R-squared:  0.04071 
## F-statistic: 43.39 on 1 and 998 DF,  p-value: 7.227e-11
```

---
### Linear Model: One Numerical Predictor

**Fitted Linear Model**:

**Coefficient Interpretation**:

**Significance of Variable:** Hypothesis Test

---
### Linear Model: One Categorical Predictor

**Example**: What about carrier?

---
### Linear Model: One Categorical Predictor

When a categorical explanatory variable has many levels, they're encoded to **dummy variables**

| carrier        | AA           | B6  | MQ| 9E | WN |
| ------------- |:----:| -----:|-----:|-----:|-----:|-----:|
| UA | 0 | 0 | 0 |  0 | 0|
| AA | 1 | 0 | 0 |  0 | 0|
| B6 | 0 | 1 | 0 |  0 | 0|
| MQ | 0 | 0 | 1 |  0 | 0|
| 9E | 0 | 0 | 0 |  1 | 0|
| WN | 0 | 0 | 0 |  0 | 1|

---
### Linear Model: One Categorical Predictor

``` r
Chicago1000$carrier = relevel(as.factor(Chicago1000$carrier), 
                              ref = "UA")

model2 = lm(arr_delay ~ carrier, data = Chicago1000)
broom::tidy(model2)
```

```
## # A tibble: 6 × 5
## term estimate std.error statistic p.value
## <chr> <dbl> <dbl> <dbl> <dbl>
## 1 (Intercept) 4.86 2.59 1.88 0.0607 
## 2 carrier9E 4.35 6.33 0.688 0.492 
## 3 carrierAA -10.7 3.76 -2.83 0.00468
## 4 carrierB6 20.8 7.58 2.74 0.00620
## 5 carrierMQ 12.7 5.30 2.40 0.0165 
## 6 carrierWN 12.6 4.24 2.98 0.00297
```

Carrier UA is the reference level and every other levels will be compared to it.

---
### Linear Model with Categorical Explanatory

Writing out the proper model:

`$$\hat{\text{y}} = 4.86+4.35*\text{9E}-10.7*\text{AA}+20.8*\text{B6}+12.7*\text{MQ}+12.6*\text{WN}$$`

**For UA (Reference level)**:

**For AA:**

---
### Linear Model with Categorical Explanatory

Interpreting p-value:

???

### Analysis of Variance (ANOVA)

**Analysis of variance**: used to analyze the differences among means.

**Question**: Does at least one of the carriers have a different mean arrival delay than the others

``` r
model3 = aov(arr_delay ~ carrier, data = Chicago1000)
summary(model3)
```

```
##              Df  Sum Sq Mean Sq F value  Pr(>F)    
## carrier       5   95301   19060   8.932 2.6e-08 ***
## Residuals   994 2121209    2134                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
```

---
### Confounding Variables

.center[
**"Correlation does not equal causation."**
]

- A common concern with observational data is that correlations are distorted by other variables

- Addressing confounding is straightforward if these variables are measured (ex. can add them to model)

???

### Confounding Variables and Simpson's Paradox

Simpson's paradox is a phenomenon where a trend appears in several groups of data but disappears or reverses when the groups are combined. 
- The paradox can be resolved when confounding variables are appropriately addressed, but can lead to problematic interpretations if they are not.

---
### Carrier as a confounder?

__Example__: Does carrier confound the relationship between arrival delay and hour?

``` r
Chicago1000 %>% ggplot(aes(x = hour, y = arr_delay)) + 
  geom_point() + 
  geom_smooth(method = "lm") + 
  aes(color = carrier)
```

---
#### Multiple Linear Regression

``` r
model4 = lm(arr_delay ~ hour + carrier, data = Chicago1000)
```

```
## 
## Call:
## lm(formula = arr_delay ~ hour + carrier, data = Chicago1000)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -64.26 -24.12 -10.68   8.32 390.10 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -19.4396     4.6483  -4.182 3.14e-05 ***
## hour          1.9055     0.3051   6.246 6.25e-10 ***
## carrier9E     2.8147     6.2161   0.453  0.65078    
## carrierAA    -9.9018     3.6933  -2.681  0.00746 ** 
## carrierB6    20.2891     7.4437   2.726  0.00653 ** 
## carrierMQ    11.9620     5.1993   2.301  0.02161 *  
## carrierWN    11.6865     4.1668   2.805  0.00513 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 45.34 on 993 degrees of freedom
## Multiple R-squared:  0.07917,	Adjusted R-squared:  0.07361 
## F-statistic: 14.23 on 6 and 993 DF,  p-value: 1.337e-15
```

???

### Overall F Test

Testing for significance of all predictor variables at once: 
- `$H_0$`:

- `$H_A$`:

---
### Individual Significance Tests

Similar to model with one predictor variable, but p-value is calculated assuming all of the other variables are in the model.

- `$H_0$`:

- `$H_A$`:

---
### Variable Selection: Adjusted `$R^2$`

How much of the variation in our response is explained by the model

The Adjusted `$R^2$` does not automatically increase for additional variables

---
### Now we can create a linear model, but is it appropriate?

Has 5 Assumptions that needs to be met:

---
### Assumption: Normality

- If the response has a non-normal distribution, the estimates will be biased

To check this assumption, we use a QQ-Plot
  - A scatterplot created by plotting two sets of quantiles against one another.
  - If data is normal, we should see a roughly straight line

``` r
plot(model4, which = 2)
```

---
### Assumption: Linearity

Look at the Residuals vs Fitted Plot
- Look to see if there are any noticeable pattern in any of the plots

``` r
plot(model4, which = 1)
```

---
### Assumption: Constant Variance

- Looking for a constant spread of the residuals
- Don't want to see a narrowing or widening of points (fan shaped)

``` r
plot(model4, which = 3)
```

---
### Assumption: No Perfect Collinearity

Collinearity - a linear relationship between two explanatory variables

Why bad?

When is this an issue?

---
### Assumption: No Perfect Collinearity

Let's look at the relationship between `air_time` and `distance`

``` r
library(corrplot)
M <- cor(Chicago1000[,c(15:16)])
corrplot(M, method = "number", type = "upper")
```

---
### Outliers

Outliers: a data point that differs significantly from other observations

When you have an outlier:

???

### Last Thing

- Non-constant variance is one of the most common model violations, however it is usually fixable by transforming the response (y) variable.
  + Common transformation is the log-transformation or square root

- If you notice a curve in your data, you can add a quadratic term
  + `$\hat{Y} = \hat{\beta_0} + \hat{\beta_1}X_1 + \hat{\beta_2}X_1^2$`
  + `lm(arr_delay ~ hour + I(hour^2), data = Chicago1000)`