28 통계 분석 결과의 활용도를 높여주는 broom 패키지
전통적인 base R에서 통계 분석 결과(통계적 모델)를 활용하는 방법은 summary() 함수를 사용하는 것이다. 하지만 이 방법은 결과를 데이터프레임으로 변환할 수 없고, 결과를 정리하기도 어렵다. 다음 예를 보자.
28.1 전통적 통계 분석과 활용법
27 장에서 보았던 귀뚜라미(crickets) 데이터셋을 사용하여 귀뚜라미의 초당 소리 횟수(rate)와 온도(temp), 아종(species) 간의 관계를 분석한 예이다.
lm() 함수를 사용하여 선형 회귀 모델을 만들었다. 이 모델은 온도와 아종이 초당 소리 횟수에 미치는 영향을 분석한다. 이렇게 만들어진 모델 객체를 interaction_fit 변수에 저장했다. 보통 이 객체를 출력해 보거나, summary() 함수를 사용하여 결과를 확인한다.
attributes(interaction_fit)$names
[1] "coefficients" "residuals" "effects" "rank"
[5] "fitted.values" "assign" "qr" "df.residual"
[9] "contrasts" "xlevels" "call" "terms"
[13] "model"
$class
[1] "lm"이 가운데 잔차(residuals)와 회귀 계수(coefficients)를 확인해 보자.`
interaction_fit$residuals 1 2 3 4 5 6
0.91074559 -1.88925441 -1.69388557 -0.29388557 0.40611443 1.40611443
7 8 9 10 11 12
-1.44706949 -0.64706949 0.25293051 1.85293051 3.09974659 3.04887825
13 14 15 16 17 18
-3.70314789 -1.30314789 -0.81108570 -1.78029355 -1.38029355 0.61970645
19 20 21 22 23 24
-0.79077057 0.70922943 3.63303691 0.02255989 0.42255989 -1.64664796
25 26 27 28 29 30
2.52890568 1.16668250 0.48699763 -1.72347938 -0.82347938 -0.12347938
31
-0.51014892 interaction_fit$coefficients (Intercept) temp speciesO. niveus
-11.0408481 3.7514472 -4.3484072
temp:speciesO. niveus
-0.2339856 또는 이 객체를 summary() 함수를 사용하여 요약할 수 있다.
result_summary <- summary(interaction_fit)
result_summary
Call:
lm(formula = rate ~ temp * species, data = crickets)
Residuals:
Min 1Q Median 3Q Max
-3.7031 -1.3417 -0.1235 0.8100 3.6330
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -11.0408 4.1515 -2.659 0.013 *
temp 3.7514 0.1601 23.429 <2e-16 ***
speciesO. niveus -4.3484 4.9617 -0.876 0.389
temp:speciesO. niveus -0.2340 0.2009 -1.165 0.254
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 1.775 on 27 degrees of freedom
Multiple R-squared: 0.9901, Adjusted R-squared: 0.989
F-statistic: 898.9 on 3 and 27 DF, p-value: < 2.2e-16이렇게 만들어진 객체에 대해 여러 속성(attributes)을 확인할 수 있다. 예를 들어 회귀 계수, 결정 계수, 조정된 결정 계수, F-통계량 등을 확인할 수 있다.
result_summary$coefficients Estimate Std. Error t value Pr(>|t|)
(Intercept) -11.0408481 4.1514800 -2.6594969 1.300079e-02
temp 3.7514472 0.1601220 23.4286850 1.780831e-19
speciesO. niveus -4.3484072 4.9616805 -0.8763981 3.885447e-01
temp:speciesO. niveus -0.2339856 0.2008622 -1.1649059 2.542464e-01result_summary$r.squared[1] 0.9900871result_summary$adj.r.squared[1] 0.9889857result_summary$fstatistic value numdf dendf
898.9095 3.0000 27.0000 베이스 R에서는 t.test(), aov(), lm() 등 통계 분석 함수를 사용하여 모델을 만들고, 그 결과를 summary() 함수를 사용하여 요약하는 것이 일반적이다. 그런데 이렇게 하는 경우는 다음 단계의 분석이나 시각화 등 후속 작업을 하기에 불편한다. 이런 불편을 해소하기 위해 broom 패키지가 만들어졌다.
28.2 broom 패키지 소개
broom 패키지는 통계 분석 결과를 깔끔하게 정리하여 데이터프레임 형태로 변환해 주는 패키지이다. 이 패키지의 주 함수는 다음 3가지이다.
-
broom::tidy(): 통계 분석 결과를 데이터프레임 형태로 변환한다. -
broom::glance(): 통계 분석 결과의 요약 정보를 데이터프레임 형태로 변환한다. -
broom::augment(): 통계 분석 결과를 원본 데이터에 추가하여 데이터프레임 형태로 변환한다.
broom이 지원하는 통계 모델은 lm(), glm(), nls(), lmer(), glmer(), gam() 등 다양한 모델을 지원한다. 패키지 비니에트 “Available methods”를 참고하면 어떤 모델이 지원되는지 확인할 수 있다.
28.3 broom 패키지로 통계 분석 결과 정리하기
broom 패키지를 사용하여 통계 분석 결과를 정리해 보자. 먼저 tidy() 함수를 사용하여 회귀 계수와 잔차를 데이터프레임 형태로 변환한다.
tidy_result <- tidy(interaction_fit)
tidy_result# A tibble: 4 × 5
term estimate std.error statistic p.value
<chr> <dbl> <dbl> <dbl> <dbl>
1 (Intercept) -11.0 4.15 -2.66 1.30e- 2
2 temp 3.75 0.160 23.4 1.78e-19
3 speciesO. niveus -4.35 4.96 -0.876 3.89e- 1
4 temp:speciesO. niveus -0.234 0.201 -1.16 2.54e- 1tidy() 함수는 통계 분석 결과를 데이터프레임 형태로 변환하여, 각 회귀 계수와 그에 대한 통계적 검정 결과를 포함한다. 이 데이터프레임은 term, estimate, std.error, statistic, p.value 등의 열을 포함한다.
이제 glance() 함수를 사용하여 모델의 요약 정보를 데이터프레임 형태로 변환하자.
glance_result <- glance(interaction_fit)
glance_result# A tibble: 1 × 12
r.squared adj.r.squared sigma statistic p.value df logLik AIC BIC
<dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 0.990 0.989 1.78 899. 3.77e-27 3 -59.6 129. 136.
# ℹ 3 more variables: deviance <dbl>, df.residual <int>, nobs <int>glance() 함수는 모델의 요약 정보를 데이터프레임 형태로 변환하여, 결정 계수, 조정된 결정 계수, F-통계량, p-value 등의 정보를 포함한다. 이 데이터프레임은 모델의 전반적인 성능을 평가하는 데 유용하다.
이제 augment() 함수를 사용하여 원본 데이터에 통계 분석 결과를 추가한다.
augment_result <- augment(interaction_fit)
augment_result# A tibble: 31 × 9
rate temp species .fitted .resid .hat .sigma .cooksd .std.resid
<dbl> <dbl> <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 67.9 20.8 O. exclamationis 67.0 0.911 0.271 1.80 0.0336 0.601
2 65.1 20.8 O. exclamationis 67.0 -1.89 0.271 1.76 0.145 -1.25
3 77.3 24 O. exclamationis 79.0 -1.69 0.0966 1.77 0.0269 -1.00
4 78.7 24 O. exclamationis 79.0 -0.294 0.0966 1.81 0.000811 -0.174
5 79.4 24 O. exclamationis 79.0 0.406 0.0966 1.81 0.00155 0.241
6 80.4 24 O. exclamationis 79.0 1.41 0.0966 1.79 0.0186 0.833
7 85.8 26.2 O. exclamationis 87.2 -1.45 0.0730 1.78 0.0141 -0.847
8 86.6 26.2 O. exclamationis 87.2 -0.647 0.0730 1.80 0.00282 -0.379
9 87.5 26.2 O. exclamationis 87.2 0.253 0.0730 1.81 0.000431 0.148
10 89.1 26.2 O. exclamationis 87.2 1.85 0.0730 1.77 0.0232 1.08
# ℹ 21 more rowsaugment() 함수는 원본 데이터에 통계 분석 결과를 추가하여, 각 관측치에 대한 예측값, 잔차, 표준화된 잔차 등을 포함하는 데이터프레임을 생성한다. 이 데이터프레임은 모델의 적합도를 평가하고, 이상치를 식별하는 데 유용하다.