@split( loadExpressionData, "*.expression_genes.tsv" )
def summarizeExpressionPerGene( infile, outfile ):

    intable = "expression_data"

    dbhandle = sqlite3.connect( PARAMS["database"] )
    columns = Database.getColumnNames( dbhandle, intable ) 
    strains = list(set([ "_".join( x.split("_")[:-1] ) for x in columns if x != "transcript_id" ]))

    total = collections.defaultdict( list )

    for strain in strains:
        E.info( "summarizing strain %s" % strain )
        per_gene = collections.defaultdict( list )
        matched_columns = [ x for x in columns if re.match( strain, x ) ]
        
        fields = ",".join(matched_columns)
        statement = '''SELECT i.gene_id, %(fields)s 
                    FROM %(intable)s as a, transcript_info AS i 
                    WHERE i.transcript_id = a.transcript_id
                    ''' % locals()
        cc = dbhandle.cursor()

        for data in cc.execute( statement ):
            per_gene[data[0]].extend( data[1:] )
            total[data[0]].extend( data[1:] )
            
        outf = open( "%s.expression_genes.tsv" % strain, "w")
        outf.write("gene_id\t%s\n" % Stats.Summary().getHeader() )
        
        for gene_id,d in per_gene.iteritems():
            outf.write( "%s\t%s\n" % (gene_id, str( Stats.Summary( d ) ) ))
            
        outf.close()

    outf = open("all.expression_genes.tsv", "w")
    outf.write("gene_id\t%s\n" % Stats.Summary().getHeader() )
    
    for gene_id,d in total.iteritems():
        outf.write( "%s\t%s\n" % (gene_id, str( Stats.Summary( d ) ) ))
    outf.close()
summarize expression data per gene by combining all transcripts
dataset/ETHPy150Open CGATOxford/cgat/obsolete/pipeline_snps.py/summarizeExpressionPerGene