SOLiD의 BAM파일의 경우에는 base quality을 ord('A') 를 한 후 -33 을 하면 됨. (0~40 까지의 범위)

QUAL: ASCII of base QUALity plus 33 (same as the quality string in the Sanger FASTQ format).
A base quality is the phred-scaled base error probability which equals 10 log10 Pr{base is wrong}.
This eld can be a `*' when quality is not stored. If not a `*', SEQ must not be a `*' and the
length of the quality string ought to equal the length of SEQ.


기본적으로 모두 옛날부터 시퀀서에서 적용되던 phred quality score의 개념을 따른다. 10 10%의 에러 확률, 20 1%의 에러 확률, 30은 0.1%의 에러 확률을 의미한다예를 들어 어떤 시퀀서가 99.99%의 정확도를 냈다고 한다면 그건 생산된 데이터(reads)의 대부분이 QV40 이상이었다는 의미가 된다장비마다 데이터를 생산하면서 각 메커니즘에 맞게 어떤 신호가 어떤 형식으로 나와서 그게 base call 또는 color call을 할 때 어느 정도의 정확성을 보이는지 미리 training 시켜서 얻은 경험(?)으로 나타낸다보통 다양한 생물종의 데이터를 준비하고 같은 기종이라도 여러 대에서 실험하면서 일종의 점수표를 만드는 것으로 안다따라서 개념은 같지만 서로 다른 기종의 QV를 그대로 비교하는 건 좀 위험하며기종에 따라 QV를 좀 더 좋게 보여주는 것이 있을 수도 있다시퀀싱을 한 후에 일차적인 평가를 하는데 중요한 단서이기는 하지만 실제 최종적인 서열의 정확도를 보여주는 것은 아니다참고로 다른 NGS들과 달리 SOLiD에서는 QV를 기반으로 한 필터링을 하지 않고 일단 모두 raw data로 생산한다복잡한 genome에서 QV가 특별하게 낮은 영역도 있을 수 있으므로, 그러한 곳에 대한 정보를 전부 잃기보다는 일단 분석 과정까지 가지고 간다는 의미가 있다. 
Posted by 옥탑방람보
,
자바에서 시스템 명령어나 외부 명령어를 시스템 상에서 사용해야할때
standard out 이나 standard error 가 지속적으로 발생하여
정상적인 프로세스 종료를 하지 못하고 프로세스가 좀비화 되는 경우가 있습니다.

이것은 자바가  standard out/error 를 직접 가져오게 되어서 발생되는 현상으로
일정 이상의 데이터가 쌓이게 되면 프로세스는 죽게 됩니다.

따라서 이때에는  standard out/error 을 파일로 떨구어 주어야 합니다.
시스템 상에서 실행되는 명령어에 직접  standard out/error 를 파일로 떨구는 부분을 추가해야 합니다.

예) bcftools -bcvg test.bam > standard.out 2> standard.error

 
Posted by 옥탑방람보
,
레퍼런스 순서와 bam 안의 염색체 순서가 일치해야함.
Posted by 옥탑방람보
,

1. download (설치된 파이썬 버젼과 일치하는 것)
   http://code.google.com/p/pygr/downloads/list
2.  egg 파일 풀기 (알집 등 이용)
3. 파이썬 라이브러리 밑에 위치 (import 가능한 위치)

 
Posted by 옥탑방람보
,
1. 다운로드 
2. 압축 풀기
3. 설치 디렉토리에 mv
4. configure
5. make
6. ~/.bashrc

#python
export PATH=/.../install/Python-2.6.7/:$PATH
export PYTHONPATH=/.../kimps/lib:/.../install/Python-2.6.7/Lib:/.../install/Python-2.6.7/Lib/site-packages 
Posted by 옥탑방람보
,
MD tag and cigar
10A5^AC6

REF:         ATCGTAGCTAATTTGGACATCGGT
READ:        ATCGTAGCTATTTTGG--ATCGGT
MD TAG:      10        A5   ^AC6
CIGAR:       16M             2D6M
READ:        atcGTAGCTATTTTGGATA..GGT (ATCGTAGCTATTTTGGATAAAGGT)
MD TAG:      17               C1TC3
CIGAR:       3S 16M             2N3M
READ:        ATCGTAGCTAATTTGGACATCGGT (ATCGTGGAGCTAATTTGGACATCGGT)
CIGAR:       5M   2I19M


MD TAG
The MD eld aims to achieve SNP/indel calling without looking at the reference. For example, a string `10A5^AC6' means from the leftmost reference base in the alignment, there are 10 matches followed by an A on the reference which is di erent from the aligned read base; the next 5 reference bases are matches followed by a 2bp deletion from the reference; the deleted sequence is AC; the last 6 bases are matches. The MD eld ought to match the CIGAR string.

CIGAR

M     alignment match (can be a sequence match or mismatch)
I     insertion to the reference
D     deletion from the reference
N     skipped region from the reference
S     soft clipping (clipped sequences present in SEQ)
H     hard clipping (clipped sequences NOT present in SEQ)
P     padding (silent deletion from padded reference)
=     sequence match
X     sequence mismatch

H can only be present as the rst and/or last operation.
S may only have H operations between them and the ends of the CIGAR string.
For mRNA-to-genome alignment, an N operation represents an intron. For other types of alignments, the interpretation of N is not de ned.
Sum of lengths of the M/I/S/=/X operations ought to equal the length of SEQ.
Posted by 옥탑방람보
,
import os

a = '../work/'
print os.path.abspath(a) 
Posted by 옥탑방람보
,
Exception in thread "main" java.lang.RuntimeException: SAM validation error: ERROR: Record 62094965, Read name ILLUMINA-A16956_100211:4:14:19403:10471#0, MAPQ should be 0 for unmapped read

이러한 에러 발생시,

VALIDATION_STRINGENCY=
LENIENT

옵션을 추가해준다.


This is a common problem, and you'll run into it with all the Picard suite of tools.

There's a setting that goes something like VALIDATION_STRINGENCY, and if you set it to LENIENT, it will complain about those reads, but it won't stop on them.

This happens with bwa, because it concatenates reference sequence, which leads to slightly odd things happening when a read aligns over the overlap. So this might be the source of your problem. Regardless of what's causing it, you can examine the problem reads, and cut them out, or change the stringency to let them go through.

 
Posted by 옥탑방람보
,
java -jar MarkDuplicates.jar INPUT=test.bam OUTPUT=test.marked.bam METRICS_FILE=test.txt TMP_DIR=. ASSUME_SORTED=true
Posted by 옥탑방람보
,
import inspec
print inspect.getfile( inspect.currentframe() )

import os
print os.path.abspath( __file__ )
 
import sys
sys._getframe().f_code.co_filename 
Posted by 옥탑방람보
,