Create whole genome contigs we want to process, only non-alts.

    Skips problem contigs like HLAs for downstream analysis.
973044
def _create_genome_regions(data):
    
    work_dir = utils.safe_makedir(os.path.join(dd.get_work_dir(data), "coverage", dd.get_sample_name(data)))
    variant_regions = os.path.join(work_dir, "target-genome.bed")
    with file_transaction(data, variant_regions) as tx_variant_regions:
        with open(tx_variant_regions, "w") as out_handle:
            for c in shared.get_noalt_contigs(data):
                out_handle.write("%s\t%s\t%s\n" % (c.name, 0, c.size))
    return variant_regions
1151992
def _prep_genome(out_dir, data):
    
    genome_name = utils.splitext_plus(os.path.basename(dd.get_ref_file(data)))[0]
    out_dir = utils.safe_makedir(os.path.join(out_dir, genome_name))
    ref_file = dd.get_ref_file(data)
    utils.symlink_plus(ref_file, os.path.join(out_dir, os.path.basename(ref_file)))
    with open(os.path.join(out_dir, "GenomeSize.xml"), "w") as out_handle:
        out_handle.write('<sequenceSizes genomeName="%s">' % genome_name)
        for c in pysam.AlignmentFile("%s.dict" % utils.splitext_plus(ref_file)[0]).header["SQ"]:
            cur_ploidy = ploidy.get_ploidy([data], region=[c["SN"]])
            out_handle.write('<chromosome fileName="%s" contigName="%s" totalBases="%s" knownBases="%s" '
                             'isCircular="false" ploidy="%s" md5="%s"/>' %
                             (os.path.basename(ref_file), c["SN"], c["LN"], c["LN"], cur_ploidy, c["M5"]))
        out_handle.write('</sequenceSizes>')
    return out_dir
966991
def main(args):

  
  sys.stderr.write("reading reference genome\n")
  ref = FASTAData(open(args.genome).read())
  sys.stderr.write("read in "+str(len(ref.keys()))+" chromosomes\n")

  """ Next make the transcriptome """
  txome = {}
  sys.stderr.write("write the transcriptome\n")
  inf = None
  if args.gpd[-3:] == '.gz':
    inf = gzip.open(args.gpd)
  else:
    inf = open(args.gpd)
  stream = GPDStream(inf)
  tof = open(args.tempdir+'/transcriptome.fa','w')
  z = 0
  for gpd in stream:
    z += 1
    if gpd.transcript_name in txome:
      sys.stderr.write("WARNING already have a transcript "+gpd.transcript_name+" ignoring line "+str(z)+" of the gpd\n")
      continue
    txome[gpd.transcript_name] = gpd.gene_name
    tof.write('>'+gpd.transcript_name+"\n"+str(gpd.get_sequence(ref))+"\n")
  tof.close()
  inf.close()
  sys.stderr.write("wrote "+str(len(txome.keys()))+" transcripts\n")

  """Build the salmon index"""
  sys.stderr.write("building a salmon index\n")
  cmd = 'salmon index -p '+str(args.numThreads)+' -t '+args.tempdir+'/transcriptome.fa -i '+args.tempdir+'/salmon_index'
  p = Popen(cmd.split())
  p.communicate()
  sys.stderr.write("finished building the index\n")

  """Use the index to quanitfy"""
  sys.stderr.write("quanitfy reads\n")
  reads = ''
  if args.rU:
    reads = '-r '+args.rU
  else:
    reads = '-1 '+args.r1+' -2 '+args.r2
  cmd = 'salmon quant -p '+str(args.numThreads)+' -i '+args.tempdir+'/salmon_index -l A '+reads+' -o '+args.tempdir+'/output_quant'
  p = Popen(cmd.split())
  p.communicate()
  sys.stderr.write("finished quanitfying\n")

  """Now parse the salmon output to add gene name"""
  salmon = {}
  with open(args.tempdir+'/output_quant/quant.sf') as inf:
    header = inf.readline()
    for line in inf:
      f = line.rstrip().split("\t")
      # by each transcript name hold a data strcture of the information
      salmon[f[0]] = {'name':f[0],'length':int(f[1]),'EffectiveLength':float(f[2]),'TPM':float(f[3]),'NumReads':float(f[4])}
  genes = {}
  for name in salmon:
    gene = txome[name]
    if gene not in genes: genes[gene] = []
    genes[gene].append(salmon[name])
  genetot = {}
  for gene in genes:
    tot = sum([x['TPM'] for x in genes[gene]])
    genetot[gene] = tot
  ordered_gene_names = sorted(genetot.keys(), key=lambda x: genetot[x],reverse=True)
      
  """Collected enough information to make output"""
  sys.stderr.write("generating output\n")
  of = sys.stdout
  if args.output != '-':
    of = open(args.output,'w')
  of.write("geneName\ttranscriptName\tlength\tEffectiveLength\ttxTPM\tNumReads\tgeneTPM\n")
  for gene in ordered_gene_names:
    txs = sorted(genes[gene],key=lambda x: x['TPM'],reverse=True)
    for tx in txs:
      of.write(gene+"\t"+tx['name']+"\t"+str(tx['length'])+"\t"+str(tx['EffectiveLength'])+"\t"+str(tx['TPM'])+"\t"+str(tx['NumReads'])+"\t"+str(genetot[gene])+"\n")
  of.close()

  # Temporary working directory step 3 of 3 - Cleanup
  if not args.specific_tempdir:
    rmtree(args.tempdir)
674824
def process_nodes_only(names, options, run_list, create_command_options) # rubocop:disable CyclomaticComplexity
      nodenames = []
      if PROVIDERS.member?(names[0])
        count = names.length == 2 ? names[1] : 1
        do_provider_members(count, nodenames, options)
      elsif names[0].start_with?("windows_")
        nodenames.push(names[1..-1])
      else # standard nodes
        nodenames.push(names)
      end
      nodenames.flatten.each do |node|
        node_names_flatten(create_command_options, node, run_list)
      end
    end
1024609
def mask(args):
    
    p = OptionParser(mask.__doc__)
    p.add_option("--hard", dest="hard", default=False, action="store_true",
            help="Hard mask the low-complexity bases [default: %default]")
    opts, args = p.parse_args(args)

    if len(args) != 1:
        sys.exit(not p.print_help())

    genomefile, = args

    # entire pipeline
    countsfile = genomefile + ".counts"
    wm_mk_counts(infile=genomefile, outfile=countsfile)

    maskedfastafile = "%s.masked%s" % op.splitext(genomefile)
    wm_mk_masks(infile=countsfile, outfile=maskedfastafile, genomefile=genomefile)

    if opts.hard:
        hardmask(maskedfastafile)
1325418
def create_inputs(data):
    
    from bcbio.pipeline import sample
    data = cwlutils.normalize_missing(data)
    aligner = tz.get_in(("config", "algorithm", "aligner"), data)
    # CRAM files must be converted to bgzipped fastq, unless not aligning.
    # Also need to prep and download remote files.
    if not ("files" in data and data["files"] and aligner and (_is_cram_input(data["files"]) or
                                                               objectstore.is_remote(data["files"][0]))):
        # skip indexing on samples without input files or not doing alignment
        if ("files" not in data or not data["files"] or data["files"][0] is None or not aligner):
            return [[data]]
    data["files_orig"] = data["files"]
    data["files"] = prep_fastq_inputs(data["files"], data)
    # preparation converts illumina into sanger format
    data["config"]["algorithm"]["quality_format"] = "standard"
    # Handle any necessary trimming
    data = utils.to_single_data(sample.trim_sample(data)[0])
    _prep_grabix_indexes(data["files"], data)
    data = _set_align_split_size(data)
    out = []
    if tz.get_in(["config", "algorithm", "align_split_size"], data):
        splits = _find_read_splits(data["files"][0], int(data["config"]["algorithm"]["align_split_size"]))
        for split in splits:
            cur_data = copy.deepcopy(data)
            cur_data["align_split"] = split
            out.append([cur_data])
    else:
        out.append([data])
    if "output_cwl_keys" in data:
        out = cwlutils.samples_to_records([utils.to_single_data(x) for x in out],
                                          ["files", "align_split", "config__algorithm__quality_format"])
    return out
748198
def instantiate(args):
    
    p = OptionParser(instantiate.__doc__)
    p.set_annot_reformat_opts()
    p.add_option("--extended_stride", default=False, action="store_true",
                 help="Toggle extended strides for gene numbering")
    opts, args = p.parse_args(args)

    if len(args) != 3:
        sys.exit(not p.print_help())

    taggedbed, blacklist, gapsbed = args
    r = NameRegister(prefix=opts.prefix, pad0=opts.pad0, uc=opts.uc)
    r.get_blacklist(blacklist)
    r.get_gaps(gapsbed)

    # Run through the bed, identify stretch of NEW ids to instantiate,
    # identify the flanking FRAMEs, interpolate!
    bed = Bed(taggedbed)
    outputbed = taggedbed.rsplit(".", 1)[0] + ".new.bed"
    fw = open(outputbed, "w")

    tagkey = lambda x: x.rsplit("|", 1)[-1]
    for chr, sbed in bed.sub_beds():
        current_chr = chr_number(chr)
        if not current_chr:
            continue

        sbed = list(sbed)

        ranks = []
        for i, s in enumerate(sbed):
            nametag = s.extra[0]
            tag = tagkey(nametag)

            if tag in (NEW, FRAME):
                ranks.append((i, nametag))

        blocks = []
        for tag, names in groupby(ranks, key=lambda x: tagkey(x[-1])):
            names = list(names)
            if tag == NEW:
                blocks.append((tag, [sbed[x[0]] for x in names]))
            else:
                start, end = names[0][-1], names[-1][-1]
                start, end = atg_name(start, retval="rank"), atg_name(end, retval="rank")
                blocks.append((tag, [start, end]))

        id_table = {}  # old to new name conversion
        for i, (tag, info) in enumerate(blocks):
            if tag != NEW:
                continue

            start_id = 0 if i == 0 else blocks[i - 1][1][-1]
            end_id = start_id + 10000 if i == len(blocks) -1 \
                        else blocks[i + 1][1][0]

            r.allocate(info, chr, start_id, end_id, id_table, extended_stride=opts.extended_stride)

        # Output new names
        for i, s in enumerate(sbed):
            nametag = s.extra[0]
            name, tag = nametag.split("|")

            if tag == NEW:
                assert name == '.'
                name = id_table[s.accn]
            elif tag == OVERLAP:
                if name in id_table:
                    name = id_table[name]

            s.extra[0] = "|".join((name, tag))
            print(s, file=fw)

    fw.close()
1198219
def clean_sequences(self):
        """"""
        logging.info('Removing plasmids and masking phages')
        plasmid_db = os.path.join(self.reffilepath, 'plasmidfinder', 'plasmid_database.fa')
        phage_db = os.path.join(self.reffilepath, 'prophages', 'combinedtargets.tfa')
        with progressbar(self.runmetadata.samples) as bar:
            for sample in bar:
                plasmid_removal = 'bbduk.sh ref={} in={} out={} overwrite'\
                    .format(plasmid_db, sample.general.combined, sample.general.combined.replace('.f', '_noplasmid.f'))
                subprocess.call(plasmid_removal, shell=True, stdout=self.devnull, stderr=self.devnull)
                phage_masking = 'bbduk.sh ref={} in={} out={} kmask=N overwrite'\
                    .format(phage_db, sample.general.combined.replace('.f', '_noplasmid.f'),
                            sample.general.combined.replace('.f', '_clean.f'))
                subprocess.call(phage_masking, shell=True, stdout=self.devnull, stderr=self.devnull)
                os.remove(sample.general.combined)
                os.rename(sample.general.combined.replace('.f', '_clean.f'), sample.general.combined)
                os.remove(sample.general.combined.replace('.f', '_noplasmid.f'))
1008397
def generate_combined_fasta(self, genome_list, genome_dir):
        '''
        '''

        fasta = []
            
        for genome in genome_list:
            
            full_path = genome_dir + genome
            handle = open(full_path, "rU")
            print 'making combined fasta for', genome
            try:
                seq_record = SeqIO.read(handle, 'genbank')
                org_accession = seq_record.name
                
            except AssertionError,e:
                print str(e), genome
                        
            for feature in seq_record.features:
                if feature.type == 'CDS':
                    try:                                  
                        prot_accession = feature.qualifiers['protein_id'][0]
                        prot_translation = feature.qualifiers['translation'][0]
                        
                        newfast = '>' + org_accession + ',' + prot_accession + \
                         '\n' + prot_translation + '\n'
                        
                        #if newfast not in fasta:
                        fasta.append(newfast)
                            
                    except AttributeError,e:
                        print "organism %s, protein %s did not have \
                        the right attributes" % (org_accession, prot_accession)
                        print str(e)
                    
                    except KeyError,e:
                        print "organism %s, protein %s did not have \
                        the right key" % (org_accession, prot_accession)
                        print str(e)
                          
            handle.close()
            
            print "%s proteins were added" % len(fasta)

            set_fasta = set(fasta)

            print "%s unique proteins were added -- dropping redundant ones" % len(set_fasta)
                                                
        faastring = "".join(set_fasta)

        write_fasta = open('combined_fasta', 'w')
        write_fasta.write(faastring)
        write_fasta.close()

        return set_fasta
1032106
def main():
    
    # Define Parser object and add to toil
    url_prefix = 'https://s3-us-west-2.amazonaws.com/cgl-pipeline-inputs/'
    parser = argparse.ArgumentParser(description=main.__doc__, formatter_class=argparse.RawTextHelpFormatter)
    parser.add_argument('--config', required=True,
                        help='Path to configuration file for samples, one per line. UUID,URL_to_bamfile. '
                             'The URL may be a standard "http://", a "file://<abs_path>", or "s3://<bucket>/<key>"')
    parser.add_argument('--gtf', help='URL to annotation GTF file',
                        default=url_prefix + 'rnaseq_cgl/gencode.v23.annotation.gtf')
    parser.add_argument('--gtf-pickle', help='Pickled GTF file',
                        default=url_prefix + 'spladder/gencode.v23.annotation.gtf.pickle')
    parser.add_argument('--gtf-m53', help='M53 preprocessing annotation table',
                        default=url_prefix + 'spladder/gencode.v23.annotation.gtf.m53')
    parser.add_argument('--positions', help='URL to SNP positions over genes file (TSV)',
                        default=url_prefix + 'spladder/positions_fixed.tsv')
    parser.add_argument('--genome', help='URL to Genome fasta',
                        default=url_prefix + 'rnaseq_cgl/hg38_no_alt.fa')
    parser.add_argument('--genome-index', help='Index file (fai) of genome',
                        default=url_prefix + 'spladder/hg38_no_alt.fa.fai')
    parser.add_argument('--ssec', default=None, help='Path to master key used for downloading encrypted files.')
    parser.add_argument('--output-s3-dir', default=None, help='S3 Directory of the form: s3://bucket/directory')
    parser.add_argument('--output-dir', default=None, help='full path where final results will be output')
    parser.add_argument('--sudo', action='store_true', default=False,
                        help='Set flag if sudo is required to run Docker.')
    parser.add_argument('--star-index', help='URL to download STAR Index built from HG38/gencodev23 annotation.',
                        default=url_prefix + 'rnaseq_cgl/starIndex_hg38_no_alt.tar.gz')
    parser.add_argument('--fwd-3pr-adapter', help="Sequence for the FWD 3' Read Adapter.", default='AGATCGGAAGAG')
    parser.add_argument('--rev-3pr-adapter', help="Sequence for the REV 3' Read Adapter.", default='AGATCGGAAGAG')
    Job.Runner.addToilOptions(parser)
    args = parser.parse_args()
    # Sanity Checks
    if args.config:
        assert os.path.isfile(args.config), 'Config not found at: {}'.format(args.config)
    if args.ssec:
        assert os.path.isfile(args.ssec), 'Encryption key not found at: {}'.format(args.config)
    if args.output_s3_dir:
        assert args.output_s3_dir.startswith('s3://'), 'Wrong format for output s3 directory'
    # Program checks
    for program in ['curl', 'docker']:
        assert which(program), 'Program "{}" must be installed on every node.'.format(program)

    Job.Runner.startToil(Job.wrapJobFn(parse_input_samples, args), args)
671876
def from_interbase_coordinates(contig, start, end=None):
        '''
        
        '''
        typechecks.require_string(contig)
        typechecks.require_integer(start)
        if end is None:
            end = start + 1
        typechecks.require_integer(end)
        contig = pyensembl.locus.normalize_chromosome(contig)
        return Locus(contig, start, end)
25431
def get_noalt_contigs(data):
    
    alts = set([])
    alt_files = [f for f in tz.get_in(["reference", "bwa", "indexes"], data, []) if f.endswith("alt")]
    if alt_files:
        for alt_file in alt_files:
            with open(alt_file) as in_handle:
                for line in in_handle:
                    if not line.startswith("@"):
                        alts.add(line.split()[0].strip())
    else:
        for contig in ref.file_contigs(dd.get_ref_file(data)):
            if ("_alt" in contig.name or "_decoy" in contig.name or
                  contig.name.startswith("HLA-") or ":" in contig.name):
                alts.add(contig.name)
    return [c for c in ref.file_contigs(dd.get_ref_file(data)) if c.name not in alts]
590422
def split(*items):
    
    out = []
    for data in [x[0] for x in items]:
        dis_orgs = data["config"]["algorithm"].get("disambiguate")
        if dis_orgs:
            if not data.get("disambiguate", None):
                data["disambiguate"] = {"genome_build": data["genome_build"],
                                        "base": True}
            out.append([data])
            # handle the instance where a single organism is disambiguated
            if isinstance(dis_orgs, six.string_types):
                dis_orgs = [dis_orgs]
            for dis_org in dis_orgs:
                dis_data = copy.deepcopy(data)
                dis_data["disambiguate"] = {"genome_build": dis_org}
                dis_data["genome_build"] = dis_org
                dis_data["config"]["algorithm"]["effects"] = False
                dis_data = run_info.add_reference_resources(dis_data)
                out.append([dis_data])
        else:
            out.append([data])
    return out
4794
def skesa_assemble(self):
        
        with progressbar(self.metadata) as bar:
            for sample in bar:
                # Initialise the assembly command
                sample.commands.assemble = str()
                try:
                    if sample.general.trimmedcorrectedfastqfiles:
                        # If the sample is a pure isolate, assemble it. Otherwise, run the pre-metagenome pipeline
                        try:
                            status = sample.run.Description
                        except AttributeError:
                            status = 'unknown'
                        if status == 'metagenome':
                            self.merge(sample)
                        else:
                            # Set the output directory
                            sample.general.assembly_output = os.path.join(sample.general.outputdirectory,
                                                                          'assembly_output')
                            make_path(sample.general.assembly_output)
                            sample.general.assemblyfile = os.path.join(sample.general.assembly_output,
                                                                       '{name}_unfiltered.fasta'
                                                                       .format(name=sample.name))
                            sample.general.bestassemblyfile = os.path.join(sample.general.assembly_output,
                                                                           '{name}.fasta'
                                                                           .format(name=sample.name))
                            fastqfiles = sample.general.trimmedcorrectedfastqfiles

                            # Set the the forward fastq files
                            sample.general.assemblyfastq = fastqfiles
                            forward = fastqfiles[0]
                            gz = True if '.gz' in forward else False
                            # If there are two fastq files
                            if len(fastqfiles) == 2:
                                # Set the reverse fastq name https://github.com/ncbi/SKESA/issues/7
                                sample.commands.assemble = 'skesa --fastq {fastqfiles} --cores {threads} ' \
                                                           '--use_paired_ends --vector_percent 1 ' \
                                                           '--contigs_out {contigs}'\
                                    .format(fastqfiles=','.join(fastqfiles),
                                            threads=self.cpus,
                                            contigs=sample.general.assemblyfile)
                            # Same as above, but use single read settings for the assembler
                            else:
                                sample.commands.assemble = 'skesa --fastq {fastqfiles} --cores {threads} ' \
                                                           '--vector_percent 1 --contigs_out {contigs}'\
                                    .format(fastqfiles=','.join(fastqfiles),
                                            threads=self.cpus,
                                            contigs=sample.general.assemblyfile)
                    # If there are no fastq files, populate the metadata appropriately
                    else:
                        sample.general.assembly_output = 'NA'
                        sample.general.assemblyfastq = 'NA'
                        sample.general.bestassemblyfile = 'NA'
                except AttributeError:
                    sample.general.assembly_output = 'NA'
                    sample.general.assemblyfastq = 'NA'
                    sample.general.trimmedcorrectedfastqfiles = 'NA'
                    sample.general.bestassemblyfile = 'NA'
                if sample.commands.assemble and not os.path.isfile(sample.general.assemblyfile):
                    # Run the assembly
                    out, err = run_subprocess(sample.commands.assemble)
                    write_to_logfile(sample.commands.assemble,
                                     sample.commands.assemble,
                                     self.logfile,
                                     sample.general.logout,
                                     sample.general.logerr,
                                     None,
                                     None)
                    write_to_logfile(out,
                                     err,
                                     self.logfile,
                                     sample.general.logout,
                                     sample.general.logerr,
                                     None,
                                     None)
387012
def _create_binary(self, binary_tgt, results_dir):
    """"""
    # Note that we rebuild a chroot from scratch, instead of using the REQUIREMENTS_PEX
    # and PYTHON_SOURCES products, because those products are already-built pexes, and there's
    # no easy way to merge them into a single pex file (for example, they each have a __main__.py,
    # metadata, and so on, which the merging code would have to handle specially).
    interpreter = self.context.products.get_data(PythonInterpreter)
    with temporary_dir() as tmpdir:
      # Create the pex_info for the binary.
      run_info_dict = self.context.run_tracker.run_info.get_as_dict()
      build_properties = PexInfo.make_build_properties()
      build_properties.update(run_info_dict)
      pex_info = binary_tgt.pexinfo.copy()
      pex_info.build_properties = build_properties

      pex_builder = PexBuilderWrapper.Factory.create(
        builder=PEXBuilder(path=tmpdir, interpreter=interpreter, pex_info=pex_info, copy=True),
        log=self.context.log)

      if binary_tgt.shebang:
        self.context.log.info('Found Python binary target {} with customized shebang, using it: {}'
          .format(binary_tgt.name, binary_tgt.shebang))
        pex_builder.set_shebang(binary_tgt.shebang)
      else:
        self.context.log.debug('No customized shebang found for {}'.format(binary_tgt.name))

      # Find which targets provide sources and which specify requirements.
      source_tgts = []
      req_tgts = []
      constraint_tgts = []
      for tgt in binary_tgt.closure(exclude_scopes=Scopes.COMPILE):
        if has_python_sources(tgt) or has_resources(tgt):
          source_tgts.append(tgt)
        elif has_python_requirements(tgt):
          req_tgts.append(tgt)
        if is_python_target(tgt):
          constraint_tgts.append(tgt)

      # Add interpreter compatibility constraints to pex info. This will first check the targets for any
      # constraints, and if they do not have any will resort to the global constraints.
      pex_builder.add_interpreter_constraints_from(constraint_tgts)

      # Dump everything into the builder's chroot.
      for tgt in source_tgts:
        pex_builder.add_sources_from(tgt)

      # We need to ensure that we are resolving for only the current platform if we are
      # including local python dist targets that have native extensions.
      self._python_native_code_settings.check_build_for_current_platform_only(self.context.targets())
      pex_builder.add_requirement_libs_from(req_tgts, platforms=binary_tgt.platforms)

      # Build the .pex file.
      pex_path = os.path.join(results_dir, '{}.pex'.format(binary_tgt.name))
      pex_builder.build(pex_path)
      return pex_path
685505
def split_genome(genome, chunk_size=10000):
    

    chunks = []
    from Bio import SeqIO
    with open(genome) as handle:
        for record in SeqIO.parse(handle, "fasta"):
            sequence = record.seq
            n = len(sequence)
            chunks += [str(sequence[i:min(i + chunk_size, n)])
                       for i in range(0, n, chunk_size)]
    return np.array(chunks)
295940
def _handle_special_inputs(inputs, variables):
    
    from bcbio import structural
    optional = [["config", "algorithm", "coverage"],
                ["config", "algorithm", "variant_regions"],
                ["config", "algorithm", "sv_regions"],
                ["config", "algorithm", "validate"],
                ["config", "algorithm", "validate_regions"]]
    all_vs = set([get_base_id(v["id"]) for v in variables])
    out = []
    for input in inputs:
        if input == ["reference", "aligner", "indexes"]:
            for v in variables:
                vid = get_base_id(v["id"]).split("__")
                if vid[0] == "reference" and vid[1] in alignment.TOOLS:
                    out.append(vid)
        elif input == ["reference", "snpeff", "genome_build"]:
            found_indexes = False
            for v in variables:
                vid = get_base_id(v["id"]).split("__")
                if vid[0] == "reference" and vid[1] == "snpeff":
                    out.append(vid)
                    found_indexes = True
            assert found_indexes, "Found no snpEff indexes in %s" % [v["id"] for v in variables]
        elif input == ["config", "algorithm", "background", "cnv_reference"]:
            for v in variables:
                vid = get_base_id(v["id"]).split("__")
                if (vid[:4] == ["config", "algorithm", "background", "cnv_reference"] and
                      structural.supports_cnv_reference(vid[4])):
                    out.append(vid)
        elif input in optional:
            if _get_string_vid(input) in all_vs:
                out.append(input)
        else:
            out.append(input)
    return out