This pipeline opens the result of ICAalamelodic.m, lets the user interactively label the components that look like neuronal activity (rather than movement artefacts or noise), sort them by label, plots a final summary for the chosen components, and save the reordered maps and time series.



In [1]:

    
import matplotlib
import numpy as np
import matplotlib.pyplot as plt
from scipy import stats
from scipy import io
%matplotlib inline 
import pylab

Open time series



In [2]:

    
import scipy.io as sio



In [3]:

    
# from http://stackoverflow.com/questions/3579568/choosing-a-file-in-python-with-simple-dialog
from Tkinter import Tk
from tkFileDialog import askopenfilename

Tk().withdraw() # we don't want a full GUI, so keep the root window from appearing
filename = askopenfilename() # show an "Open" dialog box and return the path to the selected file
print(filename)









    



/media/test7/THDDCGCaMP62/100142/100142Final/100142ss2onc250regcdFF20sMpsfkfintMB92Smith0_4_60TS.mat



In [4]:

    
Ua=sio.loadmat(filename)



In [5]:

    
DT=Ua['TSmean']



In [6]:

    
DT.shape









    Out[6]:





(11945, 92)



In [7]:

    
S1=DT.shape



In [8]:

    
DTmean=np.zeros(S1)
DTvar=np.zeros(S1)
Var=np.zeros(S1[1])



In [9]:

    
for i in range(S1[1]):
    DTmean[:,i]=DT[:,i]-np.mean(DT[:,i],0)



In [10]:

    
for i in range(S1[1]):
    Var[i]=np.sqrt(np.var(DTmean[:,i]))
    DTvar[:,i]=DTmean[:,i]/Var[i]



In [11]:

    
DTvar.shape









    Out[11]:





(11945, 92)

open maps



In [12]:

    
import nibabel as nb



In [13]:

    
# from http://stackoverflow.com/questions/3579568/choosing-a-file-in-python-with-simple-dialog
from Tkinter import Tk
from tkFileDialog import askopenfilename

Tk().withdraw() # we don't want a full GUI, so keep the root window from appearing
filename2 = askopenfilename() # show an "Open" dialog box and return the path to the selected file
print(filename2)









    



/media/test7/THDDCGCaMP62/100142/100142Final/100142ss2onc250regcdFF20sMpsfkfintMB92Smith0_4_60IC.nii



In [14]:

    
img1 = nb.load(filename2)



In [15]:

    
data = img1.get_data()



In [16]:

    
S=data.shape



In [17]:

    
S









    Out[17]:





(95, 44, 8, 92)

Zscore maps



In [18]:

    
Demean=np.zeros(S)
Dmaps=np.zeros(S)
Dvar=np.zeros(S)
Var=np.zeros(S[3])
D2=np.zeros([S[0],S[1],5,S[3]])
Tvar=np.zeros(S[3])

Transform the maps to have zero mean



In [19]:

    
for i in range(S[3]):
    Demean[:,:,:,i]=data[:,:,:,i]-np.mean(np.mean(np.mean(data[:,:,:,i],0),0),0)

Transform the maps to have unit variance and zscore



In [20]:

    
for i in range(S[3]):
    Dsq=np.reshape(Demean[:,:,:,i],S[0]*S[1]*S[2])
    Var[i]=np.sqrt(np.var(Dsq))
    Dvar=Demean[:,:,:,i]/Var[i]
    Dmaps[:,:,:,i]=Dvar-2.5
Dmaps[Dmaps<0]=0

Order ICs by variance



In [21]:

    
datao=data
Dmapso=Dmaps



In [22]:

    
plt.plot(Var)









    Out[22]:





[<matplotlib.lines.Line2D at 0x7f74bfb03c90>]

Separate maps in substacks, sort the independent components by brain regions



In [23]:

    
my_cmap=plt.cm.jet
my_cmap.set_bad(alpha=0)
Good_ICs=np.zeros(S[3])
Label_ICs=[]
pylab.rcParams['figure.figsize'] = (13, 2.5)



In [24]:

    
Dtemp=data[:,:,:,0]



In [25]:

    
%%javascript
IPython.OutputArea.auto_scroll_threshold =4000;



In [26]:

    
if S[2]>5:
    Nstack=5
    Int100=[(i+1)*100/Nstack for i in range(Nstack)]
    Percs=np.percentile(range(S[2]),Int100)
    Indices=np.split(range(S[2]),Percs)
    D1=np.zeros([S[0],S[1],Nstack])
    Dmean=Dtemp[:,:,range(Nstack)]
    for i in range(Nstack):
        Vmean=np.mean(Dtemp[:,:,Indices[i]],2)
        #Dmean[:,:,i]=np.max(Vmean,0)
        Dmean[:,:,i]=Vmean
else:
    Nstack=S[2]
    D1=np.zeros([S[0],S[1],S[2]])
    Dmean=data[:,:,range(S[2])]  
    Dmean=np.squeeze(Dtemp[:,:,:])









    



/usr/local/lib/python2.7/dist-packages/numpy/lib/shape_base.py:422: VisibleDeprecationWarning: using a non-integer number instead of an integer will result in an error in the future
  sub_arys.append(_nx.swapaxes(sary[st:end], axis, 0))



In [27]:

    
DTvar.shape









    Out[27]:





(11945, 92)



In [28]:

    
S









    Out[28]:





(95, 44, 8, 92)



In [29]:

    
# from http://stackoverflow.com/questions/3579568/choosing-a-file-in-python-with-simple-dialog
Tk().withdraw() # we don't want a full GUI, so keep the root window from appearing
filename = askopenfilename() # show an "Open" dialog box and return the path to the selected file
print(filename)









    



/media/test7/THDDCGCaMP62/100142/100142Final/Xk100142.mat



In [30]:

    
Ua=sio.loadmat(filename)
Xk=Ua['Xk']
#Xk[1,:]=Ua['Walk']



In [31]:

    
Xk.shape









    Out[31]:





(11945, 6)



In [32]:

    
# from http://stackoverflow.com/questions/3579568/choosing-a-file-in-python-with-simple-dialog
from Tkinter import Tk
from tkFileDialog import askopenfilename

Tk().withdraw() # we don't want a full GUI, so keep the root window from appearing
filenamet = askopenfilename() # show an "Open" dialog box and return the path to the selected file
print(filenamet)
nimt=nb.load(filenamet)
Dtemp=np.squeeze(nimt.get_data())
Dtemp.shape


if S[2]>5:
    Nstack=5
    Int100=[(i+1)*100/Nstack for i in range(Nstack)]
    Percs=np.percentile(range(S[2]),Int100)
    Indices=np.split(range(S[2]),Percs)
    Dmean=np.zeros([S[0],S[1],Nstack])
    #Dmean=np.squeeze(data[:,:,range(Nstack),2])
    for i in range(Nstack):
        Vmean=np.mean(Dtemp[:,:,Indices[i]],2)
        Dmean[:,:,i]=Vmean

plt.imshow(Vmean,cmap=plt.cm.gray)









    



/media/test7/THDDCGCaMP62/100142/100142Final/MAX_100142ss2onc250regcdFF20sMpsfkfintMB92Smith0_4_60IC.nii






    Out[32]:





<matplotlib.image.AxesImage at 0x7f74bf98bd90>



In [33]:

    
Xk=Xk.T



In [33]:

    
Label_ICs=[]



In [34]:

    
for j in range(S[3]):

    if S[2]>5:
        for i in range(Nstack):
            V=Dmaps[:,:,Indices[i],j]
            D1[:,:,i]=np.max(V,2)
        D2[:,:,:,j]=D1
        D1[D1==0]=np.nan
           
    else:
        for i in range(S[2]):
            V=Dmaps[:,:,i,j]
            D1[:,:,i]=V 
            

    print(j)
    for i in range(Nstack):
        plt.subplot(1,5,i+1)
        plt.imshow(Dmean[:,:,i],cmap=plt.cm.gray)
        plt.imshow(D1[:,:,i], cmap=my_cmap,interpolation='none')
        frame1 = plt.gca()
        frame1.axes.get_xaxis().set_visible(False)
        frame1.axes.get_yaxis().set_visible(False)
        
    plt.show()
    
   # plt.plot(TS_ROI[Order[j],:])
    plt.plot(DTvar[:,j])
    plt.plot(Xk[0,:]/np.std(Xk[0,:])+0.5,color=(1,0,0))   
    plt.plot(Xk[1,:]/np.std(Xk[1,:])+0.5,color=(0,1,0))
    plt.plot(Xk[2,:]/np.std(Xk[2,:])+0.5,color=(0.5,0.5,0))    
    #plt.plot(Xk[3,:]/np.std(Xk[1,:])+0.5,color=(0,0.5,1))
    
    plt.show()
    
    Label_ICs.append(raw_input())
    if Label_ICs[j]=='':
        Good_ICs[j]=0
    else:
        Good_ICs[j]=1









    



0






    



/usr/local/lib/python2.7/dist-packages/numpy/ma/core.py:4185: UserWarning: Warning: converting a masked element to nan.
  warnings.warn("Warning: converting a masked element to nan.")






    












    



/usr/local/lib/python2.7/dist-packages/ipykernel/__main__.py:29: RuntimeWarning: invalid value encountered in divide
/usr/local/lib/python2.7/dist-packages/ipykernel/__main__.py:30: RuntimeWarning: invalid value encountered in divide






    












    



1






    












    












    



2






    












    












    



3






    












    












    



4






    












    












    



5






    












    












    



6






    












    












    



7






    












    












    



8






    












    












    



9






    












    












    



alpha
10






    












    












    



11






    












    












    



12






    












    












    



13






    












    












    



alpha
14






    












    












    



15






    












    












    



alpha
16






    












    












    



alpha
17






    












    












    



18






    












    












    



19






    












    












    



alpha
20






    












    












    



21






    












    












    



gamma
22






    












    












    



23






    












    












    



24






    












    












    



25






    












    












    



26






    












    












    



27






    












    












    



28






    












    












    



beta
29






    












    












    



30






    












    












    



31






    












    












    



32






    












    












    



33






    












    












    



34






    












    












    



35






    












    












    



36






    












    












    



37






    












    












    



38






    












    












    



39






    












    












    



40






    












    












    



beta
41






    












    












    



42






    












    












    



43






    












    












    



44






    












    












    



45






    












    












    



46






    












    












    



47






    












    












    



48






    












    












    



49






    












    












    



50






    












    












    



51






    












    












    



52






    












    












    



53






    












    












    



54






    












    












    



55






    












    












    



56






    












    












    



57






    












    












    



58






    












    












    



59






    












    












    



60






    












    












    



61






    












    












    



62






    












    












    



63






    












    












    



64






    












    












    



65






    












    












    



66






    












    












    



67






    












    












    



68






    












    












    



69






    












    












    



beta
70






    












    












    



beta
71






    












    












    



72






    












    












    



73






    












    












    



74






    












    












    



75






    












    












    



76






    












    












    



77






    












    












    



78






    












    












    



79






    












    












    



80






    












    












    



81






    












    












    



82






    












    












    



83






    












    












    



84






    












    












    



85






    












    












    



86






    












    












    



87






    












    












    



88






    












    












    



89






    












    












    



90






    












    












    



91



In [36]:

    
#zip(range(S[3]),Label_ICs)



In [37]:

    
set(Label_ICs)









    Out[37]:





{'', 'alpha', 'gamma'}



In [38]:

    
#Label_ICs[94]='M'



In [ ]:



In [37]:

    
Xk=Xk.T



In [36]:

    
Xk.shape









    Out[36]:





(11945, 6)



In [38]:

    
Xksmoothed=np.zeros(Xk.shape)

Xksmoothed[0,:]=np.convolve(Xk[0,999:Xk.shape[1]-1000],np.ones(2000)/2000)

Xksmoothed[1,:]=np.convolve(Xk[1,999:Xk.shape[1]-1000],np.ones(2000)/2000)

Xksmoothed[2,:]=np.convolve(Xk[2,999:Xk.shape[1]-1000],np.ones(2000)/2000)

Xkdff=Xk-Xksmoothed



In [45]:

    
Rsq=np.zeros((1,S[3]))
Betas=np.zeros((3,S[3]))



In [46]:

    
from sklearn import linear_model



In [47]:

    
algorithm = linear_model.LinearRegression()



In [48]:

    
for j in range(S[3]):
    model = algorithm.fit(Xkdff[range(3),:].T, DT[:,j])
    Betas[:,j] = model.coef_
    Rsq[:,j] = model.score(Xkdff[range(3),:].T,DT[:,j])



In [49]:

    
max(max(Rsq))









    Out[49]:





0.73803527281687442



In [50]:

    
List1=[(Label_ICs[i],i) for i in range(S[3])]
Newlist=sorted(List1, key=lambda List1: List1[0])

Neworder=[Newlist[i][1] for i in range(S[3])]

NewDT=DTvar[:,Neworder[:]].T

for j in range(len(Neworder)):
    A=NewDT[:,j]
    V=np.sqrt(np.var(A))
    NewDT[:,j]=A/V

C1=np.zeros([16,3])
C1[0][:]=(1,0,0)
C1[1][:]=(0,1,0)
C1[2][:]=(0,0,1)
C1[3][:]=(0.8,0.8,0)
C1[4][:]=(0,1,1)
C1[5][:]=(1,0,1)
C1[6][:]=(1,0.5,0)
C1[7][:]=(0,1,0.5)
C1[8][:]=(0.5,0,1)
C1[9][:]=(0.8,0.8,0.5)
C1[10][:]=(0.5,1,1)
C1[11][:]=(1,0.5,1)
C1[12]=(0.5,0.5,0.5)
C1[13]=(0.2,0.5,0.5)
C1[14]=(0.5,0.2,0.5)
C1[15]=(0.5,0.5,0.2)
h=3

Newmaps=Dmaps[:,:,:,Neworder[:]]

L=len(set([Label_ICs[Neworder[i]] for i in range(len(Neworder))]))

Regionmaps=np.zeros([S[0],S[1],L,3])
Datasort=np.zeros([S[0],S[1],S[2],L,3])

Regionname=[]

DMapsordered=Dmapso[:,:,:,Neworder[:]]

j=0
i=0
k=Label_ICs[Neworder[0]]
m=0
Regionname.append(Label_ICs[Neworder[i]])
for i in range(len(Neworder)):
    
    #C2=C1[i%6][:]
    for l in range(3):
        M=np.max(np.squeeze(np.reshape(Newmaps[:,:,:,i],S[0]*S[1]*S[2])))
        Regionmaps[:,:,j,l]=Regionmaps[:,:,j,l]+0.6*np.max(DMapsordered[:,:,:,i],2)*C1[i%12+1][l]/M
        Datasort[:,:,:,j,l]=Datasort[:,:,:,j,l]+Dmaps[:,:,:,Neworder[i]]*C1[i%15+1][l] 
    i=i+1
    m=m+1
    if i<len(Neworder):
        k1=Label_ICs[Neworder[i]]
        
        
    if k1 != k:
        j=j+1
        k=k1
        m=0
        Regionname.append(Label_ICs[Neworder[i]])

pylab.rcParams['figure.figsize'] = (14, 5)
import scipy
from scipy import ndimage
j=0
m=0
L=0
k=Label_ICs[Neworder[0]]
for i in range(len(Neworder)):
    m=m+1
    
    
    if i<len(Neworder):
        k1=Label_ICs[Neworder[i]]
        
    if k1 != k:
        
        k=k1
        m=0
        
        plt.show()
        plt.figure(2*j+1)
        Rotated_Plot = ndimage.rotate(Regionmaps[:,:,j], -90)
        IM=plt.imshow(Rotated_Plot) 
        frame1 = plt.gca()
        frame1.axes.get_xaxis().set_visible(False)
        frame1.axes.get_yaxis().set_visible(False)
        j=j+1
        plt.figure(2*j)
        plt.plot(Xk[0,:]/np.std(Xk[0,:])+0.5,color=(1,0,0))   
        plt.plot(Xk[1,:]/np.std(Xk[1,:])+0.5,color=(0,1,0))
        plt.plot(Xk[2,:]/np.std(Xk[1,:])+0.5,color=(0.5,0.5,0))    
        #plt.plot(Xk[3,:]/np.std(Xk[1,:])+0.5,color=(0,0.5,1))
    plt.plot(NewDT[i,:]+h*m,color=C1[i%12+1][:])
    print(Neworder[i])
    print(Rsq[:,Neworder[i]])
    print(Betas[:,Neworder[i]])
plt.figure(2*j+1)
Rotated_Plot = ndimage.rotate(Regionmaps[:,:,j], -90)
IM=plt.imshow(Rotated_Plot)
frame1 = plt.gca()
frame1.axes.get_xaxis().set_visible(False)
frame1.axes.get_yaxis().set_visible(False)
print(Neworder)









    



0
[ 0.71396066]
[  7.68901928e-06   6.78309722e-06  -1.08424676e-06]
1
[ 0.6667187]
[  4.13280377e-06   2.43637153e-06  -1.44292196e-06]
2
[ 0.7356329]
[  2.25057999e-06   2.46042131e-06  -7.98514995e-08]
3
[ 0.52820035]
[  4.16977542e-06   2.46216267e-06  -1.80536057e-06]
4
[ 0.72091026]
[  4.08119806e-06   4.02259590e-06   2.33747385e-07]
5
[ 0.66948874]
[  2.12553362e-06   2.07504821e-06  -1.47424423e-07]
6
[ 0.64089588]
[  1.27589292e-06   1.19479436e-06  -6.29355164e-08]
7
[ 0.68388673]
[  7.44541518e-07   7.44140179e-07  -1.18686630e-07]
8
[ 0.11594127]
[ -1.24117304e-06  -1.34636169e-06   8.06889894e-07]
10
[ 0.42516305]
[  3.07901853e-06   1.28017036e-06  -1.22637600e-06]
11
[ 0.73803527]
[  2.80120256e-06   6.71266217e-07  -7.14453232e-07]
12
[ 0.51314828]
[  2.00525594e-06   2.57012160e-06  -8.06623647e-07]
14
[ 0.61155432]
[  2.47739734e-06   2.81479502e-06  -6.93465312e-07]
17
[ 0.31985058]
[  4.04667312e-07   4.47985387e-07  -1.03020136e-07]
18
[ 0.01034361]
[  4.89772143e-07  -5.79268123e-07  -1.34666360e-07]
20
[ 0.61676926]
[  4.43644144e-06   4.70722236e-06  -8.22084691e-07]
22
[ 0.62035323]
[  9.59100101e-07   9.52314054e-07  -3.92418489e-08]
23
[ 0.0089629]
[ -4.50951784e-07  -3.40528914e-07   4.42798929e-07]
24
[ 0.54884463]
[  2.48207532e-06   1.82705551e-06  -8.22437506e-07]
25
[ 0.60153391]
[  8.74952059e-07   8.38124220e-07  -2.10006243e-07]
26
[ 0.4553805]
[  1.06495349e-06   9.36030191e-07  -2.32437117e-07]
27
[ 0.05294751]
[ -5.24934788e-06  -4.23493796e-06   3.48303080e-06]
29
[ 0.51825464]
[ -2.91222906e-05  -2.88477826e-05   1.89204115e-06]
30
[ 0.5143058]
[ -6.78483082e-06  -4.57205800e-06   2.58819766e-06]
31
[ 0.02424544]
[ -3.03540722e-07  -2.13236740e-07   1.36604060e-07]
32
[ 0.1292088]
[ -1.40275742e-06  -1.78205470e-06   7.70453334e-07]
33
[ 0.03898803]
[ -3.65318461e-07  -3.70985530e-07   2.77248670e-07]
34
[ 0.39612532]
[ -1.74915924e-05  -1.38763778e-05   7.85518445e-06]
35
[ 0.60148521]
[  6.35938525e-07   6.17955961e-07  -7.36136936e-08]
36
[ 0.16941697]
[ -5.32535034e-07  -3.84528029e-07   2.52557839e-07]
37
[ 0.25516005]
[  1.74627519e-07   1.17865128e-06   6.31222829e-08]
38
[ 0.53695979]
[  4.21649130e-07   4.31077653e-07  -5.01972581e-08]
39
[ 0.55460713]
[  9.51596465e-07   9.08055685e-07  -6.01973625e-08]
41
[ 0.3768832]
[  1.32711427e-06   1.20604138e-06  -4.61309140e-07]
42
[ 0.51384288]
[  2.04955833e-06   2.38287742e-06  -3.89971612e-07]
43
[ 0.51175252]
[  1.78628870e-06   1.43939355e-06  -4.97429781e-07]
44
[ 0.04753798]
[ -2.75272211e-07  -3.01587607e-07   2.10128664e-07]
45
[ 0.30755747]
[ -7.34169975e-07  -6.17742788e-07   1.50200119e-07]
46
[ 0.46329214]
[  9.30661612e-07   8.63326840e-07  -2.38143729e-07]
47
[ 0.4130708]
[  6.49042640e-07   1.09563074e-06  -4.00023288e-08]
48
[ 0.31903248]
[ -1.70138451e-06  -1.06085525e-06   7.62101652e-07]
49
[ 0.21213565]
[  2.65946239e-07   3.52513812e-07   6.99862620e-08]
50
[ 0.02977441]
[ -2.37999711e-08  -3.09807262e-07   4.58958360e-08]
51
[ 0.32299241]
[ -5.53286216e-07  -5.28150164e-07   1.93290212e-07]
52
[ 0.30946836]
[ -1.37426155e-06  -1.13213053e-06   6.80439068e-07]
53
[ 0.32881097]
[  1.28509574e-06   1.35166051e-06  -6.79310217e-07]
54
[ 0.4281177]
[  3.85233078e-07   2.96046079e-07  -1.12352305e-07]
55
[ 0.4208409]
[ -3.16386937e-04  -4.80002267e-04   4.96126708e-05]
56
[ 0.02028715]
[  4.53348254e-08   3.46340747e-08   5.38848802e-08]
57
[ 0.07651597]
[  1.36200917e-06   2.06987548e-06  -7.83140139e-07]
58
[ 0.25733655]
[  4.34234168e-07   3.68642653e-07  -1.16408029e-07]
59
[ 0.38302066]
[  1.68771973e-06   1.17864552e-06  -4.37591401e-07]
60
[ 0.37412381]
[ -1.92460554e-06  -1.86532287e-06   7.06426345e-07]
61
[ 0.39493549]
[  1.49480971e-06   1.48196309e-06  -2.31360305e-08]
62
[ 0.15846807]
[ -1.83344204e-07  -2.13963224e-07   3.28744945e-08]
63
[ 0.1887232]
[ -9.59591963e-07  -2.03658434e-07   1.81083857e-07]
64
[ 0.24359928]
[ -3.99903027e-07  -8.81805541e-07  -3.18244260e-07]
65
[ 0.30862389]
[ -5.15806792e-06  -5.39172860e-06   1.04334814e-06]
66
[ 0.09480307]
[  2.38929751e-07   1.81631972e-07   4.17372812e-08]
67
[ 0.05450108]
[ -2.84918872e-07  -1.39678955e-06  -1.53939073e-07]
68
[ 0.27686129]
[ -2.77897606e-07  -3.88028752e-07  -1.01454699e-07]
71
[ 0.24604912]
[  1.86865794e-07   1.69599054e-07   9.61359397e-09]
72
[ 0.24558294]
[  6.16809231e-07   2.25861009e-07  -1.17599965e-07]
73
[ 0.30776404]
[  1.29511012e-06   1.85795066e-06  -9.26046944e-08]
74
[ 0.17762882]
[  1.55499471e-06   8.97090285e-07  -4.56916094e-07]
75
[ 0.03456288]
[ -7.03995832e-07  -1.18639084e-06   5.48734072e-07]
76
[ 0.14783182]
[  1.86913997e-07   3.16887745e-07  -3.75212133e-08]
77
[ 0.14503777]
[ -2.12128962e-06  -5.92525222e-07   4.64606900e-07]
78
[ 0.05241588]
[ -1.04540184e-06  -1.18096232e-06   5.76175074e-07]
79
[ 0.04324752]
[  4.10735245e-08  -2.71907160e-07   2.50905190e-07]
80
[ 0.19700388]
[ -2.48619695e-06  -2.11100390e-06   8.19793639e-07]
81
[ 0.15171951]
[  3.74551823e-07   7.05508113e-07   1.91286277e-07]
82
[ 0.0386612]
[ -3.59975924e-06  -4.07879981e-06   2.05125742e-06]
83
[ 0.13554778]
[  1.08638080e-06   1.34496847e-06  -5.82521678e-07]
84
[ 0.05912226]
[  4.97559567e-07   2.48156546e-07  -2.34927908e-07]
85
[ 0.0113403]
[ -3.06545988e-07  -2.31363832e-08   8.32665510e-08]
86
[ 0.01710032]
[  3.62928646e-08   1.50053560e-07   7.40371994e-08]
87
[ 0.07991076]
[ -6.11794091e-07  -7.64822944e-07   2.33439151e-07]
88
[ 0.01105098]
[ -1.77430998e-09   1.42331916e-07   1.59266574e-07]
89
[ 0.03116551]
[  5.19323775e-06   2.18479917e-06  -5.48539118e-06]
90
[ 0.00377477]
[ -2.19217211e-07  -1.97513434e-09   3.41441658e-07]
91
[ 0.00848126]
[ -1.30960947e-07  -3.47132266e-08   7.62495274e-08]






    












    



9
[ 0.02702032]
[ -1.55105383e-06  -1.29671567e-06   1.35551027e-06]
13
[ 0.07914616]
[  5.98524639e-06   5.19905776e-06  -3.91649480e-06]
15
[ 0.15704352]
[ -1.03408869e-06  -1.02123325e-06   5.52391107e-07]
16
[ 0.08540065]
[ -1.73737297e-06  -1.80289132e-06   1.04720391e-06]
19
[ 0.1205097]
[ -1.46857740e-06  -1.40683571e-06   1.15498725e-06]






    












    












    



28
[ 0.46972427]
[  6.62901675e-07   5.02223156e-07  -1.91148698e-07]
40
[ 0.00863091]
[  1.16262513e-06   2.01084514e-07  -1.45368523e-06]
69
[ 0.11593996]
[  3.07593276e-07   2.71233737e-07  -1.42185332e-07]
70
[ 0.11991367]
[ -5.16629375e-07  -6.81944985e-07   1.24355447e-07]






    












    












    



21
[ 0.02107053]
[  3.80515924e-08   2.59065605e-07   2.56259177e-07]
[0, 1, 2, 3, 4, 5, 6, 7, 8, 10, 11, 12, 14, 17, 18, 20, 22, 23, 24, 25, 26, 27, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 9, 13, 15, 16, 19, 28, 40, 69, 70, 21]



In [ ]:

    
# from http://stackoverflow.com/questions/3579568/choosing-a-file-in-python-with-simple-dialog
Tk().withdraw() # we don't want a full GUI, so keep the root window from appearing
filename = askopenfilename() # show an "Open" dialog box and return the path to the selected file
print(filename)



In [ ]:

    
Ua=sio.loadmat(filename)
Xk=Ua['Xk']



In [ ]:

    
plt.plot(Xkdff[range(3),:].T)



In [ ]: