<commit_before>#!/usr/bin/env python
# -*- coding: UTF-8 -*-

"""
Finishing pipeline, starting with a phase1/2 BAC. The pipeline ideally should
include the following components

+ BLAST against the Illumina contigs to fish out additional seqs
+ Use minimus2 to combine the contigs through overlaps
+ Map the mates to the contigs and perform scaffolding
"""

import os
import os.path as op
import sys
import logging

from collections import defaultdict

from jcvi.formats.contig import ContigFile
from jcvi.formats.fasta import Fasta, SeqIO, gaps, format, tidy
from jcvi.formats.sizes import Sizes
from jcvi.utils.cbook import depends
from jcvi.assembly.base import n50
from jcvi.apps.align import run_megablast
from jcvi.apps.base import OptionParser, ActionDispatcher, sh, mkdir, need_update


def main():

    actions = (
        ('screen', 'screen sequences against library'),
        ('dedup', 'remove duplicate contigs within assembly'),
        ('build', 'build assembly files after a set of clean-ups'),
        ('overlap', 'build larger contig set by fishing additional seqs'),
        ('overlapbatch', 'call overlap on a set of sequences'),
        ('scaffold', 'build scaffolds based on the ordering in the AGP file'),
            )
    p = ActionDispatcher(actions)
    p.dispatch(globals())


def dedup(args):
    """
    %prog dedup assembly.assembly.blast assembly.fasta

    Remove duplicate contigs within assembly.
    """
    from jcvi.formats.blast import BlastLine

    p = OptionParser(dedup.__doc__)
    p.set_align(pctid=0, pctcov=98)
    opts, args = p.parse_args(args)

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

    blastfile, fastafile = args
    cov = opts.pctcov / 100.
    sizes = Sizes(fastafile).mapping
    fp = open(blastfile)
    removed = set()
    for row in fp:
        b = BlastLine(row)
        query, subject = b.query, b.subject
        if query == subject:
            continue
        qsize, ssize = sizes[query], sizes[subject]<commit_msg>add assembly.postprocess.dust()<commit_after>#!/usr/bin/env python
# -*- coding: UTF-8 -*-

"""
Finishing pipeline, starting with a phase1/2 BAC. The pipeline ideally should
include the following components

+ BLAST against the Illumina contigs to fish out additional seqs
+ Use minimus2 to combine the contigs through overlaps
+ Map the mates to the contigs and perform scaffolding
"""

import os
import os.path as op
import sys
import logging

from collections import defaultdict

from jcvi.formats.contig import ContigFile
from jcvi.formats.fasta import Fasta, SeqIO, gaps, format, parse_fasta, tidy
from jcvi.formats.sizes import Sizes
from jcvi.utils.cbook import depends
from jcvi.assembly.base import n50
from jcvi.apps.align import run_megablast
from jcvi.apps.base import OptionParser, ActionDispatcher, sh, mkdir, need_update


def main():

    actions = (
        ('screen', 'screen sequences against library'),
        ('dedup', 'remove duplicate contigs within assembly'),
        ('dust', 'remove low-complexity contigs within assembly'),
        ('build', 'build assembly files after a set of clean-ups'),
        ('overlap', 'build larger contig set by fishing additional seqs'),
        ('overlapbatch', 'call overlap on a set of sequences'),
        ('scaffold', 'build scaffolds based on the ordering in the AGP file'),
            )
    p = ActionDispatcher(actions)
    p.dispatch(globals())


def dust(args):
    """
    %prog dust assembly.fasta

    Remove low-complexity contigs within assembly.
    """
    p = OptionParser(dust.__doc__)
    opts, args = p.parse_args(args)

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

    fastafile, = args
    dustfastafile = fastafile.rsplit(".", 1)[0] + ".dust.fasta"
    if need_update(fastafile, dustfastafile):
        cmd = "dustmasker -in {0}".format(fastafile)
        cmd += " -out {1} -outfmt fasta".format(dustfastafile)
        sh(cmd)

    for name, seq in parse_fasta(dustfastafile):
        nlow = sum(1 for x in seq if x in "acgtN")
        pctlow = nlow * 100. / len(seq)
        if pctlow < 98:
            continue
        #print "{0}\t{1:.1f}".format(name, pctlow)
        print name


def dedup(args):
    """
    %prog dedup assembly.assembly.blast assembly.fasta

    Remove duplicate contigs within assembly.
    """
    from jcvi.formats.blast import BlastLine

    p = OptionParser(dedup.__doc__)
    p.set_align(pctid=0, pctcov=98)
    opts, args = p.parse_args(args)

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

    blastfile, fastafile = args
    cov = opts.pctcov / 100.
    sizes = Sizes(fastafile).mapping
    fp = open(blastfile)
    removed = set()
    for row in fp:
        b = BlastLine(row)
        query, subject = b.query, b.subject
        if query == subject:
            continue
        qsize, ssize = sizes[query], sizes[subject]
{'fraction_non_alphanumeric': 0.10367481010059536, 'fraction_numerical': 0.008211866146581811, 'mean_word_length': 3.033112582781457, 'pattern_counts': {'":': 0, '<': 4, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '4002134', 'n_tokens_mistral': 1566, 'n_tokens_neox': 1482, 'n_words': 488}