Search in sources :

Example 16 with ProbeSet

use of uk.ac.babraham.SeqMonk.DataTypes.Probes.ProbeSet in project SeqMonk by s-andrews.

the class HeatmapMatrix method findCommonProbeListParent.

@SuppressWarnings({ "rawtypes", "unchecked" })
private ProbeList findCommonProbeListParent() {
    if (probeLists.length == 1) {
        return probeLists[0];
    } else {
        HashSet[] parents = new HashSet[probeLists.length - 1];
        for (int i = 1; i < probeLists.length; i++) {
            parents[i - 1] = new HashSet<ProbeList>();
            ProbeList currentList = probeLists[i];
            parents[i - 1].add(currentList);
            while (!(currentList instanceof ProbeSet)) {
                parents[i - 1].add(currentList.parent());
                currentList = currentList.parent();
            }
        }
        // Now we go through the first list to find the common parent
        ProbeList firstList = probeLists[0];
        boolean notFound = false;
        while (true) {
            notFound = false;
            for (int i = 0; i < parents.length; i++) {
                if (!parents[i].contains(firstList)) {
                    notFound = true;
                }
            }
            if (notFound) {
                firstList = firstList.parent();
            } else {
                return firstList;
            }
        }
    }
}
Also used : ProbeSet(uk.ac.babraham.SeqMonk.DataTypes.Probes.ProbeSet) ProbeList(uk.ac.babraham.SeqMonk.DataTypes.Probes.ProbeList) HashSet(java.util.HashSet)

Example 17 with ProbeSet

use of uk.ac.babraham.SeqMonk.DataTypes.Probes.ProbeSet in project SeqMonk by s-andrews.

the class RunningWindowProbeGenerator method run.

/* (non-Javadoc)
	 * @see java.lang.Runnable#run()
	 */
public void run() {
    ProbeSet finalSet;
    if (designWithinExisting) {
        finalSet = designPerProbe();
    } else {
        finalSet = designPerChromosome();
    }
    generationComplete(finalSet);
}
Also used : ProbeSet(uk.ac.babraham.SeqMonk.DataTypes.Probes.ProbeSet)

Example 18 with ProbeSet

use of uk.ac.babraham.SeqMonk.DataTypes.Probes.ProbeSet in project SeqMonk by s-andrews.

the class SeqMonkParser method parseProbes.

/**
 * Parses the list of probes.
 *
 * @param sections The tab split initial line from the probes section
 * @throws SeqMonkException
 * @throws IOException Signals that an I/O exception has occurred.
 */
private void parseProbes(String[] sections) throws SeqMonkException, IOException {
    if (sections.length < 2) {
        throw new SeqMonkException("Probe line didn't contain at least 2 sections");
    }
    if (!sections[0].equals("Probes")) {
        throw new SeqMonkException("Couldn't find expected probes line");
    }
    int n = Integer.parseInt(sections[1]);
    probes = new Probe[n];
    String description = "No generator description available";
    if (sections.length > 2) {
        description = sections[2];
    }
    ProbeSet probeSet = new ProbeSet(description, n);
    if (sections.length > 3) {
        if (sections[3].length() > 0) {
            probeSet.setCurrentQuantitation(sections[3]);
        }
    }
    if (sections.length > 4) {
        probeSet.setComments(sections[4].replaceAll("`", "\n"));
    }
    // We need to save the probeset to the dataset at this point so we can add the probe
    // lists as we get to them.
    application.dataCollection().setProbeSet(probeSet);
    int positionOffset;
    // We used to store chr start and end
    if (thisDataVersion < 8) {
        positionOffset = 4;
    } else // We now store chr and packed position (to give start end and strand)
    {
        positionOffset = 3;
    }
    int expectedSectionLength = 3 + dataSets.length + dataGroups.length;
    String line;
    for (int i = 0; i < n; i++) {
        line = br.readLine();
        if (line == null) {
            throw new SeqMonkException("Ran out of probe data at line " + i + " (expected " + n + " probes)");
        }
        // Since the probes section can have blank trailing sections we need
        // to not trim these, hence the -1 limit.
        sections = line.split("\\t", -1);
        if (i == 0) {
            /*
				 * Older versions of this format put down data for just
				 * datasets.  Newer versions include data for datagroups
				 * as well.  We need to figure out which one we're looking
				 * at
				 */
            if (sections.length == positionOffset + dataSets.length) {
                expectedSectionLength = positionOffset + dataSets.length;
            } else if (sections.length == positionOffset + dataSets.length + dataGroups.length) {
                expectedSectionLength = positionOffset + dataSets.length + dataGroups.length;
            }
        }
        if (sections.length != expectedSectionLength) {
            throw new SeqMonkException("Expected " + expectedSectionLength + " sections in data file for " + sections[0] + " but got " + sections.length);
        }
        if (i % 10000 == 0) {
            Enumeration<ProgressListener> e = listeners.elements();
            while (e.hasMoreElements()) {
                e.nextElement().progressUpdated("Processed data for " + i + " probes", i, n);
            }
        }
        Chromosome c = application.dataCollection().genome().getChromosome(sections[1]).chromosome();
        // Sanity check
        if (c == null) {
            throw new SeqMonkException("Couldn't find a chromosome called " + sections[1]);
        }
        Probe p;
        if (thisDataVersion < 8) {
            int start = Integer.parseInt(sections[2]);
            int end = Integer.parseInt(sections[3]);
            p = new Probe(c, start, end);
        } else {
            long packedValue = Long.parseLong(sections[2]);
            p = new Probe(c, packedValue);
        }
        if (!sections[0].equals("null")) {
            p.setName(sections[0]);
        }
        probes[i] = p;
        probeSet.addProbe(probes[i], null);
        for (int j = positionOffset; j < sections.length; j++) {
            if (sections[j].length() == 0)
                continue;
            if ((j - positionOffset) >= dataSets.length) {
                dataGroups[j - (positionOffset + dataSets.length)].setValueForProbe(p, Float.parseFloat(sections[j]));
            } else {
                dataSets[j - positionOffset].setValueForProbe(p, Float.parseFloat(sections[j]));
            }
        }
    }
    application.dataCollection().activeProbeListChanged(probeSet);
    // This rename doesn't actually change the name.  We put this in because
    // the All Probes group is drawn in the data view before probes have been
    // added to it.  This means that it's name isn't updated when the probes have
    // been added and it appears labelled with 0 probes.  This doesn't happen if
    // there are any probe lists under all probes as they cause it to be refreshed,
    // but if you only have the probe set then you need this to make the display show
    // the correct information.
    probeSet.setName("All Probes");
}
Also used : ProbeSet(uk.ac.babraham.SeqMonk.DataTypes.Probes.ProbeSet) ProgressListener(uk.ac.babraham.SeqMonk.DataTypes.ProgressListener) Chromosome(uk.ac.babraham.SeqMonk.DataTypes.Genome.Chromosome) SeqMonkException(uk.ac.babraham.SeqMonk.SeqMonkException) Probe(uk.ac.babraham.SeqMonk.DataTypes.Probes.Probe)

Example 19 with ProbeSet

use of uk.ac.babraham.SeqMonk.DataTypes.Probes.ProbeSet in project SeqMonk by s-andrews.

the class ProbeSetTreeModel method probeListRemoved.

/* (non-Javadoc)
	 * @see uk.ac.babraham.SeqMonk.DataTypes.Probes.ProbeSetChangeListener#probeListRemoved(uk.ac.babraham.SeqMonk.DataTypes.Probes.ProbeList)
	 */
public void probeListRemoved(ProbeList l) {
    TreeModelEvent me;
    if (l instanceof ProbeSet) {
        me = new TreeModelEvent(l, getPathToRoot(l.parent()), new int[] { 0 }, new ProbeList[] { l });
    } else {
        me = new TreeModelEvent(l, getPathToRoot(l.parent()), new int[] { getIndexOfChild(l.parent(), l) }, new ProbeList[] { l });
    }
    Enumeration<TreeModelListener> e = listeners.elements();
    while (e.hasMoreElements()) {
        e.nextElement().treeNodesRemoved(me);
    }
}
Also used : ProbeSet(uk.ac.babraham.SeqMonk.DataTypes.Probes.ProbeSet) TreeModelEvent(javax.swing.event.TreeModelEvent) ProbeList(uk.ac.babraham.SeqMonk.DataTypes.Probes.ProbeList) TreeModelListener(javax.swing.event.TreeModelListener)

Example 20 with ProbeSet

use of uk.ac.babraham.SeqMonk.DataTypes.Probes.ProbeSet in project SeqMonk by s-andrews.

the class ContigProbeGenerator method run.

/* (non-Javadoc)
	 * @see java.lang.Runnable#run()
	 */
public void run() {
    Chromosome[] chromosomes = collection.genome().getAllChromosomes();
    Vector<Probe> newProbes = new Vector<Probe>();
    for (int c = 0; c < chromosomes.length; c++) {
        // Time for an update
        updateGenerationProgress("Processed " + c + " chromosomes", c, chromosomes.length);
        // We'll merge together the reads for all of the selected DataStores and
        // compute a single set of probes which covers all of them.
        ReadsWithCounts[] v = new ReadsWithCounts[selectedStores.length];
        for (int s = 0; s < selectedStores.length; s++) {
            v[s] = selectedStores[s].getReadsForChromosome(chromosomes[c]);
        }
        ReadsWithCounts rawReads = new ReadsWithCounts(v);
        v = null;
        // We now want to convert this list into a non-redundant set of
        // read positions with counts.  If we don't do this then we get
        // appalling performance where we have many reads mapped at the
        // same position
        // Our default is to do all strands at once
        int[] strandsToTry = new int[] { 100 };
        if (separateStrands.isSelected()) {
            strandsToTry = new int[] { Location.FORWARD, Location.REVERSE, Location.UNKNOWN };
        }
        for (int strand = 0; strand < strandsToTry.length; strand++) {
            ReadsWithCounts reads = getNonRedundantReads(rawReads, strandsToTry[strand]);
            if (reads.totalCount() == 0) {
                // System.err.println("Skipping strand "+strandsToTry[strand]+" on chr "+chromosomes[c]);
                continue;
            }
            int strandForNewProbes = Location.UNKNOWN;
            if (strandsToTry.length > 1) {
                strandForNewProbes = strandsToTry[strand];
            }
            int start = -1;
            // We now start a process where we work out at what point we cross the
            // threshold of having more than depthCutoff reads overlapping at any
            // point
            LinkedList<SequenceReadWithCount> currentSet = new LinkedList<SequenceReadWithCount>();
            int currentSetSize = 0;
            for (int r = 0; r < reads.reads.length; r++) {
                // See if we need to quit
                if (cancel) {
                    generationCancelled();
                    return;
                }
                while (currentSetSize > 0 && SequenceRead.end(currentSet.getFirst().read) < SequenceRead.start(reads.reads[r])) {
                    SequenceReadWithCount lastRead = currentSet.removeFirst();
                    currentSetSize -= lastRead.count;
                    if (start > 0 && currentSetSize < depthCutoff) {
                        // We just got to the end of a probe
                        Probe p = new Probe(chromosomes[c], start, SequenceRead.end(lastRead.read), strandForNewProbes);
                        // Check to see if we have a previous probe against which we can check
                        Probe lastProbe = null;
                        if (!newProbes.isEmpty())
                            lastProbe = newProbes.lastElement();
                        // Can we merge?
                        if (lastProbe != null && p.chromosome() == lastProbe.chromosome() && p.strand() == lastProbe.strand() && p.start() - lastProbe.end() <= distance) {
                            // Remove the last probe from the stored set
                            newProbes.remove(newProbes.size() - 1);
                            // Expand this probe to cover the last one and add it to the stored set
                            newProbes.add(new Probe(p.chromosome(), lastProbe.start(), p.end(), strandForNewProbes));
                        } else if (lastProbe != null) {
                            // We might still remove this if it's too small
                            if (lastProbe.length() < minSize) {
                                newProbes.remove(newProbes.size() - 1);
                            }
                            // We still need to add the new probe
                            newProbes.add(p);
                        } else {
                            newProbes.add(p);
                        }
                        start = -1;
                    }
                }
                // If there's nothing there already then just add it
                if (currentSetSize == 0) {
                    currentSet.add(new SequenceReadWithCount(reads.reads[r], reads.counts[r]));
                    currentSetSize += reads.counts[r];
                } else // If there are reads in the current set then we need to add this read
                // so that the current set is ordered by the end positions of the
                // reads, with the earliest end first.  We therefore start from the back
                // and work our way to the front, as soon as we see an entry whose end
                // is lower than ours we add ourselves after that
                {
                    // Now we add this read at a position based on its end position
                    ListIterator<SequenceReadWithCount> it = currentSet.listIterator(currentSet.size());
                    while (true) {
                        // If we reach the front of the set then we add ourselves to the front
                        if (!it.hasPrevious()) {
                            currentSet.addFirst(new SequenceReadWithCount(reads.reads[r], reads.counts[r]));
                            currentSetSize += reads.counts[r];
                            break;
                        } else {
                            SequenceReadWithCount previousRead = it.previous();
                            if (SequenceRead.end(previousRead.read) < SequenceRead.end(reads.reads[r])) {
                                // We want to add ourselves after this element so backtrack
                                // by one position (which must exist because we just went
                                // past it
                                it.next();
                                it.add(new SequenceReadWithCount(reads.reads[r], reads.counts[r]));
                                currentSetSize += reads.counts[r];
                                break;
                            }
                        }
                    }
                }
                // See if we crossed the threshold for starting a new probe
                if (start < 0 && currentSetSize >= depthCutoff) {
                    start = SequenceRead.start(reads.reads[r]);
                }
            }
            // out of the so far unprocessed reads on this chromosome
            if (start > 0) {
                Probe p = new Probe(chromosomes[c], start, SequenceRead.end(currentSet.getFirst().read), strandForNewProbes);
                // Check to see if we can merge with the last probe made
                Probe lastProbe = null;
                if (!newProbes.isEmpty())
                    lastProbe = newProbes.lastElement();
                // Can we merge?
                if (lastProbe != null && p.chromosome() == lastProbe.chromosome() && p.start() - lastProbe.end() <= distance) {
                    newProbes.remove(newProbes.size() - 1);
                    newProbes.add(new Probe(p.chromosome(), lastProbe.start(), p.end(), strandForNewProbes));
                } else if (lastProbe != null) {
                    // We might still remove this if it's too small
                    if (lastProbe.length() < minSize) {
                        newProbes.remove(newProbes.size() - 1);
                    }
                    // final chance
                    if (p.length() > minSize) {
                        newProbes.add(p);
                    }
                } else {
                    // Add the remaining probe if it's big enough.
                    if (p.length() > minSize) {
                        newProbes.add(p);
                    }
                }
            }
        }
    }
    Probe[] finalList = newProbes.toArray(new Probe[0]);
    newProbes.clear();
    ProbeSet finalSet = new ProbeSet(getDescription(), finalList);
    generationComplete(finalSet);
}
Also used : Chromosome(uk.ac.babraham.SeqMonk.DataTypes.Genome.Chromosome) ReadsWithCounts(uk.ac.babraham.SeqMonk.DataTypes.Sequence.ReadsWithCounts) Probe(uk.ac.babraham.SeqMonk.DataTypes.Probes.Probe) LinkedList(java.util.LinkedList) ProbeSet(uk.ac.babraham.SeqMonk.DataTypes.Probes.ProbeSet) Vector(java.util.Vector) IntVector(uk.ac.babraham.SeqMonk.Utilities.IntVector) LongVector(uk.ac.babraham.SeqMonk.Utilities.LongVector)

Aggregations

ProbeSet (uk.ac.babraham.SeqMonk.DataTypes.Probes.ProbeSet)30 Probe (uk.ac.babraham.SeqMonk.DataTypes.Probes.Probe)26 Vector (java.util.Vector)22 Chromosome (uk.ac.babraham.SeqMonk.DataTypes.Genome.Chromosome)22 Feature (uk.ac.babraham.SeqMonk.DataTypes.Genome.Feature)9 Location (uk.ac.babraham.SeqMonk.DataTypes.Genome.Location)5 ProbeList (uk.ac.babraham.SeqMonk.DataTypes.Probes.ProbeList)5 QuantitationStrandType (uk.ac.babraham.SeqMonk.DataTypes.Sequence.QuantitationStrandType)5 SplitLocation (uk.ac.babraham.SeqMonk.DataTypes.Genome.SplitLocation)4 SeqMonkException (uk.ac.babraham.SeqMonk.SeqMonkException)4 LongVector (uk.ac.babraham.SeqMonk.Utilities.LongVector)4 BinomialDistribution (org.apache.commons.math3.distribution.BinomialDistribution)3 ArrayList (java.util.ArrayList)2 Hashtable (java.util.Hashtable)2 Random (java.util.Random)2 ProbeTTestValue (uk.ac.babraham.SeqMonk.Analysis.Statistics.ProbeTTestValue)2 ReadsWithCounts (uk.ac.babraham.SeqMonk.DataTypes.Sequence.ReadsWithCounts)2 HashMap (java.util.HashMap)1 HashSet (java.util.HashSet)1 LinkedList (java.util.LinkedList)1