WGS, RNA-seq, sRNA-seq 등 시퀀싱 테크닉에따라 그리고 샘플의 종(species)에 따라 QC기준은 달라질 수 밖에 없지만 사람을 기준으로 대표적인 QC기준을 설명하고자 한다.

크게 실험단계에서의 라이브러리 제작이 잘 되는지 QC가 있고 시퀀싱이 진행된 후에 BI 기준의 QC가 있다.

실험 QC는 잘 모르니 생략하고 BI기준의 시퀀싱 이후 단계만 요약하고자 한다.

가장 첫 번째 단계는 Illumina 시퀀서 기준으로 PF(passing filter)가 있다.

passing filter란 시퀀싱 단계에서 일어나는 cluster가 얼마나 잘 구성되었는지 여부를 의미한다. cluster들끼리 오버랩되면 에러가 생길 가능성이 커져 오버랩이 일어난 cluster는 모두 제거하는데 이 비율이 얼마나 되는지 수치가 PF값이다.

https://www.illumina.com/content/dam/illumina-marketing/documents/products/technotes/hiseq-x-percent-pf-technical-note-770-2014-043.pdf

PF값은 참고용으로만 사용할 뿐 실제 분석에서는 대부분 raw file인 fastq 파일부터 시작하기 때문에 QC는 아래 값들을 주로 사용한다.

fastq파일의 Number of reads, GCratio, Q20, Q30, N(%), Duplication rate, Adapter sequence

BAM파일의 Mapping rate, House-keeping gene expression, Insert size distribution

각각의 수치를 구하는 방법은 여러 종류의 tool이 있으나 FastQC가 fastq 수준에서는 가장 많이 사용한다.

https://www.bioinformatics.babraham.ac.uk/projects/fastqc/

BAM파일은 alignment tool에 따라서 수치를 제공하는 프로그램도 있고 samtools stats 명령어를 통해서도 확인 가능하다. 

HLA는 MHC의 한 부분이다. HLA는 세포 표면 단백질 항원으로 모든 유핵 세포 표면에 발현되어 있고 면역세포가 '정상적인 자가 세포'인지 아닌지 판단하는 지표로 사용한다.

염색체 6p21에 약 3Mbp 정도의 크기이며 polymorphic하다. 후술하겠지만 locus에 따라 수만 종류의 allele이 밝혀졌으며 여전히 밝혀지지 않은 allele들도 다수 존재하는 것으로 보인다.

HLA는 Class I(A, B, C), Class II(DP, DM, DO, DQ, DR) 그리고 Class III로 이루어져 있다.

Class I은 세포 내에서 외부로 펩타이드를 제시한다. 만약 세포가 바이러스 등 외부 요인에 의해 감염되었다면 그 조각을 세포막 외부에 제시하고 killer T-cell 등이 이를 인지하고 면역 체계를 발동시키도록 하는 역할이다. 

Class II는 세포 밖에서 항원을 제시하여 T-림프구가 이를 인지하도록 한다. T-림프구가 이를 인지하면 T-helper 세포를 자극시켜 이 항원에 대응하는 항체를 가지는 B세포의 합성을 유도한다. 

Class III는 complement system 또는 complement cascade라고하는 면역 체계를 활성화 시킨다.

HLA는 주로 장기 이식 거부 반응이나 자가면역 질환(1형 당뇨, 실리악 스프루 등)과 연관이 깊다.

HLA는 변이가 매우 잘 일어나기때문에 동일한 HLA allele을 가지는 사람이 있을 확률은 매우 낮다. 장기 기증을 신청하면 HLA genotype을 미리 조사했다가 같은 type의 수여자가 생겨나면 알려주는 식으로 진행되며 가족일 경우 확률이 높은 것은 당연하다. 

HLA의 명명법은 위와 같다. Field 1은 2digit, Field 2는 4digit 등으로 말하기도 한다.

최신 업데이트의 HLA database는 IMGT-HLA database에서 얻을 수 있다.

https://www.ebi.ac.uk/ipd/imgt/hla/

출처 -

https://en.wikipedia.org/wiki/Human_leukocyte_antigen#Functions

http://hla.alleles.org/nomenclature/naming.html

SnpEff는 기본적으로 지원하는 genome 데이터가 있지만 manual하게 빌드하여 SNP의 효과를 예상 할 수도 있다.

SnpEff 설치 후 설치 폴더 내에 생성되는 snpEff.config 파일을 수정해 주어야 한다.

/PATH/TO/INSTALL/SnpEff/snpEff.config 파일 내에 아래 내용을 추가한다.

# genome for test 
test.genome : testtest

test.genome의 test가 이 genome의 식별 코드이다. 뒤에 testtest는 어떠한 영향을 주는지는 잘 모르겠다. 다른 글을 참조했을때 보통은 종 이름을 쓰는 듯 하다.

이후에 /PATH/TO/INSTALL/SnpEff/ 폴더 내에 /data/test/ 폴더를 생성하고 빌드하고자 하는 genome 서열과 annotation 파일을 아래의 형식으로 옮긴다.

mkdir data/test/
cp original_genome.fasta data/test/sequences.fa
cp original_annotation.gtf data/test/genes.gtf

모든 준비가 완료되면 빌드를 시작한다.

java -jar snpEff.jar build test

Illumina 시퀀싱에서 약 10-12개의 염기가 균등하게 분포하지 않는 패턴을 보인다. gDNA에서는 조금 더 드물지만 mRNA-seq에서는 대부분의 데이터가 이러한 패턴을 보이는데 이유를 찾아보았다.

Illumina에서는 이러한 현상의 원인을 랜덤프라이머를 제작하여 시퀀싱을 진행하지만 랜덤 프라이머가 완전한 랜덤이 아니기때문이라고 얘기한다.  

아래 plot은 bisulfite-seq 이다. bisulfite 처리로 인해 CtoT 변화로 C의 비율은 낮고 T의 비율이 높게 나온다. 하지만 그와 별도로 여전히 10개의 염기의 비율이 특이적이다.

해당 부분은 분석에 크게 영향을 주지 않으니 무시하고 진행하여도 상관없다.

Reference -

http://seqanswers.com/forums/showthread.php?t=11843

http://nar.oxfordjournals.org/content/38/12/e131

DNA methylation은 DNA 염기 서열의 변화 없이 변화를 일으킨다. 유전자 프로모터에 위치할 때 전사를 억제하는 방향으로 작용할 때가 많으며 포유류에서의 DNA 메틸화는 정상적으로 발달할 수 있게 하는 핵심 역할을 하며 그 외에도 genomic imprinting, X-chromosome inactivation, repression of transposable elements, aging, and carcinogenesis등에 관여하는것으로 알려져 있다.

Cytosine과 adenine이 메틸화 될 수 있으며 adenine의 경우 박테리아나 식물에서 그리고 드물게 포유류에서 관찰되지만 비율이 적어 많이 연구되지는 않았으며 cytosine에는 eukaryotes와 prokaryotes 모두에서 빈번히 일어난다.

식물이나 다른 유기체에서 DNA 메틸화는 CpG, CHG or CHH (Hsms A,T 또는 C)상황에서 일어난다. 그러나 포유류에서는 대부분 CpG에서 발견된다. Non-CpG 메틸화는 배아 줄기세포에서 발견되거나 신경 발달과정에서 확인할 수 있다. 그러나 최근에는 Non-CpG 메틸화가 질병과 연관이 있을 수 있다는 연구보고가 있어 중요도가 올라가고 있다.

Reference -

https://en.wikipedia.org/wiki/DNA_methylation

Pyrosequencing, methylation-specific polymerase chain reaction (PCR), direct Sanger sequencing이 프로모터 영역이나 CpG island 등 특정 영역의 메틸화 정도를 보려고 할 때 사용하는 기술이다. 이 기술들은 유용하지만 정확도가 낮고 read length가 짧으며 수율이 낮다는 단점이 있다.

이후에 적은 비용과 많은 양의 DNA가 필요하지 않으면서도 게놈 영역 전체를 커버할 수 있는 Microarray를 사용한 기술을 사용한 메틸화 분석 기술이 나왔으나 이는 depth에 의해 결과가 영향을 받을 수 있다.

NGS기술은 단일 염기 수준의 정확도로 거의 모든 CpG 사이트를 구분할 수 있는 기술임을 보였으나 여전히 density-biased, deficient in robustness and consistency, or incapable of analyzing 5mC specifically 등의 문제가 있다.

가장 많이 사용하는 기술의 요약이다.

1. Affinity Enrichment-Based Methods

antibody(MeDIP-Seq)나 binding protein(MBD-Seq)을 사용하는 방법으로 특정 영역을 당겨서 시퀀싱을 하고 나머지는 버리는 방법이다. enriched CpG영역에 대한 분석에 용이하다.

2. Restriction Enzymes-Based Methods

MspI 제한효소를 하용하여 CCGG motifs를 절단한다. 시간과 비용이 적게 들고 DNA도 소량만 있으면 되지만 게놈 영역 전체를 골고루 커버할 수 없다. (특정 영역이 chromatin structure 등에 의해 더 잘 잘리거나, 덜 잘리거나 하는 경향이 있을 수 있음.)

3. Bisulfite Conversion-Based Methods

bisulfite 처리를 하면 보통의 C가 U라 바뀌지만 메틸화 C는 바뀌지 않는다. 따라서 bisulfite처리를 한 샘플과 게놈 sequence를 비교하여 실제 메틸화 C가 어디 있는지를 알아낼 수 있으며 따라서 이 분석에 맞는 alignment프로그램을 사용하여야 한다. 가장 많이 쓰이는 프로그램은 BISMARK.

WGBS는 가장 많은 정보를 담을 수 있지만 비싸다는 단점이 있다.

4. Oxidative Bisulfite Conversion-Based Methods

최근에 개발된 기술로 5hmC와 5caC, 5fC를 포함하여 cytosine modification을 찾아낼 수 있다. 방법으로는 Ox-BS나 TAB-Seq 등이 있다.

5. Capture-Based Methods

기존의 방법들에 비해 많은 진보가 있는 기술로서 whole genome sequencing에 비해 cost-effective 하면서도 수율도 높고 특히나 complex regions에서 강점을 보인다.  MethylCap-seq은 CpG islands 외의 질병과 관련되어있는 메틸화 영역에 대한 연구가 가능하며 reproducibility도 높음을 여러 샘플에서 보여주었다. MethylCap-seq 외에도 SeqCap Epi CpGiant 등 bisulfite-converted DNA를 사용하여 방법들이 있다. 

6. Third-Generation Sequencing

최근에는 chemical conversion이 없이도 DNA modification 분석이 가능한 기술이 가능하다. SMRT DNA sequencing은 kinetics of DNA polymerase를 인지함으로 modifed DNA를 찾아낼 수 있다. nanopore sequencer또한 DNA base의 modification 여부를 reading과정 중에 찾아낼 수 있다.

참고 -

Barros-Silva, Daniela et al. “Profiling DNA Methylation Based on Next-Generation Sequencing Approaches: New Insights and Clinical Applications.” Genes vol. 9,9 429. 23 Aug. 2018, doi:10.3390/genes9090429

KEGG Mapper 사용법

KEGG Mapper는 KEGG PATHWAY에 실제 유전자의 발현량을 색으로 입혀 시각적으로 한 눈에 알아볼 수 있게 하는데 도움을 주는 웹 기반 프로그램이다.

홈페이지 : 

https://www.genome.jp/kegg/mapper.html

위의 예시에선 유전자에 붉은 색은 높은 발현량을 의미한다.

실제 데이터를 만들어 넣는 법은 매우 간단하다.

위의 홈페이지에서 Color Pathway 항목으로 가면 데이터를 입력할 수 있는 항목들이 있다.

Select KEGG pathway map : 색을 입히고자 하는 kegg pathway의 category를 입력하면 된다.

Enter data : examples를 보면 어떤식으로 데이터를 입력해야하는지 보여주는데 아래에서 실제 데이터로 설명하도록 하겠다.

Option : 색을 직접 RGB로 입력하거나 상대적인 값을 측정해 색을 부여하는 방식을 고를 수 있다.

여기서 가장 까다로운게 입력 데이터의 유전자 ID를 KEGG ID로 맞춰야 한다는 것인데 이는 아래 포스팅을 참조하면 해결할 수 있다.

2018/11/15 - [bioinformatics] - Gene ID conversion

위의 포스팅에서 엑셀로 KEGG Gene ID를 찾은 후 아래처럼 입력한다.

저장은 반드시 텍스트 (탭으로 분리) (*.txt)로 저장해야 한다.

위와같은 결과를 얻을 수 있다. 

condtion을 3개 주었기 때문에 가장 위에 contion1,2,3 탭이 보인다. 클릭을 하면 색만 바뀌기 때문에 비교하기 쉽다.

Reference -

https://www.genome.jp/kegg/mapper.html

Gene ID conversion

하나의 유전자를 지칭하는 명칭은 ensembl, kegg, refseq 등 분석 방법에 따라 달라지고 분석 중에 gene id를 다른 방식으로 맞춰야 하는 일들이 생긴다.

R에서 biomaRt등의 라이브러리를 사용하여 스크립트 내에서 변환하는 방법도 있지만 web 기반의 tool를 사용해서 바꾸는 방법에 대해서 설명하고자 한다.

홈페이지 : 

https://biodbnet-abcc.ncifcrf.gov/db/db2db.php

ID List에 변환하고자 하는 유전자 목록을 넣었고 ID가 ensembl ID이기 때문에 input에는 Ensembl Gene ID, 결과는 Gene Symbol로 맞추었다.

Organism은 9606이 human이며 다른 종을 찾고싶다면 Taxon ID를 클릭해서 들어가면 검색이 가능하다. 이 항목은 option이기 때문에 꼭 넣어주어야 하는 것은 아니다.

입력을 많이 넣지 않았기 때문에 넣어준 ID가 하나 빼고는 다 치환된 것을 확인하였다. 

Result in Excel을 클릭하여 엑셀파일로 받으면 기존의 데이터에 덮어쓰거나 추가 열을 만드는 등 편집하기 쉽다.

Reference -

https://biodbnet-abcc.ncifcrf.gov/db/db2dbRes.php

SRA data 다운로드받기

NCBI에서 SRA data를 받을 수 있는 방법은 세 가지 이다.

1. SRA Toolkit

NCBI SRA 다운로드 페이지 :

https://www.ncbi.nlm.nih.gov/sra/docs/toolkitsoft/

압축을 풀면 바로 bin폴더가 생성되어 있고 이 중에 fastq-dump를 사용하여 받을 수 있다.

사용법은 

fastq-dump -A [accession number]

2. ascp utility

aspera 홈페이지 :

https://downloads.asperasoft.com/en/downloads/50

sh 파일을 다운르도 후 root로 진행하면 된다. (보류)

3. wget

ascp가 고속 전송을 지원하기 때문에 FTP를 사용하는 것 보다 10배는 빠르지만 별도의 설치나 key파일을 필요로 하기 때문에 초기 셋팅이 번거롭다.

wget은 가장 간단하게 사용할 수 있지만 속도가 느리다.

wget /sra/sra-instant/reads/ByRun/sra/{SRR|ERR|DRR}/<first 6 characters of accession>/<accession>/<accession>.sra

wget ftp://ftp-trace.ncbi.nih.gov/sra/sra-instant/reads/ByRun/sra/SRR/SRR304/SRR304976/SRR304976.sra

아래의 bash script를 만들고 "sh sradownload.sh SRR304976" 라고 입력하면 sra를 다운받아서 fastq까지 만들어준다.

wget ftp://ftp-trace.ncbi.nih.gov/sra/sra-instant/reads/ByRun/sra/${1:0:3}/${1:0:6}/${1}/${1}.sra

fastq-dump --split-3 ${1}.sra

Reference -

https://www.ncbi.nlm.nih.gov/sra/docs/toolkitsoft/

somatic mutation과 germline mutation

somatic mutation(체세포 돌연변이)과 germline mutation(생식세포 돌연변이)의 개념이 가장 많이 등장하는 곳은 암 유전체일 것이다. 

두 돌연변이 모두 DNA 서열상에서의 돌연변이를 가리키며 각각의 개념의 정의는 아래와 같다.

체세포 돌연변이 - 체세포에서 돌연변이가 발생. 몸 전체에서 일부 영역의 세포만이 돌연변이를 갖는다. 

생식세포 돌연변이 - 생식세포 돌연변이는 부모 세대로부터 물려받은, 배아가 형성될 때부터 이미 가지고 있는 돌연변이로서 몸 전체 어느 세포든지 같은 돌연변이를 갖는다.

두 돌연변이가 암 유전체에서 많이 등장하는 이유와 구분하는 방법은 무엇일까

대부분 암은 체세포 돌연변이로 발생한다. 사람이 살아가는 동안 ROS, UV 등으로 DNA가 계속 변이의 위험에 노출되는데 이때 세포 성장이나 세포 주기에 관여하는 유전자에 돌연변이가 생기면 암세포로 변이되는 것이 가장 일반적인 암 발생 패턴이다. 그런데 이러한 유전자에 생식세포 때부터 돌연변이가 있다면 애초에 사산될 확률이 높기 때문이다.

따라서 암 치료 시 환자 개개인의 암세포를 채취해서 생식세포 돌연변이보다는 체새포 돌연변이를 찾아내고 해당 변이에 맞는 약물을 투여하는 것이 개인 맞춤형 치료이다.

물론 BRCA 돌연변이같이 유전되는 생식세포 돌연변이도 존재한다. 이러한 암을 대략 5~20% 정도로 예상하고 있다.

두 돌연변이를 구분할 수 있는 가장 확실한 방법은 N-T pair 비교이다.

환자의 정상 세포와 암세포를 각각 채취하여 GATK에서 제공하는 Mutact2 또는 이 같은 목적의 프로그램을 돌리는 것이다.

이론적으로 염색체의 염기 서열은 Homozygous인지 Heterozygous인지에 따라 0%, 50%, 100%를 가지게 된다. 

ex) G/G or G/T or T/T 이런식으로 maternal, paternal을 가진다. 

하지만 체세포 변이가 일어나면 그 비율이 달라질 수 있다. 특정 염색체에서만 변이가 일어나거나 CNV 등이 일어나서 염기 비율이 달라진 세포와 일반적인 세포가 섞여있기 때문이다. 10개의 세포중에 7개의 정상세포가 T/T를 가지고 있었고 3개가 암세포이고 G/T로 변이되었다고 하자. G/T의 비율은 3/17가 된다. 

프로그램은 정상 세포에서 발견한 돌연변이와 암세포에서 발견한 돌연변이를 고려하여 체세포 돌연변이로 구분하게 된다.

Whole genome sequencing을 사용하여 염색체 전체에서 수행하는 것이 가장 확실하지만, 비용과 시간의 문제로 인해 Whole exome sequencing으로도 많이 진행하고 있으며 RNA-seq으로도 진행할 수 있지만 신뢰도가 그리 높지는 않다. 신뢰도가 높지 않은 이유는 RNA는 기본적으로 불안전 할 뿐만 아니라 RNA editing등으로 서열이 쉽게 변하고 이것이 변이라고 분석되어 질 수 있기 때문이다.

Reference -

http://ib.bioninja.com.au/standard-level/topic-3-genetics/33-meiosis/somatic-vs-germline-mutatio.html