FDR correction in LEfSe

이전에 differential abundance analysis 수행하는 것에 대한 논문리뷰를 한 적이 있습니다.

 

해당 논문에서 저자는 LEfSe 는 CoDa (compositional data analysis) methods 를 사용하지 않으며, FDR correction 을 지원하지 않기 때문에 사용을 지양하는 것이 좋아보인다는 언급을 마지막에 하고 있습니다.

 

그런데 LEfSe 는 K-W test -> wilcoxon test 로 두 번에 걸쳐 유의한 차이가 있는지 검정을 수행하는데 fdr correction 까지 수행할 필요가 있는지 궁금해서 관련 글을 찾아봤습니다.

 

그러다 아래의 링크에서 관련 내용을 찾을 수 있었습니다.

 

Multiple comparison correcting

 

---

 

The p-values from LEfSe do require some kind of multiple comparison control, if you’re not using any other biological substructure. But if you do include any of the multiple-class-consistency or directionality tests (i.e. subclass), those can often be more strict than standard MHT correction just by themselves.

 

결국 LEfSe 는 K-W test 로 후보 feature 를 찾고, 그 feature 에 대해 다시 Wilcoxon 을 수행하는 순서인데 Wilcoxon test 는 subgrouping parameter 를 설정했을때만 수행됩니다. 따라서 subgrouping 을 설정하지 않으면 fdr correction 이 필요해보이지만, subgroup 을 설정했다면 따로 p value 를 수정할 필요 없이 그대로 사용해도 문제가 없을 것 같다는 답변입니다.

 

경우에 따라 fdr correction 을 수행할수도 있기 때문에 R 을 사용하여 fdr correction 하는 코드를 작성해두었습니다.

 

---

library('data.table')

lefse <- fread("lefse_p.csv", header = TRUE)

p_values <- lefse$pvalues

# methods 중 선택하여 사용
fdr_p_values <- p.djust(p_values, method = "fdr")
hochberg_p_values <- p.adjust(p_values, method = "hochberg")
hommel_p_values <- p.adjust(p_values, method = "hommel")
bonferroni_p_values <- p.adjust(p_values, method = "bonferroni")
BH_p_values <- p.adjust(p_values, method = "BH")
BY_p_values <- p.adjust(p_values, method = "BY")
none_p_values <- p.adjust(p_values, method = "none")

# 기존 결과에 column 추가
lefse$fdr_p_values <- fdr_p_values
lefse$hochberg_p_values <- hochberg_p_values


fwrite(lefse, file = "lefse_q.csv")

 

( run_lefse.py 를 수행할 때 '-a 1' 파라미터를 작성해야 모든 feature 에 대한 p value 가 담긴 결과파일을 얻을 수 있습니다. 저는 .res 파일을 csv 로 변환, p value 열 이름을 pvalues 라고 수정한 뒤 그 파일을 가지고 fdr correction 수행하였습니다.)