역시 변이 탐지에서도 machine learning

역시 세상에는 내마을을 알고 있는 분들이 있다는..

논문 제목은 Systematic benchmark of state-of-the-art variant calling pipelines identifies major factors affecting accuracy of coding sequence variant discovery 입니다.

doi: https://doi.org/10.1186/s12864-022-08365-3

결론은 누가누가 변이를 잘 찾나 테스트 되겠습니다.

그래서 샘플 준비도 GIAB과 함께 non-GIAB 데이터도 준비하고, 정렬 프로그램도 BWA, BOWTIE, Isaac, Novoalign 다양하게, 변이 caller도 GATK, FreeBayes, Strelka이외에도 DeepVariant나 Clair3와 Octopus 등 다양한 변이 caller를 도입해서 테스트를 진행해 보았다고 합니다.

Clair3와  Octopus는 처음들어 봤네요 @.@ 역시 사람은 공부를 해야함

그리고 GIAB외에도 non-GIAB 데이터를 활용한것은 또 새로운 접근법인듯합니다. 아무래도 많은 연구가 이뤄진 데이터와 그러지 않은 데이터간의 차이가 있을 수 있고, 아직 T2T를 사용한 변이 탐지 벤치마킹 논문을 찾아보진 못했는데 이 또한 조만간 나오겠네요. 이미 나와 있을수도.. 

여튼 이것저것 따져 봤을 때 정렬 프로그램은 일단 Bowtie2는 성능이 저하되는 결과를 보여주기에 사용하면 안되는 방법이고, 변이 탐지에서는 DeepVariant가 좋은 성능을 보이는 결과를 얻을 수 있었다고 합니다.

그리고 GIAB와 non-GIAB 샘플을 비교해본 바 현재 사용하는 분석 방법이 GIAB에서만 좋은 결과를 보이고 non-GIAB에서는 그렇지 않다는 것, aka GIAB 샘플이나 non-GIAB 샘플이나 편견없이 잘 분석된니, 현재 분석 방법들을 활용해도 문제가 없다는 것을 확인하였다고 합니다.

그러나 아직 존재하는 몇 가지 문제점을 지적하기도 했습니다.

기계학습 방법을 사용하는 caller 같은 경우 모델에 따라 변이 탐지에 문제가 발생할 수 있기 때문에 조심해서 사용하여야 하고, 앞서 얘기했던 T2T 표준 유전체를 사용하였을 때 발생할 수 있는 기존 정렬 프로그램이나 변이 탐지 프로그램의 한계점이나 영향(부정적인)은 확인된바가 없기떄문에 검토가 필요하다는 점, GIAB이나 non-GIAB이던 일반적인 벤치마킹에 사용하고 있는 표준 샘플들의 사용 region에 대한 문제, 대게 잘 정의된 변이들은 유전자 영역에 분포하고 있기 때문에 비유전자 영역에 있는 변이 또는 시퀀싱이 잘 되지 않은 영역들에 대한 정보는 제한적일 수 있기에 이런 단점들이 보강된 표준 샘플들이 필요하는 점, NGS 데이터의 경우 지속적으로 개선되고 새로운 방법이나 도구들이 개발되고 있으나 민감도, 재현성과 같은 변이 탐지 성능 지표에 집중되어 있지 분석 시간이나 사용의 편의성, 프로그램 설치 난이도와 같은 비교적 중요하게 다루지 않는 항목들에 대해서도 관심이 필요하다는 것들에 대해서 언급을 해주었습니다. 

간만에 최신 변이 탐지 벤치마킹 논문을 한번 읽어봤는데, 저는 유익한 시간이었는데 어떠셨나요? 

이제 즐거운 주말을 보내러~ :)

ps. 근데 chatGPT한테 요약 해달라고 하니 왜 이 논문에서 Dama dama (사슴)을 왜 그리 찾는지... 처음에는 내가 url을 잘못 입력한줄 알았는데... 

출처: @ye._.vely618

2019-nCoV Tree 그려보기 -2-

2019-nCoV Tree 그려보기 1편에서 어라?
실망하셨던분들을 위해 준비한 2편 되겠습니다.

일단 일전에 말씀드린것과 같이 NCBI 홈페이지가서 다음과 같이 genbank파일을 입맛에 맞게 받으시면 되겠습니다.

refseq말고 다 받고 싶으시면 다 받으셔도 되요.

[Create File] 버튼을 클릭하시면 sequences.gb? 라는 파일로 다운로드가 될것입니다.
(사실 fasta파일로 받으셔도 상관은 없는데... genbank파일로 받으시는게 나중에 더 많은 일을 하실 수가 있으십니다.)

(20년 2월 3일 기준 coronavirus nucleotide 서열 다받으니 대략 300M정도 나왔습니다.)

그럼 이 파일가지고 어쩌라고?
라떼는 fasta파일 가지고.. 어? 마? 그랬는데?
그래서 준비했습니다. gb파일을 fasta파일로 만들어주는 바로그 스크립트!

import os,sys
from Bio import SeqIO
inFile = 'sequences.gb'
for rec in SeqIO.parse(inFile, "genbank"):
seq_id = rec.id
seq_name = rec.name
seq_desc = rec.description
seq_seq = rec.seq
print '>{}|{}\n{}'.format(seq_name,seq_desc,seq_seq)

이렇게 하면 다운바은 sequences.gb파일의 nucleotide 전체 서열을 fasta파일로 만들어주게 됩니다. 만약 나는 다 필요없는데... 하시는 분은
print 전에 if문 넣어서 특정 seq_name의 서열들만 혹은 seq_desc안에 특정 문자열.. 예를 들어 wuhan과 같은 내용이 있는 서열들만 선별해서 fasta파일로 만들 수 있습니다.
참 쉽죠?

그럼 다음 이시간에는 MEGA-X를 사용하는 방법을 얘기해보도록 하겠습니다. 제발~

출처: @sana_twice.09

2019-nCoV Tree 그려보기 -1-

2020년 정월 초하루가 얼마 지나지 않은 나날 동안 한중일이 아닌 전세계가 코로나 바이러스 때문에 홍역을 치르고 있습니다.

초기 방역로 이렇게 급속도로 전세계로 퍼져나가 감염자 및 사망자가 증가하고 있는 가운데 이 짧은 기간 동안에 다양한 연구 결과들이 나오고 있는데...

신종코로나 바이러스 논문들을 보고 있노라면 정렬을 왜그리 해대는지 그리고 무슨 나무를 왜그렇게 그려대는 걸까? 라는 궁금증을 가지시는 분들이 있을 것이라고 생각해서 그런게 무엇인지 그리고 나는 그려볼 수 없는지 하시는 분들을 위해서 MEGA를 이용한 corona virus phylogenetic tree를 한번 그려보고자 합니다.

이거대로 따로하면 당신도 tree 전문가(는 아니고 그냥 tree draw skill +1)

연구자들이 NCBI에 서열들이 올라와 있다고 얘기를 하고 서열들을 서로 비교한다는데
그럼 NCBI는 어떻게 접속하고 서열들은 어떻게 받지?

NCBI는 구글에 물어봐도 되고 네이버에 물어봐도 됩니다.
https://www.ncbi.nlm.nih.gov/

여기에 그냥 검색하듯 coronavirus를 검색하면 다음과 같이 나오는데..

이걸 어떻게 받나?
하나하나 집념으로 다운받으셔도 되고요...

다음과 같은 방법도 있습니다.

아까의 화면에서 Viruses, genomic, RefSeq 항목을 클릭해주면 해당 항목을 만족하는 서열들만 선택이 됩니다.

그리고 오른쪽 상단에 있는 Send to: 라는 항목을 클릭하면...!!
다음 이미지에 나와 있는것처럼 다운 받으시면 되겠습니다.

다운로드한 genbank파일을 parsing해서 비교해보고 싶은 부위의 서열을 득해서
MEGA라는 프로그램을 다음 url(https://www.megasoftware.net/)에서 다운로드 받아 설치한 후 입력한후 서열 정렬 및 Tree 작업을 해주면...

요렇게 이뿌게 Tree를 만들어 줍니다.

이렇게 포스팅을 끝내면 저게 모야? 라고 쌍욕을 날리실겁니다!

그래서 Tree 그릴때 필요한 서열을 수집하는 작업과 관련 스크립트에 대해서도 자세히 적어보는 글은 조만간 투척하도록 하겠습니다. :)

출처: @sana_twice.09