######################################################################################### # Plots pie charts for the autosomes reflecting the different genomes and by the color # the population to which more than 50% of the individuals with the corresponding # IBD segment belong ######################################################################################### library(ggplot2) library(doBy) library(grid) load(file = paste("DenisovaMatGNew.Rda",sep = "")) load(file = paste("NeanderMatGNew.Rda",sep = "")) load(file = paste("AncDenisovaMatGNew.Rda",sep = "")) load(file = paste("AncNeanderMatGNew.Rda",sep = "")) load(file = paste("PopulationMatG.Rda",sep = "")) real1 <- which(NeanderMatG[,"IBDlengthCorrected"] > 10) real2 <- which(NeanderMatG[,"#SNVs"] > 7) real <- intersect(real1,real2) DenisovaMatG <- DenisovaMatG[real,] NeanderMatG <- NeanderMatG[real,] AncDenisovaMatG <- AncDenisovaMatG[real,] AncNeanderMatG <- AncNeanderMatG[real,] PopulationMatG <- PopulationMatG[real,] de1 <- which(DenisovaMatG[,"AltR"] > 0.15) de2 <- which(DenisovaMatG[,"ancientRatio"] > 0.3) de <- intersect(de1,de2) ne1 <- which(NeanderMatG[,"AltR"] > 0.15) ne2 <- which(NeanderMatG[,"ancientRatio"] > 0.3) ne <- intersect(ne1,ne2) ##### Ancestor correction #### nanca <- which(AncNeanderMatG[,"ancestorRatio"] < 0.3) danca <- which(AncDenisovaMatG[,"ancestorRatio"] < 0.3) dnanca <- intersect(nanca,danca) de <- intersect(de,danca) ne <- intersect(ne,nanca) ################################################################################ ################################################################################ ################################################################################ ################################################################################ ####################### Majority Population - pie ####################### lde <- length(de) lne <- length(ne) ##### maj thresh ##### thresh <- 0.5 lll <- length(real) llv <- rep("",lll) for (i in 1:lll) { wm <- which(DenisovaMatG[i,c(25, 27:29)]/sum(DenisovaMatG[i,c(25, 27:29)]) > thresh) if (length(wm) > 0) { llv[i] <- names(wm) } } dev <- rep("",lde) for (i in 1:lde) { wm <- which(DenisovaMatG[de[i],c(25,27:29)]/sum(DenisovaMatG[de[i],c(25,27:29)]) > thresh) if (length(wm) > 0) { dev[i] <- names(wm) } } nev <- rep("",lne) for (i in 1:lne) { wm <- which(NeanderMatG[ne[i],c(25,27:29)]/sum(NeanderMatG[ne[i],c(25,27:29)]) > thresh) if (length(wm) > 0) { nev[i] <- names(wm) } } llvt <- table(llv) llvtl <- paste(names(llvt), " ", llvt, sep = "") devt <- table(dev) devtl <- paste(names(devt), " ", devt, sep = "") nevt <- table(nev) nevtl <- paste(names(nevt), " ", nevt, sep = "") llvtfrac <- llvt/lll devtfrac <- devt/lde nevtfrac <- nevt/lne col <- c("grey","#00b050", "#ffc000", "#0070c0", "yellow") name = "genomesMajPopulationPieAll_majThresh05_G" # xfig(paste(name,".fig",sep = ""),width = 10, height = 10) pie(llvt, labels = llvtl, col = col, radius = 1, cex = 1.2) dev.off() name = "genomesMajPopulationPieDenisova_majThresh05_G" # xfig(paste(name,".fig",sep = ""),width = 10, height = 10) pie(devt, labels = devtl, col = col, radius = 1, cex = 1.2) dev.off() name = "genomesMajPopulationPieNeandertal_corr_wocorrAnc_woDENE_woAMR_majThresh50_G" # xfig(paste(name,".fig",sep = ""), width = 10, height = 10) pie(nevt, labels = nevtl, col = col, radius = 1, cex = 1.2) dev.off() name = "genomesMajPopulationPieComb_majThresh05_G" # xfig(paste(name, ".fig", sep = ""), width = 10, height = 10) # pdf(paste(name, ".pdf", sep = ""), width = 30, height = 10) par(mar = c(3, 6, 3, 5)) layout(matrix(c(1, 2, 3), 1, 3, byrow = TRUE), c(10, 10, 10), c(10, 10, 10)) pie(llvt, labels = llvtl, col = col, radius = 1, cex = 1) pie(nevt, labels = nevtl, col = col, radius = 1, cex = 1) pie(devt, labels = devtl, col = col, radius = 1, cex = 1) dev.off()