MultiQC 설치 및 실행하기

로그파일의 이름과 무관하게 형식만 맞으면 다 읽는 것으로 보이며 그럼에도 불구하고 실행시간은 오래걸리지 않는다. 

로그파일들의 크기가 대체로 작긴 하지만 폴더내에 있는 파일 중에 매치가되지 않으면 바로바로 점프하는 듯 하다.

설치 방법은 가장 간단한 방법은 PyPI을 사용하는 방법이다.

프로그램이 Python 기반이기 때문에 Python을 설치하고 PyPI를 사용하도록 하자.

2017/08/16 - [programming language/python] - Python 설치 및 실행하기

pip가 설치되어 있다면 아래처럼 입력하고

권한문제가 생긴다면 --user 옵션을 넣도록 하자.

--user 옵션을 사용했다면 실행 파일은 ~/.local/bin/multiqc 에 생겼을 것이다.

설치가 정상적으로 끝났다면 프로젝트 폴더에가서

명령을 수행하면 multiqc_report.html 파일이 생성된다.

html파일을 열면 프로젝트 폴더 안에 있는 대부분의 로그파일들이 올라와있는 것을 확인할 수 있다.

당연하게도 모든 중간작업마다 로그파일을 생성했어야만 정상적으로 파일을 읽고 작동한다.

table이나 plot 복사를 하면 report를 작성하는데 매우 유용하게 사용할 수 있을 것이다.

Reference-

http://multiqc.info/docs/#using-multiqc

RPKM, FPKM and TPM의 정의

RPKM (Reads Per Kilobase of transcript per Million mapped reads)

FPKM (Fragments Per Kilobase of transcripts per Million mapped reads)

TPM (Transcripts Per Million)

Reference -

http://www.incodom.kr/Expression_profiling

https://www.rna-seqblog.com/rpkm-fpkm-and-tpm-clearly-explained/

Wagner GP, Kin K, Lynch VJ. Measurement of mRNA abundance using RNA-seq data: RPKM measure is inconsistent among samples. Theory in biosciences. 2012 Dec 1;131(4):281-5.

DESeq2 에서 multiple condition 수행하기

쭉 진행하면 나오는 결과는 test와 brain간의 비교이다. DESeq2는 default로 가장 처음과 마지막에 있는 컨디션을 비교하기 때문이다.

이를 아래처럼 수정하면 test와 UHR간의 비교 결과를 보여준다.

2, 3번째 값을 수정하면 그에 맞게 바뀌는데 당연히 이미 설정해놓은 컨디션 값을 지정해주어야하며

순서를 바꾸면 log2FoldChange값이 +에서 -로 바뀌게 된다.

log2FC값 외에 변경되는 값은 없다.

리눅스에서 프록시 설정하기

Centos 6에서 R 설치를 위한 라이브러리 설치

zlib

bzip

liblzma

pcre

curl

모든 dependency가 설치되었으면 R의 configure를 해보도록 하자.

중간에 warning이 뜨긴 했지만 일단 무시하고 진행했으며 R일 실행되는 것 까지 확인하였다.

Reference -

https://unix.stackexchange.com/questions/343452/how-to-install-r-3-3-1-in-my-own-directory

Arguments in R

argument를 더 복잡하게 쓰려면 optparse라는 라이브러리를 써도 되지만 간단하게 정리하고 싶다면 위와 같이 작성할 수 있다.

추가로 argument를 입력하지 않았을 때 간단한 설명을 넣고 싶다면 아래와 같이 하면 된다.

argument가 하나도 들어오지 않았다면 ERROR 메세지 뒤에 정해놓은 문자열을 출력하고 자동 종료된다.

Remove duplicates

source -

F-measure

precision (정확도)

recall (재현율)

2018/07/11 - [bioinformatics] - 민감도와 특이도

두 값을 모두 구했다면 F값을 계산하면 된다.

위의 값대로 계산하면 precision과 recall의 조화평균을 구할 수 있다.

precision과 recall등은 파라미터 등에 의해 조절될 수 있기 때문에 ROC 커브를 그려서 프로그램의 전반적인 성능을 테스트 할 수도 있다. 

이 부분에 대해서는 추후에 다루도록 하겠다.

source -

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

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

Gene ontology analysis - DAVID

DAVID는 Database for Annotation, Visualization and Integrated Discovery의 약자로 유전자 리스트를 입력으로 받아 각 유전자의 기능은 해석해주는 웹 제공을 기반으로하는 무료 툴 이다.

주로 유전자 기능 분류를 하거나 기능을 모를때 주석을 달기 위해 사용되며 이를 위해 현재 공개되어 있는 주요한 데이터베이스의 정보를 대부분 가져와 직접시켜 DAVID만의 데이터 베이스를 만들고 있다.

주요 기능으로는 유전자 리스트가 주어졌을 때

- 특정 기능에 대한 유전자들이 많이 포함되었는지

- 비슷한 기능을 가진 유전자들의 그룹화

- BioCarta & KEGG pathway map과의 가시화된 연결

- 2-D로 유전자와 특정 묶음간의 연관성

- 유전자와 상호작용하는 단백질 리스트

- 유전자의 질병간의 연관성 리스트

- 단백질의 기능적 도메인과 모티프

- 관련 문헌

- 유전자 ID를 다른 ID로 변환 ex) ensembl id에서 refseq id로

등등이 존재한다.

DAVID를 실제로 사용하기 위해서는 gene id list가 필요하다. gene symbol인지 특정한 데이터베이스에서 사용하는 ID인지는 ID mapping 과정을 통해 변환하는 과정이 있기 때문에 크게 중요하지 않다.

step1의 A에 유전자 ID를 직접 넣거나 한 줄에 하나씩 입력된 파일을 B에 넣고 step2의 identifier는 어떤 종류의 ID를 사용하고 있는지 넣은 후 step3의 gene list를 체크 후 submit 하면 된다.

당장 파일이 없을 땐 demolist를 눌러서 진행한다.

입력한 gene id를 기반으로 어떤 종에서 찾고싶은 것인지 선택하여야 한다. 만약 넣어준 gene list의 매칭되는 종이 있다면 위와같은 화면이 나오겠지만 없다면 mapping과정을 진행하는 창이 뜰 것이다. mapping 과정은 gene id를 통일하는 과정이다. mapping이 진행되고 나면 위와 같은 화면이 나온다.

Use 버튼을 눌러 분석을 실행했다면 위와같이 특정한 분류대로 어떤 유전자가 많이 포함되어있는지를 보여준다. 

가장 밑의 Functional Annotation Clustering은 누르면 모든 cluster에 대해서 어떤 유전자가 어떻게 묶였으며 각 p-value는 어떻게 되는지 확인할 수 있다.

위와 같이 그룹화된것을 확인하고 p-value, Benjamini, FDR 값을 토대로 특정 값 이하의 그룹을 significant하다고 정의내린 후 결과를 정리하면 된다.

통계 방법에 대한 정의는 아래 포스트를 참조하면 된다.

2018/07/11 - [bioinformatics] - Multiple Comparsion Problem

DAVID에서 제공하는 FDR값은 정상적이지 않다. 제일 마지막 줄에 있는 그룹만 봐도 FDR값이 3.3인데 FDR음 0에서 1 사이의 값을 가져야 한다.

Benjamini 값을 토대로 cutoff를 정하면 될 것이다.

source -

https://www.ncbi.nlm.nih.gov/pmc/articles/PMC2375021/

https://david.ncifcrf.gov/

_PAR_Y in Genecode annotation

htseq-count를 돌렸는데 genecode annotation을 사용하였다.

결과로 나온 gene_id 중에 _PAR_Y라는 태그가 붙어서 나오는 것을 확인하였다.

특이하여 genecode 홈페이지에 직접 들어가서 무슨 의미인지 확인해보니 염색체 Y의 "pseudoautosomal region" (PAR)에 존재하는 유전자라는 뜻으로 염색체 X와 Y에 동일한 서열부분에 존재하는것을 의미한다. 

실제로 annotation file을 열어보면 염색체 X에는 "ENSG00000124333.15_2"라는 gene_id가 염색체 Y에는 "ENSG00000124333.15_2_PAR_Y"라는 gene_id가 존재하는것을 확인할 수 있었다.

source -

https://www.gencodegenes.org/faq.html