數(shù)據(jù)過濾質(zhì)控之后
PCA分析（gcta）
###make germ
nohup /home/software/gcta_1.92.3beta3/gcta64 --bfile ld.QC.75_noinclude0-502502-geno02-maf03 --make-grm --autosome-num 26 --out ld.QC.75_noinclude0-502502-geno02-maf03.gcta &
# autosome表示只選出常染色體來運(yùn)行
###pca
nohup /home/software/gcta_1.92.3beta3/gcta64 --grm ld.QC.75_noinclude0-502502-geno02-maf03.gcta --pca 5 --out ld.QC.75_noinclude0-502502-geno02-maf03.gcta.out &
########輸入的文件就是.gcta,這個文件不是上一步的生成文件，pca后面跟的是要做幾個主成分的比較
利用R繪圖
####
##在excel中將pca.eigenvec文件的第一列改為品種ID即可
a=read.table("pca.eigenvec",header=F)
a$V1 <- factor(a$V1,levels=c("Angus","Holstein","Simmental","Shorthorn",
"Charolais","Mishima-Ushi","Hanwoo","JPBC","Kazakh",
"Mongolian","Yanbian","Wenling","Brahman")) ## 修改圖例順序
Breed=b[,1]
p = ggplot(data = a ,
aes(x = a[,3],
y = a[,4], ############x = a[,3], y = a[,4]代表第3列和第4列比較也就是pc1與pc2比較）
group = Breed,
shape = Breed,
stroke = 1.5))+
geom_point(aes(color=Breed),size=3) +
scale_shape_manual(values = seq(1,15)) + ##設(shè)置形狀
# scale_shape_manual(values = rep(12,times = 13)) +
guides(color=guide_legend(override.aes = list(size=4, stroke = 1.5))) ##改變圖例大小
p + labs(x = "PC1(12.44%)", y = "PC2(5.87%)") + ############修改x軸和y軸的坐標(biāo)名稱
scale_color_manual(values=c("#66a61e", "#66a61e","#66a61e","#66a61e","#66a61e",
"#7470b3", "#e7288a","#e728e7","#e7288a","#e7288a","#e7288a",
"#b35107","#b35107")) +
theme_classic()+ # 去除灰色背景及網(wǎng)格線
theme(panel.border = element_rect(fill=NA,color="black", size=0.5, linetype ="solid")) #添加邊框
#####圖就出來啦
