Examples with Clone com.milaboratory.mixcr.basictypes.Clone used on opensource projects

Example 1 with Clone

use of com.milaboratory.mixcr.basictypes.Clone in project mixcr by milaboratory.

the class ActionExportClones method go0.

@Override
public void go0() throws Exception {
    CloneExportParameters parameters = (CloneExportParameters) this.parameters;
    try (InputStream inputStream = IOUtil.createIS(parameters.getInputFile());
        InfoWriter<Clone> writer = new InfoWriter<>(parameters.getOutputFile())) {
        CloneSet set = CloneSetIO.readClns(inputStream, VDJCLibraryRegistry.getDefault());
        set = CloneSet.transform(set, parameters.getFilter());
        writer.attachInfoProviders((List) parameters.exporters);
        writer.ensureHeader();
        long limit = parameters.getLimit();
        for (int i = 0; i < set.size(); i++) {
            if (set.get(i).getFraction() < parameters.minFraction || set.get(i).getCount() < parameters.minCount) {
                limit = i;
                break;
            }
        }
        ExportClones exportClones = new ExportClones(set, writer, limit);
        SmartProgressReporter.startProgressReport(exportClones, System.err);
        exportClones.run();
    }
}

Example 2 with Clone

use of com.milaboratory.mixcr.basictypes.Clone in project mixcr by milaboratory.

the class FullSeqAssembler method buildClone.

/* ================================= Re-align and build final clone ====================================== */
private Clone buildClone(double count, BranchSequences targets) {
    Alignment<NucleotideSequence>[] vTopHitAlignments = new Alignment[targets.ranges.length], jTopHitAlignments = new Alignment[targets.ranges.length];
    VDJCHit hit = clone.getBestHit(Variable);
    if (hit == null)
        throw new UnsupportedOperationException("No V hit.");
    NucleotideSequence vTopReferenceSequence = hit.getGene().getFeature(hit.getAlignedFeature());
    hit = clone.getBestHit(Joining);
    if (hit == null)
        throw new UnsupportedOperationException("No J hit.");
    NucleotideSequence jTopReferenceSequence = hit.getGene().getFeature(hit.getAlignedFeature());
    // Excessive optimization
    CachedIntArray cachedIntArray = new CachedIntArray();
    int assemblingFeatureLength = targets.assemblingFeatureLength;
    for (int i = 0; i < targets.ranges.length; i++) {
        Range range = targets.ranges[i];
        NucleotideSequence sequence = targets.sequences[i].getSequence();
        // Asserts
        if (range.getFrom() < N_LEFT_DUMMIES + lengthV && range.getTo() >= N_LEFT_DUMMIES + lengthV && i != targets.assemblingFeatureTargetId)
            throw new RuntimeException();
        if (range.getTo() >= rightAssemblingFeatureBound && range.getFrom() < rightAssemblingFeatureBound && i != targets.assemblingFeatureTargetId)
            throw new RuntimeException();
        if (range.getTo() < N_LEFT_DUMMIES + lengthV) {
            boolean floatingLeftBound = !parameters.alignedRegionsOnly && i == 0 && alignerParameters.getVAlignerParameters().getParameters().isFloatingLeftBound();
            // Can be reduced to a single statement
            if (range.getFrom() < N_LEFT_DUMMIES)
                // This target contain extra non-V nucleotides on the left
                vTopHitAlignments[i] = alignLinearSeq1FromRight(((LinearGapAlignmentScoring<NucleotideSequence>) alignerParameters.getVAlignerParameters().getScoring()), vTopReferenceSequence, sequence.getSequence(), 0, range.getTo() - N_LEFT_DUMMIES, 0, sequence.size(), !floatingLeftBound, cachedIntArray);
            else if (floatingLeftBound)
                vTopHitAlignments[i] = alignLinearSeq1FromRight(((LinearGapAlignmentScoring<NucleotideSequence>) alignerParameters.getVAlignerParameters().getScoring()), vTopReferenceSequence, sequence.getSequence(), range.getFrom() - N_LEFT_DUMMIES, range.length(), 0, sequence.size(), false, cachedIntArray);
            else
                vTopHitAlignments[i] = Aligner.alignGlobal(alignerParameters.getVAlignerParameters().getScoring(), vTopReferenceSequence, sequence, range.getFrom() - N_LEFT_DUMMIES, range.length(), 0, sequence.size());
        } else if (i == targets.assemblingFeatureTargetId) {
            /*
                 *  V gene
                 */
            boolean vFloatingLeftBound = !parameters.alignedRegionsOnly && i == 0 && alignerParameters.getVAlignerParameters().getParameters().isFloatingLeftBound();
            // Can be reduced to a single statement
            if (range.getFrom() < N_LEFT_DUMMIES)
                // This target contain extra non-V nucleotides on the left
                vTopHitAlignments[i] = alignLinearSeq1FromRight(((LinearGapAlignmentScoring<NucleotideSequence>) alignerParameters.getVAlignerParameters().getScoring()), vTopReferenceSequence, sequence.getSequence(), 0, lengthV, 0, targets.assemblingFeatureOffset, !vFloatingLeftBound, cachedIntArray);
            else if (vFloatingLeftBound)
                vTopHitAlignments[i] = alignLinearSeq1FromRight(((LinearGapAlignmentScoring<NucleotideSequence>) alignerParameters.getVAlignerParameters().getScoring()), vTopReferenceSequence, sequence.getSequence(), range.getFrom() - N_LEFT_DUMMIES, lengthV - (range.getFrom() - N_LEFT_DUMMIES), 0, targets.assemblingFeatureOffset, false, cachedIntArray);
            else
                vTopHitAlignments[i] = Aligner.alignGlobal(alignerParameters.getVAlignerParameters().getScoring(), vTopReferenceSequence, sequence, range.getFrom() - N_LEFT_DUMMIES, lengthV - (range.getFrom() - N_LEFT_DUMMIES), 0, targets.assemblingFeatureOffset);
            /*
                 *  J gene
                 */
            boolean jFloatingRightBound = !parameters.alignedRegionsOnly && i == targets.ranges.length - 1 && alignerParameters.getJAlignerParameters().getParameters().isFloatingRightBound();
            if (range.getTo() >= rightAssemblingFeatureBound + jLength)
                // This target contain extra non-J nucleotides on the right
                jTopHitAlignments[i] = alignLinearSeq1FromLeft(((LinearGapAlignmentScoring<NucleotideSequence>) alignerParameters.getJAlignerParameters().getScoring()), jTopReferenceSequence, sequence.getSequence(), jOffset, jLength, targets.assemblingFeatureOffset + assemblingFeatureLength, sequence.size() - (targets.assemblingFeatureOffset + assemblingFeatureLength), !jFloatingRightBound, cachedIntArray);
            else if (jFloatingRightBound)
                jTopHitAlignments[i] = alignLinearSeq1FromLeft(((LinearGapAlignmentScoring<NucleotideSequence>) alignerParameters.getJAlignerParameters().getScoring()), jTopReferenceSequence, sequence.getSequence(), jOffset, range.getTo() - rightAssemblingFeatureBound, targets.assemblingFeatureOffset + assemblingFeatureLength, sequence.size() - (targets.assemblingFeatureOffset + assemblingFeatureLength), false, cachedIntArray);
            else
                jTopHitAlignments[i] = Aligner.alignGlobal(alignerParameters.getJAlignerParameters().getScoring(), jTopReferenceSequence, sequence, jOffset, range.getTo() - rightAssemblingFeatureBound, targets.assemblingFeatureOffset + assemblingFeatureLength, sequence.size() - (targets.assemblingFeatureOffset + assemblingFeatureLength));
        } else if (range.getFrom() > rightAssemblingFeatureBound) {
            boolean floatingRightBound = !parameters.alignedRegionsOnly && i == targets.ranges.length - 1 && alignerParameters.getJAlignerParameters().getParameters().isFloatingRightBound();
            if (range.getTo() >= rightAssemblingFeatureBound + jLength)
                // This target contain extra non-J nucleotides on the right
                jTopHitAlignments[i] = alignLinearSeq1FromLeft(((LinearGapAlignmentScoring<NucleotideSequence>) alignerParameters.getJAlignerParameters().getScoring()), jTopReferenceSequence, sequence.getSequence(), jOffset + (range.getFrom() - rightAssemblingFeatureBound), jLength - (range.getFrom() - rightAssemblingFeatureBound), 0, sequence.size(), !floatingRightBound, cachedIntArray);
            else if (floatingRightBound)
                jTopHitAlignments[i] = alignLinearSeq1FromLeft(((LinearGapAlignmentScoring<NucleotideSequence>) alignerParameters.getJAlignerParameters().getScoring()), jTopReferenceSequence, sequence.getSequence(), jOffset + (range.getFrom() - rightAssemblingFeatureBound), range.length(), 0, sequence.size(), false, cachedIntArray);
            else
                jTopHitAlignments[i] = Aligner.alignGlobal(alignerParameters.getJAlignerParameters().getScoring(), jTopReferenceSequence, sequence, jOffset + (range.getFrom() - rightAssemblingFeatureBound), range.length(), 0, sequence.size());
        } else
            throw new RuntimeException();
    }
    NSequenceWithQuality assemblingFeatureSeq = targets.sequences[targets.assemblingFeatureTargetId].getRange(targets.assemblingFeatureOffset, targets.assemblingFeatureOffset + targets.assemblingFeatureLength);
    Clone clone = cloneFactory.create(0, count, geneScores, new NSequenceWithQuality[] { assemblingFeatureSeq });
    vTopHitAlignments[targets.assemblingFeatureTargetId] = mergeTwoAlignments(vTopHitAlignments[targets.assemblingFeatureTargetId], clone.getBestHit(Variable).getAlignment(0).move(targets.assemblingFeatureOffset));
    jTopHitAlignments[targets.assemblingFeatureTargetId] = mergeTwoAlignments(clone.getBestHit(Joining).getAlignment(0).move(targets.assemblingFeatureOffset), jTopHitAlignments[targets.assemblingFeatureTargetId]);
    EnumMap<GeneType, VDJCHit[]> hits = new EnumMap<>(GeneType.class);
    for (GeneType gt : GeneType.VDJC_REFERENCE) hits.put(gt, Arrays.stream(clone.getHits(gt)).map(h -> moveHitTarget(h, targets.assemblingFeatureTargetId, targets.assemblingFeatureOffset, targets.ranges.length)).toArray(VDJCHit[]::new));
    VDJCHit[] tmp = hits.get(Variable);
    tmp[0] = substituteAlignments(tmp[0], vTopHitAlignments);
    tmp = hits.get(Joining);
    tmp[0] = substituteAlignments(tmp[0], jTopHitAlignments);
    return new Clone(targets.sequences, hits, count, 0);
}

Example 3 with Clone

use of com.milaboratory.mixcr.basictypes.Clone in project mixcr by milaboratory.

the class FieldExtractors method getFields.

public static synchronized Field[] getFields() {
    if (descriptors == null) {
        List<Field> descriptorsList = new ArrayList<>();
        // Number of targets
        descriptorsList.add(new PL_O("-targets", "Export number of targets", "Number of targets", "numberOfTargets") {

            @Override
            protected String extract(VDJCObject object) {
                return Integer.toString(object.numberOfTargets());
            }
        });
        // Best hits
        for (final GeneType type : GeneType.values()) {
            char l = type.getLetter();
            descriptorsList.add(new PL_O("-" + Character.toLowerCase(l) + "Hit", "Export best " + l + " hit", "Best " + l + " hit", "best" + l + "Hit") {

                @Override
                protected String extract(VDJCObject object) {
                    VDJCHit bestHit = object.getBestHit(type);
                    if (bestHit == null)
                        return NULL;
                    return bestHit.getGene().getName();
                }
            });
        }
        // Best gene
        for (final GeneType type : GeneType.values()) {
            char l = type.getLetter();
            descriptorsList.add(new PL_O("-" + Character.toLowerCase(l) + "Gene", "Export best " + l + " hit gene name (e.g. TRBV12-3 for TRBV12-3*00)", "Best " + l + " gene", "best" + l + "Gene") {

                @Override
                protected String extract(VDJCObject object) {
                    VDJCHit bestHit = object.getBestHit(type);
                    if (bestHit == null)
                        return NULL;
                    return bestHit.getGene().getGeneName();
                }
            });
        }
        // Best family
        for (final GeneType type : GeneType.values()) {
            char l = type.getLetter();
            descriptorsList.add(new PL_O("-" + Character.toLowerCase(l) + "Family", "Export best " + l + " hit family name (e.g. TRBV12 for TRBV12-3*00)", "Best " + l + " family", "best" + l + "Family") {

                @Override
                protected String extract(VDJCObject object) {
                    VDJCHit bestHit = object.getBestHit(type);
                    if (bestHit == null)
                        return NULL;
                    return bestHit.getGene().getFamilyName();
                }
            });
        }
        // Best hit score
        for (final GeneType type : GeneType.values()) {
            char l = type.getLetter();
            descriptorsList.add(new PL_O("-" + Character.toLowerCase(l) + "HitScore", "Export score for best " + l + " hit", "Best " + l + " hit score", "best" + l + "HitScore") {

                @Override
                protected String extract(VDJCObject object) {
                    VDJCHit bestHit = object.getBestHit(type);
                    if (bestHit == null)
                        return NULL;
                    return String.valueOf(bestHit.getScore());
                }
            });
        }
        // All hits
        for (final GeneType type : GeneType.values()) {
            char l = type.getLetter();
            descriptorsList.add(new PL_O("-" + Character.toLowerCase(l) + "HitsWithScore", "Export all " + l + " hits with score", "All " + l + " hits", "all" + l + "HitsWithScore") {

                @Override
                protected String extract(VDJCObject object) {
                    VDJCHit[] hits = object.getHits(type);
                    if (hits.length == 0)
                        return "";
                    StringBuilder sb = new StringBuilder();
                    for (int i = 0; ; i++) {
                        sb.append(hits[i].getGene().getName()).append("(").append(SCORE_FORMAT.format(hits[i].getScore())).append(")");
                        if (i == hits.length - 1)
                            break;
                        sb.append(",");
                    }
                    return sb.toString();
                }
            });
        }
        // All hits without score
        for (final GeneType type : GeneType.values()) {
            char l = type.getLetter();
            descriptorsList.add(new PL_O("-" + Character.toLowerCase(l) + "Hits", "Export all " + l + " hits", "All " + l + " Hits", "all" + l + "Hits") {

                @Override
                protected String extract(VDJCObject object) {
                    VDJCHit[] hits = object.getHits(type);
                    if (hits.length == 0)
                        return "";
                    StringBuilder sb = new StringBuilder();
                    for (int i = 0; ; i++) {
                        sb.append(hits[i].getGene().getName());
                        if (i == hits.length - 1)
                            break;
                        sb.append(",");
                    }
                    return sb.toString();
                }
            });
        }
        // All gene names
        for (final GeneType type : GeneType.values()) {
            char l = type.getLetter();
            descriptorsList.add(new StringExtractor("-" + Character.toLowerCase(l) + "Genes", "Export all " + l + " gene names (e.g. TRBV12-3 for TRBV12-3*00)", "All " + l + " genes", "all" + l + "Genes", type) {

                @Override
                String extractStringForHit(VDJCHit hit) {
                    return hit.getGene().getGeneName();
                }
            });
        }
        // All families
        for (final GeneType type : GeneType.values()) {
            char l = type.getLetter();
            descriptorsList.add(new StringExtractor("-" + Character.toLowerCase(l) + "Families", "Export all " + l + " gene family anmes (e.g. TRBV12 for TRBV12-3*00)", "All " + l + " families", "all" + l + "Families", type) {

                @Override
                String extractStringForHit(VDJCHit hit) {
                    return hit.getGene().getFamilyName();
                }
            });
        }
        // Best alignment
        for (final GeneType type : GeneType.values()) {
            char l = type.getLetter();
            descriptorsList.add(new PL_O("-" + Character.toLowerCase(l) + "Alignment", "Export best " + l + " alignment", "Best " + l + " alignment", "best" + l + "Alignment") {

                @Override
                protected String extract(VDJCObject object) {
                    VDJCHit bestHit = object.getBestHit(type);
                    if (bestHit == null)
                        return NULL;
                    StringBuilder sb = new StringBuilder();
                    for (int i = 0; ; i++) {
                        Alignment<NucleotideSequence> alignment = bestHit.getAlignment(i);
                        if (alignment == null)
                            sb.append(NULL);
                        else
                            sb.append(alignment.toCompactString());
                        if (i == object.numberOfTargets() - 1)
                            break;
                        sb.append(",");
                    }
                    return sb.toString();
                }
            });
        }
        // All alignments
        for (final GeneType type : GeneType.values()) {
            char l = type.getLetter();
            descriptorsList.add(new PL_O("-" + Character.toLowerCase(l) + "Alignments", "Export all " + l + " alignments", "All " + l + " alignments", "all" + l + "Alignments") {

                @Override
                protected String extract(VDJCObject object) {
                    VDJCHit[] hits = object.getHits(type);
                    if (hits.length == 0)
                        return "";
                    StringBuilder sb = new StringBuilder();
                    for (int j = 0; ; ++j) {
                        for (int i = 0; ; i++) {
                            Alignment<NucleotideSequence> alignment = hits[j].getAlignment(i);
                            if (alignment == null)
                                sb.append(NULL);
                            else
                                sb.append(alignment.toCompactString());
                            if (i == object.numberOfTargets() - 1)
                                break;
                            sb.append(',');
                        }
                        if (j == hits.length - 1)
                            break;
                        sb.append(';');
                    }
                    return sb.toString();
                }
            });
        }
        descriptorsList.add(new FeatureExtractors.NSeqExtractor("-nFeature", "Export nucleotide sequence of specified gene feature", "N. Seq. ", "nSeq") {

            @Override
            public String convert(NSequenceWithQuality seq) {
                return seq.getSequence().toString();
            }
        });
        descriptorsList.add(new FeatureExtractors.NSeqExtractor("-qFeature", "Export quality string of specified gene feature", "Qual. ", "qual") {

            @Override
            public String convert(NSequenceWithQuality seq) {
                return seq.getQuality().toString();
            }
        });
        descriptorsList.add(new FeatureExtractors.WithHeader("-aaFeature", "Export amino acid sequence of specified gene feature", 1, new String[] { "AA. Seq. " }, new String[] { "aaSeq" }) {

            @Override
            protected String extractValue(VDJCObject object, GeneFeature[] parameters) {
                GeneFeature geneFeature = parameters[parameters.length - 1];
                NSequenceWithQuality feature = object.getFeature(geneFeature);
                if (feature == null)
                    return NULL;
                int targetId = object.getTargetContainingFeature(geneFeature);
                TranslationParameters tr = targetId == -1 ? TranslationParameters.FromLeftWithIncompleteCodon : object.getPartitionedTarget(targetId).getPartitioning().getTranslationParameters(geneFeature);
                if (tr == null)
                    return NULL;
                return AminoAcidSequence.translate(feature.getSequence(), tr).toString();
            }
        });
        // descriptorsList.add(new FeatureExtractorDescriptor("-aaFeatureFromLeft", "Export amino acid sequence of " +
        // "specified gene feature starting from the leftmost nucleotide (differs from -aaFeature only for " +
        // "sequences which length are not multiple of 3)", "AA. Seq.", "aaSeq") {
        // @Override
        // public String convert(NSequenceWithQuality seq) {
        // return AminoAcidSequence.translate(seq.getSequence(), FromLeftWithoutIncompleteCodon).toString();
        // }
        // });
        // 
        // descriptorsList.add(new FeatureExtractorDescriptor("-aaFeatureFromRight", "Export amino acid sequence of " +
        // "specified gene feature starting from the rightmost nucleotide (differs from -aaFeature only for " +
        // "sequences which length are not multiple of 3)", "AA. Seq.", "aaSeq") {
        // @Override
        // public String convert(NSequenceWithQuality seq) {
        // return AminoAcidSequence.translate(seq.getSequence(), FromRightWithoutIncompleteCodon).toString();
        // }
        // });
        descriptorsList.add(new FeatureExtractors.NSeqExtractor("-minFeatureQuality", "Export minimal quality of specified gene feature", "Min. qual. ", "minQual") {

            @Override
            public String convert(NSequenceWithQuality seq) {
                return "" + seq.getQuality().minValue();
            }
        });
        descriptorsList.add(new FeatureExtractors.NSeqExtractor("-avrgFeatureQuality", "Export average quality of specified gene feature", "Mean. qual. ", "meanQual") {

            @Override
            public String convert(NSequenceWithQuality seq) {
                return "" + seq.getQuality().meanValue();
            }
        });
        descriptorsList.add(new FeatureExtractors.NSeqExtractor("-lengthOf", "Exports length of specified gene feature.", "Length of ", "lengthOf") {

            @Override
            public String convert(NSequenceWithQuality seq) {
                return "" + seq.size();
            }
        });
        descriptorsList.add(new FeatureExtractors.MutationsExtractor("-nMutations", "Extract nucleotide mutations for specific gene feature; relative to germline sequence.", 1, new String[] { "N. Mutations in " }, new String[] { "nMutations" }) {

            @Override
            String convert(Mutations<NucleotideSequence> mutations, NucleotideSequence seq1, NucleotideSequence seq2, TranslationParameters tr) {
                return mutations.encode(",");
            }
        });
        descriptorsList.add(new FeatureExtractors.MutationsExtractor("-nMutationsRelative", "Extract nucleotide mutations for specific gene feature relative to another feature.", 2, new String[] { "N. Mutations in ", " relative to " }, new String[] { "nMutationsIn", "Relative" }) {

            @Override
            String convert(Mutations<NucleotideSequence> mutations, NucleotideSequence seq1, NucleotideSequence seq2, TranslationParameters tr) {
                return mutations.encode(",");
            }
        });
        final class AAMutations extends FeatureExtractors.MutationsExtractor {

            AAMutations(String command, String description, int nArgs, String[] hPrefix, String[] sPrefix) {
                super(command, description, nArgs, hPrefix, sPrefix);
            }

            @Override
            String convert(Mutations<NucleotideSequence> mutations, NucleotideSequence seq1, NucleotideSequence seq2, TranslationParameters tr) {
                if (tr == null)
                    return "-";
                Mutations<AminoAcidSequence> aaMuts = MutationsUtil.nt2aa(seq1, mutations, tr);
                if (aaMuts == null)
                    return "-";
                return aaMuts.encode(",");
            }
        }
        descriptorsList.add(new AAMutations("-aaMutations", "Extract amino acid mutations for specific gene feature", 1, new String[] { "AA. Mutations in " }, new String[] { "aaMutations" }));
        descriptorsList.add(new AAMutations("-aaMutationsRelative", "Extract amino acid mutations for specific gene feature relative to another feature.", 2, new String[] { "AA. Mutations in ", " relative to " }, new String[] { "aaMutationsIn", "Relative" }));
        final class MutationsDetailed extends FeatureExtractors.MutationsExtractor {

            MutationsDetailed(String command, String description, int nArgs, String[] hPrefix, String[] sPrefix) {
                super(command, description, nArgs, hPrefix, sPrefix);
            }

            @Override
            String convert(Mutations<NucleotideSequence> mutations, NucleotideSequence seq1, NucleotideSequence seq2, TranslationParameters tr) {
                if (tr == null)
                    return "-";
                MutationsUtil.MutationNt2AADescriptor[] descriptors = MutationsUtil.nt2aaDetailed(seq1, mutations, tr, 10);
                if (descriptors == null)
                    return "-";
                StringBuilder sb = new StringBuilder();
                for (int i = 0; i < descriptors.length; i++) {
                    sb.append(descriptors[i]);
                    if (i == descriptors.length - 1)
                        break;
                    sb.append(",");
                }
                return sb.toString();
            }
        }
        String detailedMutationsFormat = "Format <nt_mutation>:<aa_mutation_individual>:<aa_mutation_cumulative>, where <aa_mutation_individual> is an expected amino acid " + "mutation given no other mutations have occurred, and <aa_mutation_cumulative> amino acid mutation is the observed amino acid " + "mutation combining effect from all other. WARNING: format may change in following versions.";
        descriptorsList.add(new MutationsDetailed("-mutationsDetailed", "Detailed list of nucleotide and corresponding amino acid mutations. " + detailedMutationsFormat, 1, new String[] { "Detailed mutations in " }, new String[] { "mutationsDetailedIn" }));
        descriptorsList.add(new MutationsDetailed("-mutationsDetailedRelative", "Detailed list of nucleotide and corresponding amino acid mutations written, positions relative to specified gene feature. " + detailedMutationsFormat, 2, new String[] { "Detailed mutations in ", " relative to " }, new String[] { "mutationsDetailedIn", "Relative" }));
        descriptorsList.add(new ExtractReferencePointPosition());
        descriptorsList.add(new ExtractDefaultReferencePointsPositions());
        descriptorsList.add(new PL_A("-readId", "Export id of read corresponding to alignment", "Read id", "readId") {

            @Override
            protected String extract(VDJCAlignments object) {
                return "" + object.getMinReadId();
            }

            @Override
            public FieldExtractor<VDJCAlignments> create(OutputMode outputMode, String[] args) {
                System.out.println("WARNING: -readId is deprecated. Use -readIds");
                return super.create(outputMode, args);
            }
        });
        descriptorsList.add(new PL_A("-readIds", "Export id of read corresponding to alignment", "Read id", "readId") {

            @Override
            protected String extract(VDJCAlignments object) {
                long[] readIds = object.getReadIds();
                StringBuilder sb = new StringBuilder();
                for (int i = 0; ; i++) {
                    sb.append(readIds[i]);
                    if (i == readIds.length - 1)
                        return sb.toString();
                    sb.append(",");
                }
            }
        });
        descriptorsList.add(new ExtractSequence(VDJCAlignments.class, "-sequence", "Export aligned sequence (initial read), or 2 sequences in case of paired-end reads", "Read(s) sequence", "readSequence"));
        descriptorsList.add(new ExtractSequenceQuality(VDJCAlignments.class, "-quality", "Export initial read quality, or 2 qualities in case of paired-end reads", "Read(s) sequence qualities", "readQuality"));
        descriptorsList.add(new PL_C("-cloneId", "Unique clone identifier", "Clone ID", "cloneId") {

            @Override
            protected String extract(Clone object) {
                return "" + object.getId();
            }
        });
        descriptorsList.add(new PL_C("-count", "Export clone count", "Clone count", "cloneCount") {

            @Override
            protected String extract(Clone object) {
                return "" + object.getCount();
            }
        });
        descriptorsList.add(new PL_C("-fraction", "Export clone fraction", "Clone fraction", "cloneFraction") {

            @Override
            protected String extract(Clone object) {
                return "" + object.getFraction();
            }
        });
        descriptorsList.add(new ExtractSequence(Clone.class, "-sequence", "Export aligned sequence (initial read), or 2 sequences in case of paired-end reads", "Clonal sequence(s)", "clonalSequence"));
        descriptorsList.add(new ExtractSequenceQuality(Clone.class, "-quality", "Export initial read quality, or 2 qualities in case of paired-end reads", "Clonal sequence quality(s)", "clonalSequenceQuality"));
        descriptorsList.add(new PL_A("-descrR1", "Export description line from initial .fasta or .fastq file " + "of the first read (only available if --save-description was used in align command)", "Description R1", "descrR1") {

            @Override
            protected String extract(VDJCAlignments object) {
                List<SequenceRead> reads = object.getOriginalReads();
                if (reads == null)
                    throw new IllegalArgumentException("Error for option \'-descrR1\':\n" + "No description available for read: either re-run align action with -OsaveOriginalReads=true option " + "or don't use \'-descrR1\' in exportAlignments");
                return reads.get(0).getRead(0).getDescription();
            }

            @Override
            public FieldExtractor<VDJCAlignments> create(OutputMode outputMode, String[] args) {
                System.out.println("WARNING: -descrR1 is deprecated. Use -descrsR1");
                return super.create(outputMode, args);
            }
        });
        descriptorsList.add(new PL_A("-descrR2", "Export description line from initial .fasta or .fastq file " + "of the second read (only available if --save-description was used in align command)", "Description R2", "descrR2") {

            @Override
            protected String extract(VDJCAlignments object) {
                List<SequenceRead> reads = object.getOriginalReads();
                if (reads == null)
                    throw new IllegalArgumentException("Error for option \'-descrR1\':\n" + "No description available for read: either re-run align action with -OsaveOriginalReads=true option " + "or don't use \'-descrR1\' in exportAlignments");
                SequenceRead read = reads.get(0);
                if (read.numberOfReads() < 2)
                    throw new IllegalArgumentException("Error for option \'-descrR2\':\n" + "No description available for second read: your input data was single-end");
                return read.getRead(1).getDescription();
            }

            @Override
            public FieldExtractor<VDJCAlignments> create(OutputMode outputMode, String[] args) {
                System.out.println("WARNING: -descrR2 is deprecated. Use -descrsR2");
                return super.create(outputMode, args);
            }
        });
        descriptorsList.add(new PL_A("-descrsR1", "Export description lines from initial .fasta or .fastq file " + "of the first reads (only available if -OsaveOriginalReads=true was used in align command)", "Descriptions R1", "descrsR1") {

            @Override
            protected String extract(VDJCAlignments object) {
                List<SequenceRead> reads = object.getOriginalReads();
                if (reads == null)
                    throw new IllegalArgumentException("Error for option \'-descrR1\':\n" + "No description available for read: either re-run align action with -OsaveOriginalReads option " + "or don't use \'-descrR1\' in exportAlignments");
                StringBuilder sb = new StringBuilder();
                for (int i = 0; ; i++) {
                    sb.append(reads.get(i).getRead(0).getDescription());
                    if (i == reads.size() - 1)
                        return sb.toString();
                    sb.append(",");
                }
            }
        });
        descriptorsList.add(new PL_A("-descrsR2", "Export description lines from initial .fasta or .fastq file " + "of the second reads (only available if -OsaveOriginalReads=true was used in align command)", "Descriptions R2", "descrsR2") {

            @Override
            protected String extract(VDJCAlignments object) {
                List<SequenceRead> reads = object.getOriginalReads();
                if (reads == null)
                    throw new IllegalArgumentException("Error for option \'-descrR1\':\n" + "No description available for read: either re-run align action with -OsaveOriginalReads option " + "or don't use \'-descrR1\' in exportAlignments");
                StringBuilder sb = new StringBuilder();
                for (int i = 0; ; i++) {
                    SequenceRead read = reads.get(i);
                    if (read.numberOfReads() < 2)
                        throw new IllegalArgumentException("Error for option \'-descrsR2\':\n" + "No description available for second read: your input data was single-end");
                    sb.append(read.getRead(1).getDescription());
                    if (i == reads.size() - 1)
                        return sb.toString();
                    sb.append(",");
                }
            }
        });
        descriptorsList.add(new PL_A("-readHistory", "Export read history", "Read history", "readHistory") {

            @Override
            protected String extract(VDJCAlignments object) {
                try {
                    return GlobalObjectMappers.toOneLine(object.getHistory());
                } catch (JsonProcessingException ex) {
                    throw new RuntimeException(ex);
                }
            }
        });
        for (final GeneType type : GeneType.values()) {
            String c = Character.toLowerCase(type.getLetter()) + "IdentityPercents";
            descriptorsList.add(new PL_O("-" + c, type.getLetter() + " alignment identity percents", type.getLetter() + " alignment identity percents", c) {

                @Override
                protected String extract(VDJCObject object) {
                    VDJCHit[] hits = object.getHits(type);
                    if (hits == null)
                        return NULL;
                    StringBuilder sb = new StringBuilder();
                    sb.append("");
                    for (int i = 0; ; i++) {
                        sb.append(hits[i].getIdentity());
                        if (i == hits.length - 1)
                            return sb.toString();
                        sb.append(",");
                    }
                }
            });
        }
        for (final GeneType type : GeneType.values()) {
            String c = Character.toLowerCase(type.getLetter()) + "BestIdentityPercent";
            descriptorsList.add(new PL_O("-" + c, type.getLetter() + "best alignment identity percent", type.getLetter() + "best alignment identity percent", c) {

                @Override
                protected String extract(VDJCObject object) {
                    VDJCHit hit = object.getBestHit(type);
                    if (hit == null)
                        return NULL;
                    return Float.toString(hit.getIdentity());
                }
            });
        }
        descriptorsList.add(new PL_O("-chains", "Chains", "Chains", "Chains") {

            @Override
            protected String extract(VDJCObject object) {
                return object.commonChains().toString();
            }
        });
        descriptorsList.add(new PL_O("-topChains", "Top chains", "Top chains", "topChains") {

            @Override
            protected String extract(VDJCObject object) {
                return object.commonTopChains().toString();
            }
        });
        descriptors = descriptorsList.toArray(new Field[descriptorsList.size()]);
    }
    return descriptors;
}

Example 4 with Clone

use of com.milaboratory.mixcr.basictypes.Clone in project mixcr by milaboratory.

the class CloneAssemblerRunnerTest method test1.

@Ignore
@Test
public void test1() throws Exception {
    String[] str = { "sequences/sample_IGH_R1.fastq", "sequences/sample_IGH_R2.fastq" };
    CloneSet cloneSet = runFullPipeline(str);
    System.out.println("\n\n");
    for (Clone clone : cloneSet) {
        System.out.println(clone);
        System.out.println(Arrays.toString(clone.getHits(GeneType.Variable)));
        System.out.println(Arrays.toString(clone.getHits(GeneType.Joining)));
        System.out.println(clone.getFeature(GeneFeature.CDR3));
        System.out.println("" + clone.getFeature(GeneFeature.VCDR3Part) + " " + clone.getFeature(GeneFeature.VJJunction) + " " + clone.getFeature(GeneFeature.JCDR3Part));
    }
}

Example 5 with Clone

use of com.milaboratory.mixcr.basictypes.Clone in project mixcr by milaboratory.

the class FullSeqAssemblerTest method testRandom1.

@Test
public void testRandom1() throws Exception {
    CloneFraction[] clones = { new CloneFraction(750, masterSeq1WT), // V: S346:G->T
    new CloneFraction(1000, masterSeq1VSub1), // J: D55:A
    new CloneFraction(1000, masterSeq1VDel1JDel1), // J: D62:C
    new CloneFraction(500, masterSeq1VDel1JDelVSub2) };
    Well19937c rand = new Well19937c();
    rand.setSeed(12345);
    RandomDataGenerator rdg = new RandomDataGenerator(rand);
    List<SequenceRead> readsOrig = new ArrayList<>();
    int readLength = 100;
    int id = -1;
    for (CloneFraction clone : clones) {
        for (int i = 0; i < clone.count; i++) {
            // Left read with CDR3
            ++id;
            readsOrig.add(new PairedRead(new SingleReadImpl(id, new NSequenceWithQuality(clone.seq.getRangeFromCDR3Begin(-rand.nextInt(readLength - clone.seq.cdr3Part), readLength)), "R1_" + id), new SingleReadImpl(id, new NSequenceWithQuality(clone.seq.getRangeFromCDR3End(rdg.nextInt(-clone.seq.cdr3Part / 2, clone.seq.jPart), readLength).getReverseComplement()), "R2_" + id)));
            ++id;
            readsOrig.add(new PairedRead(new SingleReadImpl(id, new NSequenceWithQuality(clone.seq.getRangeFromCDR3Begin(rdg.nextInt(-clone.seq.vPart, clone.seq.cdr3Part / 2 - readLength), readLength)), "R1_" + id), new SingleReadImpl(id, new NSequenceWithQuality(clone.seq.getRangeFromCDR3Begin(-rand.nextInt(readLength - clone.seq.cdr3Part), readLength)).getReverseComplement(), "R2_" + id)));
        }
    }
    // readsOrig = Arrays.asList(setReadId(0, readsOrig.get(12)), setReadId(1, readsOrig.get(13)));
    int[] perm = rdg.nextPermutation(readsOrig.size(), readsOrig.size());
    List<SequenceRead> reads = new ArrayList<>();
    for (int i = 0; i < readsOrig.size(); i++) reads.add(readsOrig.get(perm[i]));
    RunMiXCR.RunMiXCRAnalysis params = new RunMiXCR.RunMiXCRAnalysis(new SequenceReaderCloseable<SequenceRead>() {

        int counter = 0;

        @Override
        public void close() {
        }

        @Override
        public long getNumberOfReads() {
            return counter;
        }

        @Override
        public synchronized SequenceRead take() {
            if (counter == reads.size())
                return null;
            return reads.get(counter++);
        }
    }, true);
    params.alignerParameters = VDJCParametersPresets.getByName("rna-seq");
    params.alignerParameters.setSaveOriginalReads(true);
    params.alignerParameters.setVAlignmentParameters(params.alignerParameters.getVAlignerParameters().setGeneFeatureToAlign(GeneFeature.VTranscriptWithP));
    RunMiXCR.AlignResult align = RunMiXCR.align(params);
    // // TODO exception for translation
    // for (VDJCAlignments al : align.alignments) {
    // for (int i = 0; i < al.numberOfTargets(); i++) {
    // System.out.println(VDJCAlignmentsFormatter.getTargetAsMultiAlignment(al, i));
    // System.out.println();
    // }
    // System.out.println();
    // System.out.println(" ================================================ ");
    // System.out.println();
    // }
    RunMiXCR.AssembleResult assemble = RunMiXCR.assemble(align);
    Assert.assertEquals(1, assemble.cloneSet.size());
    CloneFactory cloneFactory = new CloneFactory(align.parameters.cloneAssemblerParameters.getCloneFactoryParameters(), align.parameters.cloneAssemblerParameters.getAssemblingFeatures(), align.usedGenes, align.parameters.alignerParameters.getFeaturesToAlignMap());
    FullSeqAssembler agg = new FullSeqAssembler(cloneFactory, DEFAULT_PARAMETERS, assemble.cloneSet.get(0), align.parameters.alignerParameters);
    FullSeqAssembler.RawVariantsData prep = agg.calculateRawData(() -> CUtils.asOutputPort(align.alignments.stream().filter(a -> a.getFeature(GeneFeature.CDR3) != null).collect(Collectors.toList())));
    List<Clone> clns = new ArrayList<>(new CloneSet(Arrays.asList(agg.callVariants(prep))).getClones());
    clns.sort(Comparator.comparingDouble(Clone::getCount).reversed());
    System.out.println("# Clones: " + clns.size());
    id = 0;
    for (Clone clone : clns) {
        clone = clone.setId(id++);
        System.out.println(clone.numberOfTargets());
        System.out.println(clone.getCount());
        System.out.println(clone.getFraction());
        System.out.println(clone.getBestHit(GeneType.Variable).getAlignment(0).getAbsoluteMutations());
        System.out.println(clone.getBestHit(GeneType.Joining).getAlignment(0).getAbsoluteMutations());
        System.out.println();
    // ActionExportClonesPretty.outputCompact(System.out, clone);
    }
}