working version v0.3

Author: ddecap

halvade/src/be/ugent/intec/halvade/HalvadeOptions.java

@@ -81,7 +81,7 @@ public class HalvadeOptions {
     protected boolean reuseJVM = false;
     protected boolean justAlign = false;
     protected String exomeBedFile = null;
-    protected double coverage = 50; 
+    protected double coverage = -1.0; 
     protected String halvadeBinaries;
     protected String bin;
     protected boolean combineVcf = true;
@@ -143,7 +143,8 @@ public int GetOptions(String[] args, Configuration hConf) throws IOException, UR
 
             parseDictFile(hConf);
             double inputSize = getInputSize(in, hConf);
-            coverage = DEFAULT_COVERAGE * (inputSize / DEFAULT_COVERAGE_SIZE);
+            if(coverage == -1.0)
+                coverage = Math.max(1.0, DEFAULT_COVERAGE * (inputSize / DEFAULT_COVERAGE_SIZE));
             Logger.DEBUG("Estimated coverage: " + roundOneDecimal(coverage));
             // set a minimum first where the real amount is based on
             reduces = (int) (coverage * REDUCE_TASKS_FACTOR);
@@ -185,7 +186,7 @@ protected double getInputSize(String input, Configuration conf) throws URISyntax
 
     private static final String DICT_SUFFIX = ".dict";
     private void parseDictFile(Configuration conf) {
-        be.ugent.intec.halvade.utils.Logger.DEBUG("parsing dictionary file...");
+        be.ugent.intec.halvade.utils.Logger.DEBUG("parsing dictionary " + ref + DICT_SUFFIX);
         try {
             FileSystem fs = FileSystem.get(new URI(ref + DICT_SUFFIX), conf);
             FSDataInputStream stream = fs.open(new Path(ref + DICT_SUFFIX));

halvade/src/be/ugent/intec/halvade/HalvadeResourceManager.java

@@ -1,7 +1,18 @@
 /*
- * To change this license header, choose License Headers in Project Properties.
- * To change this template file, choose Tools | Templates
- * and open the template in the editor.
+ * Copyright (C) 2014 ddecap
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program.  If not, see <http://www.gnu.org/licenses/>.
  */
 
 package be.ugent.intec.halvade;
@@ -24,13 +35,13 @@ public class HalvadeResourceManager {
     protected static final int ALL = Integer.MAX_VALUE;
     protected static final int MEM_AM = (int) (1.5*1024);
     protected static final int VCORES_AM = 1;
-    protected static final int MEM_STAR = (int) (12*1024);
+    protected static final int MEM_STAR = (int) (31*1024);
     protected static final int[][] RESOURCE_REQ = { 
-        //mapmem, mapcpu, redmem, redcpu
+        //mapmem, redmem
         {MEM_STAR,  ALL},     // RNA with shared memory pass1
         {MEM_STAR,  14*1024}, // RNA with shared memory pass2
         {MEM_STAR,  14*1024}, // RNA without shared memory
-        {14*1024,   14*1024}, // DNA
+        {16*1024,   16*1024}, // DNA
         {4*1024,    4*1024}   // combine
     };
 
@@ -65,7 +76,7 @@ public static void setJobResources(HalvadeOptions opt, Configuration conf, int t
         conf.set("mapreduce.map.memory.mb", "" + mmem); 
         conf.set("mapreduce.reduce.java.opts", "-Xmx" + (int)(0.8*rmem) + "m");
         conf.set("mapreduce.map.java.opts", "-Xmx" + (int)(0.8*mmem) + "m"); 
-        conf.set("mapreduce.job.reduce.slowstart.completedmaps", "" + 1.0);
+        conf.set("mapreduce.job.reduce.slowstart.completedmaps", "0.99");
 
         HalvadeConf.setMapThreads(conf, opt.mthreads);
         HalvadeConf.setReducerThreads(conf, opt.rthreads);     

halvade/src/be/ugent/intec/halvade/MapReduceRunner.java

@@ -20,9 +20,9 @@
 import fi.tkk.ics.hadoop.bam.SAMRecordWritable;
 import fi.tkk.ics.hadoop.bam.VariantContextWritable;
 import be.ugent.intec.halvade.hadoop.datatypes.ChromosomeRegion;
+import be.ugent.intec.halvade.hadoop.datatypes.GenomeSJ;
 import be.ugent.intec.halvade.hadoop.mapreduce.HalvadeTextInputFormat;
 import be.ugent.intec.halvade.hadoop.partitioners.*;
-import be.ugent.intec.halvade.utils.HalvadeFileUtils;
 import org.apache.hadoop.conf.Configuration;
 import org.apache.hadoop.conf.Configured;
 import org.apache.hadoop.fs.FileStatus;
@@ -133,9 +133,11 @@ protected int runPass1RNAJob(Configuration pass1Conf, String tmpOutDir) throws I
             pass1Job.setMapperClass(be.ugent.intec.halvade.hadoop.mapreduce.StarAlignPassXMapper.class);
 
             pass1Job.setInputFormatClass(HalvadeTextInputFormat.class);
-            pass1Job.setMapOutputKeyClass(LongWritable.class);
+            pass1Job.setMapOutputKeyClass(GenomeSJ.class);
             pass1Job.setMapOutputValueClass(Text.class);
 
+            pass1Job.setSortComparatorClass(GenomeSJSortComparator.class);
+            pass1Job.setGroupingComparatorClass(GenomeSJGroupingComparator.class);
             pass1Job.setNumReduceTasks(1); 
             pass1Job.setReducerClass(be.ugent.intec.halvade.hadoop.mapreduce.RebuildStarGenomeReducer.class);          
             pass1Job.setOutputKeyClass(LongWritable.class);
@@ -194,8 +196,8 @@ protected int runHalvadeJob(Configuration halvadeConf, String tmpOutDir, int job
         halvadeJob.setMapOutputValueClass(SAMRecordWritable.class);
         halvadeJob.setInputFormatClass(HalvadeTextInputFormat.class);
         halvadeJob.setPartitionerClass(ChrRgPartitioner.class);
-        halvadeJob.setSortComparatorClass(ChrRgPositionComparator.class);
-        halvadeJob.setGroupingComparatorClass(ChrRgRegionComparator.class);
+        halvadeJob.setSortComparatorClass(ChrRgSortComparator.class);
+        halvadeJob.setGroupingComparatorClass(ChrRgGroupingComparator.class);
         halvadeJob.setOutputKeyClass(Text.class);
         halvadeJob.setOutputValueClass(VariantContextWritable.class);
 

halvade/src/be/ugent/intec/halvade/hadoop/datatypes/GenomeSJ.java

@@ -0,0 +1,103 @@
+/*
+ * Copyright (C) 2014 ddecap
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program.  If not, see <http://www.gnu.org/licenses/>.
+ */
+
+package be.ugent.intec.halvade.hadoop.datatypes;
+
+import java.io.DataInput;
+import java.io.DataOutput;
+import java.io.IOException;
+import net.sf.samtools.SAMSequenceDictionary;
+import org.apache.hadoop.io.WritableComparable;
+
+/**
+ *
+ * @author ddecap
+ */
+public class GenomeSJ implements WritableComparable<GenomeSJ> {
+    protected int type; // 0 = sj string, 1 = overhang length
+    protected int chromosome;
+    protected int firstBaseIntron;
+    protected int lastBaseIntron;
+    protected int overhang;
+
+    public void setOverhang(int overhang) {
+        this.type = 1;
+        this.overhang = overhang;
+    }
+    
+    public void parseSJString(String sjString, SAMSequenceDictionary dict) {
+        String columns[] = sjString.split("\t");
+        this.type = 0;
+        this.chromosome = dict.getSequenceIndex(columns[0]);
+        this.firstBaseIntron = Integer.parseInt(columns[1]);
+        this.lastBaseIntron = Integer.parseInt(columns[2]);
+    }
+    
+    public int getType() {
+        return type;
+    }
+    public int getOverhang() {
+        return overhang;
+    }
+    public int getChromosome() {
+        return chromosome;
+    }
+    public int getFirstBaseIntron() {
+        return firstBaseIntron;
+    }
+    public int getLastBaseIntron() {
+        return lastBaseIntron;
+    }
+    
+    public GenomeSJ() {
+        this.type = 0;
+        this.chromosome = Integer.MAX_VALUE;
+        this.firstBaseIntron = -1;
+        this.lastBaseIntron = -1;
+        this.overhang = -1;
+    }
+    
+    @Override
+    public void write(DataOutput d) throws IOException {
+        d.writeInt(type);
+        d.writeInt(chromosome);
+        d.writeInt(firstBaseIntron);
+        d.writeInt(lastBaseIntron);
+        d.writeInt(overhang);
+    }
+
+    @Override
+    public void readFields(DataInput di) throws IOException {
+        type = di.readInt();
+        chromosome = di.readInt();
+        firstBaseIntron = di.readInt();
+        lastBaseIntron = di.readInt();
+        overhang = di.readInt();
+    }
+
+    @Override
+    public int compareTo(GenomeSJ o) {
+        if(chromosome == o.chromosome) {
+            if(firstBaseIntron == o.firstBaseIntron)
+                return lastBaseIntron - o.lastBaseIntron;
+            else
+                return firstBaseIntron - o.firstBaseIntron;
+        } else 
+            return chromosome - o.chromosome;
+    }
+    
+}

halvade/src/be/ugent/intec/halvade/hadoop/mapreduce/DnaGATKReducer.java

@@ -37,12 +37,14 @@ protected void processAlignments(Iterable<SAMRecordWritable> values, Context con
         ChromosomeRange r = new ChromosomeRange();
         SAMRecordIterator SAMit = new SAMRecordIterator(values.iterator(), header, r);
 
-        if(useElPrep)
+        if(useElPrep && isFirstAttempt)
             elPrepPreprocess(context, tools, SAMit, preprocess);
-        else 
-            PicardPreprocess(context, tools, SAMit, preprocess); 
+        else  {
+            if(!isFirstAttempt) Logger.DEBUG("attempt " + taskId + ", preprocessing with Picard for smaller peak memory");
+            PicardPreprocess(context, tools, SAMit, preprocess);
+        }
         region = makeRegionFile(context, r, tools, region);
-        if(region == null) return;        
+        if(region == null) return;
 
         indelRealignment(context, region, gatk, preprocess, tmpFile1);        
         baseQualityScoreRecalibration(context, region, r, tools, gatk, tmpFile1, tmpFile2);        

halvade/src/be/ugent/intec/halvade/hadoop/mapreduce/GATKReducer.java

@@ -40,7 +40,7 @@
  * @author ddecap
  */
 public abstract class GATKReducer extends HalvadeReducer {
-
+    protected boolean isFirstAttempt;
     protected boolean useBedTools;
     protected boolean useUnifiedGenotyper;
     protected double sec, scc;
@@ -68,13 +68,14 @@ protected void reduce(ChromosomeRegion key, Iterable<SAMRecordWritable> values,
             processAlignments(values, context, tools, gatk);
         } catch (URISyntaxException | QualityException | ProcessException ex) {
             Logger.EXCEPTION(ex);
-            throw new InterruptedException();
+            throw new InterruptedException(ex.getMessage());
         }
     }
 
     @Override
     protected void setup(Context context) throws IOException, InterruptedException {
         super.setup(context);
+        isFirstAttempt = taskId.endsWith("_0");
         isRNA = HalvadeConf.getIsRNA(context.getConfiguration());
         scc = HalvadeConf.getSCC(context.getConfiguration(), isRNA);
         sec = HalvadeConf.getSEC(context.getConfiguration(), isRNA);
@@ -130,17 +131,20 @@ protected void PicardPreprocess(Context context, PreprocessingTools tools, SAMRe
 
         long startTime = System.currentTimeMillis();
 
+        int count = 0;
         SAMRecord sam;
         while(input.hasNext()) {
             sam = input.next();
             writer.addAlignment(sam);
+            count++;
         }
         int reads = input.getCount();
         writer.close();
 
         context.getCounter(HalvadeCounters.IN_PREP_READS).increment(reads);
         long estimatedTime = System.currentTimeMillis() - startTime;
         context.getCounter(HalvadeCounters.TIME_HADOOP_SAMTOBAM).increment(estimatedTime);
+        Logger.DEBUG("time writing " + count + " records to disk: " + estimatedTime / 1000);
 
         //preprocess steps of iprep
         Logger.DEBUG("clean sam");
@@ -212,7 +216,7 @@ protected void baseQualityScoreRecalibration(Context context, String region, Chr
         String[] snpslocal = HalvadeFileUtils.downloadSites(context, taskId);
         String[] newKnownSites = new String[snpslocal.length];
         for(int i = 0 ; i < snpslocal.length; i++) {
-            if(useBedTools) newKnownSites[i] = tools.filterDBSnps(snpslocal[i], r); 
+            if(useBedTools) newKnownSites[i] = tools.filterDBSnps(ref.replaceAll("fasta", "dict"), snpslocal[i], r, tmpFileBase, threads); 
             else newKnownSites[i] = snpslocal[i]; 
             if(newKnownSites[i].endsWith(".gz"))
                 newKnownSites[i] = HalvadeFileUtils.Unzip(newKnownSites[i]);
@@ -282,11 +286,11 @@ protected void splitNTrim(Context context, String region, GATKTools gatk, String
     }
 
     // TODO improve annotate/filter
-    protected void filterVariants(Context context, GATKTools gatk, String input, String output) throws InterruptedException {      
+    protected void filterVariants(Context context, String region, GATKTools gatk, String input, String output) throws InterruptedException {      
         Logger.DEBUG("run VariantFiltration");
         context.setStatus("run VariantFiltration");
         context.getCounter(HalvadeCounters.TOOLS_GATK).increment(1);
-        gatk.runVariantFiltration(input, output, ref, windows, cluster, minFS, maxQD);
+        gatk.runVariantFiltration(input, output, ref, region, windows, cluster, minFS, maxQD);
 
         HalvadeFileUtils.removeLocalFile(keep, input, context, HalvadeCounters.FOUT_GATK_TMP);
     }
@@ -295,7 +299,7 @@ protected void annotateVariants(Context context, String region, GATKTools gatk,
         Logger.DEBUG("run VariantAnnotator");
         context.setStatus("run VariantAnnotator");
         context.getCounter(HalvadeCounters.TOOLS_GATK).increment(1);
-        gatk.runVariantAnnotator(input, output, ref);
+        gatk.runVariantAnnotator(input, output, ref, region);
 
         HalvadeFileUtils.removeLocalFile(keep, input, context, HalvadeCounters.FOUT_GATK_TMP);
 

halvade/src/be/ugent/intec/halvade/hadoop/mapreduce/HalvadeMapper.java

@@ -11,11 +11,8 @@
 import be.ugent.intec.halvade.utils.HalvadeConf;
 import java.io.File;
 import java.io.IOException;
-import java.io.RandomAccessFile;
 import java.net.URI;
 import java.net.URISyntaxException;
-import java.nio.channels.FileChannel;
-import java.nio.channels.FileLock;
 import java.nio.file.FileSystems;
 import java.nio.file.Files;
 import org.apache.hadoop.fs.Path;
@@ -34,10 +31,10 @@ protected void cleanup(Context context) throws IOException, InterruptedException
         super.cleanup(context);
         try {
             Logger.DEBUG(readcount + " fastq reads processed");
-            allTasksHaveStarted = HalvadeConf.allTasksCompleted(context.getConfiguration());
             Logger.DEBUG("starting cleanup: closing aligner");
             instance.closeAligner();
             Logger.DEBUG("finished cleanup");
+            allTasksHaveStarted = HalvadeConf.allTasksCompleted(context.getConfiguration());
         } catch (URISyntaxException ex) {
             Logger.EXCEPTION(ex);
         }