R(프로그래밍 언어) qtl2 #
Find similar titles

Structured data

Category: Programming

R/qtl2 #

R/qtl2는 Quantitative Trait Loci(QTL) mapping을 하기 위한 대화형 소프트웨어다. R/qtl의 확장 버전으로, R과 C++ 프로그래밍 언어로 작성된 오픈소스 소프트웨어다. Collaborative Cross and Diversity Outbred mice, heterogeneous stocks, MAGIC plant populations와 같은 여러 조상에서 유래된 다중 부모 집단을 포함한다. 유전자 발현과 단백체학을 포함하는 고차원적 분자 표현형과 현대 고밀도 유전자형을 다루기 위해 설계되었으며 QTL mapping에 필요한 모든 기초적인 기능을 제공한다.

데이터 형식 #

R/qtl과는 다른 포맷형식을 사용하며, control file을 포함한 몇 개의 파일이 필요하다. control file을 제외한 모든 file은 csv 확장자를 가진다.

control file : YAML 또는 JSON 확장자를 가지며, 기본적으로 genotype, phenotype, gmap, pmap 파일명을 포함한다.

crosstype: f2
geno: iron_geno.csv
pheno: iron_pheno.csv
phenocovar: iron_phenocovar.csv
covar: iron_covar.csv
gmap: iron_gmap.csv
alleles:
- S
- B
genotypes:
  SS: 1
  SB: 2
  BB: 3
sex:
  covar: sex
  f: female
  m: male
cross_info:
  covar: cross_direction
 (SxB)x(SxB): 0
 (BxS)x(BxS): 1
x_chr: X
na.strings:
- '-'
- NA

genotype file : 각 marker에 대한 genotype 정보
genetic map file : 각 marker에 대한 genetic position 정보(단위 : cM)
physical map file : 각 marker에 대한 physical position 정보(단위 : bp)
phenotype file : 양적형질에 대한 정보
covariation file : 공변이 정보 (성별 등)
pheno covariation file : 형질에 대한 공변이 정보(시간 등)

사용법 요약 #

1. 데이터 로드

control 파일을 불러온다.

library(qtl2)
grav2 <- read_cross2("~/my_data/iron.yaml")

2. 유전자형 확률 계산

QTL을 계산하기 위해, marker들 사이에 pseudomarker를 삽입하고 이에 대한 유전자형 확률을 계산한다.

map <- insert_pseudomarkers(iron$gmap, step=1)
pr <- calc_genoprob(iron, map, error_prob=0.002)

3. kinship matrix 계산

선형혼합모형(linear mixed model)을 이용하여 genome scan을 하기 위해 개체들간의 관계를 설명하는 kinship matrix를 계산한다. type="loco" 옵션을 넣으면 LOCO(leave one chromosome out) 방법을 사용한다.

kinship <- calc_kinship(pr)
kinship_loco <- calc_kinship(pr, type="loco")

4. genome scan 수행

선형혼합모형(linear mixed model)을 이용하여 genome scan을 수행한다. kinship matrix를 사용하지 않을 경우, Haley-Knott regression 방법을 사용하며, LOCO 방식으로 구한 kinship matrix를 넣으면 LOCO 방식으로 genome scan을 수행한다. genome scan을 통해 LOD score를 얻을 수 있으며, LOD score를 이용하여 그래프를 그릴 수 있다.

out <- scan1(pr, iron$pheno, Xcovar=Xcovar) # H-K
out_pg <- scan1(pr, iron$pheno, kinship, Xcovar=Xcovar) # LMM
out_loco <- scan1(pr, iron$pheno, kinship_loco, Xcovar=Xcovar) # LOCO

color <- c("slateblue", "violetred", "green3")
par(mar=c(4.1, 4.1, 1.6, 1.1))
ymx <- max(maxlod(out), maxlod(out_pg), maxlod(out_pg_loco))
for(i in 1:2) {
    plot(out, map, lodcolumn=i, col=color[1],
         main=colnames(iron$pheno)[i], ylim=c(0, ymx*1.02))
    plot(out_pg, map, lodcolumn=i, col=color[2], add=TRUE)
    plot(out_pg_loco, map, lodcolumn=i, col=color[3], add=TRUE, lty=2)
    legend("topleft", lwd=2, col=color, c("H-K", "LMM", "LOCO"),
           bg="gray90", lty=c(1,1,2))}

5. permutation test 수행

genome scan 결과에 대한 통계적 유의성을 검증하기 위해 permutation test를 수행한다. r/qtl의 scanone()과 동일하며, threshold보다 큰 값을 유의하다고 판단한다.

operm <- scan1perm(pr, iron$pheno, Xcovar=Xcovar, n_perm=1000)
summary(operm)

## LOD thresholds (1000 permutations)
##      liver spleen
## 0.05  3.46   3.46

6. QTL effect 추정

4번 과정에서는 오직 LOD score만 얻을 수 있기 때문에, QTL effect를 따로 추정한다. scan1coef() 함수는 linear mixed model을 사용하며, scan1blup() 함수는 best linear unbiased predictor(BLUPs)를 사용한다. 두 함수 모두 REML(restricted maximum likelihood) 방법을 default로 적용한다.

c2eff <- scan1coef(pr[,"2"], iron$pheno[,"liver"])

par(mar=c(4.1, 4.1, 1.1, 2.6), las=1)
col <- c("slateblue", "violetred", "green3")
plot(c2eff, map["2"], columns=1:3, col=col)
last_coef <- unclass(c2eff)[nrow(c2eff),]
for(i in seq(along=last_coef)){
    axis(side=4, at=last_coef[i], names(last_coef)[i],
         tick=FALSE, col.axis=col[i])}

single QTL position에서 genotype에 대한 phenotype plot을 그릴 수 있다. 아래 예제는 liver phenotype(28.6cM)에서의 염색체 2번 QTL에 대한 genotype을 나타낸 것이다.

g <- maxmarg(pr, map, chr=2, pos=28.6, return_char=TRUE)

par(mar=c(4.1, 4.1, 0.6, 0.6))
plot_pxg(g, iron$pheno[,"liver"], ylab="Liver phenotype")

더 나아가, SNP database를 활용하면 SNP association 분석도 진행할 수 있다.

Reference #

R/qtl2 user guide : https://kbroman.org/qtl2/assets/vignettes/user_guide.html
R/qtl2 sample data : https://kbroman.org/qtl2/pages/sampledata.html