Differences

This shows you the differences between two versions of the page.

--- bioinfo:itol에서_annotated_tree_만들기 [2023/08/11 17:50] – [트리 자료의 준비] hyjeong
+++ bioinfo:itol에서_annotated_tree_만들기 [2023/08/17 09:06] (current) – [BSR matrix의 heatmap 표현] hyjeong
@@ Line 68: / Line 68: @@
   ...
-labels.txt 파일을 만드는 방법은 비교적 쉽다. 실습용 zip 파일의 압축을 해제하였을 때 나오는 accessions.txt 파일을 ncbi-genome-download('-n' or '%%--$$dry-run' 옵션)에 적용한 뒤, 출력물을 파일로 리다이렉트한 다음 적절히 편집하여 DATA 라인 아래에 채우면 된다. 구분자인 콤마 뒤에 공백을 넣지 않도록 주의한다.
+labels.txt 파일을 만드는 방법은 비교적 쉽다. 실습용 zip 파일의 압축을 해제하였을 때 나오는 accessions.txt 파일을 ncbi-genome-download('-n' or '%%--%%$dry-run' 옵션)에 적용한 뒤, 출력물을 파일로 리다이렉트한 다음 적절히 편집하여 DATA 라인 아래에 채우면 된다. 구분자인 콤마 뒤에 공백을 넣지 않도록 주의한다.
   $ ncbi-genome-download -A accessions.txt -n bacteria
-실제로 이 명령을 실행해 보면 accessions.txt에 수록된 51개 genome 중 50개에 대한 정보만이 반환될 것이다(2022-08-11 기준). 이는 [[https://www.ncbi.nlm.nih.gov/datasets/genome/GCF_001663585.1/|GCF_001663585.1(//Paenibacillus polymyxa ND24//)]]이 'unverified source organism'이라는 이유로 RefSeq에서 suppress되었기 때문이다.
+실제로 이 명령을 실행해 보면 accessions.txt에 수록된 51개 genome 중 50개에 대한 정보만이 반환될 것이다(2022-08-11 기준). 이는 [[https://www.ncbi.nlm.nih.gov/datasets/genome/GCF_001663585.1/|GCF_001663585.1]](//Paenibacillus polymyxa// ND24)이 'unverified source organism'이라는 이유로 RefSeq에서 suppress되었기 때문이다. 이를 적절히 수정한 labels.txt 파일({{:bioinfo:labels.zip|zipped}})을 업로드해 두었다.
 ===== Heatmap 자료의 표현 =====
+Heamap 데이터셋에서는 각 노드에 연결된 복수의 수치 자료(필드, row 형태로 나열)가 색상으로 표현된다. ANI 혹은 BSR 등이 heatmap으로 나타내기에 매우 적합한 자료 유형이다.
+수치는 annotation file([[https://itol.embl.de/help/dataset_heatmap_template.txt|template]])에 정의된 color gradient(COLOR_MIN, COLOR_MAX)에 매핑된 결과에 의하여 색으로 나타난다. 각 필드에 해당하는 라벨은 annotation file 내에서 FIELD_LABEL로 정의되는데, FIELD_TREE에서 그 순서에 대한 정보를 별도로 제공할 수 있다. FIELD_TREE 라인이 없거나 주석화되어 있으면 FIELD_LABEL과 실제 데이터에서 각 필드의 순서를 주의 깊게 결정해야 한다.
+앞 과정에서 만들었던 ANIm_modified_labels_sorted.csv를 heatmap으로 표현하는 방법을 알아보자. 컬럼/로 라벨은 전부 assembly accession(예: GCF_001573585.1)으로 바꾸어 놓은 상태이므로, heatmap에서 윗부분에 표현될 FIELD_LABEL은 실제 균주명 정보를 활용하여 “Clostridium botulinum CDC42961” 등의 형태로 표현되게 만들어야 한다. 그러기 위해서는 전 단계에서 만든 labels.txt 파일을 이용하여 컬럼 라벨을 되돌려 놓는 것이 핵심이다. 역시 R 환경에서 named vector를 이용한 데이터 조작 작업을 해야 한다. Labels.txt 파일은 실제 DATA 라인만을 추출하여 accession-label.txt 파일에 수록해 놓되, 구분자는 쉼표라고 가정하자.
+  # 이전 작업 내용을 지운다.
+  > rm(list=ls())
+  > ani = read.table("ANIm_modified_labels_sorted.csv",sep=",",header=T,row.names=1)
+  > data = read.table("accession-label.txt",sep=",",stringsAsFactors=F)
+  > key = data[,1]
+  > names(key) = data[,2]
+  # 만약 일부의 컬럼만 취하고 싶다면 해당되는 accession number를 수록한 벡터(key2)를 만들고
+  # 컬럼 라벨을 바꾸기 전에 다음을 실행한다.
+  > key2 = c("CBAW1310","GCF_000063585.1","GCF_000017025.1","GCF_001276985.1","GCF_001020205.1","GCF_001889325.1","GCF_000020285.1","GCF_000020165.1","GCF_000204565.1","GCF_000768655.1","GCF_000171095.1")
+  > ani.subset = ani[,which(colnames(ani) %in% key2)]
+  # key2 벡터에 들어갈 값이 너무 많아서 별도의 텍스트 파일(data.txt)에 line-by-line
+  # 형태로 들어 있다면, 이를 읽어들여서 처리하면 된다.
+  # data.txt 파일이 multi-column file이라면 $V2, $V3..을 붙여넣어 지정한다.
+  > key2 = read.table(“data.txt”,sep=”\t”, stringsAsFactors)$V1
+  # 꼭 which(colnames() %in% ) 구문을 쓰지 않아도 된다. 다음과 같이 하면
+  # key2 벡터에 출현하는 순서대로 컬럼 이름이 정렬되는 효과를 갖는다.
+  > ani.subset = ani[,key2]
+  # key2의 순서에 맞게 컬럼을 정렬한다.
+  > ani.subset = ani.subset[,key2]
+  # column name substitution using named vector("key") and for loop
+  > for (i in names(ani)) {
+    temp = names(key)[key==i]
+    names(ani)[names(ani)==i] = temp
+    }
+  > for (i in names(ani.subset)) {
+    temp = names(key)[key==i]
+    names(ani.subset)[names(ani.subset)==i] = temp
+    }
+  > write.table(ani,"ANIm_with_strain_names.csv",sep=",",quote=F)
+  > write.table(ani.subset,"ANIm_subset.csv",sep=",",quote=F)
+[[https://itol.embl.de/help/dataset_heatmap_template.txthttps://itol.embl.de/help/dataset_heatmap_template.txtheatmap.txt 템플레이트 파일]]에서 SEPARATOR를 COMMA로 선언한 뒤 ANIm_with_strain_names.csv의 내용을 적절히 재배치하고 수정하면 된다. 첫줄은 “FIELD_LABELS,” 바로 뒤에 한 줄로 넣도록 하고, 나머지 라인 전체는 “DATA” 라인 바로 밑에 넣으면 된다. 아울러 FIELD_TREE 라인은 주석 처리한다.
+이상의 방법에서는 ANI 매트릭스를 미리 정렬해 두고 iTOL에서도 그대로 표현되게 하느라 상당히 공을 들였다. 행 단위의 정렬은 iTOL에 최초로 업로드한 트리 자료를 따르게 되므로, 여기에서 설명한 정렬이란 열(컬럼) 단위, 즉 수치화된 필드의 정렬을 의미한다. 그런데 열 단위의 정렬 정보를 별도의 트리 파일로 보유하고 있다면 이를 heatmap.txt 데이터 파일에 넣는 것으로 작업이 단순해진다(FIELD_TREE 항목). 단, 콤마가 섞인 newick format의 자료가 데이터 파일에 포함되므로 field separator는 콤마가 아닌 것으로 바꾸어야 하고, leaf label에는 괄호 등이 들어가지 않게 유의해야 한다. 또한 newick file의 leaf label은 field label과 같아야 한다. 다시 말해서 GCA_ 또는 GCF_(assembly accession)로 시작하는 파일명을 지닌 유전체 염기서열을 이용하여 ANI 매트릭스를 만들었고 이것으로부터 hierarchical clustering을 하여 newick file을 만들게 되면 당연히 assembly accession이 leaf label이 될 것이다. 그러나 iTOL annotated tree에서 heatmap 데이터의 필드를 assembly accession으로 표현하고 싶지는 않을 것이다. 따라서 accession-label.txt 파일을 참조하여 assembly accession을 균주명 등의 정보로 치환하는 것이 바람직하다. 그러나 괄호나 콜론 등 균주명에서는 허용될지라도 newick 파일에서는 특수한 용도로 쓰이는 문자는 다른 것으로 치환해야 하니 세심한 작업이 필요하다.
 ===== BSR matrix의 heatmap 표현 =====
+LS-BSR에서 얻어진 마커 유전자에 대한 blast score ratio 매트릭스도 마찬가지 방법으로 정리하여 iTOL에서 heatmap 데이터셋으로 표현할 수 있다. LS-BSR의 설치 및 일반적인 사용 방법은 본 매뉴얼 내의 [[bioinfo:blast_score_ratio_bsr_을_이용한_마커_유전자의_탐색|다른 곳]]에서 설명하였다.
 앞에서 설명한 ANI 매트릭스는 column/raw의 이름이 iTOL 트리 자료의 leaf name 그대로이고 대각선에 대하여 대칭 구조이므로 원하는 자료만을 선별하거나 이름을 바꾸는 것 정도의 조작만 가하면 된다. 그러나 BSR 매트릭스는 각 유전체에 대한 BSR 값이 컬럼 단위로 위치하므로 이를 트랜스포즈해야 iTOL의 데이터셋 파일로 쓸 수 있게 된다. R에서 이를 실제로 어떻게 처리하는지 다음의 코드를 통해 알아보도록 한다.
+  > bsr = read.table("RUN_toxins_bsr_matrix.txt",sep="\t",row.names=1,header=T,check=F)
+  > bsr.t = t(bsr)
+  # 자료 구조 확인용
+  > dim(bsr.t)
+  > View(bsr.t)
+  > row.names(bsr.t) = gsub("^.*_GCF_","GCF_",row.names(bsr.t))
+  > write.table(bsr.t,file="bsr_matrix_modified.csv",sep=",",quote=F)
+컬럼으로 전환된 마커 유전자의 순서는 markers.fasta(or .pep) 파일에 있던 것과 달라질 수 있다. 아마도 LS-BSR 프로그램 내부에서 FASTA 파일의 서열 ID를 추출하여 문자열 기준으로 다시 정렬을 하는 것으로 보인다. 원래 상태로의 재정렬이 필요하다면 R 환경에서 하는 것이 좋을 것이다. 바로 위의 “Heatmap 자료의 표현” 항목에서 컬럼의 일부만을 추출하는 방법을 참고하기 바란다.