User Tools

Site Tools


cog_assignment

Differences

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

Link to this comparison view

Both sides previous revisionPrevious revision
Next revision
Previous revision
Next revisionBoth sides next revision
cog_assignment [2018/02/19 08:40] – [자료 다운로드] hyjeongcog_assignment [2019/03/14 15:20] – [Query sequence file의 준비] hyjeong
Line 38: Line 38:
 각 디렉토리로 들어가서 make를 실행하면 디렉토리 이름에 해당하는 실행파일이 생긴다. gcc 4.4.7로는 잘 빌드되지만 gcc 5.3.0(Linuxbrew)로는 오류가 발생하였다. 만들어진 실행파일을 $PATH에 위치시킨다. 각 디렉토리로 들어가서 make를 실행하면 디렉토리 이름에 해당하는 실행파일이 생긴다. gcc 4.4.7로는 잘 빌드되지만 gcc 5.3.0(Linuxbrew)로는 오류가 발생하였다. 만들어진 실행파일을 $PATH에 위치시킨다.
  
-COG data 파일 중 당장 필요한 것은 다음과 같다. /data/Utilities/DB/COG/COG2014에 있는 것을 현재의 작업 디렉토리로 가지고 온다. 아니면 이 디렉토리로 query sequence file을 복사하여 거기에서 작업을 해라! fun2003-2014.tab(one-letter functional category)과 cognames2003-2014.tab(COG 번호, function code, COG 명칭) 파일은 COG 계산 작업에 직접적으로 필요하지는 않다. +COG data 파일 중 당장 필요한 것은 다음과 같다. /data/Utilities/DB/COG/COG2014에 있는 것을 현재의 작업 디렉토리로 가지고 온다. 아니면 이 디렉토리로 query sequence file을 복사하여 거기에서 작업을 해라! 이렇게 하는 것이 더 편리하다. 그 대신 이전에 작업한 흔적은 잘 지워야 한다. fun2003-2014.tab(one-letter functional category)과 cognames2003-2014.tab(COG 번호, function code, COG 명칭) 파일은 COG 계산 작업에 직접적으로 필요하지는 않다. 
-  * **prot2003-2014.fa**+  * **prot2003-2014.fa** BLASTDB로 전환된 다음에는 필요하지 않다.
   * cog2003-2014.csv => **COGs.csv**라는 이름의 심볼릭 링크를 만들든지 복사를 해라.   * cog2003-2014.csv => **COGs.csv**라는 이름의 심볼릭 링크를 만들든지 복사를 해라.
   * **COG.p2o.csv**: 다음과 같이 <protein id>,<genome id>로 만들어진 파일이다. cog2003-2014.csv 파일의 첫번째와 두번째 컬럼을 awk 명령으로 떼어내면 된다. 아래에 예를 든 awk 명령어는 별로 마음에 들지 않는다.   * **COG.p2o.csv**: 다음과 같이 <protein id>,<genome id>로 만들어진 파일이다. cog2003-2014.csv 파일의 첫번째와 두번째 컬럼을 awk 명령으로 떼어내면 된다. 아래에 예를 든 awk 명령어는 별로 마음에 들지 않는다.
Line 62: Line 62:
   $ grep '>' Gen2.fa | sed -e 's/^>//' | awk 'OFS=","{print $1,"Gen2"}' > file2   $ grep '>' Gen2.fa | sed -e 's/^>//' | awk 'OFS=","{print $1,"Gen2"}' > file2
   $ cat file1 file2 > GenQuery.p2o.csv   $ cat file1 file2 > GenQuery.p2o.csv
 +
 +파일의 구조를 간단히 그림으로 설명하면 다음과 같다. 그림 아래에서 **GenQueryquery.p2o.csv**라는 파일명은 **GenQuery.p2o.csv**로 고쳐야 한다. 그림 파일의 원본이 없어서 볼썽사납지만 그대로 둔다.
 +{{ :cog.png?600 |}}
 +
 +이상의 작업을 하는 것이 너무 성가셔서 Gen1,fa, Gen2.fa.. 파일을 인수로 주면 query.fa 및 GenQuery.p2o.csv 파일을 만드는 스크립트 **process_COG_query.pl**를 하나 작성해 보았다.
 +
 +  #!/bin/sh
 +  
 +  # usage : process_COG_query.sh query1.fa query2.fa ...
 +  
 +  for query in "$@"
 +  do
 +      cat $query >> tmp.query.fa
 +      prefix=${query%.*}
 +      grep '>' $query | sed -e 's/^>//' | awk -v gen="$prefix" 'OFS=","{print $1,gen}' >> tmp.file
 +  done
 +  
 +  sed -e 's/^lcl|//' tmp.file > GenQuery.p2o.csv
 +  
 +  awk '/^>/ {
 +            split($0, field, "|")
 +    if (field[1] == ">lcl")
 +        {
 +        print
 +        }
 +    else
 +        {
 +        gsub(">",">lcl|", $0);
 +        print 
 +        }
 +    ;next}
 +  {print}' tmp.query.fa > query.fa
 +  
 +  rm tmp.query.fa tmp.file
  
 Query protein set와 COG protein에 대해서 blast DB를 만든다. 사용한 fasta file의 이름과 DB name이 다른 것에 유의한다. 혼동을 피하기 위해서 이렇게 한 것이니 각자 자율적으로 결정해도 된다. prot2003-2014.fa에서 만들어진 DB는 일정한 곳에 보관하여 계속 재사용하는 것이 좋을 것이다. Query protein set와 COG protein에 대해서 blast DB를 만든다. 사용한 fasta file의 이름과 DB name이 다른 것에 유의한다. 혼동을 피하기 위해서 이렇게 한 것이니 각자 자율적으로 결정해도 된다. prot2003-2014.fa에서 만들어진 DB는 일정한 곳에 보관하여 계속 재사용하는 것이 좋을 것이다.
Line 76: Line 110:
 ===== COGnitor process ===== ===== COGnitor process =====
 ==== BLAST 검색 세 차례 ==== ==== BLAST 검색 세 차례 ====
-query의 자체 검색(filter OFF), 그리고 query -> COGs에 대한 검색(filter OFF and ON)을 실시한다. blastp가 아니라 psiblast를 이용한다는 점이 의외였다. Multiple processor를 쓴다면 -num_threads <NUM>을 해도 좋다. 세번째 psiblast 검색에 가장 많은 시간이 걸린다. 이를 단축할만한 획기적인 방법이 있으면 정말 좋을 것이다.+query의 자체 검색(filter OFF), 그리고 query -> COGs에 대한 검색(filter OFF and ON)을 실시한다. blastp가 아니라 psiblast를 이용한다는 점이 의외였다. Multiple processor를 쓴다면 -num_threads <NUM>을 해도 좋다. 두번째 및 세번째 psiblast 검색에 가장 많은 시간이 걸린다(-seg no로 필터링을 하지 않는 것이 좀 더 빠름. 당연하지 않은가? [[http://www.biology.wustl.edu/gcg/seg.html|SEG 프로그램]]을 실행하지 않으므로!). 이를 단축할만한 획기적인 방법이 있으면 정말 좋을 것이다. 아래 사례에서는 -num_threads int_value 옵션이 빠져 있으니 그대로 실행하면 시간이 대단히 많이 걸림에 유의하라.
  
   $ psiblast -query query.fa -db Query -show_gis -outfmt 7 -num_descriptions 10 -num_alignments 10 -dbsize 100000000 -comp_based_stats F -seg no -out BLASTss/QuerySelf.tab   $ psiblast -query query.fa -db Query -show_gis -outfmt 7 -num_descriptions 10 -num_alignments 10 -dbsize 100000000 -comp_based_stats F -seg no -out BLASTss/QuerySelf.tab
Line 99: Line 133:
   $ COGcognitor -i=./BLASTcogn -t=COGs.csv -q=GenQuery.p2o.csv -o=GenQuery.COG.csv   $ COGcognitor -i=./BLASTcogn -t=COGs.csv -q=GenQuery.p2o.csv -o=GenQuery.COG.csv
      
-모든 것이 끝났다. **GenQuery.COG.csv**에 최종적으로 할당된 COG 정보가 수록된다. 주의할 점은 하나의 Query protein에 대하여 복수의 COG가 부여될 수 있. 이러한 경우에는 5번째 필드인 cognitor-score를 가지고 1등을 선별해야 할 것이다.+모든 것이 끝났다. **GenQuery.COG.csv**에 최종적으로 할당된 COG 정보가 수록된다. 주의할 점은 하나의 Query protein에 대하여 복수의 COG가 부여될 수 있. 이러한 경우에는 5번째 필드인 cognitor-score를 가지고 1등을 선별해야 할 것이다.
  
   $ cat GenQuery.COG.csv   $ cat GenQuery.COG.csv
Line 111: Line 145:
 [[findBestFromCOGs.pl|findBestFromCOGs.pl]]이라는 스크립트를 사용하여 GenQuery.COG.csv를 처리하면 각 protein query에 대한 best COG를 출력한다. [[findBestFromCOGs.pl|findBestFromCOGs.pl]]이라는 스크립트를 사용하여 GenQuery.COG.csv를 처리하면 각 protein query에 대한 best COG를 출력한다.
  
-  $ findBestFromCOGs.pl GenQuery.COG.csv+  $ findBestFromCOGs.pl GenQuery.COG.csv > GenQuery.COG.csv.bestHit 
 +  $ cat GenQuery.COG.csv.bestHit
   AND37645.1,448,1,448,32900.1,COG0593,sinlge   AND37645.1,448,1,448,32900.1,COG0593,sinlge
   AND37646.1,378,1,378,17170.9,COG0592,sinlge   AND37646.1,378,1,378,17170.9,COG0592,sinlge
Line 123: Line 158:
 마찬가지로 Clustering data(Readme 파일 섹션 2.10)란, COG를 만드는 프로그램인 COGtriangle의 결과물이다. 바로 COG data package에 포함된 cog2003-2014.csv 파일 아니겠는가?  마찬가지로 Clustering data(Readme 파일 섹션 2.10)란, COG를 만드는 프로그램인 COGtriangle의 결과물이다. 바로 COG data package에 포함된 cog2003-2014.csv 파일 아니겠는가? 
  
 +==== Functional classification ====
 +COG의 각 functional class에 몇 개씩의 유전자가 분포하는지를 집계하고 싶다면 [[cogclass2018.pl|COGclass2018.pl]] 스크립트를 이용하라. 하나의 COG에 대하여 복수의 functional class가 있는 경우는 무작위로 하나만 선택한다.
 +
 +  $ COGclass2018.pl GenQuery.COG.csv.bestHit > COGclass
 +  ...
 +  COG5503 KV (random selection) ===> K # STDERR 출력
 +  COG5524 CT (random selection) ===> C # STDERR 출력
 +  
 +COGclass 파일의 끝부분에 다음과 같은 집계표가 나온다.
 +{{ :cog_class.png?600 |}}
 ===== 개선 아이디어 ===== ===== 개선 아이디어 =====
 Gen1.fa, Gen2.fa... 파일을 인수로 공급하면 일괄적으로 COG assignment를 실시하는 스크립트를 만들어 보자. Gen1.fa, Gen2.fa... 파일을 인수로 공급하면 일괄적으로 COG assignment를 실시하는 스크립트를 만들어 보자.
cog_assignment.txt · Last modified: 2021/08/10 12:51 by hyjeong