In [1]:

    

from __future__ import division
%pylab inline
import pandas as pd
from Bio import SeqIO
import nwalign
import itertools
from scipy import signal
import mahotas
import sklearn
from skimage.feature import local_binary_pattern,match_template,structure_tensor,structure_tensor_eigvals
from skimage.transform import rotate


    

    




Populating the interactive namespace from numpy and matplotlib

In [6]:

    

def scaleMatrix(mat):
    mMin = mat.min()
    mMax = mat.max()
    
    s = mat-mMin
    s =255*mat/(mMax-mMin)
    
    return s

def hist(ax, lbp):
    n_bins = lbp.max() + 1
    return ax.hist(lbp.ravel(), normed=True, bins=n_bins, range=(0,n_bins),facecolor='0.5')


    

In [37]:

    

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


    

In [ ]:

    

for i,mat in matD.iteritems():
    print i
    
    plt.figure(figsize=(15,4))
    subplot(121)
    imshow(mat,interpolation='none')

    numPts = 8
    histo=[]
    
    for j in range(1,400,5):
        lbp = local_binary_pattern(mat,numPts,j,'uniform')
        nBins = int(lbp.max()+1)
        
        hgram,edges=np.histogram(lbp.ravel(),normed=False,bins=nBins)
        
        histo+=[hgram.tolist()]
            
    hMatrix = np.matrix(histo)    
    subplot(122)
    imshow(hMatrix,interpolation='none',aspect=0.25)
#     subplot(122)
#     plot(histo)
    plt.show()


    

    




ERROR! Session/line number was not unique in database. History logging moved to new session 254
S1_775






    













    




S1_774






    













    




S1_478






    













    




S1_479






    













    




S1_771






    













    




S1_770






    













    




S1_773






    













    




S1_772






    













    




S1_472






    













    




S1_473






    













    




S1_470






    













    




S1_471






    













    




S1_476






    













    




S1_477






    













    




S1_474






    













    




S1_475






    













    




S1_776






    

In [ ]:

    

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 makeTemplates(minPeriod,maxPeriod,repLen,step):
    
    size = int(maxPeriod*repLen/step)
    
    adjPeriods = [ int(np.round(p*repLen/step)) for p in range(minPeriod,maxPeriod+1)]
    matSizes = [ size for aP in adjPeriods ]
    
    tList = [ getTemplateMatrix(i,s) for i,s in zip(adjPeriods,matSizes)]
    
    return tList


# refs = makeTemplates(2,30,170,34)

for p in range(2,20):

    print 'p=', p
    
    tMat = getTemplateMatrix(p,500)

    plt.figure(figsize=(15,4))
    subplot(141)
    imshow(tMat, interpolation='none',cmap=cm.gray_r)

    r1 = p-1
    np1 = 8

    r2 = p
    np2 = 8
    
    r3 = p+1
    np3 = 8

    lbp1 = local_binary_pattern(tMat,np1,r1,'uniform')
    lbp2 = local_binary_pattern(tMat,np2,r2,'uniform')
    lbp3 = local_binary_pattern(tMat,np3,r3,'uniform')

    ax1=subplot(142)
    ax2=subplot(143)
    ax3=subplot(144)

    hist(ax1,lbp1)
    hist(ax2,lbp2)
    hist(ax3,lbp3)

    plt.show()