본 글은 다음 논문에 대하여 정리한 글입니다.
Kim, M., Oh, H. S., Park, S. C., & Chun, J. (2014). Towards a taxonomic coherence between average nucleotide identity and 16S rRNA gene sequence similarity for species demarcation of prokaryotes. International journal of systematic and evolutionary microbiology, 64(2), 346-351.
3줄 요약
1. 본 연구에서는 ANI와 16S rRNA gene 서열 유사도 사이의 관계를 비교하였음.
2. 0.01%의 16S rRNA gene 서열 유사도를 간격으로 95%와 96% ANI 경계값에 대해 pairwise comparision을 진행한 결과 두 ANI 경계값 모두 98.65%를 최적의 16S rRNA gene 서열 유사도 경계값으로 갖는다는 것을 확인함.
3. 몇 가지 예외가 있는데 유전체 내의 여러 16S rRNA gene 사이의 유전적 다양성이 큰 경우 등이 있으며 이는 추가적인 연구가 필요할 것으로 보임.
지난 50년간 DDH (DNA-DNA 혼성화, DNA-DNA hybridization)는 박테리아 종을 구분하는데 명확하고 객관적인 경계값을 주며 'golden standard'로써 널리 사용되어 왔다. 그러나 DDH 실험 자체가 노동 집약적이고 오류가 발생하기 쉽기 때문에 유전자를 바탕으로 한 대체 표준에 대한 요구가 꾸준히 있어왔다.
Genome sequence 유사도는 DDH를 대신할만한 분류학적 파라미터가 될 수 있다는 것이 명백했기 때문에 digital한 Genome sequence 비교를 위한 방법이 많이 개발되었다. 이는 ANI (average nucleotide identity) [1] , GBDP (genome BLAST distance phylogeny) [2], MUMi (maximal unique matches index) [3] 등을 포함한다. 이들 중 ANI는 종 분류를 위한 차세대 golden standard로 가장 널리 사용되어 왔다. ANI는 두 유전체의 상응하는 영역에 대한 유사도를 나타내어주는 값으로 95~96%의 ANI value을 대략 70%의 DDH value로 동일시하여 종 구분에 사용하고 있다.
초기의 DDH와 16S rRNA 서열 유사도 사이의 관계는 97%의 16S rRNA 서열 유사도가 70%의 DDH value에 대응된다고 [4] 연구되어 97% 이상의 유사도를 갖는 서열에 대해서만 DDH가 요구되었지만 요즘에는 380개의 균주 비교를 통하여 얻어진 98.7%~99.0%의 새로운 경계값이 사용되기도 한다. [5] 그러나 16S rRNA 서열 유사도가 99% 이상임에도 DDH를 진행했을 때 서로 구분되는 경우도 많이 나오고 있기 때문에 16S rRNA 서열만을 통해서 두 종을 구분하는 것은 어렵다. 더 최근에는 98.2%~99.0%의 값이 571개의 균주 비교를 통하여 제안되기도 했다. [6] 이들 연구는 DDH와 16S rRNA gene sequence 유사도를 비교하는데 초점을 맞추었다. 이 연구에서는 ANI와 16S rRNA gene 서열 유사도를 비교하는데 초점을 맞추겠다.
과정
11995개의 원핵생물 genome sequence가 GenBank database에서 검색되었다. 종 분류는 표준 균주와 16S rRNA gene 서열 유사도와 ANI를 계산하여 확인하였다. Low-quality genome, >2000 contigs의 genome, draft genomes stemming from single-cell genomics using multiple displacement amplification은 낮은 genome coverage와 poor quality로 인하여 사용하지 않았다. 또한 full 16S rRNA gene 서열을 포함하지 않는 genome sequence도 사용하지 않았다. 결과적으로 6787개의 genome이 후속 연구로 진행되었다.
ANI value는 같은 family (과)에 속하는 모든 genome에 대하여 상호 계산되었다. 다른 family 간의 ANI 비교는 <60%의 낮은 ANI value를 가지기 때문에 고려되지 않았다. 모든 ANI 계산은 BLASTN과 Goris et al(2007) [7]의 알고리즘을 이용한 내부 소프트웨어를 이용하여 계산되었다. 계산은 family 내에 최소한 3개의 genome과 2개의 종이 있는 경우에만 진행하였다. ANI value의 분포에 대한 통계는 R statistical software (http://www.r-project.org) 를 이용했다. ANI value를 normalize (정규화) 하기 위하여 같은 종에서 나온 모든 균주 사이의 a mean of pairwise ANI value가 계산되었다. rRNAselector [8] 프로그램이 각 genome에서 16S rRNA gene 서열을 얻기 위하여 사용되었다. 16S rRNA gene 서열의 완전성은 [9]를 통하여 검증하였다. 16S rRNA gene 서열 사이의 Pairwise similarities는 EzTaxon server에 있는 robust global sequence alignment algorithms을 이용하였다. [10] ANI와 16S rRNA gene 서열 유사도 사이의 관계는 logarithmic regression (로그 회귀)을 이용하여 확인하였다.
주어진 ANI 경계값에 대응하는 최적의 16S rRNA gene 서열 유사도 경계값을 결정하기 위하여 precision-recall과 F score를 16S rRNA gene 서열 유사도에 대하여 0.01% 간격으로 계산하였다. F measure는 정보 검색 과정의 분류 성능을 측정하기 위하여 고안되었고 [11] 이진, 또는 다중 classifiers의 성능을 평가하는데 사용되고 있다. 최근에는 원핵생물 종 구분을 위하여 최적 서열 유사도 경계값을 결정하고 [12] 분류, 계통학적으로 일관적인 위치를 확인하는데 [13] 사용되고 있다. 본 실험에서는 pairwise ANI value가 95% 또는 96% 미만이라면 다른 종, 아니라면 같은 종으로 추정하였다. [7][14] 모든 pairwise 16S rRNA gene 서열 유사도 값은 주어진 경계값 X ANI, Y 16S rRNA 서열 유사도를 기준으로 4개의 카테고리로 분류되었다.
TP (True, Positives) 카테고리 : ANI $\geq$ X, 16S rRNA 서열 유사도 $\geq$ Y
FN (False, Negatives) 카테고리 : ANI $\geq$ X, 16S rRNA 서열 유사도 < Y
FP (False, Positives) 카테고리 : ANI < X, 16S rRNA 서열 유사도 $\geq$ Y
TN (True, Negatives) 카테고리 : ANI < X, 16S rRNA 서열 유사도 < Y
각 ANI와 16S rRNA 서열 유사도 경계값에 대하여 평가 라운드가 끝날 때마다 이중 교차검증 (twofold cross-validation)을 하기 위하여 전체 dataset은 랜덤하게 두개의 subset으로 갈라졌고 각 subset에 대하여 precision and recall value가 위의 4개의 카테고리를 이용하여 추정되었다. 최적의 경계값은 sensitivity (recall)을 최대화하는 동시에 false discovery rate (1-specificity)를 최소화하는 값으로 결정되었다. 그 후에 precision and recall values가 precision and recall의 harmonic mean임과 동시에 테스트의 정확성을 의미하는 F score를 계산하기 위하여 사용되었다. 16S rRNA 서열 유사도를 대상으로 가장 높은 F score를 가지는 값이 각 subset의 컷오프로 선정되었고 각 컷오프의 성능을 holdout subset을 대상으로 평가하였다.
결과
최종 dataset은 22개 phyla (문), 41 classes (강), 93 order (목), 202 families (과), 655 genera (속), 1738 species (종)의 총 1,044,179개의 16S rRNA gene 서열 유사도와 ANI의 pairwise value로 구성되었다. 같은 family (과)에 속하는 균주들을 대상으로 계산된 ANI value는 81.0%~96.0% 범위의 낮은 빈도의 영역과 >96% 범위의 높은 빈도 영역으로 나뉘는 등 고르지 못한 분포를 보였다. 대부분의 종내 ANI value는 96% 이상으로 [14]에서 권고되었던 종 구분 값에 해당되었다. DDH에서도 비슷한 패턴을 볼 수 있었는데, 70% 근처에서 움푹 꺼진 모양을 몰 수 있었다. [6] 반면 속 사이의 ANI value 비교에서는 62%~100% 사이의 넓은 분포를 보였다.
공개된 genome sequence 수가 분류학적 그룹에 따라 불균등하기 때문에 이런 편향을 완화시키기 위하여 종 내에서 single averaged ANI value을 얻는 것으로 sampling을 normalize (정규화)했다. 이는 결과가 정규화하기 이전과 전반적으로 비슷한 패턴을 보였으나 정규화했을 때에 96%~100% ANI value에서도 상당히 낮은 빈도를 보였다는 점에서 차이를 확인할 수 있다. 이러한 차이가 생기는 이유는 의학적으로 중요한 종들에 대해 시퀀싱이 많이 되어있기 때문인데, 특히 Escherichia coli genome만 512개가 본 연구에서 사용되었다.
가끔 서로 다른 종임에도 ANI value가 96%를 넘는 경우가 있었다. Escherichia coli-Shigella 종들, Burkholderia mallei-Burkholderia pseudomallei (주. 작성자가 직접 ggdc.dsmz.de 서버로 표준 균주에 대한 dDDH test를 진행한 결과 recommended인 identities/HSP length에서 Escherichia coli-Shigella 사이에서는 >73%, Burkholderia 사이에서는 92.70%를 보였음, golden standard인 DDH 기준으로는 동일 종으로 분류해야 함), Bordetella bronchiseptica-Bordetella parapertussis-Bordetella pertussis 가 잘 알려진 예시다. [14] 또한 본 연구에서는 추가로 Bacillus anthracis–Bacillus thuringiensis–Bacillus cereus, Yersinia pseudotuberculosis– Yersinia pestis, Brucella species, Lactobacillus casei–Lactobacillus paracasei, Mycobacterium tuberculosis–Mycobacterium bovis, Leptospira kirschneri–Leptospira interrogans, Treponema paraluiscuniculi–Treponema pallidum 의 경우를 발견하였다. 이들은 추가적인 분류학적 조사를 통하여 독립된 분류학적 상태를 유지하는 것이 가치있는지 확인해보아야 할 것이다. (주. Yersinia를 제외하고는 taxonomic status에 변동이 없다)
ANI value와 16S rRNA gene 서열 유사도 사이의 관계는 선형이 아닌 것으로 확인되었다. 이는 70개의 작은 dataset으로 진행한 [15] 연구와 동일한 경향을 보인다. 그러나 선형 관계는 로그 변환 후에 $r^2=0.805$과 P<0.001으로 찾을 수 있다. 비슷한 경향으로 기존의 연구에 의하면 DDH와 16S rRNA gene 서열 유사도는 비선형 관계를 보이지만 로그 변환, 또는 로그-로그 변환 후에는 높은 상관관계가 얻어졌다.
95% 또는 96% ANI에 대한 최적의 16S rRNA gene 서열 유사도의 최적 경계값을 찾기 위해서 precision-recall과 F score가 0.01% 간격으로 모든 16S rRNA gene 서열 유사도에 대하여 계산되었다. 본 연구에서는 이중 교차검증 (twofold cross-validation)을 통하여 98.65% 16S rRNA gene 서열 유사도가 95%와 96%의 두 ANI value에 대해서 동시에 가장 높은 F score 값을 가짐을 확인하였다. 이 값은 지난 연구에서 [6] 제안된 범위 안에 존재한다.
본 연구에서는 이 경계값을 현재의 분류를 검증하기 위하여 dataset에 적용하였다. ANI와 16S rRNA gene 서열 유사도에 대한 genome pairs의 분포는 위 Fig. 3. 와 같다. 16S rRNA gene 서열 유사도 경계값인 <98.65%을 대상으로 ANI value를 <95%, 95%~96%, $\leq$96%의 세 기준으로 나누어보았다. <95% ANI (667,705 comparision)에서는 Fusobacterium nucleatum subsp. fusiforme ATCC $51190^T$ - F. nucleatum subsp. animalis ATCC $51191^T$ (16S rRNA gene 서열 유사도/ANI : 98.51/91.80 %), F. nucleatum subsp. animalis ATCC $51191^T$ - F. nucleatum subsp. nucleatum ATCC $25586^T$ (16S rRNA gene 서열 유사도/ANI : 98.56/90.64 %)의 오직 2가지의 예외만이 존재했다. (주. 작성자가 직접 Fusobacterium nucleatum의 subspecies를 대상으로 한 ggdc.dsmz.de의 dDDH test에서는 모두 서로 <50%의 값을 보였다. golden standard에 의하면 Fusobacterium nucleatum의 subspecies는 모두 아종이 아닌 종 수준으로 격상되어야 한다. )
95~96% ANI 범위에서는 0.53% (27/5104)가 <98.65%의 16S rRNA gene 서열 유사도를 보였다. $\leq$96% 범위에서는 0.10% (356/371370)가 <98.65%의 16S rRNA gene 서율 유사도를 보였다. 이런 예외들이 16S rRNA gene 서열 유사도 컷오프가 ANI를 바탕으로 결과에 어긋남에 따라 각 케이스에 대해 더 살펴볼 필요가 있다.
절반 이상의 예외 케이스들은 상대적으로 높은 수준의 16S rRNA 유전자간 이질성을 보였다. 383개의 >98.65% 16S rRNA gene 서열 유사도와 <95% ANI인 예외 케이스들 중 150 케이스가 E.coli-Shigella 그룹에서 발견되었다. E.coli는 특히 multiple 16S rRNA gene간의 유전자간 차이 (1.10%)가 있는 것으로 알려졌다. [16] 다른 13개의 경우도 16S rRNA gene에서 높은 수준의 유전체 내 변이를 보였다. Lactobacillus rhamnosus (0–7.67 %), Caldanaerobacter subterraneus (0.03–6.23 %), Desulfitobacterium hafniense (0.06–3.73 %), Bacteroides ovatus (0.07–3.30 %), Yersinia enterocolitica (0–2.67 %)
and Desulfitobacterium dehalogenans (0–2.14 %). 이런 16S rRNA gene에서의 높은 유전체 내 변이가 최적 16S rRNA gene 서열 유사도 경계값에 대한 차이를 만들었을 것이다. 또한 16S rRNA gene 서열에서 낮은 수준의 유전체 내 변이를 보여주는 경우도 있다. Aggregatibacter actinomycetemcomitans (44개 케이스)는 낮은 16S rRNA 서열 유사도 (98.00~98.20%)를 갖지만 높은 수준의 ANI (97.43%~99.09%)를 공유한다. 왜 해당 균주에서 높은 수준의 16S rRNA 수준의 다양성이 확인되는지는 알 수 없다. 추가적인 생태학적, 계통발생학적 조사가 이런 차이의 실마리를 줄 것이다.
36% (138 케이스)는 Neisseria meningitidis alpha704 (GenBank accession no. CAJS00000000)와 다른 N.meningitidis 균주들이 97.00–98.64%의 16S rRNA gene 서열 유사도를 갖지만 96.68–98.12%의 ANI를 갖는 경우이다. 이 경우 N.meningitidis alpha704의 퀄리티가 의문스러운데 시퀀싱의 draft stage에 있으면서 42개의 contig와 단 하나의 완전한 16S rRNA gene 서열을 포함하기 때문이다. 보통 N.meningitidis 균주는 네개의 rrn operons을 갖는다. 이에 대한 원인도 탐구되어야 할 것이다.
이런 종 내의 예외들과 더불어 종간의 예외들도 발견된다. Halomicrobium katesii DSM $19301^T$ 와 Halomicrobium mukohataei DSM $12286^T$ (96.71% ANI와 95.59% 16S rRNA 서열 유사도), Thermoanaerobacter 종들 (95.03–98.16% ANI와 97.02–98.54% 16S rRNA gene 서열 유사도), Caldicellulosiruptor kronotskyensis $2002^T$ 와
Caldicellulosiruptor bescii DSM $6725^T$ (95.59% ANI와 97.99% 16S rRNA gene 서열 유사도)가 있다. (주. 작성자가 직접 Halomicrobium의 두 종을 바탕으로 한 ggdc.dsmz.de dDDH test에서는 73.50%의 값으로 species가 아닌 subspecies 정도의 관계가 되어야 하는 것으로 보이는 반면 Caldicellulosiruptor에서 한 dDDH에서는 65.30%가 나와 별개의 종으로 인정할 수 있다는 결과가 나왔다. 이 경우 97.99%라는 상대적으로 높은 ANI 값을 가지는 동시에 70% 이하인 dDDH 값을 가지는 것을 보면 ANI value와 dDDH value 사이의 관계 또한 조사될 필요가 있을 것이다.) Thermoanaerobacter와 Caldicellulosiruptor에 속하는 종들은 높은 수준의 multiple 16S rRNA gene 사이의 차이가 있다고 알려졌다. [16][17] 호염성 고세균 속인 Halomicrobium 역시 종내 높은 수준의 16S rRNA gene 변이가 보고되었다. [18] 이러한 종들의 분류는 예외적으로 낮은 16S rRNA gene 서열 유사도를 감안할 때 더 많은 주의를 요한다.
전 지구상에 수백만의 원핵생물 종이 있을 것으로 추정되지만 오직 일부만이 발견되었다. 현재 validly published된 원핵생물 이름은 오직 11,000개에 불과하다. 새로운 종을 확인하는데 사용되는 16S rRNA gene 서열 유사도 경계값의 역할은 그것이 원핵생물 분류학에 도입된 이래로 매우 중요했다. 본 연구에서는 기존의 연구보다 더 많은 샘플을 비교하였으며 최근에 출판되는 DDH value들은 70%보다 낮은 경우에만 주로 출판되기 때문에 이로 인하여 발생하는 편향을 제거하였다는 점에서 큰 의미가 있다. ANI를 DDH에 대한 대체제로 사용하면서, 우리는 98.65%의 경계값이 새로운 종을 확인하는데 있어 속도를 크게 높여줄 것으로 예상한다. 다만 이 경계값은 다음 조건에서 조심스레 사용되어야 하는데, 최소한 양쪽 가닥이 모두 시퀀싱되어야 하며 16S rRNA gene 서열은 완전해야 한다.
[1] Konstantinidis, K. T., Ramette, A., & Tiedje, J. M. (2006). The bacterial species definition in the genomic era. Philosophical Transactions of the Royal Society B: Biological Sciences, 361(1475), 1929-1940.
[2] Henz, S. R., Huson, D. H., Auch, A. F., Nieselt-Struwe, K., & Schuster, S. C. (2005). Whole-genome prokaryotic phylogeny. Bioinformatics, 21(10), 2329-2335.
[3] Deloger, M., El Karoui, M., & Petit, M. A. (2009). A genomic distance based on MUM indicates discontinuity between most bacterial species and genera. Journal of bacteriology, 191(1), 91-99.
[4] Stackebrandt, E., & GOEBEL, B. M. (1994). Taxonomic note: a place for DNA-DNA reassociation and 16S rRNA sequence analysis in the present species definition in bacteriology. International journal of systematic and evolutionary microbiology, 44(4), 846-849.
[5] Stackebrandt, E. (2006). Taxonomic parameters revisited: tarnished gold standards. Microbiol. Today, 33, 152-155.
[6] Meier-Kolthoff, J. P., Auch, A. F., Klenk, H. P., & Göker, M. (2013). Genome sequence-based species delimitation with confidence intervals and improved distance functions. BMC bioinformatics, 14(1), 1-14.
[7] Goris, J., Konstantinidis, K. T., Klappenbach, J. A., Coenye, T., Vandamme, P., & Tiedje, J. M. (2007). DNA–DNA hybridization values and their relationship to whole-genome sequence similarities. International journal of systematic and evolutionary microbiology, 57(1), 81-91.
[8] Lee, J. H., Yi, H., & Chun, J. (2011). rRNASelector: a computer program for selecting ribosomal RNA encoding sequences from metagenomic and metatranscriptomic shotgun libraries. The Journal of Microbiology, 49(4), 689-691.
[9] Kim, O. S., Cho, Y. J., Lee, K., Yoon, S. H., Kim, M., Na, H., ... & Chun, J. (2012). Introducing EzTaxon-e: a prokaryotic 16S rRNA gene sequence database with phylotypes that represent uncultured species. International journal of systematic and evolutionary microbiology, 62(Pt_3), 716-721.
[10] Chun, J., Lee, J. H., Jung, Y., Kim, M., Kim, S., Kim, B. K., & Lim, Y. W. (2007). EzTaxon: a web-based tool for the identification of prokaryotes based on 16S ribosomal RNA gene sequences. International journal of systematic and evolutionary microbiology, 57(10), 2259-2261.
[11] van Rijsbergen, C. J. (1979). Information Retrieval, 2nd edButterworths.
[12] Mende, D. R., Sunagawa, S., Zeller, G., & Bork, P. (2013). Accurate and universal delineation of prokaryotic species. Nature methods, 10(9), 881-884.
[13] McDonald, D., Price, M. N., Goodrich, J., Nawrocki, E. P., DeSantis, T. Z., Probst, A., ... & Hugenholtz, P. (2012). An improved Greengenes taxonomy with explicit ranks for ecological and evolutionary analyses of bacteria and archaea. The ISME journal, 6(3), 610-618.
[14] Richter, M., & Rosselló-Móra, R. (2009). Shifting the genomic gold standard for the prokaryotic species definition. Proceedings of the National Academy of Sciences, 106(45), 19126-19131.
[15] Konstantinidis, K. T., & Tiedje, J. M. (2005). Genomic insights that advance the species definition for prokaryotes. Proceedings of the National Academy of Sciences, 102(7), 2567-2572.
[16] Pei, A. Y., Oberdorf, W. E., Nossa, C. W., Agarwal, A., Chokshi, P., Gerz, E. A., ... & Pei, Z. (2010). Diversity of 16S rRNA genes within individual prokaryotic genomes. Applied and environmental microbiology, 76(12), 3886-3897.
[17] Acinas, S. G., Marcelino, L. A., Klepac-Ceraj, V., & Polz, M. F. (2004). Divergence and redundancy of 16S rRNA sequences in genomes with multiple rrn operons. Journal of bacteriology, 186(9), 2629-2635.
[18] Cui, H. L., Zhou, P. J., Oren, A., & Liu, S. J. (2009). Intraspecific polymorphism of 16S rRNA genes in two halophilic archaeal genera, Haloarcula and Halomicrobium. Extremophiles, 13(1), 31-37.
