From 91824f478e962ffe66ecea663aebfe325b966db2 Mon Sep 17 00:00:00 2001 From: depristo Date: Thu, 13 Jan 2011 15:16:06 +0000 Subject: [PATCH] FASTQ directory is gone git-svn-id: file:///humgen/gsa-scr1/gsa-engineering/svn_contents/trunk@4986 348d0f76-0448-11de-a6fe-93d51630548a --- R/privateMutations.R | 45 +++++++- build.xml | 6 +- .../threading/GenomeLocProcessingTracker.java | 104 ++++++++++++++++++ ...haredMemoryGenomeLocProcessingTracker.java | 31 ++++++ 4 files changed, 176 insertions(+), 10 deletions(-) create mode 100644 java/src/org/broadinstitute/sting/utils/threading/GenomeLocProcessingTracker.java create mode 100644 java/src/org/broadinstitute/sting/utils/threading/SharedMemoryGenomeLocProcessingTracker.java diff --git a/R/privateMutations.R b/R/privateMutations.R index a23e28247..2fe1f8ab3 100644 --- a/R/privateMutations.R +++ b/R/privateMutations.R @@ -16,13 +16,46 @@ f <- function(d, acs) { legend("topright", legend=lapply(acs, function(x) paste("AC =", x)), fill=cols, title="Sub-population") } -expected <- function(ps, N, eps) { - co = 2 * N / ( 1 - eps ) - v = co * ((1 - ps)/(1-eps))^(2 * N - 1) - v / sum(v) +expected <- function(maxAN, N, eps, ac1scale = F) { + scale = 10 + + f <- function(ps, N) { + co = 2 * N / ( 1 - eps ) + co * ((1 - ps)/(1-eps))^(2 * N - 1) + } + + # these are the points that we'll actually show, but we need to do the calculation + # special for the AC = 1 given the equation actually fits an infinite population + # not a discrete population with max chromosomes + ps = 1:maxAN / maxAN + v = f(ps, N) + v = v / sum(v) + + if ( ac1scale ) { + subps = seq(1, maxAN*scale) / (maxAN * scale) + #print(subps) + subv = f(subps, N) + #print(subv) + #print(v[1:10]) + pBelowAC1 = sum(subv[1:scale] / sum(subv)) + #print(list(pBelowAC1=pBelowAC1, v1=v[1])) + v[1] = v[1] + pBelowAC1 + } + + list(ps = ps, pr = v) } f(afs, c(1,2,3,5,10,50)) -x = 1:MAX_AC / 200000 -points(x, expected(x,1000,1e-8),type="l",lty=3,lwd=3) + +if ( F ) { +scale = 100 +ex1 = expected(200000, 1000, 1e-8) +ex2 = expected(200000*scale, 1000, 1e-8) +i = 1:(200000*scale) %% scale == 1 +plot(ex2$ps[i], cumsum(ex1$pr), type="l",lty=3,lwd=3, log="x", col="red") +points(ex2$ps[i], cumsum(ex2$pr)[i], type="l",lty=3,lwd=3, log="x") +} + +ex = expected(200000, 1000, 1e-8, T) +points(ex$ps, ex$pr, type="l",lty=3,lwd=3) diff --git a/build.xml b/build.xml index ee7cb2851..ee6c4debc 100644 --- a/build.xml +++ b/build.xml @@ -519,10 +519,8 @@ - + + diff --git a/java/src/org/broadinstitute/sting/utils/threading/GenomeLocProcessingTracker.java b/java/src/org/broadinstitute/sting/utils/threading/GenomeLocProcessingTracker.java new file mode 100644 index 000000000..4d3507e4c --- /dev/null +++ b/java/src/org/broadinstitute/sting/utils/threading/GenomeLocProcessingTracker.java @@ -0,0 +1,104 @@ +package org.broadinstitute.sting.utils.threading; + +import org.broadinstitute.sting.utils.GenomeLoc; + +import java.util.List; + +/** + * Created by IntelliJ IDEA. + * User: depristo + * Date: 1/13/11 + * Time: 9:38 AM + * To change this template use File | Settings | File Templates. + */ +public abstract class GenomeLocProcessingTracker { + /** + * Information about processing locations and their owners + */ + public static final class ProcessingLoc { + GenomeLoc loc; + String owner; + + /** + * Create an unowned loc + * @param loc + */ + public ProcessingLoc(GenomeLoc loc) { + this(loc, null); + + } + + /** + * Create a loc that's already owned + * @param loc + * @param owner + */ + public ProcessingLoc(GenomeLoc loc, String owner) { + this.loc = loc; + this.owner = owner; + } + + public GenomeLoc getLoc() { + return loc; + } + + public String getOwner() { + return owner; + } + + public void setOwner(String owner) { + this.owner = owner; + } + + public boolean isOwned() { return owner != null; } + } + + // -------------------------------------------------------------------------------- + // + // Code to claim intervals for processing and query for their ownership + // + // -------------------------------------------------------------------------------- + /** + * Queries the current database if a location is owned. Does not guarantee that the + * loc can be owned in a future call, though. + * + * @param loc + * @return + */ + public boolean locIsOwned(GenomeLoc loc) { + return findOwner(loc) != null; + } + + public ProcessingLoc findOwner(GenomeLoc loc) { + for ( ProcessingLoc l : getProcessingLocs() ) { + if ( l.getLoc().equals(loc) ) + return l; + } + + return null; + } + + /** + * The workhorse routine. Attempt to claim processing ownership of loc, with my name. + * This is an atomic operation -- other threads / processes will wait until this function + * returns. The return result is the ProcessingLoc object describing who owns this + * location. If the location isn't already claimed and we now own the location, the pl owner + * will be myName. Otherwise, the name of the owner can found in the pl. + * + * @param loc + * @param myName + * @return + */ + public abstract ProcessingLoc claimOwnership(GenomeLoc loc, String myName); + + /** + * Returns the list of currently owned locations, updating the database as necessary. + * DO NOT MODIFY THIS LIST! As with all parallelizing data structures, the list may be + * out of date immediately after the call returns, or may be updating on the fly. + * + * This is really useful for printing, counting, etc. operations that aren't mission critical + * + * @return + */ + protected abstract List getProcessingLocs(); +} diff --git a/java/src/org/broadinstitute/sting/utils/threading/SharedMemoryGenomeLocProcessingTracker.java b/java/src/org/broadinstitute/sting/utils/threading/SharedMemoryGenomeLocProcessingTracker.java new file mode 100644 index 000000000..a57921c4a --- /dev/null +++ b/java/src/org/broadinstitute/sting/utils/threading/SharedMemoryGenomeLocProcessingTracker.java @@ -0,0 +1,31 @@ +package org.broadinstitute.sting.utils.threading; + +import org.broadinstitute.sting.utils.GenomeLoc; + +import java.util.ArrayList; +import java.util.List; + +/** + * For algorithmic testing purposes only. Uses synchronization to keep a consistent + * processing list in shared memory. + */ +public class SharedMemoryGenomeLocProcessingTracker extends GenomeLocProcessingTracker { + List processingLocs = new ArrayList(); + + public synchronized ProcessingLoc claimOwnership(GenomeLoc loc, String myName) { + // processingLocs is a shared memory synchronized object, and this + // method is synchonized, so we can just do our processing + ProcessingLoc owner = super.findOwner(loc); + + if ( owner == null ) { // we are unowned + owner = new ProcessingLoc(loc, myName); + processingLocs.add(owner); + } + + return owner; + } + + protected synchronized List getProcessingLocs() { + return processingLocs; + } +}