Skip to content

Commit 6cba97b

Browse files
committed
Support specHLA
1 parent 6889e4d commit 6cba97b

2 files changed

Lines changed: 127 additions & 0 deletions

File tree

SequenceAnalysis/src/org/labkey/sequenceanalysis/SequenceAnalysisModule.java

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -127,6 +127,7 @@
127127
import org.labkey.sequenceanalysis.run.analysis.SawfishJointCallingHandler;
128128
import org.labkey.sequenceanalysis.run.analysis.SequenceBasedTypingAnalysis;
129129
import org.labkey.sequenceanalysis.run.analysis.SnpCountAnalysis;
130+
import org.labkey.sequenceanalysis.run.analysis.SpecHlaAnalysis;
130131
import org.labkey.sequenceanalysis.run.analysis.SubreadAnalysis;
131132
import org.labkey.sequenceanalysis.run.analysis.UnmappedReadExportHandler;
132133
import org.labkey.sequenceanalysis.run.analysis.ViralAnalysis;
@@ -343,6 +344,7 @@ public static void registerPipelineSteps()
343344
SequencePipelineService.get().registerPipelineStep(new PindelAnalysis.Provider());
344345
SequencePipelineService.get().registerPipelineStep(new PbsvAnalysis.Provider());
345346
SequencePipelineService.get().registerPipelineStep(new GenrichStep.Provider());
347+
SequencePipelineService.get().registerPipelineStep(new SpecHlaAnalysis.Provider());
346348
SequencePipelineService.get().registerPipelineStep(new SawfishAnalysis.Provider());
347349

348350
SequencePipelineService.get().registerPipelineStep(new PARalyzerAnalysis.Provider());
Lines changed: 125 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,125 @@
1+
package org.labkey.sequenceanalysis.run.analysis;
2+
3+
import htsjdk.samtools.SAMFileHeader;
4+
import org.labkey.api.pipeline.PipelineJobException;
5+
import org.labkey.api.sequenceanalysis.model.AnalysisModel;
6+
import org.labkey.api.sequenceanalysis.model.Readset;
7+
import org.labkey.api.sequenceanalysis.pipeline.AbstractAnalysisStepProvider;
8+
import org.labkey.api.sequenceanalysis.pipeline.AnalysisOutputImpl;
9+
import org.labkey.api.sequenceanalysis.pipeline.AnalysisStep;
10+
import org.labkey.api.sequenceanalysis.pipeline.PipelineContext;
11+
import org.labkey.api.sequenceanalysis.pipeline.PipelineStepProvider;
12+
import org.labkey.api.sequenceanalysis.pipeline.ReferenceGenome;
13+
import org.labkey.api.sequenceanalysis.pipeline.SamSorter;
14+
import org.labkey.api.sequenceanalysis.pipeline.SamtoolsRunner;
15+
import org.labkey.api.sequenceanalysis.pipeline.SequencePipelineService;
16+
import org.labkey.api.sequenceanalysis.run.AbstractCommandPipelineStep;
17+
import org.labkey.api.sequenceanalysis.run.AbstractCommandWrapper;
18+
import org.labkey.api.sequenceanalysis.run.SimpleScriptWrapper;
19+
import org.labkey.api.util.FileUtil;
20+
import org.labkey.api.util.Path;
21+
22+
import java.io.File;
23+
import java.util.ArrayList;
24+
import java.util.Arrays;
25+
import java.util.List;
26+
27+
public class SpecHlaAnalysis extends AbstractCommandPipelineStep<SimpleScriptWrapper> implements AnalysisStep
28+
{
29+
public SpecHlaAnalysis(PipelineStepProvider<?> provider, PipelineContext ctx)
30+
{
31+
super(provider, ctx, new SimpleScriptWrapper(ctx.getLogger()));
32+
}
33+
34+
public static class Provider extends AbstractAnalysisStepProvider<SpecHlaAnalysis>
35+
{
36+
public Provider()
37+
{
38+
super("SpecHlaStep", "SpecHLA", null, "This will run SpecHLA for HLA genotyping from WGS/WXS data. This should use a BAM aligned to a custom HLA DB, rather than aligned to the full genome", Arrays.asList(
39+
40+
), null, "https://github.com/deepomicslab/SpecHLA/");
41+
}
42+
43+
@Override
44+
public SpecHlaAnalysis create(PipelineContext ctx)
45+
{
46+
return new SpecHlaAnalysis(this, ctx);
47+
}
48+
}
49+
50+
@Override
51+
public Output performAnalysisPerSampleRemote(Readset rs, File inputBam, ReferenceGenome referenceGenome, File outputDir) throws PipelineJobException
52+
{
53+
AnalysisOutputImpl output = new AnalysisOutputImpl();
54+
55+
File subsetBam = FileUtil.appendName(outputDir, FileUtil.getBaseName(inputBam) + ".subset.bam");
56+
SamtoolsRunner sr = new SamtoolsRunner(getWrapper().getLogger());
57+
sr.execute(Arrays.asList(
58+
sr.getSamtoolsPath().getPath(),
59+
"view",
60+
"-h",
61+
"-F", "12", //This selects pairs where either mate is mapped
62+
"-o", subsetBam.getPath(),
63+
inputBam.getPath()
64+
));
65+
output.addIntermediateFile(subsetBam);
66+
67+
File queryNameSortBam = new SamSorter(getPipelineCtx().getLogger()).execute(subsetBam, FileUtil.appendName(outputDir, FileUtil.getBaseName(inputBam) + ".querySort.bam"), SAMFileHeader.SortOrder.queryname);
68+
output.addIntermediateFile(queryNameSortBam);
69+
70+
File fq1 = FileUtil.appendName(outputDir, FileUtil.getBaseName(inputBam) + ".R1.fastq.gz");
71+
File fq2 = FileUtil.appendName(outputDir, FileUtil.getBaseName(inputBam) + ".R2.fastq.gz");
72+
sr.execute(Arrays.asList(
73+
sr.getSamtoolsPath().getPath(),
74+
"fastq",
75+
"-1",
76+
fq1.getPath(),
77+
"-2",
78+
fq2.getPath(),
79+
queryNameSortBam.getPath()
80+
));
81+
output.addIntermediateFile(fq1);
82+
output.addIntermediateFile(fq2);
83+
84+
File specHlaExe = AbstractCommandWrapper.resolveFileInPath("spechla", null, true);
85+
86+
List<String> toRun = new ArrayList<>(Arrays.asList(
87+
specHlaExe.getPath(),
88+
"-n",
89+
"specHLA",
90+
"-u",
91+
"1", // 1 = exon. 0 = full-length
92+
"-1",
93+
fq1.getPath(),
94+
"-2",
95+
fq2.getPath(),
96+
"-o",
97+
outputDir.getPath()
98+
));
99+
100+
Integer maxThreads = SequencePipelineService.get().getMaxThreads(getWrapper().getLogger());
101+
if (maxThreads != null)
102+
{
103+
toRun.add("-j");
104+
toRun.add(maxThreads.toString());
105+
}
106+
107+
getWrapper().execute(toRun);
108+
109+
File outFile = FileUtil.appendPath(outputDir, Path.parse("specHLA/hla.result.txt"));
110+
if (!outFile.exists())
111+
{
112+
throw new PipelineJobException("SpecHLA result file does not exist: " + outFile.getPath());
113+
}
114+
115+
output.addSequenceOutput(outFile, FileUtil.getBaseName(inputBam) + ": HLA Typing", "specHLA Genotyping", rs.getReadsetId(), null, referenceGenome.getGenomeId(), null);
116+
117+
return output;
118+
}
119+
120+
@Override
121+
public Output performAnalysisPerSampleLocal(AnalysisModel model, File inputBam, File referenceFasta, File outDir) throws PipelineJobException
122+
{
123+
return null;
124+
}
125+
}

0 commit comments

Comments
 (0)