User Tools

Site Tools


haplotype_distance_between_within_two_populations

haplotype distance between/within two populations

  • file: matrix_HapDistance_between-within.r
  • document page: 11


code:

#----read haplotype list----
Data <- read.table("D:/Heidi/Master/R_Daten/HaplotypeDistance/ListHaplotype_betweenBsp.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_betweenBsp.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_betweenBsp.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
}

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

DistanceMatrix <- rotate270.matrix(DistanceMatrix)


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

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

 image.plot(x,y,DistanceMatrix, col=ColorRamp, main="haplotype distance matrixbetween/within populations",
             xlab="Haplotype", ylab="Haplotype", axes = FALSE)
        contour(DistanceMatrix, add = TRUE)
        axis(1, at = c(1:a), labels=Data[,1], cex.axis=0.7)
        axis(2, at = c(1:b), labels=Data[(Row:1),1], cex.axis=0.7)
        box()
        
        half <- (Row/2) + 0.5

        lines(c(0,half),c(half,half), lwd=2)
        lines(c(half,half),c(0,half), lwd=2)        
        mtext(side = 1, at =(Row/4), line = 2, text = "Population 1", cex=0.8)
        mtext(side = 1, at =(3*Row/4), line = 2, text = "Population 2", cex=0.8)
        mtext(side = 2, at =(Row/4), line = 2, text = "Population 2", cex=0.8)
        mtext(side = 2, at =(3*Row/4), line = 2, text = "Population 1", cex=0.8)

mixed data from several populations

  • file: matrix_HapDistance_between-within_complex.r
  • document page: 13


code:

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


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

 image.plot(x,y,DistanceMatrixTogether, col=ColorRamp, main="haplotype distancematrix between/within populations", 
             xlab="Haplotype", ylab="Haplotype", axes = FALSE)
        contour(DistanceMatrixTogether, add = TRUE)
        axis(1, at = c(1:a), labels=c(Pop1,Pop2), cex.axis=0.6)
        axis(2, at = c(1:b), labels=c(Pop2[NROW(Pop2):1],Pop1[NROW(Pop1):1]), cex.axis=0.6)
        box()  
        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, text = "Population 1", cex=0.8)
        mtext(side = 1, at =(3*Row/4), line = 2, text = "Population 2", cex=0.8)
        mtext(side = 2, at =(Row/4), line = 2, text = "Population 2", cex=0.8)
        mtext(side = 2, at =(3*Row/4), line = 2, text = "Population 1", cex=0.8)
haplotype_distance_between_within_two_populations.txt · Last modified: 2008/07/22 13:31 by 127.0.0.1