Convert MDmodule output to motifs
        
        Parameters
        ----------
        fo : file-like
            File object containing MDmodule output.

        Returns
        -------
        motifs : list
            List of Motif instances.
483196
def parse(self, fo):
        
        motifs = []
        nucs = {"A":0,"C":1,"G":2,"T":3}
        p = re.compile(r'(\d+)\s+(\d+\.\d+)\s+(\d+\.\d+)\s+(\d+\.\d+)\s+(\d+\.\d+)')
        pf = re.compile(r'>.+\s+[bf]\d+\s+(\w+)')

        pwm = []
        pfm = []
        align = []
        m_id = ""
        for line in fo.readlines():
            if line.startswith("Motif"):
                if m_id:
                    motifs.append(Motif())
                    motifs[-1].id = m_id
                    motifs[-1].pwm = pwm
                    motifs[-1].pfm = pfm
                    motifs[-1].align = align
                    pwm = []
                    pfm = []
                    align = []
                m_id = line.split("\t")[0]
            else: 
                m = p.search(line)
                if m:
                    pwm.append([float(m.group(x))/100 for x in [2,3,4,5]])
                m = pf.search(line)
                if m:
                    if not pfm:
                        pfm = [[0 for x in range(4)] for x in range(len(m.group(1)))]
                    for i in range(len(m.group(1))):
                        pfm[i][nucs[m.group(1)[i]]] += 1
                            
                    align.append(m.group(1))
        
        if pwm:
            motifs.append(Motif())
            motifs[-1].id = m_id
            motifs[-1].pwm = pwm
            motifs[-1].pfm = pfm
            motifs[-1].align = align

        return motifs
502816
def _extract_mod_from_feature(mf):
        """"""
        # ModificationFeature / SequenceModificationVocabulary
        mf_type = mf.getModificationType()
        if mf_type is None:
            return None
        mf_type_terms = mf_type.getTerm().toArray()
        known_mf_type = None
        for t in mf_type_terms:
            if t.startswith('MOD_RES '):
                t = t[8:]
            mf_type_indra = _mftype_dict.get(t)
            if mf_type_indra is not None:
                known_mf_type = mf_type_indra
                break
        if not known_mf_type:
            logger.debug('Skipping modification with unknown terms: %s' %
                        ', '.join(mf_type_terms))
            return None

        mod_type, residue = known_mf_type

        # getFeatureLocation returns SequenceLocation, which is the
        # generic parent class of SequenceSite and SequenceInterval.
        # Here we need to cast to SequenceSite in order to get to
        # the sequence position.
        mf_pos = mf.getFeatureLocation()
        if mf_pos is not None:
            # If it is not a SequenceSite we can't handle it
            if not mf_pos.modelInterface.getName() == \
                'org.biopax.paxtools.model.level3.SequenceSite':
                mod_pos = None
            else:
                mf_site = cast(_bp('SequenceSite'), mf_pos)
                mf_pos_status = mf_site.getPositionStatus()
                if mf_pos_status is None:
                    mod_pos = None
                elif mf_pos_status and mf_pos_status.toString() != 'EQUAL':
                    logger.debug('Modification site position is %s' %
                                mf_pos_status.toString())
                else:
                    mod_pos = mf_site.getSequencePosition()
                    mod_pos = '%s' % mod_pos
        else:
            mod_pos = None
        mc = ModCondition(mod_type, residue, mod_pos, True)
        return mc
777719
def location(args):
    
    fastafile = args.fastafile
    pwmfile = args.pwmfile

    lwidth = args.width
    if not lwidth:
        f = Fasta(fastafile)
        lwidth = len(f.items()[0][1])
        f = None

    jobs = []
    motifs = pwmfile_to_motifs(pwmfile)
    ids = [motif.id for motif in motifs]
    if args.ids:
        ids = args.ids.split(",")
    
    n_cpus = int(MotifConfig().get_default_params()["ncpus"])
    pool = Pool(processes=n_cpus, maxtasksperchild=1000) 
    for motif in motifs:
        if motif.id in ids:
            outfile = os.path.join("%s_histogram" % motif.id)
            jobs.append(
                    pool.apply_async(
                        motif_localization, 
                        (fastafile,motif,lwidth,outfile, args.cutoff)
                        ))
    
    for job in jobs:
        job.get()
847742
def generateFeatureGroups(fgiContainer, linkageGroups, matchArr, timeKey,
                          massKey, logMassKey, massScalingFactor):
    
    #Generate feature groups from the linked features
    newFgiIds = list()
    for linkageGroup in viewvalues(linkageGroups):
        fgiId = fgiContainer._getNextFgiId()
        fgi = fiGroupFromLinkageGroup(matchArr, linkageGroup, fgiId,
                                      timeKey, massKey
                                      )
        fgiContainer.container[fgiId] = fgi
        fgi.metrics = clusterMetrics(matchArr[timeKey][linkageGroup],
                                     matchArr[logMassKey][linkageGroup],
                                     massScalingFactor=massScalingFactor
                                     )
        fgi.rt = fgi.metrics['meanTime']
        fgi.mz = fgi.metrics['meanMass']
        newFgiIds.append(fgiId)
    return newFgiIds
749920
def decompile(f):
    
    co = f.__code__
    args, kwonly, varargs, varkwargs = paramnames(co)
    annotations = f.__annotations__ or {}
    defaults = list(f.__defaults__ or ())
    kw_defaults = f.__kwdefaults__ or {}

    if f.__name__ == '<lambda>':
        node = ast.Lambda
        body = pycode_to_body(co, DecompilationContext(in_lambda=True))[0]
        extra_kwargs = {}
    else:
        node = ast.FunctionDef
        body = pycode_to_body(co, DecompilationContext(in_function_block=True))
        extra_kwargs = {
            'decorator_list': [],
            'returns': annotations.get('return')
        }

    return node(
        name=f.__name__,
        args=make_function_arguments(
            args=args,
            kwonly=kwonly,
            varargs=varargs,
            varkwargs=varkwargs,
            defaults=defaults,
            kw_defaults=kw_defaults,
            annotations=annotations,
        ),
        body=body,
        **extra_kwargs
    )
1372783
def create_denovo_motif_report(inputfile, pwmfile, fgfa, background, locfa, outdir, params, stats=None):
    """"""
    logger.info("creating reports")

    motifs = read_motifs(pwmfile, fmt="pwm")
    
    # ROC plots
    create_roc_plots(pwmfile, fgfa, background, outdir)
    
    # Closest match in database
    mc = MotifComparer()
    closest_match = mc.get_closest_match(motifs)
    
    if stats is None:
        stats = {}
        for bg, bgfa in background.items():
            for m, s in calc_stats(motifs, fgfa, bgfa).items():
                if m not in stats:
                    stats[m] = {}
                stats[m][bg] = s

    stats = add_star(stats)

    if not params:
        params = {}
    cutoff_fpr = params.get('cutoff_fpr', 0.9)
    lwidth = np.median([len(seq) for seq in Fasta(locfa).seqs])

    # Location plots
    logger.debug("Creating localization plots")
    for motif in motifs:
        logger.debug("  {} {}".format(motif.id, motif))
        outfile = os.path.join(outdir, "images/{}_histogram.svg".format(motif.id))
        motif_localization(locfa, motif, lwidth, outfile, cutoff=cutoff_fpr)

    # Create reports
    _create_text_report(inputfile, motifs, closest_match, stats, outdir)
    _create_graphical_report(inputfile, pwmfile, background, closest_match, outdir, stats)
495027
def filterMotifs(memeFile,outFile, minSites):
    

    with open(memeFile, "r") as mF:
        oldMEME=mF.readlines()
        newMEME=oldMEME[:7]
        i=7
        while i < len(oldMEME):
            if oldMEME[i].split(" ")[0] == "MOTIF":
                print(oldMEME[i].split("\n")[0], int(oldMEME[i+2].split("nsites= ")[1].split(" ")[0]))
                sys.stdout.flush()
                if int(oldMEME[i+2].split("nsites= ")[1].split(" ")[0]) > minSites:
                    newMEME.append(oldMEME[i])
                    f=i+1
                    while oldMEME[f].split(" ")[0] != "MOTIF":
                        newMEME.append(oldMEME[f])
                        f=f+1
                    i=i+1
                else:
                    i=i+1
            else:
                i=i+1
    with open(outFile, "w+") as out:
        out.write("".join(newMEME) )

    return newMEME
810317
def _plot_weights_motif(self, index, plot_type="motif_raw",
                            background_probs=DEFAULT_BASE_BACKGROUND,
                            ncol=1,
                            figsize=None):
        

        w_all = self.get_weights()

        if len(w_all) == 0:
            raise Exception("Layer needs to be initialized first")
        W = w_all[0]
        if index is None:
            index = np.arange(W.shape[2])

        if isinstance(index, int):
            index = [index]
        fig = plt.figure(figsize=figsize)

        if plot_type == "motif_pwm" and plot_type in self.AVAILABLE_PLOTS:
            arr = pssm_array2pwm_array(W, background_probs)
        elif plot_type == "motif_raw" and plot_type in self.AVAILABLE_PLOTS:
            arr = W
        elif plot_type == "motif_pwm_info" and plot_type in self.AVAILABLE_PLOTS:
            quasi_pwm = pssm_array2pwm_array(W, background_probs)
            arr = _pwm2pwm_info(quasi_pwm)
        else:
            raise ValueError("plot_type needs to be from {0}".format(self.AVAILABLE_PLOTS))

        fig = seqlogo_fig(arr, vocab=self.VOCAB_name, figsize=figsize, ncol=ncol, plot_name="filter: ")
        # fig.show()
        return fig
512234
def parse(self, fo):
        
        motifs = []
        p = re.compile(r'\d+\s+@\s+\d+\.\d+\s+sd\s+\d+\.\d+\s+(\w+)$')

        line = fo.readline()
        while line and line.find("Color") == -1:
            m = p.search(line)
            if m:
                pwm_data = {}
                for i in range(4):
                    vals = [x.strip() for x in fo.readline().strip().split(" ") if x]
                    pwm_data[vals[0].upper()] = vals[1:]
                pwm = []
                for i in range(len(pwm_data["A"])):
                    pwm.append([float(pwm_data[x][i]) for x in ["A","C","G","T"]])
                motifs.append(Motif(pwm))
                motifs[-1].id = "%s_%s" % (self.name, m.group(1))
            line = fo.readline()
        
        return motifs
1130575
def get_structure_from_dict(d):
    

    species = []
    frac_coords = []
    masses = []
    if 'points' in d:
        for p in d['points']:
            species.append(p['symbol'])
            frac_coords.append(p['coordinates'])
            masses.append(p['mass'])
    elif 'atoms' in d:
        for p in d['atoms']:
            species.append(p['symbol'])
            frac_coords.append(p['position'])
            masses.append(p['mass'])
    else:
        raise ValueError('The dict does not contain structural information')

    return Structure(d['lattice'], species, frac_coords,
                     site_properties={"phonopy_masses": masses})
290885
def parse(self, fo):
        
        motifs = []
        m = [[float(x) for x in fo.readline().strip().split(" ")] for i in range(4)]
        matrix = [[m[0][i], m[1][i],m[2][i],m[3][i]] for i in range(len(m[0]))]
        motifs = [Motif(matrix)]
        motifs[-1].id = self.name
        
        return motifs
1213053
def to_python(fname, converter):
    
    with hdf5.File(fname, 'r') as f:
        source_model = f['/']
    for sg in source_model:
        for src in sg:
            if hasattr(src, 'mfd'):
                # multipoint source
                # src.tom = converter.tom
                kwargs = getattr(src.mfd, 'kwargs', {})
                if 'bin_width' not in kwargs:
                    kwargs['bin_width'] = [converter.width_of_mfd_bin]
    return source_model
628418
def parse_out(self, fo):
        
        motifs = []
        nucs = {"A":0,"C":1,"G":2,"T":3}
        pseudo = 0.0 # Should be 1/sqrt(# of seqs)
        aligns = {}
        for line in fo.readlines():
            if line.startswith("#"):
                pass
            elif len(line) > 1:
                vals = line.strip().split("\t")
                m_id, site = [x.strip().split(" ")[1].replace('"',"") for x in vals[8].split(";") if x]
                #if vals[6] == "+":
                if site.upper().find("N") == -1:
                    aligns.setdefault(m_id, []).append(site)
                #else:
                #    print site, rc(site)
                #    aligns.setdefault(id, []).append(rc(site))
                        
        for m_id, align in aligns.items():
            #print id, len(align)

            width = len(align[0])
            pfm =  [[0 for x in range(4)] for x in range(width)]
            for row in align:
                for i in range(len(row)):
                    pfm[i][nucs[row[i]]] += 1
            total = float(len(align))
            pwm = [[(x + pseudo/4)/total+(pseudo) for x in row] for row in pfm]
            m = Motif()
            m.align = align[:]
            m.pwm = pwm[:]
            m.pfm = pfm[:]
            m.id = m_id
            motifs.append(m)    
        return motifs
763273
def filesfile_string(self):
        """"""
        lines = []
        app = lines.append

        app(self.input_file.path)          # 1) Path of the input file
        app(self.output_file.path)         # 2) Path of the output file
        app(self.ddb_filepath)             # 3) Input derivative database e.g. t13.ddb.in
        app(self.md_filepath)              # 4) Output molecular dynamics e.g. t13.md
        app(self.gkk_filepath)             # 5) Input elphon matrix elements  (GKK file)
        app(self.outdir.path_join("out"))  # 6) Base name for elphon output files e.g. t13
        app(self.ddk_filepath)             # 7) File containing ddk filenames for elphon/transport.

        return "\n".join(lines)
1005767
def to_dict(_dict, onts, mdb_safe=False):
    
    res = {}
    if isinstance(_dict, SamlBase):
        res["__class__"] = "%s&%s" % (_dict.c_namespace, _dict.c_tag)
        for key in _dict.keyswv():
            if key in IMP_SKIP:
                continue
            val = getattr(_dict, key)
            if key == "extension_elements":
                _eel = extension_elements_to_elements(val, onts)
                _val = [_eval(_v, onts, mdb_safe) for _v in _eel]
            elif key == "extension_attributes":
                if mdb_safe:
                    _val = dict([(k.replace(".", "__"), v) for k, v in
                                 val.items()])
                    #_val = {k.replace(".", "__"): v for k, v in val.items()}
                else:
                    _val = val
            else:
                _val = _eval(val, onts, mdb_safe)

            if _val:
                if mdb_safe:
                    key = key.replace(".", "__")
                res[key] = _val
    else:
        for key, val in _dict.items():
            _val = _eval(val, onts, mdb_safe)
            if _val:
                if mdb_safe and "." in key:
                    key = key.replace(".", "__")
                res[key] = _val
    return res
237253
def create_roc_plots(pwmfile, fgfa, background, outdir):
    """"""
    motifs = read_motifs(pwmfile, fmt="pwm", as_dict=True)
    ncpus = int(MotifConfig().get_default_params()['ncpus'])
    pool = Pool(processes=ncpus)
    jobs = {}
    for bg,fname in background.items():
        for m_id, m in motifs.items():

            k = "{}_{}".format(str(m), bg)
            jobs[k] = pool.apply_async(
                                            get_roc_values,
                                            (motifs[m_id], fgfa, fname,)
                                            )
    imgdir = os.path.join(outdir, "images")
    if not os.path.exists(imgdir):
        os.mkdir(imgdir)
    
    roc_img_file = os.path.join(outdir, "images", "{}_roc.{}.png")

    for motif in motifs.values():
        for bg in background:
            k = "{}_{}".format(str(motif), bg)
            error, x, y = jobs[k].get()
            if error:
                logger.error("Error in thread: %s", error)
                logger.error("Motif: %s", motif)
                sys.exit(1)
            roc_plot(roc_img_file.format(motif.id, bg), x, y)
1013855
def get_sequences(options, fasta_file_contents):
    ''' 
    '''
    errors = []
    fasta_files_str = ", ".join(fasta_file_contents.keys())
    fasta_records = None
    reverse_mapping = {}

    try:
        fasta_records, reverse_mapping = parse_FASTA_files(options, fasta_file_contents)
        if not fasta_records:
            errors.append("No protein sequences found in the FASTA file(s) %s." % fasta_files_str)
    except Exception, e:
        e = '\n'.join([l for l in traceback.format_exc(), str('e') if l.strip()])
        errors.append("Error parsing FASTA file(s) %s:\n%s" % (fasta_files_str, str(e)))

    if not fasta_records:
        return None, {}, errors

    colorprinter.message('Found %d protein sequence(s).' % len(fasta_records))
    return fasta_records, reverse_mapping, errors