Cross-ID mapping 문제 해결? BioGSP

Web address : http://biogps.gnf.org/

Bioinformatics 연구에서 다양한 소스의 데이터를 통합해야 하는 경우가
빈번하다. 유전자 이름만 해도 NCBI 내에서만 Gene ID, Entrez ID, Genbank ID
등등 다양하고, 기관 마다 차이가 존재하기 때문에, 이러한 다양한 기관의
ID 들 간의 cross-mapping 은 필수적이다.

그러나 매번 같은 소소의 정보를 이용하지 않는 경우가 많다보니, 매번 서로다른
기관들의 ID 매핑이 어려운 과정은 아니지만 반복되어 짜증을 유발하는 경우가
많다.



이번에 Genome Biology 에 소개된 ( http://genomebiology.com/2009/10/11/R130 )
BioGPS 는 이러한 cross-mapping 문제 해결을 위한 web-server 다.
BioGPS 에서 커버하는 ID 소스는

* Gene symbol ( from various institution )
* GO
* Interpro
* Affymetrix ID

등이고, primary ID mapping 파일을 다운로드 받을 수 있게 제공하고 있다.

DEG selection 방법 어떤 걸 써야 할까? [ Comparison among nonparametric DEG selection methods ]

초창기 Microarray 연구에서 가장 큰 관심사 중 하나는 어떻게 하면 robust한 DEG ( Differentially Expressed Gene ) 을 찾아낼 수 있는가 였다.

가장 먼저 시도되었던 기본적인 방법들이 fold change 비교와 t-test 였는데, 이는 지금까지도 microarray data 분석의 가장 기본적인 방법으로 널리 이용되고 있다. 또한 MAQC( MicroArray Quality Control) 그룹의 보고에 따르면, 이 두가지 방법을 느슨한 stringency 를 주어 DEG를 뽑을 때, cross-platform, cross-laboratory microarray data에서 가장 재연성있고 안정정으로 DEG 셋을 얻을 수 있다고 하니, 가장 기본적이면서 또한 가장 중요한 DEG selection 방법이라고 할 수 있다.

이렇듯 가장 기본적인 DEG selection 방법인 t-test에도 약점이 있었으니, 그것은 t-test가 normal distribution을 가정한 parametric test라는 것이다. microarray data는 normal distribution을 따르지 않는 경우가 대부분이라고 할 수 있다. 자연스럽게 t-test의 이런 약점을 보완할 수 있는 방법들이 DEG selection 방법으로 적용되었는데, 대표적으로 non-parametric t-test와 rank sum test 를 들 수 있다.

오늘 이 posting 을 맘먹게 한 논문 ,Nonparametric identification methods for differentially expressed genes, bioinformatics, 18, 1454, 은 현재 Princeton 대학에서 교수로 재직중인 Olga G. Troyanskaya 가 Stanford 에서 박사과정 중 출판한 논문으로, 앞서 언급한 microarray data의 DEG selection methods 들 중, nonparametric method 3가지를 비교 분석한 논문이다.

논문에서 비교 분석 대상이 된 Non-parametric method 3가지는 non-parametric t-test, Wilcoxon rank sum test, Ideal discriminator method 되겠다. 이 3가지 방법을 simulated data 셋에 적용하여 TPR, FPR 을 측정하여 상대적인 비교우위에 대한 분석을 논하고 있는데, 결론부터 얘기하자면, Rank sum test 를 쓰는 것이 다른 방법들에 비해 stringent 한 결과를 주고, multiple permutation 을 통한 significance 를 측정할 필요가 없기 때문에 일반적인 상황에서는 rank sum test를 쓰는 것이 이 세가지 방법 중 가장 computationally efficient 하고, 상대적으로 보수적이지만 정확하게 DEG set 을 뽑아낼 수 있다는 것. 덧붙여 p-value 0.1 정도에서 sensitivity와 specificity 의 trade-off 에서 가장 만족스러운 결과를 얻을 수 있었다는 결과가 리포팅 되어 있다.

차후에 SAM,GSEA 등을 비롯한 다양한 DEG selection 방법이 소개되었으나, 기본적인 single gene based DEG selection analysis 는 기본적으로 필요한 과정이니, microarray 분석을 할 때 기 이를 염두에 두고 nonparametric test 를 적용해보는 것이 좋을 듯 하다.

이 논문이 출판되고 4년 후, MAQC 그룹이 출판한 Nature 논문 에 따르면 느슨한 t-test p-value 0.1 와 fold change 1.2~1.5 정도를 함께 사용하여 DEG set 을 선정하였을 때, 가장 만족스러운 결과를 얻을 수 있었다고 보고하고 있으니, 이 결과도 참고해둘 만 하다.

* Non-parametric t-test : normal distribution 을 가정한 상태로 p-value를 구하지 않고, permutation test에 의해 t-test statistics 의 distribution 을 만들고 이와 원 t-test statistics 값을 비교하여 p-value significance 를 구한다.

Random permutation generator

Bioinformatics 분야에서 널리 쓰이는 통계적 테스트 방법으로 Permutation test가 있다.

샘플 숫자가 통계적 테스트를 할 만큼 크지 않은 경우, 주어진 샘플을 permutation 하여
인공적으로 샘플 숫자를 늘려, 전체 population 을 통계적 테스트가 가능한 크기만큼 키운 후,
원래의 주어진 샘플에 대한 수치 ( e.g. mean, standard deviation , etc ) 가 전체 셋과 비교하여
얼마나 유의한지 ( significant ) 테스트 하는 방식이다.

permutation 할 수 있는 개수는 원래의 샘플의 개수에 따라 이미 정해져 있지만, 샘플 개수가
10개만 넘어가도 가능한 전체 숫자 만큼의 permutation set을 모두 고려하기 힘들다. 따라서
permutation set의 generation은 random 하게 통계 테스트가 가능할 만큼만 하면 충분하다.
보통 그 숫자는 1000개, 10,000 개 정도다.

최근 이 Permutation test 를 코딩해야할 일이 있어 코딩하던 중, CPAN의 random permutation
generator를 써봤는데, 이 모듈 사용이 직관적이지 않고, 쓸데없는 정보를 덧붙여 클래스 형태로 데이터를 format해야 해서 쓰기가 번거롭다.

그래서 random permutation 모듈을 하나 만들어봤다. 그냥 간단한 서브루틴 하나다.

sub rand_permutation{
   my @array=@_;
   my @permuted_array;
  
   while( @array ){
      my $rand_pos=int(rand(@array));
      push @permuted_array, $array[$rand_pos];
      delete $array[$rand_pos];
      @array=grep /\w/, @array;
   }
   return @permuted_array;

}

1,2,3,4,5 로 구성된 어레이를 넘기면 서브루틴이 random 하게 이 순서를 뒤섞어 준다.
그래서 결과로는 5,1,3,2,4 와 같이 random하게 permutation 된 셋을 리턴한다.

Dealing with gzip file in Perl

There are several ways to process gzip compressed files in Perl.
I tried to show some of them according to the chronological order
that I experienced.


Year 2005

Before I joined Jong Park's lab, I have no chance to deal with gzip compressed file processing. But the situation changed. The first data set in gzip format is PDB data files, which contain the information for 3-D structure of protein.

I solve the problem just by 'uncompressed' the gzip files !!
Simple and easy, hur?


Year 2006

One day, my lab member let me know there's linux command 'zcat', which 'cat' the gzip compressed file. So I just used that command in my Perl script.

 my @file=`zcat xx.gz`;
 while( @file ){ .. } 

Year 2007

I was not CPAN lover before 2007. So I didn't tried to find any solution from CPAN. But after I participated 'Agile programming education', I changed my mind and started to love CPAN as I used it.

The module I found at that time for this problem is 'PerlIO::gzip'.

use PerilIO::gzip;
open $fh,'<:gzip',"xx.gz";
while(<$fh>){.. }

Year 2009

Recently, I set up my new computer with Ubuntu 8.04 and 8.10 version. But there was a problem to install PerlIO::gzip module. It seemed that there's conflict between OS set-up and module's basic configuration. So I tried to find if there's another module to process gzip files in CPAN.

And there is as TMTOWTDI !

use IO::Zlib;
my $fh=IO::Zlib->new("xx.gz",'r');
while(<$fh>){ .. }

Symbolic reference 사용시 유의점

$x=10;
$sym='x';
$$sym=20;

print "$x\n";

위의 코드를 실행하면 $x 값으로 20이 출력되는 것을 확인할 수 있다.
Symbolic reference 가 적용되어 $sym = 'x'에 대해
$ { $sym } = $ { x} = $x 로 해석되어
$$sym=20 은 $x=20과 같이 해석되기 때문이다.

그럼 아래 코드는 어떨까?

my $x=10;
my $sym='x';
$$sym=20;

print "$x\n";

첫번째 코드의 결과와 같은 값이 출력될까?
정답은 NO!
이 코드는 그대로 $x에 대해 첫 라인의 정의대로 10이 출력된다.

해답은 'my' 과 symbolic reference의 동작 원리에 담겨있다.
my로 선언된 lexical 변수는 symbol table에 저장되지 않는다.
symbolic reference는 오직 lexical 변수를 제외한 package 변수에만 적용된다.


따라서 두번째 코드에서 변수 선언을 my 대신 local이나 our로 해주면 첫번째 코드와 같은 결과를 얻을 수 있다.

my $x=10;
$sym='x';
$$sym=20;

prnt "$x\n";

이 코드는 어떨까?
$$sym = ${x}= $x 로 해석되어 최종적으로 symbolic reference
의 적용을 받는 변수는 $x다. 그런데 $x가 my 로 선언되어 있어
symbol table에 저장되어 있지 않다.

따라서 $$sym=20 은 lexical 영역의 $x에 영향을 주지 못하고
이 코드는 결과적으로 $x=10 을 출력하게 된다.

그렇다면 $$sym=20 는 완전히 아무런 의미도 없는 코드일까?

그렇지 않다. $$sym 은 원래의 역할 그대로의 역할을 수행했다.
즉, symbol table에 존재하는 package global $x 에 20을
저장한 것이다.

이를 확인하기 위해서는 package 변수 $x를 출력해 보면 된다.

my $x=10;
$sym='x';
$$sym=20;

print "Lexical 변수 : $x \n";     # 10을 출력
print "Package 변수 : $::x \n ";  # 20을 출력

BioTool::NCBIfetch

Nowadays, I've tried to build my own perl pacakges selections on routine bioinformatics tasks.
Althrough BioPerl is already out there, it usually doesn't have modules for my own purpose.
For example, a package for learning PAML packages from A to Z, not requiring specific formatting to run the BioPerl. Well, so I started to build my own perl packages.

I named the base package name 'BioTool' and added packages whenever I confronted any situation
that I found some tasks might be used repeatly and routinely.

Recenlty, I made 'BioTool::NCBIfetch' package, which can extract sequences and all the related references db information including Gene symbol, Gene description, Chromosome location, Ensembl, Unigene, Uniprot, KEGG and GO for a given gene query.

Since the package fetches the information from result of NCBI eutils query, which results with up-to-date information, the analysis result for a given gene is up-to-dated. So users don't need to worried about whether the reference data is out-dated when they use a program which works the same taks based on localized data. In short, users could free themselves from updating all related indenpent databases day to day.

I'll show some example usage of this package

use BioTool::NCBIfetch;
my $ncbi=BioTool::NCBIfetch->new;
$fetch->set_query(gene,780,xml);  # Input : search DB, ID(geneID, gID), result type, result format type)

$fetch->get_result;
my $symbol=$fetch->get_symbol;
my @TreEMBL=$fetch->get_uniprotTreEMBL;

print "$symbol\t@TreEMBL\n"; 

This package might be very helpful for ones who want to cross mapping between major biology databases based on NCBI gene identifier. I actually made it for cross-mapping genes in GEO platform since GEO platform annotation is frequently incomplete and inaccurate.

If anyone interests in this pacakge and want to use, request!

Economic crisis & public data deposition in GEO

Does stunning global economic crisis starting from the late 2007 influence each gorvernment to cut down the budget for scientific research? If you're an independent researcher, you've already experienced how serious it is. If you are not, you don't know the actual influence of it on science research.

Starting with just interest, I tried to answer to above question by analyze the trend of data deposition to the largest public gene expression database, GEO. Since the number of gene expression data deposition has been increasing from the creation of GEO, I hypothesize that the data deposition might be not continuously increasing in year 2008 if the global economic crisis was really a matter for cutting down the budget for the science research.




Look at the table above! The number of data deposition to GEO in 2008 is slightly decreased as I expected. The number of data deposition has continuously increased from 2001 through 2007 but this trend stop in the year 2008.



This trend is repeately found for the top 5 countries of budget for science. By the way, USA is really a huge player in the basic science. Look at the size of data! US itself deposited several fold more data than the sum of all the data the other four countries deposited. Anyhow, back to the main point, I think this trend represent the effect of economic crisis on shirinking science budget in each gorvernment.




Except the USA from the chart, another interesting point shows up. That's UK. Deposition of data from UK was decreased from 2007 not 2008. Their contribution to public data was peak in 2006 but the number is less than a half in 2007. What happened in UK at that time? I don't have any evidence on this matter. Just I expect that UK changed their research plan for biology using microarray since microarray generally lack reproducibility and consistency. Or there might be some big scientific project to producing microarray data in 2006 only.

Among the 5 countries I mentioned here, Germany showed least decreasing nubmer of deposited
data. They deposited only 21 data set less compared to previous year ( 14% ) while the other four countries decreased that number around 50%.

The small contributing countries with less than 50 deposition in each year, some countries seem that they are not affected by economic crisis. South Korea is one of the exam. Their deposition increased from 16 in 2007 to 19 in 2008. As I'm researching in Korea, I think it's not because Korean gorvernment did not cut budget for science but most of researchers who producing microarray data don't know GEO at all or don't consider to deposit their data. Why? most of them are just wet-lab biologist. So they usually don't use public data for their own resarch. Even when they want to contribute, depositing data to GEO needs some computational skills for easy or automatic depositing. Without any programming skill, the process to deposit the data is very tedious and boring. So three more data sets of deposited data does not show the actual trend of Korean science budget related to economic crisis. I think the minor contributing countries might have similar situation.

That's all I want to say today! In one sentence, global economic crisis seems actually to force each gorvenment to cut their science budget.