Section B.2 Fitting a Simple Regression Model
In order to draw a regression line (or in fact any line in a Cartesian coordinate system), we need to know two things:
-
Where the line begins: what is the value of Y (our dependent variable) when X (our independent variable) is 0, so that we have a point from which to start drawing the line. The technical name for this point is the intercept.
-
What is the slope of that line, that is, how inclined the line is, the angle of the line.
If you recall from elementary algebra, the equation for any straight line is: \(y = mx + b\text{.}\) In statistics we use a slightly different notation, although the equation remains the same: \(y = \beta_0 + \beta_1 x\text{.}\)
We need the origin of the line (\(\beta_0\)) and the slope of the line (\(\beta_1\)). How does R get the intercept and the slope for the line? How does R know where to draw this line? We need to estimate these parameters (or coefficients) from the data. How? We don’t have the time to get into these more mathematical details now, but will suggest further reading at the end of this appendix for the curious. For now, suffice to say that for linear regression models like the one we cover here, when drawing the line, R tries to minimise the distance from every point in the scatterplot to the regression line using a method called least squares estimation.
In order to fit the model, we use the
lm() function using the formula specification \((Y \sim X)\text{.}\) Typically you want to store your regression model in a “variable”, let’s call it fit_1:
fit_1 <- lm(HR90 ~ RD90, data = ncovr)
You will see in your RStudio global environment space that there is a new object called
fit_1 with 12 elements on it. We can get a sense for what this object is and includes using the functions we introduced in previous weeks:
class(fit_1)
attributes(fit_1)
## [1] "lm"
## $names ## [1] "coefficients" "residuals" "effects" ## [4] "rank" "fitted.values" "assign" ## [7] "qr" "df.residual" "xlevels" ## [10] "call" "terms" "model" ## ## $class ## [1] "lm"
R is telling us that this is an object of class
lm and that it includes a number of attributes. One of the beauties of R is that you are producing all the results from running the model, putting them in an object, and then giving you the opportunity for using them later on. If you want to simply see the basic results from running the model, you can use the summary() function.
summary(fit_1)
## ## Call: ## lm(formula = HR90 ~ RD90, data = ncovr) ## ## Residuals: ## Min 1Q Median 3Q Max ## -17.80 -3.42 -0.72 2.54 67.10 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) 6.1829 0.0984 62.8 <2e-16 *** ## RD90 3.7712 0.0985 38.3 <2e-16 *** ## --- ## Signif. codes: ## 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Residual standard error: 5.47 on 3083 degrees of freedom ## Multiple R-squared: 0.322, Adjusted R-squared: 0.322 ## F-statistic: 1.47e+03 on 1 and 3083 DF, p-value: <2e-16
Or if you prefer more parsimonious presentation, you could use the
display() function of the arm package:
library(arm)
display(fit_1)
## lm(formula = HR90 ~ RD90, data = ncovr) ## coef.est coef.se ## (Intercept) 6.18 0.10 ## RD90 3.77 0.10 ## --- ## n = 3085, k = 2 ## residual sd = 5.47, R-Squared = 0.32
For now I just want you to focus on the numbers in the “Estimate” (or coef.est) column. The value of 6.18 estimated for the intercept is the “predicted” value for Y when X equals zero. This is the predicted value of the homicide score when resource deprivation has a value of zero.
summary(ncovr$RD90)
## Min. 1st Qu. Median Mean 3rd Qu. Max. ## -2.410 -0.667 -0.202 0.000 0.439 5.583
RD90 is a variable that has been centered at 0. It has been created by the researchers in such a way that it has a mean value of 0. Since we only have one explanatory variable in the model, this corresponds to the mean of the homicide rate, 6.18. In many other contexts the intercept has less of a meaning.
We then need the \(\beta_1\) regression coefficient for our independent variable, the value that will shape the slope in this scenario. This value is 3.77. This estimated regression coefficient for our independent variable has a convenient interpretation. When the value is positive, it tells us that for every one unit increase in X, there is a \(\beta_1\) increase on Y. If the coefficient is negative, then it represents a decrease on Y. Here, we can read it as “for every one unit increase in the resource deprivation score, there is a 3.77 unit increase in the homicide rate.”
Knowing these two parameters not only allows us to draw the line, but we can also solve for any given value of X. Let’s go back to our guess-the-homicide-rate game. Imagine I tell you the level of resource deprivation is 1. What would be your best bet now? We can simply go back to our regression line equation and insert the estimated parameters:
\begin{equation*}
y = b_0 + b_1 x
\end{equation*}
\begin{equation*}
y = 6.18 + 3.77 \times 1
\end{equation*}
\begin{equation*}
y = 9.95
\end{equation*}
Or if you don’t want to do the calculation yourself, you can use the
predict() function (differences are due to rounding error). First, you name your stored model and then you identify the new data (which has to be in a dataframe format and with a variable name matching the one in the original data set).
predict(fit_1, data.frame(RD90 = c(1)))
## 1 ## 9.954
This is the expected value of Y, homicide rate, when X, resource deprivation is 1 according to our model (according to our simplification of the real world, our simplification of the whole cloud of points into just one straight line). Look back at the scatterplot we produced earlier. Does it look as if the red line when X is 1 corresponds to a value of Y of 9.95?
