親緣關(guān)系矩陣
? ? ? ? 親緣關(guān)系矩陣被廣泛用于全基因組關(guān)聯(lián)分析(GWAS)和全基因組選擇(GS)中放祟,當(dāng)然主要是混合線性模型。在GWAS中親緣關(guān)系矩陣主要用于校正遺傳背景鹅髓,因?yàn)镚WAS一般采用單點(diǎn)掃描的模型近顷,這種模型極易形成假陽性鳞贷。而遺傳背景的本質(zhì)是除當(dāng)前掃描標(biāo)記外其馏,其他所有標(biāo)記的效應(yīng)(具體的我們后面再談)剧蹂。親緣關(guān)系矩陣可分為A陣和G陣,A陣是通過系譜進(jìn)行計(jì)算所得醒串,G陣是通過基因標(biāo)記計(jì)算所得执桌。G陣又被稱作realized relationship matrix(現(xiàn)實(shí)關(guān)系矩陣)。為什么G陣更realized呢芜赌?因?yàn)楦鶕?jù)系譜推導(dǎo)除的A陣是無法區(qū)分全同胞的仰挣,全同胞的關(guān)系是完全一樣的。但是我們知道即使是同卵雙胞胎也會(huì)有不同缠沈,而通過基因標(biāo)記就可以找到他們之間的不同的SNP膘壶。因此通過基因標(biāo)記計(jì)算的親緣關(guān)系矩陣更真實(shí)错蝴。
? ? ? ? 關(guān)于G矩陣被應(yīng)用最多的是Vanraden于2008在奶業(yè)科學(xué)上發(fā)表的文章。該文論述了幾種G矩陣的計(jì)算方法颓芭,但目前主流用的是前兩種顷锰。兩種方法的區(qū)別在于校正標(biāo)記期望方差的方法不同(具體我下次寫一篇文章專門說公式推導(dǎo))。`Vanraden`自己的模擬結(jié)果是第一種方法最好(使用所有標(biāo)記期望方差的均值來校正)亡问,但是楊劍老師主推的是第二種方法(每個(gè)標(biāo)記校正自己的期望方差)官紫。目前在GS中主要用的是第一種方法來計(jì)算kinship。其實(shí)兩種方法差別并不大州藕,因?yàn)闃?biāo)記本來就是0束世,1,2三種分型床玻,方差的區(qū)別也大不到哪去毁涉。
? ? ? ? 以下正式講解GCTA計(jì)算G矩陣的方法。
? ? ? ? 今天我們來講講如何利用GCTA來計(jì)算親緣關(guān)系矩陣笨枯。目前GCTA已經(jīng)支持windows了薪丁,但作者仍建議我們在linux下使用。安裝CCTA比較簡單馅精,直接下載解壓,運(yùn)行可執(zhí)行文件即可粱檀。最好軟鏈接到系統(tǒng)環(huán)境目錄下(如~/bin)洲敢,就可在任意目錄下輸入gcta64運(yùn)行。
? ? ? ? 進(jìn)入GCTA官網(wǎng)茄蚯,可以看到使用 ```--make-grm``` 就可以計(jì)算G矩陣:
? ? 使用```--make-grm-alg``` 可以選擇不同的算法压彭,因?yàn)闂顒蠋煴容^認(rèn)同第二種算法,所以默認(rèn)是第二種算法啦渗常。
這里我們使用最常用的第一種G矩陣計(jì)算方法(設(shè)置make-grm-alg 1即可)壮不。想用第二種算法的可以使用make-grm-alg 0。完整代碼如下:
gcta --bfile test --make-grm --make-grm-alg 1 --out kinship
bfile參數(shù)用于指定待運(yùn)行的plink二進(jìn)制文件皱碘,這里用的數(shù)據(jù)是GCTA自帶的測試數(shù)據(jù)询一。運(yùn)行結(jié)果如下:
? ? ? ? 運(yùn)行后,會(huì)得到如下四個(gè)文件:
? ? ? ? 1. test.grm.bin? 含G陣下三角元素癌椿,是二進(jìn)制文件
? ? ? ? 2. test.grm.N.bin 記錄計(jì)算G陣的SNP個(gè)數(shù)健蕊,是二進(jìn)制文件
? ? ? ? 3. test.grm.id 記錄個(gè)體的family號(hào)和id號(hào),即plink fam文件的前兩列
? ? ? ? 4. kinship.log 日志文件踢俄。
使用以上代碼生成的是二進(jìn)制文件的親緣關(guān)系矩陣缩功,不利于進(jìn)行其他分析。如果僅僅為了計(jì)算親緣關(guān)系矩陣都办,并需要用其他軟件進(jìn)一步分析嫡锌,推薦使用以下代碼:
gcta64 --bfile test --make-grm-gz --make-grm-alg 1 --out kinship
運(yùn)行后虑稼,會(huì)得到如下四個(gè)文件:
1. kinship.grm.gz? 個(gè)體間親緣關(guān)系的gz壓縮文件,解壓后可獲得四列數(shù)據(jù)結(jié)果势木,分別是 id1蛛倦, id2, 用于計(jì)算id1和id2親緣關(guān)系的SNP總數(shù)跟压, id1和id2的親緣關(guān)系系數(shù)胰蝠。
2. kinship.grm.id 同上
3.kinship.log 同上