Search in sources :

Example 1 with PairedReadMergingResult

use of com.milaboratory.core.merger.PairedReadMergingResult in project mixcr by milaboratory.

the class TargetMerger method merge.

/**
 * @param targetLeft       left sequence
 * @param targetRight      right sequence
 * @param trySequenceMerge whether to try merging using sequence overlap (if alignment overlap failed)
 */
public TargetMergingResult merge(AlignedTarget targetLeft, AlignedTarget targetRight, boolean trySequenceMerge) {
    for (GeneType geneType : GeneType.VJC_REFERENCE) {
        if (!hasAlignments(targetLeft, geneType) || !hasAlignments(targetRight, geneType))
            continue;
        List<HitMappingRecord> als = extractSortedHits(targetLeft, targetRight, geneType);
        Alignment<NucleotideSequence>[] topHits = als.get(0).alignments;
        if (topHits[0] != null && topHits[1] != null) {
            final Alignment<NucleotideSequence> left = topHits[0];
            final Alignment<NucleotideSequence> right = topHits[1];
            final int from = Math.max(left.getSequence1Range().getFrom(), right.getSequence1Range().getFrom());
            final int to = Math.min(left.getSequence1Range().getTo(), right.getSequence1Range().getTo());
            if (to <= from)
                continue;
            int delta = left.convertToSeq2Position(from) - right.convertToSeq2Position(from);
            if (delta != left.convertToSeq2Position(to) - right.convertToSeq2Position(to))
                continue;
            int seq1Offset = delta > 0 ? delta : 0;
            int seq2Offset = delta > 0 ? 0 : -delta;
            int overlap = Math.min(targetLeft.getTarget().size() - seq1Offset, targetRight.getTarget().size() - seq2Offset);
            int mismatches = SequencesUtils.mismatchCount(targetLeft.getTarget().getSequence(), seq1Offset, targetRight.getTarget().getSequence(), seq2Offset, overlap);
            double identity = MismatchOnlyPairedReadMerger.identity(identityType, targetLeft.getTarget(), seq1Offset, targetRight.getTarget(), seq2Offset, overlap);
            if (identity < minimalAlignmentMergeIdentity)
                return new TargetMergingResult(geneType);
            final AlignedTarget merge = merge(targetLeft, targetRight, delta, OverlapType.AlignmentOverlap, mismatches);
            return new TargetMergingResult(true, null, merge, PairedReadMergingResult.MATCH_SCORE * (overlap - mismatches) + PairedReadMergingResult.MISMATCH_SCORE * mismatches, overlap, mismatches, delta);
        }
    }
    if (!trySequenceMerge)
        return new TargetMergingResult();
    final PairedReadMergingResult merge = merger.merge(targetLeft.getTarget(), targetRight.getTarget());
    if (!merge.isSuccessful())
        return new TargetMergingResult();
    return new TargetMergingResult(false, null, merge(targetLeft, targetRight, merge.getOffset(), OverlapType.SequenceOverlap, merge.getErrors()), merge.score(), merge.getOverlap(), merge.getErrors(), merge.getOffset());
}
Also used : Alignment(com.milaboratory.core.alignment.Alignment) NucleotideSequence(com.milaboratory.core.sequence.NucleotideSequence) PairedReadMergingResult(com.milaboratory.core.merger.PairedReadMergingResult) GeneType(io.repseq.core.GeneType)

Example 2 with PairedReadMergingResult

use of com.milaboratory.core.merger.PairedReadMergingResult in project mixcr by milaboratory.

the class VDJCAlignerWithMerge method process0.

@Override
protected VDJCAlignmentResult<PairedRead> process0(final PairedRead read) {
    PairedReadMergingResult merged = merger.process(read);
    if (merged.isSuccessful()) {
        VDJCAlignments alignment = singleAligner.process(new SingleReadImpl(read.getId(), merged.getOverlappedSequence(), "")).alignment;
        if (listener != null)
            listener.onSuccessfulSequenceOverlap(read, alignment);
        if (alignment != null) {
            boolean isRC = ((SequenceHistory.RawSequence) alignment.getHistory(0)).index.isReverseComplement;
            alignment = alignment.setHistory(new SequenceHistory[] { new SequenceHistory.Merge(SequenceHistory.OverlapType.SequenceOverlap, new SequenceHistory.RawSequence(read.getId(), (byte) (isRC ? 1 : 0), false, read.getR1().getData().size()), new SequenceHistory.RawSequence(read.getId(), (byte) (isRC ? 0 : 1), merged.isReversed(), read.getR2().getData().size()), merged.getOffset(), merged.getErrors()) }, new SequenceRead[] { read });
        }
        return new VDJCAlignmentResult<>(read, alignment);
    } else
        return pairedAligner.process(read);
}
Also used : SingleReadImpl(com.milaboratory.core.io.sequence.SingleReadImpl) PairedReadMergingResult(com.milaboratory.core.merger.PairedReadMergingResult) SequenceRead(com.milaboratory.core.io.sequence.SequenceRead) SequenceHistory(com.milaboratory.mixcr.basictypes.SequenceHistory) VDJCAlignments(com.milaboratory.mixcr.basictypes.VDJCAlignments)

Aggregations

PairedReadMergingResult (com.milaboratory.core.merger.PairedReadMergingResult)2 Alignment (com.milaboratory.core.alignment.Alignment)1 SequenceRead (com.milaboratory.core.io.sequence.SequenceRead)1 SingleReadImpl (com.milaboratory.core.io.sequence.SingleReadImpl)1 NucleotideSequence (com.milaboratory.core.sequence.NucleotideSequence)1 SequenceHistory (com.milaboratory.mixcr.basictypes.SequenceHistory)1 VDJCAlignments (com.milaboratory.mixcr.basictypes.VDJCAlignments)1 GeneType (io.repseq.core.GeneType)1