From 77caa6a460d9ac6d104d9c85c973bae74064e812 Mon Sep 17 00:00:00 2001
From: hanna julienne <hanna.julienne@pasteur.fr>
Date: Tue, 29 Oct 2019 14:46:28 +0100
Subject: [PATCH] corrected sumZ and added pleiotropy index

---
 jass/models/stats.py     |  1 -
 jass/models/worktable.py | 14 +++++++++++++-
 2 files changed, 13 insertions(+), 2 deletions(-)

diff --git a/jass/models/stats.py b/jass/models/stats.py
index d64a581c..c25bb57c 100755
--- a/jass/models/stats.py
+++ b/jass/models/stats.py
@@ -108,7 +108,6 @@ def omnibus_stat(z, cov, invcov):
 omnibus_stat.can_use_pattern = True
 
 
-
 def fisher_test(z, cov, invcov):
     """
     fisher combination of univariate p-values corrected with bonferoni
diff --git a/jass/models/worktable.py b/jass/models/worktable.py
index e1b5fb78..3f1feac6 100755
--- a/jass/models/worktable.py
+++ b/jass/models/worktable.py
@@ -122,6 +122,15 @@ def post_computation_filtering(worktable_chunk, significant_treshold = 5*10**-8)
 
     return worktable_chunk
 
+def compute_pleiotropy_index(W,significance_treshold):
+
+    N_significatif = (2.0 * spst.norm.sf(W.fillna(0, inplace=False).abs())<significance_treshold).sum(1)
+    N_pheno = (~W.isnull()).sum(1)
+    #pleiotropy index is not meaningful for too few phenotype
+    S = N_significatif/N_pheno
+    S.loc[N_pheno < 4] = np.nan
+    return S
+
 def create_worktable_file(
     phenotype_ids: List[str],
     init_file_path: str,
@@ -255,8 +264,11 @@ def create_worktable_file(
         sum_stat_tab["UNIVARIATE_MIN_QVAL"] =  sum_stat_tab["UNIVARIATE_MIN_PVAL"]*(1-np.isnan(sum_stat_tab[phenotype_ids]).astype(int)).sum(1)
         sum_stat_tab.loc[sum_stat_tab.UNIVARIATE_MIN_QVAL>1 , "UNIVARIATE_MIN_QVAL"] = 1
 
+        #Computing pleiotropy
+        sum_stat_tab["PLEIOTROPY_INDEX"] = compute_pleiotropy_index(sum_stat_tab[phenotype_ids], significance_treshold)
+
         sum_stat_tab = sum_stat_tab[
-            ["Region", "CHR", "snp_ids", "position", 'Ref_allele', 'Alt_allele',  "MiddlePosition", "JASS_PVAL", "UNIVARIATE_MIN_PVAL", "UNIVARIATE_MIN_QVAL"]
+            ["Region", "CHR", "snp_ids", "position", 'Ref_allele', 'Alt_allele',  "MiddlePosition", "JASS_PVAL", "UNIVARIATE_MIN_PVAL", "UNIVARIATE_MIN_QVAL","PLEIOTROPY_INDEX" ]
             + phenotype_ids]
 
         if post_filtering:
-- 
GitLab