박테리아 분류의 새로운 바람

그동안 Bacillus라는 genus는 최신의 잣대로는 별로 관계가 없는 박테리아까지 한데 품어주는 넉넉한 집 역할을 해 왔었다. 작년이었던가, 유전체학 시대에 걸맞는 분석 결과를 근거로 이렇게 많은 식구를 거느리던 바실러스 가족은 새로운 genus로 하나 둘 독립해 나가가기 시작하였다. 가장 큰 집안이었던 Bacillus subtilis와 Bacillus cereus clade는 변화가 없지만, 나머지 것들은 생소한 genus 명을 갖게 된 것이다.

이러한 혁명적(?) 수준의 변화를 가져다 준 2020년의 주요 논문을 직접 찾아서 읽어보는 것이 올바른 자세임에는 틀림이 없다. 그러나 그것이 성가시다면, 중요한 사실을 요약하여 작성한 가벼운 글을 인터넷에서 찾아서 읽는 것도 나쁘지 않을 것이다. 다음의 글은 찰스 리버 연구소에서 근무하는 Christine Farrance라는 사람이 지난 7월에 게재한 글이다.

Reclassification of Bacteria Happens 

이 글이 실린 Eureka: a dose of science라는 웹사이트는 찰스 리버 연구소의 과학 블로그라고 한다. 블로그를 방문해 보면 히스패닉 문화나 흑인 과학자에게 경의를 표하는 등 과학 그 이상의 정신을 담고 있었다. 별 영양가도 없는 유튜브에서 이상한 뉴스, 해괴한 동영상을 보느니 이런 곳에 올라온 자료를 보는 것이 훨씬 나을 것이다.

새롭게 제안된 분류 체계에서 B. velezensis나 B. siamensis는 어떻게 자리를 잡고 있는지 아직 상세하게 살펴보지는 못하였다. Farrance가 소개한 논문(Patel and Gupta, Int. J. Syst. Evol. Microbiol. 2020;70:406–438 링크)의 그림을 보면 이들은 전부 B. subtilis clade에 속하므로 그 분류학상 위치가 달라지지는 않을 것이다.

Clostridium botulinum도 더하면 더했지 바실러스보다 나은 상황이 아니다. 내가 몇 가지 분석을 해 보았을 때 genus 수준에서 재분류를 해야 함이 명백하다. 그런데 누군가 이를 논문으로 정식으로 제안을 하여 미생물 분류학계에서 인정을 받는다면, 질병관리청과 같이 이 세균을 법률로써 다루어야 하는 공식 기관에서는 많은 어려움을 겪게 된다. 예를 들어서 보툴리눔이 4개의 새로운 종으로 나뉜다고 가정해 보자. 몇 가지 법에서 정의한 고위험성병원체 체계를 이에 맞추어 고쳐야 한다.

대장균(Escherichia coli)는 세균성 이질을 일으키는 Shigella와 유전체학적 기준으로는 다를 바가 없다. 그러나 이를 만약 한 종으로 합쳐 버린다면? 식품에는 허용하는 대장균의 수가 있지만 내가 알기로 이질균은 있어서는 안된다. 그러나 과학적으로 상당한 근거가 있다 하더라도 그 경계를 허물거나 혹은 새롭게 만드는 것의 파급 효과는 상당히 크다. 미생물 분류의 문제가 과학만으로 끝나는 것이 아니라 국민의 안전 및 보건을 다루는 영역까지 건드린다는 것이 엄연한 사실이다.

'아버지를 아버지라 부르지 못하는' 안타까운 일이 미생물 학계에서도 벌어지고 있다. 

EMBOSS "msbar" 명령어는 좀 이상하다

$ wossname msbar
Find programs by keywords in their short description
SEARCH FOR 'MSBAR'
msbar  Mutate a sequence

EMBOSS의 mabar 명령으로 박테리아 유전체 염기서열에 인위적으로 SNP를 도입한 뒤 몇 가지 계산 작업을 하였다. SNP의 수는 10개에서 200까지 총 20 단위로 하였고, 30번씩을 반복하였다. msbar는 도입한 염기 변이의 위치를 별도로 제공하지는 않는다. 

SNP 검출 프로그램을 사용하여 내가 과연 제대로 일을 하고 있는지 궁금하여 각 변이 유전체에 대하여 최종 점검을 해 보았다. 결과는 너무나 이상하였다. msbar를 통하여 도입한 수에 비하여 약 75%에 불과한 숫자가 나오는 것이 아닌가. Sensitivity나 PPV(positive predictive value) 전부 99%를 훨씬 넘는 프로그램을 썼는데 이럴 수는 없다.

T. Seemann의 snp-dists를 사용하여 원본 유전체와 돌연변이 유전체를 직접 비교해 보았다. insertion/deletion은 전혀 없이 염기 치완 돌연변이만 도입하였으므로, 두 파일을 합쳐 놓으면 aligned FASTA 파일과 다를 이유가 없기 때문이었다. 결과는 역시 이상하게 나왔다. 모든 샘플을 다 건드리지는 않았지만 도입한 변이 수의 80%만이 실제로 존재하는 것으로 확인되었다. 이런 결과를 이용할 수는 없다.

CFSAN SNP Pipeline(PeerJ Computer Science 2015년 논문 링크)에서 성능평가용으로 개발한 SNP Mutator 프로그램을 써 보기로 하였다. 설정할 파라미터가 많지만 무엇보다도 돌연변이를 도입한 위치와 염기 정보를 고스란히 제공하기 때문에 결과물을 믿을 수 있다. 이 프로그램으로만든 돌연변이 유전체를 이용하여 필요한 계산을 하고 확인을 위해 SNP를 검출해 보았다. 놀랍게도 도입한 숫자가 거의 그대로 검출되었다. 결국 성능이 의심스러운 돌연변이 도입 프로그램을 쓰는 바람에 일주일을 허비한 셈이 되었다. 하지만 모든 과정을 일괄적으로 실행하는 스크립트와 결과 분석용 R 코드를 그런대로 잘 정비해 두었기에 SNP Mutator로부터 새로 시작하는 것은 그렇게 어렵지 않았다.

SNP Mutator의 중요한 특징은 다음과 같다.

• Mutations can be any number of single-base substitutions, insertions, and deletions at randomly chosen positions, uniformly distributed across the genome.
• Mutations can be chosen from a subset (pool) of all possible positions.
• Replicates can be partitioned into multiple groups with each group sharing a pool of eligible positions.

의심의 눈초리를 거두지 말아야 한다. 언제 어떤 실수를 할지 모르므로.

[하루에 한 R] for loop를 이용하여 데이터프레임 생성하기

Column/Row를 구성하는 벡터가 전부 준비된 상태에서 데이터프레임을 만드는 것은 비교적 쉽다. 기존에 존재하는 데이터프레임을 cbind() 함수로 합치는 것도 마찬가지이다. 그런데 합쳐야 할 벡터가 수십 개 이상 된다고 하자. 이를 일일이 타이핑하기는 곤란하다. 반복문을 쓰지 않고서는 매우 어렵다. 

EMBOSS의 msbar 명령을 이용하여 박테리아 유전체 서열에 염기치환돌연변이를 일정하게 도입하였다. 그 각각이 원본 유전체 서열('ancestor'에 해당)에 대하여 얼마나 멀어졌는지를 Mash로 측정해 보았다. 목표로 하는 돌연변이 수에 대하여 msbar는 총 30회를 실시하였고, 돌연변이의 수는 10, 20, 30..200으로 하였다. 결과물은 전부 하나의 파일(dist.tab)에 저장하였다. 

$ head -n 5 dist.tab 
GCF_000008445.1.010.00.c.fnm	GCF_000008445.1.fna	1.61302e-06	0	9999/10000
GCF_000008445.1.010.01.c.fnm	GCF_000008445.1.fna	4.8398e-06	0	9997/10000
GCF_000008445.1.010.02.c.fnm	GCF_000008445.1.fna	0	0	10000/10000
GCF_000008445.1.010.03.c.fnm	GCF_000008445.1.fna	1.61302e-06	0	9999/10000
GCF_000008445.1.010.04.c.fnm	GCF_000008445.1.fna	9.68178e-06	0	9994/10000

첫번째 컬럼이 돌연변이 유전체 FASTA 파일의 이름이다. fna라는 흔한 확장자를 fnm('m' for mutation)으로 고쳤다. 파일명 중간쯤에 '010'이라는 숫자가 바로 돌연변이 수에 해당한다. 만약 10, 20...으로 표기를 했다면, 100을 넘어가는 결과와 함께 다룰 때 sort에 주의해야 한다. 그래서 약간 귀찮지만 이런 파일명을 쓴 것이다. 

하나의 돌연변이 숫자에는 총 30개의 Mash distance(위 자료에서 세번째 컬럼)가 있다. 이것을 뽑아내어 각기 다른 컬럼으로 만들고자 하는 것이 오늘의 목표이다. 돌연변이의 수는 총 10에서 200까지 10 단위로 증가하므로 20개가 되고, 이것이 목표로 하는 데이터프레임의 전체 컬럼 수가 된다. Row의 수는 30개가 된다.

일반적인 슬라이싱을 통해서 원본 데이터프레임의 일부를 잘라낸 뒤 새 데이터프레임에 한 컬럼씩 붙여 넣는 방법이 있을 것이라 생각했는데 잘 되지 않았다. 검색을 거듭한 결과 리스트를 써야 함을 알게 되었다. 참조했던 웹사이트가 어디였는지는 기억을 하지 못하겠다. do.call() 함수를 쓰는 방법을 잘 알아야 될 것이다.

> header = sprintf(seq(10,200,10), fmt='%03d')

> header
 [1] "010" "020" "030" "040" "050" "060" "070" "080" "090" "100" "110" "120"
[13] "130" "140" "150" "160" "170" "180" "190" "200"
> df = read.table("dist.tab", sep="\t", row.names=1)
> datalist = list()
> for(i in header) {
+     pattern = paste(".", i, ".", sep="")
+     datalist[[i]] = df[grepl(pattern, row.names(df)), 2]
+ }
> df.2 = do.call(cbind, datalist)
> View(df.2)

최적화된 R code라고 확신을 할 수는 없다. '.010.' 형태의 패턴을 추출할 때, 앞뒤의 점은 정규표현식에서 임의의 문자 하나에 대응하는 점이 아니라 문자 그대로의 '.'임이 중요한데, 이것이 grepl() 함수에서 정확히 인식되고 있는 것인지를 잘 모르겠다. grep()은 인덱스를 반환하지만 grepl()은 TRUE 또는 FALSE를 반환한다고 한다. 어쨌거나 View(df.2)로 확인을 하면 다음과 같이 내가 의도한 데이터 변환이 이루어졌음을 알 수 있다.

나의 '데이터 주무르기'는 그야말로 마음 내키는 대로이다. 모든 조작을 Perl에서 하던 시절도 있었고, 요즘은 bash script와 awk/sed 등을 적절히 섞어서 쓰는 빈도가 늘어났다. 그러다가도 예외 사항이 별로 없이 매트릭스 형태로 잘 짜여진 데이터 파일을 보면 R에서 작업을 하고 싶어진다. 사용했던 코드를 노트 파일과 컴퓨터 작업 공간에 잘 저장해 두는 것도 좋은 버릇이라고 자평한다.

점점 까다로와지는 공공기관의 전산 보안

주말에 집에서 일을 하기 위해 사무실에서 쓰던 노트북 컴퓨터를 들고 퇴근하였다. 월요일(어제) 외부에서 있을 과제 중간 점검 발표를 준비하기 위해 파워포인트 자료를 마지막으로 점검하려는데 보안 정책에 위배된다면서 와이파이 접속이 차단되는 것이었다. 아, 드디어 이런 날이 오고 말았구나! 집에서 이더넷 케이블을 직접 연결해 보지는 않았다. 와이파이 접속이 안되는 노트북이 무슨 쓸모가 있는가?

컴퓨터 전원을 껐다가 다시 켰더니 이제는 얼굴 인식을 위해 작동하던 내장 카메라도 제 기능을 하지 못하여 PIN을 입력해야만 했다. 노트북이 연구소 전산망 내에 있지 않은 상황이라면 이를 잠재적인 위험 환경으로 간주하고, 입출력이 이루어질 수 있는 주변 기기를 통제하는 보안 강화 모드를 강제하는 것으로 여겨진다. 그것이 이제는 너무나 당연하게 여겨지는 와이파이라 해도 말이다.

연구소 안에서 상용메일 접속이 차단된 것은 꽤 오래 전의 일이다. 아마존웹서비스와 같은 클라우드도 당연히 되지 않는다. USB 매체 사용도 매우 까다로와서 아예 쓰지 않는다. 리눅스를 설치한 컴퓨터는 상황이 약간 낫지만, 이런 사실을 함부로 발설(?)하다가는 그나마도 쓰지 못하게 될까 두렵다. 이러다가 USB 마이크로폰도 작동이 안되는 날이 오는 것은 아닐까? 사무실 데스크탑에 연결한 웹캠은 다행스럽게도 아직 작동을 잘 한다.

공공기관에서 접속하면 안되는 '나쁜' 웹사이트를 국정원에서 직접 관리하는 것이 정말 옳은 일인지 잘 모르겠다.

공공기관의 전산망 보안은 외부의 해킹 시도나 악성 코드 침투를 방지하는 것으로부터 한층 더 강화되어 근무자가 내부의 자료를 외부로 가져가지 못하게 하는 것으로 바뀌어 가는 추세이다. 물론 일반 회사는 이보더 더 강도 높은 보안 정책을 유지한다. 가장 대표적인 것은 문서에 자물쇠를 거는 DRM이 있고(보안 해제 권한은 몇몇 높은 사람들만 갖고 있음), 프린터도 철저하게 통제한다. 휴대폰은 당연히 보안 프로그램이 깔려 있어서 출근을 하여 회사 공간에 들어가면 카메라가 작동되지 않는다.

정부 조직은 내가 근무하는 출연연보다는 불편하고, 일반 기업보다는 약간 느슨한 정책을 쓰는 것으로 알고 있다.

아마도 노트북 컴퓨터의 경우 외부 반출을 위해 공식적으로 거치는 절차가 있을 것이다. 문제는 연구비로 구입한 노트북이든 개인이 구입한 노트북이든 연구소 내에서 전산망에 접속하려면 보안 솔루션을 깔아야 한다는 것이다. 아주 엄격하게 말하자면, 개인 노트북을 연구소에 가져와서 연구소 전산망에 접속하면 안 된다. 회의실 근처에 있는 공개형 와이파이 정도만 써야 한다. 그러나 이것으로는 내부 리눅스 서버에 접속을 할 수 없다. 코로나 시대를 맞아서 재택 근무를 권장하면서도 실제 환경은 이를 따라가지 못한다.

민감한 정보를 다루는 정보 기관이나 민간 기업이 아니라면 약간은 유연한 정책을 실행하면 안 되는 것일까? 모처럼 연구비로 성능이 좋은 노트북을 구입하였는데 외부에서 쓰지 못한다면 여간 안타까운 일이 아니다. 

공식 반출 허가를 받지 않고 노트북을 외부에 너무 들고 다니는 것도 문제가 될 소지는 있다. 아예 집에다 가져다 놓고 가정용(=개인용)으로 쓰거나, 심지어는 퇴직 이후에도 반납을 하지 않는 철면피 같은 사람도 있으니 말이다.

정말 마시기 싫은 음료

새벽 4시에 일어나서 벌컥벌컥 들이마시는 쿨프렙 산.정말 먹기 싫은 맛이다. 누군가 평하기를 이온음료에 조미료를 탄 것 같은 시큼하고 느끼한 맛이라고 했는데, 그 말이 딱 맞다. 1리터를 조제해서 1시간 안에 마시고, 이어서 맹물 1리터를 또 1시간 안에. 이 과정을 어제 저녁 일곱 시에도 했고, 새벽 네 시에도 한 차례를 더 해야 한다. 

으아! 이 찝찌름한 액체를 아직 250ml 더 마셔야 한다. 그러고 나서 또 화장실을 들락거리게 될 것이다. 아침에는 일찍 병원 건강검진 센터로 가야 하니 그 사이에 잠을 제대로 자기는 틀렸다.

폴리에틸렌 클리콜과 아스코르빈산을 주성분으로 하는 장세척 용도의 의약품이다. 아마도 엄청난 삼투압에 의햐여 설사를 유발하는 것일 게다. 안전성은 충분히 입증된 의약품이지만, 장에 큰 부담을 주는 일은 아닐까? 장내 미생물 살림터를 홀랑 뒤집어 엎는 일이니 말이다. 

장이 깨끗한 편이라 이제 겨우 생애 두 번째 대장 내시경 검사를 받는다. 최초 검사는 5년 전이었다. 만약 이것을 매년 받으라면? 정말 내키지 않는다. 어차피 의식하 진정(소위 '수면') 내시경이라서 검사 과정 자체는 괴롭지 않은데, 전날 밤부터 당일 새벽까지 의식처럼 치러야 하는 '쿨프렙 대잔치'는 여간 고역이 아니다. 아픔을 유발하는 것은 아니지만 결코 유쾌하지 못한 자극을 통한 괴로움이란...

맹물이 가장 맛있게 느껴질 때는? 격하게 땀을 흘리고 나서 마실 때가 아니다. 쿨프렙을 1리터 마시고 나서 바로 뒤에 맹물을 마실 때이다.

마지막 남은 쿨프렙을 다 마셔버렸다. 뱃속이 다시금 요동을 친다.

[독서 기록] 다정한 무관심 및 내 친구 압둘와합을 소개합니다

다정한 무관심(함께 살기 위한 개인주의 연습) 지은이: 한승혜

개인주의는 이기주의와 동의어가 아니다. 개인주의의 반대발이 전체주의라는 것을 떠올린다면 이런 오해는 쉽게 풀릴 것이다. SNS에 짧게 쓴 글들을 정리하여 책으로 묶었다. 나는 블로그를 제외하면 남들이 다 한다는 SNS를 거의 하지 않는 사람이라서, 이런 공간에서 의견을 올리고 공감을 얻어가 토론을 한다는 것 자체를 중요하게 생각하지 않는다. 나에게는 '없는 세상'이라고나 할까... '메타버스'도 그런 것이 아닐까한다. 여담이지만 세상에 가장 쓸모없는 것은 공개 게시판이나 덧들이라고 생각한다. 어쨌거나 이것도 책을 만드는 한 방식이 될 수도 있겠다.

내 친구 압둘와합을 소개합니다(어느 수줍은 국어 교사의 특별한 시리아 친구) 지은이: 김혜진

시리아의 좋은 집안에서 태어난 청년이 자국을 방문했던 한국인 여행자와 겪은 일을 계기로 국교도 없던 한국으로 법학 공부를 하러 내한하게 된다. 그 사이에 시리아에서는 민주화 혁명과 내전이 일어나고 가족들은 난민 신세가 된다. 우연한 기회에 압둘와합을 알게 된 저자는 얼떨결에 시리아를 돕기 위한 단체를 만들게 되고, 어려움 속에서도 난민 어린이들을 가르칠 수 있는 학교까지 설립을 하게 된다. 

중간 중간에 시리아를 소개하기 위해 압둘와합이 직접 쓴 글이 삽입되어 있다. 만약 내가 한국을 소개하는 글을 써 달라는 부탁을 받는다면, 정확하고도 공정한 글을 쓰기 위해 얼마나 고민을 하게 될까? 어쩌면 소심한 나는 부담스럽다고 거절을 할지도 모른다.

현재진행형인 아프가니스탄 사태에 밀려서 시리아에서 일어났던 일에 대해 한국인들이 관심을 덜 갖게 된 것 같다. 한국과 시리아는 여전히 미수교 상태이고, 압둘와합은 한국인으로 귀화하였다.

교양으로 읽는 세계 종교사 지은이: 시마다 히로미, 옮긴이 김성순

이리저리 건너뛰어가면서 읽다가 대여 기간이 다 되어서 부득이하게 손을 놓게 된 책이다. 최근 불교에 대한 관심을 갖게 되고, 이슬람 세계에서 벌어지는 일을 이해해 보고자 빌리게 되었다. 여름날 산사에 들르면 백중절(음력 7월 15일)과 관련한 행사로 분주한 것을 볼 수 있다. 사후세계나 지옥, 윤회와 관련한 개념은 원래 석가모니가 처음 가르침을 시작하던 당시에는 없는 개념이 아니었을까? 삼국시대때 우리나라에 전래된 불교는 중국을 거치면서 원래의 모습과는 다르게 된 것이 아니었을까? 이런 궁금증을 해결해 보려고 책을 읽기 시작했으나 아직 답을 얻지는 못했다.

똑똑한 다중서열정렬(multiple sequence alignment, MSA) 트리밍 도구, ClipKIT

이발소를 연상하게 하는 ClipKIT의 로고부터 살펴보자. 논문은 2020년 PLoS Biology에 "ClipKIT: A multiple sequence alignment trimming software for accurate phylogenomic inference"라는 제목으로 발표되었다.

계통수를 그리려면 Clustal 시리즈의 소프트웨어, 혹은 muscle, mafft, prank 등의 프로그램이 만든 다중서열정렬(MSA) 결과물을 이용해야 한다. MSA 내부에 존재하는 무의미한 gap 등을 처리하는데 보통 'trimAl -automated1'을 즐겨 사용하였다. 문제는 십만 개가 넘어가는 바이러스 게놈의 MSA를 trimAl에 밀어넣으니 도무지 끝이 날 조짐이 보이지 않는 것이었다. 이 소프트웨어는 다중 쓰레드를 지원하지도 않으므로 그냥 기다리는 것 외에는 방법이 없었다.

대안을 찾기 위하여 검색을 거듭한 결과 몇 달 전에 ClipKIT이라는 소프트웨어를 알게 되었다. 트리밍 전략이 매우 단순하고, 상식적인 시간 안에 결과를 내어놓기 때문에 즐겨 사용하게 되었다. 

MSA를 트리밍하려면 각 컬럼(즉 염기 위치)에 갭이 몇 개인지, 서로 다른 염기가 있는지, 혹은 동일한 염기가 있는지 등을 수치화해야 한다. MSA의 품질을 나쁘게 만드는 가장 중요한 요인은 갭이므로(물론 진화적인 의미가 있는 갭도 있겠지만), 이것이 MSA 트리밍 도구에 제공하는 옵션에서 매우 중요한 역할을 한다.

갭과 관련한 수치는 그 정의를 명확하게 이해하는 것이 중요하다. 이것을 제대로 알지 못해서 초창기에 많은 혼란을 겪었다. 예를 들어서 100개의 염기서열을 정렬하였더니 어느 위치에서 15개는 염기가 있고 85개는 갭이 있다고 하자. trimAl에서 사용하는 gap score는 0.15가 된다. 이 프로그램에서 사용하는 gap score의 정의는 얄궂게도 갭이 없는 것들의 비율이다.

Sg = 1.0 - (number of sequences with a gap / n)

정의를 이렇게 해 놓으니 갭이 많을수록 gap score는 작아진다. 상식적으로 생각해도 이해가 잘 안된다. 갭이 많으면 gap score도 크게 나타나는 것이 자연스럽지 않을까? 어쩌면 'score가 커야 좋은 것'이라는 인식에 맞추어서 이렇게 정의를 내렸는지도 모른다. 이를 트리밍에 적용하려면 'trimal -gt 0.190'과 같이 적용한다. 즉, gap score가 최소한 0.19가 되지 않는 컬럼은 MSA에서 잘라(트림)내라는 뜻이다. 100개의 시퀀스 중에서 19개 혹은 그 이상이 염기를 제시하고 있어야 트리밍을 당하지 않는다는 것이다. '-nogaps'는 '-gt 1'과 동등하며, 어느 한 시퀀스라도 갭이 있으면 그 컬럼을 제거한다는 뜻이다. 반대로 '-noallgaps'는 모든 염기서열이 갭으로만 구성된 컬럼을 제거하게 된다.

CipKIT에서는 gappyness라는 지표를 쓴다. gapyness의 정의는 다음과 같다.

Gappyness = number of gaps / entries in alignment

갭이 많으면 gappyness도 커진다. Gap score보다는 직관적으로 이해하기 더 좋은 지표이다. 그리고 gap score + gappyness = 1.0이 됨을 쉽게 알 수 있다. ClipKIT에서 gappyness를 기준으로 트리밍을 하려면 'clipkit -m gappy 0.8'이라고 옵션을 제공하면 된다. 기본 gappyness는 0.9이다. 

정해진 threshold 값을 옵션으로 주지 않고 프로그램이 스스로 결정하게 만들 수도 있다. CipKIT에서는 smart-gap이라는 방법이 이에 해당한다. 이를 사용하면 일정한 값을 주어서 실행한 것보다 트리밍한 뒤 남는 MSA의 길이가 더 길다고 한다.

ClipKIT의 개발자에게 경의를 표하는 바이다.