@transform(buildExpectedGenomeCoverage, regex("(\S+)/(\S+).coverage"), add_inputs(buildCoverageOverGenomes), r"genome_coverage.dir/\2_observed.tsv")
def mergeExpectedAndObservedGenomeCoverage(infiles, outfile):

    expected = open(infiles[0])
    expected_header = expected.readline()
    observed = open(infiles[1])
    observed_header = observed.readline()

    expected_data = {}
    E.info("reading expected coverage over genomes")
    for line in expected.readlines():
        data = line[:-1].split("\t")
        gi, coverage = data[0], data[1]
        expected_data[gi] = coverage

    outf = open(outfile, "w")
    E.info("writing results")
    outf.write("track\tgi\tspecies\tobserved\texpected\n")
    for line in observed.readlines():
        data = line[:-1].split("\t")
        track, gi, species, coverage = data[0], data[
            1], "_".join(data[2].split("_")[5:7]), data[3]
        outf.write("%s\t%s\t%s\t%s\t%s\n" %
                   (track, gi, species, coverage, expected_data[gi]))
    outf.close()

###################################################
###################################################
###################################################
merge the expected and actual estimates
dataset/ETHPy150Open CGATOxford/cgat/obsolete/pipeline_metagenomebenchmark.py/mergeExpectedAndObservedGenomeCoverage