[Patchwork-commits] r216 - .git .git/logs .git/logs/refs/heads .git/logs/refs/remotes/origin .git/refs/heads .git/refs/remotes/origin pkg/TAPS/R
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Thu Jun 18 10:25:41 CEST 2015
Author: sebastian_d
Date: 2015-06-18 10:25:41 +0200 (Thu, 18 Jun 2015)
New Revision: 216
Modified:
.git/index
.git/logs/HEAD
.git/logs/refs/heads/master
.git/logs/refs/remotes/origin/master
.git/refs/heads/master
.git/refs/remotes/origin/master
pkg/TAPS/R/TAPS.r
Log:
markus update taps
Modified: .git/index
===================================================================
(Binary files differ)
Modified: .git/logs/HEAD
===================================================================
--- .git/logs/HEAD 2014-09-29 08:35:21 UTC (rev 215)
+++ .git/logs/HEAD 2015-06-18 08:25:41 UTC (rev 216)
@@ -98,3 +98,4 @@
6734fab24795a8c823c6a97e0131fcfe5a235107 cefeee2662b45c7daf3b8916b967173f68b6ddcf Sebastian DiLorenzo <dilorenzo.sebastian at gmail.com> 1409562736 +0200 pull: Fast-forward
cefeee2662b45c7daf3b8916b967173f68b6ddcf e8cc7370da527f1d46bb81720f663ceed281e07d Sebastian DiLorenzo <dilorenzo.sebastian at gmail.com> 1409832387 +0200 pull: Fast-forward
e8cc7370da527f1d46bb81720f663ceed281e07d b5798bb1874d9dc8e50fa5083961a8836d79f76c Sebastian DiLorenzo <dilorenzo.sebastian at gmail.com> 1411979666 +0200 pull: Fast-forward
+b5798bb1874d9dc8e50fa5083961a8836d79f76c 52402622419759bc3642b6f49dd0f4abd2003d38 Sebastian DiLorenzo <dilorenzo.sebastian at gmail.com> 1434615906 +0200 pull: Fast-forward
Modified: .git/logs/refs/heads/master
===================================================================
--- .git/logs/refs/heads/master 2014-09-29 08:35:21 UTC (rev 215)
+++ .git/logs/refs/heads/master 2015-06-18 08:25:41 UTC (rev 216)
@@ -98,3 +98,4 @@
6734fab24795a8c823c6a97e0131fcfe5a235107 cefeee2662b45c7daf3b8916b967173f68b6ddcf Sebastian DiLorenzo <dilorenzo.sebastian at gmail.com> 1409562736 +0200 pull: Fast-forward
cefeee2662b45c7daf3b8916b967173f68b6ddcf e8cc7370da527f1d46bb81720f663ceed281e07d Sebastian DiLorenzo <dilorenzo.sebastian at gmail.com> 1409832387 +0200 pull: Fast-forward
e8cc7370da527f1d46bb81720f663ceed281e07d b5798bb1874d9dc8e50fa5083961a8836d79f76c Sebastian DiLorenzo <dilorenzo.sebastian at gmail.com> 1411979666 +0200 pull: Fast-forward
+b5798bb1874d9dc8e50fa5083961a8836d79f76c 52402622419759bc3642b6f49dd0f4abd2003d38 Sebastian DiLorenzo <dilorenzo.sebastian at gmail.com> 1434615906 +0200 pull: Fast-forward
Modified: .git/logs/refs/remotes/origin/master
===================================================================
--- .git/logs/refs/remotes/origin/master 2014-09-29 08:35:21 UTC (rev 215)
+++ .git/logs/refs/remotes/origin/master 2015-06-18 08:25:41 UTC (rev 216)
@@ -95,3 +95,4 @@
6734fab24795a8c823c6a97e0131fcfe5a235107 cefeee2662b45c7daf3b8916b967173f68b6ddcf Sebastian DiLorenzo <S_D at Sebastians-iMac.local> 1409562735 +0200 pull: fast-forward
cefeee2662b45c7daf3b8916b967173f68b6ddcf e8cc7370da527f1d46bb81720f663ceed281e07d Sebastian DiLorenzo <S_D at array-47-13.medsci.uu.se> 1409832387 +0200 pull: fast-forward
e8cc7370da527f1d46bb81720f663ceed281e07d b5798bb1874d9dc8e50fa5083961a8836d79f76c Sebastian DiLorenzo <S_D at array-47-13.medsci.uu.se> 1411979666 +0200 pull: fast-forward
+b5798bb1874d9dc8e50fa5083961a8836d79f76c 52402622419759bc3642b6f49dd0f4abd2003d38 Sebastian DiLorenzo <dilorenzo.sebastian at gmail.com> 1434615905 +0200 pull: fast-forward
Modified: .git/refs/heads/master
===================================================================
--- .git/refs/heads/master 2014-09-29 08:35:21 UTC (rev 215)
+++ .git/refs/heads/master 2015-06-18 08:25:41 UTC (rev 216)
@@ -1 +1 @@
-b5798bb1874d9dc8e50fa5083961a8836d79f76c
+52402622419759bc3642b6f49dd0f4abd2003d38
Modified: .git/refs/remotes/origin/master
===================================================================
--- .git/refs/remotes/origin/master 2014-09-29 08:35:21 UTC (rev 215)
+++ .git/refs/remotes/origin/master 2015-06-18 08:25:41 UTC (rev 216)
@@ -1 +1 @@
-b5798bb1874d9dc8e50fa5083961a8836d79f76c
+52402622419759bc3642b6f49dd0f4abd2003d38
Modified: pkg/TAPS/R/TAPS.r
===================================================================
--- pkg/TAPS/R/TAPS.r 2014-09-29 08:35:21 UTC (rev 215)
+++ pkg/TAPS/R/TAPS.r 2015-06-18 08:25:41 UTC (rev 216)
@@ -14,7 +14,7 @@
TAPS_plot <- function(#samples='all',
directory=NULL,autoEstimate=FALSE,
- bin=400,cores=1,matched=FALSE,allelePeaks=FALSE) {
+ bin=250,cores=1,matched=FALSE,allelePeaks=FALSE) {
#Automatically check, and if needed install, packages stats and fields
#Load stats. It should be in all, at least semi-new, R distributions so we dont need to install.package it or
@@ -141,13 +141,17 @@
alf$Start = as.integer(as.character(alf$Start))
alf$End = as.integer(as.character(alf$End))
alf$Value = as.numeric(as.character(alf$Value))
- }
- else
- {
+ } else if (file.exists('rawcopy.Rdata')) {
+ load('rawcopy.Rdata')
+ Log2=probes.txt[,2:5]
+ alf=snps.txt[,2:5]
+ } else {
Log2 <- readLog2() ## Log-R
alf <- readAlf() ## Allele Frequency
}
+ #browser()
+
Log2=Log2[!is.nan(Log2$Value),]
Log2=Log2[!is.na(Log2$Value),]
Log2 <- Log2[which(Log2$Value != -Inf & Log2$Value != +Inf ),]
@@ -156,7 +160,7 @@
alf$Value[alf$Value<0]=0; alf$Value[alf$Value>1]=1
alf <- alf[which(alf$Value != -Inf & alf$Value != +Inf ),]
- segments <- readSegments() ## segments if available (CBS recommended)
+ segments <- readSegments() ## segments if available (CBS recommended)
# cat(' ..processing')
if (is.null(segments)) { ## segmentation using DNA-copy if needed (must then be installed)
@@ -171,10 +175,12 @@
segments$Value <- segments$Value-median(Log2$Value,na.rm=T) ## Median-centering
Log2$Value <- Log2$Value-median(Log2$Value,na.rm=T) ## Median-centering
+ #browser()
allRegions=NULL; #if ('allRegions.Rdata' %in% dir()) load('allRegions.Rdata')
if (is.null(allRegions)) allRegions <- makeRegions(Log2, alf, segments,matched=matched,min=30,allelePeaks=allelePeaks) ## Calculates necessary data for segments (all functions are in this file)
save(allRegions,file='allRegions.Rdata')
regs=NULL;# if ('shortRegions.Rdata' %in% dir()) load('shortRegions.Rdata')
+ #browser()
if (is.null(regs)) {
regs <- regsFromSegs(Log2,alf,segments,bin=bin,min=30,matched=matched,allelePeaks=allelePeaks) ## Calculates the same data for shortened segments
save(regs,file='shortRegions.Rdata')
@@ -211,6 +217,7 @@
}
# cat('..plotting.\n')
+ #browser()
OverviewPlot(regs$chr,regs$start,regs$end,regs$logs,regs$scores,hg18=hg18,
as.character(Log2$Chromosome),Log2$Start,Log2$Value,as.character(alf$Chromosome),alf$Start,alf$Value,
name=name,MAPD=sampleData$MAPD[i],MHOF=sampleData$MHOF[i])
@@ -391,69 +398,99 @@
}
###
-regsFromSegs <- function (Log2,alf, segments, bin=200,min=1,matched=F,allelePeaks=FALSE) {
+segMatch <- function(segments.txt,probes.txt) {
+ ## Returns a data frame with probes.startIx probes.endIx
+ #verry fast
+ colnames(segments.txt)[colnames(segments.txt)=='chr']='Chromosome'
+ colnames(segments.txt)[colnames(segments.txt)=='start']='Start'
+ colnames(segments.txt)[colnames(segments.txt)=='end']='End'
+
+ segStarts=data.frame(n.seg=1:nrow(segments.txt),n.pos=NA,chr=chrom_n(segments.txt$Chromosome),pos=segments.txt$Start,type='start')
+ segEnds=data.frame(n.seg=1:nrow(segments.txt),n.pos=NA,chr=chrom_n(segments.txt$Chromosome),pos=segments.txt$End,type='end')
+ lpos=data.frame(n.seg=NA,n.pos=1:nrow(probes.txt),chr=chrom_n(probes.txt$Chromosome),pos=(probes.txt$Start+probes.txt$End)/2,type='logr')
+ dummy1=data.frame(n.seg=NA,n.pos=NA,chr=-Inf,pos=-Inf,type='logr')
+ dummy2=data.frame(n.seg=NA,n.pos=NA,chr=Inf,pos=Inf,type='logr')
+
+ table=rbind(segStarts,segEnds,lpos,dummy1,dummy2)
+ table=table[order(table$chr,table$pos),]
+
+ start=(table$type=='start')
+ end=(table$type=='end')
+ #l=(table$type=='logr')
+
+ # omitting segment ends, let segment start "n.pos" be that of the next marker
+ table$n.pos[!end][start[!end]] <- table$n.pos[!end][which(start[!end])+1]
+
+ # omitting segment starts, let segment end "n.pos" be that of the previous marker
+ table$n.pos[!start][end[!start]] <- table$n.pos[!start][which(end[!start])-1]
+
+ # re-create a data frame with startIx and endIx (in probes.txt) per segment
+ indices=data.frame(startIx=table$n.pos[start],endIx=table$n.pos[end])
+
+ na=is.na(indices$startIx+indices$endIx)
+ na[!na] = indices$startIx[!na]>indices$endIx[!na]
+ indices$startIx[na] <- indices$endIx[na] <- NA
+
+ return(indices)
+}
+regsFromSegs <- function (Log2,alf, segments, bin=150,min=1,matched=F,allelePeaks=FALSE) {
## This function builds short segments and calcualtes their average Log-R and Allelic Imbalance Ratio.
rownames(Log2)=1:nrows(Log2)
rownames(alf)=1:nrows(alf)
- regs=list('chr'=NULL,'start'=NULL,'end'=NULL,'logs'=NULL,'scores'=NULL,'probes'=NULL,'snps'=NULL)
+ regs=NULL #list('chr'=NULL,'start'=NULL,'end'=NULL,'logs'=NULL,'scores'=NULL,'probes'=NULL,'snps'=NULL)
#,'key1'=rep(NA,nrow(Log2)),'key2'=rep(NA,nrow(alf)))
n=nrow(segments)
s_check=NULL
+
+ pMatch=segMatch(segments,Log2)
+ sMatch=segMatch(segments,alf)
+
for (c in 1:n) { ## for every segment
- tlog=Log2[Log2$Chromosome==segments$Chromosome[c],] ## Log-R on this chromosome
- talf=alf[alf$Chromosome==segments$Chromosome[c],] ## Allele Freq on this chromosome
- tlog=tlog[(tlog$Start>=segments$Start[c])&(tlog$Start<segments$End[c]),] ## Log-R on this segment
- talf=talf[(talf$Start>=segments$Start[c])&(talf$Start<segments$End[c]),] ## Allele Freq on this segment
+ #tlog=Log2[Log2$Chromosome==segments$Chromosome[c],] ## Log-R on this chromosome
+ #talf=alf[alf$Chromosome==segments$Chromosome[c],] ## Allele Freq on this chromosome
+ #tlog=tlog[(tlog$Start>=segments$Start[c])&(tlog$Start<segments$End[c]),] ## Log-R on this segment
+ #talf=talf[(talf$Start>=segments$Start[c])&(talf$Start<segments$End[c]),] ## Allele Freq on this segment
- tsnps=nrow(talf) ## number of snps and probes in this segment
- tprobes=nrow(tlog)
+ tsnps=sMatch[c,2]-sMatch[c,1]+1; if (is.na(tsnps)) tsnps=0 #nrow(talf) ## number of snps and probes in this segment
+ tprobes=pMatch[c,2]-pMatch[c,1]+1; if (is.na(tsnps)) tsnps=0 #nrow(tlog)
tnregs=max(trunc(tsnps/bin),1) ## Further split into how many (a minimum of 1) subsegments
tcuts=segments$Start[c] ## The first start pos
tlength=segments$End[c]-segments$Start[c] ## Length of this whole segment
for (i in 1:tnregs) tcuts = c(tcuts, tcuts[1]+round(i/tnregs*tlength)) ## Break the segment at these positions
- for (r in 1:(tnregs)) { ## build the subsegments
- regs$chr=c(regs$chr,as.character(segments$Chromosome[c])) ## Chromosome
- s_=tcuts[r] ## Start
- e_=tcuts[r+1] ## End
- thisalf=talf[(talf$Start>=s_)&(talf$Start<=e_),] ## get the Log-R values
- thislog=tlog[(tlog$Start>=s_)&(tlog$Start<=e_),] ## and the allele frequency
- #regs$key1[as.integer(rownames(thislog))]=length(regs$chr) ## store their positions for fast access during plotting
- #regs$key2[as.integer(rownames(thisalf))]=length(regs$chr) ## --"--
- regs$logs=c( regs$logs, mean(thislog$Value) ) ## store average log ratio of this segment
- regs$probes=c(regs$probes,nrow(thislog)) ## store number of probes
- regs$snps=c(regs$snps,nrow(thisalf)) ## store number of bi-allelic probes (SNPs)
- #regs$or_seg=c(regs$or_seg,c)
- regs$start=c(regs$start,s_) ## store start and end positions
- regs$end=c(regs$end,e_)
+ tab=data.frame(chr=segments$Chromosome[c],start=tcuts[-length(tcuts)],end=tcuts[-1],logs=NA,scores=NA,probes=0,snps=0)
+ if (is.null(regs)) regs=tab else regs=rbind(regs,tab)
+ }
+
+ pMatch=segMatch(regs,Log2)
+ sMatch=segMatch(regs,alf)
+
+ for (i in 1:nrow(regs)) { ## build the subsegments
+ #regs$chr=c(regs$chr,as.character(segments$Chromosome[c])) ## Chromosome
+ #s_=tcuts[r] ## Start
+ #e_=tcuts[r+1] ## End
+ #thisalf=talf[(talf$Start>=s_)&(talf$Start<=e_),] ## get the Log-R values
+ #thislog=tlog[(tlog$Start>=s_)&(tlog$Start<=e_),] ## and the allele frequency
+ #regs$key1[as.integer(rownames(thislog))]=length(regs$chr) ## store their positions for fast access during plotting
+ #regs$key2[as.integer(rownames(thisalf))]=length(regs$chr) ## --"--
+ if (!is.na(pMatch[i,1])) {
+ regs$logs[i]=median(Log2$Value[pMatch[i,1]:pMatch[i,2]],na.rm=T) ## store average log ratio of this segment
+ regs$probes[i]=pMatch[i,2]-pMatch[i,1]+1 ## store number of probes
+ }
+ if (!is.na(sMatch[i,1])) {
temp = NA
- try(temp <- allelicImbalance(thisalf$Value,min,matched,allelePeaks),silent=T)
- regs$scores =c(regs$scores,temp)
- # if (nrow(thisalf)>min) { ## Time to calculate Allelic Imbalance Ratio (if enough SNPs)
- # t1=sort( abs(thisalf$Value-0.5) ) ## distance from middle (het) in the allele freq pattern, ascending
- # if (length(unique(t1))>3) { ## do not attempt clustering with too few snps
- # xx=NULL
- # try(xx <- kmeans(t1, 2),silent=T) ## Attempt k-means (Hartigan-Wong: has proven very stable)
- # if (!is.null(xx)) if (min(xx$size) > 0.05*max(xx$size)) { ## On some occations data quality is poor, requiring 5%+ heterozygous SNPs avoids most such cases.
- # xx=xx$centers
- # } else xx=NA
- # } else xx=NA
- # } else xx=NA
- ##try (if (is.na(xx)) xx=0:1, silent=T)
- # try (if (length(xx)==0) xx=0:1, silent=T)
- # regs$scores=c(regs$scores, min(xx)/max(xx) ) ## Allelic Imbalance Ratio = inner / outer cluster.
- # regs$het=c(regs$het, min(xx)) ## $het and $hom are no longer in use.
- # regs$hom=c(regs$hom, max(xx))
- # if(matched == T | is.na(regs$scores[length(regs$scores)])) {
- # regs$scores[length(regs$scores)] <- 2*median(abs(thisalf$Value-.5),na.rm=T)
- # }
+ try(temp <- allelicImbalance(alf$Value[sMatch[i,1]:sMatch[i,2]],min,matched,allelePeaks),silent=T)
+ regs$scores[i]=temp
+ regs$snps[i]=sMatch[i,2]-sMatch[i,1]+1
}
}
- regs=as.data.frame(regs)
+
+ #regs=as.data.frame(regs)
regs=regs[!is.na(regs$logs),] ### MODDAT MARKUS MAJ 2013
return (regs)
}
allelicImbalance <- function (data,min=30,matched=F,allelePeaks=F) {
+ #browser()
if(matched == T ) {
return(2*median(abs(data-.5),na.rm=T))
}
@@ -604,8 +641,19 @@
return (segments)
}
###
+chrom_n <- function(data) {
+ out=rep(Inf,length(data))
+ for (c in c(1:22)) {
+ out[data==paste('chr',c, sep='')]=c
+ }
+ out[data=='chrX']=23
+ out[data=='chrY']=24
+ out[data=='chrM']=25
+ return(out)
+}
makeRegions <- function(Log2, alf, segments,dataType='Nexus',min=30,matched=FALSE,allelePeaks=F) {
## makeRegions is similar to "regsfromsegs" except regions are not subdivided before calculation of mean Log-R and Allelic Imbalance Ratio.
+ #browser()
regions=segments
regions$Chromosome=as.character(segments$Chromosome) ## Chromosome
regions$lengthMB=round((regions$End-regions$Start)/1000000,3) ## length in megabases
@@ -616,19 +664,31 @@
regionIx=NULL ## Not currently used
regionIx$Log2 <- list()
regionIx$alf <- list()
+
+ pMatch=segMatch(regions,Log2)
+ na=is.na(pMatch[,1])
+ regions=regions[!na,]
+ pMatch=pMatch[!na,]
+ sMatch=segMatch(regions,alf)
+
for (i in 1:nrows(regions)) {
- log2temp=which(equals(Log2$Chromosome,regions$Chromosome[i])) ## index of Log-R (current chrom)
- alftemp=which(equals(alf$Chromosome,regions$Chromosome[i])) ## index of Allele frequency (current chrom)
+ #log2temp=which(equals(Log2$Chromosome,regions$Chromosome[i])) ## index of Log-R (current chrom)
+ #alftemp=which(equals(alf$Chromosome,regions$Chromosome[i])) ## index of Allele frequency (current chrom)
- log2temp=log2temp [Log2$Start[log2temp]>=regions$Start[i] & Log2$Start[log2temp]<regions$End[i]] ## index of Log-R (current segment)
- alftemp=alftemp [alf$Start[alftemp]>=regions$Start[i] & alf$Start[alftemp]<regions$End[i]] ## index of Allele frequency (current segment)
+ #log2temp=log2temp [Log2$Start[log2temp]>=regions$Start[i] & Log2$Start[log2temp]<regions$End[i]] ## index of Log-R (current segment)
+ #alftemp=alftemp [alf$Start[alftemp]>=regions$Start[i] & alf$Start[alftemp]<regions$End[i]] ## index of Allele frequency (current segment)
+
+ log2temp=Log2$Value[pMatch[i,1]:pMatch[i,2]] ## index of Log-R (current segment)
+ alftemp=numeric(0)
+ if (!is.na(sMatch[i,1])) alftemp=alf$Value[sMatch[i,1]:sMatch[i,2]] ## index of Allele frequency (current segment)
+
regions$probes[i]=length(log2temp)
regions$snps[i]=length(alftemp)
- regionIx$Log2[[i]]=log2temp ## indexes of Log-R and Allele Frequency are saved for future use (plot color coding)
- regionIx$alf[[i]]=alftemp
- log2temp=Log2$Value[log2temp]
- regions$log2[i]=median(log2temp)
- alftemp=alf$Value[alftemp]
+ #regionIx$Log2[[i]]=log2temp ## indexes of Log-R and Allele Frequency are saved for future use (plot color coding)
+ #regionIx$alf[[i]]=alftemp
+ #log2temp=Log2$Value[log2temp]
+ regions$log2[i]=median(log2temp, na.rm=T)
+ #alftemp=alf$Value[alftemp]
temp = NA
try(temp <- allelicImbalance(alftemp,min,matched=matched,allelePeaks=allelePeaks),silent=T)
regions$imba[i] <- temp
@@ -1243,6 +1303,7 @@
## Pileups
pileup <- function(regions) { #, nsamples) {
+ if (nrow(regions)==0) return(data.frame(Chromosome=integer(0),Start=integer(0),End=integer(0),Count=integer(0)))
chrs <- unique(as.character(regions$Chromosome))
if (length(chrs)>1) { # if more than one chrom, run them separately.
pile=NULL
@@ -1268,6 +1329,7 @@
pileup_dif <- function(regions1,regions2) {
chrs <- unique(c(as.character(regions1$Chr),as.character(regions2$Chr)))
+ if (length(chrs)==0) return(data.frame(Chromosome=integer(0),Start=integer(0),End=integer(0),Count1=integer(0),Count2=integer(0)))
if (length(chrs)>1) { # if more than one chrom, run them separately.
pile=NULL
for (i in 1:length(chrs)) pile <- rbind(pile,pileup_dif(regions1[as.character(regions1$Chr)==chrs[i],], regions2[as.character(regions2$Chr)==chrs[i],]))
@@ -1308,6 +1370,9 @@
pvals <- function(regions,tot1,tot2) {
+ if (nrow(regions)==0) return(data.frame(Chromosome=integer(0),Start=integer(0),End=integer(0),Count1=integer(0),Count2=integer(0),
+ percent1=integer(0),precent2=integer(0),precent=integer(0),p.value=integer(0),
+ odds.ratio=integer(0),conf_low=integer(0),conf_high=integer(0)))
p <- rep(1,nrow(regions))
or <- rep(1,nrow(regions))
conf_low <- rep(NA,nrow(regions))
@@ -1575,18 +1640,19 @@
#mtext(text="Allelic imbalance",side=2,line=4,las=3,cex=2)
#mtext(text="Coverage, all segments",side=1,line=4,cex=2,las=1)
- data <- pile
- data$e <- 0 -> data$s
- for (i in 1:23) data$s[data$Chromosome==as.character(chroms$Chr[i])] <- chroms$before[i] -> data$e[data$Chromosome==as.character(chroms$Chr[i])]
- data$s <- data$s+data$Start
- data$e <- data$e+data$End
- rect(
- xleft=data$s, xright=data$e,
- ybottom=0,ytop=data$Percent,
- col=color,
- border=NA,
- )
-
+ if (nrow(pile)>0) {
+ data <- pile
+ data$e <- 0 -> data$s
+ for (i in 1:23) data$s[data$Chromosome==as.character(chroms$Chr[i])] <- chroms$before[i] -> data$e[data$Chromosome==as.character(chroms$Chr[i])]
+ data$s <- data$s+data$Start
+ data$e <- data$e+data$End
+ rect(
+ xleft=data$s, xright=data$e,
+ ybottom=0,ytop=data$Percent,
+ col=color,
+ border=NA,
+ )
+ }
## The chrom labels
text(
x=chroms$before+chroms$length/2,
@@ -1755,71 +1821,75 @@
)
#mtext(text="Allelic imbalance",side=2,line=4,las=3,cex=2)
#mtext(text="Coverage, all segments",side=1,line=4,cex=2,las=1)
+ mhcolor=colorRampPalette(c("#FFFFFF",color),space="rgb")(5)[3]
- ## Group 1
- data <- pile1
- data$e <- 0 -> data$s
- for (i in 1:23) data$s[data$Chromosome==as.character(chroms$Chr[i])] <- chroms$before[i] -> data$e[data$Chromosome==as.character(chroms$Chr[i])]
- data$s <- data$s+data$Start
- data$e <- data$e+data$End
- rect(
- xleft=data$s, xright=data$e,
- ybottom=0,ytop=data$Percent,
- col=color,
- border=NA,
- )
-
- ## Group 2
- data <- pile2
- data$e <- 0 -> data$s
- for (i in 1:23) data$s[data$Chromosome==as.character(chroms$Chr[i])] <- chroms$before[i] -> data$e[data$Chromosome==as.character(chroms$Chr[i])]
- data$s <- data$s+data$Start
- data$e <- data$e+data$End
- rect(
- xleft=data$s, xright=data$e,
- ybottom=0,ytop= -data$Percent,
- col=color,
- border=NA,
- )
-
- ### The difference
- data <- dif
- data$e <- 0 -> data$s
- for (i in 1:23) data$s[data$Chromosome==as.character(chroms$Chr[i])] <- chroms$before[i] -> data$e[data$Chromosome==as.character(chroms$Chr[i])]
- data$s <- data$s+data$Start
- data$e <- data$e+data$End
- data <- data[order(data$percent,decreasing=T),]
-
- #difcolor <- colorRampPalette(c("#FFFFFF",color),space="rgb")(4)[2]
-
- rect(
- xleft=data$s, xright=data$e,
- ybottom=0,ytop=data$percent,
- #ybottom=data$conf_low,ytop=data$conf_high,
- col=difcolor,
- border=NA
- )
-
- ## The significant
- mhcolor=colorRampPalette(c("#FFFFFF",color),space="rgb")(5)[3]
- data=data[data$p.value<p_cutoff & abs(data$percent)>freq_cutoff,]
- if (nrow(data)>0) {
+ if (nrow(pile1)>0) {
+ ## Group 1
+ data <- pile1
+ data$e <- 0 -> data$s
+ for (i in 1:23) data$s[data$Chromosome==as.character(chroms$Chr[i])] <- chroms$before[i] -> data$e[data$Chromosome==as.character(chroms$Chr[i])]
+ data$s <- data$s+data$Start
+ data$e <- data$e+data$End
rect(
xleft=data$s, xright=data$e,
- ybottom= -100,ytop=-95,
- col='#000000',
+ ybottom=0,ytop=data$Percent,
+ col=color,
+ border=NA,
+ )
+ }
+ if (nrow(pile2)>0) {
+ ## Group 2
+ data <- pile2
+ data$e <- 0 -> data$s
+ for (i in 1:23) data$s[data$Chromosome==as.character(chroms$Chr[i])] <- chroms$before[i] -> data$e[data$Chromosome==as.character(chroms$Chr[i])]
+ data$s <- data$s+data$Start
+ data$e <- data$e+data$End
+ rect(
+ xleft=data$s, xright=data$e,
+ ybottom=0,ytop= -data$Percent,
+ col=color,
+ border=NA,
+ )
+ }
+ if (nrow(dif)>0) {
+ ### The difference
+ data <- dif
+ data$e <- 0 -> data$s
+ for (i in 1:23) data$s[data$Chromosome==as.character(chroms$Chr[i])] <- chroms$before[i] -> data$e[data$Chromosome==as.character(chroms$Chr[i])]
+ data$s <- data$s+data$Start
+ data$e <- data$e+data$End
+ data <- data[order(data$percent,decreasing=T),]
+
+ #difcolor <- colorRampPalette(c("#FFFFFF",color),space="rgb")(4)[2]
+
+ rect(
+ xleft=data$s, xright=data$e,
+ ybottom=0,ytop=data$percent,
+ #ybottom=data$conf_low,ytop=data$conf_high,
+ col=difcolor,
border=NA
)
+
+ ## The significant
+ data=data[data$p.value<p_cutoff & abs(data$percent)>freq_cutoff,]
+ if (nrow(data)>0) {
+ rect(
+ xleft=data$s, xright=data$e,
+ ybottom= -100,ytop=-95,
+ col='#000000',
+ border=NA
+ )
+ }
+ ## The OR
+ #for (i in 1:23) {
+ # temp=data[data$Chr==as.character(chroms$chr[i]),]
+ # if (nrow(temp)>5) {
+ # temp$odds.ratio[temp$odds.ratio<1 & temp$odds.ratio>0]=1/temp$odds.ratio[temp$odds.ratio<1 & temp$odds.ratio>0]
+ # temp=temp[order(temp$odds.ratio,decreasing=T),][1,]
+ # text(x=chroms$before[i]+mean(temp$Start,temp$End), y=-105, cex=1.3, label=round(temp$odds.ratio,1))
+ # }
+ #}
}
- ## The OR
- #for (i in 1:23) {
- # temp=data[data$Chr==as.character(chroms$chr[i]),]
- # if (nrow(temp)>5) {
- # temp$odds.ratio[temp$odds.ratio<1 & temp$odds.ratio>0]=1/temp$odds.ratio[temp$odds.ratio<1 & temp$odds.ratio>0]
- # temp=temp[order(temp$odds.ratio,decreasing=T),][1,]
- # text(x=chroms$before[i]+mean(temp$Start,temp$End), y=-105, cex=1.3, label=round(temp$odds.ratio,1))
- # }
- #}
axis(2,at=c(-97),#,-105),
labels=c('sign.'),#, 'BestOR'),
pos=0,
@@ -1855,6 +1925,12 @@
col='#00000070',
lwd=1,lty=3
)
+ segments(
+ x0=0,x1=sum(chroms$length),
+ y0=0,y1=0,
+ col='#00000090',
+ lwd=1,lty=1
+ )
## Y axis ticks
axis(2,at=seq(-80,80,20),
@@ -1890,7 +1966,6 @@
OverviewPlot <- function(chr,start,end,int,ai,hg18,mchr,mpos,mval,schr,spos,sval,name='',xlim=c(-1,1.5),ylim=c(0,1),MAPD,MHOF)
{
-
if(hg18==T)
{
chroms=chroms_hg18
@@ -3887,7 +3962,7 @@
write.table(x=newSampleData[newSampleData$Sample == sample,],file=paste(root,'backup.csv',sep='/'),sep='\t',row.names=F,append=T,col.names=F)
return(newSampleData[newSampleData$Sample == sample,])
} else if(answer == 29) { #Get a browser()
- browser()
+ #browser()
}
} else { #If pre-existing values for the sample exists, grab the values.
More information about the Patchwork-commits
mailing list