This is an automated email from the git hooks/post-receive script. bob.dybian-guest pushed a commit to branch master in repository beagle.
commit 1f3da147909be8f61b3a2df9880d50ed2b3726e2 Author: Dylan Aïssi <[email protected]> Date: Tue Feb 7 22:57:07 2017 +0100 Imported Upstream version 4.1~170121-6cc+dfsg --- ibd/HaploidIbd.java | 31 +++++++++++++++++++++++-------- ibd/IbsHapSegments.java | 11 +++++------ main/Main.java | 6 +++--- main/MainHelper.java | 20 ++++++++++---------- main/Par.java | 13 ++++++++++++- main/WindowWriter.java | 8 ++++---- vcf/VcfWindow.java | 3 +-- 7 files changed, 58 insertions(+), 34 deletions(-) diff --git a/ibd/HaploidIbd.java b/ibd/HaploidIbd.java index e4575e2..06f606b 100644 --- a/ibd/HaploidIbd.java +++ b/ibd/HaploidIbd.java @@ -34,6 +34,7 @@ import java.util.concurrent.ConcurrentMap; import java.util.concurrent.ExecutorService; import java.util.concurrent.Executors; import java.util.concurrent.TimeUnit; +import main.GeneticMap; import vcf.GL; /** @@ -52,32 +53,46 @@ import vcf.GL; */ public final class HaploidIbd { + private final GeneticMap genMap; private final int ibdTrim; private final float minIbdLod; - private final float minIbsLength; // positions from Dag.posArray() - private final float minFreqLod; // for shared haplotype + private final float minFreqLod; + private final float minCm; /** * Constructs a new {@code HaploidIbd} instance from the specified data. + * @param genMap the genetic map * @param ibdTrim the number of markers to trim from an IBS segment * when computing the IBD versus non-IBD likelihood ratio * @param minIbdLod the minimum IBD LOD score of reported IBD segments + * @param minCm the minimum cM length of reported IBD segments * * @throws IllegalArgumentException if {@code ibdTrim < 0 } * @throws IllegalArgumentException if * {@code ibdLod <= 0.0f || Float.isFinite(ibdLod) == false} + * @throws IllegalArgumentException if + * {@code minCm <= 0.0f || Float.isFinite(minCm) == false} + * @throws NullPointerException if {@code genMap == null} */ - public HaploidIbd(int ibdTrim, float minIbdLod) { + public HaploidIbd(GeneticMap genMap, int ibdTrim, float minIbdLod, + float minCm) { + if (genMap==null) { + throw new IllegalArgumentException(GeneticMap.class.toString()); + } if (ibdTrim < 0) { - throw new IllegalArgumentException("trim: " + ibdTrim); + throw new IllegalArgumentException(String.valueOf(ibdTrim)); } if (minIbdLod <= 0.0 || Float.isFinite(minIbdLod) == false) { - throw new IllegalArgumentException("minIbdlod: " + minIbdLod); + throw new IllegalArgumentException(String.valueOf(minIbdLod)); + } + if (minCm <= 0.0 || Float.isFinite(minCm) == false) { + throw new IllegalArgumentException(String.valueOf(minCm)); } + this.genMap = genMap; this.ibdTrim = ibdTrim; this.minIbdLod = minIbdLod; - this.minIbsLength = 0.8f*minIbdLod; this.minFreqLod = minIbdLod; + this.minCm = minCm; } /** @@ -105,8 +120,8 @@ public final class HaploidIbd { public Map<IntPair, List<IbdSegment>> run(GL gl, Dag dag, SampleHapPairs haps, final int nThreads) { checkParameters(gl, dag, haps); - double[] pos = dag.posArray(); - IbsHapSegments ibsSegments = new IbsHapSegments(haps, pos, minIbsLength); + double[] pos = genMap.genPos(dag.markers()); + IbsHapSegments ibsSegments = new IbsHapSegments(haps, pos, minCm); ConcurrentMap<IntPair, List<IbdSegment>> ibdMap = new ConcurrentHashMap<>(); diff --git a/ibd/IbsHapSegments.java b/ibd/IbsHapSegments.java index 8579264..234b37b 100644 --- a/ibd/IbsHapSegments.java +++ b/ibd/IbsHapSegments.java @@ -64,9 +64,8 @@ public final class IbsHapSegments { * * @throws IllegalArgumentException if * {@code haps.nMarkers() != pos.length} - * @throws IllegalArgumentException if {@code pos[0] < 0}, or if - * {@code pos[j] < pos[j-1]} for any {@code j} satisfing - * {@code (0 < j && j < pos.length)} + * @throws IllegalArgumentException if {@code pos[j] < pos[j-1]} + * for any {@code j} satisfing {@code (0 < j && j < pos.length)} * @throws IllegalArgumentException if * {@code (Double.isNaN(pos[j])==true || Double.isInfinite(pos[j]) == true)} * for any {@code j} satisfying {@code (0 <= j && j < pos.length)} @@ -110,12 +109,12 @@ public final class IbsHapSegments { if (haps.nMarkers()!= pos.length) { throw new IllegalArgumentException("haps.nMarkers()!= pos.length"); } - if (pos[0]<0 || Double.isNaN(pos[0]) || Double.isInfinite(pos[0]) ) { - throw new IllegalArgumentException("pos=" + pos[0]); + if (Double.isNaN(pos[0]) || Double.isInfinite(pos[0]) ) { + throw new IllegalArgumentException(String.valueOf(pos[0])); } for (int j=1; j<pos.length; ++j) { if (Double.isNaN(pos[j]) || Double.isInfinite(pos[j]) ) { - throw new IllegalArgumentException("pos=" + pos[j]); + throw new IllegalArgumentException(String.valueOf(pos[j])); } if (pos[j] < pos[j-1]) { String s = "positions are not non-decreasing"; diff --git a/main/Main.java b/main/Main.java index 30c0eb3..64f7d71 100644 --- a/main/Main.java +++ b/main/Main.java @@ -65,8 +65,8 @@ public class Main { /** * The program name and version. */ - public static final String program = "beagle.27Jul16.86a.jar (version 4.1)"; - public static final String command = "java -jar beagle.27Jul16.86a.jar"; + public static final String program = "beagle.21Jan17.6cc.jar (version 4.1)"; + public static final String command = "java -jar beagle.21Jan17.6cc.jar"; /** * The copyright string. @@ -78,7 +78,7 @@ public class Main { */ public static final String shortHelp = Main.program + Const.nl + Main.copyright - + Const.nl + "Enter \"java -jar beagle.27Jul16.86a.jar\" for a " + + Const.nl + "Enter \"java -jar beagle.21Jan17.6cc.jar\" for a " + "summary of command line " + "arguments."; private final Par par; diff --git a/main/MainHelper.java b/main/MainHelper.java index dd68eb4..b1ace1e 100644 --- a/main/MainHelper.java +++ b/main/MainHelper.java @@ -58,15 +58,20 @@ public class MainHelper { /** * Constructs a new {@code MainHelper} instance. * @param par the command line parameters - * @param genMap the genetic map + * @param genMap the genetic map or {@code null} if no genetic map + * was specified * @param runStats the class for collecting and printing run-time statistics * @throws NullPointerException - * if {@code (par == null || genMap == null || runStarts == null)} + * if {@code (par == null || runStarts == null)} */ MainHelper(Par par, GeneticMap genMap, RunStats runStats) { if (runStats==null) { throw new NullPointerException("runStats==null"); } + if (genMap==null) { + double scaleFactor = 1e-6; + genMap = new PositionMap(scaleFactor); + } this.par = par; this.hapSampler = new HapPairSampler(par, runStats); this.recombSampler = new RecombHapPairSampler(par, runStats); @@ -188,7 +193,8 @@ public class MainHelper { float scale = par.adjustedIbdScale(nSamples); Dag dag = ibdDag(cd, targetHapPairs, scale); - HaploidIbd hapIbd = new HaploidIbd(par.ibdtrim(), par.ibdlod()); + HaploidIbd hapIbd = new HaploidIbd(genMap, par.ibdtrim(), + par.ibdlod(), par.ibdcm()); GL ibdGL = new NoPhaseGL(cd.targetGL()); Map<IntPair, List<IbdSegment>> ibdMap = @@ -244,14 +250,8 @@ public class MainHelper { return new SampleHapPairAlleleProbs(shp); } long t0 = System.nanoTime(); - GeneticMap imputationMap = genMap; - if (par.map()==null) { - double scaleFactor = 1e-6; - imputationMap = new PositionMap(scaleFactor); - } - LiAndStephensHapSampler recombHapSampler = - new LiAndStephensHapSampler(par, imputationMap); + new LiAndStephensHapSampler(par, genMap); BasicAlleleProbs alProbs = recombHapSampler.sample(cd, shp); runStats.imputationNanos(System.nanoTime() - t0); diff --git a/main/Par.java b/main/Par.java index bd72cef..eb0abc4 100644 --- a/main/Par.java +++ b/main/Par.java @@ -66,6 +66,7 @@ public final class Par { // ibd parameters private final boolean ibd; private final float ibdlod; + private final float ibdcm; private final float ibdscale; private final int ibdtrim; @@ -137,6 +138,7 @@ public final class Par { // ibd parameters ibd = Validate.booleanArg("ibd", argsMap, false, false); ibdlod = Validate.floatArg("ibdlod", argsMap, false, 3.0f, FMIN, FMAX); + ibdcm = Validate.floatArg("ibdcm", argsMap, false, 1.5f, FMIN, FMAX); ibdscale = Validate.floatArg("ibdscale", argsMap, false, 0.0f, 0.0f, FMAX); ibdtrim = Validate.intArg("ibdtrim", argsMap, false, 40, 0, IMAX); @@ -199,7 +201,8 @@ public final class Par { + "IBD parameters ..." + nl + " ibd=<perform IBD detection (true/false)> (default=false)" + nl - + " ibdlod=<min LOD score for reporting IBD> (default=3.0)" + nl + + " ibdlod=<min LOD score of reported IBD segments> (default=3.0)" + nl + + " ibdcm=<min cM length of reported IBD segments> (default=1.5)" + nl + " ibdscale=<model scale factor for Refined IBD> (default: data-dependent)" + nl + " ibdtrim=<markers at each segment end> (default=40)" + nl; } @@ -473,6 +476,14 @@ public final class Par { } /** + * Returns the ibdcm parameter. + * @return the ibdcm parameter + */ + public float ibdcm() { + return ibdcm; + } + + /** * Returns the ibdscale parameter. * @return the ibdscale parameter */ diff --git a/main/WindowWriter.java b/main/WindowWriter.java index 2a7d6a8..1edf341 100644 --- a/main/WindowWriter.java +++ b/main/WindowWriter.java @@ -85,8 +85,8 @@ public class WindowWriter implements Closeable { this.samples = samples; this.outPrefix = outPrefix; this.vcfOutFile = new File(outPrefix + ".vcf.gz"); - this.ibdOutFile = new File(outPrefix + ".ibd"); - this.hbdOutFile = new File(outPrefix + ".hbd"); + this.ibdOutFile = new File(outPrefix + ".ibd.gz"); + this.hbdOutFile = new File(outPrefix + ".hbd.gz"); ByteArrayOutputStream baos = new ByteArrayOutputStream(); try (PrintWriter vcfOut=new PrintWriter( @@ -294,8 +294,8 @@ public class WindowWriter implements Closeable { int nextOverlap, int nextSplice, int nMarkers) { Map<IntPair, IbdSegment> lastBuffer = new HashMap<>(ibdBuffer); ibdBuffer.clear(); - try (PrintWriter ibdOut = FileUtil.printWriter(ibdOutFile, appendIbd); - PrintWriter hbdOut = FileUtil.printWriter(hbdOutFile, appendIbd)) { + try (PrintWriter ibdOut = FileUtil.bgzipPrintWriter(ibdOutFile, appendIbd); + PrintWriter hbdOut = FileUtil.bgzipPrintWriter(hbdOutFile, appendIbd)) { Iterator<IntPair> keyIt = ibd.keySet().iterator(); while (keyIt.hasNext()) { IntPair key = keyIt.next(); diff --git a/vcf/VcfWindow.java b/vcf/VcfWindow.java index 7f61c7e..3f8b155 100644 --- a/vcf/VcfWindow.java +++ b/vcf/VcfWindow.java @@ -24,7 +24,6 @@ import java.io.Closeable; import java.io.File; import java.util.ArrayList; import java.util.List; -import main.GeneticMap; /** * <p>Class {@code VcfWindow} represents a sliding window of VCF records. @@ -51,7 +50,7 @@ public class VcfWindow implements Closeable { */ public VcfWindow(SampleFileIt<? extends VcfEmission> it) { if (it.hasNext()==false) { - throw new IllegalArgumentException("it.hasNext()==false"); + throw new IllegalArgumentException("No VCF records after filtering"); } this.it = it; this.overlap = 0; -- Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-med/beagle.git _______________________________________________ debian-med-commit mailing list [email protected] http://lists.alioth.debian.org/cgi-bin/mailman/listinfo/debian-med-commit
