Validation_motion.py 3.0 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132
  1. #!/usr/bin/env python
  2. # coding: utf-8
  3. database_path = '/media/andrey/My Passport/GIN/Anesthesia_CA1/meta_data/meta_recordings_transition_state.xlsx'
  4. # ### Select the range of recordings for the analysis (see "Number" row in the meta data file)
  5. # In[4]:
  6. rec = [x for x in range(0,198)] #198
  7. # In[1]:
  8. import numpy as np
  9. import numpy.ma as ma
  10. import matplotlib.pyplot as plt
  11. import matplotlib.ticker as ticker
  12. import pandas as pd
  13. import seaborn as sns
  14. import pickle
  15. import os
  16. sns.set()
  17. sns.set_style("whitegrid")
  18. from scipy.signal import medfilt
  19. from scipy.stats import skew, kurtosis, zscore
  20. from scipy import signal
  21. from sklearn.linear_model import LinearRegression, TheilSenRegressor
  22. plt.rcParams['figure.figsize'] = [8, 8]
  23. # In[2]:
  24. from capipeline import *
  25. motion_index = np.zeros((len(rec),100000),dtype='float')
  26. recording_length = np.zeros((len(rec)),dtype='int')
  27. start_quiet_period = np.zeros((len(rec)),dtype='int')
  28. stop_quiet_period = np.zeros((len(rec)),dtype='int')
  29. for i, r in enumerate(rec):
  30. print("Reconrding # ", r)
  31. animal = get_animal_from_recording(r, database_path)
  32. #condition = get_condition(r, database_path)
  33. #print("#" + str(r) + " " + str(animal) + " " + str(condition) + " ")
  34. meta_data = pd.read_excel(database_path)
  35. path_excel_rec = str(meta_data['Folder'][r]) + str(meta_data['Subfolder'][r]) + str(meta_data['Recording idx'][r]) + '/suite2p/'
  36. stat = np.load(path_excel_rec+ '/plane0/ops.npy', allow_pickle=True)
  37. #plt.plot(stat.item(0)['yoff'],alpha=0.5)
  38. #plt.plot(stat.item(0)['xoff'],alpha=0.5)
  39. motion_index[i,:len(stat.item(0)['yoff'])] = np.sqrt(stat.item(0)['yoff']**2 + stat.item(0)['xoff']**2)
  40. print(min(motion_index[i,:]))
  41. print(max(motion_index[i,:]))
  42. recording_length[i] = len(stat.item(0)['yoff'])
  43. print(meta_data['Quiet periods'][r].split(','))
  44. start_quiet_period[i] = int(meta_data['Quiet periods'][r].split(',')[0])
  45. stop_quiet_period[i] = int(meta_data['Quiet periods'][r].split(',')[1])
  46. #np.save("./xy-motion.npy",motion_index)
  47. mi = motion_index
  48. #mi = np.load("./xy-motion.npy")
  49. print(np.max(mi))
  50. print(np.min(mi))
  51. mi_av = np.mean(mi[:,:30000].reshape(mi.shape[0],100,300), axis=2)
  52. plt.rcParams["axes.grid"] = False
  53. plt.figure(figsize = (10,15))
  54. plt.imshow(mi_av,cmap='Purples',vmin = 0, vmax = 10) #,aspect='equal'
  55. line = plt.scatter(start_quiet_period[0:len(rec)]/300,np.arange(len(rec)),marker='>',color='k',s=10,label='start of quiet period')
  56. line.set_clip_on(False)
  57. line = plt.scatter(stop_quiet_period[0:len(rec)]/300,np.arange(len(rec)),marker='<',color='k',s=10,label='end of quiet period')
  58. line.set_clip_on(False)
  59. line = plt.scatter(recording_length[0:len(rec)]/300, np.arange(len(rec)) ,marker='|',color='k',s=10,label='end of the recording')
  60. line.set_clip_on(False)
  61. plt.legend()
  62. plt.xlabel('x 300 frames')
  63. plt.ylabel('recording')
  64. plt.title('Transition state dataset motion validation')
  65. #plt.gcf().set_facecolor("white")
  66. plt.xlim([0,100])
  67. plt.savefig("Validation_motion.png")
  68. plt.savefig("Validation_motion.svg")
  69. plt.show()