Differences

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

--- to_be_renamed [2019/01/25 11:27] – hyjeong
+++ to_be_renamed [2022/06/18 13:12] (current) – [Roary] hyjeong
@@ Line 61: / Line 61: @@
 또한 RAST server에서 export한 GFF3 파일도 Roary에서 그대로 쓰일 수가 없다. 왜냐하면 염기서열이 뒷부분에 있지 않기 때문이다. 뿐만 아니라 gene 없이 cds feature만 있다는 것도 문제가 된다. 몇 가지를 테스트해 본 경험으로 가장 바람직한 것은, RAST에서 export한 GFF3 file에서 CDS feature만 골라낸 것 + '##FASTA' 라인 + contig sequence 파일(유전자 염기서열 파일이 아님!)을 합쳐서 새로 만든 GFF3 파일을 사용하는 것이 좋다. 'Name='도 'Product='로 바꾸는 것을 강력 권장한다. 왜냐하면 이것이 gene id로 쓰이게 되면 중간에 공백이 들어 있어서 나중에 매우 불편해지기 때문이다.
+Roary를 실행하면 십중팔구는 다음과 같은 메시지와 함께 GFF 파일을 수정하여 fixed_input_files 디렉토리에 복사하게 된다. 나중에 query_pan_genome 스크립트로 GFF 파일을 대상으로 하는 작업을 할 때에는 수정된 것을 써야 한다.
+/06/30 15:53:28 Input file contains duplicate gene IDs, attempting to fix by adding a unique suffix, new GFF in the fixed_input_files directory:
 === Output files ===
 {{:output.png?600|test}}
@@ Line 80: / Line 82: @@
 === Command line tool의 사용법 ===
-테스트를 해 보았는데 그 동작이 완벽한 것 같지는 않다.
+원본 GFF 파일이 아니라 roary가 수정한 것(fixed_input_files/ 디렉토리에 있는 것)을 사용해야 한다. 인수로 주어질 수 있는 것은 *gff 또는 단일 GFF 파일이다.
   $ query_pan_genome -h # for help
   $ query_pan_genome -a union *.gff # 결과물: pan_genome_results
-  $ query_pan_genome -a intersection *.gff # 결과물: pan_genome_results (텅 비었음. why?)
+  $ query_pan_genome -a intersection *.gff # 결과물: pan_genome_results
   # difference
   # complement
@@ Line 95: / Line 97: @@
 {{ :svg.png?400 |}}
+=== Strain-specific gene 찾기 ===
+R에서 gene_presence_absence.Rtab과 gene_presence_absence.csv 두 파일을 다루면 된다. 다음의 사례에서는 Lb_1-46 균주에서 특이적인 유전자의 id를 추출하는 사례를 보여준다. Lb_1-46은 데이터프레임으로 읽어들이면 Lb_1.46으로 바뀌는 것에 유의해라. 구글링을 잘 하면 이를 원래 이름 그대로 유지하는 방법이 있다([[https://blog.genoglobe.com/2018/12/r-readtable.html|하루에 한 R 관련글]]).
+  > dat = read.table("gene_presence_absence.Rtab",sep="\t",header=T,row.names=1)
+  > dat$Lb_1.46
+  > dat.s = subset(dat, rowSums(as.matrix(dat))==1)
+  > Lb.s = dat.s[which(dat.s$Lb_1.46==1),]
+  > Lb.s.genes = rownames(Lb.s)
+  > dat.2 = read.table("gene_presence_absence.csv",sep=",",header=T,row.names=1)
+  > dat.2$Lb_1.46
+  > dat.f = dat.2[which(rownames(dat.2) %in% Lb.s.genes),]
+  > dat.f$Lb_1.46
+  > write.table(dat.f$Lb_1.46,"out.txt",sep="\t",quote=F,col.names=F)
 === 기타 해결할 문제 ===
 결과 파일을 열어보면 일부 단백질의 ID가 변형되어 쓰인 것을 알 수 있다. 즉 원본 annotation file에서 MT_RS20470라는 locus tag을 갖는 유전자가 "MT_RS20470.p01_16540"으로 바뀐 것이다. .p01은 유전자에서 단백질로 번역됨을 나타내는 것이겠지만 "_16540"은 무엇인가? (원래 '_'가 3회 반복된 것이나 DokuWiki에서 그렇게 타이핑하면 긴 가로선이 나오기 때문에 부득이하게 하나만 타이핑하였다)
+Prokka에서 만든 gff3 파일을 사용하였더니 tRNA gene을 결과물 중에 포함시키는 현상이 가끔 관찰된다.
 === Roary 이후 개발된 프로그램 ===
@@ Line 142: / Line 158: @@
 LS-BSR의 간단한 사용법을 알아보자. genome sequence 파일은 확장자가 .fasta가 아니면 작동을 하지 않는다. markers.fasta는 항상 유전자
  염기서열 파일이어야 한다. 검색용 프로그램은 기본이 tblastn이며(단백질로 번역을 먼저 함), 만약 blastn을 쓰고 싶으면 -b blastn을 더해야 한다. -x test라고 지정하면 test라는 디렉토리를 임시로 생성하여 중간 계산 파일을 넣은 뒤, 최종
- 결과 파일은 'test_'로 시작하는 접두사를 갖고 만들어진다. 분석이 끝나면 test 디렉토리는 없어진다. /home/hyjeong/github/LS-BSR을 사용해라. --- //[[hyjeong@kribb.re.kr|Haeyoung Jeong]] 2018/04/26 15:52//
+ 결과 파일은 'test_'로 시작하는 접두사를 갖고 만들어진다. 분석이 끝나면 test 디렉토리는 없어진다. /usr/local/apps/LS-BSR/ls_bsr.py을 사용해라. 비록 conda를 이용하여 LS-BSR을 설치했다 하여도 이는 environment만 구성할 뿐, 실제 작동 스크립트 등은 git로 깔았기 때문이다. --- //[[hyjeong@kribb.re.kr|Haeyoung Jeong]] 2019/06/27 11:22//
   # gene screen method (peptides in interest are ready)
-  # /data/apps/LS-BSR/ls_bsr.py <= 스크립트 위치
   $ python path_to_LS-BSR/ls_bsr.py -d genomes -g markers.fasta -x test
   # de novo gene prediction method
   $ python path_to_LS-BSR/ls_bsr.py -d genomes -c usearch -x test
-gene prediction method에서는 prodigal로 유전자를 예측한 뒤 단백질로 번역한 다음 usearch로 클러스터링(default: -i 0.9)하여 사용한다. 클러스터링 방법은 usearch/vsearch/cd-hit이 가능하다($PATH에 있어야 함).
+계산이 끝나면 현 디렉토리에 markers.fasta가 번역된 genes.pep 파일이 생긴다. gene prediction method에서는 prodigal로 유전자를 예측한 뒤 단백질로 번역한 다음 usearch로 클러스터링(default: -i 0.9)하여 사용한다. 클러스터링 방법은 usearch/vsearch/cd-hit이 가능하다($PATH에 있어야 함).