diff --git a/main.py b/main.py index 34f26ec..e2770d5 100644 --- a/main.py +++ b/main.py @@ -1,6 +1,7 @@ import sys from get_data import get_data_homology,get_data_genome +from process_data import create_data_homology_ls arg=sys.argv arg=arg[1:] @@ -18,10 +19,24 @@ ldg=[] dir_g="data" ld,ldg,a,d=get_data_genome(arg,dir_g,a,d,ld,ldg) -#print(a[0][0:10],"\n",a[2][0:10],"\n",d,"\n",ld[0][0:10],"\n",ld[2][0:10]) +#print(a[0][0:10],"\n",d,"\n",ld[0][0:10],"\n") dir_hom="data_homology" a_h,d_h=get_data_homology(arg,dir_hom,a_h,d_h) #print(a_h[0][0:10],"\n",d_h) -x=input() +n=2 #no. of numbers neighbors + +lsy,lcmap=create_data_homology_ls(a_h,d_h,n,a,d,ld,ldg) +""" +print(lsy,"\n",lcmap) +lt=ldg[0] +for x in lsy: + + xr=lsy[x]['f'] + xl=lsy[x]['b'] + for g in range(len(xl)-1,-1,-1): + print(a[0].iloc[lt[xl[g]],[3,4]]) + print("------------------\n",a[0].iloc[lt[x],[3,4]],"\n-------------------") + for g in xr: + print(a[0].iloc[lt[g],[3,4]])""" diff --git a/process_data.py b/process_data.py index ee8202d..a9dd9b0 100644 --- a/process_data.py +++ b/process_data.py @@ -1,6 +1,50 @@ import pandas import gc +import numpy as np +def get_nearest_neighbors(g,gs,n,a,d,ld,ldg): + ne=[] + nr=[] + gi=d[gs.capitalize()] + sldf=a[gi] + sld=ld[gi] + sldg=ldg[gi] + if g not in sldg: + return ne + i=sldg[g] + #get the -n neighbors + start=int(sldf.iloc[i,[3]]) + for j in range(n): + #select the column + itemp=0 + end=sldf.iloc[:,4] + end=np.array(end) + assert(len(end)==len(sld)) + end=end-start + end_s=np.argsort(end) + for k in end_s: + if end[k]<0 and end[k+1]>=0: + itemp=k + break + ne.append(sld[itemp]) + start=int(sldf.iloc[itemp,[3]]) + #print(start) + #get the +n neighbors + end=int(sldf.iloc[i,[4]]) + for j in range(n): + itemp=0 + start=sldf.iloc[:,3] + start=np.array(start) + start=start-end + start_s=np.argsort(start) + for k in start_s: + if start[k]>0: + itemp=k + break + nr.append(sld[itemp]) + end=int(sldf.iloc[itemp,[4]]) + + return ne,nr ls=[] ld=[] @@ -13,3 +57,27 @@ def list_dict_genomes(a,n): ls.append(uc) ld.append(ldg) return ls,ld + +def create_data_homology_ls(a_h,d_h,n,a,d,ld,ldg): + lsy={} #dictionary which stores +/- n genes of the given gene by id. Each key is a gene id which corresponds to the one in center. + lcmap={} #dictionary which stores the gene pairs already considered + t=0 + for df in a_h: + for index,row in df.iterrows(): + x=row["gene_stable_id"] + y=row["homology_gene_stable_id"] + xs=row["species"] + ys=row["homology_species"] + if x+y in lcmap or y+x in lcmap: + continue + if x not in lsy: + xl,xr=get_nearest_neighbors(x,xs,n,a,d,ld,ldg) + lsy[x]=dict(b=xl,f=xr) + if y not in lsy: + yarr=[] + yl,yr=get_nearest_neighbors(y,ys,n,a,d,ld,ldg) + lsy[y]=dict(b=yl,f=yr) + lcmap[x+y]=1 + + + return lsy,lcmap