[Linux] Human Genome 데이터 예시

* 이번 예시에서 사용한 sequence data는 https://hgdownload.cse.ucsc.edu/goldenPath/hg38/chromosomes/에서 찾을 수 있다. 해당 사이트에 업로드되어 있는 데이터는 hg38(GRCh38 Genome Reference Consortium Human Reference 38)의 sequence data이다.

Index of /goldenPath/hg38/chromosomes

 

1. 파일 다운 받아 압축 풀고 새 이름의 파일로 저장

$ curl http://hgdownload.cse.ucsc.edu/goldenPath/hg38/chromosomes/chr22.fa.gz | gunzip -c > chr22.fa

 

- 22번 염색체 파일 사용

 

 

2. 파일 내용 보기

$ cat chr22.fa | head
>chr22
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN

 

- N = unknown base

 

 

3. N의 개수

$ cat chr22.fa | grep -o N | head
N
N
N
N
N
N
N
N
N
N

$ cat chr22.fa | grep -o N | wc -l
11658691

 

→ 11,658,691 bases are unknown.

(-o: only matching)

 

 

4. 22번 염색체 전체 base 수

- 3가지 방법으로 추출 가능

 

(1) seqkit

$ seqkit stat chr22.fa
file      format  type  num_seqs     sum_len     min_len     avg_len     max_len
chr22.fa  FASTA   DNA          1  50,818,468  50,818,468  50,818,468  50,818,468

 

→ 50,818,468 bases

 

(2) infoseq

$ infoseq chr22.fa
Display basic information about sequences
USA                      Database  Name           Accession      Type Length %GC    Organism            Description
fasta::chr22.fa:chr22    -              chr22          -              N    50818468 47.69

 

→ 50,818,468 bases

 

(3) bioawk

$ cat chr22.fa | bioawk -c fastx '{ print length($seq) }'
50818468

 

→ 50,818,468 bases

 

=> 22번 염색체의 약 23% 'unknown' (sequenced, but not yet localized = missing)

---

Reference

The Biostar Handbook: 2nd Edition - István Albert