hem rbc sex age diag
1 14.8 5.27 female 12.7 appendicitis
2 15.7 5.26 male 14.1 no appendicitis
3 11.4 3.98 female 14.1 no appendicitis
4 13.6 4.64 female 16.4 no appendicitis
5 12.6 4.44 female 11.1 appendicitis
6 12.5 4.96 male 11.1 no appendicitisSTAT 516 Lec 02
Simple linear regression
Karl Gregory
2026-10-26
Hemoglobin versus RBC count example
These data are a subset of the dataset Marcinkevičs et al. (2023)
Some outliers and missing values were removed.
Hemoglobin versus RBC
Hemoglobin level vs red blood cell count for \(n = 762\) children.
Simple linear regression
For \((x_1,Y_1),\dots,(x_n,Y_n)\), the simple linear regression model is \[ Y_i = \beta_0 + \beta_1 x_i + \varepsilon_i, \quad i = 1,\dots,n \] where
- \(x_1,\dots,x_n\) are the covariate or predictor values.
- \(Y_1,\dots,Y_n\) are the response values.
- \(\beta_0\) and \(\beta_1\) are the intercept and slope parameters, respectively.
- \(\varepsilon_1,\dots,\varepsilon_n\) are independent \(\text{Normal}(0,\sigma^2)\) error terms.
- \(\sigma^2\) is the error term variance.
Goals in simple linear regression
We will learn how to:
- Estimate the intercept and slope parameters \(\beta_0\) and \(\beta_1\).
- Estimate the error term variance \(\sigma^2\).
- Perform inference on \(\beta_1\).
- Build a confidence interval for \(\beta_0 + \beta_1 x_{\operatorname{new}}\) at any \(x_{\operatorname{new}}\).
- Build a prediction interval for \(Y\) at any \(x_{\operatorname{new}}\).
- Decompose the variation in \(Y\) into sums of squares.
- Check whether the model assumptions are satisfied.
- Identify outliers and understand their effects.
Least-squares estimation of slope and intercept
The squared error criterion given by the sum \[ Q(b_0,b_1) = \sum_{i=1}^n(Y_i - (b_0 + b_1 x_i))^2 \] of squared vertical distances of \(Y_i\) from the line \(y = b_0 + b_1x\).
We find \(Q(b_0,b_1)\) is minimized at \((b_0,b_1) = (\hat \beta_0,\hat\beta_1)\), where
- \(\hat \beta_0 = \bar Y_n - \hat \beta_1 \bar x_n\)
- \(\displaystyle \hat \beta_1 = \frac{\sum_{i=1}^n (x_i - \bar x_n)(Y_i-\bar Y_n) }{\sum_{i=1}^n(x_i - \bar x_n)^2}\)
provided \(\sum_{i=1}^n(x_i - \bar x_n)^2 > 0\).
The least-squares line or fitted line is the line \(y = \hat \beta_0 + \hat \beta_1x\).
Pearson’s correlation coefficient
Given \((x_1,Y_n),\dots,(x_n,Y_n)\), the quantity \[ r_{xY} = \frac{\sum_{i=1}^n (x_i - \bar x_n)(Y_i-\bar Y_n)}{\sqrt{\sum_{i=1}^n(x_i - \bar x_n)^2 \sum_{i=1}^n(Y_i - \bar Y_n)^2}} \] is called Pearson’s correlation coefficient.
- Describes strength and direction of linear relationships.
- Must satisfy \(r_{xY} \in [-1,1]\).
- Values close to zero indicate a weak linear relationship.
- Is related to \(\hat \beta_1\) by \[ \hat \beta_1 = r_{xY} \frac{S_{Y}}{S_{x}}. \] where \(S_Y\) and \(S_x\) are the sample std devs of the \(Y\) and \(x\) values.
Hemoglobin versus RBC
Find the least-squares line on the hemoglobin data.
Hemoglobin versus RBC
Estimating the error term variance
After obtaining \(\hat \beta_0\) and \(\hat \beta_1\), define the
- fitted values as \(\hat Y_i = \hat \beta_0 + \hat\beta_1 x_i\)
- residuals as \(\hat \varepsilon_i = Y_i - \hat Y_i\)
for \(i = 1,\dots,n\).
Then an unbiased estimator of \(\sigma^2\) is given by \[ \hat \sigma^2 = \frac{1}{n-2} \sum_{i=1}^n \hat \varepsilon_i^2. \]
Confidence interval for the slope parameter
Assume \(\varepsilon_1,\dots,\varepsilon_n \overset{\operatorname{ind}}{\sim}\text{Normal}(0,\sigma^2)\) and set \(S_{xx} = \sum_{i=1}^n(x_i - \bar x_n)^2\).
- The slope estimator \(\hat \beta_1\) is distributed as \[ \hat \beta_1 \sim \text{Normal}(\beta_1,\sigma^2/S_{xx}). \]
- “Studentizing” the above gives \[ \frac{\hat \beta_1 - \beta_1}{\hat \sigma/\sqrt{S_{xx}}} \sim t_{n-2}. \]
- So a \((1-\alpha)100\%\) confidence interval for \(\beta_1\) is \[ \hat \beta_1 \pm t_{n-2,\alpha/2} \hat \sigma / \sqrt{S_{xx}}. \]
- We often write \(\widehat{\text{s.e.}}(\hat \beta_1) = \hat \sigma / \sqrt{S_{xx}}\), s.e. for standard error.
Exercise
Justify the CI using the sampling distribution result.
Hemoglobin versus RBC
Obtain an estimate of the error term variance.
We obtain \(\hat \sigma^2 = 0.626\).
Now construct a 95% confidence interval for \(\beta_1\).
The 95% CI is \((2.011,2.314)\).
Tests of hypotheses about the slope
We most often test hypotheses about \(\beta_1\) of the form \[ \begin{array}{lclcl} H_0\text{: } \beta_1 \geq 0 &\text{ or } & H_0\text{: } \beta_1 = 0 &\text{ or } & H_0\text{: } \beta_1 \leq 0\\ H_1\text{: } \beta_1 < 0 && H_1\text{: } \beta_1 \neq 0 && H_1\text{: } \beta_1 > 0. \end{array} \] Reject or fail to reject \(H_0\) based on the value of the test statistic \[ T_{\operatorname{test}} = \frac{\hat \beta_1 }{\hat\sigma / \sqrt{S_{xx}}}. \] Rejection rules for the above at significance level \(\alpha\) are \[ \begin{array}{lclcl} T_{\operatorname{test}} < -t_{n-2,\alpha} &\text{ or } & |T_{\operatorname{test}}| > t_{n-2,\alpha/2}&\text{ or } & T_{\operatorname{test}} > t_{n-2,\alpha}. \end{array} \]
The corresponding p-values are, with \(T \sim t_{n-2}\), the probabilities \[ \begin{array}{lclcl} P(T < T_{\operatorname{test}}) &\text{ or } & 2 \times P(T > |T_{\operatorname{test}}|) & \text{ or } & P(T > T_{\operatorname{test}}). \end{array} \]
Hemoglobin versus RBC
Test the hypotheses \(H_0\): \(\beta_1 = 0\) vs \(H_1\): \(\beta_1 \neq 0\) at \(\alpha = 0.05\).
We get \(T_{\operatorname{test}} = 28.021\), \(t_{n-2,\alpha/2} = 1.963\), and \(p\)-value 0; we reject \(H_0\).
Test the hypotheses \(H_0\): \(\beta_1 \leq 2\) vs \(H_1\): \(\beta_1 > 2\).
We get \(T_{\operatorname{test}} = 2.107\), \(t_{n-2,\alpha} = 1.647\), and \(p\)-value 0.018; we reject \(H_0\).
Confidence interval for the height of the line
Assume \(\varepsilon_1,\dots,\varepsilon_n \overset{\operatorname{ind}}{\sim}\text{Normal}(0,\sigma^2)\). Then:
The estimator \(\hat \beta_0 + \hat \beta_1 x_{\operatorname{new}}\) is distributed as \[ \hat \beta_0 + \hat \beta_1 x_{\operatorname{new}} \sim \text{Normal}\left(\beta_0 + \beta_1 x_{\operatorname{new}},\sigma^2\left[\frac{1}{n} + \frac{(x_{\operatorname{new}} - \bar x_n)^2}{S_{xx}}\right]\right). \]
“Studentizing” the above gives \[ \frac{\hat \beta_0 + \hat \beta_1 x_{\operatorname{new}} - (\beta_0 + \beta_1 x_{\operatorname{new}})}{\hat \sigma \sqrt{ \frac{1}{n} + \frac{(x_{\operatorname{new}} - \bar x_n)^2}{S_{xx}}}} \sim t_{n-2}. \]
So a \((1-\alpha)100\%\) confidence interval for \(\beta_0 + \beta_1 x_{\operatorname{new}}\) is \[ \hat \beta_0 + \hat \beta_1 x_{\operatorname{new}} \pm t_{n-2,\alpha/2}\hat \sigma \sqrt{ \frac{1}{n} + \frac{(x_{\operatorname{new}} - \bar x_n)^2}{S_{xx}}}. \]
Hemoglobin versus RBC
Give a 95% CI for the mean hemoglobin level of individuals with RBC count 5.5.
We are \(95\%\) confident that the mean hemoglobin level of individuals with RBC count 5.5 lies in the interval \((14.767,15.011)\).
Prediction interval for a new value of the response
Assume \(\varepsilon_1,\dots,\varepsilon_n \overset{\operatorname{ind}}{\sim}\text{Normal}(0,\sigma^2)\). Then:
The new residual \(Y_{\operatorname{new}} - \hat \beta_0 + \hat \beta_1 x_{\operatorname{new}}\) is distributed as \[ Y_{\operatorname{new}} - \hat \beta_0 + \hat \beta_1 x_{\operatorname{new}} \sim \text{Normal}\left(0,\sigma^2\left[1 + \frac{1}{n} + \frac{(x_{\operatorname{new}} - \bar x_n)^2}{S_{xx}}\right]\right). \]
“Studentizing” the above gives \[ \frac{Y_{\operatorname{new}} - \hat \beta_0 + \hat \beta_1 x_{\operatorname{new}}}{\hat \sigma \sqrt{ 1 + \frac{1}{n} + \frac{(x_{\operatorname{new}} - \bar x_n)^2}{S_{xx}}}} \sim t_{n-2}. \]
So a \((1-\alpha)100\%\) prediction interval interval for \(Y_{\operatorname{new}}\) is \[ \hat \beta_0 + \hat \beta_1 x_{\operatorname{new}} \pm t_{n-2,\alpha/2}\hat \sigma \sqrt{1 + \frac{1}{n} + \frac{(x_{\operatorname{new}} - \bar x_n)^2}{S_{xx}}}. \]
Hemoglobin versus RBC
Give a 95% prediction interval for the hemoglobin level when the RBC count is 5.5.
We are \(95\%\) confident that an individual with RBC count 5.5 will have a hemoglobin level in the interval \((13.331,16.447)\).
Hemoglobin versus RBC
Plot confidence and prediction limits over the range of RBC counts.
alpha <- 0.05
ta2 <- qt(1-alpha/2,n-2)
xseq <- seq(min(x),max(x),length = 500)
xseq_se <- sqrt(sgsqhat)*sqrt(1/n+(xseq-xbar)^2/Sxx)
loci <- b0hat + b1hat * xseq - ta2 * xseq_se
upci <- b0hat + b1hat * xseq + ta2 * xseq_se
xseq_pse <- sqrt(sgsqhat)*sqrt(1+1/n+(xseq-xbar)^2/Sxx)
lopi <- b0hat + b1hat * xseq - ta2 * xseq_pse
uppi <- b0hat + b1hat * xseq + ta2 * xseq_pseHemoglobin versus RBC
The predict() function in R
Can obtain CI for \(\beta_0 + \beta_1x_\text{new}\) and PI for \(Y_{\text{new}}\) with predict() function.
Sums of squares in simple linear regression
We decompose the variation in \(Y_1,\dots,Y_n\) by defining the:
- Total sum of squares: \(\operatorname{SS}_{\operatorname{Tot}} = \sum_{i=1}^n(Y_i - \bar Y_n)^2\)
- Regression sum of squares: \(\operatorname{SS}_{\operatorname{Reg}} = \sum_{i=1}^n(\hat Y_i - \bar Y_n)^2\)
- Error sum of squares: \(\operatorname{SS}_{\operatorname{Error}} = \sum_{i=1}^n(Y_i - \hat Y_i)^2\)
We have \(\operatorname{SS}_{\operatorname{Tot}} = \operatorname{SS}_{\operatorname{Reg}} + \operatorname{SS}_{\operatorname{Error}}\).
The coefficient of determination is defined as \(\displaystyle R^2 = \frac{\operatorname{SS}_{\operatorname{Reg}}}{\operatorname{SS}_{\operatorname{Tot}}}\).
- \(R^2 \in [0,1]\)
- Proportion of variation in \(Y\) “explained” by the covariate \(x\).
- In simple linear regression we have \(R^2 = r_{xY}^2\).
The mean squares in simple linear regression
The SS, appropriately scaled, follow chi-square distributions:
- \(\displaystyle \frac{\operatorname{SS}_{\operatorname{Tot}}}{\sigma^2} \sim \chi^2_{n-1}(\phi_{\operatorname{Total}})\)
- \(\displaystyle \frac{\operatorname{SS}_{\operatorname{Reg}}}{\sigma^2} \sim \chi^2_1(\phi_{\operatorname{Reg}})\)
- \(\displaystyle \frac{\operatorname{SS}_{\operatorname{Error}}}{\sigma^2} \sim \chi^2_{n-2}\)
The quantities \(\phi_{\operatorname{Tot}}\) and \(\phi_{\operatorname{Reg}}\) are called noncentrality parameters.
Dividing \(\operatorname{SS}_{\operatorname{Reg}}\) and \(\operatorname{SS}_{\operatorname{Error}}\) by their dfs, we define:
- Regression mean square: \(\displaystyle \operatorname{MS}_{\operatorname{Reg}} = \frac{\operatorname{SS}_{\operatorname{Reg}}}{1}\)
- Error mean square: \(\displaystyle \operatorname{MS}_{\operatorname{Error}} = \frac{\operatorname{SS}_{\operatorname{Error}}}{n-2}\)
The Analysis of Variance (ANOVA) table
We often present the SS, df, and MS values in a table like this:
| Source | Df | SS | MS | F value | p-value |
|---|---|---|---|---|---|
| Regression | 1 | \(\operatorname{SS}_{\operatorname{Reg}}\) | \(\operatorname{MS}_{\operatorname{Reg}}\) | \(F_{\operatorname{test}}\) | \(P(F > F_{\operatorname{test}})\) |
| Error | \(n-2\) | \(\operatorname{SS}_{\operatorname{Error}}\) | \(\operatorname{MS}_{\operatorname{Error}}\) | ||
| Total | \(n-1\) | \(\operatorname{SS}_{\operatorname{Tot}}\) |
This is an example of an ANOVA table.
Overall F test
In addition to the SS, df, and MS value, the ANOVA table presents
- \(\displaystyle F_{\operatorname{test}} = \frac{\operatorname{MS}_{\operatorname{Reg}}}{\operatorname{MS}_{\operatorname{Error}}}\)
- \(P(F > F_{\operatorname{test}})\), where this is computed under \(F \sim F_{1,n-2}\)
These are the test statistic and p-value of the overall F test.
In simple linear regression this p-value is the same as the one for testing \(H_0\): \(\beta_1 = 0\) versus \(H_1\): \(\beta_1 \neq 0\) with the t test; moreover \(F_{\operatorname{test}} = T_{\operatorname{test}}^2\).
For what it’s worth, one can show \(\displaystyle F_{\operatorname{test}} = \frac{(n-2)r_{xY}^2}{1-r_{xY}^2}\) in SLR.
We will discuss the overall F test in greater detail later.
Building the ANOVA table
| Source | Df | SS | MS | F value | p-value |
|---|---|---|---|---|---|
| x | 1 | 491.37 | 491.37 | 785.15 | 0 |
| Error | 760 | 475.63 | 0.63 | ||
| Total | 761 | 967 |
Hemoglobin versus RBC
Demonstrate the lm(), summary(), and anova() functions in R
Hemoglobin versus RBC
Call:
lm(formula = Y ~ x)
Residuals:
Min 1Q Median 3Q Max
-5.970 -0.423 0.007 0.464 2.379
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2.9945 0.3706 8.08 2.6e-15 ***
x 2.1626 0.0772 28.02 < 2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.791 on 760 degrees of freedom
Multiple R-squared: 0.508, Adjusted R-squared: 0.507
F-statistic: 785 on 1 and 760 DF, p-value: <2e-16Hemoglobin versus RBC
Invariance of tests to shifting/scaling of data values
Our test statistics and p-values do not depend on the units in which the data are recorded:
Hemoglobin measured in grams per deciliter. What if we convert to grams per liter?
RBC measured in millions of cells per microliter. What if we convert to thousands of cells per microliter?
And what if we center the RBC values around their mean?
Hemoglobin versus RBC
Call:
lm(formula = hem2 ~ rbc2)
Residuals:
Min 1Q Median 3Q Max
-59.70 -4.23 0.07 4.64 23.79
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 1.33e+02 2.87e-01 466 <2e-16 ***
rbc2 2.16e+04 7.72e+02 28 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 7.91 on 760 degrees of freedom
Multiple R-squared: 0.508, Adjusted R-squared: 0.507
F-statistic: 785 on 1 and 760 DF, p-value: <2e-16Checking model assumptions
Validity of the foregoing analyses depends on these assumptions:
The responses are normally distributed around the regression line (Check QQ plot of residuals). If \(n\) is large this only matters when making prediction intervals.
The response has the same variance for all values of the covariate (Check residuals vs fitted values plot).
The covariate and the response are linearly related (Check residuals vs fitted values plot).
The response values are independent of each other (No way to check; must trust experimental design).
Hemoglobin versus RBC
Can generating diagnostic plots from lm() with plot()
Hemoglobin versus RBC
The abalone
Photo on the left by Sharktopus - Own work, CC BY-SA 3.0, https://commons.wikimedia.org/w/index.php?curid=14082271
Abalone data example
Predict shucked weight of an abalone by its length.
There are \(n = 4176\) records. Data come from Nash and Ford (1995).
Abalone
Abalone
Try transforming x:
Abalone
Abalone
Could transform both Y and x:
Abalone
Abalone
Abalone
Call:
lm(formula = logY ~ logx)
Residuals:
Min 1Q Median 3Q Max
-0.7775 -0.1244 -0.0109 0.1094 2.7749
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.80387 0.00827 97.2 <2e-16 ***
logx 3.08666 0.01132 272.8 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.198 on 4174 degrees of freedom
Multiple R-squared: 0.947, Adjusted R-squared: 0.947
F-statistic: 7.44e+04 on 1 and 4174 DF, p-value: <2e-16Interpreting the slope with a log response
Suppose the fitted simple linear regression model gives \[ \widehat{\log Y_i} = \hat \beta_0 + \hat \beta_1 x_i, \quad i=1,\dots,n. \]
Let \(\hat Y_i = \exp(\widehat{\log Y_i})\) be the back-transformed fitted values.
If \(x_2 = x_1+1\) then \(\log \hat Y_2 - \log \hat Y_1 = \hat \beta_1\), which gives \[ \hat Y_2 = e^{\hat \beta_1} \hat Y_1. \]
So an increase in the covariate by one unit changes the back-transformed fitted value by a factor of \(e^{\hat \beta_1}\).
Interpreting the slope with a log response and log covariate
Suppose the fitted simple linear regression model gives \[ \widehat{\log Y_i} = \hat \beta_0 + \hat \beta_1 \log x_i, \quad i=1,\dots,n. \]
Let \(\hat Y_i = \exp(\widehat{\log Y_i})\) be the back-transformed fitted values.
If \(x_2 = x_1 c\) for some \(c>0\), then \[ \log \hat Y_2 - \log \hat Y_1 = \hat \beta_1 (\log (x_1c) - \log x_1) \iff \hat Y_2 = c^{\hat \beta_1}\hat Y_1. \]
So a change in the covariate by a factor of \(c>0\) changes the back-transformed fitted value by a factor of \(c^{\hat \beta_1}\).
Further if \(x_2 = x_1(1+\delta)\) with \(\delta\) close to zero, we have \[ \hat Y_2 \approx (1 + \hat \beta_1 \delta)\hat Y_1 \] which we can obtain by using \(\log(x) \approx x - 1\) for \(x\) close to \(1\).
Thus a 1% change (\(\delta = 1/100\)) in the covariate leads to a \(\hat \beta_1\%\) change (\(\hat \beta_1\delta = \hat \beta_1 / 100\)) in the back-transformed fitted value.
Interpreting the slope with a log covariate
Suppose the fitted simple linear regression model gives \[ \hat Y_i = \hat \beta_0 + \hat \beta_1 \log x_i, \quad i=1,\dots,n. \]
If \(x_2 = x_1 (1+\delta)\) for some small \(\delta\), then \[ \hat Y_2 - \hat Y_1 = \hat \beta_1 \log(1+\delta) \approx \hat \beta_1 \delta \] using \(\log(x) = x - 1\) for \(x\) close to \(1\).
Rearranging the above gives \(\hat Y_2 = \hat Y_1 + \hat \beta_1 \delta\).
So a change in the covariate by \(1\%\) changes the back-transformed fitted value by the amount \(\hat \beta_1/100\).
Abalone
For the abalone data we the fitted model is \[ \log (\texttt{Shucked\_Wt}) = 0.804 + 3.087 \log (\texttt{Diameter}) \]
According to the fitted model:
Doubling the diameter is expected to change the shucked weight by a factor of \(2^{3.087} = 8.495\).
Increasing the diameter by \(1%\) is expected to increase the shucked weight by \(3.087\%\).
Abalone
We must back-transform prediction intervals if we have transformed \(Y\).
Abalone
plot(Y~x,col="gray"); logx <- seq(min(logx),max(logx),length=500)
newdata <- data.frame(logx = logx)
logy_hat <- predict(lm3,newdata = newdata,int = "pred")
lines(exp(logy_hat[,1]) ~ exp(logx), col = "red")
lines(exp(logy_hat[,2]) ~ exp(logx), col = "red", lty = 3,lwd=1.5)
lines(exp(logy_hat[,3]) ~ exp(logx), col = "red", lty = 3,lwd=1.5)
Outliers in simple linear regression
Outlying data points can have a large influence on the estimated regression function.
Let’s generate some data and then add an outlier:
The red data point appears to exert an undue influence over the fit.
Leverage
The leverage of a point \((x_i,Y_i)\) among \((x_1,Y_1),\dots,(x_n,Y_n)\) is \[ \text{lev}_i = \frac{1}{n} + \frac{(x_i - \bar x_n)^2}{S_{xx}} \]
Leverage only shows outlying-ness in the \(x\) direction.
Least-squares line must pass through the point \((\bar x_n,\bar Y_n)\).
Greater leverage means greater influence on the least-squares line.
Cook’s distance
Cook’s Distance measures how much each data point changes the fit: \[ D_i = \frac{1}{2 \hat \sigma^2}\sum_{j = 1 }^n(\hat Y_j - \hat Y_{j(i)})^2 \quad \text{ for } i = 1,\dots,n, \] where \(\hat Y_{j(i)}\) is the \(j\)th fitted value from the model fitted without obs \(i\).
Can also write \(\displaystyle D_i = \frac{\hat e_i^2}{2\hat \sigma^2}\frac{\text{lev}_i}{(1 - \text{lev}_i)^2}\) for \(i = 1,\dots,n\).
Make a plot of the Cook’s distances
Code to compute Cook’s distances
n <- length(Y)
xbar <- mean(x)
Ybar <- mean(Y)
rxY <- cor(x,Y) # Pearson's correlation coefficient
b1hat <- rxY * sd(Y)/sd(x)
b0hat <- Ybar - b1hat * xbar
Sxx <- sum((x - xbar)^2)
lev <- 1/n + (x - xbar)^2/Sxx
Yhat <- b0hat + b1hat * x
ehat <- Y - Yhat
sgsqhat <- sum(ehat^2)/(n-2)
cooksD <- ehat^2 / (2*sgsqhat) * lev / (1 - lev)^2Two-sample t-test by simple linear regression
Let \(Y_{ij} \overset{\operatorname{ind}}{\sim}\text{Normal}(\mu_i,\sigma^2)\), \(j = 1,\dots,n_i\), \(i=1,2\) and consider \[ \text{$H_0$: $\mu_2 - \mu_1 = 0$ versus $H_1$: $\mu_2 - \mu_1 \neq 0$.} \] The (equal-variances) two-sample t-test uses the test statistic \[ T_{\operatorname{test}}= \frac{\bar Y_2 - \bar Y_1}{S_{\operatorname{pooled}}\sqrt{\frac{1}{n_1} + \frac{1}{n_2}}}, \]
where \(\bar y_i =n_i^{-1}\sum_{j=1}^{n_i}Y_{ij}\), \(i=1,2\) and \[ S_{\operatorname{pooled}}^2 = \frac{(n_1 - 1)S_1^2 + (n_2 - 1)S_2^2}{n_1 + n_2 - 2}, \quad S_i = \frac{1}{n_i - 1}\sum_{j=1}^{n_i}(Y_{ij} - \bar Y_i)^2. \] We reject \(H_0\) at significance level \(\alpha\) if \(|T_{\operatorname{test}}| > t_{n-2,\alpha/2}\).
Appendicitis example
Look again at the data from Marcinkevičs et al. (2023).
hem rbc sex age diag
1 14.8 5.27 female 12.7 appendicitis
2 15.7 5.26 male 14.1 no appendicitis
3 11.4 3.98 female 14.1 no appendicitis
4 13.6 4.64 female 16.4 no appendicitis
5 12.6 4.44 female 11.1 appendicitis
6 12.5 4.96 male 11.1 no appendicitisIs the mean hemaglobin level the same in children with and without appendicitis (ignoring rbc, age, and sex)?
Appendicitis
Appendicitis
Two Sample t-test
data: data$hem by data$diag
t = -0.5, df = 760, p-value = 0.6
alternative hypothesis: true difference in means between group appendicitis and group no appendicitis is not equal to 0
95 percent confidence interval:
-0.204 0.122
sample estimates:
mean in group appendicitis mean in group no appendicitis
13.3 13.4 Appendicitis
Let the \(Y_i\) be the hemaglobin values and define an indicator variable as \[ x_i = \left\{\begin{array}{ll} 0& \text{if no appendicitis}\\ 1& \text{if appendicitis} \end{array}\right.\quad \text{ for } i = 1,\dots,n. \] Then in the SLR model \(Y_i = \beta_0 + \beta_1 x_i + \varepsilon_i\) we have
- \(\beta_0 = \mu_{\operatorname{no~app}}\)
- \(\beta_0 + \beta_1 = \mu_{\operatorname{app}}\)
- \(\beta_1 = \mu_{\operatorname{app}} - \mu_{\operatorname{no~app}}\)
The t test in the simple linear regression setup of \[ \text{$H_0$: $\beta_1 = 0$ versus $H_1$: $\beta_1 \neq 0$} \] will give the same p value as the equal-variances two-sample t test of \[ \text{$H_0$: $\mu_{\operatorname{app}} - \mu_{\operatorname{no~app}} = 0$ versus $H_1$: $\mu_{\operatorname{app}} - \mu_{\operatorname{no~app}} \neq 0$.} \quad \text{ Cool!} \]
Exercise
Show that in the above setup we have
\[ \frac{\bar Y_2 - \bar Y_1}{S_{\operatorname{pooled}}\sqrt{\frac{1}{n_1} + \frac{1}{n_2}}} = \frac{\hat \beta_1}{\hat \sigma /\sqrt{S_{xx}}}. \] Do it in steps, showing:
- \(\hat \beta_0 = \bar Y_1\)
- \(\hat \beta_1 = \bar Y_2 - \bar Y_1\)
- \(\hat \sigma = S_{\operatorname{pooled}}\)
- \(1 / \sqrt{S_{xx}} = \sqrt{\frac{1}{n_1} + \frac{1}{n_2}}\)
Appendicitis
Prepare the data:
Appendicitis
Call:
lm(formula = Y ~ x)
Residuals:
Min 1Q Median 3Q Max
-5.132 -0.732 -0.032 0.668 4.168
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 13.3732 0.0637 209.78 <2e-16 ***
x -0.0409 0.0830 -0.49 0.62
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 1.13 on 760 degrees of freedom
Multiple R-squared: 0.000319, Adjusted R-squared: -0.000997
F-statistic: 0.242 on 1 and 760 DF, p-value: 0.623Appendicitis
Automatic if we designate the predictor as a “factor” (but watch sign!).
Call:
lm(formula = data$hem ~ x)
Residuals:
Min 1Q Median 3Q Max
-5.132 -0.732 -0.032 0.668 4.168
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 13.3323 0.0532 250.49 <2e-16 ***
xno appendicitis 0.0409 0.0830 0.49 0.62
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 1.13 on 760 degrees of freedom
Multiple R-squared: 0.000319, Adjusted R-squared: -0.000997
F-statistic: 0.242 on 1 and 760 DF, p-value: 0.623Importance (or unimportance) of Normal error terms
Generate some data with a right-skewed error term distribution:
Importance (or unimportance) of Normal error terms
Check the Normal quantile-quantile plot of the residuals.
Importance (or unimportance) of Normal error terms
Generate a large number of such data sets and obtain \(\hat \beta_1\) for each one.
The importance (or unimportance) of Normal error terms
Check if the \(\hat \beta_1\) values have a Normal distribution.
References
Marcinkevičs, Ričards, Patricia Reis Wolfertstetter, Ugne Klimiene, Ece Ozkan, Kieran Chin-Cheong, Alyssia Paschke, Julia Zerres, et al. 2023. “Regensburg Pediatric Appendicitis Dataset.” Zenodo. https://doi.org/10.5281/zenodo.7711412.
Nash, Sellers, Warwick, and Wes Ford. 1995. “Abalone.” UCI Machine Learning Repository. https://doi.org/10.24432/C55C7W.