SILVAngs

QIIME2이후 새로운 녀석이 나타났다!!

silvangs음.. 왠지 

MG-RAST 대항마로 나온듯

MG-RAST야 분석 파이프라인만 제공하는것에 반해

silvangs의 경우 막강한 16S rRNA DB를 탑재하고 있으니...

그리고 silvangs 뒤에는 de.NBI (German Network for Bioinformatics Infrastructure)가 든든히 버티고 있는중 좋은 피드백이면 반영해주지 않을까? 하는 생각이..

그리고 silvangs도 일단은 credit으로 운영중 가입하면 300k주고 시작하는데
과금 방법이 분석횟수인지 분석에 필요한 스텝기준인지 아직
안사용해봐서 잘 모르겠다능.. :)

여튼 단점은 rDNA용 커뮤니티 분석이라는것이 한계

WGS가 플랫폼으로 나오면 대박일듯...  물론 WGS이 나오기에는 아직
넘어야할 산이 좀 있으니...
불안불안한 MG-RAST에 이어 web 기반 분석 플랫폼이 나와준것에 박수를..

HaploMerger2

HaploMerger2

HaploMerger2 논문

HaploMerger2 다운로드 사이트 , HM2 메뉴얼


2012년 HaploMerger1 개발 이후 2014년에 성능이 향상된 HaploMerger2로 2015년 11월자 파일이 최신인 프로그램으로 간단히 정리하면 polymorphism으로 인해 assembly가 제대로 되지 않는 genome을 assembly하여 고 퀄러티의 genome을 얻는데 도움을 주는 도구 되겠습니다. :)


기본적인 작업 프로세스는 아래의 순서도를 따릅니다.

HM2에서 제공되는 batch 스크립트를 보시면 위에 모식도를 그대로 구현해 놓은것이라고 보시면 됩니다.

HM2 설치는 다른 것들과 달리 컴파일하는거 없고 그냥 다운받아서 압축 파일만 풀면 되는데 run_all.batch파일이 속을 썩일수도 있습니다.

run_all.batch파일안에  치환해야 하는 파일 이름들 중 mp_ref/alt_ss이나 mp_ref/alt_ss_re이 것을
mp_ref_ss와 mp_alt_ss로 분리해주시면 중간에 에러나서 실행되지 않는 step들이 해결될것입니다. 그리고 추가적으로 우리 kent옹께서 만드신 faToNib도 필요하고 여기저기 숨어있는 풀어야할 문제들이 약간씩은 있을 수 있으나 log파일 확인하시면 다 나오는것이니 너무 걱정하지 마시고 일단 run 해보시면됩니다. libraries안에 gapCloser.cfg와 sspace_libraries.list는 GapCloser와 SSPACE할때 필요한 라이브러리 설정 파일입니다. 그 안에 형식에 맞게 fastq파일 넣어주시면됩니다.


마지막으로 HM2가 많이 궁금하신분은 제가 아는 선에는 HM2를 가장 잘 다루는 사람은 (주)테라젠이텍스의 박신기대리님이신듯!! 궁금한것 박대리님께!! 물어보십시요!! ㅋ :)

HGAP 버전 비교

PacBio에서 assembly를 담당하고 있는 HGAP과 Falcon을 주로 사용하고 있는것으로 알고 있는데 오늘 그 중에서 HGAP을 잠시 확인 하도록 하겠습니다.
왜 HGAP만 확인하느냐?
제가 HGAP을 많이 쓸 예정이거든요 ㅋ
Falcon궁금하시면 직접 하시면됩니다!! :)

우선 HGAP 버전을 한번 살펴보겠습니다.

HGAP버전을 보시면 3개가 존재합니다. 그러나 HGAPv1의 경우 이미 이번 SMRT버전에서 퇴출당했습니다. 당시에는 long read로 어셈블리해서 획기적으로 계산량을 줄이긴 했으나 Celera Assembler 자체가 시간이 많이 걸리기 때문에  (CA자체라기 보다는 특정 모듈들이겠죠?) 자체 프로그램을 고안해서 HGAPv2외에 추가로 HGAPv3까지 내놓은 상황입니다.

근데 보통 New버전 나오면 예전거 갈아치우는데(한국이 그런 경향이 많긴하지만.. 모 특정 위치에 있는 분들은 상당히 잘 안바뀌긴하죠..) HGAPv2와 HGAPv3를 계속 사용하고 있습니다. 두개의 큰차이는 consensus 부분인데 이것에 따라 속도와 결과물 차이가 나기때문에 두개 결과중에 좋은거 사용해라하는 느낌이 있습니다.

그래서 현재 HGAPv2와 HGAPv3에 대해서 테스트를 계획중이고
각 단계별로 parameter를 바꿔주면 어떻게 바뀌는지 비교해볼 계획입니다.
PacBio를 사용하시는 업자분들께서는 어떤 버전을 최적화해서 사용하는지는 잘 모르겠으나 아마 각각의 장단점이 있을겁니다. 그래서 직접 테스트를 해보겠다능~ :)

그럼 HGAP 비교글은 내년에 좀 더 볼만한 내용을 찾아뵙도록 하겠습니다. :)

getopts.pl 이 없다면!!!

perl script를 실행시킬때 Can't locate getopts.pl 이라는 에러가 났을때

당신의 OS가 우분투라면?

apt-get install libperl4-corelibs-perl

로 해결가능합니다. orz

16S rRNA와 시퀀싱플랫폼

블로그에 쓰는 내용이 16S rRNA에 많이 집중되고 있긴하죠? ㅎㅎ
하는일이 이거다보니.. :)

여튼 오늘은 16S rRNA와 시퀀싱플랫폼에 대해서 잠시 이야기 하도록 하겠습니다.

최근 16S rRNA sequencing의 최강자 454가 서비스를 bye bye한 관계로
많은 연구자들이 MiSeq체제로 변환하고 있는데 (물론 Ion도 있고, PacBio도 있습니다.)
기존에 454를 사용했을 때와 다른 V region을 사용하고 다소 다른 결과들을 보이는 것들이
있을 것 입니다. 그래서 이리저리 검색하다가 걸린 논문 두 개를 가지고 잠시 얘기해보고자 합니다.

Nucl. Acids Res. (2010) 38 (22): e200.

PeerJ (2016) 4:e1869

논문을 찾아보기 시작한 이유는 454가 막을 내린 후 다른 시퀀싱 플랫폼에서는 왜 다른 region을 target하고 있고 왜 diversity에 차이를 보이는지... (내가 분석을 잘못했나.. ㅎㄷㄷㄷ)

우선 V region, 454는 V1-V2였는데 MiSeq은 V4 region을, Ion은 4개 region? 7개 region? 을 동시에target하고 있다는..

모 논문보면 아시겠지만
2010년 논문은 454 vs MiSeq을 비교했습니다. 두 플랫폼이 차이가 날까? 어떤것이 차이가 날까?

2010년도 논문을 한장 figure 요약하자면 이거죠

상단은 phylum abundance/ 하단은 genus abundance 그리고 좌측에 있는 V4는 이전 연구에서 사용되었던 기준이라고 생각되는 참고용 분포입니다. 실험 결과가 V4 region이 저 분포를 나타내면 실험이 잘됐다고 확인하는 용도로 사용됩니다.

사실 phylum은 크게 차이없죠, 차이가 있기도 쉽지않습니다.
관심사는 genus되겠습니다. 일루미나 데이터의 대부분은 unclassified입니다. 논문에서는 error때문에 이렇게 나왔으니 error좀 낮아지면 일루미나가 output이 많으니 sequence error 문제점을 개선된다면 미생물 분석에 적합한 킹왕짱 시퀀서가 될거라고 하는데 개인적인 생각으로는 sequence error(454가 일루미나가 한테 sequence error ㅋㅋ 좀 웃겼다능..
systematic error는 눈에 안뵈냐라고 한다면 눼눼, 하긴 이 논문이 2010년이란 것을 감안한다면 무리는 아니긴 합니다.)보다는 db선택이 unclassified문제는 보정할 수 있지 않을까 합니다.

그렇다면 이제 16s rRNA 입문한지 얼마안된 님께서 그렇게 느끼는 느낌적인 이유는
무엇인가?

바로 2016년 논문되겠습니다.

2016년 논문은 PacBio로 시퀀싱한 것을 db에 따라 분석 결과가 달라지는냐에 대한 내용으로
다음 한 장으로 요약 할 수 있겠습니다.




동일한 시료를 가지고 시퀀싱한 데이터를 16S rRNA db에 따라 classification되는 정도를 확인해본 그래프 입니다. 느낌 딹오시죠?
다만 나중에 뒤통수 맞았다는 느낌 안 받게 하나 말씀드리자면 RDP와 Silva의 경우 genus까지만 확인했고 gg는 genus와 genus+species 두가지로 확인한 것 입니다.

이제 PacBio의 Sequel에서 CCS로 생산된 16S rRNA 서열가지고 연구해도 나쁘지 않을 것 같다는..

아.. 이제 돈만 많으면 되는건가!!! 

MrBayes 설치

MCMC방법으로 Phylogeny Tree를 그려주는 MrBayes!!
(컴파일해서 나오는 실행 파일명이 그닥 맘에 들지 않는 tool중 하나 ㅋㅋ)

자세한건 물어보시지 마시고 요기서는 설치만!
(그렇다고 다음에는 자세한거 말씀드리지는 않는다능)

자 설치를 위해서는 요기에서 말하는 Library는 미리미리 설치하면 암에 안 걸리니 먼저 잘 설치해주시기 바랍니다.


CentOS

yum install make,automake,gcc,gcc-c++,kernel-devel,git,autoconf,automake,libtool,subversion,pkgconfig,java-1.6.0-openjdk-devel,openmpi,openmpi-devel

CentOS에서 Autoconf 버전으로 인해 문제가 발생했을 때는 그냥 컴파일 하세요. :)

$ wget http://ftp.gnu.org/gnu/autoconf/autoconf-2.69.tar.gz
$ tar xvfvz autoconf-2.69.tar.gz
$ cd autoconf-2.69
$ ./configure
$ make
$ sudo make install

Ubuntu

apt-get install build-essential,autoconf,automake,libtool,subversion,pkg-config,openjdk-6-jdk,git, openmpi-bin,openmpi-doc,libopenmpi-dev

그리고 추가적으로 beagle를 설치하셔야 할겁니다.

Beagle 설치방법은 요기를 참고하시고 대략 밑에처럼 하시면 됩니다.

git clone --depth=1 https://github.com/beagle-dev/beagle-lib.git
cd beagle-lib
./autogen.sh
./configure --prefix=$HOME
make install

마지막으로
MrBayes 소스파일을 다운로드 받으시고 설치해 주시면 되겠습니다.

$ autoconf
$ ./configure --prefix=/install/to/path --enable-mpi=yes
혹은
$ ./configure --prefix=/install/to/path --enable-mpi=yes --with-beagle=/install/to/path
$ make

QIIME2

QIIME2 출시(된지는 한달여정도... 지났지만..)

설치만 해두려고 qiime2 설치하다가 환경 설정이 꼬였는지
쓰고있던 qiime 1이 작동안해서 멘붕이 왔다는건 비밀;;; 빼고는
아직까지는 호감인 녀석!

여튼 이번 qiime2 출시하면서 변화를 준게 완전히 qiime 플랫폼을 만든듯..
통상적으로 사용하였던 파일 포맷이 아닌 아예 새롭게
qza/qzv를 만들어서 사용한다는... 음..
qiime2를 설치하지않은 사용자들을 위한 qza/qzv view사이트까지
(모하자는거지... 유료화로 바꿀건가)

그리고 CLI기반이 아닌 q2studio라는 GUI 인터페이스를 제공한다는
기쁜소식!! (인지 아닌지는 까봐야..)
GUI는 자동화하기 귀찮은데... 제길..

qiime1보다 report가 얼마나 세련되졌는 확인할 필요가 있고
그에 따라 유료화를 점쳐볼수도...
그러므로 나는 일단 소스 다운로드;;;;

16S rRNA Variable regions

Variable regions of the 16S ribosomal RNA.

출처: Nat Rev Microbiol. 2014 Sep;12(9):635-45. doi: 10.1038/nrmicro3330.

PacBio SMRT Portal 설치

우선 서버사양은 Dell Workstation인데
CPU대비 메모리가 풍족해서 아직 메모리로 문제가 된적은 없는
Ubuntu 16.04.1 LTS 로 생정보하는데 필요한 프로그램이 설치되어 있고
관련 라이브러리가 설치된 조건하에서
SMRT Portal 설치를 시도하였고 전혀 Error없이 설치!! :)

SMRT Analysis Download 페이지에서 하단의 run 파일 두개 다운로드
무조건 2개 모두 다운로드!!

설치필수파일: smrtanalysis_2.3.0.140936.run

패치파일: smrtanalysis-patch_2.3.0.140936.p5.run

다운로드 후에
smrtanalysis 유저와 동일한 이름의 유저그룹도 생성한 후 설치 진행하시면되겠습니다.

$su -l smrtanalysis
$bash smrtanalysis_2.3.0.140936.run -p smrtanalysis-patch_2.3.0.140936.p4.run --rootdir
 /path/to/smrtanalysis

설치에 필요한 프로그램들 번들로 설치하고 기존에 mysql 서버 사용중이라도 영향을 미치지않게 포트도 별도로 사용합니다. tomcat의 경우는 기존에 설치되어 있는 서비스가 있다면 영향을 미칠지 명확하지는 않습니다.(아마 포트를 변경해서 영향이 최소화되도록 진행될것 같습니다.)

설치 후 localhost가 아닌경우 방화벽 포트 열어주시면 웹페이지로 접속하여
멋있는 SMRT portal home 화면을 보실수 있으실겁니다. :)

ufw allow [설치할때 보였던 port들]

참고문헌 [Installation SMRT Portal]

rRNA database 간단 비교

근래 microbiome데이터를 다루면서
16s rRNA db 사용에 대해서 확인을 하고자
데이터가 꽤나 많은 그리고 non-academic의 경우 fee도 지불하는 silva를 qiime의 default db인 greengene과 한번 비교를 해봤습니다.

rdp말고 왜 silva를 비교했느냐?
rdp는 qiime에서 사용가능한 format으로 데이터를 제공안해서
그냥 테스트로 하는건데 format 변환하는거마저 신경을 쓰기에는
좀 귀찮아서 qiime용 format을 제공하는 silva를 사용하였습니다.

모 별거 없습니다. ㅋ

	GreenGene	Silva 119	Silva 123
Kindom	2	3	3
Phylum	90	95	113
Class	249	480	602
Order	405	1,468	1,686
Family	514	3,123	3,670
Genus	1,816	8,056	9,835
Species	1,651	31,337	42,974
Sum	99,322	173,838	251,764

위의 GreenGene는 qiime 설치시 함께 따라오는 기본 파일을 사용했고, silva는 silva측에서 qiime format으로 제공하는 데이터중 SSU 파일입니다.
리눅스에서 cut 명령어로 -d, -f 옵션을 사용해서 unique한것을 세어본것입니다.

GreenGene과 silva의 차이점은 GG는 archaea, bacteria만 있고, silva는 Eukaryota가 함께 있다는.. silva는 엄연히 따지자면 16s/18s를 모아놓은 것으로 생각하시면되겠습니다. :)

그리고 마지막에 Sum의 의미는 이 rRNA db를 이루는 classified/unclassified 되는 species(?)의 unique 개수입니다.
이 unique 개수가 많으면 좋을 것 같은 느낌이 본능적으로 딱 오죠?

근데 말입니다.

숫자가 크다고 다 좋은것일까요?

GG와 silva에서 큰 차이점이 GG의 경우 unclassified로 분류되는 항목은 딱 하나씩 포함되어 있습니다. (그 이유는 숙제~!!)
대신 silva에는 unclassified와 같은 부류인 아직 의미가 없는 즉 정보가 없는 항목들이 다수 포함되어 있다는것입니다. unclassified, uncultured, unknow....

무엇이 우월하다 좋다 말하기는 어렵습니다.
분석 툴에 따라 궁합이 맞는 db가 존재하기 마련이고 또 목적에 따라 db는 변화해야하기에
GG는 더 curation을 했으니 좋은것이다 혹은 silva가 더 seq정보가 많으니 좋은거다라는
것은 지양하고 결국 각자가 알아서 잘 사용하면 된다능.. :)

seqtk

Site: seqtk

최근에 접한 fasta/fastq파일 핸들링 툴로 

가벼움과 사용이 쉬운것에 무한한 감사를...

사이트를 방문하시면 나와있듯이

$git clone https://github.com/lh3/seqtk.git

$cd seqtk

$make

하면 설치가 성공!! 이제 당신은 seqtk 마스터!!

잘쓰세요~ :)

ps. 전 subsampling하느라 사용했다능.. :)

16S rRNA Workflow

요즘하는게 요 16S rRNA작업인데
workflow와 해당 작업에서 사용할만한 tool을 잘 정리한글이 있어서..

H3abionet.org

덤으로 GOBLET이라는 사이트가 있는것도 알게 되었네요..
GOBLET이 뭔지 궁금해요?
그럼 링크에 들어가시면됩니다. :)

pie charts 어디까지 그려봤니? Krona

요즘 논문 낼때 고퀄의 figure들은 이제 당연한시대에 살고 있는데
나도 그런 고퀄 figure만들어보고 싶은데 잘 안된다.
네 저도 힘듭니다. 근데 다들 하니 저도 하는거죠
그 중에 독보적인 것은 circos이죠 ㅋ
그리고 pie chart를 좀 이쁘게 그리고 싶을때 사용하는것이 이것
Krona (OMIC-TOOLS에도 소개 되고있습니다.) 되겠습니다.

설치는 모 그닥 어렵지 않습니다.

Krona github가셔서 다운로드 받고 압축풀면 그걸로 끝입니다.
다만 Krona의 경우 metagenomic의 visualization을 위해서 만들어져 있는 관계로

taxonomy정보를 만들어줘야 합니다.
물론 Krona에서 해당 작업을 해주는 스크립트가 존재하니 걱정하지 마시기 바랍니다.
거기다 taxonomy 정보는 unique하고 찾기 힘든 data가 아닌 ncbi taxonomy정보입니다. 걱정마세요 :)

사용법은 제가 BLAST output을 사용한 pie charts 밖에 그려보지 않아서
Krona중에 ktImportBLAST밖에 사용해보지 못해서 말씀드릴게 별로 없네요 ㅎㅎ
-설명서를 잘 보시면 EC 카테고리별로도 pie chart를 작성할수도 있습니다.


다음은 제가 그냥 편할려고 만든 bash입니다.

#!/bin/bash

BLASTN=/path/to/ncbi-blast-plus/bin/blastn
KRONA=/path/to/Krona/bin/ktImportBLAST

ARGS=2

if [ $# -eq "$ARGS" ];
then

        $BLASTN -query $1 -db blast_database -outfmt 6 -num_threads 8 -evalue 1e-5 -out $2.outfmt6.blastn
        $KRONA -o $2.html $2.outfmt6.blastn
else
        echo "krona_BLAST.sh query.fa prefix_name"

fi
exit 0

snpEff에서 사용하는 genome annotation은 어떻게 추가할 수 있는가?

snp annotation에 빈번히 사용되는
snpEff (snpSift는 난 모르겠고)

snpEff에서도 나름 최신의 genome 정보를 제공하고 있지만
내가 de novo진행한것은?? (제가 곰팽이 de novo들을 많이하다보니.. ㅋㅋ)
어떻하라는 말인가...

고갱님 genome과 gene정보를 NCBI에 등록하고 NCBI gff가 공개되서 snpEff팀에서 지원해주는 시점에 사용하시면 됩니다.

근데 우리 고갱님들 그때까지 기다리시면 암걸리시죠?

그래서 간단한 가내수공업만 할 줄 아시면 곧바로 작업 가능합니다.

일단 위에 소개된 snpEff 사이트에서 snpEff 다운받으시고
압축푸시면 되겠습니다.

그리고 이번에 새로 조립하신 complete든 draft든 genome의 서열과 gff(version 2/3 택일)파일을 /path/to/snpEff/data/ 폴더 밑에 genome 이름으로 폴더 만드시고 그 밑에 복사하시면 되겠습니다. 대신 genome과 gff파일 이름은 genes.gff, sequences.fa로 바꿔주시는 센스!!

ex) 새로 조립한 genome이름이 Lee girwon이라면 Lee_girwon이라고 만드시고 그 밑에 파일을 복사해주시면 되겠습니다.

그리고 추가적으로 하나더 해야 하는 작업은 snpEff.config파일 수정
/path/to/snpEff/snpEff.config파일 끝에 새롭게 추가할 genome을 추가해 줍니다.

ex) vi snpEff.config
Lee_girwon.genome : Lee_girwon
--다음 라인은 선택사항 입니다.--
[TAB]Lee_girwon.chromosomes : AAA0001.1, AAA0002.1
[TAB]Lee_girwon.AAA0001.1.codonTable : Standard
[TAB]Lee_girwon.AAA0002.1.codonTable : Invertebrate_Mitochondrial

자 snpEff.config에 필요한 정보를 추가하였다면 이제는 database를 만들어주는 시간입니다.

java -jar /path/to/snpEff/snpEff.jar build -gff3 -v Lee_girwon

하시면 snpEff database가 뚝딱 만들어 집니다.

참 쉽죠?

위에 까만글씨로 욕은 아닌데 욕먹은 느낌이 나서 귀찮다?
그럼 작업가능한 서버계정과 흡족 할 사례비주시면 대행해드립 ㅋㅋ

CUTADAPT 사용설명서

cutadapt는 현재 1.9.1이 최신판으로

cutadapt guide에 설치에서부터 사용법
A to Z까지 너무나도 자세히 잘  소개되어 있습니다.
2011년에 논문도 있고

여튼..

cutadapt --help하시면 잘나오고
이 URL을 보시면 옵션에 따라 어느 위치의 adapter를 제거해주는지
어떤 방식으로 제거해 주는지 잘 나와있습니다.

-a,-g,-b 옵션이 adapter를 확인하고 remove하는 위치를 선정하는 중요한 옵션입니다.

single-end인 경우에는 cutadapt 사용시 큰 문제는 없지만
paired-end인 경우 cutadapt 사용시 다음과 같이 사용해야 합니다.

$cutadapt -a adapt_seq -A adapt_seq -o out_1.fq -p out_2.fq in_1.fq in_2.fq
-g, -b를 paired-end에서 사용하려면 -g -G, -b -B를 사용하시면됩니다.
-u는 일정한 길이로 cutting하는 옵션인데 second에서도 적용되게 하려면 -U를 사용하시면된다고 합니다.

SQLite3에서 import/export 하기

SQLite3를 사용할때 csv로 import, export하는 방법에 대해서
알아보겠습니다.

이게 왜 필요한거냐구요?
저한테 필요한거라서 포스팅하는겁니다.
그리고 알아두시면 쓸일이 있을지 누가 아나요? :)

csv를 특정 테이블에 import하기

#sqlite3 <database>
sqlite> .headers on
sqlite> .separator ","
sqlite> .import <in.csv> <table.name>
이때 in.csv는 쉼표 구분자로 되어 있어야 합니다.
다른 구분자로 되어 있으면 아마 제대로 import가 되지 않겠죠? :)

특정 테이블을 csv로 export하기

#sqlite3 <database>
sqlite> .headers on
sqlite> .mode csv
sqlite> .output <output.txt>
sqlite> select * from <table.name>;
이때 output.txt는 export할때 내용을 받는 파일 이름입니다.
그리고 tab구분자를 사용하려면 아마 .mode tabs 하시면될겁니다.

그리고 sqlite3 데이터베이스안에서 작업하지 않고 다음과 같이 명령어 한줄로도 끝낼수 있습니다.

sqlite3 -header -separator " " ./home/data.db "select * from table;" > out.txt      
sqlite3 -header -separator "," ./home/data.db "select * from table;" > out.txt        
sqlite3 -header -separator quot;\t" ./home/data.db "select * from table;" > out.txt 

모 그렇다고 합니다. :)

다양한 DNA-Protein Interactions Seq

출처: core-genomics blog

위의 블로그에 잘 정리된 Chip Seq들
참 별거 다 있다는..  :)

Long-read sequence assembly of the gorilla genome

서부저지고릴라(Western Lowland Gorilla, Gorilla gorilla gorilla) 중 하나인 Susie의 genome이 PacBio를 이용해서 좀더 high resolution으로 만들어졌다는 아름다운 논문입니다.

Long-read sequence assembly of the gorilla genome

왼쪽 녹색은 Susie, 오른쪽 아이보리색은 gorGor3의 contig size

각각 전체 Genome에서 10% 서열을 나타내는데 사용되는 contig 개수를 보여주는 그림으로 PacBio로 시퀀싱하여 어셈블리한 Susie가 short read assembly로 하는것보다 월등함을 확인시켜주고 있다(300M를 보여주는데 susie는 contig개수가 10개 남짓이면 되는것에 비해 gorGor3는 세어보시길;;;).

그리고 첫장 Table1에서 기존에 short read assembly한 결과보다 이번 결과가 더 월등하다는것을 여실히 보여주고 있는데 스캣폴드 개수가 554개 무슨 곰팡이 contig 개수인줄..
contig 최대 길이 서열은 36M bp, scaffold 최대 길이 서열은 110M bp. orz

그밖에 논문에서
기존 genome에서보다 gap 더 줄였구요
기존에 짧게밖에 못봤던 mobile element들 거의 full length로 확인할수있었구요
수kb에 달하는 insertion 확인해서 유전자 없는것도 확인할수 있었습니다라는
다양한 잘난척을 시전해 주고 계시는데...


결국 사용한 SMRT cell이 236개라는...
이거 PacBio 시퀀싱가격만... ㄷㄷㄷㄷ

이 논문보시고 우리도 genome 향상시킬수있어!! 라고 핑크빛 바램을 가지고 있으시는분들..  여러분들도 원래부터 좋은 genome가지고 연구할수 있었습니다.
다만 연구비가 귀여워서 못한것 뿐이고 그리고 모든 동물에 대해서 이렇게
드라마틱하게 genome 품질이 향상되지는 않습니다.
척추동물정도면 이정도 연구비 때려부으면 가능하지만 그 이하에서는 아직
해결해야할 것들이 좀 있습니다.

다 아시는 분들께서 모르시는척 하시기는... :)

그리고 PacBio의 Sequel 출시로 기존에 RSII로 했을때 보다는 반값에 가능하지 않을까합니다.
일단 SMRT cell개수를 줄일수 있으니... ㅋㅋ
그거 노리고 일단 RSII기준으로 시퀀싱비용 비싼듯보이게하고 Sequel로 하면 싼것처럼 느끼게 하려는 고도의 노림수인가;;;

여튼... 잘 따져보시고 시퀀싱하시기 바랍니다.

너도나도 앞다투어 시퀀싱하면 거지꼴 못면합니다.

Widespread Polycistronic Transcripts in Fungi Revealed by Single-Molecule mRNA Sequencing

Widespread Polycistronic Transcripts in Fungi Revealed by Single-Molecule mRNA Sequencing

간만에 읽은 저널 한편...
어느 지구정복을 꿈꾸시는 과학자분께서 운영하시는 것에 비하면 그냥 트윗터 수준임을
미리 알려드립니다.

자세히 안파해칩니다.

필요한것만 읽습니다.(제목만 보고 대충 때려맞추겠다는 심본데? 정답!!)

곰팡이중 basidiomycete fungi를 Iso-Seq을 이용하여 transcriptome 분석을 수행하였고, 우리는 기존의 short read가지고 깨작대던 님들이 못찾는거 찾는 기승전시퀀싱자랑하는 논문되겠습니다.

일단 이 논문은 JGI와 Pac이 손잡고 만들었습니다.
이 말은 곧 SMRT 비용 신경안쓰고(는 아니고 다른 연구자들보다 적게 신경썪을..) 분석에 사용하기 좋은 고 퀄러티 read들을 넉넉히 생산했을 것이기에 이런거 하고 싶다고 그냥 무작정 논문에 나온 SMRT cell 만큼 시퀀싱하시면 거지꼴 못면합니다.
suppl보시면 아시겠지만 SMRT cell 두자릿수 입니다. 앞자리가 10이 아닌건 안비밀 Orz..

복잡하고 어려운거 직접 보시면되니깐 쉬운거 말씀드리고 끝내겠습니다.

기존에 분석한 basidiomycete fungi중 Plicaturopsis crispa를 집중으로 파해쳤는데 이전까지 알고 있는 isoform 비율 10%가 아닌 한 20%정도 된다. 그리고 곰팡이도 isoform 3개 이상짜리도 엄청 많이 있음. 우리 곰팡이 무시하지 마셈.

그리고 비교셋으로 일루미나 숏-리드도 시퀀싱해서 ToFU (Transcript isOforms: Full-length and Unassemble,의 약자로 iso-seq을 분석하는 파이프라인? 시퀀싱 전략? iso-seq 결과물? 논문 보시면 아시겠지만 다양하게 사용되는것을 알수 있음)와 비교해봤는데
기승전 풉 짧은 것들은 안됨.
(하..... 지금까지 짧은것가지고 한것도 서러운데.....  ㅠ.ㅜㅋ)

여튼.. 모 검증은 해봐야하는거고 이 논문에는 RT-PCR해서 polycistronic 검증을 하긴했는데 좀더 확인해봐야 할것 같고..
새로운것 찾았다고 하고싶은 분들은 Iso-Seq 관심가지고 해보시는것도 나쁘지 않을것 같습니다.

근데 왜 너는 PacBio 관계자도 아니면서 PacBio로 실험한 논문 소개 하냐?

제 비록 제가 몸담고 있는 곳에는 PacBio가 주력이 아니지만
저는 언제나 연구자분들이 좀더 멋진 연구를 하실수 있도록 아낌없이
조언을 해드리고자 고심하는 연구자아니깐요 (캬~ 멋있다)

InterProScan Error를 대하는 방법

근래 interproscan을 돌리는데
어떤 옵션은 문제없이 작업이 진행되는데
어떤 옵션은 작업중에 주구장창 에러만 뱉어내서
내가 무엇을 잘못하는지 어리둥절하고 있었는데 은근 쉬운 문제였다는...

우선 제가 하는 작업은

1) interproscan.sh -i input.fa -appl pfam -f tsv
2) interproscan.sh -i input.fa -goterms -f tsv

1)번 command는 문제업이 계속 잘 진행
2)번 command는 미친듯이 자바에러를... ㅠ.ㅜ

그런데 역시나 그렇듯이
에러를 자세히 보면 해결할 수 있다는 :)
(문제속에 답이 있다는 드립은 아니고;; )

제가 겪은 오류의 문제는 hostname 설정만 해주면 된다는 사실
이유는 -iprlookup 때문인데.. Host를 찾을 수 없어서 쉰나게 에러를(아이 씬나~ )
그래서 /etc/hosts에 127.0.0.1에 그냥 알아서 잘 추가해주면 된답니다.

interproscan 돌리실때 에러난다고 속상해하지 마시고
에러를 잘 확인하시면 그 속에 길이 있을수 있습니다.
ps. 세그먼트 폴트가 난다면 그냥 때려치우세요..

Kmer 분석 관련 논문

간만에 뒤적이다가 Kmer 관련해서 설명된 2013년도 논문이 있어서 아실분은 아시겠지만
공유차 끄적끄적 합니다.

Estimation of genomic characteristics by analyzing kmer frequency in de novo genome projects

논문요약을 하자면

NGS의 발전으로 re-seq말고도 de novo를 하는데 short read라서 우리가 genome assembly를 잘했는지 "아 몰라요" 하지 않기위해 서열의 content frequency를 사용해서 예측하는 k-mer 분석을 해서 대략적으로라도 genome size를 확인해서 우리가 "assembly 쩔 잘했어요"를 외치는데 도움이 되는 테스트도 해봤고 소스도 있습니다. 우리 잘했죠?

모 그런 내용 되겠습니다. :)

가장맘에 들었던 figure는 error와 heterozygous에따른 Kmer 분포 였습니다.


NGS 특히 일루미나 데이터를 가지고 de novo 해보셨던 분들께서는
아마 다음과 같은 그림을 보셨을겁니다.

맨날 보는 그림이고 이그림은 시퀀싱 error를 포함하고 있는 그림입니다(1-2 depth에 peak가 무제한인게 가당키냐 한 얘기냐고요 >_< ).

다음이 이상적인 시퀀싱이 된 결과 에서 얻어질수 있는 원래 kmer 그래프 입니다(이것도 가당키나 가능한 얘기냐??).

그리고 하단의 그림은 위의 그림들을 merge하여 잘 표현해놓은 그림입니다.

그리고 마지막으로 Kmer분석하면 진짜 genome size 예측 제대로 되는거 맞어?
라는 의구심을 가지시는 분들께서는 논문 p25에 Table2. b 테이블보시면 되겠습니다.

현재 genome size가 보고된 종들에 대해서 Kmer분석을 해봤고, 어떻게 하면 그나마 실제 genome size에 가까운 Kmer 분석 결과를 얻을 수 있는지 아실수 있습니다. 그리고 보시다 싶이 genome size큰거는 (사람은 150Mb 차이가 나는데) 예측 잘 안되는거 아니냐 라고 할수 있지만 %로 보면 오차가 5%내외이니 크게 걱정하지 않으셔도 됩니다(아마추어 같이 왜그러십니까).

그리고 논문 마지막에 K-mer 분석하려면 (인간적으로) genome size에 30x 이상은 시퀀싱하자 라고 촉구하고 있습니다.

사족: 논문 맹신하지 말아라. 겪어보면 알겠지만 케바케다.