New upstream version 1.1.2

  - conda install --yes nomkl

  - conda install --yes -c conda-forge -c bioconda bcbio-nextgen

  - conda install --yes -c conda-forge -c bioconda bcbio-nextgen-vm

  - conda install --yes -c conda-forge -c bioconda cwltool toil rabix-bunny

  # Clean up space with external tools we don't need for tests

  - conda clean --yes --tarballs --index-cache

  - conda remove --yes --force qt

  # Update to latest bcbio-nextgen code within the container

  - bcbio_vm.py devel setup_install -i quay.io/bcbio/bcbio-vc

  # -- Standard bcbio variant tests

  - docker run -v `pwd`:`pwd` quay.io/bcbio/bcbio-vc bash -c "cd `pwd` && /usr/local/share/bcbio-nextgen/anaconda/bin/py.test -p no:cacheprovider tests/unit --cov=bcbio"

  - py.test -p no:cacheprovider tests/bcbio_vm -v -m docker_multicore

  - docker run -v `pwd`:`pwd` quay.io/bcbio/bcbio-vc bash -c "cd `pwd` && /usr/local/share/bcbio-nextgen/anaconda/bin/py.test -p no:cacheprovider -p no:stepwise tests/unit --cov=bcbio"

  - py.test -p no:cacheprovider -p no:stepwise tests/bcbio_vm -v -m docker_multicore

  # -- bcbio variant CWL tests

  - py.test -p no:cacheprovider tests/bcbio_vm -v -s -m cwl_docker_joint

  - py.test -p no:cacheprovider tests/bcbio_vm -v -s -m cwl_docker_somatic

  # XXX Currently commented out joint test, taking too long and causing Travis timeouts

  # - py.test -p no:cacheprovider -p no:stepwise tests/bcbio_vm -v -s -m cwl_docker_joint

  - py.test -p no:cacheprovider -p no:stepwise tests/bcbio_vm -v -s -m cwl_docker_somatic

  # -- platform integration

  - sudo mkdir -p /etc/pki/tls/certs && sudo ln -s /etc/ssl/certs/ca-certificates.crt /etc/pki/tls/certs/ca-bundle.crt

  - py.test -p no:cacheprovider tests/bcbio_vm -v -s -m cwl_arvados

  # XXX Arvados tests failing with 404 when contacting from Travis

  #- py.test -p no:cacheprovider -p no:stepwise tests/bcbio_vm -v -s -m cwl_arvados

  # -- Cleanup variant docker image

  - docker ps -a -q | xargs --no-run-if-empty docker rm

  - docker rmi -f quay.io/bcbio/bcbio-vc

  - docker images

  - df -h

  - bcbio_vm.py devel setup_install -i quay.io/bcbio/bcbio-rnaseq

  - py.test -p no:cacheprovider tests/bcbio_vm -v -s -m cwl_docker_rnaseq

  - py.test -p no:cacheprovider -p no:stepwise tests/bcbio_vm -v -s -m cwl_docker_rnaseq

  # -- Cleanup RNA-seq docker image

  - docker ps -a -q | xargs --no-run-if-empty docker rm

  - docker rmi -f quay.io/bcbio/bcbio-rnaseq

## 1.1.2 (12 December 2018)

- VarDict low frequency somatic filters: generalize strand and mismatch based

  filter based on cross-validation to avoid over filtering on high depth panels.

- strelka2 joint calling: switch to improved gvcfgenotyper approach for calling

  from gVCFs.

- VarDict RNA-seq variant calling: avoid structural variants with recent vardict-java.

- RNA-seq variation: filter RNA-seq variants close to splice junctions,

  supporting STAR and hisat2.

- RNA-seq variation: add snpEff effects to output variant calls. Thanks to Manasa Surakala.

- RNA-seq: gzip/bgzip FASTQ files in `work/fastq` instead of the original directory.

- use biobambam2 BAM to FASTQ conversion instead of Picard in all cases.

- Trimming: add built-in support for adapters from the SMARTer Universal Low Input RNA Kit

(truseq2) and the Illumina NEXTera DNA prep kit from NEB (nextera2).

- ChIP/ATAC-seq: allow skipping duplicate marking.

- joint calling: ensure correct upload to final directory when no annotations present

- Logging: fix logging in parallel runs with new joblib loky backend. Thanks to

  Ben Liesfeld and Roland Ewald.

## 1.1.1 (6 November 2018)

- single-cell RNA-seq: add built-in support for 10x_v2.

- Fix UMI support for small RNA. Compatible with Qiagen UMI small RNA protocol.

- Ignore .Renviron when running Rscript to head-off PATH conflicts.

- Support SRR ids to download samples with bcbio_prepare_samples script.

- tumor-only prioritization: do not apply LowPriority filter by default, instead

  annotate with external databases. Use `tumoronly_germline_filter` to re-enable

  previous behavior.

- UMIs: apply default filtering based on de-duplicated read depth. Uses

`--min-reads 2` with raw de-duplicated coverage of 800 or more  or `--min-reads 1`

  otherwise. Allows error correction with UMIs for higher depth samples.

- gemini: databases no longer created by default. Use `tools_on: [gemini]` or

  `tools_on: [gemini_orig]` to create a database. We now use a reduced database

  for build 37 to match build 38 and make this forward compatible with CWL.

- vcfanno: run gemini and somatic annotations by default, producing annotated

  VCFs with external information.

- alignment preparation: support a list of split files from multiple sequencing

  lanes, merging into a single fastq

- variant: support octopus variant caller for germline and somatic samples.

- peddy: fix bug where not all files uploaded on first pipeline run

- peddy: For somatic analyses use separate germline calls for tumor/normal, if

  available, or extracted germline calls from supported callers, instead of

  somatic variants.

- GATK: support ploidy specification during joint calling.

- GATK BQSR: bin qualities into static groups (10, 20, 30) to match GATK4

  recommendations. Thanks to Severine Catreux.

- GATK: support 4.0.10.0 which does not use UCSC 2bit references for Spark tools

- variant calling: support bcftools 1.9 which is more strict about duplicated

  key names in INFO and FORMAT.

- seq2c: Upload global calls, coverage and read_mapping files to project directory.

- RNA-seq variant calling: Apply annotations after joint calling for GATK to

  avoid import errors with GenomicsDB. Thanks to Komal Rathi.

- CWL: add `--cwl` target to bcbio_nextgen.py upgrade to add and maintain bcbio-vm.

- CWL: use standard null instead of string "null" for representing None values.

- CWL: support for heterogeneity and structural variant callers that make

  use of variant inputs.

- CWL: support ensemble calling for combining multiple variant callers.

- ensemble: remove no-ALT ref calls that contribute to incorrect ensemble outputs

- RNA-seq: output a matrix of un-deduped UMI counts when doing single-cell/DGE

  for quality control purposes. This is called `tagcounts-dupes.mtx` in the

  final directory.

- single-cell RNA-seq: allow pre-transformed FASTQ files as input to DGE/single-cell

  pipeline.

- single-cell RNA-seq: only create one index per specified genome instead of per

  sample

- fgbio: back compatibility for older quality setting `--min-consensus-base-quality`

- RNA-seq: fix for `fusion_caller` getting interpreted as a path, leading to

  memoization/upload issues.

- RNA-seq: memoize rRNA quality calculations, speeding up reruns.

- RNA-seq: prefix `description` with an X if it starts with a number, for R

  compatibility.

  Thanks to Avinash Reddy and Dan Stetson at AstraZeneca.

- single-cell RNA-seq: respect `--positional` flag with the new tag counting. Thanks to

  Babak Alaei at AstraZeneca.

- RNA-seq: turn on `--seqBias` flag by default for Salmon as early-version overfitting

  issues have been fixed.

- RNA-seq: report insert size from Salmon fragment distribution, not samtools stats.

- RNA-seq: when processing explant samples, produce a combined tx2gene.csv file from

  all organisms processed.

## 1.1.0 (11 July 2018)

- Germline calls: rename outputs to `samplename-germline` to provide easier

  to understand outputs in final directory.

- Add bcbioRNASeq object creation and automatic quality report generation

  with `tools_on: [bcbiornaseq]`

- CWL: Support germline/somatic calling for tumor samples.

- CNVkit: improve whole genome runs. Better speed in normalize_sv_coverage

  through parallelization and avoiding logging. Avoid memory errors in segmentation.

- UMI: upload prepared UMI bam file (pre-consensus) to final output directory

- Add support for bbmap as an aligner

- RNA-seq variant calling: parallelize GATK HaplotypeCaller over regions to

  avoid memory and timeout issues.

- Support joint calling with GATK using pre-prepared gVCF inputs.

- RNA-seq variant calling: allow annotation of output variants with vcfanno

- Support hg38 builds with peddy QC

- QC: support VerifyBamID2 for contamination detection

- CWL: adjust defaults for align_split_size and nomap_split_targets to match

  different parallelization and overhead for these runs

- CWL: support for Cromwell runner

- custom genomes: Unzip GTF file prior to installation.

- Avoid making variant_regions required during processing (by filling with

  coverage) to differentiate targeted and non analyses downstream.

- Avoid attempts to download pre-installed S3 genomes, providing better

  errors with missing genome installs.

- Trimming: add explicit `polyg` option for removing 3' G stretches in NovaSeq

  and NextSeq data. Now defaults to no polyG trimming unless turned on.

- Chip-seq: Add RiP calculation for chip-seq data.

- DeepVariant and Strelka2 support for customized targeted/genome calling models

  per region to handle heterogeneous inputs.

- STAR: enable passing custom options for alignment.

- Add `tools_off: [coverage_qc]` option to skip calculating coverage stats (samtools-stats and picard).

- Adding BAM file for each sample in small-RNAseq pipeline, samtools

  and qualimap qc metrics to multiqc report.

- Allow arbitrary genomes for ChIP-seq. Thanks to @evchambers for pointing out the issue.

## 1.0.9 (10 April 2018)

- Use smoove for lumpy variant calling and genotyping, replacing custom lumpyexpress

include *.txt

include *.md

include bcbio/data/umis/*.txt

include bcbio/data/umis/*.txt.gz

include bcbio/data/umis/*.json

include config/*.yaml

include config/*.ini

    """Determine if a BAM file is empty

    """

    bam_file = objectstore.cl_input(bam_file)

    sambamba = config_utils.get_program("sambamba", {})

    cmd = ("set -o pipefail; "

           "{sambamba} view {bam_file} | head -1 | wc -l")

           "samtools view {bam_file} | head -1 | wc -l")

    p = subprocess.Popen(cmd.format(**locals()), shell=True,

                         executable=do.find_bash(),

                         stdout=subprocess.PIPE, stderr=subprocess.PIPE,

    http://stackoverflow.com/a/12451083/252589

    """

    bam_file = objectstore.cl_input(bam_file)

    sambamba = config_utils.get_program("sambamba", {})

    cmd = ("set -o pipefail; "

           "{sambamba} view {bam_file} | head -50000 | "

           "{sambamba} view -S -F paired /dev/stdin  | head -1 | wc -l")

           "samtools view -h {bam_file} | head -50000 | "

           "samtools view -S -f 1 /dev/stdin  | head -1 | wc -l")

    p = subprocess.Popen(cmd.format(**locals()), shell=True,

                         executable=do.find_bash(),

                         stdout=subprocess.PIPE, stderr=subprocess.PIPE,

def get_downsample_pct(in_bam, target_counts, data):

    """Retrieve percentage of file to downsample to get to target counts.

    Avoids minimal downsample which is not especially useful for

    improving QC times; 90& or more of reads.

    """

    total = sum(x.aligned for x in idxstats(in_bam, data))

    with pysam.Samfile(in_bam, "rb") as work_bam:

        n_rgs = max(1, len(work_bam.header.get("RG", [])))

    rg_target = n_rgs * target_counts

    if total > rg_target:

        return float(rg_target) / float(total)

        pct = float(rg_target) / float(total)

        if pct < 0.9:

            return pct

def get_aligned_reads(in_bam, data):

    index(in_bam, data["config"], check_timestamp=False)

        sort_base = sort_base.split(suffix)[0]

    return sort_base + SUFFIXES[order]

def aligner_from_header(in_bam):

    """Identify aligner from the BAM header; handling pre-aligned inputs.

    """

    from bcbio.pipeline.alignment import TOOLS

    with pysam.Samfile(in_bam, "rb") as bamfile:

        for pg in bamfile.header.get("PG", []):

            for ka in TOOLS.keys():

                if pg.get("PN", "").lower().find(ka) >= 0:

                    return ka

def sample_name(in_bam):

    """Get sample name from BAM file.

    """

        return 0

    return int(numpy.median(lengths))

def filter_stream_cmd(bam_file, data, filter_flag):

    """

    return a command to keep only alignments matching the filter flag

    see https://github.com/lomereiter/sambamba/wiki/%5Bsambamba-view%5D-Filter-expression-syntax for examples

    """

    sambamba = config_utils.get_program("sambamba", data["config"])

    num_cores = dd.get_num_cores(data)

    cmd = ('{sambamba} view -t {num_cores} -f bam -F "{filter_flag}" {bam_file}')

    return cmd.format(**locals())

def filter_primary_stream_cmd(bam_file, data):

    return filter_stream_cmd(bam_file, data, "not secondary_alignment")

def filter_primary(bam_file, data):

    """Filter reads to primary only BAM.

    Removes:

      - not primary alignment (0x100) 256

      - supplementary alignment (0x800) 2048

    """

    stem, ext = os.path.splitext(bam_file)

    out_file = stem + ".primary" + ext

    if utils.file_exists(out_file):

        return out_file

    if not utils.file_exists(out_file):

        with file_transaction(data, out_file) as tx_out_file:

        cmd = filter_primary_stream_cmd(bam_file, data)

        cmd += "> {tx_out_file}"

            cores = dd.get_num_cores(data)

            cmd = ("samtools view -@ {cores} -F 2304 -b {bam_file} > {tx_out_file}")

            do.run(cmd.format(**locals()), ("Filtering primary alignments in %s." %

                                            os.path.basename(bam_file)))

    return out_file

from bcbio.pipeline import datadict as dd

from bcbio.variation import coverage

from bcbio.variation import multi as vmulti

from bcbio.structural import regions

def sample_callable_bed(bam_file, ref_file, data):

        return r.name == "CALLABLE" and (not noalt_calling or chromhacks.is_nonalt(r.chrom))

    out_file = "%s-callable_sample.bed" % os.path.splitext(bam_file)[0]

    with shared.bedtools_tmpdir(data):

        callable_bed, depth_files = coverage.calculate(bam_file, data)

        sv_bed = regions.get_sv_bed(data)

        callable_bed, depth_files = coverage.calculate(bam_file, data, sv_bed)

        input_regions_bed = dd.get_variant_regions(data)

        if not utils.file_uptodate(out_file, callable_bed):

            with file_transaction(data, out_file) as tx_out_file:

    producing a global set of callable regions.

    """

    samples = utils.unpack_worlds(samples)

    samples = [cwlutils.unpack_tarballs(x, x) for x in samples]

    samples = cwlutils.unpack_tarballs(samples, samples[0])

    # back compatibility -- global file for entire sample set

    global_analysis_file = os.path.join(samples[0]["dirs"]["work"], "analysis_blocks.bed")

    if utils.file_exists(global_analysis_file) and not _needs_region_update(global_analysis_file, samples):

        return analysis_file, no_analysis_file

    else:

        return None, None

def get_split_regions(bed_file, data):

    """Retrieve a set of split regions using the input BED for callable regions.

    Provides a less inclusive hook for parallelizing over multiple regions.

    """

    out_file = "%s-analysis_blocks.bed" % utils.splitext_plus(bed_file)[0]

    with shared.bedtools_tmpdir(data):

        if not utils.file_uptodate(out_file, bed_file):

            ref_regions = get_ref_bedtool(dd.get_ref_file(data), data["config"])

            nblock_regions = ref_regions.subtract(pybedtools.BedTool(bed_file)).saveas()

            min_n_size = int(tz.get_in(["config", "algorithm", "nomap_split_size"], data, 250))

            block_filter = NBlockRegionPicker(ref_regions, data["config"], min_n_size)

            final_nblock_regions = nblock_regions.filter(

                block_filter.include_block).saveas().each(block_filter.expand_block).saveas()

            with file_transaction(data, out_file) as tx_out_file:

                final_regions = ref_regions.subtract(final_nblock_regions, nonamecheck=True).\

                                saveas().merge(d=min_n_size).saveas(tx_out_file)

        chroms = set([])

        with shared.bedtools_tmpdir(data):

            for r in pybedtools.BedTool(bed_file):

                chroms.add(r.chrom)

        out = []

        for r in pybedtools.BedTool(out_file):

            if r.chrom in chroms:

                out.append((r.chrom, r.start, r.stop))

        return out