Search in sources :

Example 41 with Range

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

the class SequenceBaseTest method test4.

@Test
public void test4() throws Exception {
    SequenceBase base = new SequenceBase();
    base.put("A1", 18, new NucleotideSequence("CAC"));
    base.put("A1", 10, new NucleotideSequence("ATTAGACACACAC"));
    assertEquals(new NucleotideSequence("ATTAGACACAC"), base.get("A1", new Range(10, 21)));
}
Also used : NucleotideSequence(com.milaboratory.core.sequence.NucleotideSequence) Range(com.milaboratory.core.Range) Test(org.junit.Test)

Example 42 with Range

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

the class SequenceBaseTest method test2m.

@Test
public void test2m() throws Exception {
    SequenceBase base = new SequenceBase();
    base.put("A1", 100, new NucleotideSequence("ATTAGACACACAC"));
    base.put("A1", 10, new NucleotideSequence("ATTAGACACACAC"));
    base.put("A1", 20, new NucleotideSequence("CACATA"));
    assertEquals(new NucleotideSequence("ACACA"), base.get("A1", new Range(19, 24)));
    assertEquals(new NucleotideSequence("TTAG"), base.get("A1", new Range(101, 105)));
}
Also used : NucleotideSequence(com.milaboratory.core.sequence.NucleotideSequence) Range(com.milaboratory.core.Range) Test(org.junit.Test)

Example 43 with Range

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

the class VDJCAlignmentsFormatter method getTargetAsMultiAlignment.

public static MultiAlignmentHelper getTargetAsMultiAlignment(VDJCObject vdjcObject, int targetId, boolean addHitScore, boolean addReads) {
    if (addReads && !(vdjcObject instanceof VDJCAlignments))
        throw new IllegalArgumentException("Read alignments supported only for VDJCAlignments.");
    NSequenceWithQuality target = vdjcObject.getTarget(targetId);
    NucleotideSequence targetSeq = target.getSequence();
    SequencePartitioning partitioning = vdjcObject.getPartitionedTarget(targetId).getPartitioning();
    List<Alignment<NucleotideSequence>> alignments = new ArrayList<>();
    List<String> alignmentLeftComments = new ArrayList<>();
    List<String> alignmentRightComments = new ArrayList<>();
    for (GeneType gt : GeneType.values()) for (VDJCHit hit : vdjcObject.getHits(gt)) {
        Alignment<NucleotideSequence> alignment = hit.getAlignment(targetId);
        if (alignment == null)
            continue;
        alignment = alignment.invert(targetSeq);
        alignments.add(alignment);
        alignmentLeftComments.add(hit.getGene().getName());
        alignmentRightComments.add(" " + (int) (hit.getAlignment(targetId).getScore()) + (addHitScore ? " (" + (int) (hit.getScore()) + ")" : ""));
    }
    // Adding read information
    if (addReads) {
        VDJCAlignments vdjcAlignments = (VDJCAlignments) vdjcObject;
        SequenceHistory history = vdjcAlignments.getHistory(targetId);
        List<SequenceHistory.RawSequence> reads = history.rawReads();
        for (SequenceHistory.RawSequence read : reads) {
            NucleotideSequence seq = vdjcAlignments.getOriginalSequence(read.index).getSequence();
            int offset = history.offset(read.index);
            Alignment<NucleotideSequence> alignment = Aligner.alignOnlySubstitutions(targetSeq, seq, offset, seq.size(), 0, seq.size(), AffineGapAlignmentScoring.IGBLAST_NUCLEOTIDE_SCORING);
            alignments.add(alignment);
            alignmentLeftComments.add(read.index.toString());
            alignmentRightComments.add("");
        }
    }
    MultiAlignmentHelper helper = MultiAlignmentHelper.build(MultiAlignmentHelper.DEFAULT_SETTINGS, new Range(0, target.size()), targetSeq, alignments.toArray(new Alignment[alignments.size()]));
    if (!alignments.isEmpty())
        drawPoints(helper, partitioning, POINTS_FOR_REARRANGED);
    drawAASequence(helper, partitioning, targetSeq);
    helper.addSubjectQuality("Quality", target.getQuality());
    helper.setSubjectLeftTitle("Target" + targetId);
    helper.setSubjectRightTitle(" Score" + (addHitScore ? " (hit score)" : ""));
    for (int i = 0; i < alignmentLeftComments.size(); i++) {
        helper.setQueryLeftTitle(i, alignmentLeftComments.get(i));
        helper.setQueryRightTitle(i, alignmentRightComments.get(i));
    }
    return helper;
}
Also used : MultiAlignmentHelper(com.milaboratory.core.alignment.MultiAlignmentHelper) ArrayList(java.util.ArrayList) SequencePartitioning(io.repseq.core.SequencePartitioning) Range(com.milaboratory.core.Range) ReferencePoint(io.repseq.core.ReferencePoint) Alignment(com.milaboratory.core.alignment.Alignment) NSequenceWithQuality(com.milaboratory.core.sequence.NSequenceWithQuality) NucleotideSequence(com.milaboratory.core.sequence.NucleotideSequence) GeneType(io.repseq.core.GeneType)

Example 44 with Range

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

the class CloneFactory method create.

public Clone create(int id, double count, EnumMap<GeneType, TObjectFloatHashMap<VDJCGeneId>> geneScores, NSequenceWithQuality[] targets) {
    EnumMap<GeneType, VDJCHit[]> hits = new EnumMap<>(GeneType.class);
    for (GeneType geneType : GeneType.VJC_REFERENCE) {
        VJCClonalAlignerParameters vjcParameters = parameters.getVJCParameters(geneType);
        if (vjcParameters == null)
            continue;
        GeneFeature featureToAlign = featuresToAlign.get(geneType);
        TObjectFloatHashMap<VDJCGeneId> gtGeneScores = geneScores.get(geneType);
        if (gtGeneScores == null)
            continue;
        GeneFeature[] intersectingFeatures = new GeneFeature[assemblingFeatures.length];
        for (int i = 0; i < assemblingFeatures.length; ++i) {
            intersectingFeatures[i] = GeneFeature.intersection(featureToAlign, assemblingFeatures[i]);
            if (intersectingFeatures[i] != null)
                switch(geneType) {
                    case Variable:
                        if (!intersectingFeatures[i].getFirstPoint().equals(assemblingFeatures[i].getFirstPoint()))
                            throw new IllegalArgumentException();
                        break;
                    case Joining:
                        if (!intersectingFeatures[i].getLastPoint().equals(assemblingFeatures[i].getLastPoint()))
                            throw new IllegalArgumentException();
                        break;
                }
        }
        VDJCHit[] result = new VDJCHit[gtGeneScores.size()];
        int pointer = 0;
        TObjectFloatIterator<VDJCGeneId> iterator = gtGeneScores.iterator();
        while (iterator.hasNext()) {
            iterator.advance();
            VDJCGene gene = usedGenes.get(iterator.key());
            Alignment<NucleotideSequence>[] alignments = new Alignment[assemblingFeatures.length];
            for (int i = 0; i < assemblingFeatures.length; ++i) {
                if (intersectingFeatures[i] == null)
                    continue;
                NucleotideSequence referenceSequence = gene.getFeature(featureToAlign);
                Range rangeInReference = gene.getPartitioning().getRelativeRange(featureToAlign, intersectingFeatures[i]);
                if (rangeInReference == null || referenceSequence == null)
                    continue;
                Boolean leftSide;
                switch(geneType) {
                    case Variable:
                        leftSide = intersectingFeatures[i].getLastPoint().isTrimmable() ? true : null;
                        break;
                    case Joining:
                        leftSide = intersectingFeatures[i].getFirstPoint().isTrimmable() ? false : null;
                        break;
                    case Constant:
                        leftSide = null;
                        break;
                    default:
                        throw new RuntimeException();
                }
                BandedAlignerParameters<NucleotideSequence> alignmentParameters = vjcParameters.getAlignmentParameters();
                int referenceLength = rangeInReference.length();
                NucleotideSequence target = targets[i].getSequence();
                if (alignmentParameters.getScoring() instanceof LinearGapAlignmentScoring) {
                    if (leftSide == null) {
                        alignments[i] = BandedLinearAligner.align((LinearGapAlignmentScoring<NucleotideSequence>) alignmentParameters.getScoring(), referenceSequence, target, rangeInReference.getFrom(), referenceLength, 0, target.size(), alignmentParameters.getWidth());
                    } else if (leftSide) {
                        assert rangeInReference.getFrom() + referenceLength == referenceSequence.size();
                        alignments[i] = BandedLinearAligner.alignSemiLocalLeft((LinearGapAlignmentScoring<NucleotideSequence>) alignmentParameters.getScoring(), referenceSequence, target, rangeInReference.getFrom(), referenceLength, 0, target.size(), alignmentParameters.getWidth(), alignmentParameters.getStopPenalty());
                    } else {
                        assert rangeInReference.getFrom() == 0;
                        // int offset2 = Math.max(0, target.size() - referenceLength);
                        alignments[i] = BandedLinearAligner.alignSemiLocalRight((LinearGapAlignmentScoring<NucleotideSequence>) alignmentParameters.getScoring(), referenceSequence, target, rangeInReference.getFrom(), referenceLength, 0, target.size(), alignmentParameters.getWidth(), alignmentParameters.getStopPenalty());
                    }
                } else {
                    if (leftSide == null) {
                        alignments[i] = BandedAffineAligner.align((AffineGapAlignmentScoring<NucleotideSequence>) alignmentParameters.getScoring(), referenceSequence, target, rangeInReference.getFrom(), referenceLength, 0, target.size(), alignmentParameters.getWidth());
                    } else if (leftSide) {
                        assert rangeInReference.getFrom() + referenceLength == referenceSequence.size();
                        alignments[i] = BandedAffineAligner.semiLocalRight((AffineGapAlignmentScoring<NucleotideSequence>) alignmentParameters.getScoring(), referenceSequence, target, rangeInReference.getFrom(), referenceLength, 0, target.size(), alignmentParameters.getWidth());
                    } else {
                        assert rangeInReference.getFrom() == 0;
                        // int offset2 = Math.max(0, target.size() - referenceLength);
                        alignments[i] = BandedAffineAligner.semiLocalLeft((AffineGapAlignmentScoring<NucleotideSequence>) alignmentParameters.getScoring(), referenceSequence, target, rangeInReference.getFrom(), referenceLength, 0, target.size(), alignmentParameters.getWidth());
                    }
                }
            }
            result[pointer++] = new VDJCHit(gene, alignments, featureToAlign, iterator.value());
        }
        Arrays.sort(result, 0, pointer);
        hits.put(geneType, pointer < result.length ? Arrays.copyOf(result, pointer) : result);
    }
    // D
    NucleotideSequence sequenceToAlign = targets[indexOfAssemblingFeatureWithD].getSequence();
    int from = 0;
    int to = sequenceToAlign.size();
    VDJCHit[] hs = hits.get(GeneType.Variable);
    if (hs != null && hs.length > 0) {
        int p = hs[0].getPartitioningForTarget(indexOfAssemblingFeatureWithD).getPosition(ReferencePoint.VEndTrimmed);
        if (p != -1) {
            if (p < 0)
                p = -2 - p;
            from = p;
        }
    }
    hs = hits.get(GeneType.Joining);
    if (hs != null && hs.length > 0) {
        int p = hs[0].getPartitioningForTarget(indexOfAssemblingFeatureWithD).getPosition(ReferencePoint.JBeginTrimmed);
        if (p != -1) {
            if (p < 0)
                p = -2 - p;
            to = p;
        }
    }
    if (from < to)
        hits.put(GeneType.Diversity, dAligner.align(sequenceToAlign, VDJCAligner.getPossibleDLoci(hits.get(GeneType.Variable), hits.get(GeneType.Joining)), from, to, indexOfAssemblingFeatureWithD, assemblingFeatures.length));
    else
        hits.put(GeneType.Diversity, new VDJCHit[0]);
    return new Clone(targets, hits, count, id);
}
Also used : Range(com.milaboratory.core.Range) NucleotideSequence(com.milaboratory.core.sequence.NucleotideSequence) VDJCHit(com.milaboratory.mixcr.basictypes.VDJCHit) Clone(com.milaboratory.mixcr.basictypes.Clone)

Example 45 with Range

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

the class RangeSetTest method testIntersection1.

@Test
public void testIntersection1() throws Exception {
    RangeSet s10 = createUnsafe(10, 20, 30, 40, 50, 60, 70, 80);
    RangeSet s12 = s10.intersection(new Range(19, 55));
    Assert.assertEquals(s12, createUnsafe(19, 20, 30, 40, 50, 55));
}
Also used : Range(com.milaboratory.core.Range) 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