How to get enriched pathways in the data using continous statistic measure?

Question

I was doing pathway enrichment analysis using the below code

p_value = scipy.stats.hypergeom.sf(len(intersection) - 1, geneset_size, len(pathway), len(hits))
        if p_value < 0.05:
            for gene in list(pathway):
                results = results.append({"Gene": gene, "Model": model, "Enriched Pathway": Pathway_name}, ignore_index=True)

where,

hits = df["Gene"][np.logical_and(df["Model"] == model, df["p_value_anova"] < 0.05)]
intersection = pathway.intersection(hits)
Model = model name (like A, B, C) in input file
Gene = Gene from the input gene
Enriched pathway=Pathways that are enriched from the Pathway input file provided

Now the problem is that earlier this checking pathway enrichment on the basis of p_value in hits was working but now using cutoff (0.05) is resulting in loss of all pathways. I am looking to test the enrichment using the statistical test that will work on continous measure than binarizing the P-value.

Answer 1

Just to mention I personally think there is a higher loop above the loop stated above.

I've not looked at each bit of code in detail but there's a bug is on line 4:

results = results.append({"Gene": gene, "Model": model, "Enriched Pathway": Pathway_name}, ignore_index=True)

Firstly, this is fashionable coding because you are dumping straight into a pandas dataframe (described below). The correction is simply ... declare pandas in advance

import pandas as pd
results = pd.DataFrame()
p_value = # code missing
if p_value < 0.05:
    for gene in list(pathway):
        results.append({"Gene": gene, "Model": model, "Enriched Pathway": Pathway_name}, ignore_index=True)

The thing is if there is a higher loop above list(pathway) then results = results.append(.... will over-ride the last set of results. This would explain why results 'disappear' because if the last data in the loop is a blank - or very little output, all previous stuff gets deleted.

I personally would use a dictionary here. Thus

# declare the dictionary before all loops     
results = {"Gene": [], "Model": [], "Enriched Pathway": []} 
p_value = # code missing
if p_value < 0.05:
    for gene in list(pathway):
       genedict = {"Gene": gene, "Model": model, "Enriched Pathway": Pathway_name}
       for k,v in mydicti.values():
           results[v].append(v)

... you could then transfer this to pandas once the loops are complete via,

df = pd.DataFrame(genedict)

I do suspect there's other issues, but as a general rule results = results.append(... is better declaring results first and then appending to results, unless you've a very specific reason (e.g. writing each loop into a separate file). The pandas thing is kinda clever for big data sets, because pandas is fast, whilst dictionaries are slow - for big data (like >>500 000 big). The interchangeability between dictionaries and pandas means it's kinda easier conceptually (and code-wise) to switch at the coders convenience in my personal opinion.