Set uniform bar widths for ADMIXTURE barplots in R
0
0
Entering edit mode
8.7 years ago
devenvyas ▴ 740

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!

admixture R • 3.1k views
ADD COMMENT

Login before adding your answer.

Traffic: 1849 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6