How to use this document. Run the code along with me. Comments start with #. Each section mirrors the live lab; code chunks show the command and its output together.

Today’s plan

  1. Clean environment + setup
  2. Load & inspect the Prestige data
  3. Visualize distributions
  4. Linear regression: prestige ~ education
  5. Confidence intervals for the coefficients
  6. Quick practice
  7. Wrap-up

1 Clean Environment + Setup

rm(list = ls())   # clean environment

# Load the packages
library(ggplot2)
library(dplyr)
library(carData)  # contains the Prestige dataset

2 Load & Inspect the Prestige Data

data("Prestige")

# In RStudio you would run View(Prestige); here we preview the first rows.
knitr::kable(head(Prestige), caption = "First rows of the Prestige data")
First rows of the Prestige data
education income women prestige census type
gov.administrators 13.11 12351 11.16 68.8 1113 prof
general.managers 12.26 25879 4.02 69.1 1130 prof
accountants 12.77 9271 15.70 63.4 1171 prof
purchasing.officers 11.42 8865 9.11 56.8 1175 prof
chemists 14.62 8403 11.68 73.5 2111 prof
physicists 15.64 11030 5.13 77.6 2113 prof
summary(Prestige)
#>    education          income          women           prestige    
#>  Min.   : 6.380   Min.   :  611   Min.   : 0.000   Min.   :14.80  
#>  1st Qu.: 8.445   1st Qu.: 4106   1st Qu.: 3.592   1st Qu.:35.23  
#>  Median :10.540   Median : 5930   Median :13.600   Median :43.60  
#>  Mean   :10.738   Mean   : 6798   Mean   :28.979   Mean   :46.83  
#>  3rd Qu.:12.648   3rd Qu.: 8187   3rd Qu.:52.203   3rd Qu.:59.27  
#>  Max.   :15.970   Max.   :25879   Max.   :97.510   Max.   :87.20  
#>      census       type   
#>  Min.   :1113   bc  :44  
#>  1st Qu.:3120   prof:31  
#>  Median :5135   wc  :23  
#>  Mean   :5402   NA's: 4  
#>  3rd Qu.:8312            
#>  Max.   :9517

Variable notes

  • income = average income (in Canadian dollars)
  • education = average years of education
  • women = % of women in the occupation
  • prestige = social prestige score (our dependent variable)

3 Visualize Distributions

Histogram of prestige (default binning).

ggplot(Prestige, aes(x = prestige)) +
  geom_histogram(color = "pink") +
  labs(title = "Distribution of Occupational Prestige",
       x = "Prestige", y = "Frequency")

Histogram of prestige (fixed bin width).

ggplot(Prestige, aes(x = prestige)) +
  geom_histogram(binwidth = 5, fill = "skyblue", color = "black") +
  labs(title = "Distribution of Occupational Prestige",
       x = "Prestige", y = "Frequency")

Histogram of education.

ggplot(Prestige, aes(x = education)) +
  geom_histogram(binwidth = 1, fill = "lightgreen", color = "black") +
  labs(title = "Distribution of Education Years",
       x = "Education (Years)", y = "Frequency")

4 Linear Regression: prestige ~ education

Predict the social prestige score from years of education.

model <- lm(prestige ~ education, data = Prestige)
summary(model)
#> 
#> Call:
#> lm(formula = prestige ~ education, data = Prestige)
#> 
#> Residuals:
#>      Min       1Q   Median       3Q      Max 
#> -26.0397  -6.5228   0.6611   6.7430  18.1636 
#> 
#> Coefficients:
#>             Estimate Std. Error t value Pr(>|t|)    
#> (Intercept)  -10.732      3.677  -2.919  0.00434 ** 
#> education      5.361      0.332  16.148  < 2e-16 ***
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> Residual standard error: 9.103 on 100 degrees of freedom
#> Multiple R-squared:  0.7228, Adjusted R-squared:   0.72 
#> F-statistic: 260.8 on 1 and 100 DF,  p-value: < 2.2e-16
  • Intercept = predicted prestige when education = 0.
  • education = expected change in prestige for a 1-year increase in education.

5 Confidence Intervals (95%)

confint(model, level = 0.95)
#>                  2.5 %    97.5 %
#> (Intercept) -18.027220 -3.436744
#> education     4.702223  6.019533
# Lower and upper bounds for each coefficient

6 Quick Practice

P1) Add income as another predictor. Include multiple predictors by connecting them with +. Compare coefficients and p-values.

model2 <- lm(prestige ~ education + income, data = Prestige)
summary(model2)
#> 
#> Call:
#> lm(formula = prestige ~ education + income, data = Prestige)
#> 
#> Residuals:
#>      Min       1Q   Median       3Q      Max 
#> -19.4040  -5.3308   0.0154   4.9803  17.6889 
#> 
#> Coefficients:
#>               Estimate Std. Error t value Pr(>|t|)    
#> (Intercept) -6.8477787  3.2189771  -2.127   0.0359 *  
#> education    4.1374444  0.3489120  11.858  < 2e-16 ***
#> income       0.0013612  0.0002242   6.071 2.36e-08 ***
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> Residual standard error: 7.81 on 99 degrees of freedom
#> Multiple R-squared:  0.798,  Adjusted R-squared:  0.7939 
#> F-statistic: 195.6 on 2 and 99 DF,  p-value: < 2.2e-16

P2) Using confint, which predictors have CIs that include 0?

confint(model2, level = 0.95)
#>                     2.5 %       97.5 %
#> (Intercept) -1.323493e+01 -0.460629799
#> education    3.445127e+00  4.829761535
#> income       9.162805e-04  0.001806051

P3) Try predicting prestige using the % of women.

model_women <- lm(prestige ~ women, data = Prestige)
summary(model_women)
#> 
#> Call:
#> lm(formula = prestige ~ women, data = Prestige)
#> 
#> Residuals:
#>     Min      1Q  Median      3Q     Max 
#> -33.444 -12.391  -4.126  13.034  39.185 
#> 
#> Coefficients:
#>             Estimate Std. Error t value Pr(>|t|)    
#> (Intercept) 48.69300    2.30760  21.101   <2e-16 ***
#> women       -0.06417    0.05385  -1.192    0.236    
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> Residual standard error: 17.17 on 100 degrees of freedom
#> Multiple R-squared:  0.014,  Adjusted R-squared:  0.004143 
#> F-statistic:  1.42 on 1 and 100 DF,  p-value: 0.2362
confint(model_women, level = 0.95)
#>                 2.5 %      97.5 %
#> (Intercept) 44.114782 53.27121657
#> women       -0.171008  0.04266233

7 Wrap-up

Key takeaways

  • confint() gives the CI for regression coefficients.
  • p-values tell us whether relationships are statistically significant.
  • Education and income often predict social prestige — a core social science insight.
  • A CI that includes 0 means the predictor is not statistically distinguishable from “no effect.”
  • Next week: one-sample t-test for group means.
LS0tCnRpdGxlOiAiTGFiIDYg4oCUIENvbmZpZGVuY2UgSW50ZXJ2YWxzICYgUmVncmVzc2lvbiBCYXNpY3MiCnN1YnRpdGxlOiAiUXVhbnRpdGF0aXZlIFJlYXNvbmluZyDCtyBMQUIgNDEzIgphdXRob3I6ICJJbnN0cnVjdG9yOiBTdWJpbiBOYSIKZGF0ZTogIk9jdG9iZXIgMTUsIDIwMjUiCm91dHB1dDoKICBodG1sX2RvY3VtZW50OgogICAgdGhlbWU6IGZsYXRseQogICAgaGlnaGxpZ2h0OiB0YW5nbwogICAgdG9jOiB0cnVlCiAgICB0b2NfZmxvYXQ6IHRydWUKICAgIHRvY19kZXB0aDogMgogICAgbnVtYmVyX3NlY3Rpb25zOiB0cnVlCiAgICBkZl9wcmludDogcGFnZWQKICAgIGNvZGVfZG93bmxvYWQ6IHRydWUKLS0tCgpgYGB7ciBzZXR1cCwgaW5jbHVkZT1GQUxTRX0Ka25pdHI6Om9wdHNfY2h1bmskc2V0KAogIGVjaG8gPSBUUlVFLCBtZXNzYWdlID0gRkFMU0UsIHdhcm5pbmcgPSBGQUxTRSwKICBmaWcuYWxpZ24gPSAiY2VudGVyIiwgZmlnLndpZHRoID0gNywgZmlnLmhlaWdodCA9IDQuMiwKICBjb21tZW50ID0gIiM+IgopCmBgYAoKPiAqKkhvdyB0byB1c2UgdGhpcyBkb2N1bWVudC4qKiBSdW4gdGhlIGNvZGUgYWxvbmcgd2l0aCBtZS4gQ29tbWVudHMgc3RhcnQgd2l0aCBgI2AuCj4gRWFjaCBzZWN0aW9uIG1pcnJvcnMgdGhlIGxpdmUgbGFiOyBjb2RlIGNodW5rcyBzaG93IHRoZSBjb21tYW5kIGFuZCBpdHMgb3V0cHV0IHRvZ2V0aGVyLgoKKipUb2RheSdzIHBsYW4qKgoKMS4gQ2xlYW4gZW52aXJvbm1lbnQgKyBzZXR1cAoyLiBMb2FkICYgaW5zcGVjdCB0aGUgYFByZXN0aWdlYCBkYXRhCjMuIFZpc3VhbGl6ZSBkaXN0cmlidXRpb25zCjQuIExpbmVhciByZWdyZXNzaW9uOiBgcHJlc3RpZ2UgfiBlZHVjYXRpb25gCjUuIENvbmZpZGVuY2UgaW50ZXJ2YWxzIGZvciB0aGUgY29lZmZpY2llbnRzCjYuIFF1aWNrIHByYWN0aWNlCjcuIFdyYXAtdXAKCi0tLQoKIyBDbGVhbiBFbnZpcm9ubWVudCArIFNldHVwCgpgYGB7ciBjbGVhbi1zZXR1cH0Kcm0obGlzdCA9IGxzKCkpICAgIyBjbGVhbiBlbnZpcm9ubWVudAoKIyBMb2FkIHRoZSBwYWNrYWdlcwpsaWJyYXJ5KGdncGxvdDIpCmxpYnJhcnkoZHBseXIpCmxpYnJhcnkoY2FyRGF0YSkgICMgY29udGFpbnMgdGhlIFByZXN0aWdlIGRhdGFzZXQKYGBgCgojIExvYWQgJiBJbnNwZWN0IHRoZSBQcmVzdGlnZSBEYXRhCgpgYGB7ciBpbXBvcnR9CmRhdGEoIlByZXN0aWdlIikKCiMgSW4gUlN0dWRpbyB5b3Ugd291bGQgcnVuIFZpZXcoUHJlc3RpZ2UpOyBoZXJlIHdlIHByZXZpZXcgdGhlIGZpcnN0IHJvd3MuCmtuaXRyOjprYWJsZShoZWFkKFByZXN0aWdlKSwgY2FwdGlvbiA9ICJGaXJzdCByb3dzIG9mIHRoZSBQcmVzdGlnZSBkYXRhIikKCnN1bW1hcnkoUHJlc3RpZ2UpCmBgYAoKPiAqKlZhcmlhYmxlIG5vdGVzKioKPgo+IC0gYGluY29tZWAgPSBhdmVyYWdlIGluY29tZSAoaW4gQ2FuYWRpYW4gZG9sbGFycykKPiAtIGBlZHVjYXRpb25gID0gYXZlcmFnZSB5ZWFycyBvZiBlZHVjYXRpb24KPiAtIGB3b21lbmAgPSAlIG9mIHdvbWVuIGluIHRoZSBvY2N1cGF0aW9uCj4gLSBgcHJlc3RpZ2VgID0gc29jaWFsIHByZXN0aWdlIHNjb3JlIChvdXIgZGVwZW5kZW50IHZhcmlhYmxlKQoKIyBWaXN1YWxpemUgRGlzdHJpYnV0aW9ucwoKKipIaXN0b2dyYW0gb2YgcHJlc3RpZ2UgKGRlZmF1bHQgYmlubmluZykuKioKCmBgYHtyIGhpc3QtcHJlc3RpZ2UtZGVmYXVsdH0KZ2dwbG90KFByZXN0aWdlLCBhZXMoeCA9IHByZXN0aWdlKSkgKwogIGdlb21faGlzdG9ncmFtKGNvbG9yID0gInBpbmsiKSArCiAgbGFicyh0aXRsZSA9ICJEaXN0cmlidXRpb24gb2YgT2NjdXBhdGlvbmFsIFByZXN0aWdlIiwKICAgICAgIHggPSAiUHJlc3RpZ2UiLCB5ID0gIkZyZXF1ZW5jeSIpCmBgYAoKKipIaXN0b2dyYW0gb2YgcHJlc3RpZ2UgKGZpeGVkIGJpbiB3aWR0aCkuKioKCmBgYHtyIGhpc3QtcHJlc3RpZ2V9CmdncGxvdChQcmVzdGlnZSwgYWVzKHggPSBwcmVzdGlnZSkpICsKICBnZW9tX2hpc3RvZ3JhbShiaW53aWR0aCA9IDUsIGZpbGwgPSAic2t5Ymx1ZSIsIGNvbG9yID0gImJsYWNrIikgKwogIGxhYnModGl0bGUgPSAiRGlzdHJpYnV0aW9uIG9mIE9jY3VwYXRpb25hbCBQcmVzdGlnZSIsCiAgICAgICB4ID0gIlByZXN0aWdlIiwgeSA9ICJGcmVxdWVuY3kiKQpgYGAKCioqSGlzdG9ncmFtIG9mIGVkdWNhdGlvbi4qKgoKYGBge3IgaGlzdC1lZHVjYXRpb259CmdncGxvdChQcmVzdGlnZSwgYWVzKHggPSBlZHVjYXRpb24pKSArCiAgZ2VvbV9oaXN0b2dyYW0oYmlud2lkdGggPSAxLCBmaWxsID0gImxpZ2h0Z3JlZW4iLCBjb2xvciA9ICJibGFjayIpICsKICBsYWJzKHRpdGxlID0gIkRpc3RyaWJ1dGlvbiBvZiBFZHVjYXRpb24gWWVhcnMiLAogICAgICAgeCA9ICJFZHVjYXRpb24gKFllYXJzKSIsIHkgPSAiRnJlcXVlbmN5IikKYGBgCgojIExpbmVhciBSZWdyZXNzaW9uOiBgcHJlc3RpZ2UgfiBlZHVjYXRpb25gCgpQcmVkaWN0IHRoZSBzb2NpYWwgcHJlc3RpZ2Ugc2NvcmUgZnJvbSB5ZWFycyBvZiBlZHVjYXRpb24uCgpgYGB7ciBtb2RlbH0KbW9kZWwgPC0gbG0ocHJlc3RpZ2UgfiBlZHVjYXRpb24sIGRhdGEgPSBQcmVzdGlnZSkKc3VtbWFyeShtb2RlbCkKYGBgCgo+IC0gKipJbnRlcmNlcHQqKiA9IHByZWRpY3RlZCBwcmVzdGlnZSB3aGVuIGBlZHVjYXRpb25gID0gMC4KPiAtICoqZWR1Y2F0aW9uKiogPSBleHBlY3RlZCBjaGFuZ2UgaW4gcHJlc3RpZ2UgZm9yIGEgMS15ZWFyIGluY3JlYXNlIGluIGVkdWNhdGlvbi4KCiMgQ29uZmlkZW5jZSBJbnRlcnZhbHMgKDk1JSkKCmBgYHtyIGNvbmZpbnR9CmNvbmZpbnQobW9kZWwsIGxldmVsID0gMC45NSkKIyBMb3dlciBhbmQgdXBwZXIgYm91bmRzIGZvciBlYWNoIGNvZWZmaWNpZW50CmBgYAoKIyBRdWljayBQcmFjdGljZQoKKipQMSkgQWRkIGBpbmNvbWVgIGFzIGFub3RoZXIgcHJlZGljdG9yLioqCkluY2x1ZGUgbXVsdGlwbGUgcHJlZGljdG9ycyBieSBjb25uZWN0aW5nIHRoZW0gd2l0aCBgK2AuIENvbXBhcmUgY29lZmZpY2llbnRzIGFuZCAqcCotdmFsdWVzLgoKYGBge3IgcDF9Cm1vZGVsMiA8LSBsbShwcmVzdGlnZSB+IGVkdWNhdGlvbiArIGluY29tZSwgZGF0YSA9IFByZXN0aWdlKQpzdW1tYXJ5KG1vZGVsMikKYGBgCgoqKlAyKSBVc2luZyBgY29uZmludGAsIHdoaWNoIHByZWRpY3RvcnMgaGF2ZSBDSXMgdGhhdCBpbmNsdWRlIDA/KioKCmBgYHtyIHAyfQpjb25maW50KG1vZGVsMiwgbGV2ZWwgPSAwLjk1KQpgYGAKCioqUDMpIFRyeSBwcmVkaWN0aW5nIHByZXN0aWdlIHVzaW5nIHRoZSAlIG9mIHdvbWVuLioqCgpgYGB7ciBwM30KbW9kZWxfd29tZW4gPC0gbG0ocHJlc3RpZ2UgfiB3b21lbiwgZGF0YSA9IFByZXN0aWdlKQpzdW1tYXJ5KG1vZGVsX3dvbWVuKQpjb25maW50KG1vZGVsX3dvbWVuLCBsZXZlbCA9IDAuOTUpCmBgYAoKIyBXcmFwLXVwCgoqKktleSB0YWtlYXdheXMqKgoKLSBgY29uZmludCgpYCBnaXZlcyB0aGUgQ0kgZm9yIHJlZ3Jlc3Npb24gY29lZmZpY2llbnRzLgotICpwKi12YWx1ZXMgdGVsbCB1cyB3aGV0aGVyIHJlbGF0aW9uc2hpcHMgYXJlIHN0YXRpc3RpY2FsbHkgc2lnbmlmaWNhbnQuCi0gRWR1Y2F0aW9uIGFuZCBpbmNvbWUgb2Z0ZW4gcHJlZGljdCBzb2NpYWwgcHJlc3RpZ2Ug4oCUIGEgY29yZSBzb2NpYWwgc2NpZW5jZSBpbnNpZ2h0LgotIEEgQ0kgdGhhdCBpbmNsdWRlcyAwIG1lYW5zIHRoZSBwcmVkaWN0b3IgaXMgbm90IHN0YXRpc3RpY2FsbHkgZGlzdGluZ3Vpc2hhYmxlIGZyb20gIm5vIGVmZmVjdC4iCi0gTmV4dCB3ZWVrOiBvbmUtc2FtcGxlICp0Ki10ZXN0IGZvciBncm91cCBtZWFucy4K