#!/usr/bin/env python3

import sys
import os
import glob
import re
from scipy.stats import hypergeom

infilename = sys.argv[1]
popsize = int(sys.argv[2])
outfilename = sys.argv[3]

# open input file and read file names to process
infile = open(infilename, "r")
filelines = infile.readlines()
infile.close()

for line in filelines:
    line = line.strip()
    filelist = line.split("\t")
    # get file names from list
    filename1 = filelist[0]
    filename2 = filelist[1]
    m = re.search('(.*)-miso', filename1)
    samplepath = m.group(1)
    dirlist = samplepath.split("/")
    samplename = dirlist[-1]
    
    # open the files and calculate overlap and p-value
    file1 = open(filename1, "r")
    file2 = open(filename2, "r")
    genelist1 = file1.readlines()
    genelist2 = file2.readlines()
    file1.close()
    file2.close()

    # make a set of all the genes in genelist2 separating comma-delimeted gene multiplets
    geneset2 = set()
    for gene in genelist2:
        gene = gene.strip()
        genemult = gene.split(",")
        for element in genemult:
            geneset2.add(element)

    # loop through genes in genelist1 and make a set of all genes in genelist2
    commonset = set()
    for gene in genelist1:
        gene = gene.strip()
        genemult = gene.split(",")
        for element in genemult:
            if element in geneset2:
                commonset.add(element)
                break

    commongenes = list(commonset)
    numgenes1 = len(genelist1)
    numgenes2 = len(genelist2)
    numcommon = len(commongenes)

    # Calculate p-value of overlap
    pval = hypergeom.sf(numcommon-1, popsize, numgenes1, numgenes2)
    print("%s\t%d\t%d\t%d\t%g" % (samplename, numgenes1, numgenes2, numcommon, pval))

    # Output overlap in overlap file
    outputdir = "/".join(dirlist)
    outputdir = outputdir + "-miso/"
    outfname = outputdir + outfilename
    outfile = open(outfname, "w")
    for gene in commongenes:
        outfile.write("%s\n" %(gene))
    outfile.close()
    
