• file: matrix_multiplePlots.r
• document page: 25

code:

#----read data------------------------------------------------------------------
#----read haplotype list----
Data <- read.table("D:/Heidi/Master/R_Daten/HaplotypeDistance/ListHaplotype_betweenPop.txt", skip=1)

Row <- nrow(Data)
Columns <- Row


#----read data row by row:----
x <- 0
n <- 1

DistanceMatrix <- as.matrix(scan("D:/Heidi/Master/R_Daten/HaplotypeDistance/
                        HapDistanceMatrix_betweenPop.txt", what=double(0), skip=x, nlines=1, nmax=n))
DistanceMatrix <- cbind(DistanceMatrix, matrix(NA, ncol=Columns, nrow=1))
DistanceMatrix <- DistanceMatrix[,1:Columns]

n <- n + 1
x <- x + 1

while(n<(Row+1)){
  nextrow <- as.matrix(scan("D:/Heidi/Master/R_Daten/HaplotypeDistance/
                       HapDistanceMatrix_betweenPop.txt", what=double(0), skip=x, nlines=1, nmax=n))
  nextrow <- cbind(t(nextrow), matrix(NA, ncol=Columns, nrow=1))
  nextrow <- nextrow[,1:Columns]

  DistanceMatrix <- rbind(DistanceMatrix, nextrow)

  n <- n + 1
  x <- x + 1
}


#----select the data of DistanceMatrix------------------------------------------
dimnames(DistanceMatrix) <- list(Data[,1], Data[,1])

#--population1:
Pop1 <- read.table("D:/Heidi/Master/R_Daten/HaplotypeDistance/HapDistanceMatrix_withinPop1.txt", skip=5)
Pop1 <- as.character(Pop1[,1])
DistanceMatrixPop1 <- DistanceMatrix[Pop1, Pop1]

#--population2:
Pop2 <- read.table("D:/Heidi/Master/R_Daten/HaplotypeDistance/HapDistanceMatrix_withinPop2.txt", skip=5)
Pop2 <- as.character(Pop2[,1])
DistanceMatrixPop2 <- DistanceMatrix[Pop2, Pop2]

#--between population1/2:
#whole DistanceMatrix
wholeDistanceMatrix <- DistanceMatrix
x <- 1
while(x <= ncol(wholeDistanceMatrix)){

  twholeDistanceMatrix <- t(wholeDistanceMatrix) 
  wholeDistanceMatrix[x,(x:ncol(wholeDistanceMatrix))] <- twholeDistanceMatrix[x,(x:ncol(wholeDistanceMatrix))]

  x <- x + 1
}                                                                  
wholeDistanceMatrixBetween <- wholeDistanceMatrix[Pop2, Pop1]

 
#----together:-----
DistanceMatrixUp <- cbind(DistanceMatrixPop1, matrix(NA, ncol(DistanceMatrixPop1), nrow(DistanceMatrixPop2)))
DistanceMatrixLo <- cbind(wholeDistanceMatrixBetween, DistanceMatrixPop2) 
DistanceMatrixTogether <- rbind(DistanceMatrixUp, DistanceMatrixLo)


#----whole DistanceMatrix both Populations together:----
Pop_together <- c(Pop1, Pop2)

wholeDistanceMatrixTogether <- wholeDistanceMatrix[Pop_together, Pop_together[c(length(Pop_together)):1]]


#----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))
}

DistanceMatrixTogether <- rotate270.matrix(DistanceMatrixTogether)


#----draw matrix----------------------------------------------------------------
windows()  #new graphic window  
library(fields)
 ColorRamp <- rgb( seq(1,0,length=256),  # Red
                   seq(1,0,length=256),  # Green
                   seq(1,1,length=256))  # Blue
                   
 a <- ncol(DistanceMatrixTogether)
 b <- nrow(DistanceMatrixTogether)

 x <- c(1:a)
 y <- c(1:b)
 
png("D:/Heidi/Master/R_Graphiken/matrix_HapDistance_between-within Pop and Groups.png", width = 680, height = 660) 
op <- par(mfrow=c(2,2), oma=c(0, 0 ,3, 0))

#----matrix one---- 
op1 <- par(mfg=c(1,1),mar=c(0.2, 5.2, 4.2, 0.2))
  image(x,y,DistanceMatrixTogether, col=ColorRamp, ylab="Group 1", font.lab=2, axes = FALSE)
     contour(DistanceMatrixTogether, add = TRUE)
     axis(2, at = c(1:b), labels=c(Pop2[NROW(Pop2):1],Pop1[NROW(Pop1):1]), cex.axis=0.6, las=2)
     box(bty="l")
        
     lines(c(0,NROW(Pop1)+0.5),c(NROW(Pop2)+0.5,NROW(Pop2)+0.5), lwd=2)
     lines(c(NROW(Pop1)+0.5,NROW(Pop1)+0.5),c(0,NROW(Pop2)+0.5), lwd=2)
        
     mtext(side = 2, at =(3.2*Row/4), line = 2.7, text = "Population 1", cex=0.7)
     mtext(side = 2, at =(Row/4), line = 2.7, text = "Population 2", cex=0.7)
par(op1)

#----matrix two----
op2 <- par(mfg=c(2,1), mar=c(4.2, 5.2, 0.2, 0.2))
  image(x,y, wholeDistanceMatrixTogether, col=ColorRamp, 
             xlab="Group 1", ylab="Group 2", font.lab=2, axes = FALSE)
     contour(DistanceMatrixTogether, add = TRUE)
     axis(1, at = c(1:a), labels=c(Pop_together), cex.axis=0.6, las=2)
     axis(2, at = c(1:b), labels=c(Pop2[NROW(Pop2):1],Pop1[NROW(Pop1):1]), cex.axis=0.6, las=2)
     box()
        
     lines(c(0,NROW(Pop1)+ NROW(Pop2)+0.5),c(NROW(Pop2)+0.5,NROW(Pop2)+0.5), lwd=2)
     lines(c(NROW(Pop1)+0.5,NROW(Pop1)+0.5),c(0,NROW(Pop1)+ NROW(Pop2)+0.5), lwd=2)
        
     mtext(side = 1, at =(Row/4), line = 2.5, text = "Population 1", cex=0.7)
     mtext(side = 1, at =(3.2*Row/4), line = 2.5, text = "Population 2", cex=0.7)
     mtext(side = 2, at =(3.2*Row/4), line = 2.7, text = "Population 3", cex=0.7)
     mtext(side = 2, at =(Row/4), line = 2.7, text = "Population 4", cex=0.7)
par(op2)

#----matrix three----
op3 <- par(mfg=c(2,2), mar=c(4.2, 0.2, 0.2, 6.2))
  image(x,y,DistanceMatrixTogether, col=ColorRamp, xlab="Group 2", font.lab=2, axes = FALSE)
     contour(DistanceMatrixTogether, add = TRUE)
     axis(1, at = c(1:a), labels=c(Pop1,Pop2), cex.axis=0.6, las=2)
     box(bty="l")
        
     lines(c(0,NROW(Pop1)+0.5),c(NROW(Pop2)+0.5,NROW(Pop2)+0.5), lwd=2)
     lines(c(NROW(Pop1)+0.5,NROW(Pop1)+0.5),c(0,NROW(Pop2)+0.5), lwd=2)
        
     mtext(side = 1, at =(Row/4), line = 2.5, text = "Population 3", cex=0.7)
     mtext(side = 1, at =(3.2*Row/4), line = 2.5, text = "Population 4", cex=0.7)
par(op3)

#----legend---- 
maximum <- max(DistanceMatrixTogether, wholeDistanceMatrixTogether, DistanceMatrixTogether , na.rm = TRUE)

Legend <- seq(from=0, to=maximum, length=100)
Legend <- as.matrix(Legend)

op4 <- par(mfg=c(1,2), mar=c(0, 15.6, 3.2, 6.2) ) 
  image(y=Legend, t(Legend), col=ColorRamp, axes=FALSE)
     axis(side=4, las=2)
     box()
     
  title("haplotype distance matrix between/within populations and groups", line=0, outer=TRUE)
par(op4)

#----reset par parameter----
par(op)
dev.off()