# ---
# jupyter:
#   jupytext:
#     text_representation:
#       extension: .py
#       format_name: light
#       format_version: '1.5'
#       jupytext_version: 1.14.4
#   kernelspec:
#     display_name: Python 2
#     language: python
#     name: python2
# ---

# +
# Set thread count for network propagation functions
import mkl
mkl.set_num_threads(4)

# Import other required packages
import os
import pandas as pd
# -

# # Propagating GWAS Hits from WTCCC

# Load functions for performing network GWAS from python module file
# %run './network_GWAS_functions.py'

# Load network for propagation
network = load_network('./Data/PCNet.txt', delimiter='\t')

# Load example GWAS-based gene-disease associations from file as pandas Series
wd = './Data/'
diseases, GWAS_Associations = [], {}
for fn in os.listdir(wd):
    if fn.endswith('GWAS_Associations.csv'):
        disease_name = fn.split('_')[0]
        diseases.append(disease_name)
        GWAS_Association_table = pd.read_csv(wd+fn, index_col=0, header=-1)
        GWAS_Associations[disease_name] = GWAS_Association_table[1]
        GWAS_Associations[disease_name].name = 'GWAS P'

p_thresh = '1e-5'

# Perform random walk propagation for each disease
save_dir = './Results/'
for disease in diseases:
    propagation_results_table = random_walk(network, GWAS_Associations[disease], p_thresh=1e-5)
    propagation_results_table.to_csv(save_dir+disease+'_p'+p_thresh+'_prop_results.csv')
    print

# # Evaluation of Propagation Results vs Disease Gold Standard

# Load gold standard gene sets
f = open('./Data/WTCCC_DisGeNET_Gold_Standards.txt')
lines = f.read().splitlines()
gold_std_disease_geneset = {}
for line in lines:
    line_split = line.split('\t')
    gold_std_disease_geneset[line_split[0]] = line_split[1:]

# Load functions for evaluating network GWAS results
# %run './network_GWAS_results_analysis_functions.py'

# +
# Load propagation results and evaluate it against disease gold standard
GWAS_AUPRCs, prop_AUPRCs, adjp_AUPRCs = {}, {}, {}

for disease in diseases:
    # Load propagation results table
    results_table = pd.read_csv(save_dir+disease+'_p'+p_thresh+'_prop_results.csv', index_col=0)
    # Perform precision/recall calculations
    GWAS_AUPRC, prop_AUPRC, adjp_AUPRC = PRC_plots(results_table, gold_std_disease_geneset[disease], p_thresh=p_thresh, 
                                                   outdir=save_dir, plot_title=disease+' (Propagation Threshold: p <= '+p_thresh+')',
                                                   file_prefix = disease+'_p'+p_thresh)
    # Record AUPRC results
    GWAS_AUPRCs[disease] = GWAS_AUPRC
    prop_AUPRCs[disease] = prop_AUPRC
    adjp_AUPRCs[disease] = adjp_AUPRC
# -

weak_p_thresh='1e-4'

# Load propagation results and visualize propagation result distributions
for disease in diseases:
    # Load propagation results table
    results_table = pd.read_csv(save_dir+disease+'_p'+p_thresh+'_prop_results.csv', index_col=0)
    # Perform precision/recall calculations
    prop_score_distributions(results_table, gold_std_disease_geneset[disease], p_thresh=p_thresh, weak_p_thresh=weak_p_thresh,
                             outdir=save_dir, plot_title=disease+' (Propagation Threshold: p <= '+p_thresh+')',
                             file_prefix = disease+'_p'+p_thresh)

StarcoderdataJupyterScriptsDedupFiltered