"""
Module.
Examples
--------
"""
__author__ = "Yeremia Gunawan Adhisantoso"
__email__ = "adhisant@tnt.uni-hannover.de"
__license__ = "Clear BSD"
__version__ = "1.0.0"
import numpy
from scipy.sparse.linalg import eigsh
[docs]
def get_Laplacian(M):
"""
Function get_Laplacian.
Parameters
----------
Returns
-------
Examples
--------
Notes
-----
"""
S=M.sum(1)
i_nz=numpy.where(S>0)[0]
S=S[i_nz]
M=(M[i_nz].T)[i_nz].T
S=1/numpy.sqrt(S)
M=S*M
M=(S*M.T).T
n=numpy.size(S)
M=numpy.identity(int(n))-M
M=(M+M.T)/2
return M
[docs]
def evec_distance(v1,v2):
"""
Function evec_distance.
Parameters
----------
Returns
-------
Examples
--------
Notes
-----
"""
d1=numpy.dot(v1-v2,v1-v2)
d2=numpy.dot(v1+v2,v1+v2)
if d1<d2:
d=d1
else:
d=d2
return numpy.sqrt(d)
[docs]
def get_ipr(evec):
"""
Function get_ipr.
Parameters
----------
Returns
-------
Examples
--------
Notes
-----
"""
ipr=1.0/(evec*evec*evec*evec).sum()
return ipr
[docs]
def get_score(M1,M2,num_evec):
"""
Function get_score.
Parameters
----------
Returns
-------
Examples
--------
Notes
-----
"""
k1=numpy.sign(M1).sum(1)
k2=numpy.sign(M2).sum(1)
d1=numpy.diag(M1)
d2=numpy.diag(M2)
kd1=~((k1==1)*(d1>0))
kd2=~((k2==1)*(d2>0))
iz=numpy.nonzero((k1+k2>0)*(kd1>0)*(kd2>0))[0]
M1b=(M1[iz].T)[iz].T
M2b=(M2[iz].T)[iz].T
i_nz1=numpy.where(M1b.sum(1)>0)[0]
i_nz2=numpy.where(M2b.sum(1)>0)[0]
i_z1=numpy.where(M1b.sum(1)==0)[0]
i_z2=numpy.where(M2b.sum(1)==0)[0]
M1b_L=get_Laplacian(M1b)
M2b_L=get_Laplacian(M2b)
a1, b1=eigsh(M1b_L,k=num_evec,which="SM")
a2, b2=eigsh(M2b_L,k=num_evec,which="SM")
b1_extend=numpy.zeros((numpy.size(M1b,0),num_evec))
b2_extend=numpy.zeros((numpy.size(M2b,0),num_evec))
for i in range(num_evec):
b1_extend[i_nz1,i]=b1[:,i]
b2_extend[i_nz2,i]=b2[:,i]
b1_extend_eff=b1_extend
b2_extend_eff=b2_extend
num_evec_eff=num_evec
evd=numpy.zeros(num_evec_eff)
for i in range(num_evec_eff):
evd[i]=evec_distance(b1_extend_eff[:,i],b2_extend_eff[:,i])
Sd=evd.sum()
l=numpy.sqrt(2)
evs=abs(l-Sd/num_evec_eff)/l
# logger.info("reproducibility score: %6.3f " %(evs))
return round(evs, 4)
[docs]
def get_hic_spector_train(M1,M2):
"""
Function get_hic_spector_train.
Parameters
----------
Returns
-------
Examples
--------
Notes
-----
"""
num_evec = 10 # number of eigenvectors
m1 = M1[0, 0].cpu().numpy()
m2 = M2[0, 0].cpu().numpy()
score = get_score(m1, m2,num_evec)
return score
[docs]
def get_hic_spector_test(M1,M2):
"""
Function get_hic_spector_test.
Parameters
----------
Returns
-------
Examples
--------
Notes
-----
"""
num_evec = 20 # number of eigenvectors
m1 = M1.cpu().numpy()
m2 = M2.cpu().numpy()
score = get_score(m1[0], m2[0],num_evec)
return score