jgrewe/gp_neurobio
1
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity

Merge branch 'master' of https://whale.am28.uni-tuebingen.de/git/jgrewe/gp_neurobio

Browse Source
This commit is contained in:
efish 2018-11-16 10:10:08 +01:00
commit 96ae72ff4f
2 changed files with 48 additions and 46 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

79
code/analysis_rs.py
View File

@ -1,58 +1,55 @@
import numpy as np
import matplotlib.pyplot as plt
from read_baseline_data import *
from IPython import embed
from NixFrame import *
from utility import *
from IPython import embed
# plot and data values
inch_factor = 2.54
data_dir = '../data'
dataset = '2018-11-09-ad-invivo-1'
# read eod and time of baseline
time, eod = read_baseline_eod(os.path.join(data_dir, dataset))
#fig = plt.figure(figsize=(12/inch_factor, 8/inch_factor))
#ax = fig.add_subplot(111)
#ax.plot(time[:1000], eod[:1000])
#ax.set_xlabel('time [ms]', fontsize=12)
#ax.set_ylabel('voltage [mV]', fontsize=12)
#plt.xticks(fontsize = 8)
#plt.yticks(fontsize = 8)
#fig.tight_layout()
fig, ax = plt.subplots(figsize=(12/inch_factor, 8/inch_factor))
ax.plot(time[:1000], eod[:1000])
ax.set_xlabel('time [ms]', fontsize=12)
ax.set_ylabel('voltage [mV]', fontsize=12)
plt.xticks(fontsize=8)
plt.yticks(fontsize=8)
fig.tight_layout()
#plt.savefig('eod.pdf')
plt.show()
#interspikeintervalhistogram, windowsize = 1 ms
#plt.hist
#coefficient of variation
#embed()
#exit()
# read spikes during baseline activity
spikes = read_baseline_spikes(os.path.join(data_dir, dataset))
# calculate interpike intervals and plot them
interspikeintervals = np.diff(spikes)
#fig = plt.figure()
#plt.hist(interspikeintervals, bins=np.arange(0, np.max(interspikeintervals), 0.0001))
#plt.show()
fig, ax = plt.subplots(figsize=(12/inch_factor, 8/inch_factor))
plt.hist(interspikeintervals, bins=np.arange(0, np.max(interspikeintervals), 0.0001))
plt.show()
# calculate coefficient of variation
mu = np.mean(interspikeintervals)
sigma = np.std(interspikeintervals)
cv = sigma/mu
#print(cv)
# calculate zero crossings of the eod
# plot mean of eod circles
# plot std of eod circles
# plot psth into the same plot
# calculate vector strength
print(cv)
threshold = 0;
# calculate eod times and indices by zero crossings
threshold = 0
shift_eod = np.roll(eod, 1)
eod_times = time[(eod >= threshold) & (shift_eod < threshold)]
sampling_rate = 40000.0
eod_idx = eod_times*sampling_rate
# align eods and spikes to eods
max_cut = int(np.max(np.diff(eod_idx)))
eod_cuts = np.zeros([len(eod_idx)-1, max_cut])
# eods 15 + 16 are to short
relative_times = []
spike_times = []
eod_durations = []
for i, idx in enumerate(eod_idx[:-1]):
eod_cut = eod[int(idx):int(eod_idx[i+1])]
@ -60,34 +57,38 @@ for i, idx in enumerate(eod_idx[:-1]):
eod_cuts[i, len(eod_cut):] = np.nan
time_cut = time[int(idx):int(eod_idx[i+1])]
spike_cut = spikes[(spikes > time_cut[0]) & (spikes < time_cut[-1])]
relative_time = spike_cut - time_cut[0]
if len(relative_time) > 0:
relative_times.append(relative_time[:][0]*1000)
spike_time = spike_cut - time_cut[0]
if len(spike_time) > 0:
spike_times.append(spike_time[:][0]*1000)
eod_durations.append(len(eod_cut)/sampling_rate*1000)
# calculate vector strength
vs = vector_strength(spike_times, eod_durations)
# determine means and stds of eod for plot
# determine time axis
mu_eod = np.nanmean(eod_cuts, axis=0)
std_eod = np.nanstd(eod_cuts, axis=0)*3
time_axis = np.arange(max_cut)/sampling_rate*1000
#fig = plt.figure(figsize=(12/inch_factor, 8/inch_factor))
# plot eod form and spike histogram
fig, ax1 = plt.subplots(figsize=(12/inch_factor, 8/inch_factor))
ax1.hist(relative_times, color='crimson')
ax1.hist(spike_times, color='crimson')
ax1.set_xlabel('time [ms]', fontsize=12)
ax1.set_ylabel('number', fontsize=12)
ax1.tick_params(axis='y', labelcolor='crimson')
plt.yticks(fontsize = 8)
plt.yticks(fontsize=8)
ax1.spines['top'].set_visible(False)
ax2 = ax1.twinx()
ax2.fill_between(time_axis, mu_eod+std_eod, mu_eod-std_eod, color='dodgerblue', alpha=0.5)
ax2.plot(time_axis, mu_eod, color='black', lw=2)
ax2.set_ylabel('voltage [mV]', fontsize=12)
ax2.tick_params(axis='y', labelcolor='dodgerblue')
plt.xticks(fontsize = 8)
plt.yticks(fontsize = 8)
plt.xticks(fontsize=8)
plt.yticks(fontsize=8)
fig.tight_layout()
plt.show()
#NixToFrame(data_dir)

15
code/utility.py
View File

@ -1,8 +1,8 @@
import numpy as np
def zero_crossing(eod,time):
threshold = 0;
def zero_crossing(eod, time):
threshold = 0
shift_eod = np.roll(eod, 1)
eod_times = time[(eod >= threshold) & (shift_eod < threshold)]
sampling_rate = 40000.0
@ -10,9 +10,10 @@ def zero_crossing(eod,time):
return eod_idx
def vector_strength(spike_times, eod_durations)
alphas = spike_times/ eod_durations
cs = (1/len(spike_times))*np.sum(np.cos(alphas))^2
sn = (1/len(spike_times))*np.sum(np.sin(alphas))^2
vs = np.sprt(cs+sn)
def vector_strength(spike_times, eod_durations):
n = len(spike_times)
phase_times = np.zeros(len(spike_times))
for i, idx in enumerate(spike_times):
phase_times[i] = (spike_times[i] / eod_durations[i]) * 2 * np.pi
vs = np.sqrt((1/n*sum(np.cos(phase_times)))**2 + (1/n*sum(np.sin(phase_times)))**2)
return vs