18.05 S22 Class 26: Linear Regression

Linear regression Class 26, 18.05 1 Learning Goals 1. Be able to use the method of least squares to fit a line to bivariate data. 2. Be able to give a formula for the total squared error when fitting any type of curve to data. 3. Be able to say the words homoscedasticity and heteroscedasticity. 2 Introduction Suppose we have collected bivariate data (𝑥𝑖 , 𝑦𝑖 ), 𝑖 = 1, … , 𝑛. The goal of linear regression is to model the relationship between 𝑥 and 𝑦 by finding a function 𝑦 = 𝑓(𝑥) that is a close fit to the data. The modeling assumptions we will use are that 𝑥𝑖 is not random and that 𝑦𝑖 is a function of 𝑥𝑖 plus some random noise. With these assumptions 𝑥 is called the independent or predictor variable and 𝑦 is called the dependent or response variable. Here is a series of examples showing the results of linear regression. We will discuss the details of how to do linear regression in the next section. Example 1. 0.05 (1963) 0.22 (1985) 0.39 (2006) 0.49 (2017) The cost of a 0.06 (1968) 0.25 (1988) 0.41 (2007) 0.50 (2018) first class stamp in cents over time is given in the following list. 0.08 (1971) 0.10 (1974) 0.13 (1975) 0.15 (1978) 0.20 (1981) 0.29 (1991) 0.32 (1995) 0.33 (1999) 0.34 (2001) 0.37 (2002) 0.42 (2008) 0.44 (2009) 0.45 (2012) 0.46 (2013) 0.49 (2015) 0.55 (2019) Using the R function lm we found the ‘least squares fit’ for a line to this data is 𝑦 = −0.21390 + 0.88203𝑥, where 𝑥 is the number of years since 1960 and 𝑦 is in cents. Using this result we ‘predict’ that in 2021 (𝑥 = 61) the cost of a stamp will be 53.6 cents (since −0.21390 + 0.88203 ⋅ 61 = 53.6). 1 2 30 10 20 y 40 50 18.05 Class 26, Linear regression, Spring 2022 0 10 20 30 40 50 60 x Stamp cost (cents) vs. time (years since 1960). Orange dot is predicted cost in 2021. Note that none of the data points actually lie on the line. Rather this line has the ‘best fit’ with respect to all the data, with a small error for each data point. (Note, the actual cost of a stamp dropped in January 2021 was 55 cents. See https: //en.wikipedia.org/wiki/History_of_United_States_postage_rates) Example 2. Suppose we have 𝑛 pairs of fathers and adult sons. Let 𝑥𝑖 and 𝑦𝑖 be the heights of the 𝑖th father and son, respectively. The least squares line for this data could be used to predict the adult height of a young boy from that of his father. Example 3. We are not limited to best fit lines. For all positive 𝑑, the method of least squares may be used to find a polynomial of degree 𝑑 with the ‘best fit’ to the data. Here’s a figure showing the least squares fit of a parabola (𝑑 = 2). Fitting a parabola, 𝑎𝑥2 + 𝑏𝑥 + 𝑐, to data Example 4. In fact, we can use linear regression to fit many other types of curves to bivariate data. 18.05 Class 26, Linear regression, Spring 2022 3 3 Fitting a line using least squares Suppose we have data (𝑥𝑖 , 𝑦𝑖 ) as above. Our first goal is to find the line 𝑦 = 𝑎𝑥 + 𝑏 that ‘best fits’ the data. Our model says that each 𝑦𝑖 is predicted by 𝑥𝑖 up to some error 𝜖𝑖 : 𝑦𝑖 = 𝑎𝑥𝑖 + 𝑏 + 𝜖𝑖 . So 𝜖𝑖 = 𝑦𝑖 − 𝑎𝑥𝑖 − 𝑏. The method of least squares finds the values 𝑎̂ and 𝑏̂ of 𝑎 and 𝑏 that minimize the sum of the squared errors: 𝑆(𝑎, 𝑏) = ∑ 𝜖2𝑖 = ∑(𝑦𝑖 − 𝑎𝑥𝑖 − 𝑏)2 . 𝑖 Using calculus or linear algebra (details in the appendix), we find 𝑎̂ = 𝑠𝑥𝑦 𝑠𝑥𝑥 𝑏̂ = 𝑦 ̄ − 𝑎̂ 𝑥̄ (1) where 𝑥̄ = 1 ∑ 𝑥𝑖 , 𝑛 𝑦̄ = 1 ∑ 𝑦𝑖 , 𝑛 𝑠𝑥𝑥 = 1 ∑(𝑥𝑖 −𝑥)̄ 2 , (𝑛 − 1) 𝑠𝑥𝑦 = 1 ∑(𝑥𝑖 −𝑥)(𝑦 ̄ 𝑖 −𝑦). ̄ (𝑛 − 1) Here 𝑥̄ is the sample mean of 𝑥, 𝑦̄ is the sample mean of 𝑦, 𝑠𝑥𝑥 is the sample variance of 𝑥, and 𝑠𝑥𝑦 is the sample covariance of 𝑥 and 𝑦. Example 5. Use least squares to fit a line to the following data: (0,1), (2,1), (3,4). Solution: In our case, (𝑥1 , 𝑦1 ) = (0, 1), (𝑥2 , 𝑦2 ) = (2, 1) and (𝑥3 , 𝑦3 ) = (3, 4). So 5 𝑥̄ = , 𝑦̄ = 2, 3 𝑠𝑥𝑥 = 7 , 3 𝑠𝑥𝑦 = 2 Using the above formulas we get 6 4 , 𝑏̂ = . 7 7 6 4 So the least squares line has equation 𝑦 = 𝑥 + . This is shown as the orange line in the 7 7 following figure. We will discuss the blue parabola soon. 𝑎̂ = 𝑦 4 1 𝑥 1 2 3 Least squares fit of a line (orange) and a parabola (blue) 18.05 Class 26, Linear regression, Spring 2022 4 Simple linear regression: It’s a little confusing, but the word linear in ‘linear regression’ does not refer to fitting a line. We will explain its meaning below. However, the most common curve to fit is a line. When we fit a line to bivariate data it is called simple linear regression. 3.1 Residuals For a line the model is 𝑦𝑖 = 𝑎𝑥 ̂ + 𝑏̂ + 𝜖 𝑖 . 20 We think of 𝑎𝑥 ̂ 𝑖 +𝑏̂ as predicting or explaining 𝑦𝑖 . The left-over term 𝜖𝑖 is called the residual, which we think of as random noise or measurement error. A useful visual check of the linear regression model is to plot the residuals. The data points should hover near the regression line. The residuals should look about the same across the range of 𝑥. ● 4 ● ● ● ● 15 ● ● ● ● ● ● ● ●● ● ●● ● ●● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ●● ● 2 4 6 x 8 10 ● ● ● ● ● ● ● ● 0 ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● −3 0 ● ● ●● ● ● ● ● ● ● ● ●● ● ● ● 0 ● ●● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● −1 ● ●● ● ● ● ● ● ● ● ● ● ● ●●● ● −2 5 ● ●● ● ● ● ● e 1 ● ● ● 0 y 10 ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ●● ● ● ● 2 ● ●● ● ●● ● ● ● ● 3 ● ● ● ●● ● 2 4 6 8 10 x Data with regression line (left) and residuals (right). Note the homoscedasticity. 3.2 Homoscedasticity An important assumption of the linear regression model is that the residuals 𝜖𝑖 have the same variance for all 𝑖. This is called homoscedasticity. You can see this is the case for both figures above. The data hovers in the band of fixed width around the regression line and at every 𝑥 the residuals have about the same vertical spread. Below is a figure showing heteroscedastic data. The vertical spread of the data increases as 𝑥 increases. Before using least squares on this data we would need to transform the data to be homoscedastic. 5 20 18.05 Class 26, Linear regression, Spring 2022 ● ● ● ● ● ● ● ● 15 ● ●● ● ● ● ● ● ● y 10 ● ●● ●● ● ● ● ● ● ● ● ● ● ● ● ● ●●● ● ● ● ● ●● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ●● ● ●● ●● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ●●● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ●● ● ● ●● ● ● ● ● ●● ● ● ● ● ●● ● ● ● ●● ● 0 5 ● ● ● ● ● 0 2 4 6 8 10 x Heteroscedastic Data 4 Linear regression for fitting polynomials When we fit a line to data it is called simple linear regression. We can also use linear regression to fit polynomials to data. The use of the word linear in both cases may seem confusing. This is because the word ‘linear’ in linear regression does not refer to fitting a line. Rather it refers to the linear algebraic equations for the unknown parameters. Example 6. Take the same data as in Example 5 and use least squares to find the best fitting parabola to the data. Solution: A parabola has the formula 𝑦 = 𝑎𝑥2 + 𝑏𝑥 + 𝑐. The squared error is 𝑆(𝑎, 𝑏, 𝑐) = ∑(𝑦𝑖 − (𝑎𝑥2𝑖 + 𝑏𝑥𝑖 + 𝑐))2 . After substituting the given values for each 𝑥𝑖 and 𝑦𝑖 , we can use calculus to find the triple (𝑎, 𝑏, 𝑐) that minimizes 𝑆. With this data, we find that the least squares parabola has equation 𝑦 = 𝑥2 − 2𝑥 + 1. Note that for 3 points the quadratic fit is perfect. 𝑦 4 1 𝑥 1 2 3 Least squares fit of a line (orange) and a parabola (blue) Example 7. The pairs (𝑥𝑖 , 𝑦𝑖 ) may give the age and vocabulary size of a 𝑛 children. Since we expect that young children acquire new words at an accelerating pace, we might guess that a higher order polynomial would best fit the data. 18.05 Class 26, Linear regression, Spring 2022 6 Example 8. (Transforming the data) Sometimes it is necessary to transform the data before using linear regression. For example, let’s suppose the relationship is exponential, i.e. 𝑦 = 𝑐𝑒𝑎𝑥 . Then ln(𝑦) = 𝑎𝑥 + ln(𝑐). So we can use simple linear regression on the data (𝑥𝑖 , ln(𝑦𝑖 )) to obtain a model ln(𝑦) = 𝑎𝑥 ̂ + 𝑏̂ and then exponentiate to obtain the exponential model ̂ 𝑦 = 𝑒𝑏 𝑒𝑎𝑥̂ . 4.1 Overfitting You can always achieve a better fit by using a higher order polynomial. For instance, given 6 data points (with distinct 𝑥𝑖 ) one can always find a fifth order polynomial that goes through all of them. This can result in what’s called overfitting. That is, fitting the noise as well as the true relationship between 𝑥 and 𝑦. An overfit model will fit the original data better but perform less well on predicting 𝑦 for new values of 𝑥. Indeed, a primary challenge of statistical modeling is balancing model fit against model complexity. Example 9. In the plot below, we fit polynomials of degree 1, 2, and 9 to bivariate data consisting of 10 data points. The degree 2 model (maroon) gives a significantly better fit than the degree 1 model (blue). The degree 10 model (orange) gives fits the data exactly, but at a glance we would guess it is overfit. That is, we don’t expect it to do a good job fitting the next data point we see. 5 0 y 10 In fact, we generated this data using a quadratic model, so the degree 2 model will tend to perform best fitting new data points. 0 2 4 6 8 10 x 4.2 R function lm As you would expect we don’t actually do linear regression by hand. Computationally, linear regression reduces to solving simultaneous equations, i.e. to matrix calculations. The R function lm can be used to fit any order polynomial to data. (lm stands for linear model). We will explore this in the next studio class. In fact lm can fit many types of functions besides polynomials, as you can explore using R help or google. 18.05 Class 26, Linear regression, Spring 2022 5 7 Multiple linear regression Data is not always bivariate. It can be trivariate or even of some higher dimension. Suppose we have data in the form of tuples (𝑦𝑖 , 𝑥1,𝑖 , 𝑥2,𝑖 , … 𝑥𝑚,𝑖 ) We can analyze this in a manner very similar to linear regression on bivariate data. That is, we can use least squares to fit the model 𝑦 = 𝛽0 + 𝛽1 𝑥1 + 𝛽2 𝑥2 + … + 𝛽𝑚 𝑥𝑚 . Here each 𝑥𝑗 is a predictor variable and 𝑦 is the response variable. For example, we might be interested in how a fish population varies with measured levels of several pollutants, or we might want to predict the adult height of a son based on the height of the mother and the height of the father. We don’t have time in 18.05 to study multiple linear regression, but we wanted you to see the name. 6 Least squares as a statistical model The linear regression model for fitting a line says that the value 𝑦𝑖 in the pair (𝑥𝑖 , 𝑦𝑖 ) is drawn from a random variable 𝑌𝑖 = 𝑎𝑥𝑖 + 𝑏 + 𝜀𝑖 where the ‘error’ terms 𝜀𝑖 are independent random variables with mean 0 and standard deviation 𝜎. The standard assumption is that the 𝜀𝑖 are i.i.d. with distribution 𝑁 (0, 𝜎2 ). So, the mean of 𝑌𝑖 is given by: 𝐸[𝑌𝑖 ] = 𝑎𝑥𝑖 + 𝑏 + 𝐸[𝜀𝑖 ] = 𝑎𝑥𝑖 + 𝑏. From this perspective, the least squares method chooses the values of 𝑎 and 𝑏 which minimize the sample variance about the line. In fact, under the assumption that 𝜀𝑖 ∼ 𝑁 (0, 𝜎2 ), the least square estimate (𝑎,̂ 𝑏)̂ coincides with the maximum likelihood estimate for the parameters (𝑎, 𝑏); that is, among all possible coefficients, (𝑎,̂ 𝑏)̂ are the ones that make the observed data most probable. 7 Regression to the mean The reason for the term ‘regression’ is that the predicted response variable 𝑦 will tend to be ‘closer’ to (i.e., regress to) its mean than the predictor variable 𝑥 is to its mean. Here closer is in quotes because we have to control for the scale (i.e. standard deviation) of each variable. The way we control for scale is to first standardize each variable. 𝑥 − 𝑥̄ 𝑢𝑖 = √𝑖 , 𝑠𝑥𝑥 𝑦𝑖 − 𝑦̄ . 𝑣𝑖 = √ 𝑠𝑦𝑦 18.05 Class 26, Linear regression, Spring 2022 8 Standardization changes the mean to 0 and variance to 1: 𝑢̄ = 𝑣̄ = 0, 𝑠𝑢𝑢 = 𝑠𝑣𝑣 = 1. The algebraic properties of covariance show 𝑠𝑢𝑣 = √ 𝑠𝑥𝑦 = 𝜌, 𝑠𝑥𝑥 𝑠𝑦𝑦 the correlation coefficient. Thus the least squares fit to 𝑣 = 𝑎𝑢 + 𝑏 has 𝑎̂ = 𝑠𝑢𝑣 =𝜌 𝑠𝑢𝑢 and 𝑏̂ = 𝑣 ̄ − 𝑎𝑢 ̂ ̄ = 0. So the least squares line is 𝑣 = 𝜌𝑢. Since 𝜌 is the correlation coefficient, it is between -1 and 1. Let’s assume it is positive and less than 1 (i.e., 𝑥 and 𝑦 are positively but not perfectly correlated). Then the formula 𝑣 = 𝜌𝑢 means that if 𝑢 is positive then the predicted value of 𝑣 is less than 𝑢. That is, 𝑣 is closer to 0 than 𝑢. Equivalently, 𝑦 − 𝑦̄ 𝑥 − 𝑥̄ < √ √ 𝑠𝑥𝑥 𝑠𝑦𝑦 i.e., 𝑦 regresses to 𝑦.̄ Notice how the standardization takes care of controlling the scale. Consider the extreme case of 0 correlation between 𝑥 and 𝑦. Then, no matter what the 𝑥 value, the predicted value of 𝑦 is always 𝑦.̄ That is, 𝑦 has regressed all the way to its mean. Note also that the regression line always goes through the point (𝑥,̄ 𝑦)̄. Example 10. Regression to the mean is important in longitudinal studies. Rice (Mathematical Statistics and Data Analysis) gives the following example. Suppose children are given an IQ test at age 4 and another at age 5 we expect the results will be positively correlated. The above analysis says that, on average, those kids who do poorly on the first test will tend to show improvement (i.e. regress to the mean) on the second test. Thus, a useless intervention might be misinterpreted as useful since it seems to improve scores. Example 11. Another example with practical consequences is reward and punishment. Imagine a school where high performance on an exam is rewarded and low performance is punished. Regression to the mean tells us that (on average) the high performing students will do slightly worse on the next exam and the low performing students will do slightly better. An unsophisticated view of the data will make it seem that punishment improved performance and reward actually hurt performance. There are real consequences if those in authority act on this idea. 8 Appendix We collect in this appendix a few things you might find interesting. You will not be asked to know these things for exams. 18.05 Class 26, Linear regression, Spring 2022 8.1 9 Proof of the formula for least square fit of a line The most straightforward proof is to use calculus. The sum of the squared errors is 𝑛 𝑆(𝑏, 𝑎) = ∑(𝑦𝑖 − 𝑎𝑥𝑖 − 𝑏)2 . 𝑖=1 Taking partial derivatives (and remembering that 𝑥𝑖 and 𝑦𝑖 are the data, hence constant) 𝑛 𝜕𝑆 = ∑ −2(𝑦𝑖 − 𝑎𝑥𝑖 − 𝑏) = 0 𝜕𝑏 𝑖=1 𝑛 𝜕𝑆 = ∑ −2𝑥𝑖 (𝑦𝑖 − 𝑎𝑥𝑖 − 𝑏) = 0 𝜕𝑎 𝑖=1 Summing this up we get two linear equations in the unknowns 𝑏 and 𝑎: (∑ 𝑥𝑖 ) 𝑎 + 𝑛𝑏 = ∑ 𝑦𝑖 (∑ 𝑥2𝑖 ) 𝑎 + (∑ 𝑥𝑖 ) 𝑏 = ∑ 𝑥𝑖 𝑦𝑖 Solving for 𝑎 and 𝑏 gives the formulas in Equation (1). A sneakier approach which avoids calculus is to standardize the data, find the best fit line, and then unstandardize. We omit the details. For a slew of applications across disciplines see: https://en.wikipedia.org/wiki/Linear_regression#Applications_of_linear_regression 8.2 Measuring the fit Once one computes the regression coefficients, it is important to check how well the regression model fits the data (i.e., how closely the best fit line tracks the data). A common but crude ‘goodness of fit’ measure is the coefficient of determination, denoted 𝑅2 . We’ll need some notation to define it. The total sum of squares is given by: TSS = ∑(𝑦𝑖 − 𝑦)̄2 . The residual sum of squares is given by the sum of the squares of the residuals. When fitting a line, this is: RSS = ∑(𝑦𝑖 − 𝑎̂ 𝑥𝑖 − 𝑏)̂2 . The RSS is the “unexplained” portion of the total sum of squares, i.e. unexplained by the regression equation. The difference TSS − RSS is the “explained” portion of the total sum of squares. The coefficient of determination 𝑅2 is the ratio of the “explained” portion to the total sum of squares: TSS − RSS 𝑅2 = . TSS In other words, 𝑅2 measures the proportion of the variability of the data that is accounted for by the regression model. A value close to 1 indicates a good fit, while a value close to 0 18.05 Class 26, Linear regression, Spring 2022 10 indicates a poor fit. In the case of simple linear regression, 𝑅2 is simply the square of the correlation coefficient between the observed values 𝑦𝑖 and the predicted values 𝑎𝑥𝑖 + 𝑏. Example 12. In the overfitting example (9), the values of 𝑅2 are: degree 1 2 9 𝑅2 0.3968 0.9455 1.0000 Notice the goodness of fit measure increases as 𝑛 increases. The fit is better, but the model also becomes more complex, since it takes more coefficients to describe higher order polynomials.