[R] 다중 선형 회귀 분석 (Multiple linear regression) - lm(), factor()

[R] 다중 선형 회귀 분석 (Multiple linear regression) - lm(), factor()

 

 지난 포스팅에서 연속형 변수의 단순 선형 회귀 분석 (Simple linear regression: 2022.12.22 - [선형 회귀 분석/R] - [R] 단순 선형 회귀 분석 (Simple linear regression) - lm()), 범주형 변수의 단순 선형 회귀 분석 (Simple linear regression with categorical variables:2022.12.22 - [선형 회귀 분석/R] - [R] 범주형 변수의 선형 회귀 분석 (Simple linear regression with categorical variable) - lm(), factor())에 대해 알아보았다.  위 분석들은 모두 독립 변수가 1개인 단순 선형 회귀 분석이다. 하지만 세상 일은 그렇게 단순하지 않다. 

 

교란 변수 (Confounder, confounding variable)

 우리는 앞서  연속형 변수의 단순 선형 회귀 분석 (Simple linear regression: 2022.12.22 - [선형 회귀 분석/R] - [R] 단순 선형 회귀 분석 (Simple linear regression) - lm()) 포스팅에서 수축기 혈압과 심혈관 질환의 관련성을 찾아보았다. 그런데, 비만인 사람은 혈압이 대체적으로 높은 것으로 알려져 있다. 더구나 비만인 사람은 심혈관 질환의 위험도도 높아지는 것으로 알려져 있다. 그렇다면, "수축기 혈압이 높았던 사람에게서 높은 심혈관 질환 위험 점수가 발견되었던 것은 그냥 그 사람이 비만도가 높아 당연히 수축기 혈압도 높고, 심혈관 질환 위험 점수도 높았던 것이 아닐까?"라는 의문이 든다. 즉, "우리가 이전 단순 선형 회귀 분석에서 발견한 과학적 현상(수축기 혈압이 높으면 심혈관 질환 위험 점수가 높아진다.)은 사실은 비만 때문이 아니었을까?"라는 의문이 들게 된다. 이렇게 독립변수와 종속변수 모두의 원인이 되는 변수를 교란 변수라고 하며, 이런 현상(교란)은 통계 분석을 할 때 반드시 교정을 해야 한다. (사실 독립변수는 인과적 관계가 아니라 연관 관계이기만 해도 충분하다.)

교란 현상의 해결 방법: 다중 회귀 분석

 교란 현상을 해결하는 방법은 보정(adjustment)과 층화(stratification)가 있다. 층화도 좋은 방법이고, 때로는 보정을 할 수 없어 층화를 해야만 할 때도 있지만, 보통 표본의 수가 부족해지는 문제가 있기 때문에, 층화가 아니라 보정을 시행하곤 한다. 보정을 시행할 수 있는 방법으로는 이전에 다룬 편상관 분석(2022.12.18 - [상관분석/R] - [R] 편상관 계수, 부분 상관 계수 (Partial correlation coefficient) - pcor.test())이 있다. 하지만 이는 상관 계수만 구할 수 있어 "예측할 수 있는 정도"만 표현할 수 있다. 그 대신 비만도의 영향을 보정했을 때, 수축기 혈압이 1단위 증가할 때 심혈관 질환 위험 점수가 얼마나 증가하는지(기울기)를 알 수 있는 방법으로는 본 포스팅에서 다룰 다중 선형 회귀 분석이 있다.

 

*실습용 데이터는 아래 링크를 클릭하면 다운로드할 수 있습니다.

2022.08.04 - [공지사항 및 소개] - 분석용 데이터 (update 22.12.18)

분석용 데이터 (update 22.12.18)

 

코드를 보여드리기에 앞서 워킹 디렉토리부터 지정하겠다.

워킹 디렉토리에 관한 설명은 다음 링크된 포스트에서 볼 수 있다.

2022.08.05 - [통계 프로그램 사용 방법/R] - [R] 작업 디렉토리 (Working Directory) 지정 - getwd(), setwd()

setwd("C:/Users/user/Documents/Tistory_blog")

 

데이터를 불러와 df에 객체로 저장하겠다.

데이터 불러오는 방법은 다음 링크에서 볼 수 있다.

2022.08.05 - [통계 프로그램 사용 방법/R] - [R] 데이터 불러오기 : EXCEL - read_excel(), read.xlsx()

2022.08.10 - [통계 프로그램 사용 방법/R] - [R] 데이터 저장하기 : CSV 파일 - write.csv(), write_csv()

2022.08.10 - [통계 프로그램 사용 방법/R] - [R] 데이터 불러오기 : SAS file (.sas7bdat) - read.sas7bdat(), read_sas()

install.packages("readr")
library("readr")
df<-read_csv("Data.csv")

 

목표:  BMI의 영향을 배제 혹은 교정했을 때 수축기 혈압(SBP)이 1단위 증가할 때 심혈관 질환 위험 점수(CVD_RISK)는 얼마나 증가하는가?

 

단순 선형 회귀 분석과 같이 lm()함수를 사용하며, 코드는 다음과 같다.

코드

LR1<-lm(CVD_RISK~SBP+BMI, data=df)
summary(LR1)

LR1<-lm(CVD_RISK~SBP+BMI, data=df) : df 데이터를 사용하여 선형 회귀 분석을 시행한다. 종속 변수는 CVD_RISK로, 독립 변수에는 SBP와 BMI를 사용하라. 그리고 그 결과를 LR1에 저장하라. SBP와 BMI의 순서가 바뀌는 것은 아무 상관이 없다. 즉, BMI를 보정했을 때의 SBP의 영향과 SBP를 보정했을 때의 BMI의 영향을 동시에 파악하는 것이다.
summary(LR1) : LR1을 요약하여 보여달라.

 

결과

Call:
lm(formula = CVD_RISK ~ SBP + BMI, data = df)

Residuals:
     Min       1Q   Median       3Q      Max 
-17.3967  -3.7096   0.0512   3.4044  15.9551 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 14.36106    1.83287   7.835 1.19e-14 ***
SBP          0.99389    0.01931  51.478  < 2e-16 ***
BMI          1.08243    0.17174   6.303 4.38e-10 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 5.079 on 997 degrees of freedom
Multiple R-squared:  0.9393,	Adjusted R-squared:  0.9392 
F-statistic:  7712 on 2 and 997 DF,  p-value: < 2.2e-16

해석의 방법은 단순 선형 회귀 분석에서와 같다. (2022.12.22 - [선형 회귀 분석/R] - [R] 단순 선형 회귀 분석 (Simple linear regression) - lm())

 

Coefficients:
                    Estimate   Std. Error t value     Pr(>|t|)    
(Intercept)  14.36106    1.83287   7.835  1.19e-14 ***
SBP              0.99389    0.01931  51.478   < 2e-16 ***
BMI               1.08243     0.17174   6.303 4.38e-10 ***

 

 BMI의 영향을 보정하더라도, SBP가 1단위 증가할 때 CVD_RISK는 0.99389만큼 증가하며 이에 대한 p-value는 $2 \times 10^{-16}$보다 작아 매우 유의함을 알 수 있다. 

 또한, SBP의 영향을 보정하더라도, BMI가 1단위 증가할 때 CVD_RISK는 1.08243만큼 증가하며 이에 대한 p-value는 $4.38 \times 10^{-10}$보다 작아 매우 유의함을 알 수 있다. 

 

즉, 단순 회귀 분석에서 나타난 수축기 혈압과 심혈관 질환 위험 점수 사이의 유의한 관련성은 BMI의 영향을 보정하더라도 유의함을 알 수 있다. 또한, 단순 회귀 분석에서의 회귀 계수  1.101935보다 작으므로, 일부분의 설명력은 BMI에게 빼았겼음을 알 수 있다.

 

 

보정 변수가 두 개 이상일 때 코드

 BMI외에 흡연 상태로도 보정하고 싶을 수 있다. 이럴 땐 어떻게 해야 할까? 답은, 단순히 BMI 뒤에 추가적인 보정변수를 '+'연산자를 이용하여 나열하면 된다. 그런데, 선형 회귀 분석 전에 연속형 변수와 범주형 변수의 전처리가 다르다는 것을 일전에 소개한 적이 있다. 범주형 변수는 분석 전에 범주형 자료임을 선언해주어야 한다. 내용은 다음 링크를 확인하면 된다.

연속형 변수: 2022.12.22 - [선형 회귀 분석/R] - [R] 단순 선형 회귀 분석 (Simple linear regression) - lm()

범주형 변수: 2022.12.22 - [선형 회귀 분석/R] - [R] 범주형 변수의 선형 회귀 분석 (Simple linear regression with categorical variable) - lm(), factor()

 그런데, 흡연 상태 (SMOK)는 범주형 자료다. 따라서 선형 회귀 분석을 시행하기 전에 범주형임을 선언해주어야 한다.

df$SMOK<-factor(df$SMOK)
LR2<-lm(CVD_RISK~SBP+BMI+SMOK, data=df)
summary(LR2)

 

결과

Call:
lm(formula = CVD_RISK ~ SBP + BMI + SMOK, data = df)

Residuals:
    Min      1Q  Median      3Q     Max 
-17.280  -3.663   0.133   3.454  16.141 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  15.4269     1.9747   7.812 1.42e-14 ***
SBP           0.9842     0.0205  48.009  < 2e-16 ***
BMI           1.0832     0.1717   6.310 4.21e-10 ***
SMOK1         0.1361     0.4189   0.325    0.745    
SMOK2         0.7824     0.5036   1.554    0.121    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 5.077 on 995 degrees of freedom
Multiple R-squared:  0.9394,	Adjusted R-squared:  0.9392 
F-statistic:  3859 on 4 and 995 DF,  p-value: < 2.2e-16

SBP, BMI는 유의한 결과를 보이고 있지만, SMOK는 1(과거 흡연자), 2(현재 흡연자) 모두 유의한 결과를 보이고 있지 않다. 

"BMI와 흡연 상태로 보정을 하였을 때, 수축기 혈압 1 단위 증가할 때, CVD_RISK는 0.9842점 증가한다"라고 할 수 있지만

"SBP와 BMI로 보정을 하였을 때, 비흡연자와 과거 흡연자 혹은 비흡연자와 현재 흡연자 사이의 CVD_RISK의 차이는 존재하지만 유의하다고 할 수 없다."라고 결론 내리게 된다.

 

 코드 정리

##워킹 디렉토리 지정
setwd("C:/Users/user/Documents/Tistory_blog")

##데이터 불러오기
install.packages("readr")
library("readr")
df<-read_csv("Data.csv")

##다중 선형 회귀 분석 1
LR1<-lm(CVD_RISK~SBP+BMI, data=df)
summary(LR1)

##다중 선형 회귀 분석 2
df$SMOK<-factor(df$SMOK)
LR2<-lm(CVD_RISK~SBP+BMI+SMOK, data=df)
summary(LR2)

 

[R] 다중 선형 회귀 분석 (Multiple linear regression) 정복 완료!

작성일: 2022.12.23.

최종 수정일: 2022.12.23.

이용 프로그램: R 4.2.2

RStudio v2022.07.2

RStudio 2022.07.2+576 "Spotted Wakerobin" Release

운영체제: Windows 10, Mac OS 12.6.1