Anesthesia_CA1/validation/calcium_imaging/Validation_motion.py

#!/usr/bin/env python
# coding: utf-8

database_path = '/media/andrey/My Passport/GIN_new/Anesthesia_CA1/meta_data/meta_recordings - anesthesia.xlsx'


# ### Select the range of recordings for the analysis (see "Number" row in the meta data file)

# In[4]:
rec = [x for x in range(0,188)]


# In[1]:

import numpy as np

import numpy.ma as ma

import matplotlib.pyplot as plt
import matplotlib.ticker as ticker

import pandas as pd
import seaborn as sns
import pickle
import os

sns.set()
sns.set_style("whitegrid")

from scipy.signal import medfilt 

from scipy.stats import skew, kurtosis, zscore

from scipy import signal

from sklearn.linear_model import LinearRegression, TheilSenRegressor


plt.rcParams['figure.figsize'] = [8, 8]


# In[2]:


from capipeline import *

motion_index = np.zeros((len(rec),100000),dtype='float')

recording_length = np.zeros((len(rec)),dtype='int')

number_of_quiet_periods = 6

start_quiet_period = np.zeros((len(rec),number_of_quiet_periods),dtype='int')

stop_quiet_period = np.zeros((len(rec),number_of_quiet_periods),dtype='int')

n_quiet_periods_in_recording = np.zeros((len(rec)),dtype='int')


for i, r in enumerate(rec):

    print("Reconrding # ", r)
    
    animal = get_animal_from_recording(r, database_path)
        
    #condition = get_condition(r, database_path)
    
    #print("#" +  str(r) + " " + str(animal) + " " + str(condition) + " ")

    meta_data = pd.read_excel(database_path)
 
    path_excel_rec = str(meta_data['Folder'][r]) + str(meta_data['Subfolder'][r])  + str(meta_data['Recording idx'][r]) + '/suite2p/'

    stat = np.load(path_excel_rec+ '/plane0/ops.npy', allow_pickle=True)

    #plt.plot(stat.item(0)['yoff'],alpha=0.5)
    #plt.plot(stat.item(0)['xoff'],alpha=0.5)

    motion_index[i,:len(stat.item(0)['yoff'])] = np.sqrt(stat.item(0)['yoff']**2 + stat.item(0)['xoff']**2) 

    print("Min motion index:", min(motion_index[i,:]))
    print("Max motion index:", max(motion_index[i,:]))

    recording_length[i] = len(stat.item(0)['yoff'])

    print("Recording length:", recording_length[i])

 

    if (type(meta_data['Quiet periods'][r]) == type('str')):
        print("Quiet periods: ", meta_data['Quiet periods'][r].split(','))

        n_quiet_periods_in_recording[i] = int(len(meta_data['Quiet periods'][r].split(','))/2)

        for k in range(n_quiet_periods_in_recording[i]):
            start_quiet_period[i,k] = int(meta_data['Quiet periods'][r].split(',')[k*2])
            stop_quiet_period[i,k] = int(meta_data['Quiet periods'][r].split(',')[k*2+1])  
            print(start_quiet_period[i,k])
            print(stop_quiet_period[i,k] )
    else:
        print("Nan for recording:", r)
        start_quiet_period[i,0] = 0
        stop_quiet_period[i,0] = recording_length[i]

#np.save("./xy-motion.npy",motion_index)

#mi = motion_index

mi = np.load("./xy-motion.npy")

print(np.max(mi))

print(np.min(mi))

mi_av = np.mean(mi[:,:18000].reshape(mi.shape[0],60,300), axis=2)
plt.rcParams["axes.grid"] = False

plt.figure(figsize = (10,15))

pos = plt.imshow(mi_av,cmap='Purples',vmin = 0, vmax = 10) #,aspect='equal'

plt.colorbar(pos)


for i in range(len(rec)):
    for k in range(int(n_quiet_periods_in_recording[i])):
        line = plt.scatter(start_quiet_period[i,k]/300,i,marker='>',color='k',s = 10,label='start of quiet period')
        line.set_clip_on(False)
        line = plt.scatter(stop_quiet_period[i,k]/300,i,marker='<',color='k',s = 10, label='end of quiet period')
        line.set_clip_on(False)

line = plt.scatter(recording_length[0:len(rec)]/300, np.arange(len(rec)) ,marker='|',color='k', s = 10, label='end of the recording')
line.set_clip_on(False)


#plt.legend(loc='upper right')

plt.xlabel('x 300 frames')
plt.ylabel('recording')

plt.title('Anesthesia dataset motion validation')

#plt.gcf().set_facecolor("white")

plt.xlim([0,60])

line.set_clip_on(False)

plt.savefig("Validation_motion.png")
plt.savefig("Validation_motion.svg")



plt.show()