use of com.milaboratory.mixcr.vdjaligners.VDJCAlignerParameters in project mixcr by milaboratory.
the class ActionAssemble method go.
@Override
public void go(ActionHelper helper) throws Exception {
// Saving initial timestamp
long beginTimestamp = System.currentTimeMillis();
// Checking consistency between actionParameters.doWriteClnA() value and file extension
if ((actionParameters.getOutputFileName().toLowerCase().endsWith(".clna") && !actionParameters.doWriteClnA()) || (actionParameters.getOutputFileName().toLowerCase().endsWith(".clns") && actionParameters.doWriteClnA()))
System.out.println("WARNING: Unexpected file extension, use .clns extension for clones-only (normal) output and\n" + ".clna if -a / --write-alignments options specified.");
// Extracting V/D/J/C gene list from input vdjca file
final List<VDJCGene> genes;
final VDJCAlignerParameters alignerParameters;
try (VDJCAlignmentsReader reader = new VDJCAlignmentsReader(actionParameters.getInputFileName(), VDJCLibraryRegistry.getDefault())) {
genes = reader.getUsedGenes();
// Saving aligner parameters to correct assembler parameters
alignerParameters = reader.getParameters();
}
AlignmentsProvider alignmentsProvider = AlignmentsProvider.Util.createProvider(actionParameters.getInputFileName(), VDJCLibraryRegistry.getDefault());
CloneAssemblerParameters assemblerParameters = actionParameters.getCloneAssemblerParameters();
// set aligner parameters
assemblerParameters.updateFrom(alignerParameters);
// Overriding JSON parameters
if (!actionParameters.overrides.isEmpty()) {
assemblerParameters = JsonOverrider.override(assemblerParameters, CloneAssemblerParameters.class, actionParameters.overrides);
if (assemblerParameters == null) {
System.err.println("Failed to override some parameter.");
System.exit(1);
}
}
// Performing assembly
try (CloneAssembler assembler = new CloneAssembler(assemblerParameters, actionParameters.doWriteClnA() || actionParameters.events != null, genes, alignerParameters.getFeaturesToAlignMap())) {
// Creating event listener to collect run statistics
CloneAssemblerReport report = new CloneAssemblerReport();
report.setStartMillis(beginTimestamp);
report.setInputFiles(new String[] { actionParameters.getInputFileName() });
report.setOutputFiles(new String[] { actionParameters.getOutputFileName() });
report.setCommandLine(helper.getCommandLineArguments());
assembler.setListener(report);
// Running assembler
CloneAssemblerRunner assemblerRunner = new CloneAssemblerRunner(alignmentsProvider, assembler, actionParameters.threads);
SmartProgressReporter.startProgressReport(assemblerRunner);
assemblerRunner.run();
// Getting results
final CloneSet cloneSet = assemblerRunner.getCloneSet(alignerParameters);
// Passing final cloneset to assemble last pieces of statistics for report
report.onClonesetFinished(cloneSet);
// Writing results
if (actionParameters.doWriteClnA())
try (ClnAWriter writer = new ClnAWriter(actionParameters.getOutputFileName())) {
// writer will supply current stage and completion percent to the progress reporter
SmartProgressReporter.startProgressReport(writer);
// Writing clone block
writer.writeClones(cloneSet);
// Pre-soring alignments
try (AlignmentsMappingMerger merged = new AlignmentsMappingMerger(alignmentsProvider.create(), assembler.getAssembledReadsPort())) {
writer.sortAlignments(merged, assembler.getAlignmentsCount());
}
writer.writeAlignmentsAndIndex();
}
else
try (CloneSetIO.CloneSetWriter writer = new CloneSetIO.CloneSetWriter(cloneSet, actionParameters.getOutputFileName())) {
SmartProgressReporter.startProgressReport(writer);
writer.write();
}
// Writing report
report.setFinishMillis(System.currentTimeMillis());
assert cloneSet.getClones().size() == report.getCloneCount();
report.setTotalReads(alignmentsProvider.getTotalNumberOfReads());
// Writing report to stout
System.out.println("============= Report ==============");
Util.writeReportToStdout(report);
if (actionParameters.report != null)
Util.writeReport(actionParameters.report, report);
if (actionParameters.jsonReport != null)
Util.writeJsonReport(actionParameters.jsonReport, report);
// Writing raw events (not documented feature)
if (actionParameters.events != null)
try (PipeWriter<ReadToCloneMapping> writer = new PipeWriter<>(actionParameters.events)) {
CUtils.drain(assembler.getAssembledReadsPort(), writer);
}
}
}
use of com.milaboratory.mixcr.vdjaligners.VDJCAlignerParameters in project mixcr by milaboratory.
the class PartialAlignmentsAssemblerTest method createTestData.
public static InputTestData createTestData(long seed) throws Exception {
EnumMap<GeneType, String> geneNames = new EnumMap<GeneType, String>(GeneType.class) {
{
put(Variable, "TRBV20-1*00");
put(Diversity, "TRBD2*00");
put(Joining, "TRBJ2-6*00");
put(Constant, "TRBC2*00");
}
};
// config
RandomGenerator rnd = RandomUtil.getThreadLocalRandom();
rnd.setSeed(seed);
final VDJCAlignerParameters defaultFeatures = VDJCParametersPresets.getByName("default");
defaultFeatures.getVAlignerParameters().setGeneFeatureToAlign(VRegion);
defaultFeatures.getDAlignerParameters().setGeneFeatureToAlign(DRegion);
defaultFeatures.getJAlignerParameters().setGeneFeatureToAlign(JRegion);
// used alleles
EnumMap<GeneType, VDJCGene> genes = new EnumMap<>(GeneType.class);
// germline parts of sequences
EnumMap<GeneType, NucleotideSequence> germlineRegions = gtMap();
// left, right cut of germline
EnumMap<GeneType, int[]> germlineCuts = gtMap();
// begin, end positions in assembled sequence
EnumMap<GeneType, int[]> refPositions = gtMap();
// single assembled sequence
SequenceBuilder<NucleotideSequence> referenceBuilder = NucleotideSequence.ALPHABET.createBuilder();
NucleotideSequence VDJunction = TestUtil.randomSequence(NucleotideSequence.ALPHABET, 3, 10);
NucleotideSequence DJJunction = TestUtil.randomSequence(NucleotideSequence.ALPHABET, 3, 10);
for (GeneType gt : GeneType.VDJC_REFERENCE) {
VDJCGene gene = VDJCLibraryRegistry.getDefault().getLibrary("default", "hs").get(geneNames.get(gt));
NucleotideSequence seq = gene.getFeature(defaultFeatures.getFeatureToAlign(gt));
int[] cuts = null;
switch(gt) {
case Variable:
cuts = new int[] { 0, rnd.nextInt(gene.getFeature(GermlineVCDR3Part).size() - 5) };
break;
case Diversity:
cuts = new int[] { rnd.nextInt(seq.size() / 3), rnd.nextInt(seq.size() / 3) };
break;
case Joining:
cuts = new int[] { rnd.nextInt(gene.getFeature(GermlineJCDR3Part).size() - 5), 0 };
break;
case Constant:
cuts = new int[] { 0, rnd.nextInt(seq.size() / 2) };
break;
}
NucleotideSequence gSeq = seq.getRange(cuts[0], seq.size() - cuts[1]);
int[] positions = new int[2];
positions[0] = referenceBuilder.size();
referenceBuilder.append(gSeq);
positions[1] = referenceBuilder.size();
if (gt == Variable)
referenceBuilder.append(VDJunction);
if (gt == Diversity)
referenceBuilder.append(DJJunction);
genes.put(gt, gene);
germlineCuts.put(gt, cuts);
germlineRegions.put(gt, gSeq);
refPositions.put(gt, positions);
}
NucleotideSequence VJJunction = NucleotideSequence.ALPHABET.createBuilder().append(VDJunction).append(germlineRegions.get(Diversity)).append(DJJunction).createAndDestroy();
NucleotideSequence reference = referenceBuilder.createAndDestroy();
return new InputTestData(genes, germlineRegions, germlineCuts, refPositions, VDJunction, DJJunction, reference, VJJunction);
}
Aggregations