• file: matrix_PairwiseDifferences.r
• document page: 33

code:

#----read data------------------------------------------------------------------
Data <- read.table("D:/Heidi/Master/R_Daten/SummaryStatistics/pairwise_differences.txt",
                    skip=10, row.names=1)
DataMatrix <- as.matrix(Data)

#----UnderMatrix----
n <- 2
x <- 1

UnderMatrix <- matrix(NA, ncol=ncol(DataMatrix), nrow=1)

while(n<=nrow(DataMatrix)){
  nextrow <- DataMatrix[n,1:x]
  nextrow <- cbind(t(nextrow), matrix(NA, ncol=ncol(DataMatrix)-x, nrow=1))
  
  UnderMatrix <- rbind(UnderMatrix, nextrow)
  
  n <- n+1
  x <- x+1
}


#----UpperMatrix----
n <- 2
x <- 1

UpperMatrix <- matrix(NA, ncol=x, nrow=1)
UpperMatrix <- cbind(UpperMatrix, t(DataMatrix[x,n:ncol(DataMatrix)]))

n <- n+1
x <- x+1
while(n<=nrow(DataMatrix)){
  nextrow <- matrix(NA, ncol=x, nrow=1)
  nextrow <- cbind(nextrow, t(DataMatrix[x,n:ncol(DataMatrix)]))
  
  UpperMatrix <- rbind(UpperMatrix, nextrow)
  
  n <- n+1
  x <- x+1
}
UpperMatrix <- rbind(UpperMatrix, matrix(NA, ncol=ncol(DataMatrix), nrow=1))


#----DiagonalMatrix----
n <- 1

DiagonalMatrix <- DataMatrix[n,n]
DiagonalMatrix <- cbind(DiagonalMatrix, matrix(NA, ncol=ncol(DataMatrix)-n, nrow=1))

n <- n+1
x <- 1
while(n<=nrow(DataMatrix)){
  nextrow <- matrix(NA, ncol=x, nrow=1)
  nextrow <- cbind(nextrow, t(DataMatrix[n,n]))
  nextrow <- cbind(nextrow, matrix(NA, ncol=ncol(DataMatrix)-n, nrow=1))
  
  DiagonalMatrix <- rbind(DiagonalMatrix, nextrow)
  
  n <- n+1
  x <- x+1
}


#----plot-----------------------------------------------------------------------
a <- ncol(DataMatrix)
b <- nrow(DataMatrix)

x <- c(1:a)
y <- c(1:b)

ColorRamp <- rgb( seq(1,0,length=256),  # Red
                   seq(1,0,length=256),  # Green
                   seq(1,1,length=256))  # Blue
                   
ColorRamp2 <- rgb( seq(1,0,length=256),  # Red
                   seq(1,0.6,length=256),  # Green
                   seq(1,0,length=256))  # Blue

ColorRamp3 <- rgb( seq(1,1,length=256),  # Red
                   seq(1,0.5,length=256),  # Green
                   seq(1,0,length=256))  # Blue

#----Mirror matrix (left-right)----
mirror.matrix <- function(x) {
  xx <- as.data.frame(x);
  xx <- rev(xx);
  xx <- as.matrix(xx);
  xx;
}

#----Rotate matrix 270 clockworks----
rotate270.matrix <- function(x) {
  mirror.matrix(t(x))
}

UnderMatrix <- rotate270.matrix(UnderMatrix)
UpperMatrix <- rotate270.matrix(UpperMatrix)
DiagonalMatrix <- rotate270.matrix(DiagonalMatrix)


#-----draw matrix---------------------------------------------------------------
#png("D:/Heidi/Master/R_Graphiken/matrix_PairwiseDifferences.png", width = 680, height = 660)

#----devide plot region in 4 parts----
def.par <- par(no.readonly = TRUE) # save default, for resetting...
  layout(rbind(c(1,2), c(1,3), c(1,4)),
         heights=rbind(c(2,1), c(2,1), c(2,1)),
         respect=rbind(c(0,1), c(0,0), c(0,0)))

#-----draw matrixe plot----       
  op <- par(mar=c(7.1, 5.1, 9.1, 2.1))
    image(x,y, UnderMatrix, col=ColorRamp, xlab="", ylab="", axes = FALSE)
    image(x,y, UpperMatrix, col=ColorRamp2, xlab="", ylab="", axes = FALSE, add=TRUE)
    image(x,y, DiagonalMatrix, col=ColorRamp3, xlab="", ylab="", axes = FALSE, add=TRUE)
       axis(1, at = c(1:a))
       axis(2, at = c(1:b), labels=c(b:1))
       box()
       mtext(text="Population", side=1, line=2.5)
       mtext(text="Population", side=2, line=2.5)
       mtext(text="Population average pairwise difference", line=4.5, cex=1.2, font=2)
  par(op)


#----draw legends----
#--upper legend--
  op2 <- par(mar=c(0.6, 0, 3.6, 7.6))
    UpperLegend <- seq(from=0, to=max(UpperMatrix, na.rm=TRUE), length=100)
    UpperLegend <- as.matrix(UpperLegend)

    image(y=UpperLegend, t(UpperLegend), col=ColorRamp2, axes=FALSE)
       axis(side=4, las=2)
       mtext(text="Average number of pairwise differences
between populations (PiXY)", side=4, line=4.5, cex=0.75)
       box()
  par(op2)

#--diagonal legend--
  op3 <- par(mar=c(2.1, 0, 2.1, 7.6))
    DiagonalLegend <- seq(from=0, to=max(DiagonalMatrix, na.rm=TRUE), length=100)
    DiagonalLegend <- as.matrix(DiagonalLegend)

    image(y=DiagonalLegend, t(DiagonalLegend), col=ColorRamp3, axes=FALSE)
       axis(side=4, las=2)
       mtext(text="Average number of pairwise
differences within population (PiX)", side=4, line=4.5, cex=0.75)
       box()
  par(op3)

#--under legend--
  op4 <- par(mar=c(3.6, 0, 0.6, 7.6))
    UnderLegend <- seq(from=0, to=max(UnderMatrix, na.rm=TRUE), length=100)
    UnderLegend <- as.matrix(UnderLegend)

    image(y=UnderLegend, t(UnderLegend), col=ColorRamp, axes=FALSE)
       axis(side=4, las=2)
       mtext(text="Corrected average pairwise
difference(PiXY-(PiX+PiY)/2)", side=4, line=4.5, cex=0.75)
       box()
  par(op4)

#dev.off()
par(def.par)  #reset to default