plink群体遗传学GWAS

8. GWAS:亲缘关系——TASSEL&GCTA

2021-08-24  本文已影响0人  Wei_Sun

• 指在非家系群体或系谱不明确的群体中两特定材料之间的遗传相似度与任意材料之间的遗传相似度的相对值;

• 材料间不平衡的血缘关系是导致标记出现非连锁相关的另一个重要原因,小家系的存在会使关联分析结果出现假阳性。为了避免这种情况的产生,往往会把亲缘关系矩阵作为随机效应协变量矩阵(K矩阵)加入GWAS模型;

• 常用软件:TASSEL、GCTA、LDAK。

1.TASSEL

1.1下载及安装

1.1.1 下载地址

https://tassel.bitbucket.io/
这里下载Linux版本TASSEL 5.0

1.1.2 安装

cd /gss1/home/fzhang/software
sh TASSEL_5_unix.sh
# 一直 enter 填写安装路径

1.2 TASSEL计算亲缘结构

1.2.1 计算IBS亲缘关系矩阵

设置内存:-Xmx10g -Xms512m
排序:-SortGenotypeFilePlugin
输入文件全名:-inputFile
输出文件前缀:-outputFile

cd /gss1/home/fzhang/sunwei/TASSEL5
#对vcf文件进行排序
perl run_pipeline.pl -Xmx10g -Xms512m -SortGenotypeFilePlugin -inputFile root.id.vcf -outputFile Troot -fileType VCF
#排序后文件计算IBS亲缘关系矩阵
perl run_pipeline.pl -Xmx10g -Xms512m -importGuess Troot.vcf -KinshipPlugin -method Centered_IBS -endPlugin -export tassel_kinship_Troot.txt -exportType SqrMatrix

1.2.2 在R中绘制热图,进行展示

> setwd("D:/研究生/数据/GWAS/亲缘关系")
> library("pheatmap")
> kinship<-read.table("tassel_kinship_Troot.txt",header = > F,row.names = 1,skip = 3)
> colnames(kinship) <- row.names(kinship)
> kinship[kinship < 0] <- 0
> B <- kinship[kinship == 0]
> diag(kinship) <- NA
> pdf("Histogram of Kinship.pdf",width=10,height=8)
> hist_data <- hist(as.matrix(kinship), xlab = "Kinship", col = "grey", main = "Histogram of Kinship")
> pdf("heapmap of kinship.pdf",width=8,height=8)
> pheatmap(kinship, fontsize_row = 0.3, fontsize_col = 0.3)
> dev.off()

2.GCTA

2.1下载及安装

详情见:

2.2 GCTA计算亲缘关系

这里需要将vcf文件在plink中转换为.bed .fam .bim格式。

./plink --vcf root.id.vcf --allow-extra-chr --make-bed --out root.id --autosome-num 27

2.2.1 Linux下计算亲缘关系

cd /gss1/home/fzhang/sunwei/gcta_1.92.0beta3
./gcta64 --make-grm-gz --out root.gcta --bfile root.id --autosome-num 27

生成两个文件,root.gcta.grm.id和root.gcta.grm.gz

2.2.1 在R中整理结果

> setwd("D:/研究生/数据/GWAS/亲缘关系")
> library(reshape2)
> tmp <- read.table(gzfile("root.gcta.grm.gz"), header = F, stringsAsFactors = F)
> ids <- read.table("root.gcta.grm.id", header = F, stringsAsFactors = F)
> tmp <- tmp[,c(1,2,4)]
> result_matrix <- acast(tmp, V1~V2, value.var="V4", drop = F)
> makeSymm <- function(m) {
  m[upper.tri(m)] <- t(m)[upper.tri(m)]
  return(m)
}
> result_full <- makeSymm(result_matrix)
> diag(result_full) <- 2
> result_df <- as.data.frame(result_full)
> row.names(result_df) <- ids$V2
> colnames(result_df) <- ids$V2
> write.table(result_df, file = "gcta.kinship.txt", row.names = T, col.names = NA, sep = "\t", quote = F)

2.2.2 在R中可视化结果

> gcta<- read.table("gcta.kinship.txt")
> library("pheatmap")
> pdf("heapmap of gcta_kinship.pdf",width=8,height=8)
> pheatmap(gcta, fontsize_row = 0.3, fontsize_col = 0.3)
> dev.off()

群体结构和亲缘关系的使用

动植物的关联分析中,群体结构(或PCA)作为固定效应,亲缘关系作为随机效应加入混合线性模型中。

引用转载请注明出处,如有错误敬请指出。

上一篇下一篇

猜你喜欢

热点阅读