Merge pull request #19 from sdentro/dev

Merging v2.2.2

Author: sdentro

DESCRIPTION

@@ -3,7 +3,7 @@ Maintainer: Stefan Dentro <[email protected]>
 License: GPL-3
 Type: Package
 Title: Battenberg subclonal copy number caller
-Version: 2.2.1
+Version: 2.2.2
 Authors@R: c(person("David", "Wedge", role=c("aut"), email="[email protected]"),
     person("Peter", "Van Loo", role=c("aut")),
     person("Stefan","Dentro", email="[email protected]", role=c("aut", "cre")),

R/fitcopynumber.R

@@ -273,7 +273,7 @@ callSubclones = function(sample.name, baf.segmented.file, logr.file, rho.psi.fil
     #if no points to plot, skip
     if (length(pos)==0) { next }
 
-    if (!is.null(sv_breakpoints_file)) {
+    if (!is.null(sv_breakpoints_file) & !sv_breakpoints_file=="NA") {
       svs_pos = svs[svs$chromosome==chr,]$position / 1000000
     } else {
       svs_pos = NULL
@@ -315,7 +315,7 @@ callSubclones = function(sample.name, baf.segmented.file, logr.file, rho.psi.fil
   is_subclonal_min[is.na(is_subclonal_min)] = F
   segment_states_min = subclones$nMin1_A * ifelse(is_subclonal_min, subclones$frac1_A, 1)  + ifelse(is_subclonal_min, subclones$nMin2_A, 0) * ifelse(is_subclonal_min, subclones$frac2_A, 0) 
   segment_states_maj = subclones$nMaj1_A * ifelse(is_subclonal_maj, subclones$frac1_A, 1)  + ifelse(is_subclonal_maj, subclones$nMaj2_A, 0) * ifelse(is_subclonal_maj, subclones$frac2_A, 0) 
-  ploidy = sum((segment_states_min+segment_states_maj) * seg_length) / sum(seg_length)
+  ploidy = sum((segment_states_min+segment_states_maj) * seg_length, na.rm=T) / sum(seg_length, na.rm=T)
 
   # Plot genome wide figures
   plot.gw.subclonal.cn(subclones=subclones, BAFvals=BAFvals, rho=rho, ploidy=ploidy, goodness=goodness, output.gw.figures.prefix=output.gw.figures.prefix, chr.names=chr_names)

R/prepare_wgs.R

@@ -288,8 +288,23 @@ gc.correct.wgs = function(Tumour_LogR_file, outfile, correlations_outfile, gc_co
 
   GCcorrected = Tumor_LogR[!flag_NA,]
   GCcorrected[,3] = model$residuals
-    
+  rm(Tumor_LogR)
+  
   corr = data.frame(windowsize=names(corr), correlation=corr)
-  write.table(corr, file=correlations_outfile, sep="\t", quote=F, row.names=F)
+  write.table(corr, file=gsub(".txt", "_beforeCorrection.txt", correlations_outfile), sep="\t", quote=F, row.names=F)
   write.table(GCcorrected, file=outfile, sep="\t", quote=F, row.names=F)
+  
+  # Recalculate the correlations to see how much there is left
+  snps_gccorrected = paste(GCcorrected[,1], GCcorrected[,2], sep="_")
+  snps_gcdata = paste(GC_data[,1], GC_data[,2], sep="_")
+  snps_intersect = intersect(snps_gccorrected, snps_gcdata)
+  GCcorrected = GCcorrected[snps_gccorrected %in% snps_intersect, ]
+  GC_data = GC_data[snps_gcdata %in% snps_intersect, ]
+  
+  flag_nona = is.finite(GCcorrected[,3])
+  corr = cor(GC_data[flag_nona, 3:ncol(GC_data)], GCcorrected[flag_nona,3], use="complete.obs")
+  length = nrow(GCcorrected)
+  corr = apply(corr, 1, function(x) sum(abs(x*length))/sum(length))
+  corr = data.frame(windowsize=names(corr), correlation=corr)
+  write.table(corr, file=gsub(".txt", "_afterCorrection.txt", correlations_outfile), sep="\t", quote=F, row.names=F)
 }

man/callSubclones.Rd

@@ -6,7 +6,7 @@
 \usage{
 callSubclones(sample.name, baf.segmented.file, logr.file, rho.psi.file,
   output.file, output.figures.prefix, output.gw.figures.prefix, chr_names,
-  masking_output_file, max_allowed_state = 100, sv_breakpoints_file = NULL,
+  masking_output_file, max_allowed_state = 250, sv_breakpoints_file = NULL,
   gamma = 1, segmentation.gamma = NA, siglevel = 0.05, maxdist = 0.01,
   noperms = 1000, seed = as.integer(Sys.time()))
 }

man/determine_copynumber.Rd

@@ -0,0 +1,30 @@
+% Generated by roxygen2 (4.1.1): do not edit by hand
+% Please edit documentation in R/fitcopynumber.R
+\name{determine_copynumber}
+\alias{determine_copynumber}
+\title{Given all the determined values make a copy number call for each segment}
+\usage{
+determine_copynumber(BAFvals, LogRvals, rho, psi, gamma, ctrans, ctrans.logR,
+  maxdist, siglevel, noperms)
+}
+\arguments{
+\item{BAFvals}{BAFsegmented data.frame with 5 columns}
+
+\item{LogRvals}{Raw logR values in data.frame with 3 columns}
+
+\item{rho}{Optimal rho value, the choosen cellularity}
+
+\item{psi}{Optimal psi value, the choosen ploidy}
+
+\item{gamma}{Platform gamma parameter}
+}
+\value{
+A data.frame with copy number determined for each segment
+}
+\description{
+Given all the determined values make a copy number call for each segment
+}
+\author{
+dw9
+}
+

man/mask_high_cn_segments.Rd

@@ -0,0 +1,25 @@
+% Generated by roxygen2 (4.1.1): do not edit by hand
+% Please edit documentation in R/fitcopynumber.R
+\name{mask_high_cn_segments}
+\alias{mask_high_cn_segments}
+\title{Mask segments that have a too high CN state}
+\usage{
+mask_high_cn_segments(subclones, bafsegmented, max_allowed_state)
+}
+\arguments{
+\item{subclones}{}
+
+\item{bafsegmented}{}
+
+\item{max_allowed_state}{}
+}
+\value{
+A list with the masked subclones, bafsegmented and the number of segments masked and their total genome size
+}
+\description{
+Mask segments that have a too high CN state
+}
+\author{
+sd11
+}
+