I am having a problem generating a barplot in R where some bars are failing to visualize and other bars are taking up double or triple the width they should. I am trying to plot some ADMIXTURE results from ADMIXTURE. I got R some code from one of my adviser's collaborators to plot my results, but I am having serious problems.
I have noticed that R will randomly make certain individual bars double width and other bars zero width (i.e., some samples will fail to appear and simultaneously others will be double). In my dataset, there are 2033 individuals in the dataset from 166 populations (=165 separating bars). The double bars will often crossover population separating bars. I have been trying to tinker with the output resolution, but nothing seems to fix the issue.
Is there anyway to set it so that no bar can be wider than one pixel wide?
If each bar is locked to one pixel wide, the barplot itself should only be 2198 pixels wide.
The code
options(stringsAsFactors=F)
library(MASS)
library(RColorBrewer)
############################
structure.plot <- function(x, cex.axis=0.5) {
qcols <- sapply(names(x), function(y) {length(grep("Q", y)) > 0})
q <- x[,qcols]
n <- ncol(q)
cols <- c("#6BC367","#BF59CB","#CC5735","#8BB5D0","#563C6A","#CDAF4F","#CA4A7D","#64322D","#5A6630","#CD9B9B","#7D7CCE","#8FD0AD","#455E5D","#9BDA3E")
barplot(t(as.matrix(q)), space=0, border=NA, col=cols, axes=F, axisnames=F)
ppos <- c()
for (p in levels(x$pop)) {
ppos <- c(ppos, mean(which(x$pop == p)))
}
axis(1, at=ppos, labels=levels(qdata$pop), tick=F, las=2, cex.axis=cex.axis)
pmin <- c()
for (p in levels(x$pop)) {
pmin <- c(pmin, min(which(x$pop == p)))
}
pmin <- pmin - 0.5
pmin <- pmin[pmin > 1]
abline(v=pmin)
}
# reorder a Q-dataframe by the component with the largest values
mq_reorder <- function(x) {
qmeans <- apply(x, 2, mean)
qmax <- which(qmeans == max(qmeans))
x[order(x[,qmax], decreasing=T),]
}
##############################################################
qdata <- read.delim("Qdata.txt")
popinfo <- read.delim("order.txt")
##############################################################
# merge order and region information from popinfo,
# while also dropping the initial region column from qdata
qdata <- merge(subset(qdata, select=setdiff(names(qdata), c("region"))), popinfo,
by.x='pop', by.y='pop')
qdata$pop <- factor(qdata$pop, levels=popinfo$pop, ordered=T)
qdata <- qdata[order(qdata$order),]
##############################################################
qcols <- sapply(names(qdata), function(y) {length(grep("Q", y)) > 0})
for (p in levels(qdata$pop)) {
prows <- (qdata$pop == p)
qdata[prows, qcols] <- mq_reorder(qdata[prows, qcols])
}
tiff(file="structure-plot-A.tiff", width=2500, height=500, antialias="default")
structure.plot(qdata, cex.axis=.8)
dev.off()
The Qdata.txt
file takes the format of:
iid Q1 Q2 Q13 Q14 pop region bigregion
Sample1 0 0 0 0 popname regionname bigregionname
Sample2 0 0 0 0 popname regionname bigregionname
Sample3 0 0 0 0 popname regionname bigregionname
Sample4 0 0 0 0 popname regionname bigregionname
Sample5 0 0 0 0 popname regionname bigregionname
Sample6 0 0 0 0 popname regionname bigregionname
Sample7 0 0 0 0 popname regionname bigregionname
Sample8 0 0 0 0 popname regionname bigregionname
Sample9 0 0 0 0 popname regionname bigregionname
Sample10 0 0 0 0 popname regionname bigregionname
All those 0s are stand-ins for values between 0 and 1 that add up to 1. The columns for Q3-12 are omitted for brevity
The order.txt files takes the format of:
order region pop
1 EasternEurope Albanian
2 EasternEurope Belarusian
3 EasternEurope Bulgarian
4 EasternEurope Chuvash
Anyone have an idea of how to fix this code, so each bar is locked at one pixel?
Thanks!