Differences

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

--- bioinfo:pan-genome_analysis [2023/06/27 13:25] – [기타 pan-genome 분석 도구] hyjeong
+++ bioinfo:pan-genome_analysis [2023/06/28 13:07] (current) – [Scoary: pan-genome-wide association study] hyjeong
@@ Line 5: / Line 5: @@
 ===== Roary 활용법 =====
-기존에 작성한 문서: [[bioinfo:roary|Roary를 이용한 미생물 유전체의 비교 분석]] - 겹치는 내용이 많음!
+기존에 작성한 문서: [[bioinfo:roary|Roary를 이용한 미생물 유전체의 비교 분석]] <- 겹치는 내용이 많음!
 기본적인 pan-genome analysis는 [[https://github.com/sanger-pathogens/Roary|Roary]] 파이프라인을 통하여 진행할 수 있다. Roary는 주석화된 유전체 서열 정보로부터 단백질 서열을 추출하여 [[https://sites.google.com/view/cd-hit|CD-HIT]]로 사전 클러스터링을 실시한 뒤, BLASTP로 all-against-all comparison을 거친 [[https://micans.org/mcl/|MCL algorithm(Markov Cluster Algorithm)]]으로 클러스터를 만드는 과정으로 이어진다. '-e [%%--%%mafft]' 옵션을 주면 core gene에 대한 alignment까지 실시하므로 작업 시간이 오래 소요된다. 이러한 기본 흐름 외에도 여러 부속 유틸리티를 포함한다.
@@ Line 56: / Line 56: @@
 Roary가 생성하는 파일은 [[https://sanger-pathogens.github.io/Roary/|Roary 웹사이트]]의 Output files에서 상세히 설명하였다. * 표시가 있는 3개의 파일은 '-e [%%--%%mafft]' 옵션과 함께 roary를 실행할 때에만 생성되는 것이다.
+{{ :bioinfo:roary_test_outdir_20230627.zip |}} <- 21개의 genome을 이용한 roary 결과 파일 묶음. 분석에 사용한 미생물 유전체의 RefSeq assembly accession은 내부의 00_accessions.txt 파일에서 확인하라.
   accessory.header.embl             *core_gene_alignment.aln
@@ Line 81: / Line 83: @@
   dnaA: a1a82563115b75beda45b943a4982aba_3      c89560a585264e8dfaed7ccb9edd3b24_19403  b468e412484598827c4d804e35e27cd8_7      ...(생략)
-clustered_proteins 파일의 총 라인 수는 summary_statistics.txt 파일의 'Total genes' 숫자와 동일하다. 클러스터 식별자(ID)는 **genes_presence_absence.csv(or .Rtab)**의 'Gene' 컬럼 및 pan_genome_reference.fa의 description field에서도 공통적으로 쓰인다. 각 클러스터를 구성하는 서열에서 가장 많이 나타나는 gene name을 가져와서 클러스터 식별자로 쓰는 것이 원칙이지만 여의치 않은 경우 cluster_XXX의 일반명을 사용한다. 클러스터는 단일 멤버 유전자도 전부 포함하기 때문에 **summary_statistics.txt** 파일의 total genes와 그 숫자가 일치한다. Clustered_proteins 파일은 나중에 설명할 query_pan_genome 명령에도 GFF(fixed) 파일과 함께 쓰이게 된다.
+clustered_proteins 파일의 총 라인 수는 summary_statistics.txt 파일의 'Total genes' 숫자와 동일하다. 클러스터 식별자(ID)는 **genes_presence_absence.csv(or .Rtab)**의 'Gene' 컬럼 및 pan_genome_reference.fa의 description field에서도 공통적으로 쓰인다. 각 클러스터를 구성하는 서열에서 가장 많이 나타나는 gene name을 가져와서 클러스터 식별자로 쓰는 것이 원칙이지만 여의치 않은 경우 cluster_XXX의 일반명을 사용한다. 클러스터는 단일 멤버 유전자도 전부 포함하기 때문에 **summary_statistics.txt** 파일의 total genes와 그 숫자가 일치한다. Clustered_proteins 파일은 나중에 설명할 query_pan_genome 명령에도 GFF(fixed) 파일과 함께 쓰이게 된다. genes_presence_absence.csv 파일에는 각 genome(strain) 컬럼에 개별 유전자의 ID가 수록되어 있으나, genes_presence_absence.Rtab은 유전자의 유무에 따라 0 또는 1의 값이으로 바뀌어 있고 유전자 기능 등 보조적인 정보를 담는 컬럼은 전부 삭제되어 있다. 즉 "Gene" 컬럼과 각 균주 컬럼으로만 구성된다.
 Prank를 이용하여 느리지만 정확한 core gene의 codon-aware alignment를 하려면 -e 옵션을, mafft를 사용한 빠른 alignment를 하려면 -e %%--%%mafft 옵션을 사용한다. Core gene alignment를 생략하면 전체 계산 과정이 매우 빠르게 끝나지만 pan_genome_reference.fa 파일이 생성되지 않음에 유의하라. **pan_genome_reference.fa**는 pan genome(core 및 accessory) 클러스터 각각에 대한 대표 서열(nucleotide)을 수록하고 있다. clustered_proteins 파일에서 가장 앞에 위치한 유전자가 각 클러스터의 대표 서열에 해당한다. 염기서열의 ID는 수정된 것을 따르며, 바로 뒤의 description field에는 클러스터 번호정보가 나온다.
@@ Line 180: / Line 182: @@
 oary를 이용하여 특정 균주의 유전체로부터 core와 accessory gene을 서로 구분해 냈다고 가정하자. 이어지는 질문은 accessory gene에 어떠한 생물학적 기능 혹은 경로가 더 많은지를 알아보는 것이다. 또는 분석에 이용한 모든 균주의 pan genome으로부터 특정 균주 그룹만이 공유하는 유전자 세트를 발굴하였다고 가정하자. 이러한 유전자 서브셋이 갖고 있는 enriched function or pathway에는 무엇이 있을까? Over-Representation Analysis(ORA)는 가장 오랫동안 쓰여온 pathway analysis 도구로서, GO나 KEGG Pathway 등 기능 정보가 부가된 전체 유전자 집합에 대하여 이로부터 유래한 유전자 서브셋이 어떠한 기능을 더 많이 갖고 있는지를 통계적 테스트와 함께 추출해 준다. 여기에서는 [[https://webgestalt.org|WebGestalt(WEB-based GEne SeT AnaLysis Toolkit)]]의 사용 방법을 알아보고자 한다. 총 다섯 개의 균주 유전체를 비교하는 것으로 가정한다. 목표는 LMT16-62가 보유한 전체 유전자에 대하여 특이적인 유전자에 더욱 많이 존재하는 기능을 찾는 것이다.
 ===== Panseq을 이용한 균주 특이적 염기서열 추출 =====
+<color #ed1c24>Panseq는 설치 및 사용 방법이 다소 난해하여 이를 사용하는 것을 적극 권장하지는 않는다. 그저 이런 프로그램도 있다는 정도로 이해하면 좋을 것이다.</color>
-===== 유전체의 1:1 비교를 통해서 특이적 서열을 영역 정보화 함께 찾는 방법 =====
+비교유전체 분석을 통해서 특정 균주의 유전체에만 존재하는 염기서열을 추출할 필요가 종종 있다. 앞서 소개한 roary나 orthoMCL에서는 유전자 단위에서 모든 유전체에 공통으로 존재하는 core gene과 나머지에 해당하는 accessory gene을 추출해 주므로 후자를 활용하면 각 균주에만 존재하는 strain-specific gene을 찾을 수 있다. 그러나 이러한 특이적 영역을 반드시 유전자 단위로만 추출할 필요는 없다. 유전체 단위로 특이적인 영역을 찾아낸다면 qPCR 프라이머 등을 설계할 때 더욱 유용하게 활용할 수 있다. [[https://github.com/chadlaing/Panseq|Panseq]] 프로그램을 novel region finder 모드로 실행하면 이러한 목적에 맞는 염기서열을 얻게 된다.
+Panseq의 실행에 필요한 파라미터는 settings.txt(파일 이름은 중요하지 않음)에 미리 설정해 두어야 한다. Panseq는 query 파일을 일정 길이(fragmentationSize 파라미터로 변경 가능)로 잘라서 MUMmer 3 패키지의 nucmer를 사용하여 reference 염기서열과 비교한다.
+유전체 서열 파일이 여러 contig로 구성된 경우에는 이것이 같은 균주에 속함을 명시적으로 나타내기 위하여 다음과 같이 서열 ID를 lcl|unique_identifier|sequence_id 포맷으로 바꾸어야 한다. 단일 서열로 구성된 경우는 수정할 필요는 없다.
+  >lcl|strain1|contig1
+  ATGCTTT…
+  >lcl|strain1|contig2
+  ATGCTTT…
+그러나 실제 상황에서는 서열 ID 뒤에 description 정보가 있으면 제대로 결과가 나오지 않는 것을 발견하였다. 따라서 서열 ID 라인을 일괄적으로 정리하는 작업이 필요하다. 염기서열 파일이 field1_field2_field3.fna일 때 field1을 새로 만들어질 서열 ID의 unique_identifier로 사용함과 동시에 description을 삭제하여 최종적으로 field1.fasta 파일에 저장하려면 다음과 같이 실행한다.
+  $ ls *fna | while read f
+  > do
+  > DATA=$(cut -d '_' -f 1,2 <<<$f)
+  > awk -v -s="lcl|${DATA}" '$1~/^>/{sub (/>/, ""); printf ">%s|%s\n", s, $1; next}{print}' $f > ${DATA}.fasta
+  > done
+  $ rm *fna
+Panseq를 실행하려면 텍스트 파일(settings.txt)에 파라미터 값들을 먼저 저장해 두어야 한다. 다음은 novel region finding 모드로 실행할 때 필요한 최소한도의 settings.txt 파일 사례이다. 이를 참조하되 /usr/local/apps/Panseq/settings.txt 파일을 이용하여 적절히 변경하도록 한다. 특이적 영역을 추출할 유전체 염기서열 파일은 queryDirectory에, 비교를 할 나머지 유전체 파일은 referenceDirectory에 넣어둔다. MUMmer와 BLAST+ 및 Muscle은 PATH 환경변수에 저장되어 있어야 한다.
+  # 모든 디렉토리는 '/'로 끝나야 한다.
+  queryDirectory	<query sequence가 있는 디렉토리>
+  referenceDirectory	<reference sequence가 있는 디렉토리>
+  baseDirectory	<출력 파일을 저장할 디렉토리>
+  numberOfCores	8
+  minimumNovelRegionSize	500
+  novelRegionFinderMode	no_duplicates # or unique
+  fragmentationSize	500
+  percentIdentityCutoff	95
+  runMode	novel
+  overwrite	1
+Novel region finder는 “no_duplicates”와 “unique”의 두 가지 모드로 작동한다. no_duplicates 모드에서는 query로부터 reference에는 없는 영역을 찾아서 multi-fasta format으로 출력하는 반면, unique 모드에서는 각 query에는 고유하되 reference에는 없는 영역을 출력한다. 실행 방법은 다음과 같으며, 결과 파일은 refereneceDirectory/novelRegions.fasta이다.
+  $ /usr/local/apps/Panseq/lib/panseq settings.txt
+Core/Accessary analysis 또는 loci selector 기능에 대해서는 매뉴얼을 참조하라.
 ===== 기타 pan-genome 분석 도구 =====
@@ Line 205: / Line 246: @@
 link-to-server.py를 실행하면 data/TestSet 디렉토리의 사본이 ~/pan-genome-visualization/public/dataset에 만들어지므로 npm을 종료한 뒤 재시작해도 그대로 접속 가능하다.
+뒤에서 소개할 LS-BSR도 pan genome의 분석이 가능하다.
 ===== Scoary: pan-genome-wide association study =====
 Genome-wide association study(GWAS)는 주로 인간의 SNP와 형질의 연관성을 알아보기 위한 대규모 집단 유전학 및 유전체학의 응용 분야로 알려져 있다. 박테리아의 GWAS, 즉 주로 임상적으로 중요한 형질(병원성이나 항생제 내성 등)과 SNP의 연관성을 탐구하는 연구는 매우 최근에 들어서 시작되었다.
@@ Line 211: / Line 254: @@
 Scoary의 입력물로는 gene_presence_absence.csv 파일과 더불어 trait table 파일(traits.csv)이 필요하다. Trait의 유무는 0과 1로 표현하되, 분석할 특성이 여러 개라면 콤마로 구분한다. 첫 컬럼은 균주명이고, 첫 줄은 제목에 해당한다. 다음은 trait table의 구조와 실행 사례이다.
+===== 읽어야 할 논문 =====
+[[https://genomebiology.biomedcentral.com/articles/10.1186/s13059-021-02473-1|Pandora: nucleotide-resolution bacterial pan-genomics with reference graphs]] Genome Biology (2021)