def identical(args):
    """
    %prog identical *.fasta

    Given multiple fasta files, find all the exactly identical records
    based on the computed md5 hexdigest or GCG checksum of each sequence.

    Output is an N + 1 column file (where N = number of input fasta files).
    If there are duplicates within a given fasta file, they will all be
    listed out in the same row separated by a comma.

    Example output:
    ---------------------------
	       tta1.fsa    tta2.fsa
	t0         2131          na
	t1         3420          na
	t2    3836,3847         852
	t3          148         890
	t4          584         614
	t5          623         684
	t6         1281         470
	t7         3367          na
    """
    from jcvi.utils.cbook import AutoVivification

    allowed_checksum = ["MD5", "GCG"]

    p = OptionParser(identical.__doc__)
    p.add_option("--ignore_case", default=False, action="store_true",
            help="ignore case when comparing sequences [default: %default]")
    p.add_option("--ignore_N", default=False, action="store_true",
            help="ignore N and X's when comparing sequences [default: %default]")
    p.add_option("--ignore_stop", default=False, action="store_true",
            help="ignore stop codon when comparing sequences [default: %default]")
    p.add_option("--output_uniq", default=False, action="store_true",
            help="output uniq sequences in FASTA format" + \
                 " [default: %default]")
    p.add_option("--checksum", default="MD5", choices=allowed_checksum,
            help="specify checksum method [default: %default]")
    p.set_outfile()

    opts, args = p.parse_args(args)

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

    d = AutoVivification()
    files = []
    for fastafile in args:
        f = Fasta(fastafile)
        pf = fastafile.rsplit(".", 1)[0]
        files.append(pf)

        logging.debug("Hashing individual elements of {0}".format(fastafile))
        for name, rec in f.iteritems_ordered():
            seq = re.sub(' ', '', str(rec.seq))
            hashed = hash_fasta(seq, ignore_case=opts.ignore_case, ignore_N=opts.ignore_N, \
                ignore_stop=opts.ignore_stop, checksum=opts.checksum)
            if not d[hashed]:
                d[hashed]['seq'] = seq
                d[hashed]['count'] = 0
            if not d[hashed]['names'][pf]:
                d[hashed]['names'][pf] = set()
            d[hashed]['names'][pf].add(name)

    fw = must_open(opts.outfile, "w")
    if opts.output_uniq:
        uniqfile = "_".join(files) + ".uniq.fasta"
        uniqfw = must_open(uniqfile, "w")

    header = "\t".join(str(x) for x in (args))
    print >> fw, "\t".join(str(x) for x in ("", header))
    for idx, hashed in enumerate(d.keys()):
        line = []
        line.append("t{0}".format(idx))
        for fastafile in files:
            if fastafile in d[hashed]['names'].keys():
                line.append(",".join(d[hashed]['names'][fastafile]))
                if opts.output_uniq:
                    d[hashed]['count'] += len(d[hashed]['names'][fastafile])
            else:
                line.append("na")
        print >> fw, "\t".join(line)

        if opts.output_uniq:
            seqid = "\t".join(str(x) for x in ("t{0}".format(idx), d[hashed]['count']))
            rec = SeqRecord(Seq(d[hashed]['seq']), id=seqid, description="")
            SeqIO.write([rec], uniqfw, "fasta")

    fw.close()
    if opts.output_uniq:
        logging.debug("Uniq sequences written to `{0}`".format(uniqfile))
        uniqfw.close()
dataset/ETHPy150Open tanghaibao/jcvi/formats/fasta.py/identical/original