Perform Simple Linear Regression

To estimate the linear relationship between the variables explored above, we’ll use lm() to fit a regression model using Ordinary Least Squares (OLS). The lm() function has two required arguments:

1. a formula that indicates which variable is the dependent variable and which is the independent variable;
2. a data argument that identifies the data frame for the variables in (i).

First let’s consider the relationship between crab size and latitude. We assume there is a linear population relationship, specified as:

\[\text{crab size} = \beta_0 + \beta_1 \text{latitude} + \epsilon\]

Exercise:

1. Use lm() to estimate \(\hat\beta_0\) and \(\hat\beta_1\) using this sample of data. Use summary(lm()) or gt() to visualize the regression results.

2. Interpret your two coefficients, paying careful attention to units.

3. Is your estimate of \(\beta_0\) directly interpretable? What might you be concerned about when interpreting \(\beta_0\)?

Answers:

# (i) latitude is associated with the size of crabs:
lm(size ~ latitude, data = pie_crab) %>% 
 tidy %>% 
  gt()

#alternatively, can write: summary(lm(size ~ latitude, data = pie_crab))

The coefficient on latitude tells us that for each degree increase in latitude, we expect crab size to increase by 0.49mm on average.

The intercept estimate tells us that on average, crabs living on the equator, i.e. at latitude=0, have a size of -3.62mm. Note that this doesn’t make sense! Do these crabs even live at the equator? This is a case where you want to be careful about extrapolating. \(\hat\beta_0\) is helping us fit the data well in the range of our sample, but you shouldn’t go making predictions about crabs in the equatorial region using this model estimated from the Atlantic coast of the US.

Exercises, continued

4. Does your model suggest anything about the relationship between crab size and water temperature? Is this evidence that water temperature is driving the Bergmann’s rule (i.e. the colder climate is causing the selection of larger organisms)?

Because crab size is increasing in latitude and water temperature is decreasing in latitude in the sample region, we can infer that as water temperature decreases, crab size increases. This is a statistical correlation we could estimate, but it is not causal evidence. You might think that colder water temperature is selecting for larger organisms by favoring larger bodies as they have lower surface-area-to-volume ratio (the heat conservation hypothesis), suggesting a causal link. However, for this hypothesis to be true, we should not see this relationship in “cold-blooded” animals, but many of them also follow this relationship.

“Bergmann originally hypothesized that the organisms he studied, birds, were larger in the colder, higher latitudes due to heat-conservation. But the heat-conservation hypothesis relies on internal regulation of body temperature and therefore does not apply to ectotherms,33 Commonly referred to as “cold-blooded animals”, ectotherms are animal that depends on external sources of body heat; such as fishes, amphibians, reptiles, and invertebrates. some of which also follow Bergmann’s rule. There is likely no universal mechanism underpinning Bergmann’s rule, regardless of ecto- or endothermy. As a result, other mechanisms have been proposed to explain Bergmann’s rule, including the starvation-resistant hypothesis, the diet-quality hypothesis, the enemy hypothesis, the resource rule, seasonality hypothesis, and the temperature–size rule (Johnson et al., 2019).”

5. [If we have time] Repeat the above steps for the relationship between latitude and water temperature.

Here are analogous results for the relationship between latitude and water temperature:

# (ii) latitude is associated with water temperature: 
lm(water_temp ~ latitude, data = pie_crab) %>% 
  tidy %>% 
  gt()

For each additional degree increase in latitude, we expect water temperature to fall by 0.77 degrees Celsius on average.

Visualize Regression [if we have time]

Here, we visualize our simple linear regression models above using the scatter plots from Section 1.

Exercise:

Use the geom_smooth() function with method argument set to lm, which stands for “linear model”, to add the best fit OLS line to your scatter plots from Section 1. Use the option se = FALSE to turn off estimates of uncertainty. Don’t worry, we’ll dig into those soon.

Answers:

# (i) latitude is associated with the size of crabs:
ggplot(data = pie_crab, aes(x = latitude, y = size)) +
  geom_point() +
  geom_smooth(method = "lm", se = FALSE) +
  labs(x = "Latitude", 
       y = "Size of crab carapace (mm)") + theme_classic()

## `geom_smooth()` using formula = 'y ~ x'

# (ii) latitude is associated with water temperature:
ggplot(data = pie_crab, aes(y = water_temp, x = latitude)) +
  geom_point() + 
  geom_smooth(method = "lm", formula= y~x, se=FALSE) + 
  labs(x = "Latitude",
       y = "Water temperature (Celsius)") + theme_classic()

Residuals, Sum of Squared Errors

The difference between what was observed in the data and what was predicted from the regression line is what we call an “error” or a “residual.” Observations that lie above the regression line exceeded their predicted value and have a positive residual. Observations that lie below the regression line are less than their predicted value and have a negative residual.

Recall that one of our assumptions needed in order to prove that OLS has the lowest variance of all unbiased estimators was that the population error \(u\) is: normally distributed, has a mean of zero, and has constant variance. While these assumptions are not directly testable because we can never observe \(u\), we can use our regression residuals \(e\) to assess the plausibility of these assumptions in the population.

Exercise:

Using the regression of crab size on latitude:

1. Add to the pie_crab dataframe a column of predicted total lengths \(\hat y_i\) and a column of residuals \(e_i = y_i-\hat y_i\).44 Hint: In dplyr you can use add_predictions(mod) where mod is a regression model to generate a new variable containing your model’s predicted values.

# regression
model_1 <- lm(size ~ latitude, data = pie_crab)

# create predictions and residuals
predictions <- pie_crab %>% add_predictions(model_1) %>%
  mutate(residuals = size-pred)

# take a look to understand what just happened
head(predictions)

## # A tibble: 6 × 11
##   date       latitude site   size air_temp air_tem…¹ water…² water…³ name   pred
##   <date>        <dbl> <chr> <dbl>    <dbl>     <dbl>   <dbl>   <dbl> <chr> <dbl>
## 1 2016-07-24       30 GTM    12.4     21.8      6.39    24.5    6.12 Guan…  10.9
## 2 2016-07-24       30 GTM    14.2     21.8      6.39    24.5    6.12 Guan…  10.9
## 3 2016-07-24       30 GTM    14.5     21.8      6.39    24.5    6.12 Guan…  10.9
## 4 2016-07-24       30 GTM    12.9     21.8      6.39    24.5    6.12 Guan…  10.9
## 5 2016-07-24       30 GTM    12.4     21.8      6.39    24.5    6.12 Guan…  10.9
## 6 2016-07-24       30 GTM    13.0     21.8      6.39    24.5    6.12 Guan…  10.9
## # … with 1 more variable: residuals <dbl>, and abbreviated variable names
## #   ¹​air_temp_sd, ²​water_temp, ³​water_temp_sd

2. Make a histogram of your residuals. Do they look approximately normally distributed?

# histogram
ggplot(data=predictions) + geom_histogram(aes(residuals), bins=25)+ theme_classic()

Pretty okay! Probably would merit a Q-Q plot to confirm, given the small sample size.

3. Compute the mean of your residuals. Is it basically zero?55 Note: OLS gives you this one for free. The OLS regression line will always go through the middle of your data, as it aims to minimize squared errors. So, you’ll always end up with sample errors that are mean zero. This doesn’t necessarily mean that in your population \(u\) is mean zero!

# mean
mean(predictions$residuals)

## [1] -0.00000000000003905975

4. Plot residuals against your independent variable, latitude. Does the variance look approximately constant across all values of latitude?

# variance in residuals against latitude
ggplot(predictions) + geom_point(aes(x=latitude, y=residuals)) + theme_classic()

Roughly, yes. Doesn’t look like there’s a systematic increase or decrease in variance, except maybe at the higher end of the latitude distribution. Overall, this looks like pretty strong evidence to support the assumption that the population error is normally distributed, mean zero, and constant variance. Since we are using our sample to justify assumptions we need to make about the population, we can never be fully confident. But this is the best we can do without a census!