Kchip의 AX-ID -> RS-ID로 변경

 

SNP의 명명은 어떤 회사의 chip을 사용했는지에 따라 다릅니다.

 

Kchip은 axiom 사로 알고 있고, 그래서 SNP의 ID가 "AX-" 로 시작되게 됩니다.

 

이것을 일반적으로 dbSNP에서 사용하는 rsID로 변경해야 합니다.

 

보통은 chip 홈페이지에 들어가면, 관련된 csv 파일등이 있는데

 

Kchip은 아무리 뒤져도 나오지 않아 고생을 했습니다.

 

홈페이지 보고 메일을 보내도 답장도 없구요 ^^

 

 

우리는 ANNOVAR를 이용할 겁니다.

 

ANNOVAR에 대해서는 구글 검색해서 찾아보시면 되고

 

대부분 gene, SNP 등을 annotation 할 때 가장 많이 사용합니다.

 

사용하기 전에 개발자에게 사용하고 싶다고 등록을 하면 금방 답이 오고 다운받으실 수 있습니다.

 

굉장히 좋은 툴이니 꼭 받아둡시다.

 

 

# 우선 Kchip의 BIM 파일을 불러옵니다

fread("~/driveD/Gene_data/KoGES/Kchip_72292/Kchip.bim", header = F) -> Kchip

 

GWAS 데이터는 보통 binary인 경우 3가지 파일이 있습니다.

 

bed, bim, fam. 이에 대해 얼추 알고 계실 거라 생각(?)하고..

 

모르시면 plink 홈페이지에 자세히 나와 있어서 참조하시거나

 

다른 블로그 참조하시면 되겠습니다

 

이에 대해서는 다른 블로그에서 자세히 다룬 편입니다.

 

bim 파일에 chr, position, ref, alt 등이 포함되어 있는데

 

# 파일 정리

Kchip %>% 
  dplyr::rename(
    CHR = V1,
    AXID = V2,
    BP1 = V4,
    A1 = V5,
    A2 = V6
  ) %>% 
  mutate(
    BP2 = BP1,
    ID = paste(CHR, BP1, sep = "_")
    ) %>% 
  filter(CHR %in% seq(1:22)) %>% 
  dplyr::select(-V3) -> Kchip

 

이렇게 적당히 알기 쉽게 정리를 하고

 

# 이렇게 ANNOVAR를 할 수 있게 저장합니다.

write_tsv(Kchip %>% dplyr::select(CHR, BP1, BP2, A2, A1), 
          "~/driveD/work/rsid_update/Kchip.tsv", col_names = F)

 

ANNOVAR 로 저장하기 위해서는

CHR, BP, BP, REF, ALT 이렇게 5개 칼럼이 필요합니다.

 

position이 2번 반복되어야 합니다.

 

이렇게 저장한 뒤에

 

# ANNOVAR 시행.

# ANNOVAR는 linux command line에서 시행해야 합니다.

annotate_variation.pl \
Kchip.tsv \
~/driveD/ref_gene/ANNOVAR_dbSNP/humandb \
-filter -build hg19 -dbtype avsnp147 \
-out Kchip

 

humandb 받는 것은 ANNOVAR 홈페이지 tutorial에 자세히 써있습니다.

저희는 hg19로 avsnp147 (dbSNP을 개발자가 사용하기 편하게 조금 수정한 버전)

을 기준으로, chromosome, positioon, A1 (alternative, effect allele 의미), 

A2 (reference, non-effect allele을 의미)에 맞춰 자동으로 rsID를 찾아줍니다

 

이렇게 하시면 _dropped 파일과 _filtered 파일이 생기고, 

invalid_input 이라는 파일이 생깁니다

invalid_input에 있는 SNP들은 dbSNP에 뒤져봐도 맞는 것이 거의 없어보입니다.

 

dropped는 matching이 되는 SNP들이고

filtered는 matching이 안되어 filter된 SNP들 입니다.

보시면 A1 A2가 뒤바뀐 경우가 많습니다.

 

그래서 

 

# filter된 SNP을 불러와서 정리하고.

fread("~/driveD/work/rsid_update/Kchip.hg19_avsnp147_filtered") -> filtered

 

 

# A1, A2를 바꿔줍니다.  (즉, 이전에 A1이었던 SNP은 A2로, A2는 A1으로)
filtered %>% 
  dplyr::rename(
    CHR = V1,
    BP1 = V2,
    BP2 = V3,
    A2 = V4,
    A1 = V5
  ) %>% 
  arrange(CHR, BP1) %>% 
  mutate(ID = paste(CHR, BP1, sep = "_")) -> filtered

 

# 저장하고 다시 ANNOVAR를 합니다.

write_tsv(Kchip %>% filter(ID %in% filtered$ID) %>%  dplyr::select(CHR, BP1, BP2, A1, A2), 
          "~/driveD/work/rsid_update/Kchip_unmatched.tsv", col_names = F)

 

# linux command line에서 ANNOVAR 시행.

annotate_variation.pl \
Kchip_unmatched.tsv \
~/driveD/ref_gene/ANNOVAR_dbSNP/humandb \
-filter -build hg19 -dbtype avsnp147 \
-out Kchip_unmatched

 

 

# 그리고 matching된 SNP들을 불러옵니다.

fread("~/driveD/work/rsid_update/Kchip.hg19_avsnp147_dropped") -> dropped

fread("~/driveD/work/rsid_update/Kchip_unmatched.hg19_avsnp147_dropped") -> rev_dropped

 

 

# dropped는 첫 번째 ANNOVAR 결과이고

# rev_dropped는 unmatched된 것들로 다시 ANNOVAR해서 맞춘 것

# 파일을 하나로 합치고 저장합니다.

rbind(
  dropped %>% 
    dplyr::select(V3, V2, V4) %>% 
    dplyr::rename(
      RSID = V2,
      CHR = V3,
      BP = V4
    ),
  rev_dropped %>% 
    dplyr::select(V3, V2, V4) %>% 
    dplyr::rename(
      RSID = V2,
      CHR = V3,
      BP = V4
    )
  ) %>% 
  arrange(CHR, BP) %>% 
  mutate(ID = paste(CHR, BP, sep = "_"))-> rsid_form

 

# Kchip.bim 데이터와 합쳐서 axid, rsid만 뽑아내고
left_join(
  Kchip,
  rsid_form %>% dplyr::select(RSID, ID),
  by = "ID") -> rsid_form

 

# 저장합니다. write_delim 쓴 이유는, 구분이 space 1개라서...

write_delim(rsid_form %>% dplyr::select(AXID, RSID),
            "~/driveD/work/rsid_update/snp_update.txt", col_names = F)

 

>> less snp_update.txt

AX-37361829 rs117957074
AX-83277131 rs12095200
AX-64175821 rs117997018
AX-32416607 rs4422948
AX-32429231 rs78238606

 

이렇게 space 1개 txt 파일로 만든 후에

 

# 아래와 같이 plink로 SNP name을 update 해줍니다.

plink --bfile ~/driveD/Gene_data/KoGES/Kchip_72292/Kchip \
--update-name snp_update.txt \
--make-bed --out Kchip_RSID

 

 

이제 이 파일로 QC와 imputation 등을 하시면 됩니다

아니면 QC및 phenotype, individuals 정리 등을 먼저 하신 후에

rsid update를 하시고, 그 다음에 imputation 하셔도 됩니다

 

imputation 위해 RSID를 쓰기 때문에 필요한 과정입니다.