susceptibility1/data_overview_mod.py
2024-03-15 17:16:02 +01:00

571 lines
25 KiB
Python

#sys.path.insert(0, '..')
#from plt_RAM import plt_RAM_overview_nice
#from utils_susept import
from IPython import embed
from matplotlib import gridspec as gridspec, pyplot as plt
import numpy as np
from utils_all_down import default_settings
from utils_suseptibility import colors_overview
from utils_suseptibility import default_figsize, NLI_scorename2,pearson_label, exclude_nans_for_corr, kernel_scatter, \
plt_burst_modulation_hists, \
version_final
from utils_all import update_cell_names, load_overview_susept, make_log_ticks, p_units_to_show, save_visualization, setting_overview_score
from scipy import stats
try:
from plotstyle import plot_style, spines_params
except:
a = 5
def data_overview3():
# calcdf_RAM_overview()
save_name = 'calc_RAM_overview-noise_data8_nfft1sec_original__LocalEOD_CutatBeginning_0.05_s_NeurDelay_0.005_s_burst_corr'
save_name = 'calc_RAM_overview-noise_data8_nfft1sec_original__LocalEOD_CutatBeginning_0.05_s_NeurDelay_0.005_s'
save_name = 'calc_RAM_overview-noise_data9_nfft1sec_original__StimPreSaved4__CutatBeginning_0.05_s_NeurDelay_0.005_s'
save_name = 'calc_RAM_overview-noise_data9_nfft1sec_original__StimPreSaved4__mean5__CutatBeginning_0.05_s_NeurDelay_0.005_s'
col = 4
row = 2 # sharex=True,
plot_style()
default_figsize(column=2, length=6.2) #65.5
#fig, ax = plt.subplots(4, 2) # , figsize=(14, 7.5) constrained_layout=True,
var_it = 'Response Modulation [Hz]'
var_it2 = ''
if var_it == '':
ws = 0.35
else:
ws = 0.65
if var_it2 != '':
right = 0.9
else:
right = 0.98
right = 0.85
ws = 0.75
print(right)
grid0 = gridspec.GridSpec(3, 2, wspace=ws, bottom=0.07,
hspace=0.45, left=0.1, right=right, top=0.95)
###################################
###############################
# Das ist der Finale Score
scoreall = 'perc99/med_diagonal_proj'
scoreall = 'max(diag5Hz)/med_diagonal_proj_fr'#'max(diag5Hz)/med_diagonal_proj_fr_base_w_burstcorr'#'perc99/med_diagonal_proj'
#'max(diag5Hz)/med_diagonal_proj_fr','max(diag5Hz)/med_diagonal_proj_fr_base_w_burstcorr',
###################################
#scores = [scoreall+'_diagonal_proj']
##########################
# Auswahl: wir nehmen den mean um nicht Stimulus abhängigen Noise rauszumitteln
#save_names = []
save_names = ['calc_RAM_overview-_simplified_'+version_final()]#'calc_RAM_overview-_simplified_noise_data9_nfft1sec_original__StimPreSaved4__mean5__CutatBeginning_0.05_s_NeurDelay_0.005_s__burstIndividual_','calc_RAM_overview-noise_data9_nfft1sec_original__StimPreSaved4__mean5__CutatBeginning_0.05_s_NeurDelay_0.005_s__burstIndividual_',
save_names = ['calc_RAM_overview-_simplified_noise_data12_nfft0.5sec_original__StimPreSaved4__abs__direct_']
save_names = ['calc_RAM_overview-_simplified_noise_data12_nfft0.5sec_original__StimPreSaved4__direct_']
#save_names = ['calc_RAM_overview-_simplified_noise_data12_nfft0.5sec_original__StimPreSaved4__abs_']
#####################################################
#grid_lower_lower = gridspec.GridSpecFromSubplotSpec(1, 2, grid0[1], wspace = 0.5, hspace=0.55)#, height_ratios = [1,3]
cell_types = [' P-unit',' Ampullary', ]#, ' P-unit',]#' P-unit',
cell_types_name = ['P-units','Ampullary cells',]
species = ' Apteronotus leptorhynchus'
burst_corr_reset = 'response_modulation'
burst_fraction = [1000, 1000]# 50, # ,1,1]
burst_fraction = [1, 1] # ,1,1]
burst_corr_reset = 'burst_fraction_burst_corr_individual_stim'
redo = False
#embed()
counter = 0
tags = []
frame_load_sp = load_overview_susept(save_names[0], redo=redo, redo_class=redo)
scores = ['max(diag5Hz)/med_diagonal_proj_fr','max(diag5Hz)/med_diagonal_proj_fr',
'max(diag5Hz)/med_diagonal_proj_fr_base_w_burstcorr'] # + '_diagonal_proj'
scores = ['max(diag5Hz)/med_diagonal_proj_fr','max(diag5Hz)/med_diagonal_proj_fr',
] # + '_diagonal_proj'
max_xs = [[[],[],[]],[[],[],[]]]
for c, cell_type_here in enumerate(cell_types):
frame_file = setting_overview_score(frame_load_sp, cell_type_here, min_amp='range', species=species)
#embed()
test = False
# ok das schließe ich aus weil da irgendwas in der Detektion ist, das betrifft jetzt genau 3 Zellen, also nicht so schlimm
#63 2018-08-14-af-invivo-1
#241 2018-09-05-aj-invivo-1
#252 2022-01-08-ah-invivo-1
frame_file = frame_file[frame_file.cv_stim <5]
if test:
frame_file[frame_file.cv_base > 3].cell
frame_file[frame_file.cv_stim > 3].cv_stim
#
frame_file.groupby('cell').count()
frame_file.groupby('cell').groups.keys()
frame_file.group_by('cell')
len(frame_file.cell.unique())
##############################################
# modulatoin comparison for both cell_types
################################
# Modulation, cell type comparison
# todo: hier die diff werte über die zellen
#ax_here = []
#axd = plt.subplot(grid_lower_lower[0, c])
# embed()
#kernel_histogram(axk, colors[str(cell_type_here)], np.array(x_axis), norm=True, step=0.03, alpha=0.5)
#embed()
#axk.show_spines('lb')
#axs = plt.subplot(grid0[6+c])
colorbar = False
#if colorbar:
x_axis = ['cv_base','cv_stim','response_modulation']#,'fr_base']#
var_item_names = [var_it,var_it,var_it2]#,var_it2]#['Response Modulation [Hz]',]
var_types = ['response_modulation','response_modulation','']#,'']#'response_modulation'
max_x = max_xs[c]
x_axis_names = ['CV$_{Base}$','CV$_{stim}$','Response Modulation [Hz]']#$_{Base}$,'Fr$_{Base}$',]
#score = scores[0]
score_n = ['Perc99/Med', 'Perc99/Med', 'Perc99/Med']
score = scores[c]
scores_here = [score,score,score]#,score]
score_name = ['max(diag5Hz)/med_diagonal_proj_fr','max(diag5Hz)/med_diagonal_proj_fr']#,'max(diag5Hz)/med_diagonal_proj_fr']#'Perc99/Med'
score_name = ['Fr/Med', 'Fr/Med']#'Fr/Med'] # 'Perc99/Med'
score_name = [NLI_scorename2(), NLI_scorename2(), NLI_scorename2()]#NLI_scorename()] # 'Fr/Med''Perc99/Med'
ax_j = []
axls = []
axss = []
#embed()
#frame_max = frame_file[frame_file[score]>5]
cv_name = "cv_base"
max_val = 1.5
log = ''#'logall'#''#'logy','logall'True#False
marker = ['o']
for v, var_type in enumerate(var_types):
# ax = plt.subplot(grid0[1+v])#grid_lower[0, v]
axx, axy, axs, axls, axss, ax_j = get_grid_4(ax_j, axls, axss, grid0[v,counter])
if log == 'logy':
ymin = 'no'
else:
ymin = 0
xmin = 0
if (' P-unit' in cell_type_here) & ('cv' in x_axis[v]):
xlimk = [0, 2]
else:
xlimk = None
xlimk = None
labelpad = 0.5#-1
cmap, _, y_axis = plt_burst_modulation_hists(axx, axy, var_item_names[v], axs, cell_type_here, x_axis[v],
frame_file, scores_here[v], ymin=ymin, xmin=xmin,
burst_fraction_reset=burst_corr_reset, var_item=var_type,
max_x=max_x[v], xlim=xlimk, x_pos=1, labelpad = labelpad,
burst_fraction=burst_fraction[c], ha='right')
print(cell_type_here + ' median '+scores_here[v]+''+str(np.nanmedian(frame_file[scores_here[v]])))
print(cell_type_here + ' max ' + x_axis[v] + '' + str(np.nanmax(frame_file[x_axis[v]])))
if v == 0:
colors = colors_overview()
axx.set_title(cell_types_name[c], color = colors[cell_type_here])
axx.show_spines('')
axy.show_spines('')
axs.set_ylabel(score_name[v])
axs.set_xlabel(x_axis_names[v], labelpad = labelpad)
extra_lim = False
if extra_lim:
if (' P-unit' in cell_type_here) & ('cv' in x_axis[v]):
axs.set_xlim(xlimk)
axx.set_xlim(xlimk)
#embed()
#remove_yticks(axl)
if log == 'logy':
axy.set_yscale('log')
axs.set_yscale('log')
make_log_ticks([axs])
axy.minorticks_off()
elif log == 'logall':
axy.set_yscale('log')
axs.set_yscale('log')
make_log_ticks([axs])
axy.minorticks_off()
axx.set_xscale('log')
axs.set_xscale('log')
make_log_ticks([axs])
axx.minorticks_off()
axy.set_yticks_blank()
plt_specific_cells(axs, cell_type_here, x_axis[v], frame_file, scores_here[v], marker = ['o',"s"])
tags.append(axx)
counter += 1
#plt.show()
############################
# 1 Print scores
print('\n')
speciess = [' Apteronotus leptorhynchus', ' Eigenmannia virescens']
#for species in speciess:
for c, cell_type_here in enumerate(cell_types):
frame_file = setting_overview_score(frame_load_sp, cell_type_here, min_amp='min', species=species)
print(cell_type_here + str(len(frame_file.cell.unique())))
#embed()
cells_amp = ['2011-09-21-ab', '2010-06-21-am', '2012-05-15-ac', '2012-04-26-ae', '2012-05-07-ac',
'2010-06-21-ac']
cells_amp = update_cell_names(cells_amp)
frame_amps = frame_file[frame_file.cell.isin(cells_amp)]
frame_amps.cell.unique()
#['2012-05-07-ac-invivo-1', '2012-04-26-ae','2012-05-15-ac-invivo-1', '2010-06-21-ac-invivo-1','2010-06-21-am-invivo-1', '2011-09-21-ab-invivo-1']
#embed()
#if len(frame_amps)
#################################################################
# print the corrs
score = scores[c]#'ser_base',
score_print = [score]
score_print = [scores[c],scores[c],scores[c],scores[c],scores[c],'burst_fraction_burst_corr_individual_base']
corr_vars = ['cv_base','response_modulation','fr_base','ser_first_base','burst_fraction_burst_corr_individual_stim','cv_base']
#for score in score_print:
for c_nr, corr_var in enumerate(corr_vars):#,'ser_sum_corr'
score = score_print[c_nr]
x = frame_file['fr_base']
y = frame_file['ser_first']
c = frame_file['cv_base']
try:
c_axis, x_axis, y_axis, exclude_here = exclude_nans_for_corr(frame_file, corr_var, cv_name=corr_var , score=score)
except:
print('frame file lost')
embed()
corr, p_value = stats.pearsonr(x_axis, y_axis)
pears_l = pearson_label(corr, p_value, y_axis, n=True)
print(start_name(cell_type_here, species) + ' ' + corr_var + ' to '+str(score) + pears_l)
#if corr_var == 'ser_first_base':#
# embed()
print('\n')
################################
# cv to fr correlation
c_axis, x_axis, y_axis, exclude_here = exclude_nans_for_corr(frame_file, 'fr_base', cv_name='fr_base',
score='cv_base')
corr, p_value = stats.pearsonr(x_axis, y_axis)
pears_l = pearson_label(corr, p_value, y_axis, n=True)
print(start_name(cell_type_here, species) + ' ' + 'fr_base' + ' to ' + 'cv_base' + pears_l)
################################
# cv to fr correlation
frame_file = setting_overview_score(frame_load_sp, cell_type_here, min_amp='range', species=species)
c_axis, x_axis, y_axis, exclude_here = exclude_nans_for_corr(frame_file, 'burst_fraction_burst_corr_individual_base', cv_name='burst_fraction_burst_corr_individual_base',
score='cv_base')
corr, p_value = stats.pearsonr(x_axis, y_axis)
pears_l = pearson_label(corr, p_value, y_axis, n=True)
print(start_name(cell_type_here, species) + ' amprange: ' + 'burst_fraction_individual_base' + ' to ' + 'cv_base' + pears_l)
###############################
# fr to nonline but for both modulations
c_axis, x_axis, y_axis, exclude_here = exclude_nans_for_corr(frame_file, 'fr_base', cv_name='fr_base',
score=score)
corr, p_value = stats.pearsonr(x_axis, y_axis)
pears_l = pearson_label(corr, p_value, y_axis, n=True)
print(start_name(cell_type_here, species) + ' amprange: ' + ' ' + 'fr_base' + ' to score' + pears_l)
##################################
print('\n')
#embed()
#todo: hier noch die Werte für die Methodik printen
test = False
#embed()
if test:
plt.hist(frame_file['lim_individual'])
#print(frame_file['lim_individual'])
############################################
############################
# 2 Print names
# printe um welchen Zeitraum es sich handelt
frame_file = setting_overview_score(frame_load_sp, cell_type_here = '', f_exclude=False, snippet=None, min_amp='min',
species='')
print('\n')
print('From ' +str(np.sort(frame_file.cell)[0])+' to '+str(np.sort(frame_file.cell)[-1])+' nr '+str(len(frame_file.cell)))
#embed()
#embed()
#####################################################################
############################
# 3 Print , EODF
print('\n')
speciess = [' Apteronotus leptorhynchus', ' Eigenmannia virescens']
for species in speciess:
fish = []
for c, cell_type_here in enumerate(cell_types):
frame_file = setting_overview_score(frame_load_sp, cell_type_here, min_amp='min', species=species)
if len(frame_file)> 0:
# ja der fisch hatte halt eine ziemlich niedriege EODf
frame_here = frame_file[frame_file.cell != '2022-01-05-af-invivo-1']
try:
print(start_name(cell_type_here,
species) + ' eodf_min ' + str(int(np.round(np.nanmin(frame_here.eod_fr)))) + ' eodf_max ' + str(int(np.round(np.nanmax(frame_here.eod_fr))))+' n '+str(len(frame_here)))
except:
print('eodf thing')
embed()
my_list = np.unique(frame_here.cell)
new_list = [item[0:10] for item in my_list]
fish.extend(new_list)
print(species[0:7]+' n_fish '+str(len(np.unique(fish))))
#####################################################################
# Burst corr limits
print('\n')
for species in speciess:
fish = []
for c, cell_type_here in enumerate(cell_types):
frame_file = setting_overview_score(frame_load_sp, cell_type_here, min_amp='min', species=species)
lims = np.unique(frame_file['lim_individual'])
lims_name = ''
for lim in lims:
lims_name = lims_name + str(lim) + ': ' + str(np.sum(frame_file['lim_individual'] == lim)) + ','
print(start_name(cell_type_here, species) + ' ' + lims_name)
#####################################################################
############################
# 3 Print , EODF
# Amplituden
print('\n')
speciess = [' Apteronotus leptorhynchus', ' Eigenmannia virescens']
for species in speciess:
for c, cell_type_here in enumerate(cell_types):
frame_file = setting_overview_score(frame_load_sp, cell_type_here,min_amp = '',species=species)
if len(frame_file)> 0:
#embed()
print(start_name(cell_type_here, species) + ' amp_min ' + str(
np.nanmin(frame_file.amp)) + ' amp_max ' + str(np.nanmax(frame_file[~np.isinf(frame_file.amp)].amp)))
test = False
#embed()
if test:
frame_file[['cv_stim','cv_base']]
frame_file['cv_stim']
embed()
############################
# CVs min max
for species in speciess:
for c, cell_type_here in enumerate(cell_types):
frame_file = setting_overview_score(frame_load_sp, cell_type_here, min_amp='', species=species)
if len(frame_file) > 0:
try:
print(start_name(cell_type_here, species) + ' CV_min ' + str(np.round(np.nanmin(frame_file.cv_base), 2)) + ' CV max ' + str(np.round(np.nanmax(frame_file.cv_base), 2)) + ' FR max ' + str(np.round(np.nanmin(frame_file.fr_base))) + ' FR max ' + str(np.round(np.nanmax(frame_file.fr_base))))
except:
print('min here')
embed()
############################################
#######
# N
print('\n N')
speciess = [' Apteronotus leptorhynchus', ' Eigenmannia virescens']
for species in speciess:
for c, cell_type_here in enumerate(cell_types):
##################
# printe wie viele Zellen es am Anfang gab
frame_file = setting_overview_score(frame_load_sp, cell_type_here, f_exclude = False, snippet = None, min_amp='min', species=species)
print(start_name(cell_type_here, species) + ' nr ' + str(len(frame_file)))
#embed()
#if c in [0,2]:
#cmap, _, y_axis = plt_modulation_overview(axs, c, cell_type_here,
# cv_name, frame_file, max_val, score,
# species)
#axs.set_ylabel(score)
#embed()#frame_file[(frame_file.cv_base < 0.65) & (frame_file.response_modulation > 200)].response_modulation
#axs.set_xlabel(cv_name)
#axs.get_shared_x_axes().join(*[axs, axd])
# elif species == ' Apteronotus albifrons':
# plt_albi(ax[4, 1], cell_type_here, colors, max_val, species, x_axis, y_axis)
#ax[1,cv_n].set_xlim(0, max_val)
#set_same_ylim(np.concatenate(ax[1::, :]))
#set_same_ylim(np.concatenate(ax[1::, :]),ylim_type ='xlim')
#set_same_ylim(ax[0, :], ylim_type='xlim')
#set_ylim_same()
#ax[1, 1].get_shared_y_axes().join(*ax[1, 1::])
# embed()
#counter += 1
#embed()
############################################
# printe um welchen Zeitraum es sich handelt
frame_file = setting_overview_score(frame_load_sp, cell_type_here = '', f_exclude=False, snippet=None, min_amp='min',
species='')
print('\n sampling'+str(frame_file.sampling.unique()))
print('P-units Fr: 50-450 Hz CV: 0.15 - 1.35')
print('Ampullary Fr: 80-200 Hz, CV: 0.08 - 0.22')
#plt.show()
show = False#True
fig = plt.gcf()
#embed()
tags_final = np.concatenate([tags[0::3],tags[1::3],tags[2::3]])#,tags[2::3]
fig.tag(tags_final, xoffs=-4.2, yoffs=1.12)
save_visualization(pdf = True, show = show)
def start_name(cell_type_here, species):
return species[0:7] + ' ' + cell_type_here[0:7]
def plt_specific_cells(axs, cell_type_here, cv_name, frame_file, score, marker = []):
######################################################
# hier kommen die kontrast Punkte dazu
# für die Zellen spielt Burst correctin ja keine Rolle
# if cv_n == 0:
if cell_type_here == ' P-unit':
cells_plot2 = p_units_to_show(type_here='contrasts')[0:2]
else:
cells_plot2 = [p_units_to_show(type_here='amp')[0]]
# for cell_plt in cells_plot2:
cells_extra = frame_file[frame_file['cell'].isin(cells_plot2)].index
# ax = plt.subplot(grid[1, cv_n])
# todo: hier nur noch die kleinste und größte Amplitude nehmen
# embed()
if not marker:
axs.scatter(frame_file[cv_name].loc[cells_extra], frame_file[score].loc[cells_extra],
s=9, facecolor="None", edgecolor='black', alpha=0.7, clip_on=False) # colors[str(cell_type_here)]
else:
#embed()
axs.scatter(frame_file[cv_name].loc[cells_extra][0:2], frame_file[score].loc[cells_extra][0:2],
s=9, facecolor="None", marker = marker[1], edgecolor='black', alpha=0.7, clip_on=False) # colors[str(cell_type_here)]
axs.scatter(frame_file[cv_name].loc[cells_extra][2:4], frame_file[score].loc[cells_extra][2:4],
s=9, facecolor="None", marker = marker[0], edgecolor='black', alpha=0.7, clip_on=False) # colors[str(cell_type_here)]
def plt_var_axis(ax_j, axls, axss,score_name, burst_fraction, cell_type_here, counter, cv_name, frame_file, grid0, max_val, score,
scores_here, var_item_names, var_types, x_axis, x_axis_names, log = False):
for v, var_type in enumerate(var_types):
# ax = plt.subplot(grid0[1+v])#grid_lower[0, v]
axk, axl, axs, axls, axss, ax_j = get_grid_4(ax_j, axls, axss, grid0[counter])
counter += 1
cmap, _, y_axis = plt_burst_modulation_hists(axk, axl, var_item_names[v], axs, cell_type_here, x_axis[v],
frame_file, scores_here[v], var_item=var_type,
burst_fraction=burst_fraction[v])
axs.set_ylabel(score_name[v])
axs.set_xlabel(x_axis_names[v])
if v in [0, 1]:
if log:
axl.set_yscale('log')
axs.set_yscale('log')
axl.set_yticks_blank()
# remove_yticks(axl)
axl.minorticks_off()
if v == 0:
############################
# extra Zellen Scatter
# todo: diese Zellen müssen noch runter konvertiert werden
# todo: extra funktion für Zellen über 9 Snippets schreiben und die nochmal extra machen
cells_plot2 = p_units_to_show(type_here='bursts')
# for cell_plt in cells_plot2:
cells_extra = frame_file[frame_file['cell'].isin(cells_plot2)].index
# ax = plt.subplot(grid[1, cv_n])
axs.scatter(frame_file[x_axis[v]].loc[cells_extra], frame_file[score].loc[cells_extra],
s=5, color='white', edgecolor='black', alpha=0.5,
clip_on=False) # colors[str(cell_type_here)]
return ax_j, axls, axs, axss, counter
def species_with_both_cells(grid0, cell_types,ax_j, axls, axss, colors, cv_name_title, save_names, scores, species_all, x_axis, log = True):
for cv_n, cv_name in enumerate(x_axis):
if cv_n == 0:
redo = False
else:
redo = False
redo = False
frame_load_sp = load_overview_susept(save_names[0], redo=redo, redo_class=redo)
# frame_file = setting_overview_score(cell_type_here, frame_load_sp, min_amp=True, species=species)
# print(np.isnan(species))
score = scores[0]
max_val = 1.5
for c, cell_type_here in enumerate(cell_types):
#if c == 1:
# embed()
species = species_all[cv_n]
frame_file = setting_overview_score(frame_load_sp, cell_type_here, min_amp='min', species=species)
##################################
# modulation and species comparison
# x_axis, y_axis = get_axis(cv_name, frame_file, score)
# if cv_n == 0:
###############################
#######################
# Kernel Histogram
# plot the histograms of the values above the according vals
# grid = gridspec.GridSpecFromSubplotSpec(1, 1, grid0[0],hspace=0, wspace = 0.15)#grid[0, cv_n]
#
if c == 0:
axk, axl, axs, axls, axss, ax_j = get_grid_4(ax_j, axls, axss, grid0)
if c in [0, 2]:
axk.set_title(species)
# embed()
if len(frame_file) > 0:
axs, x_axis = kernel_scatter(axl, cell_types, axk, axs, ax_j, c, cell_type_here, colors, cv_n, cv_name,
frame_file, grid0, max_val,
score, log = log)
if log:
axl.set_yscale('log')
axl.set_yticks_blank()
axl.minorticks_off()
axs.set_xlabel(cv_name_title[cv_n])
if cv_n == 0:
axs.set_ylabel('Perc(99)/Median')
if cv_n == 0:
axm = [axs]
return ax_j, axls, axss, cell_types, frame_load_sp, max_val, score
def get_grid_4(ax_j, axls, axss, grid0):
grid_k = gridspec.GridSpecFromSubplotSpec(2, 2, grid0,
hspace=0.1, wspace=0.1, height_ratios=[0.35, 3], width_ratios=[3, 0.5])
axk = plt.subplot(grid_k[0, 0])
ax_j.append(axk)
axs = plt.subplot(grid_k[1, 0])
axss.append(axs)
axl = plt.subplot(grid_k[1, 1])
axls.append(axl)
return axk, axl, axs, axls, axss, ax_j
if __name__ == '__main__':
data_overview3()