From b24f3fcc0ff5ee9c5379e90c888816778b219af8 Mon Sep 17 00:00:00 2001
From: HarshitGupta11 <50410275+HarshitGupta11@users.noreply.github.com>
Date: Mon, 22 Jul 2019 18:30:27 +0530
Subject: [PATCH] Prediction FIx

---
 __pycache__/select.cpython-36.pyc | Bin 0 -> 2785 bytes
 prediction_thread_multiple.py     | 169 ++++++++++++++++++++++++++++++
 test_prepare_functions.py         | 128 +++++++++++-----------
 3 files changed, 233 insertions(+), 64 deletions(-)
 create mode 100644 __pycache__/select.cpython-36.pyc
 create mode 100644 prediction_thread_multiple.py

diff --git a/__pycache__/select.cpython-36.pyc b/__pycache__/select.cpython-36.pyc
new file mode 100644
index 0000000000000000000000000000000000000000..ba0586111f0a92fc6f2beb1ec3850ec196ffbb6d
GIT binary patch
literal 2785
zcmbtW%Z}Vc6s;;h`#GJV8JHxHgm*(S<e_ZPXo8e1Sj{FvNY+G_ak+cC(``FtPskv9
zH9@dt3E~G>vSiDSU+@aVs(*k*ScG%S{g??tEO4uCRoyC=tL{1HR^RG$n!o>a`+?JD
z><@P85omvmlHG%mO!6u7xUPj~NFj}<%rm7aE!38@rGwg*u56%oWK*_KyRt1is2h^s
zXM^q>E3y5MQ-#Mwu~6=!WP31UR&o^a<Z8}xUa%8>!tb#M7fYis3iJzOo9(i?E<|DM
zvPWB=F^q_z!1{cDASz1*@hGZ{X%r7k&j^&N_{<Y(R`KzVOxgI%v)kcE2cKq<$__pW
zRCYMd4@OxM9)!s(4ODPVMp1fg8U-pI$D@NRnnq#1lO9)_p^Ab$@@GNntLQP_kLU;)
zBF7jv@!RBW-s8$eKcth<<DgL6K*<zL&I+?I`Z!%_e$9qlaw(p=g}6r^fzvUGb53xo
z6TZz#t1xp>SSS1;Q<L17m=mk+W42vb6DrQcEo_N(zY^Pwet65e&!|X-$e8P(6RfoB
zmuNeEm=Qkspk4BB5Xf!zt$1YKufwEcg@}rs>lJ@8;FTR`DbB9_ijMWHxjP+aImK<S
zOY70%Hv96*D9U{q<bj{6B;DDaCSfqm_Rv_)_IiiOESV;wV?Rrya2#cWMrDlS+%w17
zI93#_O2lbp15TODs#Y4Q*-;+k<0J;SNmyCqSVm7O`*AQm!s|{Y<WRNrds$kU^ySKk
zldLjnCtkN+$k%U0S$#OJU(jc%I%hUbm}Wf~21qv-z@F!(P}}f@bb5Md22n~=_!l%Z
ze28OVJp>e`CTC%5NSvHNbpV@NxDX5@u*wG7hSb27O#lm$7Yzx}Ppq=F=r<SrHo;%C
ziuMT9C_5mARdi1HtE@B&1|oEb2p9pyCXGoGqgp~`7o#l{qDl7*O*Zt38%0O2h~0xM
z&VJQdp_3n{(VXZzC)LjBtfTM0qvzK<^S_=ye?Dj{g1&MB1re!S$EOmhtn4s}!yvD*
zG$rs*o0zGHd1`|UMNgfB8T1rUv9f2FKL)5-8WYn9PW1|nxY&_To6%#zECQ?`Q`I2r
z(}!um9|o#&P{&YjWy@jhK&>4P8v3o~nhaQH=&PF3bh?@hmljh0JeUS?7)j0M46LRE
z0ou}@vWqYbdT;Y4N?Y_qk9X1P!gn54UB`GxRHo-lw)5avQbmfID##1ayv_?E=Myn8
za==PCuw1kXfgG?)AOTE(JFzEDVE|CE)cr#5#=h4WC83{%iHb)0^DqNfbx>)%P60q0
z)F`as8JHRmwM9d8R?k!r&WlmnsB>d2b$oAtnTu>lQK16C+<<o_7(#m`m_N~3M&UmN
z<NSZY)N<hzz0+c8AmiE)oKbc$qkF;^GrCJTt%X*Pli(n+=A^Zf*ap^YqmaZl&Pc4M
zCDtpt5Sc;Z)FN~D)p?v=T_8h5Qm?_Rq_t@Ld`G<w_n=c_c8TUu{!y37Tp^>;dxKhQ
zaZ$ui0jw2vu154t8hVS&N?Yn}>b^tfT{7>HAt=@RWIiBsmCQ95T=7F6P^%sCd`L#C
zW-Xn||0W%1^pH~iO}v!R7XGOnSL=stUuy_=Da{UPb$W=tbh+L|$$o&z@xP?teI~^P
zx;^bx-_U)B{zcF~1>Cr{G;xbMrF+JQdlJ5ePPFwr2V5Xuos5@0NBio8y!5%+S108r
zrRBmc8YA4)3yd543+erKXHGZcoUrn?z8L3+<9K=REi;PON%F%;t-9xy;qsF>x;YEt
zV{dciqHgW$t$Z{MWnk4Jr`)51x8@3HL|Ik~TD?GomuYJju1i0ys|0`QhJJdZu5Rio
zse8$yck4yy?K&TPg0(VYBLUnq`n)etcpoiw8{R?Fvr}B_L8kYH?;Xw3<C@}xZcTB0
o0L^5gfPV4R_tt!KcV~$_nnrtcL1lCox!e`5>$coYx8ZL64F@4<_y7O^

literal 0
HcmV?d00001

diff --git a/prediction_thread_multiple.py b/prediction_thread_multiple.py
new file mode 100644
index 0000000..d39b203
--- /dev/null
+++ b/prediction_thread_multiple.py
@@ -0,0 +1,169 @@
+import json
+import gc
+import pandas as pd 
+import numpy as np 
+import pickle
+import tensorflow as tf
+import sys
+import time
+from test_prepare_functions import create_data_homology_ls,read_gene_sequences,create_tree_data,read_data,update_rest
+from process_data import create_map_list
+from threading import Thread,Lock
+import traceback
+from threads import Procerssrunner
+
+def main():
+    arg=sys.argv
+    fname=arg[-6]
+    model_name=arg[-5]
+    n_of_t=int(arg[-4])
+    st=int(arg[-3])
+    end=int(arg[-2])
+    name=arg[-1]
+
+    df=pd.read_csv(fname,sep="\t",header=None)
+    label_dict=dict(ortholog_one2one="Orthologs",
+                    other_paralog="Paralogs",
+                    ortholog_one2many="Orthologs",
+                    ortholog_many2many="Orthologs",
+                    within_species_paralog="Paralogs",
+                    gene_split="Gene Split",
+                    non_homolog="non_homolog")
+    label_dict_2=dict(ortholog_one2one=1,other_paralog=0,non_homolog=2,ortholog_one2many=1,ortholog_many2many=1,within_species_paralog=0,gene_split=4)
+    label_2=[]
+    label_1=[]
+    for _,row in df.iterrows():
+        label_1.append(label_dict_2[row[7]])
+        label_2.append(label_dict[row[7]])
+    df=df.assign(label=label_1)
+    df[7]=label_2
+    try:
+        if end<len(df):
+            if st<end:
+                df=df.loc[df.index.values[st:end]]
+            else:
+                raise ValueError()
+    except:
+        print("Making Predictions for the complete dataframe:)")
+    print(len(df))
+    
+    data={}    
+    with open("genome_maps","rb") as file:
+        data=pickle.load(file)
+    cmap=data["cmap"]
+    cimap=data["cimap"]
+    ld=data["ld"]
+    ldg=data["ldg"]
+    a=data["a"]
+    d=data["d"]
+    lsy=create_data_homology_ls(df,3,a,d,ld,ldg,cmap,cimap)
+    print(len(lsy))
+    data={}
+    a=[]
+    d=[]
+    ldg=[]
+    ld=[]
+    cmap=[]
+    cimap=[]
+    gc.collect()
+    
+    gene_sequences=read_gene_sequences(df,lsy,"geneseq","gene_sequences",name)
+    print("Gene Sequences Loaded")
+
+    print("Going to update not found sequences:")
+    gene_sequences=update_rest(gene_sequences,name)
+    
+    gc.collect()
+    part=len(df)//n_of_t
+    """
+    lsy={}
+    gene_sequences={}
+    part=1"""
+    return lsy,gene_sequences,df,model_name,fname,n_of_t,part,name
+
+
+if __name__=='__main__':
+    n=3
+    lsy,gene_sequences,df,model_name,fname,n_of_t,part,name=main()
+    pr=Procerssrunner()
+    pr.start_processes(n_of_t,df,gene_sequences,lsy,part,n,name)
+    smg,sml,indexes=read_data(n_of_t,name)
+    sml=np.array(sml)
+    smg=np.array(smg)
+    indexes=np.array(indexes)
+    print(indexes.shape)
+    df_temp=df.loc[indexes]
+    bls,blhs,dis,dps,dphs=create_tree_data("species_tree.tree",df_temp)
+
+    assert(len(bls)==len(indexes))
+    assert(len(blhs)==len(smg))
+    assert(len(df_temp)==len(dps))
+
+    print("Data Prepared.")
+    index_dict=create_map_list(indexes)    
+    preds=np.zeros((len(smg),3))
+    w=[0.86,0.8,0.06,0.0]
+    for i in range(1,4):
+        try:
+            model=tf.train.import_meta_graph(model_name+'_v'+str(i)+'/model.ckpt.meta')
+        except:
+            print("Something wrong with the model.")
+            continue
+        with tf.Session() as sess:
+            try:
+                model.restore(sess,model_name+'_v'+str(i)+"/model.ckpt")
+                graph = tf.get_default_graph()
+                synmgt,synmlt,blst,blhst,dpst,dphst,dist,lrt,yt=graph.get_collection("input_nodes")
+                predictions=graph.get_tensor_by_name("Predictions/BiasAdd:0")
+                print("Model Loaded Successfully :)")
+            except:
+                print(":(")
+                sys.exit()
+
+            fd={synmgt:smg,
+                synmlt:sml,
+                blst:bls,
+                blhst:blhs,
+                dpst:dps.reshape((len(blhs),1)),
+                dist:dis.reshape((len(blhs),1)),
+                dphst:dphs.reshape((len(blhs),1))}
+            preds_t_1=sess.run([predictions],feed_dict=fd)
+            preds_t_1=np.array(preds_t_1)[0]
+            fd={synmgt:smg.transpose((0,2,1,3)),
+                synmlt:sml.transpose((0,2,1,3)),
+                blst:blhs,
+                blhst:bls,
+                dpst:dphs.reshape((len(blhs),1)),
+                dist:dis.reshape((len(blhs),1)),
+                dphst:dps.reshape((len(blhs),1))}
+            preds_t_2=sess.run([predictions],feed_dict=fd)
+            preds_t_2=np.array(preds_t_2)[0]
+            preds=preds+w[i-1]*(preds_t_1+preds_t_2)/2
+        tf.reset_default_graph()
+    
+    #preds=preds/6
+    preds=np.argmax(preds,axis=1)
+    print(preds.shape)
+        
+    with open("prediction_"+fname+"_"+model_name+"_"+name+"_multiple.txt","w") as file:
+        for index,row in df.iterrows():
+            file.write(str(row[0]))
+            file.write("\t")
+            file.write(row[1])
+            file.write("\t")
+            file.write(row[3])
+            file.write("\t")
+            file.write(str(row["label"]))
+            file.write("\t")
+            if index in index_dict:
+                file.write(str(preds[index_dict[index]]))
+                file.write("\t")
+                if preds[index_dict[index]]==row["label"]:
+                    file.write(str(1))
+                else:
+                    file.write(str(0))
+            else:
+                file.write("Error")
+                file.write("\t")
+                file.write("NaN")
+            file.write("\n")
\ No newline at end of file
diff --git a/test_prepare_functions.py b/test_prepare_functions.py
index 80d72c2..b58ec4a 100644
--- a/test_prepare_functions.py
+++ b/test_prepare_functions.py
@@ -3,12 +3,13 @@ import progressbar
 import json
 import numpy as np 
 import time
-from multiprocessing import Process 
-from threading import Thread
+import pickle
 from ete3 import Tree
 from read_get_gene_seq import read_gene_seq,create_dict
 from process_data import get_nearest_neighbors
-from create_synteny_matrix import create_synteny_matrix_mul
+from create_synteny_matrix import create_synteny_matrix_mul,update
+import requests
+import json
 
 def create_data_homology_ls(df,n,a,d,ld,ldg,cmap,cimap):
     buf=open("not_found.txt","w")
@@ -54,7 +55,7 @@ def group_seq_by_species(df,g_to_sp):
     create_dict(ghsp,sph,g_to_sp)
     return g_to_sp
 
-def read_gene_sequences(df,lsy,data_dir,fname):
+def read_gene_sequences(df,lsy,data_dir,fname,name):
 
     """The basic idea here is to create a list/dictionary of all the genes by their species.
        Once the mapping is done, all the respective fasta sequence files are read by Species
@@ -101,10 +102,8 @@ def read_gene_sequences(df,lsy,data_dir,fname):
             except:
                 not_found[gene]=1
 
-    with open("processed/not_found.json","w") as file:
+    with open("processed/not_found_"+name+"_.json","w") as file:
         json.dump(not_found,file)
-    with open("processed/"+fname+".json","w") as file:#save the data
-        json.dump(data,file)
     return data
 
 def intermediate_process(gene_seq,x,y,n,index,sl,sg,ind):
@@ -115,63 +114,6 @@ def intermediate_process(gene_seq,x,y,n,index,sl,sg,ind):
     sl.append(smltemp)
     ind.append(index)
 
-def synteny_matrix(gene_seq,hdf,lsy,n,enable_break,sg,sl,ind):
-    #sg=[]
-    #sl=[]
-    t=0
-    #ind=[]
-    
-    for index,row in progressbar.progressbar(hdf.iterrows()):
-        g1=str(row[1])
-        g2=str(row[3])
-        x=[]
-        y=[]
-        t+=1
-        try:
-            temp=lsy[g1]
-        except:
-            continue
-        try:
-            temp=lsy[g2]
-        except:
-            continue
-        for i in range(len(lsy[g1]['b'])-1,-1,-1):
-            x.append(lsy[g1]['b'][i])
-        x.append(g1)
-        for k in lsy[g1]['f']:
-            x.append(k)
-
-        for i in range(len(lsy[g2]['b'])-1,-1,-1):
-            y.append(lsy[g2]['b'][i])
-        y.append(g2)
-        for k in lsy[g2]['f']:
-            y.append(k)
-        
-        assert(len(x)==len(y))
-        assert(len(x)==(2*n+1))
-        smgtemp,smltemp=create_synteny_matrix_mul(gene_seq,x,y,2*n+1)
-        if np.all(smgtemp==0):
-            continue
-        sg.append(smgtemp)
-        sl.append(smltemp)
-        ind.append(index)
-        """try:
-            th=Thread(target=intermediate_process,name="TimeOutDetector",args=(gene_seq,x,y,2*n+1,index,sl,sg,ind,))
-            th.start()
-            th.join(30)
-            if th.is_alive():
-                raise SystemError("Long Time")
-                print(row)
-                th.join()
-        except Exception as e:
-            print(e)"""
-        if t==5 and enable_break==1:
-            break
-    #print("Time Taken:",end-start)
-    #print("Average Time:",(end-start)/len(sg))
-    print(t)
-    return np.array(sg),np.array(sl),np.array(ind)
-
 
 def create_branch_length_padding(bl):
     maxlen=29
@@ -215,3 +157,61 @@ def create_tree_data(treename,df):
     create_branch_length_padding(branch_lengths_s)
     create_branch_length_padding(branch_lengths_hs)
     return np.array(branch_lengths_s),np.array(branch_lengths_hs),np.array(dist),np.array(ns),np.array(nhs)
+
+def read_data(nop,name):
+    smg=[]
+    sml=[]
+    indexes=[]
+    for i in range(nop):
+        try:
+            with open("temp_"+name+"/thread_"+str(i+1)+"_smg.temp","rb") as file:
+                smg=smg+pickle.load(file)
+            with open("temp_"+name+"/thread_"+str(i+1)+"_sml.temp","rb") as file:
+                sml=sml+pickle.load(file)
+            with open("temp_"+name+"/thread_"+str(i+1)+"_indexes.temp","rb") as file:
+                indexes=indexes+pickle.load(file)
+        except:
+            continue
+    print(len(indexes))
+    return smg,sml,indexes
+
+def update_rest(data,name):
+    gids={}
+    with open("processed/not_found_"+name+"_.json","r") as file:
+        gids=dict(json.load(file))
+
+    gids=list(gids.keys())
+
+    geneseq={}
+
+    server = "https://rest.ensembl.org"
+    ext = "/sequence/id?type=cds"
+    headers={ "Content-Type" : "application/json", "Accept" : "application/json"}
+
+    for i in progressbar.progressbar(range(0,len(gids)-50,50)):
+        ids=dict(ids=list(gids[i:i+50]))
+        while(1):
+            try:
+                r = requests.post(server+ext, headers=headers, data=str(json.dumps(ids)))
+                if not r.ok:
+                    r.raise_for_status()
+                gs=r.json()
+                tgs={}
+                for g in gs:
+                    tgs[g["query"]]=g["seq"]
+                geneseq.update(tgs)
+                break
+            except Exception as e:
+                print("Error:",e)
+                continue
+
+    data.update(geneseq)
+    for genes in gids:
+        try:
+            _=data[genes]
+        except:
+            print(genes)
+            update(data,genes)
+
+    print("Gene Sequences Updated Successfully")
+    return data