In [64]:

    

from __future__ import division
%pylab inline
import networkx as nx
from scipy.sparse import lil_matrix, kron, identity
from scipy.sparse.linalg import lsqr
from skimage.filters import threshold_otsu

def getNormAdjMatrix(graph):
    
    adjMat = nx.to_numpy_matrix(graph)
    norm = adjMat.sum(axis=0)
    norm[norm==0] = 1
    return adjMat/norm

def compare (G1,G2,alpha):
    
    normAdj1 = getNormAdjMatrix(G1)
    normAdj2 = getNormAdjMatrix(G2)
    
    wProd = kron(lil_matrix(normAdj1),lil_matrix(normAdj2))
    
    startP = np.ones(wProd.shape[0])/(wProd.shape[0])
    stopP = startP
    A = identity(wProd.shape[0]) - (wProd * alpha)
    x = lsqr(A,startP)
        
    return x[0].sum()
    

def getTemplateMatrix(period,size):
    tMat = np.random.random((size,size))
    np.fill_diagonal(tMat,val=2)
    
    for k in range(period,size,period):
        tMat += np.diag(np.ones(size-k),k)
        tMat += np.diag(np.ones(size-k),-k)
    
    return tMat

def binarizeMatrix(mat):
    thresh = threshold_otsu(mat)
    return mat > thresh


    

    




Populating the interactive namespace from numpy and matplotlib






    




WARNING: pylab import has clobbered these variables: ['identity', 'kron']
`%matplotlib` prevents importing * from pylab and numpy

In [107]:

    

matD=np.load('../simulate_sequence/len.gte.8k/ref.depth30.seed0.len.gte.8k_000.npz')

for i,mat in matD.iteritems():
    
    M=binarizeMatrix(mat)
    np.fill_diagonal(M,0)
    M1 = nx.to_networkx_graph(M,create_using=nx.Graph())
#     print i, '\t', nx.number_connected_components(M1)
    imshow(M,interpolation='none',cmap=cm.gray)
#     nx.draw(M1,M1.edges(),edge_color='b')
    plt.show()


    

    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    




---------------------------------------------------------------------------
NetworkXError                             Traceback (most recent call last)
<ipython-input-107-6fb5ea616a93> in <module>()
      5     M=binarizeMatrix(mat)
      6     np.fill_diagonal(M,0)
----> 7     M1 = nx.to_networkx_graph(M,create_using=nx.Graph())
      8 #     print i, '\t', nx.number_connected_components(M1)
      9     imshow(M,interpolation='none',cmap=cm.gray)

/Users/skasinat/anaconda/lib/python2.7/site-packages/networkx/convert.pyc in to_networkx_graph(data, create_using, multigraph_input)
    148             except:
    149                 raise nx.NetworkXError(\
--> 150                   "Input is not a correct numpy matrix or array.")
    151     except ImportError:
    152         warnings.warn('numpy not found, skipping conversion test.',

NetworkXError: Input is not a correct numpy matrix or array.

In [116]:

    

nx.number_connected_components(nx.to_networkx_graph(binarizeMatrix(getTemplateMatrix(6,10)),create_using=nx.Graph()))


    

    
Out[116]:





6

In [9]:

    

pcg?