R을 사용한 MANOVA(다변량 분산 분석)

MANOVA(다변량 분산 분석)란 무엇입니까?

  • MANOVA는 독립 변수의  여러 그룹 간의 차이를 분석하기 위해  두 개 이상의 종속 변수 를  포함하는 일  변량 ANOVA 의 확장  입니다.
  • 개별 그룹 평균을 비교하는 ANOVA와 달리 MANOVA는 각 종속 변수의 그룹 평균을 포함하는 벡터를 비교합니다.
  • MANOVA는 옴니버스 Wilk’s Lambda, Pillai’s Trace(가정에서 벗어나는 데 가장 강력함), Roy’s Largest Root 또는 Hotelling-Lawley의 검정을 사용하며, 이는 나중에  그룹 차이의 중요성을 평가하기 위해 F 통계량으로 변환됩니다  . Pillai’s Trace는 통계적 검정력 이 가장 높습니다  .
  • MANOVA는 다중 반응 변수의 최상의 선형 조합을 기반으로 그룹 내보다 그룹 간 차별을 최대화합니다.

MANOVA가 유용한 이유는 무엇입니까?

  • MANOVA는 여러 종속변수를 고려하여 독립변수에서 그룹 간의 차이를 분석합니다. 이것은  각 종속 변수에 대해 별도의 일변량 ANOVA를 수행하여 부풀릴 수 있는 제1종 오류 를 줄  입니다.
  • MANOVA는 또한 데이터 세트의 여러 종속 변수 간의 상호 상관 관계를 제어합니다.
  • 일변량 ANOVA와 비교할 때 MANOVA는 종속 변수의 더 많은 정보를 사용합니다. 즉, MANOVA는 여러 종속 변수의 결합된 정보를 기반으로 그룹 간의 차이를 찾을 수 있습니다.

그룹 간의 차이를 수량화하는 데 관심이 있는 종속 변수가 하나만 있는 경우 MANOVA가 필요하지 않습니다.

MANOVA의 가정

MANOVA  는 관측치의 독립성과 분산의 동질성에 대해 ANOVA에서 와 유사한 가정을 따릅니다. 

또한 MANOVA는 다음 가정을 충족해야 합니다.

  • 다변량 정규성: 데이터 또는 잔차는 독립 및 종속 변수의 각 조합에 대해 다변량 정규 분포를 가져야 합니다(단변량 정규성에 대해서는 Shapiro-Wilk 검정으로 확인하고 다변량 정규성에 대해서는 Mardia의 왜도 및 첨도를 확인함).
  • 분산-공분산 행렬의 동질성: 데이터는 독립 변수의 각 그룹에 의해 형성된 각 조합에 대해 등분산-공분산 행렬을 가져야 합니다. 이것은 일 변량 ANOVA 에서 확인되는 분산의 동질성의 다변량 버전입니다  . Box의 M 테스트를 사용하여 테스트할 수 있습니다. Box의 M-검정은 검정력이 거의 없으며  유의성에 대한 p 값 을 평가하기 위해 0.001과 같은 더 낮은 알파 수준을 사용합니다  .
  • 다중 공선성 :  종속변수 간에 다중공선성 이 없어야 함  (매우 높은 상관관계, 즉 > 0.9)
  • 선형성: 종속 변수는 독립 변수의 각 그룹에 대해 선형적으로 관련되어야 합니다. 종속변수가 2개 이상인 경우 종속변수 쌍은 선형으로 관련되어야 합니다.
  • 종속변수(다변량 이상값)에 이상값이 없어야 함(Mahalanobis 거리로 확인)
  • 종속 변수는 연속적이어야 합니다.

MANOVA 가설

  • 귀무 가설 : 그룹 평균 벡터는 모든 그룹에 대해 동일합니다.
  • 대립 가설 : 그룹 평균 벡터는 모든 그룹에 대해 동일하지 않습니다.

R의 일원(일 요인) MANOVA

일원 분산 분석에는 두 개 이상의 그룹과 두 개 이상의 종속 변수가 있는 하나의 독립 변수가 있습니다.

예시 데이터세트

다양한 식물 품종(

plant_var

) 및 식물 높이에 대한 관련 표현형 측정(

) 및 캐노피 볼륨(

캐노피_볼

). MANOVA를 사용하여 식물 높이와 캐노피 부피가 다른 식물 품종과 연관되어 있는지 확인하고 싶습니다.

MANOVA의 경우 데이터 세트는 종속 변수의 수보다 독립 변수의 그룹당 더 많은 관측치를 포함해야 합니다. Box의 M-검정을 사용하여 분산-공분산 행렬의 동질성을 테스트하는 데 특히 중요합니다.

데이터 세트 로드,

라이브러리 ( 타이디버스 )
머리 ( df, 2 )
# 출력
plant_var 높이 canopy_vol
< chr > < dbl > < dbl >
1 20 0.7
2 22 0.8

데이터 세트의 요약 통계 및 시각화

각 종속 변수를 기반으로 요약 통계를 가져오고,

# 종속변수 높이에 대한 요약 통계
df % > % group_by ( plant_var ) % > % 요약 ( n = n () , 평균 = 평균 ( 높이 ) , sd = sd ( 높이 ))
# 출력
plant_var n 평균 SD
< chr > < 정수 > < dbl > < dbl >
1 10 18.9 2.92
2 B 10 16.5 1.92
3 C 10 3.05 1.04
4 D 10 9.35 2.11
# 종속변수 canopy_vol에 대한 요약 통계
df % > % group_by ( plant_var ) % > % 요약 ( n = n () , 평균 = 평균 ( canopy_vol ) , sd = sd ( canopy_vol ))
# 출력
plant_var n 평균 SD
< chr > < 정수 > < dbl > < dbl >
1 10 0.784 0.121
(2) B 10 0.608 0.0968
(3) C 10 0.272 0.143
4 D 10 0.474 0.0945

데이터세트 시각화,

라이브러리 ( gridExtra )
p1 < ggplot ( df, aes ( x = plant_var, y = 높이, 채우기 = plant_var )) + geom_boxplot ( outlier. shape = NA ) + geom_jitter ( 너비 = 0.2 ) + 테마 ( 범례. 위치 = “top” )
p2 < ggplot ( df, aes ( x = plant_var, y = canopy_vol, fill = plant_var )) + geom_boxplot ( outlier. shape = NA ) + geom_jitter ( 너비 = 0.2 ) + 테마 ( 범례. 위치 = “top” )
그리드. 정렬 ( p1, p2, ncol= 2 )

단방향 MANOVA 수행

dep_vars < cbind ( df$height, df$canopy_vol )
적합 < – 마 노바 ( dep_vars ~ plant_var, 데이터 = df )
요약 ( 적합 )
# 출력
Df Pillai approx F num Df den Df Pr (> F )
plant_var 3 1.0365 12.909 6 72 7.575e-10 ***
잔차 36
# 효과 크기 가져오기
라이브러리 ( 효과크기 )
효과 크기:: eta_squared ( 맞춤 )
# 출력
매개변수 | Eta2 ( 부분 ) | 95 % CI
——————————————
plant_var | 0.52 | [ 0.36 , 1.00 ]

Pillai의 Trace 테스트 통계는 통계적으로 유의미하며[Pillai’s Trace = 1.03, F(6, 72) = 12.90, p < 0.001], 식물 품종이 결합된 식물 높이 및 캐노피 부피 모두와 통계적으로 유의한 연관성이 있음을 나타냅니다.

효과 크기의 측정값(부분적 Eta 제곱, η p 2 )은 0.52이며 식물 높이와 수관 부피 모두에 식물 품종의 큰 영향이 있음을 시사합니다.

사후 테스트

MANOVA 결과는  식물 품종 간에 통계적으로 유의한( p < 0.001) 차이 가 있음을 시사 하지만 어떤 그룹이 서로 다른지는 알 수 없습니다. 어떤 그룹이 크게 다른지 알기 위해서는 사후 테스트를 수행해야 합니다.

그룹 간 차이를 테스트하기 위해 각 종속 변수에 대해 일변량 ANOVA를 수행할 수 있지만 이는 적절하지 않으며 여러 변수에서 함께 얻을 수 있는 정보를 잃게 됩니다.

여기서는 선형 판별 분석(LDA)을 수행하여 각 그룹 간의 차이점을 확인합니다. LDA는 두 종속 변수의 정보를 사용하여 그룹을 구별합니다.

도서관 ( 질량 )
post_hoc < lda ( df$plant_var ~ dep_vars, CV=F )
사후
# 출력
부르다:
lda ( df$plant_var ~ dep_vars, CV = F )
이전 확률 그룹 :
ABCD
0.25 0.25 0.25 0.25
그룹 의미:
dep_vars1 dep_vars2
A 18.90 0.784
16.54 0.608
C 3.05 0.272
D 9.35 0.474
계수 선형 판별 식 :
LD1 LD2
dep_vars1 -0.4388374 -0.2751091
dep_vars2 -1.3949158 9.3256280
비율 추적 :
LD1 LD2
0.9855 0.0145
# 구성
plot_lda < – 데이터. 프레임 ( df [ , “plant_var” ] , lda = 예측 ( post_hoc ) $x )
ggplot ( plot_lda ) + geom_point ( aes ( x = lda. LD1 , y = lda. LD2 , color = plant_var ) , 크기 = 4 )

LDA 산점도는 두 개의 종속 변수를 기반으로 여러 식물 품종을 구별합니다. C 및 D 식물 품종은 A 및 B에 비해 상당한 차이(잘 분리됨)가 있습니다. A 및 B 식물 품종은 서로 더 유사합니다. 전반적으로 LDA는 여러 식물 품종을 구별했습니다.

MANOVA 가정 테스트

다변량 정규성의 가정

이 테스트는 모든 통계 소프트웨어 패키지에서 사용할 수 없기 때문에 테스트하기 쉽지 않을 수 있습니다. 독립 변수와 종속 변수의 각 조합에 대해 일변량 정규성을 초기에 확인할 수 있습니다. 이 테스트가 통과하지 못하면(상당한  pv 값) 다변수 정규성이 위반될 수 있습니다.

참고: 다변량 중심 극한 정리에 따라 독립 변수와 종속 변수의 각 조합에 대해 표본 크기가 큰 경우(예: n > 20) 다변량 정규성을 가정할 수 있습니다.

라이브러리 ( rstatix )
df % > % group_by ( plant_var ) % > % shapiro_test ( 높이, canopy_vol )
plant_var 변수 통계 p
< chr > < chr > < dbl > < dbl >
1 A 캐노피_vol 0.968 0.876
2 높이 0.980 0.964
(3) B의 canopy_vol 0.882 0.137
4 B 높이 0.939 0.540
(5) C의 canopy_vol 0.917 0.333
6 C 높이 0.895 0.194
7 D 캐노피_vol 0.873 0.109
8 D 높이 0.902 0.231

애즈  태양 광 ALUE 비 상당한이다 ( P  독립 변수와 종속 변수의 각 조합> 0.05), 우리는 귀무 가설을 기각하고 데이터 정상 변량을 따라 체결 못한다.

표본 크기가 크면(예: n > 50) QQ-plot 및 histogram과 같은 시각적 접근 방식이 정규성 가정을 평가하는 데 더 좋습니다.

이제 Mardia의 왜도 및 첨도 검정을 사용하여 다변량 정규성을 확인하겠습니다.

라이브러리 ( mvnormalTest )
mardia ( df [ , c ( “높이” , “canopy_vol” )]) $mv. 시험
# 출력
검정 통계량 p-값 결과
1 왜도 2.8598 0.5815
2 첨도 -0.9326 0.351
3 MV 정규성 < NA > < NA >

 Mardia의 왜도 및 첨도 검정 에서  p-값 이 유의하지 않기 때문에( p > 0.05) 귀무 가설을 기각하지 못하고 데이터가 다변량 정규성을 따른다는 결론을 내립니다.

여기에서 왜도와 첨도 p  값은 모두  다변량 정규성을 결론짓기 위해 > 0.05여야 합니다.

분산-공분산 행렬의 동질성

우리는 분산-공분산 행렬의 동질성을 평가하기 위해 Box의 M 테스트를 사용할 것입니다. 귀무 가설 : 독립 변수의 각 그룹에 의해 형성된 각 조합에 대해 분산-공분산 행렬이 동일합니다.

도서관 ( heplots )
boxM ( Y = df [ , c ( “높이” , “canopy_vol” )] , 그룹 = df$plant_var )
# 출력
공분산 행렬의 동질성 대한 상자 M-검정
데이터: df [ , c ( “높이” , “canopy_vol” )]
Chi- 제곱은 ( 약. ) = 21.048 , DF = 9 , P 값 = 0.01244

 Box의 M 검정 에서  pv 값이 유의하지 않기 때문에( p > 0.001) 귀무 가설을 기각하지 못하고 분산-공분산 행렬이 독립 변수의 각 그룹에 의해 형성된 종속 변수의 각 조합에 대해 동일하다는 결론을 내립니다.

이 가정이 실패하면 Bartlett 또는 Levene의 검정을 사용하여 분산 가정의 동질성을 확인하여 등분산에서 실패한 변수를 식별하는 것이 좋습니다.

다변수 이상값

MANOVA는 이상값에 매우 민감하며 유형 I 또는 II 오류를 생성할 수 있습니다. Mahalanobis 거리 테스트를 사용하여 다변수 이상값을 감지할 수 있습니다. Mahalanobis 거리가 클수록 이상값일 가능성이 커집니다.

라이브러리 ( rstatix )
# 거리를 구하다
mahalanobis_distance ( data = df [ , c ( “height” , “canopy_vol” )]) $입니다. 국외자
# 출력
[ 1 ] 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓
[ 18 ] 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓 거짓
[ 35 ] 거짓 거짓 거짓 거짓 거짓 거짓 거짓

결과  에서 데이터 세트에 다변수 이상값(모두 is.outlier = FALSE 또는 p > 0.001) 이 없습니다  . is.outlier = TRUE이면 데이터 세트에 다변수 이상치가 있음을 의미합니다.

선형성 가정

선형성 가정은 각 그룹에 대한 종속 변수에 대한 쌍별 산점도를 시각화하여 확인할 수 있습니다. 선형성 가정을 충족하려면 데이터 포인트가 직선에 있어야 합니다. 선형성 가정을 위반하면 통계적 검정력이 감소합니다  .

라이브러리 ( gridExtra )
p1 < – df % > % group_by ( plant_var ) % > % filter ( plant_var == “A” ) % > % ggplot ( aes ( x = 높이, y = canopy_vol )) + geom_point () + ggtitle ( “다양성: A ” )
p2 < – df % > % group_by ( plant_var ) % > % filter ( plant_var == “B” ) % > % ggplot ( aes ( x = 높이, y = canopy_vol )) + geom_point () + ggtitle ( “종류: B ” )
p3 < – df % > % group_by ( plant_var ) % > % filter ( plant_var == “C” ) % > % ggplot ( aes ( x = 높이, y = canopy_vol )) + geom_point () + ggtitle ( “다양성: C ” )
p4 < – df % > % group_by ( plant_var ) % > % filter ( plant_var == “D” ) % > % ggplot ( aes ( x = 높이, y = canopy_vol )) + geom_point () + ggtitle ( “다양성: D ” )
그리드. 정렬 ( p1, p2, p3, p4, ncol= 2 )

산점도는 종속 변수가 독립 변수의 각 그룹에 대해 선형 관계를 갖는다는 것을 나타냅니다.

다중공선성 가정

다중  공선성은 종속변수 간의 상관관계를 통해 확인할 수 있다. 종속 변수가 3개 이상인 경우 상관 행렬 또는 분산 팽창 요인을 사용하여 다중 공선성을 평가할 수 있습니다.

종속변수 간의 상관관계는 0.9보다 크거나 너무 낮아서는 안 됩니다. 상관 관계가 너무 낮으면 각 종속 변수에 대해 별도의 일변량 ANOVA를 수행할 수 있습니다.

코. 테스트 ( x = df$height, y = df$canopy_vol, 방법 = pearson ) $estimate
# 출력
코르
0.8652548

종속변수 간의 상관계수가 0.9 미만이므로 다중공선성이 없다.

Leave a Comment