draw heatmap from rqtl datatype

czheluo · Sep 29, 2019 · 8988272 · 8988272
1 parent 2828fd3
commit 8988272
Showing 1 changed file with 126 additions and 0 deletions.
diff --git a/R/drawmap-lod.R b/R/drawmap-lod.R
@@ -0,0 +1,126 @@
+library('getopt');
+options(bitmapType='cairo')
+spec = matrix(c(
+	'mark','m',1,'character',
+	'out','o',1,'character',
+	'pop','p',1,'character',
+	'help','h',0,'logical'
+	), byrow=TRUE, ncol=4)
+opt = getopt(spec)
+print_usage <- function(spec=NULL){
+	cat(getopt(spec, usage=TRUE));
+	cat("Usage example: \n")
+	cat("
+Usage example:
+	Rscript emap.R --mark --out --pop
+
+Usage:
+	--mark	map file
+	--out	out dir
+	--pop	pop type
+	--help		usage
+\n")
+	q(status=1);
+}
+times<-Sys.time()
+if ( !is.null(opt$help) ) { print_usage(spec) }
+if ( is.null(opt$mark) ) { print_usage(spec) }
+if ( is.null(opt$pop) ) { print_usage(spec) }
+
+if ( is.null(opt$out) ) { opt$out="./";}
+if(!dir.exists(opt$out)){dir.create(opt$out)}
+library('qtl');
+library('ASMap');
+
+
+opt$pop=tolower(opt$pop)
+if(opt$pop == "cp"){
+	d<-read.cross(genfile=paste(opt$mark,"qtl",sep="."),phefile=paste(opt$mark,"trt",sep="."),mapfile=paste(opt$mark,"sexAver.map",sep="."),format="mapqtl")
+	d1<-read.cross(genfile=paste(opt$mark,"qtl",sep="."),phefile=paste(opt$mark,"trt",sep="."),mapfile=paste(opt$mark,"male.map",sep="."),format="mapqtl")
+	d2<-read.cross(genfile=paste(opt$mark,"qtl",sep="."),phefile=paste(opt$mark,"trt",sep="."),mapfile=paste(opt$mark,"female.map",sep="."),format="mapqtl")
+	setwd(opt$out)
+	d<-jittermap(d)
+	d1<-jittermap(d1)
+	d2<-jittermap(d2)
+	pdf("total.sexAver.pdf");
+	plotMap(d,shift=TRUE,alternate.chrid=TRUE)
+	dev.off()
+	png("total.sexAver.png");
+	plotMap(d,shift=TRUE,alternate.chrid=TRUE)
+	dev.off()
+	pdf("total.male.pdf");
+	plotMap(d1,shift=TRUE,alternate.chrid=TRUE)
+	dev.off()
+	png("total.male.png");
+	plotMap(d1,shift=TRUE,alternate.chrid=TRUE)
+	dev.off()
+	pdf("total.female.pdf");
+	plotMap(d2,shift=TRUE,alternate.chrid=TRUE)
+	dev.off()
+	png("total.female.png");
+	plotMap(d2,shift=TRUE,alternate.chrid=TRUE)
+	dev.off()
+	chrname<-chrnames(d);
+	for(i in c(1:length(chrname))){
+		pdf(paste(chrname[i],".heatMap.sexAver.pdf",sep=""),width=18,height=10)
+		#plot.rf(d,chr=i,lmax=50,"both")
+		#plotRF(d,chr=i,lmax=50,"both")
+		heatMap(d,chr=i,lmax=50,what="lod")
+		dev.off()
+		png(paste(chrname[i],".heatMap.sexAver.png",sep=""),width=1150,height=800)
+		#plot.rf(d,chr=i,lmax=50,"both")
+		#plotRF(d,chr=i,lmax=50,"both")
+		heatMap(d,chr=i,lmax=50,what="lod")
+		dev.off()
+	}
+	chrname<-chrnames(d1);
+	for(i in c(1:length(chrname))){
+		pdf(paste(chrname[i],".heatMap.male.pdf",sep=""),width=18,height=10)
+		#plot.rf(d1,chr=i,lmax=50,"both")
+		#plotRF(d1,chr=i,lmax=50,"both")
+		heatMap(d,chr=i,lmax=50,what="lod")
+		dev.off()
+		png(paste(chrname[i],".heatMap.male.png",sep=""),width=1150,height=800)
+		#plot.rf(d1,chr=i,lmax=50,"both")
+		#plotRF(d1,chr=i,lmax=50,"both")
+		heatMap(d,chr=i,lmax=50,what="lod")
+		dev.off()
+	}
+	chrname<-chrnames(d2);
+	for(i in c(1:length(chrname))){
+		pdf(paste(chrname[i],".heatMap.female.pdf",sep=""),width=18,height=10)
+		#plot.rf(d2,chr=i,lmax=50,"both")
+		#plotRF(d2,chr=i,lmax=50,"both")
+		heatMap(d,chr=i,lmax=50,what="lod")
+		dev.off()
+		png(paste(chrname[i],".heatMap.female.png",sep=""),width=1150,height=800)
+		#plot.rf(d2,chr=i,lmax=50,"both")
+		#plotRF(d2,chr=i,lmax=50,"both")
+		heatMap(d,chr=i,lmax=50,what="lod")
+		dev.off()
+	}
+}else{
+	d<-read.cross(file=paste(opt$mark,"csv",sep="."),format="csvr",crosstype=opt$pop)
+	setwd(opt$out)
+	d<-jittermap(d)
+	d<-est.rf(d)
+	pdf("total.lg.pdf");
+	plotMap(d,shift=TRUE,alternate.chrid=TRUE)
+	dev.off()
+	png("total.lg.png");
+	plotMap(d,shift=TRUE,alternate.chrid=TRUE)
+	dev.off()
+	chrname<-chrnames(d);
+	for(i in c(1:length(chrname))){
+		pdf(paste(i,".heatMap.pdf",sep=""))
+		plotRF(d,chr=i,what='both')
+		dev.off()
+		png(paste(i,".heatMap.png",sep=""))
+		plotRF(d,chr=i,what='both')
+		dev.off()
+	}
+}
+
+escaptime=Sys.time()-times;
+print("Done!\n")
+print(escaptime)