Model-Level Inference

The F-Distribution

In theoretical statistics the F-distribution is the ratio of two chi-squared statistics,

$$F=\frac{{\chi }_1^2/{\mathit{df}}_1}{{\chi }_2^2/{\mathit{df}}_2}$$

where ${\mathit{df}}_1$ and ${\mathit{df}}_2$ are the degrees of freedom associated with each of the chi-squared statistics, respectively. For our purposes, we don’t need to pay much attention to this other than to the fact that an F-distribution is defined using TWO parameters: ${\mathit{df}}_1$ and ${\mathit{df}}_2$. Knowing these two values completely parameterize the F-distribution (they give the shape, expected value, and variation).

In regression analysis, the F-distribution associated with model-level inference is based on the following degrees of freedom:

$$\begin{array}{rl}{\mathit{df}}_1& =p\\ {\mathit{df}}_2& ={\mathit{df}}_{\mathrm{Total}}-p\end{array}$$

where p is the number of predictors used in the model and $\mathrm{Total}$ is the total degrees of freedom in the data used in the regression model ($\mathrm{Total}=n-1$). In our example, ${\mathit{df}}_1=1$ and ${\mathit{df}}_2=99-1=98$. Using these values, we have defined the $F(1,98)$-distribution.

The F-distribution is the sampling distribution of F-values (not $R^2$-values). But, it turns out that we can easily convert an $R^2$-value to an F-value using,

$$F=\frac{R^2}{1-R^2}\times \frac{{\mathit{df}}_2}{{\mathit{df}}_1}$$

In our example,

$$\begin{array}{rl}F& =\frac{0.107}{1-0.107}\times \frac{98}{1}\\ & =0.1198\times 98\\ & =11.74\end{array}$$

Thus, our observed F-value is: $F(1,98)=11.74$. To evaluate this under the null hypothesis, we find the area under the $F(1,98)$ density curve that corresponds to F-values at least as extreme as our observed F-value of 11.74.

Figure 3: Plot of the probability curve for the F(1,98) distribution. The shaded area under the curve represents the p-value for a test evaluating whether the population rho-squared is zero using an observed F-value of 11.74.

This area (which is one-sided in the $F$-distribution) corresponds to the $p$-value. In our case this $p$-value is 0.000885. The probability of observing an $F$-value at least as extreme as we the one we observed ($F=11.74$) under the assumption that the null hypothesis is true is 0.000885. This suggests that the empirical data are inconsistent with the hypothesis that ${\rho }^2=0$, and it is unlikely that the model explains no variation in students’ GPAs.

Using the F-distribution in Practice

In practice, all of this information is provided in the output of the glance() function.

glance(lm.a)

# A tibble: 1 × 12
  r.squared adj.r.squared sigma statistic  p.value    df logLik   AIC   BIC
      <dbl>         <dbl> <dbl>     <dbl>    <dbl> <dbl>  <dbl> <dbl> <dbl>
1     0.107        0.0981  7.24      11.8 0.000885     1  -339.  684.  691.
  deviance df.residual  nobs
     <dbl>       <int> <int>
1    5136.          98   100

The observed F-value is given in the statistic column and the associated degrees of freedom are provided in the df and df.residual columns. Lastly, the p-value is given in the p.value column. When we report results from an F-test, we need to report the values for both degrees of freedom, the F-value, and the p-value.

The model-level test suggested that the empirical data are not consistent with the null hypothesis that the model explains no variation in GPAs; $F(1,98)=11.8$, $p<0.001$.

ANOVA Decomposition

We can also get the model-level inferential information from the anova() output. This gives us the ANOVA decomposition for the model.

anova(lm.a)

Analysis of Variance Table

Response: gpa
          Df Sum Sq Mean Sq F value    Pr(>F)    
homework   1  616.5  616.54  11.763 0.0008854 ***
Residuals 98 5136.4   52.41                      
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Note that the two df values for the model-level F-statistic correspond to the df in each row of the ANOVA table. The first df (in this case, 1) is the model degrees-of-freedom, and the second df (in this case, 98) is the residual degrees-of-freedom. Note the p-value is the same as that from the glance() function.

This ANOVA decomposition also breaks out the sum of squared values into the variation explained by the model (616.5) and that which is unexplained by the model (residual variation; 5136.4). Summing these two values will give the total amount of variation which can be used to compute $R^2$; $R^2={\mathrm{SS}}_{\mathrm{Model}}/{\mathrm{SS}}_{\mathrm{Total}}$.

This decomposition also gives us another way to consider the F-statistic. Recall that the F-statistic had a direct relationship to $R^2$

$$F=\frac{R^2}{1-R^2}\times \frac{{\mathit{df}}_2}{{\mathit{df}}_1}$$

Since $R^2={\mathrm{SS}}_{\mathrm{Model}}/{\mathrm{SS}}_{\mathrm{Total}}$ we can rewrite this as:

$$F=\frac{\frac{{\mathrm{SS}}_{\mathrm{Model}}}{{\mathrm{SS}}_{\mathrm{Total}}}}{1-\frac{{\mathrm{SS}}_{\mathrm{Model}}}{{\mathrm{SS}}_{\mathrm{Total}}}}\times \frac{{\mathit{df}}_2}{{\mathit{df}}_1}$$

Using algebra,

$$\begin{array}{rl}F& =\frac{\frac{{\mathrm{SS}}_{\mathrm{Model}}}{{\mathrm{SS}}_{\mathrm{Total}}}}{\frac{{\mathrm{SS}}_{\mathrm{Total}}}{{\mathrm{SS}}_{\mathrm{Total}}}-\frac{{\mathrm{SS}}_{\mathrm{Model}}}{{\mathrm{SS}}_{\mathrm{Total}}}}\times \frac{{\mathit{df}}_2}{{\mathit{df}}_1}\\ & =\frac{\frac{{\mathrm{SS}}_{\mathrm{Model}}}{{\mathrm{SS}}_{\mathrm{Total}}}}{\frac{{\mathrm{SS}}_{\mathrm{Total}}-{\mathrm{SS}}_{\mathrm{Model}}}{{\mathrm{SS}}_{\mathrm{Total}}}}\times \frac{{\mathit{df}}_2}{{\mathit{df}}_1}\\ & =\frac{{\mathrm{SS}}_{\mathrm{Model}}}{{\mathrm{SS}}_{\mathrm{Total}}-{\mathrm{SS}}_{\mathrm{Model}}}\times \frac{{\mathit{df}}_2}{{\mathit{df}}_1}\\ & =\frac{{\mathrm{SS}}_{\mathrm{Model}}}{{\mathrm{SS}}_{\mathrm{Error}}}\times \frac{{\mathit{df}}_2}{{\mathit{df}}_1}\end{array}$$

This expression of the F-statistic helps us see that the F-statistic is proportional to the ratio of the explained and unexplained variation. So long as the degrees of freedom remain the same, if the model explains more variation, the numerator of the F-statistic gets larger and the denominator will be smaller. Thus, larger F-values are associated with more explained variation by the model. We could also have seen this in the earlier expression of the F-statistic using $R^2$.

The F-Statistic as the Ratio of Two Variance Estimates

In statistical theory, a sum of squares divided by a degrees of freedom is referred to as a mean squared value—the average amount of variation per degree of freedom. Since ${\mathit{df}}_1$ is the model degrees of freedom and ${\mathit{df}}_2$ is the residual (or error) degrees of freedom we could also express the F-statistic as:

$$\begin{array}{rl}F& =\frac{{\mathrm{SS}}_{\mathrm{Model}}}{{\mathrm{SS}}_{\mathrm{Error}}}\times \frac{{\mathit{df}}_{\mathrm{Error}}}{{\mathit{df}}_{\mathrm{Model}}}\\ & =\frac{\frac{{\mathrm{SS}}_{\mathrm{Model}}}{{\mathit{df}}_{\mathrm{Model}}}}{\frac{{\mathrm{SS}}_{\mathrm{Error}}}{{\mathit{df}}_{\mathrm{Error}}}}\\ & =\frac{{\mathrm{MS}}_{\mathrm{Model}}}{{\mathrm{MS}}_{\mathrm{Error}}}\end{array}$$

Thus the F-value is the ratio of the average variation explained by the model and the average variation that remains unexplained. In our example

$$\begin{array}{rl}{\mathrm{MS}}_{\mathrm{Model}}& =\frac{616.5}{1}=616.5\\ {\mathrm{MS}}_{\mathrm{Error}}& =\frac{5136.4}{98}=52.41\end{array}$$

These values are also printed in the anova() output.

anova(lm.a)

Analysis of Variance Table

Response: gpa
          Df Sum Sq Mean Sq F value    Pr(>F)    
homework   1  616.5  616.54  11.763 0.0008854 ***
Residuals 98 5136.4   52.41                      
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

$$F=\frac{616.5}{52.41}=11.76$$

The observed F-value of 11.76 indicates that the average explained variation is 11.76 times that of the average unexplained variation. There is an awful lot more explained variation than unexplained variation, on average. Another name for a mean squared value is a variance estimate. A variance estimate is literally the average amount of variation (in the squared metric) per degree of freedom. For example, go back to the introductory statistics formula for using sample data to estimate a variance:

$${\hat{\sigma }}_Y^2=\frac{\sum (Y_i-\overline{Y}{)}^2}{n-1}$$

This numerator is a sum of squares; namely the ${\mathrm{SS}}_{\mathrm{Total}}$. The denominator is the total degrees of freedom. We could have also referred to this as a mean square

$$\begin{array}{rl}{\hat{\sigma }}_Y^2& =\frac{{\mathrm{SS}}_{\mathrm{Total}}}{{\mathit{df}}_{\mathrm{Total}}}\\ & ={\mathrm{MS}}_{\mathrm{Total}}\end{array}$$

Note that the ${\mathrm{MS}}_{\mathrm{Total}}$ is not printed in the anova() output. However, it can be computed from the values that are printed. The ${\mathrm{SS}}_{\mathrm{Total}}$ is just the sum of the printed sum of squares, and likewise the $${\mathit{df}}_{\mathrm{Total}}$$ is the sum of the df values.

$$\begin{array}{rl}{\mathrm{SS}}_{\mathrm{Total}}& =616.5+5136.4=5752.9\\ {\mathit{df}}_{\mathrm{Total}}& =1+98=99\end{array}$$

Then the ${\mathrm{MS}}_{\mathrm{Total}}$ is the ratio of these values,

$${\mathrm{MS}}_{\mathrm{Total}}=\frac{5752.9}{99}=58.11$$

Since this is a estimate of the outcome variable’s variance, we could also have computes the sample variance of the outcome variable, gpa, using the var() function.

keith %>%
  summarize(V_gpa = var(gpa))

# A tibble: 1 × 1
  V_gpa
  <dbl>
1  58.1

The total mean square, or variance estimate, is also the mean square estimate of the residuals from the intercept-only model.

# Fit intercept-only model
lm.0 = lm(gpa ~ 1, data = keith)

# ANOVA decomposition
anova(lm.0)

Analysis of Variance Table

Response: gpa
          Df Sum Sq Mean Sq F value Pr(>F)
Residuals 99 5752.9   58.11               

Remember the sum of squared residuals is $(Y_i-\widehat{Y_i}{)}^2$, but in the intercept-only model $\widehat{Y_i}$ is the marginal mean, i.e., $\widehat{Y_i}=\overline{Y}$. This is the numerator of the sample variance estimate and is why the mean square error from the intercept-only model and the sample variance for GPA are equivalent!

Relationship Between Coefficient-Level and Model-Level Inference

Lastly, we point out that in simple regression models (models with only one predictor), the results of the model-level inference (i.e., the p-value) is exactly the same as that for the coefficient-level inference for the slope.

# Model-level inference
glance(lm.a)

# A tibble: 1 × 12
  r.squared adj.r.squared sigma statistic  p.value    df logLik   AIC   BIC
      <dbl>         <dbl> <dbl>     <dbl>    <dbl> <dbl>  <dbl> <dbl> <dbl>
1     0.107        0.0981  7.24      11.8 0.000885     1  -339.  684.  691.
  deviance df.residual  nobs
     <dbl>       <int> <int>
1    5136.          98   100

# Coefficient-level inference
tidy(lm.a)

# A tibble: 2 × 5
  term        estimate std.error statistic  p.value
  <chr>          <dbl>     <dbl>     <dbl>    <dbl>
1 (Intercept)    74.3      1.94      38.3  1.01e-60
2 homework        1.21     0.354      3.43 8.85e- 4

That is because the model is composed of a single predictor, so asking whether the model accounts for variation in GPA is the same as asking whether GPA is different, on average, for students who spend a one-hour difference in time on homework. Once we have multiple predictors in the model, the model-level results and predictor-level results will not be the same.