Running a Multiple Linear Regression
Jacob Barden
2024-06-12
1 Loading Libraries
library(psych) # for the describe() command
library(car) # for the vif() command## Loading required package: carData##
## Attaching package: 'car'## The following object is masked from 'package:psych':
##
## logitlibrary(sjPlot) # to visualize our results2 Importing Data
# import the dataset you cleaned previously
# this will be the dataset you'll use throughout the rest of the semester
# use EAMMi2 data
d <- read.csv(file="data/mydata.csv", header=T)3 State Your Hypothesis
We hypothesize that Extraversion, Agreeableness (measured by the BFI), and self-esteem (measured by the RSE) will significantly predict one’s rating of their perceived anxiety in relation to the pandemic.
4 Check Your Variables
# you only need to check the variables you're using in the current analysis
# although you checked them previously, it's always a good idea to look them over again and be sure that everything is correct
str(d)## 'data.frame': 1009 obs. of 6 variables:
## $ trans : chr "no" "no" "no" "no" ...
## $ age : chr "1 under 18" "1 under 18" "4 between 36 and 45" "4 between 36 and 45" ...
## $ big5_agr : num 4.33 6.67 4.67 6.67 5.33 ...
## $ big5_ext : num 1.67 6 5 5.67 4 ...
## $ rse : num 1.6 3.9 1.7 2.4 1.8 3.5 3 3.5 2.5 3.4 ...
## $ pas_covid: num 4.56 3.33 4.22 3.56 4.56 ...cont <- na.omit(subset(d, select=c(big5_agr, big5_ext, rse, pas_covid)))
cont$big5_agr <- scale(cont$big5_agr, center=T, scale=T)
cont$big5_ext <- scale(cont$big5_ext, center=T, scale=T)
cont$rse <- scale(cont$rse, center=T, scale=T)
# you can use the describe() command on an entire dataframe (d) or just on a single variable
describe(cont)## vars n mean sd median trimmed mad min max range skew
## big5_agr 1 1009 0.00 1.00 -0.01 0.04 0.88 -3.57 1.77 5.34 -0.39
## big5_ext 2 1009 0.00 1.00 -0.01 0.03 1.02 -2.30 1.83 4.14 -0.23
## rse 3 1009 0.00 1.00 0.07 0.02 1.04 -2.32 1.89 4.21 -0.20
## pas_covid 4 1009 3.24 0.68 3.22 3.26 0.66 1.00 5.00 4.00 -0.23
## kurtosis se
## big5_agr -0.03 0.03
## big5_ext -0.77 0.03
## rse -0.74 0.03
## pas_covid 0.12 0.02# also use histograms to examine your continuous variables
hist(cont$big5_agr)
hist(cont$big5_ext)
hist(cont$rse)
hist(cont$pas_covid)
# last, use scatterplots to examine your continuous variables together
plot(cont$pas_covid, cont$big5_agr)
plot(cont$big5_agr, cont$big5_ext)
plot(cont$big5_agr, cont$rse)
plot(cont$pas_covid, cont$big5_ext)
plot(cont$pas_covid, cont$rse)
plot(cont$big5_ext, cont$rse)
5 View Your Correlations
corr_output_m <- corr.test(cont)
corr_output_m## Call:corr.test(x = cont)
## Correlation matrix
## big5_agr big5_ext rse pas_covid
## big5_agr 1.00 0.08 0.23 0.02
## big5_ext 0.08 1.00 0.38 -0.08
## rse 0.23 0.38 1.00 -0.30
## pas_covid 0.02 -0.08 -0.30 1.00
## Sample Size
## [1] 1009
## Probability values (Entries above the diagonal are adjusted for multiple tests.)
## big5_agr big5_ext rse pas_covid
## big5_agr 0.00 0.03 0 0.50
## big5_ext 0.01 0.00 0 0.03
## rse 0.00 0.00 0 0.00
## pas_covid 0.50 0.01 0 0.00
##
## To see confidence intervals of the correlations, print with the short=FALSE option6 Run a Multiple Linear Regression
# use this commented out section only if you need to remove outliers
# to drop a single outlier, remove the # at the beginning of the line and use this code:
# d <- subset(d, row_id!=c(1108))
# to drop multiple outliers, remove the # at the beginning of the line and use this code:
# d <- subset(d, row_id!=c(1108) & row_id!=c(602))
# use the lm() command to run the regression
# dependent/outcome variable on the left, independent/predictor variables on the right
reg_model <- lm(pas_covid ~ big5_agr + big5_ext + rse, data = cont)7 Check Your Assumptions
7.1 Multiple Linear Regression Assumptions
- Observations should be independent
- Number of cases should be adequate (N ≥ 80 + 8m, where m is the number of IVs)
- Independent variables should not be too correlated (aka multicollinearity)
- Variables should be continuous and normally distributed
- Outliers should be identified and removed
- Relationship between the variables should be linear
- Residuals should be normal and have constant variance
needed <- 80 + 8*3
nrow(cont) >= needed## [1] TRUE7.2 Check multicollinearity
- Higher values indicate more multicollinearity
- Cutoff is usually 5
For your homework, you will need to discuss multicollinearity and any high values, but you don’t have to drop any variables.
vif(reg_model)## big5_agr big5_ext rse
## 1.057542 1.168509 1.2276377.3 Check linearity with Residuals vs Fitted plot
The plot of Residuals vs Fitted fits well to what would be considered a “good” graph. The red line is close to a slope of zero right down the center.
plot(reg_model, 1)
7.4 Check for outliers using Cook’s distance and a Residuals vs Leverage plot
The graphs below show Cook’s distance and residuals vs leverage. As seen in the graphs there are no significant outliers that would impact our results in a way that would require us to account for them.
# Cook's distance
plot(reg_model, 4)
# Residuals vs Leverage
plot(reg_model, 5)
7.5 Check homogeneity of variance in a Scale-Location plot
The scale location graph below shows a nearly flat line showing that the variance in the residuals doesn’t change as a function of x. There are no outliers skewing the slope of the line.
plot(reg_model, 3)
7.6 Issues with My Data
Before interpreting our results, we assessed our variables to see if they met the assumptions for a multiple linear regression. We analyzed a Scale-Location plot and there were no outliers that significantly affected our plots. The Cook’s distance plot and Residuals vs fitted also looked as we would expect for a data set that meets the needed requirements.
8 View Test Output
summary(reg_model)##
## Call:
## lm(formula = pas_covid ~ big5_agr + big5_ext + rse, data = cont)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.38224 -0.41286 0.02443 0.41259 1.84929
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.23885 0.02041 158.694 < 2e-16 ***
## big5_agr 0.06710 0.02100 3.196 0.00144 **
## big5_ext 0.03064 0.02207 1.388 0.16541
## rse -0.23562 0.02262 -10.414 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.6483 on 1005 degrees of freedom
## Multiple R-squared: 0.1036, Adjusted R-squared: 0.1009
## F-statistic: 38.72 on 3 and 1005 DF, p-value: < 2.2e-16# note for section below: to type lowercase Beta below (ß) you need to hold down Alt key and type 225 on numeric keypad. If that doesn't work you should be able to copy/paste it from somewhere else9 Write Up Results
To test our hypothesis that Agreeableness (measured by the BFI), Extraversion (measured by the BFI), and Self-Esteem (measured by RSE) would significantly predict self ratings of Covid related Anxiety (measured by the PAS), we used a multiple linear regression to model the relationship between the variables. Our data met the assumptions of a linear regression and had no issues with homogeneity of variance or linearity.
Our model was statistically significant, Adj. R2 = .104, F(3,1005) = 38.72, p < .001. The relationship between Extraversion and pandemic related anxiety was trivial (p=.165). Agreeableness was positively and significantly correlated with Anxiety(P=.001) and Self-esteem had a significant negative correlation with Anxiety (P<.001)
| Subjective Wellbeing (SWLS) | ||||
|---|---|---|---|---|
| Predictors | Estimates | SE | CI | p |
| Intercept | 3.24 | 0.02 | 3.20 – 3.28 | <0.001 |
| Agreeableness (BFI) | 0.07 | 0.02 | 0.03 – 0.11 | 0.001 |
| Extraversion (BFI) | 0.03 | 0.02 | -0.01 – 0.07 | 0.165 |
| Self-Esteem (RSE) | -0.24 | 0.02 | -0.28 – -0.19 | <0.001 |
| Observations | 1009 | |||
| R2 / R2 adjusted | 0.104 / 0.101 | |||
References
Cohen J. (1988). Statistical Power Analysis for the Behavioral Sciences. New York, NY: Routledge Academic.