다중 비교 검정 및 중앙 집중식 분산 분석 (ANOVA)

독립한 군이 3군 이상 있을 때, 어느 군과 어느 군의 평균치에 유의차가 있는지, 평가를 하고 싶어지는 경우가, 연구의 현장에서는 자주 있는 이야기라고 생각합니다.
그러나 일반적으로 3군 이상 있는 경우에 각각의 군끼리를 다중 비교를 이용하여 2군에서의 검정을 몇번이나 반복하는 다중 비교 검정은 좋지 않은 수단으로 되어 있습니다. 왜냐하면, 유의 수준이 가정의 설정치 이상으로 상승해 버리기 때문입니다. (다중성 문제) (왜 유의 수준이 상승하는지에 대해서는 내 실력으로는 설명하기 어려운 ... 울음) 그래서 일원 배치 분산 분석법(ANOVA)에 의해 군간끼리 차이가 있는가 평가하는 것이 일반적이라고합니다. 그러나 ANOVA에서는 어떤 요인과 어떤 요인이 서로 차이가 있는지를 고려할 수 없습니다.
다중성의 문제를 해결하고, 다중 비교 검정을 실시하기 위한 방법으로서, 사전에 유의 수준을 좁히고, 후에 전체의 유의 수준을 5%로 조정하는 수법이 있다. 대표적인 것이 Tukey‒Kramer법입니다. Tukey‒Kramer법에서는 각 그룹의 데이터수(n)가 일치하지 않아도 문제 없기 때문에 사용하기 편리한 특징이 있습니다.
　

그러면 실제로 검정을 해 나갑니다

필요한 라이브러리와 사용할 데이터 로드

#multcompのインストール
#install.packages("multcomp", repos="http://cran.ism.ac.jp/")
library(multcomp)

setwd("パスの指定")
data<-read.table("boxplot_SG_CL.txt",header = T)
head(data)

# 　sample_no sample_class heatmap_class        SG        CL
# 1         6          10a             3 0.5318429 0.8482768
# 2         7          10b             3 0.5245027 0.8566048
# 3         8          10c             3 0.5535782 0.8916458
# 4         9          10d             3 0.5012658 0.8107002
# 5        10          10e             1 0.5120243 0.6768474
# 6        11          10f             3 0.5340269 0.8464090

●sample_no : 샘플 번호
●sample_class : 도입한 유전자(복사 기호로 나타낸다)
heatmap_class : 특정 유전자 특성에 따라 group1 ~ 5로 분류 (제어 샘플로서 WT : wild type, Vc : Vector control을 채용)
SG : 화학 분석 결과 (S-lignin/G-lignin 비율)
CL : 화학 분석 결과 (Carbohydrate/total lignin 비율)



이번에 사용한 데이터는 유전자 변형된 애기장대를 이용한 데이터가 됩니다. 유전자 변형된 식물이 어떻게 화학 성분이 변화하는지를 측정한 결과가 되고 있습니다.

그룹별 SG비 비교

SG_data<-data.frame(group=factor(data$heatmap_class),SG=data$SG)
plot(SG~group,d=SG_data) #boxplotを作成

ANOVA

SG_res1<-aov(SG~group,d=SG_data)
summary(SG_res1)

#              Df  Sum Sq  Mean Sq F value  Pr(>F)   
# group         6 0.02518 0.004197    3.72 0.00188 **
# Residuals   133 0.15005 0.001128                   
# ---
# Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

ANOVA는 1% 유의한 결과를 보여준다. 그러나 box plot을 보는 한 큰 차이가 보이지 않는다.

다중 비교 검정(Tukey‒Kramer법)

SG_res2<-glht(SG_res1,linfct=mcp(group="Tukey"))
summary(SG_res2)

# 
#         Simultaneous Tests for General Linear Hypotheses
# 
# Multiple Comparisons of Means: Tukey Contrasts
# 
# 
# Fit: aov(formula = SG ~ group, data = SG_data)
# 
# Linear Hypotheses:
#   Estimate Std. Error t value Pr(>|t|)   
# 2 - 1 == 0    0.005649   0.012937   0.437  0.99943   
# 3 - 1 == 0    0.009094   0.011589   0.785  0.98547   
# 4 - 1 == 0    0.004076   0.012482   0.327  0.99989   
# 5 - 1 == 0   -0.024931   0.010430  -2.390  0.20517   
# VC - 1 == 0  -0.008298   0.011135  -0.745  0.98891   
# WT - 1 == 0   0.007968   0.013533   0.589  0.99690   
# 3 - 2 == 0    0.003444   0.011854   0.291  0.99995   
# 4 - 2 == 0   -0.001574   0.012728  -0.124  1.00000   
# 5 - 2 == 0   -0.030580   0.010723  -2.852  0.07043 . 
# VC - 2 == 0  -0.013947   0.011410  -1.222  0.88036   
# WT - 2 == 0   0.002319   0.013761   0.169  1.00000   
# 4 - 3 == 0   -0.005018   0.011355  -0.442  0.99939   
# 5 - 3 == 0   -0.034025   0.009052  -3.759  0.00426 **
# VC - 3 == 0  -0.017392   0.009855  -1.765  0.56564   
# WT - 3 == 0  -0.001126   0.012502  -0.090  1.00000   
# 5 - 4 == 0   -0.029007   0.010170  -2.852  0.07038 . 
# VC - 4 == 0  -0.012374   0.010891  -1.136  0.91272   
# WT - 4 == 0   0.003893   0.013333   0.292  0.99995   
# VC - 5 == 0   0.016633   0.008462   1.966  0.43311   
# WT - 5 == 0   0.032899   0.011436   2.877  0.06611 . 
# WT - VC == 0  0.016266   0.012081   1.346  0.82341   
# ---
#   Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# (Adjusted p values reported -- single-step method)

cld(SG_res2)

#  1    2    3    4    5   VC   WT 
# "ab" "ab"  "b" "ab"  "a" "ab" "ab"

결과의 견해로서는, 값이 낮은 순서로 a, b, c... 라고 문자를 붙이고 있어 다른 문자가 쓰여져 있는 시료간(예를 들어, group3 「b」와 group5 「a」)로 0.5 % 유의차가 있음을 나타냅니다.

그룹별 CL 비율 비교

CL_data<-data.frame(group=factor(data$heatmap_class),CL=data$CL)
plot(CL~group,d=CL_data)

CL_res1<-aov(CL~group,d=CL_data)
summary(CL_res1)

#             Df Sum Sq Mean Sq F value   Pr(>F)    
# group         6  1.804 0.30069   18.46 1.56e-15 ***
# Residuals   133  2.166 0.01628                     
# ---
# Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

0.1% 유의한 결과를 나타낸다. 방금 전의 S/G데이터보다 변동이 큰 경향에 있는 것을 읽을 수 있다.

# Simultaneous Tests for General Linear Hypotheses
# 
# Multiple Comparisons of Means: Tukey Contrasts
# 
# 
# Fit: aov(formula = CL ~ group, data = CL_data)
# 
# Linear Hypotheses:
#   Estimate Std. Error t value Pr(>|t|)    
# 2 - 1 == 0   -1.896e-05  4.915e-02   0.000  1.00000    
# 3 - 1 == 0    2.542e-01  4.403e-02   5.773  < 0.001 ***
# 4 - 1 == 0    4.047e-01  4.742e-02   8.535  < 0.001 ***
# 5 - 1 == 0    1.367e-01  3.963e-02   3.450  0.01247 *  
# VC - 1 == 0   1.527e-01  4.230e-02   3.610  0.00735 ** 
# WT - 1 == 0   1.386e-01  5.142e-02   2.695  0.10438    
# 3 - 2 == 0    2.542e-01  4.503e-02   5.645  < 0.001 ***
# 4 - 2 == 0    4.048e-01  4.836e-02   8.371  < 0.001 ***
# 5 - 2 == 0    1.367e-01  4.074e-02   3.356  0.01704 *  
# VC - 2 == 0   1.527e-01  4.335e-02   3.523  0.00986 ** 
# WT - 2 == 0   1.386e-01  5.228e-02   2.651  0.11545    
# 4 - 3 == 0    1.506e-01  4.314e-02   3.490  0.01088 *  
# 5 - 3 == 0   -1.175e-01  3.439e-02  -3.416  0.01388 *  
# VC - 3 == 0  -1.015e-01  3.744e-02  -2.710  0.10031    
# WT - 3 == 0  -1.156e-01  4.750e-02  -2.434  0.18782    
# 5 - 4 == 0   -2.680e-01  3.864e-02  -6.937  < 0.001 ***
# VC - 4 == 0  -2.520e-01  4.138e-02  -6.091  < 0.001 ***
# WT - 4 == 0  -2.662e-01  5.066e-02  -5.254  < 0.001 ***
# VC - 5 == 0   1.599e-02  3.215e-02   0.498  0.99879    
# WT - 5 == 0   1.868e-03  4.345e-02   0.043  1.00000    
# WT - VC == 0 -1.413e-02  4.590e-02  -0.308  0.99993    
# ---
#   Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# (Adjusted p values reported -- single-step method)

cld(CL_res2)

#    1     2     3     4     5    VC    WT 
#  "a"   "a"   "c"   "d"   "b"  "bc" "abc" 

요약

* 도면은 Rstudio가 아니라 Excel에서 만들고 있습니다.
*WT:wild type의 데이터는 생략



・S/G비에 대해서는, 어느 그룹에 있어서도 큰 차이는 아니지만, 그룹 3, 5간에 5% 유의차가 있다.
・탄화수소와 총 리그닌의 비에 대해서는, 약간의 차이가 보였습니다.
또한 컨트롤과 Group3, 4에서 상승이 보였고 Group1,2에서 감소가 보였다.

되돌아가다

ANOVA・다중 비교 검정용으로 설계한 실험이 아니기 때문에, 결과의 무리한 감이 있습니다만...Rstudio에 의해, 양자의 검정을 실제로 실시할 수 있었습니다.
특히, 다중 비교 검정의 결과를 a, b, c...로 나타내는 방법은 유전자 발현량을 비교한 논문에서 자주 보입니다.
마지막으로 Rstudio를 사용하지 않고 Excel로 작도한 이유는, 사용하기 편리한 그래프를 할 수 없었기 때문입니다. 의외 Excel이 더 깨끗한 그래프를 만들 수 있습니다.