@follows(mkdir("genes.dir"))
@transform(glob.glob(os.path.join(PARAMS.get("assemblydir"), "ray.dir/*.fasta.gz")),
           regex("(\S+)/(\S+).fasta.gz"),
           r"genes.dir/\2.start_stops")
def countCompleteGenes(infile, outfile):
    start = "ATG"
    stop = ["TAG", "TAA", "TGA"]

    ntotal = 0
    nstart = 0
    nstop = 0
    nstart_nstop = 0
    for fasta in FastaIterator.iterate(IOTools.openFile(infile)):
        ntotal += 1
        if fasta.sequence.startswith(start):
            nstart += 1
        if fasta.sequence[-3:len(fasta.sequence)] in stop:
            nstop += 1
        if fasta.sequence.startswith(start) and fasta.sequence[-3:len(fasta.sequence)] in stop:
            nstart_nstop += 1
    outf = open(outfile, "w")
    outf.write("total_genes\tpstart\tpstop\tpstart_stop\n")
    outf.write("\t".join(map(str,[ntotal, float(nstart)/ntotal, float(nstop)/ntotal, float(nstart_nstop)/ntotal])) + "\n")
    outf.close()

###################################################
###################################################
###################################################









# AH: default value for results_resultsdir
count the number of genes that are classed
dataset/ETHPy150Open CGATOxford/cgat/obsolete/pipeline_metagenomebenchmark.py/countCompleteGenes