@transform(loadEstimatedTaxonomicRelativeAbundances, suffix(".load"), add_inputs(loadTrueTaxonomicAbundances), ".comparison.pdf")
def compareAbundanceOfFalsePositiveSpecies(infiles, outfile):
    tablename_estimate = P.toTable(infiles[0])

    track = P.snip(
        os.path.basename(infiles[0]).replace("metaphlan_", ""), ".load")
    tablename_true = [P.toTable(x) for x in infiles[1:] if P.snip(
        os.path.basename(x), ".load") == track][0]
    dbh = sqlite3.connect("csvdb")
    cc = dbh.cursor()
    tmp = P.getTempFile(".")
    tmp.write("taxa\tabundance\tstatus\n")
    estimate = {}
    true = set()
    for data in cc.execute("""SELECT taxon, rel_abundance FROM %s WHERE taxon_level == 'species'""" % tablename_estimate).fetchall():
        estimate[data[0]] = data[1]
    for data in cc.execute("""SELECT taxa FROM %s WHERE level == 'species'""" % tablename_true).fetchall():
        true.add(data[0])

    for taxa, abundance in estimate.iteritems():
        if taxa in true:
            tmp.write("%s\t%f\ttp\n" % (taxa, abundance))
        else:
            tmp.write("%s\t%f\tfp\n" % (taxa, abundance))
    tmp.close()

    inf = tmp.name
    if track.find("15M") != -1:
        col = "cadetblue"
    elif track.find("30M") != -1:
        col = "lightblue"
    elif track.find("50M") != -1:
        col = "slategray"

    R('''dat <- read.csv("%s", header = T, stringsAsFactors = F, sep = "\t")''' %
      inf)
    R('''library(ggplot2)''')
    R('''ggplot(dat, aes(x = status, y = log2(abundance))) + geom_boxplot(colour = "%s") + geom_hline(yintersect=0, linetype="dashed")''' %
      col)
    R('''ggsave("%s")''' % outfile)
    os.unlink(inf)

###################################################
###################################################
###################################################
boxplot the relative abundance of false positive
dataset/ETHPy150Open CGATOxford/cgat/obsolete/pipeline_metagenomebenchmark.py/compareAbundanceOfFalsePositiveSpecies