第二種方法按照-1, 0, 1, 即-1是AA, 表示major基因型, 0表示雜合, 1表示aa(minor).
下面我們模擬一下數(shù)據(jù), 看一下兩者計算的區(qū)別.
dat = data.frame(snpid = paste0("snp",1:7),type = c("GG","GG","AG","GG","AG","GG","AA"))dat| snpid | type |
|---|---|
| snp1 | GG |
| snp2 | GG |
| snp3 | AG |
| snp4 | GG |
| snp5 | AG |
| snp6 | GG |
| snp7 | AA |
GG個數(shù): 4
AG個數(shù): 2
AA個數(shù): 1
計算方法, 所有G的個數(shù), 除以總配子的個數(shù), 即為G的基因頻率
p_g = (4*2+2*1)/(7*2);p_g0.714285714285714
p_a = (2*1+1*2)/(7*2);p_a0.285714285714286
dat$type012 = 0
dat[dat$type=="GG",]$type012 =0
dat[dat$type %in% c("AG","GA"),]$type012 =1
dat[dat$type=="AA",]$type012 =2dat| snpid | type | type012 |
|---|---|---|
| snp1 | GG | 0 |
| snp2 | GG | 0 |
| snp3 | AG | 1 |
| snp4 | GG | 0 |
| snp5 | AG | 1 |
| snp6 | GG | 0 |
| snp7 | AA | 2 |
計算方法:
將type012所有列之和相加, 為minor的基因個數(shù), 除以2*n, 即為minor的基因頻率
sum(dat$type012)/(length(dat$type012)*2)0.285714285714286
dat$type_101 = dat$type012 -1dat| snpid | type | type012 | type_101 |
|---|---|---|---|
| snp1 | GG | 0 | -1 |
| snp2 | GG | 0 | -1 |
| snp3 | AG | 1 | 0 |
| snp4 | GG | 0 | -1 |
| snp5 | AG | 1 | 0 |
| snp6 | GG | 0 | -1 |
| snp7 | AA | 2 | 1 |
計算方法:
分子為: 列之和 + 行數(shù)
分母為: 行數(shù) * 2
(sum(dat$type_101) + length(dat$type_101))/(2*length(dat$type_101))0.285714285714286
結(jié)論:
1, 如果轉(zhuǎn)化為0, 1, 2, 那么他們之和即為minor的個數(shù), 除以2n即為minor的頻率. 因為major的值為0, 1為雜合, 同時表示有一個minor, 2為純合, 表示2個minor, 計算比較方便.
2, 如果基因型轉(zhuǎn)化為-1, 0, 1的話, 那么將其進(jìn)行求和, -1和1抵消. 推斷過程如下:
假設(shè)A的頻率為(1-q), a的頻率為q, n為個數(shù).
AA的個數(shù)為x
Aa的個數(shù)為y
aa的個數(shù)為z
x + y + z =n
q = (2*z + y)/(2*n)
-1*x +0*y + 1*z = z-x
2*z +y = 2*z + (n - x -z) = n + z - x #分子
因此:
q = (n + z -x)/(2*n)
這種方法, 計算量比較小, 因為-1 和1 進(jìn)行了抵消.
這種方法不太直觀, 需要用筆畫一下即可理解.
所以, 我們計算G矩陣時, 將SNP分型轉(zhuǎn)化為-1, 0, 1的形式.
如何利用系譜計算近交系數(shù)和親緣關(guān)系系數(shù)
