Be great

Bowtie2에서 mismatch 허용하기.

bowtie2에서는 bowtie와 차이점은 단순히 read length에 따른 최적화나, gap 허용 외에도 

bowtie에서는 최대 3개까지의 mismatch만을 허용하는데 반해 bowtie2에서는 mismatch 또는 indel의 각각의 페널티 점수를 입력하여 read length의 일정 비율만큼 mismatch 허용할 수 있다는 차이점이 있다.

bowtie2에 대한 자세한 옵션은 http://bowtie-bio.sourceforge.net/bowtie2/manual.shtml 를 참고하면 된다. 

이 포스팅에서는 mismatch에 대해서만 다루고자 한다.

For an alignment to be considered "valid" (i.e. "good enough") by Bowtie 2, it must have an alignment score no less than the minimum score threshold. The threshold is configurable and is expressed as a function of the read length. In end-to-end alignment mode, the default minimum score threshold is -0.6 + -0.6 * L, where L is the read length. In local alignment mode, the default minimum score threshold is 20 + 8.0 * ln(L), where L is the read length. This can be configured with the --score-min option. For details on how to set options like --score-min that correspond to functions, see the section on setting function options.

Scoring options

--ma <int>

--mp MX,MN

--np <int>

--rdg <int1>,<int2>

--rfg <int1>,<int2>

--score-min <func>

bowtie2 manual 홈페이지에서 가져온 score 계산 방법이다.

우선 위의 end-to-end score threshold는 -0.6 + -0.6 * readlength 라고 적혀있다.

read length를 100으로 놓으면 비율로 생각하기 쉬우니 그렇게 계산해보면 -60.6 보다 점수가 낮으면 align하지 않겠다는 뜻으로 해석 가능하다.

각 base가 match, mismatch, gap일 때의 score를 살펴보면 마찬가지로 end-to-end일때 mismatch score는 read의 base quality에 따라 점수가 다르게 측정 되는 것으로 보인다. (maximum mismatch penalty와 minimum mismatch penalty가 존재)

preprocess과정에서 low quality base를 자르기도 하고 요즘 sequencing을 하면 대체적으로 high quality read가 많으니 일단 quality가 좋다고 가정할 때의 score인 6으로 계산한다.

(read의 quality가 나쁘면 penalty는 작다. 예를들어 read의 base가 A이며 quliaty score가 낮다면 실제로 이 base는 T이고 mismatch가 아니라 match일 수도 있기 때문이다.)

default 설정의 경우

- mismatch가 10번 생기면 -60점이니 read length의 10%라고 보면 된다.

- gap같은 경우는 open과 extend가 각각 다르게 적용되니 gap이 하나만 있다면 18개 까지 생길 수 있다.

사용할 때는 dafult 값인 mp와 rdg등은 가급적이면 안건드리고 --score-min L,-0.6,-0.3 등으로 바꿔서 (5%의 mismatch를 허용) 해보는 것을 추천 하지만 data 특성에 따라서 mismatch적절한 값을 주고 사용해야 한다.

Python에서 rpy 모듈을 사용하다가 

cannot mkdir R_TempDir 에러메세지 발생했다.

/tmp/ 폴더를 확인해 봐야 한다.

권한이 없거나, tmp 폴더에 설정한 용량을 채우면 더 진행되지 않는다.

확인필요.

miRNA의 명명법은 규칙이 있다.

이 내용을 mirbase에서 소개하고 있으며 요약하고자 한다.

원문은 http://www.mirbase.org/help/nomenclature.shtml 참조 하면 된다.

요약하자면

1. hsa-mir-121이라는 miRNA가 있을 때 학명-mir-숫자 의 형식을 따른다.

2. 숫자는 발견된 순서로서 마지막으로 121이라는 miRNA가 있었다면 이후에 발견되는 miRNA는 122부터 시작한다.

3. genome 상의 다른 영역에서 같은 mature miRNA를 가지는 precursor miRNA가 있다면 이름은 hsa-mir-121-1, hsa-mir-121-2를 가진다.

4. genome 상의 다른 영역에서 유사한 mature miRNA를 가지고 있다면 hsa-mir-121a, hsa-mir-121b를 가진다.

5. mature miRNA는 위치에 따라 precursor miRNA 이름 뒤에 -5p , -3p를 가진다. ex) hsa-mir-121-5p, hsa-mir-121-3p

6. 이 규칙은 항상 적용되는 것은 아니며 예외가 있을 수 있다.

마지막으로 이름은 아주 일부의 정보만 가지고 있기 때문에 miRNA의 정확한 정보를 알기 위해서는 database를 검색해야지 이름에 의존하면 안된다고 하며 설명을 마무리 하고 있다.

리눅스에서 문서형식으로 된 데이터 파일들을 처리하다 보면 문서 편집 기능이 필요할 때가 많다.

일일이 스크립트를 짜서 하기엔 너무 번거로우니 간단하게 한 줄로 처리하고 싶을 때 자주 사용하게 되는 명령어들을 정리해 보았다.

• grep

grep은 사용하는 법도 간단하지만, 여러 가지 조건을 걸고 보기에는 무리가 있어서 다른 유닉스 명령어와 혼합하여 사용하는 경우가 많다.

• sed

file을 열지 않고 data를 편집가능하게 해주는 명령어 이다. 

stdout으로 바로 출력되기 때문에 원본 파일은 변화가 없으며 파일로 새로 쓰고 싶으면 stdout을 바로 파일로 만들어주면 된다.

예시.

패턴 뒤에 나오는 d, p, g 등이 어떤 기능을 수행할 지 결정하며 그 외의 특별히 기억해야 할 명령어는 없다.

• awk

익숙해지면 따로 스크립트를 작성하지 않아도 간단하고 빠르게 문서 편집이 가능하다.

특히나 리눅스에서는 문서 편집 명령어들을 쓸 일이 많기 때문에 쓰다보면 자연스럽게 익숙해질 것이다.

Windows 10에서 Redstone 패치부터 리눅스에서 사용하는 오픈소스 명령어 툴인 Bash shell을 지원합니다.

기존의 윈도우 환경에서 리눅스를 사용하기 위해서는 가상 머신이나 컨테이너 등을 설치하여야 했으나 이제는 그냥 Bash shell 자체를 윈도우에서 수행할 수 있도록 바꿨습니다.

윈도우 데스크탑과 리눅스 서버를 같이 사용하면서 맥에서는 할 수 있었지만 윈도우에서는 쉽게 하지 못했던 아쉬운 기능들이 있었는데 이번 패치 이후로는 이제는 훨씬 편하게 할 수 있을 거 같아서 매우 기분이 좋습니다.

윈도우에서 bash 사용을 위해서는 간단한 설치 과정이 필요합니다. 

설치 방법

• 개발자 모드 설정

설정 -> 업데이트 및 복구 -> 개발자용 -> 개발자 모드에 체크

• windows 기능 켜기/끄기

검색창에서 windows 기능 켜기/끄기 검색 후 실행

windows 기능 창에서 Linux용 windows 하위 시스템(베타) 에 체크

검색창에서 명령 프롬프트 실행 후 bash 입력

설치과정이 끝나고 나면 이후부터는 그냥 검색창에 bash라고 쳐도 실행되는 걸 보실 수 있습니다.

이 기능을 어떻게 응용할지는 이후 포스팅에서 계속 설명하겠습니다.