2018년 2월 7일 수요일

NCBI PGAP이 생성한 공개 전 GenBank file을 Roary에 사용할 수 없는 이유

Pan genome 분석 파이프라인인 Roary는 유전체 서열 및 유전자(translated) 서열을 내부에 갖고 있는 GFF3 파일을 입력물로 받아들인다. 이를 위해서는 Prokka로 predictiongks GFF3 파일이나 NCBI에서 다운로드한 GenBank 파일을 bp_genbank2gff3.pl 스크립트(BioPerl에 포함)로 전환하여 사용하면 된다.

요즘 Paenibacillus polymyxa 및 이와 연관된 미생물 균주의 유전치 비교 분석에 힘을 쏟고 있다. 우리 연구그룹에서 유전체를 해독하여 NCBI의 PGAP으로 주석화를 하였지만 아직 공개가 되지 않은 것이 두 건 있어서 미공개 상태의 GenBank file을 genome submission portal에서 내려받은 후 bp_genbank2gff3.pl로 처리하여 GFF3 파일을 만들었다. 이를 다른 공개된 유전체 정보에서 전환한 GFF3와 한데 섞어서 Roary를 실행하였다.

그런데 core gene이 하나도 생성되지 않는 것이 아닌가. 처음에는 적절하지 않은 균주의 유전체 자료가 섞여 있어서 그러한 것이라 생각했었다. 나머지 것들과 너무나 거리가 먼 균주가 단 하나라도 섞여 있다면, 이들 전체를 아우르는 core gene이 나오지 않을 것임이 너무나 명백하다. 실제로 그러한지 확인을 하기 위하여 roary2svg.pl 스크립트로 gene_presence_absence.csv 파일을 처리하여 SVG(scalable vector graphics) 파일을 만든 다음 웹브라우저에서 열어보았다.


그게 아니었다. 비공개 GenBank file에서 생성한 GFF3가 제대로 쓰이지 못하고 있었다. 이 두 개의 유전체는 여러개의 contig로 구성된 것이다. 그게 Roary 작동에 문제를 초래하지는 않는다. GFF3 파일을 열어보니 contig 서열 ID가 있는 곳이 전부 '>noname'으로 표시된 것을 발견하였다. Contig가 몇 개이든 관계없이 전부 동일한 ID를 달고 있는 것이었다. 왜 그런 것일까? GenBank file의 앞부분으로 가 보았다. 하나의 contig에 대해서 다음의 정보가 나타나고(feature와 sequence를 포함), '//'를 구분자로 하여 다음 contig에 대해서도 LOCUS로 시작하는 새로운 정보가 이어진다.
LOCUS       <생략>
DEFINITION  Paenibacillus sp. F4, whole genome shotgun sequence.
ACCESSION
VERSION

DBLINK      BioProject: PRJNA429493
            BioSample: SAMN08348193
ACCESSION과 VERSION 난이 비어있다는 것이 원인일 것 같았다. 간단한 Perl 스크립트를 작성하여 이 부분을 LOCUS 값으로 채워 보았다. 비로소 정상적인 결과가 나왔다. 가까운 몇 개의 species에 속하는 균주의 데이터를 50개 가량 투입하여 분석한 것이라 core gene은 1070개 정도로 그 수가 조금 적다. 사소한 문제이지만 해결을 하니 속이 다 시원하다.

2018년 2월 3일 토요일

원핵생물의 subspecies는 어떻게 제안되는가

미생물의 학명에는 Bacillus subtils와 같이 속+종으로 구성된 학술적인 명칭 이외에도 var(variety), serovar, pathovar, subspecies(subsp.) 등의 정보를 뒤에 붙이고 있는 것이 많다. 이들 중에서 공식적인 명칭의 일부로 인정을 받는 것은 오직 subspecies 뿐이다. 그래서 subsp. 뒤에 오는 이름은 이탤릭체로 쓰는 것이 옳다. Subspecies는 동일한 종 내에서 유전적, 생리적, 또는 생화학적 특성에서 차이가 나는 세균의 하위 분류 단위를 의미한다.
Bacillus subtilis subsp. subtilisBacillus thuringiensis serovar israelensis
원래 ~var.는 끝에 약자임을 알리는 점(.)을 찍어야 하는데 잘 지켜지지 않는 것 같다. Variety는 변종을 의미한다. 원래 식물학, 원예학 등의 분야에서 이러한 용어의 정의가 만들어져 왔을 것으로 생각된다.  미국 아이오와 주립대 사이트에서 Cultivar versus variety라는 제목의 글을 발견하여 링크를 남긴다. 이 개념을 미생물 분류학 세계에 그대로 적용하는 것은 쉽지 않을 것이다.

미생물 분류학에 관련된 웹 자료는 이상하게 찾기가 힘들다. 오히려 수의(獸醫)미생물학과 관련한 곳에서(Nomenclature of bacteria) 읽기 쉬운 자료를 찾아내었다. 다른 종으로 보기에는 너무 유사하고, 동일 종으로 보기에는 너무 먼(가깝고도 먼?) 미생물에 대하여 학명 뒤에 subspecies epithet를 뒤에 붙이게 된다. Epithet란 별명 또는 형용어구를 의미한다. Genome sequence를 사용하여 두 균주가 동일 종인지 아닌지를 비교적 객관적으로, 그리고 재현성있게 판별할 수 있게 됨에 따라서 과거에는 어느 종의 subspecies로 여겨졌던 것이 완전히 다른 종으로 독립해 나가는 경우도 흔하다.
When you divide a species into several subspecies, the original species always gets the same subspecies epithet as the species epithet.
이는 무슨 말인가? 어떤 종을 몇 개의 subspecies로 나누게 될 경우, 원래의 종은 종소명을 species epithet로 갖게 된다. 예를 들어 설명하자. Bacillus subtilis에 subsp. spizizenii가 새로 제안되었다고 하자. 그러면 기존의 Bacillus subtilis는 Bacillus subtilis subsp. subtilis가 된다는 뜻이다.

여기서 잠시 종소명(種小名; specific name 또는 species epithet)에 대해 설명을 해야 되겠다. 린네가 제안한 이명법에서는 속명과 종명을 이탤릭체로 나란히 배열한다고 누구나 알고 있다. 그러면 Escherichia coli 중 Escherichi = genus, coli = species라고 생각하기 쉽다. 하지만 정확히 말하자면 coli는 종소명이다. 종(species)명은 언제는 Escherichia coli인 것이다.

FZB42(T)가 Bacillus amyloliquefaciens subsp. plantarum가 된 스토리

이틀 전 작성한 미생물 분류학 토막글 - (later) heterotypic synonym이란 무엇인가에서 한때 Bacillus amyloliquefaciens subsp. plantarum으로 불리던 FZB42 균주가 B. velezensis로 편입된 과정을 설명하였다. 지금은 더 이상 유효하지 않지만, FZB42 균주에게 plantarum이라는 subspecies epithet를 처음으로 제안하게 되었던 2011년 Rainer Borriss의 논문을 살펴보고자 한다.
Relationship of Bacillus amyloliquefaciens clades associated with strains DSM 7T and FZB42T: a proposal for Bacillus amyloliquefaciens subsp. amyloliquefaciens subsp. nov. and Bacillus amyloliquefaciens subsp. plantarum subsp. nov. based on complete genome sequence comparisons. [PubMed]
독일 브란덴부르크의 사탕무 토양(식물병이 창궐했던)에서 FZB13, FZB24, 그리고 FZB42라는 균주가 발견되어 학계에 보고된 것은 1998년이었다. 이 균주들이 식물의 성장을 촉진하는 것이 알려지면서 주목을 끌게 되었다. 이러한 성질은 plant-associated rhizobacteria가 나타내는 전형적인 특성 중 하나이지만 Bacillus subtilis DSM 10과 Bacillus amyloliquefaciens DSM 7(둘 다 type strain임)은 그렇지 않았다. 파지 감수성 실험에 의하면 B. amyloliquefaciens와 가깝다는 것이 확인되었지만 16S rRNA genes sequence, ribotype pattern, SfiI macrorestriction profile 등에서는 B. amyloliquefaciens DSM7(T)와 구별되는 특징을 보였다. 특히 FZB42는 전체 유전체의 8.5%이 다양한 secondary metabolite를 합성하는데 관여하고 있음이 밝혀졌는데, 이는 DSM7의 두 배를 능가하는 것이었다. DDH 실험에서도 63.7-72.1%의 DNA-DNA related value를 보여서 두 균주를 완벽하게 다른 종으로 구별하거나 한 종으로 간주하기에도 애매한 상황이다(2017년 논문에서 DSM7과 FZB42의 ANI 값은 94.61%로 발표함). 그래서 새로운 B. amyloliquefaciens subsp. plantarum을 제안하고 FZB42를 type strain에 위치시킨 것이다. 
"우리가 아는한 Bacillus 속 안에서 전통적인 미생물 분류법과 whole genome comparison을 이용하여 새로운 subspecies를 제시한 것은 이번이 처음이다."

2018년 2월 1일 목요일

미생물 분류학 토막글 - (later) heterotypic synonym이란 무엇인가

Bacillus amyloliquefaciens subsp. plantarum FZB42라는 균주가 있었다. 왜 '있었다'라는 과거형으로 표현하였는가? 지금은 전세계 어디에도 이 균주를 보관하는 culture collection이 없거나, 지구상에서 멸종을 했다는 뜻인가? 그런 문제는 다행히도 아니고, 공식 명칭이 바뀐 것이다.

분류학은 생명과학의 여러 분야 중에서 매우 독특한 학문이다. 생물체의 특성(해부학적, 생화학적, 유전체학적 등)을 주의깊게 관찰하고 체계화하는 것은 일반적인 과학의 단계와 같다. 분류학도 과학이니 가설을 입증하기 위한 실험을 하게 된다. 분류학에서 가설이란 이 생명체가 새로운 종인지, 또는 기존에 알려진 종과 얼마나 다른지에 대한 주장일 것이고, 따라서 분석 실험을 많이 진행하게 된다. 분류학이 다른 과학과 크게 다른 것은 매우 엄격한 규칙에 의존한다는 것이다. 물론 이 규칙은 새로운 지식이 발전하면서 계속 변한다.

나는 분자생물학이 국내에 맹렬하게 도입되던 시기에 공부를 하였기에 분류학에 대한 체계적인 교육을 받지 못하였다. 분자생물학·생화학은 모든 생명체를 아우르는 중심 원리에 집착한다. 반면 분류학은 거시적인 관점을 중요시하고, 동물·식물·미생물(아, 너무나 초보적이고 편견이 가득한 생명체의 분류 기준이여!)에 대해서 적용할 수 있는 방법이 천차만별이다. 현미경으로 겨우 보이는 세균에 손발이나 얼굴, 머리카락이 있는 것이 아니기 때문이다.

요즘 환경보호 이슈뿐만 아니라 각 국가가 고유하게 보유한 생물종의 전략자원화 추세에 따라 다시 분류학이 중요한 분야로 떠오르는 것 같다. 예전에는 정말 돈이 안되는 연구였다고 해도 과언이 아니니까 말이다. 어느 지역에 자생하는 어떤 식물에 혁신적인 신약 후보 물질이 있다고 하면, 그 식물을 제대로 구별해낼 수 있는 사람의 몸값이 올라갈 것이다. 학문적인 중요성이 아니라 경제적인 판단 기준에 의해서만 인정을 받으면 좀 아쉽긴 하다. 분류학 자체는 당연히 중요한 학문이다.

오늘 아침, Paenibacillus 종의 유전체 비교 분석을 하기 위한 바깥쪽 기준으로서 Bacillus amyloliquefaciens subsp. plantarum FZB42을 선택하여 그 서열 정보를 NCBI에서 찾아보았다. 이 균주는 식물체에 붙어서 사는(plant-associated bacterium을 이렇게밖에 표현하지 못하겠다) 그람 양성 세균으로서 식물의 성장을 촉진하며 항생물질로 활용할 수 있는 다양한 이차대사산물을 만들어낸다. 이미 다양한 미생물 비료 제품으로 만들어져 시판되고 있으며, 전체 유전체가 완전히 해독되어 2007년에 Nature지에 보고가 되었었다(논문 링크). 이 연구를 주도한 사람은 독일 훔볼트 대학에 근무하는 Rainer Borriss(ResearchGate)로서 현재 Bacillus와 Paenibacillus의 분류 및 comparative genomics 연구에 집중하고 있다. 연구 관심사가 나와 매우 유사하다. Bacillus는 16S rRNA gene seuqencing가 워낙 비슷해서 전통적인 방법으로 종을 구분하기가 까다롭기로 유명하다. 그래서 논문들을 추적하다보면 수시로 균주들이 이 종에서 저 종으로 옮겨다니게 만드는 제안을 많이 접하게 된다.

그런데 NCBI genome page에서 Bacillus amyloliquefaciens subsp. plantarum FZB42가 보이지 않는 것이었다. Bacillus amyloliquefaciens에 속하는 균주의 유전체 목록에 FZB42가 나타나지 않았다. 이게 도대체 어찌된 일인가? 몇년 전에 분명히 다운로드하여 잘 활용했었는데? 검색을 거듭한 결과 Bacillus velezensis로 명칭이 바뀐 것을 알 수 있었다. 그래서 이를 대신하여 B. amyloliquefaciens의 type strain인 DSM 7(free-living soil bacterium)의 유전체를 내려받아 사용하였다.

균주의 재분류 혹은 신종 보고 등에 관한 논문은 주로 International Journal of Systematics and Evolutionary Microbiology(줄여서 IJSEM이라 부른다)에 실린다. 여기에 실려서 공인을 받아야 비로소 어떤 균주의 명칭을 '공식적'으로 쓸 수 있고, 예전에 명명한 이름을 바꾸는 것도 이를 통해야 한다. 박테리아의 공인된 명칭 목록(approved list of prokaryotes)는 독일의 DSMZLSPN bacterio.net에서 확인 가능하다. 공인되기 전까지는 "Bacillus gaemokensis"처럼 따옴표로 둘러쳐야 한다. 정확하게 표기하려면 이탤릭으로 표시해야 함은 물론이다. 여담이지만 요즘은 균주, 유전자, 유전자형, 표현형, 학명 등을 쓸 때 예전만큼 엄격하게 규칙을 지키지 못하는 것 같다.

FZB42가 다른 종 명칭으로 갈아타게 만든 근거가 된 논문은 2016년도에 IJSEM에 실린 다음의 논문이었다.

[1] Bacillus velezensis is not a later heterotypic synonym of Bacillus amyloliquefaciens; Bacillus methylotrophicus, Bacillus amyloliquefaciens subsp. plantarum and ‘Bacillus oryzicola’ are later heterotypic synonyms of Bacillus velezensis based on phylogenomics. (논문 링크)

분류학에 익숙하지 않은 사람(나를 포함하여)은 무슨 뜻인지 알기 어려운 제목이다. Bacillus velezensis는 B. amyloliquefaciens의 later heterotypic synomym이 아니라고 한다. 그러면 그 과거에는 later heterotypic synonym이었던 적이 있었나보다. 정말 그렇다!

[2] Bacillus velezensis is a later heterotypic synonym of Bacillus amyloliquefaciens. (2008년 논문 링크)

그리고 이 복잡한 관계를 정립하는 것으로 보이는 최근 논문이 2017년도에 나왔다. 역시 Borriss가 참여하였다.

[3] Bacillus amyloliquefaciens, Bacillus velezensis, and Bacillus siamensis Form an “Operational Group B. amyloliquefaciens” within the B. subtilis Species Complex. (논문 링크; PubMed가 아니라 PMC가 좋은 점은 논문 원문을 읽을 수 있다는 것이다)

자, 그러면 Bacillus를 둘러싼 복잡한 내력은 조금 뒤에 알아보는 것으로 하고, (later) heterotypic synonym이라는 가슴에 잘 와닿지 않는 용어를 공부해 보기로 하였다. 서울대학교 식물분류학연구실 장진석 교수님의 웹페이지에는 식물 분류와 관계된 자료를 많이 서비스하고 있다. 여기에서 용어 설명 페이지(링크)를 찾아가 보았다.

  • synonym 이명. [미정의] 같은 분류군에 적용되는 복수의 학명 중 1개(이형이명, 동형이명 참고)
  • homotypic synonym 동형이명(명명법상 이명). 학명이 다르지만 같은 기준표본에 근거한 학명(제14.4조); 국제동물명명규약과 국제세균명명규약(세균규약)에서는 각각 객관적 이명(objective synonym)으로 불림(제14.1조 각주)
  • heterotypic synonym (taxonomic synonym) 이형이명(분류학적 이명) 같은 분류군을 나타내는 학명과는 서로 다른 기준표본에 근거한 학명(제14.4조). 국제동물명명규약과 국제세균명명규약(세균규약)에서는 각각 주관적 이명(subjective synonym)으로 불림(제14.4조 각주)
실제 조항은 조류, 균류와 식물에 대한 국제명명규약(멜버른 규약 2012 번역 링크)을 참조하면 된다. 박테리아에 관한 규약은 어디에 있을까? NCBI Bookshelf에 있다(링크).

국문으로 된 자료를 접했다고 해서 '이형이명'이 쉽게 이해되는가? 영 그렇지 못하다. 자연계에서 순수하게 분리한 박테리아 표본이 있다고 가정하자. 이를 서로 다른 연구자가 별도의 과정을 통해서 나름대로 이름을 붙인 경우 그 중 하나가 동형이명이 되는 것으로 생각된다. 그렇다면 이명이형은? 그리고 'later' 이형이명은? 


가장 최근 논문인 [3]의 서론 부분을 읽어보았다. 2004년 Reva 등이 식물 혹은 토양에서 일곱가지의 B. subtilis group 균주를 분리하였다. 이들은 근권에서 잘 정착하여 자라면서(rhizosphere colonization) B. amyloliquefaciens의 표준 균주인 DSM7과는 다른 특성이 있었다. 여기에는 GB03과 FZB42가 있었다. GB03은 2014년 우리 연구그룹에서 유전체를 해독하여 짧은 보고를 한 적이 있다(링크). 새로운 종인 Bacillus velezensis는 1년 뒤인 2005년 신종으로 보고되었다(링크).

즉 FZB42는 B. amyloliquefaciens에 속하는 하나의 균주로 학계에 보고되었고 1년 뒤에 별도의 과정을 거쳐서 B. velezensis라는 종이 보고되었다는 것이다. B. amyloliquefaciens species가 제안된 것은 1943년인데, 실제로 공인된 것은 훨씬 뒤인 1987년이었다. 그런데 추후에 두 종의 type strain을 가지고 실험을 해 보니 "high DNA-DNA similarity and lack of phenotypic distinctive characteristics"에 의해서 B. velezensis를 B. amyloliquefaciens의 "later heterotypic synonym"이 되어야 한다고 주장한 것이 논문 [2]이다. B. velezensis BCRC 17467(T)와 B. amyloliquefaciens BCRC 11601(T) DNA 유사도(DDH 실험 결과)는 74%여서 다른 종으로 보기 어려웠다는 중요한 근거 중 하나가 되었다.

만약 B. amyloliquefaciens의 기준 표본(type strain과 동의어로 보아도 무방할 것이다)을 가지고 연구하던 사람이 이를 B. velezensis라고 명명했다면 아마도 homotypic synonym이라고 부르게 되지 않을까? 나는 최소한 이렇게 이해하였다.

그러면 다음의 어느 표현이 맞을까? 둘 다 맞을까, 아니면 어느 하나만 맞을까?

  1. Bacillus velezensis (syn. Bacillus amyloliquefaciens subsp. plantarum) (syn. Bacillus methylotrophicus) (syn. Bacillus methylotrophicus subsp. plantarum) (syn. Bacillus oryzae) 
  2. Bacillus amyloliquefaciens (syn. Bacillus velezensis) 
NCBI taxonomy 사이트에는 오직 (1)의 방식으로만 표기되어 있음을 알 수 있다(링크). 만약 (1)과 같이 표기한다면 먼저 공인된 종 명칭인 B. amyloliquefaciens가 매우 섭섭할 것이다. 그러나 아직 두 표기 형태 중 어느 것이 올바른지는 확신이 서지 않는다. 심지어 B. velezensis는 syn. B. oryzicola, 그리고 syn. B. methlotrophicus이기도 하다. 이는 논문[1]의 내용 중 하나이기도 하다.

왜 논문[1]에서는 B. velezensis가 B. amyloliquefaciens의 later heterotypic synonym이 아니라고 했을까? 저자들은 DDH 실험이 아니라 genome sequencing을 그 건거로 삼았다. B. amyloliquefaciens DSM7T에 대해서는 dDDH(digital DDH, 즉 genome sequence에 의한 계산; Genome-toGenome Distance Calculator 링크)가 55.5%에 불과한 반면, B. amyloliquefaciens (subsp. plantarum) FZB42T에 대해서는 85.8%였다. Sequencing에 의한 계산이 실제 DDH 실험을 이긴 셈이다.

어차피 B. amyloliquefaciens DSM7과 B. amyloliquefaciens (subsp. plantarum) FZB42T의 dDDH 계산수치는 56.2% 수준이었다. 이 논문으로 말미암아 B. velezensis는 B. amyloliquefaciens의 synonym을 벗어났고, 대신 B. amyloliquefaciens subsp. plantarum은 Bacillus velezensis의 later heterotypic synonym이 된 것이다. 또한 B. methylotrophicus와 'B. oryzicola' 역시 B. velezensis의 later heterotypic synonym이 되었다.

이제 논문[3]을 음미하면서 오랜만에 길게 쓴 글을 마무리하고자 한다. DSM7과 유사한 66개의 genome을 내려받아 비교한 결과(rpoB 유전자와 core genome sequence 사용) 이들은 (1) B. amyloliquefaciens, (2) B. siamensis, (3) a conspecific group(B. amyloliquefaciens subsp. platarum, B. velezensis, 그리고 B. methylotrophicus의 type strain을 전부 포함)의 세 개 그룹으로 묶을 수 있었다(근거로는 TETRA와 AAI). Monophyletic clade에 속하는 이 세 개의 그룹은 너무나 유사하여 종의 약간 상위에 존재하는 "operational group B. amyloliquefaciens"를 제안하였다. 비록 이들이 요즘 점점 중요성이 강조되고 있는 분석 기법인 ANI나 dDDH로는 한 종의 범주에 넣기에는 다소 곤란하다고 해도 말이다. 또한 이그룹은 식물과 더불어 살아가는 라이프스타일을 계쏙 발전시켜나가고 있는 중이며, 그러한 과정이 유전체에 남아있는 것이다.

2018년 1월 31일 수요일

[하루에 한 R] PhyloSift 결과 종합하기(merge 함수로 tab-delimited CSV 파일의 일괄 입력 및 병합 등)

Species boundary를 살짝 넘어가는 박테리아 유전체 50개 정도를 가지고 whole genome phylogenetic tree를 그려보고자 하였다. Harvest suite로 분석을 하면 MUMi distance cutoff에 미치지 못하는(default = 0.1) 몇 개가 빠져나가고, Roary로는 아예 core gene이 잡히질 않았다(뭘 잘못했나?). 상위 genus에서 폭넓게 분포하는 어떤 마커들의 분포(LS-BSR로 자료 생성)를 시각화하는 것이 목표라서 tree는 아주 정확할 필요가 없다. 다만 모든 샘플을 아우를 정도로만 윤곽을 잡아주면 된다.

16S rRNA gene sequence로는 해상도가 좋은 트리가 나오기 어려울 것이 자명하다. 그래서 택한 도구는 PhyloSift였다. 이 도구는 원래 shotgun metagenomic read로부터 phylogenetic distribution을 알아보기 위하여 개발된 것이지만, 단일 유전체 어셈블리로부터 universal marker(목록; 원래 40개였으나 실제로 쓰이는 것은 37개)를 찾아낼 수 있으므로 미생물 유전체 조립의 완성도를 평가하는데 쓰이기도 한다. 더욱 좋은 것은 alignment까지 된 상태의 서열이 생성되므로 이를 잘 이용하면 tree를 쉽게 그릴 수 있다. PhyloSift가 metagenome sample로부터 만들어내는 tree와는 다른 것이니 주의가 필요하다.

Bacillus_subtilis_168.fa라는 유전체 파일을 이용하여 PhyloSift를 실행했다고 가정하자(search 및 align mode). 그러면 PS_temp/Bacillus_subtilis_168.fa/ 아래에 alignDir 및 blastDir가 생기고, 그 하위에 각 마커에 대한 파일이 쌓인다. marker_summary.txt는 다음과 같이 생겼다.


각 유전체 서열마다 별도로 작성된 이 파일을 병합하여 정리하면 모든 마커가 전부 확인되었는지의 여부를 쉽게 알 수 있을 것이다. 오늘의 고행은 여기서부터 시작되었다. Perl을 사용하면 어떻게 해서든 결과를 얻을 수 있었겠지만, 왠지 R을 쓰고 싶었다. R을 이용하여 코드를 짜 놓으면 재사용성이 높을 것 같았기 때문이었다.

marker_summary.text file의 정리


각 유전체 서열에 대해서 별도의 디렉토리에 나뉜 marker_summary.txt 파일을 한 곳으로 모으되 균주명이 파일의 제목이 되게 바꾸었다. Shell script를 써서 다음과 같이 구현하였다(짧은 코드인데도 공부할 것이 많아서 2019년 3월 8일 설명과 함께 업데이트). 이 스크립트에서는 완벽하지 않다. 털어내야 할 FASTA 파일의 확장자(여기에서는 fasta)를 실제 상황에 맞게 바꾸어야 하고, 출력파일을 저장할 OUTDIR 환경변수를 미리 선언하고 해당 디렉토리도 만들어 두어야 한다. 그리고 for block의 처음에 나오는 <<<도 이해하기가 조금 어렵다. f 변수에는 aaa/bbb/ccc...라는 형태의 값이 저장되는데, 이를 슬래쉬로 구분한 다음 두번째 것인 bbb를 택하겠다는 의미이다. 만약 <<<쓰지 않으면 cut은 f 변수가 가리키는 파일의 내용을 대상으로 작업을 한다. 이 명령행은 FILE=$(echo $f | cut -d'/' -f2)로 바꾸 실행해도 완전히 동일한 결과를 낸다. 16s_reps_bac과 DNGNGWU* 마커에 대한 결과를 같이 출력하고 싶다면 egrep "16s_reps_bac|DNGNGWU" 명령을 써라.
FILES=`find $1 -name marker_summary.txt -type f | grep blastDir`
for f in $FILES
do
        FILE=$(cut -d'/' -f2 <<<$f)
        FILE=${FILE%.fasta}
#        cp $f ${OUTDIR}/${FILE}
        grep DNGNGWU $f > ${OUTDIR}/${FILE}
done
find로 찾아낸 파일명(경로 포함)은 다음과 같다. 빨강색 문자열만을 뽑아내는 것이 포인트이다. cut 명령을 사용하여 슬래쉬(/)를 기준으로 분리한 두번째 문자열을 취한 다음 뒷부분의 .fasta를 제거해야 한다.

PS_temp/Paenibacillus_sp._UNCCL52_GCF_000686825.1.fasta/blastDir/marker_summary.txt

오늘의 시도에서는 파일 내부를 grep으로 뒤져서 DNGNGWU...로 시작하는 표준 마커에 대한 결과만을 모았다. match가 없는 마커의 경우 0이 표시되는 것이 아니라 아예 summary file에 존재하지 않는다. 이것이 나중에 R에서 병합 작업을 하는 것을 약간 어렵게 하였다. 만약 모든 마커에 대한 결과가 다 있었다면, 마커의 순서가 동일한지만 확인한 뒤에 cbind()를 하는 것으로 쉽게 끝났을 것이다.

[R] for loop를 이용한 텍스트 파일의 일괄 입력 및 outer join


이러한 용도의 R 패키지가 이미 개발된 것이 있을지도 모르지만 되도록 R의 기본 기능을 충실하에 이용하여 구현하고자 하였다. 먼저 정리된 파일이 모여있는 디렉토리로 이동하여 R을 실행하였다. 다른 곳에서 R을 실행한 뒤 setwd() 함수를 실행하여 파일 위치를 지정하여도 좋다. 파일의 일괄 입력에서는 Import multiple files to R의 예제를 거의 그대로 따랐다.
list.filenames = list.files()
# 특정 확장자를 지닌 파일만 읽어들이려면 list.files(pattern="*.txt")
# create an empty list
list.data=list()
length(list.filenames)
# create a loop to read in your data
for (i in 1:length(list.filenames))
{
# 이중 대괄호 안에 인덱스를 넣어야 데이터프레임을 부르는 것이 된다.
list.data[[i]] = read.table(list.filenames[i],sep=",",header=F)
colnames(list.data[[i]]) = c("markers", list.filenames[i])
}
# 다음의 두 명령은 결과가 약간 다르다.
# 궁금하면 class() 함수를 써 보라.
list.data[1]
list.data[[1]]
# add the names of your data to the list
names(list.data) = list.filenames
# name으로 사용된 파일명에서 확장자를 제거하려면 방법을 따로 알아보아라.
다음으로는 데이터 프레임을 join해야 한다. 기준이 되는 컬럼은 markers이다. outer join을 하였으므로 값이 없는 경우에는 NA가 채워지고, 이는 나중에 0으로 채우면 된다. R에서의 join 혹은 merge에 대한 개념은 여기를 참고하였다. 국문 페이지로는 "[R Friend] R 데이터프레임 결합: rbind(), cbind(), merge()"가 매우 유용하다. merge 함수는 두 개의 데이터프레임을 대상으로 한다. 만약 인덱스로 접근 가능한 데이터프레임이 있는데 총 50개가 있다면 이를 어떻게 합치겠는가? 1번과 2번을 merge하여 임시 데이터프레임을 거쳐서 2번에 저장한다(한꺼번에 저장하는 것이 더 간단할지도 모른다). 그 다음 인덱스를 1 증기시키고 2번과 3번을 merge하여 3번에 저장한다. 이를 인덱스=(전체 데이터프레임 수 - 1)이 될 때까지 for loop를 돌린 다음 종료하면 된다.

# full outer join
for (i in 1:(length(list.filenames)-1))
{
df = merge(x=list.data[[i]],y=list.data[[i+1]],by="markers",all=TRUE)
list.data[[i+1]] = df
}
df = list.data[[length(list.filenames)]]
df[is.na(df)] = 0
rownames(df) = df[,1]
df = df[,-1]
colnames(df) = list.filenames
# 최종 확인
dim(df)
View(df)

View(df)의 결과 화면(일부).
첫번째 컬럼(markers)를 row name으로 전환하여 새로운 데이터프레임 df.2에 저장하자.
> df.2 = data.frame(df[,-1], row.names=df[,1])
> dim(df.2)  # 컬럼 수가 하나 줄어들었다.

결과 점검하기


단 하나라도 마커가 발견되지 않은 유전체는 어느 것일까? 또 어떤 마커가 발견되지 않았을까? 마커의 총 수가 37개라고 해서 컬럼의 합을 가지고 판단의 기준을 삼을 수는 없다. 왜냐하면 주어진 마커가 2회 발견되는 경우도 있을 수 있기 때문이다. apply()함수를 사용하면 이를 매우 쉽게 찾아낼 수 있다. 

각 row에 대해서 어느 하나라도 0이 존재하는지를 알아보자. 빨강색으로 표시한 마커가 TRUE를 나타내고 있다. 바로 저 row name만을 반환할 수 없을까?


> apply(df.2 == 0, 1, any) # row에 대하여
> apply(df.2 == 0, 2, any) # column에 대하여
위 명령어를 입력한 다음 눈으로 찾을 수는 없는 노릇이다. 정답은 다음과 같다.
> rownames(df.2)[apply(df.2 == 0, 1, any)]
[1] "DNGNGWU00040"
> colnames(df.2)[apply(df.2 == 0, 2, any)]
[1] "Paenibacillus_polymyxa_TD94_GCF_000520775.1"
오늘 하루가 다 지나갔다! 데이터프레임에서 특정 조건을 만족하는 서브셋을 추출하는 것도 중요하지만, 그 조건에 만족하는 column 혹은 row의 이름을 알아내는 것도 중요하다. 

2018년 1월 30일 화요일

또 진공관 앰프?

<성장 없는 번영>이라는 책을 다 읽고 나니 '소비'라는 행위를 뒤돌아보게 되었다. 우리 주위에는 '제발 나를 쇼핑해 달라고' 외치는 갖은 유혹이 얼마나 많은가? 이제는 필요에 의해서 구입을 하는 것이 아니라 이것을 가짐으로 인해서 내 삶이 과시적으로 업그레이드된다는 욕구를 채우기 위하여 구입을 한다. 그 과정 중에 생태계는 자꾸 무너진다. 별도의 독서 기록에서 다루겠지만, 모든 지구 상의 나라가 현재의 미국인처럼 소비할 수 있는 날은 오지 않을 것이다. 이에 대한 대안을 찾자는 것이 이 책의 주제이다.

그럼에도 불구하고 마우스 포인터는 이미 온라인 쇼핑 사이트를 뒤적거리고 있다. '무엇을 살까?'는 고통스런 고민을 하면서 계속 서치를 하는데 들어가는 비용은 어쩌면 일단 무엇인가를 지르는 비용보다 더 클지도 모른다. 그러니 너무 오래 고민을 반복하지 말고 일단 구매를 하자. 이러한 경제적인 자기 합리화를 해 가면서 뭔가를 샀다고 치자.

또 얼마간의 시간이 지나면 커피 찌꺼기에서 곰팡이가 피어나듯(요즘 집에서 늘 겪는 문제) 새로운 욕구는 또 생겨날 것이다.

그래... 웹 서칭이 무슨 큰 잘못이랴. 근무 중에 지나치게 다른 사이트를 들락거리는 것은 물론 옳지 않다... 일단 조사만 해 두자. 중국제 저가 앰프(자작용 키트를 포함하여)를 둘러보면 어떤 유행과 같은 것이 보인다. 일년 전에 잔뜩 쏟아져나온 물건은 어느덧 사라지고 또 새로운 트렌드의 물건이 목록을 차지한다. Class D 앰프가 그랬고, JLH 1969년 class A 앰프가 그러했다.

진공관만으로 이루어진 채널당 3 W 남짓의 그저 그런 앰프는 내가 보유한 저능률 스피커 시스템을 만족스럽게 울리기 어려울 것이다. 그렇다고 하여 출력이 약간 큰 저가의 앰프를 고르면 집에서 2014년 이래 아주 잘 쓰고 있는 12DT8 + 14GW8(=PCL86) 앰프(링크)보다 결코 나을 수가 없을 것이다. 그래서 진공관 프리앰프 + IC를 이용한 이용한 hybrid amplifier를 찾아보았다. 하이브리드 앰프라 하면 요즘 Nobsound의 MS-10D가 단연 인기이다. '진공관은 단지 장식용이다! 아니다, 진공관을 빼면 소리가 안난다! 트랜스 커버를 벗기면 아무것도 없다!' 등 여러 논란에 휘말린 앰프인데 일단 모양새는 그럴싸하다. 게다가 아날로그식 레벨 메터까지 있으니 보는 재미는 좋을 것이다. 처남이 얼마전에 주문을 했다는데 배송이 완료되었는지 궁금하다.



아래에 소개한 앰프 키트는 전부 6J1 진공관(교체관 정보)과 LM1876T를 사용한 것이다. 위의 것은 헤드폰 단자도 갖추고 있다.

인터넷에서 '무단으로' 그림을 퍼 온 뒤에 글씨를 달았다. 그림의 URL을 연결하는 방식으로 블로그에 소개하고 싶었지만(이것도 무단사용에 해당할까?) ebay 등의 상품 주소는 판매기간이 종료되면 더 이상 유효하지 않은 것 같다. 굵은 글씨로 표시한 키워드를 이용하여 검색하면 어떻게든 찾을 수 있다.

아아, 집에도 앰프가 세 개, 사무실에도 세 개나 있는데 뭘 또 만들고 싶어서 이런단 말인가?

2018년 1월 27일 토요일

독서 기록 - <나를 믿어주는 한 사람의 힘> 외 네 권

이번에는 대출 기한을 일주일 연장했음에도 불구하고 대출한 책을 다 읽지 못하였다. <처음 읽는 일본사(전국역사교사모음 지음)>, <불평등의 역사(발터 샤이델 저|조미현 역)>, 그리고 <스트레스는 어떻게 삶을 이롭게 하는가(우르스 빌만 저|장혜경 역)>는 다음 기회에 꼭 읽도록 한다.

사진 퀄리티가 영 좋지 못하다. 거실 바닥에 놓고 찍으면 조명이 늘 문제다

나를 믿어주는 한 사람의 힘

  • 부제: 문화치유전문가 박상미의 인생특강
  • 박상미 지음(더공감 마음학교 링크)
나를 믿고 지지해주는 한 사람의 힘으로 인하여 풍성한 삶을 이어나가는 10 인을 심층 인터뷰하고 쓴 책. '시련을 극복한~', '성공한~' 이라는 수식어를 붙이는 것은 적당하지 않다. 이 척박한 세상에 향기를 드리우며 다른 사람에게 좋은 영향력을 미치는 풍성한 삶을 누리는 사람이라고 하면 딱 좋을 것이다. 참여한 사람은 김혜자, 박동규, 이병복, 표재순, 신경림, 인순이, 황현산, 조벽, 김현영, 그리고 섀넌 두나 화이트. 맨 마지막은 저자 박상미가 자신을 인터뷰하는 형식으로 마무리하였다. 문학평론가 황현산이 현명한 어른이 되는 길에 대하여 말한 것을 인용해 보았다.
나이가 들수록 듣는 연습을 해야 하고, 토론을 해야 합니다... 배우기를 멈추지 말고 참신하게 생각하도록 노력해야 합니다. 가장 중요한 게 책을 읽는 것입니다... 나보다 어린 애들이 쓴 글이라고 생각하는 순간 문제가 시작됩니다... 나이의 향기를 풍길 수 없는 사람들이 나이의 권력을 탐하는 것 같다(박상미) 
교수법 전문가 조벽 교수는 이렇게 말했다
(부부가 이혼에 이르는 가장 큰 문제는) 바로 '말'입니다... 멀어지는 말이 상대의 마음을 아프게 해서 마음 문을 아예 닫아버리게 하죠... 그 말에 네 가지 독이 존재해요. '비난, 경멸, 방어, 담 쌓기'라는 독이에요... 이 독을 마음에 품는 순간 얼굴에 다 나타난대요.

어떻게 바꿀 것인가

  • 부제: 비정상 정치의 정상화를 위한 첫 질문
  • 강원택(서울대 정치외교학부 교수) 지음
긴 고난의 터널을 뚫고 대한민국 민중은 드디어 5년 단임 대통령제를 핵심으로 하는 1987년 체제를 탄생시켰다. 그러나 이제는 강력한 리더십을 지닌 슈퍼맨 대통령에게 강력한 힘을 위임하는 것으로 더 이상 국가의 변혁이 어려움을 설파하고, 제도적 대안으로서 내각제를 주장한 책이다,

미래중독자(The Invention of Tomorrow)

  • 부제: 오늘을 버리고 내일만 사는 별종, 사피엔스
  • 다니엘 S. 밀로 지음
  • 양영란 옮김
"내일 보자!" - 인류가 떠올린 가장 혁명적인 문장.
인류는 '과잉'의 시대에 살고 있다. 저자 역시 인간 사회에 만연한 과잉에 대한 '과잉된 강박'을 가지고 있다고 한다. 정속 주행은 퇴보로 여기고 '가속도'를 추구하는 것이 인간의 본성이 아닐까? 왜 5만 8천년 전 호모 사피엔스는 별안간 아프리카를 떠났을까? 이것은 바로 인간이 발명한 '내일' 때문이라는 것이다. 유발 하라리가 주장했듯이 추상적 개념·상상력이 인류의 협동을 이끌어내어 현재의 문명을 만들어냈다는 것과도 일맥상통한다. 저자의 말에 의하면 이 '빌어먹고도' 찬란한 내일에 중독된 것이 인류이다. 내일을 위해 기꺼이 오늘을 포기하는 유일한 생명체가 인간이다. 기대는 곧 불안·와 동전의 양면과 같은 관계가 아닐까?

현대미술은 처음인데요

  • 부제: 큐레이터가 들려주는 친절한 미술 이야기
  • 안휘경, 제시카 체라시 지음
  • 조경실 지음
예술가의 똥(Artist's Shit; Merda d'artista; 피에로 만초니, 1961)

어떠한 작품이 현대미술인지 아닌지를 구분하는 기준은 단순한 시간이 아니라 작품의 기준이 무엇이냐에 따라 달려 있다. 현재라는 순간을 해석하는 작업은 좀처롬 쉽지 않고, 현대라는 사회를 이해하는 일 역시 혁신적인 기술과 새로운 접근방식을 요구한다(30 쪽). 

화장실에서 흔히 볼 수 있는 기성품 변기를 가져다가 작가의 사인만을 추가한 마르셀 뒤샹의 샘('Fountain', 위키피디아)을 현대미술의 시작으로 보는 견해도 있다. 하지만 이는 '세계 3대 기타리스트'라는 식으로 지나친 단순화라고 생각한다. 현대미술의 이단아 마르셀 뒤샹의 "샘"이야기(링크)를 읽어볼 것을 권한다.

작품은 사라지고 알프레드 스티글리츠가 찍은 사진과 복제품만 남았다.

때로는 어떤 작품은 미적 가치보다 아이디어, 정치적 관심, 감정의 자극으로 주목을 받을 수도 있다(31쪽).
작가 혼자서 공들여 만들여야만 작품이 되는 것이 아니다. 현대 미술은 이미 여러 스태프에 의해 구현되는 것이 일반적이며, 심지어 개념만으로도 작품이 되며, 퍼포먼스 아트도 이미 흔하다. 다음에 현대미술관을 방문하게 되면 좀 더 열린 마음으로 편안하게 감상을 할 수 있으리라 기대해 본다. 현대미술은 항상 변화하고 있는 현대를 보여주는 것이다. 중요한 것은 작품이 아니라 작가이다.

99%를 위한 경제학

  • 부제: 낮은 곳을 향하는 주류 경제학 이야기
  • 저자: 김재수(미국 인디애나퍼듀대학 경제학과 교수) 허핑턴포스트의 김재수 블로그 글목록 링크
보이지 않는 손, 비교우위, 기회비용, 비용-편익 분석, 경제적 인센티브에 반응하여 행동하는 인간... 우리가 알고 있는 주류 경제학은 단지 상위 1%를 위한 것이 아니었을까? 시장은 '기업'과 동의어가 아니다. 시장이 올바르게 작동하려면 모든 경제적 판단을 위한 정보가 전부 투명하게 공개된 상태에서 공정한 경쟁이 이루어져야 한다. 
"우리는 인센티브를 믿사오니, 공공부문과 정부규제에 고난을 받으사 복지 예산에 못박혀 죽으시고, 민영화와 규제완화를 통해 다시 살아나시며, 산 기업과 죽은 기업을 심판하러 오시리라. 성과급을 믿사오며 부자감세와 기업 프렌들리 정책의 효과가 서로 교통하는 것과 영원히 경제성장이 이루어지는 것을 믿사옵나이다. 아멘."(129쪽)
자유시장경제에 대한 옹호·세계화의 파도에 휩쓸려 경제적 약자인 을('乙')의 위치는 날이 갈수록 열악해져만 간다. 효율성을 추구하기 위해 어쩔 수 없다는 논리와 함께 이들은 계속 나락으로 떨어져만 간다. 뭔가 잘못되었다고, 권력과 결탁한 거대 기업의 편이 되어서는 안된다고 주장하고 싶지만, 아직도 주변에는 "박정희 대통령 덕분에 우리가 이정도 먹고 사는거야"라고 소리치는 사람들이 많다.
선택이 야기하는 기회비용을 본다는 것은 경제학도가 반드시 지녀야 할 반골 정신을 의미합니다. 권력과 권위, 또는 관습에 의해 당연하게 여겨지는 것들에게도 잃고 있는 것을 찾아내야 합니다.(11-12쪽)
이번에 읽은 다섯 권의 책 중에서 가장 울림이 큰 책이었다. 

소화전함이 수령한 택배 배달 물건

집에 사람이 없을 때 아주 가끔 소화전함에 물건을 넣어두었다는 문자 메시지를 받을 때가 있다. 바쁜 택배 기사들의 입장을 생각하면 이해하지 못할 일은 아니지만 물건이 없어진 경우 책임 소재의 문제가 따른다. 만약 나쁜 마음을 먹은 사람이 아파트에 들어와서 고층부터 소화전함을 훑으면서 한번씩 다 열어본다면 수십 대 일의 가능성으로 주인을 기다리는 물건을 훔쳐갈 수도 있는 것 아니겠는가? 수령 예정자가 소화전함에 넣어 두는 것을 동의했다고 해서 배송 기사에게 아무런 책임이 없을까?

배송기사: OO택배입니다. OOO님 지금 댁에 안계신가요?
수령인: 네, 지금 집에 없으니 소화전함에 넣어주세요.
배송기사: 그건 규정 상 안됩니다. 부재중이시면 경비실에 맡기거나 내일 다시 오겠습니다.
수령인: 에이, 뭐 별 일 있겠어요. 없어지면 제가 책임질테니 그냥 소화전함에 넣어주세요.

만약 이런 대화가 오갔다면 배송기사의 책임이 없을지도 모르겠다. 법을 다루는 사람들의 해석은 어떨지 모르겠지만 말이다.

오늘 아침, 오후 5~7시 사이에 주문한 케이블이 배송될 예정이라는 메시지를 받았다. 잠시 외출하였다가 집에 돌아온 것은 3시 정도. 침대에서 졸다가 5시쯤 문자 메시지를 받았는데, "소화전"께서 인수를 했다는 내용이었다. 다른 전화 연락은 일절 받은 것이 없었고 심지어 초인종도 울리지 않았다.


우리 아파트는 지은지 오래되어서 인터폰이 고장난 곳이 대다수이다. 인터폰을 걸어보다가 받지 않았다면 최소한 집 전화나 휴대폰으로 연락을 하면 되고, 집까지 와서 문 옆의 소화전에 물건을 넣을 시간이 있었다면 차라리 초인종을 한 번만이라도 눌러보는 것이 더욱 자연스럽지 않겠는가.

인터넷으로 찾아보면 이런 경우가 워낙 많아서 특별히 불만을 표시할만한 일도 아닌 것 같고, 그냥 문 앞에 방치하고 간 것에 비하면 훨씬 나은 상황인 것은 분명하다. 그렇다 해도 약간은 아쉬움이 남는다. 아래 사진은 오늘 받은 물건이다. 사무실에서 USB DAC <-> (진공관 프리앰프) <-> TDA7265 앰프 사이를 연결하기 위한 15 cm짜리 짧은 케이블이다.


낮 시간에는 아내가 필요한 홈패션 부자재를 사러 시내 재래시장에 나갔었다. 시장 안에 있는 분식집에서 호떡을 주문하고 만들어지기를 기다리고 있었는데, 주인 아주머니께서 허리춤으로 와플을 담은 바구니를 건드려서 이를 바닥에 떨어뜨렸다. 난 떨어뜨렸으니 이건 못팔겠다고 할 줄로만 알았다. 그렇게 말하고 나서 우리가 떠난 뒤 이를 모르는 다음 손님에게 판매한다고 해도 어쩔 도리는 없다고 생각했다. 그러나 아주머니께서는 아주 태연하게 와플을 다시 바구니에 담았고, 새로 구운 와플과 함께 진열대에 올려놓는 것 아닌가.

글쎄, 이 광경을 전부 목도하였으니 투철한 소비자 정신에 입각해서 뭐라고 항의했어야 할까. 가장 소극적인 항의는 다음부터 이 분식점에 가지 않는 것이겠지만, 난 특별히 뭐라 하지는 않고 묵묵히 호떡을 들고 그곳을 나왔다.

약간의 아쉬움이 묻어나는 하루였다.