分类： 08 可视化 Visualization

$ head 41588_2015_BFng3406_MOESM37_ESM.csv
Trait1,Trait2,rg,se,z,p
ADHD,Age at Menarche,-0.153,0.08218,-1.858,0.063
ADHD,Age at Smoking,-0.036,0.2427,-0.147,0.883
ADHD,Alzheimer's,-0.055,0.2191,-0.249,0.803
ADHD,Anorexia,0.192,0.1162,1.649,0.099
ADHD,Autism Spectrum,-0.164,0.1438,-1.144,0.253
ADHD,BMI,0.287,0.08913,3.222,0.001
ADHD,BMI 2010,0.258,0.08606,2.993,0.003
ADHD,Bipolar,0.265,0.1539,1.722,0.085
ADHD,Birth Length,-0.055,0.1679,-0.329,0.742

install.packages("corrplot")
library(corrplot)
#读取数据
ldsc = read.csv(file = '41588_2015_BFng3406_MOESM37_ESM.csv',header=TRUE)

#提取需要画的表型
to_plot_col<-c('Age at Menarche','Alzheimer\'s','Anorexia','Autism Spectrum','BMI','Bipolar','Birth Length','Birth Weight','Childhood Obesity','Coronary Artery Disease','Crohn\'s Disease','Depression','Ever/Never Smoked','Fasting Glucose','HDL','Height','Infant Head Circumference','LDL','Rheumatoid Arthritis','Schizophrenia','T2D','Triglycerides','Ulcerative Colitis','Years of Education')

#构建一个rg矩阵，和一个p值矩阵,并命名row和col
gcm_rg<-matrix(1, nrow = length(to_plot_col), ncol = length(to_plot_col))
gcm_p<-matrix(1, nrow = length(to_plot_col), ncol = length(to_plot_col))
rownames(gcm_rg) <- to_plot_col
colnames(gcm_rg) <- to_plot_col
rownames(gcm_p) <-to_plot_col
colnames(gcm_p) <- to_plot_col

#将数据填进矩阵
for (row in 1:nrow(ldsc)) {
    if(ldsc[row, "Trait1"]%in%to_plot_col && ldsc[row, "Trait2"]%in%to_plot_col){
    Trait1 <- ldsc[row, "Trait1"]
    Trait2  <- ldsc[row, "Trait2"]
    gcm_rg[Trait1,Trait2]<-ldsc[row, "rg"]
    gcm_p[Trait1,Trait2]<-p<-ldsc[row, "p"]
    gcm_rg[Trait2,Trait1]<-ldsc[row, "rg"]
    gcm_p[Trait2,Trait1]<-p<-ldsc[row, "p"]
    }
    }

#使用corrplot绘图：（方法选择square）
#sig.level = 0.05/300 为显著水平（这里是多重检验调整的）
#gcm_rg为相关系数矩阵
#gcm_p为p值矩阵
#order="hclust" 为通过分层聚类对矩阵进行排序
#其他选项均为细节调整 可以参考 https://cran.r-project.org/web/packages/corrplot/vignettes/corrplot-intro.html 

corrplot(gcm_rg, p.mat=gcm_p, method="square",order="hclust",sig = "pch",sig.level = 0.05/300 ,pch="*",pch.cex=2,
         tl.col="black",tl.srt=45,addgrid.col="grey90")

install.packages("dplyr")
devtools::install_github("mkanai/corrplot")

#读取数据
ldsc = read.csv(file = '41588_2015_BFng3406_MOESM37_ESM.csv',header=TRUE)

#计算FDR调整后p值
ldsc$fdr<-p.adjust(ldsc$p, method ="fdr")

to_plot_col<-c('Age at Menarche','Alzheimer\'s','Anorexia','Autism Spectrum','BMI','Bipolar','Birth Length','Birth Weight','Childhood Obesity','Coronary Artery Disease','Crohn\'s Disease','Depression','Ever/Never Smoked','Fasting Glucose','HDL','Height','Infant Head Circumference','LDL','Rheumatoid Arthritis','Schizophrenia','T2D','Triglycerides','Ulcerative Colitis','Years of Education')

gcm_rg<-matrix(NA, nrow = length(to_plot_col), ncol = length(to_plot_col))
gcm_p<-matrix(NA, nrow = length(to_plot_col), ncol = length(to_plot_col))
rownames(gcm_rg) <- to_plot_col
colnames(gcm_rg) <- to_plot_col
rownames(gcm_p) <-to_plot_col
colnames(gcm_p) <- to_plot_col

#p值矩阵改为fdr矩阵
for (row in 1:nrow(ldsc)) {
    if(ldsc[row, "Trait1"]%in%to_plot_col && ldsc[row, "Trait2"]%in%to_plot_col){
    Trait1 <- ldsc[row, "Trait1"]
    Trait2  <- ldsc[row, "Trait2"]
    gcm_rg[Trait1,Trait2]<-ldsc[row, "rg"]
    gcm_p[Trait1,Trait2]<-ldsc[row, "fdr"]
    gcm_rg[Trait2,Trait1]<-ldsc[row, "rg"]
    gcm_p[Trait2,Trait1]<-ldsc[row, "fdr"]
    }
    }

#事先对矩阵进行排序（此修改的版本直接使用order="hclust"时有bug）
gcm_rg.hclust<-corrMatOrder(gcm_rg,order ="hclust")
gcm_rg_hclust<-gcm_rg[gcm_rg.hclust,gcm_rg.hclust]
gcm_p_hclust<-gcm_p[gcm_rg.hclust,gcm_rg.hclust]

options(repr.plot.width=15, repr.plot.height=15)
#使用修改后的corrplot绘图：（注意！！方法选择 psquare）
corrplot(gcm_rg_hclust ,p.mat=gcm_p_hclust, diag=FALSE,method="psquare",sig="pch",insig="label_sig",sig.level = 0.01 ,pch.cex=2,tl.col="black",tl.srt=45,addgrid.col="grey90")

1. 命令行版本简介与软件与参考数据准备（**只看网页版的话可以跳过 1**）
2. 输入文件格式 
3. 指定绘制的区间范围
4. 指定LD文件，族裔群体，参考基因组版本
	**-》以上内容对应的网页版本**
5. 批量模式 ！ Batch mode（命令行版本的优势所在）
	**-》批量模式对应的网页版本**
6. 自制LD文件
7 额外功能
	7.1指定多个参考SNP
	7.2 使用GWAS catalog进行注释
	7.3 精确定位的可信集
	7.4 对多个SNP的注释
	7.5 绘制额外的BED文件轨道
参考url

cd <directory where you want to place locuszoom> 
tar zxf /path/to/locuszoom.tgz
ln -s bin/locuszoom /usr/local/bin/locuszoom #根据自己的环境变量制定

locuszoom/
	bin/
		locuszoom (this is the locuszoom "executable")
		locuszoom.R (the R script which is used by locuszoom for creating the plots)
		dbmeister.py (script for creating custom user databases)
		lzupdate.py (script for creating an updated copy of the provided locuszoom database)
	conf/ (configuration file located here)
	data/
		database/ (SQLite file located here)
		hapmap/ (hapmap genotype files)
		1000G/ (1000G genotype files)
	src/ (source code for locuszoom)

locuzoom在输入时的选项：
--metal 指定metal格式的输入文件名
--delim （可选）可以指定分隔符，默认为tab
--markercol  指定输入文件表示marker的列名，默认为"MarkerName"
--pvalcol 指定输入文件表示p值的列名，默认为"P-value"
--no-transform 当p已经转换为-log10（p）时，可以使用此选项跳过log转换

#2.1指定SNP与两侧的范围
--refsnp <your snp> --flank 500kb

#2.2指定SNP与区间的起止位点
--refsnp <your snp> --chr # --start <base position> --end <base position>

#2.3指定SNP与基因（绘制基因所在的范围）
--refsnp <your snp> --refgene <your gene>

#2.4 指定基因与两侧范围
--refgene <your gene> --flank 250kb

#2.5 指定基因与区间的起止位点
--refgene <your gene> --chr # --start <base position> --end <base position>

#2.6 区间的起止位点
--chr # --start <base position> --end <base position>

Genotype files available for: 
--source hapmap
  --build hg18
    --pop YRI
    --pop CEU
    --pop JPT+CHB

--source 1000G_March2012
  --build hg19
    --pop AMR
    --pop ASN
    --pop AFR
    --pop EUR

--source 1000G_June2010
  --build hg18
    --pop YRI
    --pop CEU
    --pop JPT+CHB

--source 1000G_Nov2014
  --build hg19
    --pop AMR
    --pop ASN
    --pop AFR
    --pop EUR
  --build hg38
    --pop SAS
    --pop EAS
    --pop AMR
    --pop AFR
    --pop EUR

--pop ASN --build hg19 --source 1000G_March2012
#指定1000G_March2012数据，以hg19为参考基因坐标，来计算ASN族裔的LD