Search in sources :

Example 11 with Range

use of com.milaboratory.core.Range in project mixcr by milaboratory.

the class PartialAlignmentsAssembler method searchOverlaps.

public void searchOverlaps(VDJCAlignmentsReader reader) {
    final VDJCAlignerParameters alignerParameters = reader.getParameters();
    PartialAlignmentsAssemblerAligner aligner = new PartialAlignmentsAssemblerAligner(alignerParameters);
    targetMerger.setAlignerParameters(alignerParameters);
    for (VDJCGene gene : reader.getUsedGenes()) aligner.addGene(gene);
    for (final VDJCAlignments alignment : CUtils.it(reader)) {
        total.incrementAndGet();
        if (alignment.getFeature(GeneFeature.CDR3) != null) {
            containsCDR3.incrementAndGet();
            if (!overlappedOnly) {
                totalWritten.incrementAndGet();
                writer.write(alignment);
            }
            continue;
        }
        if (alreadyMergedIds.contains(alignment.getAlignmentsIndex()))
            continue;
        final OverlapSearchResult searchResult = searchOverlaps(alignment, alignerParameters.isAllowChimeras());
        // Common procedure to cancel processing of current input alignment if it fails to pass some good
        // overlap filtering criterion
        final Runnable cancelCurrentResult = new Runnable() {

            @Override
            public void run() {
                if (searchResult != null)
                    searchResult.cancel();
                if (writePartial && !overlappedOnly && notInLeftIndexIds.contains(alignment.getAlignmentsIndex())) {
                    totalWritten.incrementAndGet();
                    partialAsIs.incrementAndGet();
                    writer.write(alignment);
                }
            }
        };
        if (searchResult == null) {
            cancelCurrentResult.run();
            continue;
        }
        List<AlignedTarget> mergedTargets = searchResult.result;
        VDJCMultiRead mRead = new VDJCMultiRead(mergedTargets);
        final VDJCAlignments mAlignment = aligner.process(mRead).alignment;
        // Checking number of overlapped non-template (NRegion) letters
        int overlapTargetId = -1;
        Range overlapRange = null;
        for (int i = 0; i < mergedTargets.size(); i++) {
            overlapRange = AlignedTarget.getOverlapRange(mergedTargets.get(i));
            if (overlapRange != null) {
                overlapTargetId = i;
                break;
            }
        }
        if (overlapTargetId == -1) {
            // No alignments for Best V Hit and Best J Hit in central (overlapped) target
            cancelCurrentResult.run();
            continue;
        }
        int targetLength = mergedTargets.get(overlapTargetId).getTarget().size();
        VDJCHit bestVHit = mAlignment.getBestHit(GeneType.Variable), bestJHit = mAlignment.getBestHit(GeneType.Joining);
        if (bestVHit == null || bestJHit == null || bestVHit.getAlignment(overlapTargetId) == null || bestJHit.getAlignment(overlapTargetId) == null) {
            cancelCurrentResult.run();
            continue;
        }
        int ndnRegionBegin = 0;
        int ndnRegionEnd = targetLength;
        Alignment<NucleotideSequence> vAlignment = bestVHit.getAlignment(overlapTargetId);
        if (vAlignment != null)
            ndnRegionBegin = vAlignment.getSequence2Range().getTo();
        Alignment<NucleotideSequence> jAlignment = bestJHit.getAlignment(overlapTargetId);
        if (jAlignment != null)
            ndnRegionEnd = jAlignment.getSequence2Range().getFrom();
        RangeSet nRegion = ndnRegionBegin >= ndnRegionEnd ? RangeSet.EMPTY : RangeSet.create(ndnRegionBegin, ndnRegionEnd);
        Range dRange = mAlignment.getPartitionedTarget(overlapTargetId).getPartitioning().getRange(GeneFeature.DRegionTrimmed);
        if (dRange != null)
            nRegion = nRegion.subtract(dRange);
        RangeSet nRegionInOverlap = nRegion.intersection(overlapRange);
        int actualNRegionLength = nRegion.totalLength();
        int minimalN = Math.min(minimalNOverlap, actualNRegionLength);
        if (nRegionInOverlap.totalLength() < minimalN) {
            droppedSmallOverlapNRegion.incrementAndGet();
            cancelCurrentResult.run();
            continue;
        }
        // DDDDDDDDDJJJJJJJJJJJJJJJ
        if (minimalN == 0 && (!overlapRange.contains(ndnRegionBegin - 1) || !overlapRange.contains(ndnRegionEnd))) {
            droppedNoNRegion.incrementAndGet();
            cancelCurrentResult.run();
            continue;
        }
        overlapped.incrementAndGet();
        totalWritten.incrementAndGet();
        writer.write(mAlignment);
        // Saving alignment that where merge to prevent it's use as left part
        alreadyMergedIds.add(alignment.getAlignmentsIndex());
        alreadyMergedIds.add(searchResult.KMerInfo.alignments.getAlignmentsIndex());
    }
    if (writePartial && !overlappedOnly)
        for (List<KMerInfo> kMerInfos : kToIndexLeft.valueCollection()) for (KMerInfo kMerInfo : kMerInfos) if (alreadyMergedIds.add(kMerInfo.alignments.getAlignmentsIndex())) {
            totalWritten.incrementAndGet();
            partialAsIs.incrementAndGet();
            writer.write(kMerInfo.getAlignments());
        }
    writer.setNumberOfProcessedReads(reader.getNumberOfReads() - overlapped.get());
}
Also used : VDJCAlignerParameters(com.milaboratory.mixcr.vdjaligners.VDJCAlignerParameters) Range(com.milaboratory.core.Range) NucleotideSequence(com.milaboratory.core.sequence.NucleotideSequence) ArrayList(java.util.ArrayList) List(java.util.List)

Example 12 with Range

use of com.milaboratory.core.Range in project repseqio by repseqio.

the class RangeTranslationParametersTest method testIncompleteCodon.

@Test
public void testIncompleteCodon() throws Exception {
    assertEquals(new Range(9, 10), new RangeTranslationParameters(ReferencePoint.CDR3Begin, ReferencePoint.VEndTrimmed, new Range(0, 10)).rightIncompleteCodonRange());
    assertEquals(new Range(9, 11), new RangeTranslationParameters(ReferencePoint.CDR3Begin, ReferencePoint.VEndTrimmed, new Range(0, 11)).rightIncompleteCodonRange());
    assertEquals(null, new RangeTranslationParameters(ReferencePoint.CDR3Begin, ReferencePoint.VEndTrimmed, new Range(0, 12)).rightIncompleteCodonRange());
    assertEquals(null, new RangeTranslationParameters(ReferencePoint.CDR3Begin, ReferencePoint.VEndTrimmed, new Range(0, 10)).leftIncompleteCodonRange());
    assertEquals(new Range(0, 1), new RangeTranslationParameters(ReferencePoint.JBeginTrimmed, ReferencePoint.CDR3End, new Range(0, 10)).leftIncompleteCodonRange());
    assertEquals(new Range(0, 2), new RangeTranslationParameters(ReferencePoint.JBeginTrimmed, ReferencePoint.CDR3End, new Range(0, 11)).leftIncompleteCodonRange());
    assertEquals(null, new RangeTranslationParameters(ReferencePoint.JBeginTrimmed, ReferencePoint.CDR3End, new Range(0, 12)).leftIncompleteCodonRange());
    assertEquals(null, new RangeTranslationParameters(ReferencePoint.JBeginTrimmed, ReferencePoint.CDR3End, new Range(0, 10)).rightIncompleteCodonRange());
}
Also used : Range(com.milaboratory.core.Range) Test(org.junit.Test)

Example 13 with Range

use of com.milaboratory.core.Range in project repseqio by repseqio.

the class ReferencePointsTest method test9.

@Test
public void test9() throws Exception {
    ReferencePoints points = new ReferencePoints(0, new int[] { 2, 52, 63, 84, 155, 455, 645, 1255, 2142 });
    GeneFeature feature = GeneFeatureTest.create(2, 5);
    GeneFeature minor = GeneFeatureTest.create(3, 4);
    Assert.assertEquals(new Range(84 - 63, 155 - 63), points.getRelativeRange(feature, minor));
    feature = GeneFeatureTest.create(2, 5, 6, 8);
    minor = GeneFeatureTest.create(6, 8);
    Assert.assertEquals(new Range(455 - 63, 2142 - 645 + 455 - 63), points.getRelativeRange(feature, minor));
}
Also used : Range(com.milaboratory.core.Range) Test(org.junit.Test)

Example 14 with Range

use of com.milaboratory.core.Range in project repseqio by repseqio.

the class SequenceResolverTest method rawHttpTest2.

@Test
public void rawHttpTest2() throws Exception {
    Path dir = TempFileManager.getTempDir().toPath().toAbsolutePath();
    Path work = dir.resolve("work");
    Files.createDirectories(work);
    Path cache = dir.resolve("cache");
    Files.createDirectories(cache);
    SequenceResolvers.initDefaultResolver(cache);
    SequenceResolver defaultResolver = SequenceResolvers.getDefault();
    NucleotideSequence seq = defaultResolver.resolve(new SequenceAddress("http://files.milaboratory.com/test-data/test.fa.gz#testRecord")).getRegion(new Range(10, 30).inverse());
    Assert.assertEquals(new NucleotideSequence("GCTCCACCACAAGACACTCT"), seq);
}
Also used : Path(java.nio.file.Path) NucleotideSequence(com.milaboratory.core.sequence.NucleotideSequence) Range(com.milaboratory.core.Range) SingleFastqReaderTest(com.milaboratory.core.io.sequence.fastq.SingleFastqReaderTest) Test(org.junit.Test)

Example 15 with Range

use of com.milaboratory.core.Range in project repseqio by repseqio.

the class SequenceResolverTest method test1.

@Test
public void test1() throws Exception {
    Path dir = TempFileManager.getTempDir().toPath().toAbsolutePath();
    Path work = dir.resolve("work");
    Files.createDirectories(work);
    Path cache = dir.resolve("cache");
    Files.createDirectories(cache);
    SequenceResolvers.initDefaultResolver(cache);
    SequenceResolver defaultResolver = SequenceResolvers.getDefault();
    NucleotideSequence seq = defaultResolver.resolve(new SequenceAddress("nuccore://EU877942.1")).getRegion(new Range(10, 30).inverse());
    Assert.assertEquals(new NucleotideSequence("gcgagagcaagcactatggc").getReverseComplement(), seq);
    Path someFasta = work.resolve("some_fasta.fasta");
    FileUtils.copyToFile(SingleFastqReaderTest.class.getClassLoader().getResourceAsStream("sequences/some_fasta.fasta"), someFasta.toFile());
    // Path rootPath = new File(SingleFastqReaderTest.class.getClassLoader().getResource("sequences/some_fasta.fasta")
    // .toURI()).toPath().toAbsolutePath().normalize().getParent();
    seq = defaultResolver.resolve(new SequenceAddress(work, "file://some_fasta.fasta#24.6jsd21.Tut")).getRegion(new Range(10, 30));
    Assert.assertEquals(new NucleotideSequence("ATCCTGGCTTAGAACTAACG"), seq);
}
Also used : Path(java.nio.file.Path) NucleotideSequence(com.milaboratory.core.sequence.NucleotideSequence) Range(com.milaboratory.core.Range) SingleFastqReaderTest(com.milaboratory.core.io.sequence.fastq.SingleFastqReaderTest) Test(org.junit.Test)

Aggregations

Range (com.milaboratory.core.Range)45 Test (org.junit.Test)23 NucleotideSequence (com.milaboratory.core.sequence.NucleotideSequence)19 VDJCHit (com.milaboratory.mixcr.basictypes.VDJCHit)4 ArrayList (java.util.ArrayList)4 SingleFastqReaderTest (com.milaboratory.core.io.sequence.fastq.SingleFastqReaderTest)3 MutationsBuilder (com.milaboratory.core.mutations.MutationsBuilder)3 ReferencePoint (io.repseq.core.ReferencePoint)3 Path (java.nio.file.Path)3 NSequenceWithQuality (com.milaboratory.core.sequence.NSequenceWithQuality)2 Clone (com.milaboratory.mixcr.basictypes.Clone)2 VDJCPartitionedSequence (com.milaboratory.mixcr.basictypes.VDJCPartitionedSequence)2 VDJCAlignerParameters (com.milaboratory.mixcr.vdjaligners.VDJCAlignerParameters)2 CUtils (cc.redberry.pipe.CUtils)1 OutputPort (cc.redberry.pipe.OutputPort)1 com.milaboratory.core.alignment (com.milaboratory.core.alignment)1 Alignment (com.milaboratory.core.alignment.Alignment)1 LinearGapAlignmentScoring (com.milaboratory.core.alignment.LinearGapAlignmentScoring)1 MultiAlignmentHelper (com.milaboratory.core.alignment.MultiAlignmentHelper)1 com.milaboratory.core.sequence (com.milaboratory.core.sequence)1