ipyrad/assemble/cluster_within.py
def declone_3rad(data, sample):
    """
    3rad uses random adapters to identify pcr duplicates. We will
    remove pcr dupes here. Basically append the radom adapter to
    each sequence, do a regular old vsearch derep, then trim
    off the adapter, and push it down the pipeline. This will
    remove all identical seqs with identical random i5 adapters.
    """

    LOGGER.info("Entering declone_3rad - {}".format(sample.name))

    ## Append i5 adapter to the head of each read. Merged file is input, and
    ## still has fq qual score so also have to append several qscores for the
    ## adapter bases. Open the merge file, get quarts, go through each read
    ## and append the necessary stuff.

    adapter_seqs_file = tempfile.NamedTemporaryFile(mode='wb',
                                        delete=False,
                                        dir=data.dirs.edits,
                                        suffix="_append_adapters_.fastq")

    try:
        with open(sample.files.edits[0][0]) as infile:
            quarts = itertools.izip(*[iter(infile)]*4)

            ## a list to store until writing
            writing = []
            counts = 0

            while 1:
                try:
                    read = quarts.next()
                except StopIteration:
                    break

                ## Split on +, get [1], split on "_" (can be either _r1 or
                ## _m1 if merged reads) and get [0] for the i5
                ## prepend "EEEEEEEE" as qscore for the adapters
                i5 = read[0].split("+")[1].split("_")[0]

                ## If any non ACGT in the i5 then drop this sequence
                if 'N' in i5:
                    continue
                writing.append("\n".join([
                                read[0].strip(),
                                i5 + read[1].strip(),
                                read[2].strip(),
                                "E"*8 + read[3].strip()]
                            ))

                ## Write the data in chunks
                counts += 1
                if not counts % 1000:
                    adapter_seqs_file.write("\n".join(writing)+"\n")
                    writing = []
            if writing:
                adapter_seqs_file.write("\n".join(writing))
                adapter_seqs_file.close()

        tmp_outfile = tempfile.NamedTemporaryFile(mode='wb',
                                        delete=False,
                                        dir=data.dirs.edits,
                                        suffix="_decloned_w_adapters_.fastq")

        ## Close the tmp file bcz vsearch will write to it by name, then
        ## we will want to reopen it to read from it.
        tmp_outfile.close()
        ## Derep the data (adapters+seq)
        derep_and_sort(data, adapter_seqs_file.name,
                       os.path.join(data.dirs.edits, tmp_outfile.name), 2)

        ## Remove adapters from head of sequence and write out
        ## tmp_outfile is now the input file for the next step
        ## first vsearch derep discards the qscore so we iterate
        ## by pairs
        with open(tmp_outfile.name) as infile:
            with open(os.path.join(data.dirs.edits, sample.name+"_declone.fastq"),\
                                'wb') as outfile:
                duo = itertools.izip(*[iter(infile)]*2)

                ## a list to store until writing
                writing = []
                counts2 = 0

                while 1:
                    try:
                        read = duo.next()
                    except StopIteration:
                        break

                    ## Peel off the adapters. There's probably a faster
                    ## way of doing this.
                    writing.append("\n".join([
                                    read[0].strip(),
                                    read[1].strip()[8:]]
                                ))

                    ## Write the data in chunks
                    counts2 += 1
                    if not counts2 % 1000:
                        outfile.write("\n".join(writing)+"\n")
                        writing = []
                if writing:
                    outfile.write("\n".join(writing))
                    outfile.close()

        LOGGER.info("Removed pcr duplicates from {} - {}".format(sample.name, counts-counts2))

    except Exception as inst:
        raise IPyradError("    Caught error while decloning "\
                                + "3rad data - {}".format(inst))

    finally:
        ## failed samples will cause tmp file removal to raise.
        ## just ignore it.
        try:
            ## Clean up temp files
            if os.path.exists(adapter_seqs_file.name):
                os.remove(adapter_seqs_file.name)
            if os.path.exists(tmp_outfile.name):
                os.remove(tmp_outfile.name)
        except Exception as inst:
            pass
def declone_3rad(data, sample):
    """
    3rad uses random adapters to identify pcr duplicates. We will
    remove pcr dupes here. Basically append the radom adapter to
    each sequence, do a regular old vsearch derep, then trim
    off the adapter, and push it down the pipeline. This will
    remove all identical seqs with identical random i5 adapters.
    """

    LOGGER.info("Entering declone_3rad - {}".format(sample.name))

    ## Append i5 adapter to the head of each read. Merged file is input, and
    ## still has fq qual score so also have to append several qscores for the
    ## adapter bases. Open the merge file, get quarts, go through each read
    ## and append the necessary stuff.

    adapter_seqs_file = tempfile.NamedTemporaryFile(mode='wb',
                                        delete=False,
                                        dir=data.dirs.edits,
                                        suffix="_append_adapters_.fastq")

    try:
        with open(sample.files.edits[0][0]) as infile:
            quarts = itertools.izip(*[iter(infile)]*4)

            ## a list to store until writing
            writing = []
            counts = 0

            while 1:
                try:
                    read = quarts.next()
                except StopIteration:
                    break

                ## Split on +, get [1], split on "_" (can be either _r1 or
                ## _m1 if merged reads) and get [0] for the i5
                ## prepend "EEEEEEEE" as qscore for the adapters
                i5 = read[0].split("+")[1].split("_")[0]

                ## If any non ACGT in the i5 then drop this sequence
                if 'N' in i5:
                    continue
                writing.append("\n".join([
                                read[0].strip(),
                                i5 + read[1].strip(),
                                read[2].strip(),
                                "E"*8 + read[3].strip()]
                            ))

                ## Write the data in chunks
                counts += 1
                if not counts % 1000:
                    adapter_seqs_file.write("\n".join(writing)+"\n")
                    writing = []
            if writing:
                adapter_seqs_file.write("\n".join(writing))
                adapter_seqs_file.close()

        tmp_outfile = tempfile.NamedTemporaryFile(mode='wb',
                                        delete=False,
                                        dir=data.dirs.edits,
                                        suffix="_decloned_w_adapters_.fastq")

        ## Close the tmp file bcz vsearch will write to it by name, then
        ## we will want to reopen it to read from it.
        tmp_outfile.close()
        ## Derep the data (adapters+seq)
        derep_and_sort(data, adapter_seqs_file.name,
                       os.path.join(data.dirs.edits, tmp_outfile.name), 2)

        ## Remove adapters from head of sequence and write out
        ## tmp_outfile is now the input file for the next step
        ## first vsearch derep discards the qscore so we iterate
        ## by pairs
        with open(tmp_outfile.name) as infile:
            with open(os.path.join(data.dirs.edits, sample.name+"_declone.fastq"),\
                                'wb') as outfile:
                duo = itertools.izip(*[iter(infile)]*2)

                ## a list to store until writing
                writing = []
                counts2 = 0

                while 1:
                    try:
                        read = duo.next()
                    except StopIteration:
                        break

                    ## Peel off the adapters. There's probably a faster
                    ## way of doing this.
                    writing.append("\n".join([
                                    read[0].strip(),
                                    read[1].strip()[8:]]
                                ))

                    ## Write the data in chunks
                    counts2 += 1
                    if not counts2 % 1000:
                        outfile.write("\n".join(writing)+"\n")
                        writing = []
                if writing:
                    outfile.write("\n".join(writing))
                    outfile.close()

        LOGGER.info("Removed pcr duplicates from {} - {}".format(sample.name, counts-counts2))

    except Exception as inst:
        raise IPyradError("    Caught error while decloning "\
                                + "3rad data - {}".format(inst))

    finally:
        ## failed samples will cause tmp file removal to raise.
        ## just ignore it.
        try:
            ## Clean up temp files
            if os.path.exists(adapter_seqs_file.name):
                os.remove(adapter_seqs_file.name)
            if os.path.exists(tmp_outfile.name):
                os.remove(tmp_outfile.name)
        except Exception as inst:
            pass
def declone_3rad ( data , sample ) : LOGGER . info ( "Entering declone_3rad - {}" . format ( sample . name ) ) ## Append i5 adapter to the head of each read. Merged file is input, and ## still has fq qual score so also have to append several qscores for the ## adapter bases. Open the merge file, get quarts, go through each read ## and append the necessary stuff. adapter_seqs_file = tempfile . NamedTemporaryFile ( mode = 'wb' , delete = False , dir = data . dirs . edits , suffix = "_append_adapters_.fastq" ) try : with open ( sample . files . edits [ 0 ] [ 0 ] ) as infile : quarts = itertools . izip ( * [ iter ( infile ) ] * 4 ) ## a list to store until writing writing = [ ] counts = 0 while 1 : try : read = quarts . next ( ) except StopIteration : break ## Split on +, get [1], split on "_" (can be either _r1 or ## _m1 if merged reads) and get [0] for the i5 ## prepend "EEEEEEEE" as qscore for the adapters i5 = read [ 0 ] . split ( "+" ) [ 1 ] . split ( "_" ) [ 0 ] ## If any non ACGT in the i5 then drop this sequence if 'N' in i5 : continue writing . append ( "\n" . join ( [ read [ 0 ] . strip ( ) , i5 + read [ 1 ] . strip ( ) , read [ 2 ] . strip ( ) , "E" * 8 + read [ 3 ] . strip ( ) ] ) ) ## Write the data in chunks counts += 1 if not counts % 1000 : adapter_seqs_file . write ( "\n" . join ( writing ) + "\n" ) writing = [ ] if writing : adapter_seqs_file . write ( "\n" . join ( writing ) ) adapter_seqs_file . close ( ) tmp_outfile = tempfile . NamedTemporaryFile ( mode = 'wb' , delete = False , dir = data . dirs . edits , suffix = "_decloned_w_adapters_.fastq" ) ## Close the tmp file bcz vsearch will write to it by name, then ## we will want to reopen it to read from it. tmp_outfile . close ( ) ## Derep the data (adapters+seq) derep_and_sort ( data , adapter_seqs_file . name , os . path . join ( data . dirs . edits , tmp_outfile . name ) , 2 ) ## Remove adapters from head of sequence and write out ## tmp_outfile is now the input file for the next step ## first vsearch derep discards the qscore so we iterate ## by pairs with open ( tmp_outfile . name ) as infile : with open ( os . path . join ( data . dirs . edits , sample . name + "_declone.fastq" ) , 'wb' ) as outfile : duo = itertools . izip ( * [ iter ( infile ) ] * 2 ) ## a list to store until writing writing = [ ] counts2 = 0 while 1 : try : read = duo . next ( ) except StopIteration : break ## Peel off the adapters. There's probably a faster ## way of doing this. writing . append ( "\n" . join ( [ read [ 0 ] . strip ( ) , read [ 1 ] . strip ( ) [ 8 : ] ] ) ) ## Write the data in chunks counts2 += 1 if not counts2 % 1000 : outfile . write ( "\n" . join ( writing ) + "\n" ) writing = [ ] if writing : outfile . write ( "\n" . join ( writing ) ) outfile . close ( ) LOGGER . info ( "Removed pcr duplicates from {} - {}" . format ( sample . name , counts - counts2 ) ) except Exception as inst : raise IPyradError ( "    Caught error while decloning " + "3rad data - {}" . format ( inst ) ) finally : ## failed samples will cause tmp file removal to raise. ## just ignore it. try : ## Clean up temp files if os . path . exists ( adapter_seqs_file . name ) : os . remove ( adapter_seqs_file . name ) if os . path . exists ( tmp_outfile . name ) : os . remove ( tmp_outfile . name ) except Exception as inst : pass
3rad uses random adapters to identify pcr duplicates. We will
    remove pcr dupes here. Basically append the radom adapter to
    each sequence, do a regular old vsearch derep, then trim
    off the adapter, and push it down the pipeline. This will
    remove all identical seqs with identical random i5 adapters.
3rad uses random adapters to identify pcr duplicates . We will remove pcr dupes here . Basically append the radom adapter to each sequence do a regular old vsearch derep then trim off the adapter and push it down the pipeline . This will remove all identical seqs with identical random i5 adapters .