susceptibility1/model_and_data.py
2024-06-13 13:33:40 +02:00

369 lines
20 KiB
Python

##from update_project import **#model_and_data
import os
import numpy as np
import pandas as pd
from IPython import embed
from matplotlib import gridspec, pyplot as plt
from threefish.plot_subplots import plt_model_flowcharts
from threefish.values import ypos_x_modelanddata
try:
from plotstyle import plot_style, spines_params
except:
print('plotstyle not installed')
from threefish.RAM.plot_subplots import plt_data_susept, plt_single_square_modl
from threefish.RAM.values import overlap_cells, perc_model_full
from threefish.load import resave_small_files, save_visualization
from threefish.RAM.reformat_matrix import load_model_susept
from threefish.core import find_folder_name
from threefish.RAM.plot_labels import label_noise_name, nonlin_title, remove_xticks, remove_yticks, set_xlabel_arrow, \
set_ylabel_arrow, title_find_cell_add, xlabel_xpos_y_modelanddata
import itertools as it
from threefish.defaults import default_figsize, default_settings
from threefish.plot.limits import set_clim_same, set_same_ylim
#from plt_RAM import model_and_data, model_and_data_sheme, model_and_data_vertical2
def table_printen(table):
print(table.keys())
for l in range(len(table)):
list_here = np.array(table.iloc[l])
l1 = "& ".join(list_here)
print(l1)
def model_and_data2(eod_metrice=False, width=0.005, nffts=['whole'], powers=[1], cells=["2013-01-08-aa-invivo-1"],
contrasts=[0], noises_added=[''], D_extraction_method=['additiv_cv_adapt_factor_scaled'],
internal_noise=['RAM'], external_noise=['RAM'], level_extraction=[''], receiver_contrast=[1],
dendrids=[''], ref_types=[''], adapt_types=[''], c_noises=[0.1], c_signal=[0.9], cut_offs1=[300]): # ['eRAM']
stimulus_length = 1 # 20#550 # 30 # 15#45#0.5#1.5 15 45 100
trials_nrs = [1] # [100, 500, 1000, 3000, 10000, 100000, 1000000] # 500
rep = 1000000 # 500000#0
good_data, remaining = overlap_cells()
cells_all = [good_data[0]]
plot_style()
default_figsize(column=2, length=3.1) #2.9.254.75 0.75# bottom=0.07, top=0.94,
grid = gridspec.GridSpec(2, 5, wspace=0.95, bottom=0.13, hspace=0.60, top=0.88,
width_ratios=[2, 0, 2, 2, 2], left=0.09, right=0.93, ) #bottom=0.09, hspace=0.25, top=0.9,
a = 0
maxs = []
mins = []
mats = []
ims = []
iternames = [D_extraction_method, external_noise,
internal_noise, powers, nffts, dendrids, cut_offs1, trials_nrs, c_signal,
c_noises,
ref_types, adapt_types, noises_added, level_extraction, receiver_contrast, contrasts, ]
nr = '2'
# cell_contrasts = ["2013-01-08-aa-invivo-1"]
# cells_triangl_contrast = np.concatenate([cells_all,cell_contrasts])
# cells_triangl_contrast = 1
# cell_contrasts = 1
rows = len(cells_all) # len(good_data)+len(cell_contrasts)
perc = 'perc'
lp = 2
label_model = r'Nonlinearity $\frac{1}{S}$'
for all in it.product(*iternames):
var_type, stim_type_afe, stim_type_noise, power, nfft, dendrid, cut_off1, trial_nrs, c_sig, c_noise, ref_type, adapt_type, noise_added, extract, a_fr, a_fe = all
fig = plt.figure()
hs = 0.45
#################################
# data cells
grid_data = gridspec.GridSpecFromSubplotSpec(1, 1, grid[0, 0],
hspace=hs)
fr_print = False
nr = 1
ax_data, stack_spikes_all, eod_frs = plt_data_susept(fig, grid_data, cells_all, cell_type='p-unit', width=width,
cbar_label=True, fr_print=fr_print,
eod_metrice=eod_metrice, nr=nr, amp_given=1, xlabel=False,
lp=lp, title=True)
for ax_external in ax_data:
ax_external.set_xticks_delta(100)
set_ylabel_arrow(ax_external, xpos=xlabel_xpos_y_modelanddata(), ypos=0.87)
set_xlabel_arrow(ax_external, ypos=ypos_x_modelanddata())
#embed()
#plt.show()
##################################
# model part
trial_nr = 100000
cell = '2013-01-08-aa-invivo-1'
cell = '2012-07-03-ak-invivo-1'
cells_given = [cell]
save_name_rev = find_folder_name(
'calc_model') + '/' + 'calc_RAM_model-2__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_' + str(
trial_nr) + '_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV_revQuadrant_'
# for trial in trials:#.009
trial_nr = 1000000 #1000000
save_names = [
'calc_RAM_model-2__nfft_whole_power_1_afe_0.025_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-2__nfft_whole_power_1_afe_0.025_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_' + str(
trial_nr) + '_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-2__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-2__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_' + str(
trial_nr) + '_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
]
##########
# Erklärung
# Ich habe hier 0.009 und nicht 0.25 weil das Modell einen Fehler hat
# den Stimulus in den Daten habe ich überprüft der tatsächliche stimulus ist 2.3 Prozent
# sollte aber 2.5 Prozent sein
# Im Fall von 2.5 Prozent wäre das ein Fehler von 0.36 sonst von 0.39
# Hier werde ich nun mit dem Fehler von 0.36 verfahren
# das bedeutet aber das sich den Stimulus zwar mit 0.009 ins Modell reintue später für die
# Susceptiblitätsberechnung sollte ich ihn aber um den Faktor 0.36 teilen
# oben habe ich einen bias factor weil die Zelle zu sensitiv gefittet ist, also passe ich das an dass die den
# gleichen CV und feurrate hat, wie die Zelle in der Stimulation, deswegen ist dieser Bias faktor nur oben!
#
#bias_factors = [1, 1, 1, 1] # 0.3
bias_factors = [0.36, 0.36, 1, 1] # 0.36
#new
save_names = [
'calc_RAM_model-3__nfft_whole_power_1_afe_0.025_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV_old_fit_',
'calc_RAM_model-2__nfft_whole_power_1_afe_0.009_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_' + str(
trial_nr) + '_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-2__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV_old_fit_',
'calc_RAM_model-2__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_' + str(
trial_nr) + '_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
]#'calc_RAM_model-2__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
save_names = [
'calc_RAM_model-3__nfft_whole_power_1_afe_0.025_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV_old_fit_',
'calc_RAM_model-3__nfft_whole_power_1_afe_0.025_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_'+str(trial_nr) + '_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV_old_fit_',
'calc_RAM_model-3__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV_old_fit_',
'calc_RAM_model-2__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_100000_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV_old_fit_',
]
save_names = [
'calc_RAM_model-2__nfft_whole_power_1_afe_0.023_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-2__nfft_whole_power_1_afe_0.023_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_' + str(
trial_nr) + '_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-2__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-2__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_' + str(
trial_nr) + '_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
] #calc_RAM_model-2__nfft_whole_power_1_afe_2.6_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV
save_names = [
'calc_RAM_model-3__nfft_whole_power_1_afe_0.025_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-3__nfft_whole_power_1_afe_0.025_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_'+str(trial_nr) + '_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-3__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-3__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_'+str(trial_nr) + '_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
]
save_names = [
'calc_RAM_model-2__nfft_whole_power_1_afe_0.023_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-2__nfft_whole_power_1_afe_0.023_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_' + str(
trial_nr) + '_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-2__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-2__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_' + str(
trial_nr) + '_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
] #calc_RAM_model-2__nfft_whole_power_1_afe_2.6_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV
bias_factors = [0.36, 0.36, 1, 1]#0.36
save_names = [
'calc_RAM_model-2__nfft_whole_power_1_afe_0.009_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-2__nfft_whole_power_1_afe_0.009_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_' + str(
trial_nr) + '_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-2__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-2__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_' + str(
trial_nr) + '_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
]
save_names = [
'calc_RAM_model-2__nfft_whole_power_1_afe_0.009_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-2__nfft_whole_power_1_afe_0.009_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_' + str(
trial_nr) + '_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-2__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
'calc_RAM_model-2__nfft_whole_power_1_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_' + str(
trial_nr) + '_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
]
#
# 'calc_RAM_model-2__nfft_whole_power_1_afe_0.009_RAM_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_' + str(
#trial_nr) + '_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
#'calc_RAM_model-2__nfft_whole_power_1_afe_0.009_RAM_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_11_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
#'calc_RAM_model-2__nfft_whole_power_1_afe_0.009_RAM_RAM_additiv_cv_adapt_factor_scaled_cNoise_0.1_cSig_0.9_cutoff1_300_cutoff2_300no_sinz_length1_TrialsStim_500000_a_fr_1__trans1s__TrialsNr_1_fft_o_forward_fft_i_forward_Hz_mV',
nrs_s = [3, 4, 8, 9] #, 10, 11
#embed()
tr_name = trial_nr / 1000000
if tr_name == 1:
tr_name = 1 #'$c=1\,\%$','$c=0\,\%$'
c = 2.5
cs = ['$c=%.1f$' % c + '$\,\%$', '$c=0\,\%$']
titles = ['Model\n$N=11$', 'Model\n' + '$N=10^6$',
'Model\,(' + label_noise_name().lower() + ')' + '\n' + '$N=11$',
'Model\,(' + label_noise_name().lower() + ')' + '\n' + '$N=10^6$'
] #%#%s$' % (tr_name) + '\,million'
#'Model\,('+noise_name().lower()+')' + '\n' + '$N=11$\n $c=1\,\%$',$N=%s $' % (tr_name) +'\,million'
# 'Model\,('+noise_name().lower()+')' + '\n' + '$N=%s$' % (tr_name) + '\,million\n $c=1\,\%$ '
ax_model = []
for s, sav_name in enumerate(save_names):
try:
ax_external = plt.subplot(grid[nrs_s[s]])
except:
print('vers something')
embed()
ax_model.append(ax_external)
save_name = find_folder_name('calc_model') + '/' + sav_name
cell_add, cells_save = title_find_cell_add(cells_given)
perc = 'perc'
path = save_name + '.pkl' # '../'+
# model = get_stack_one_quadrant(cell, cell_add, cells_save, path, save_name)
# full_matrix = create_full_matrix2(np.array(model), np.array(stack_rev))
# stack_final = get_axis_on_full_matrix(full_matrix, model)
# im = plt_RAM_perc(ax, perc, np.abs(model))
model = load_model_susept(path, cells_save, save_name.split(r'/')[-1] + cell_add)
#embed()
if len(model) > 0:
add_nonlin_title, cbar, fig, stack_plot, im = plt_single_square_modl(ax_external, cell, model, perc,
titles[s],
width, eod_metrice=eod_metrice,
titles_plot=True,
resize=True,
bias_factor=bias_factors[s],
fr_print=fr_print, nr=nr)
# if s in [1,3,5]:
#embed()
ims.append(im)
mats.append(stack_plot)
maxs.append(np.max(np.array(stack_plot)))
mins.append(np.min(np.array(stack_plot)))
col = 2
row = 2
ax_external.set_xticks_delta(100)
ax_external.set_yticks_delta(100)
# cbar[0].set_label(nonlin_title(add_nonlin_title)) # , labelpad=100
cbar.set_label(nonlin_title(' [' + add_nonlin_title), labelpad=lp) # rotation=270,
if s in np.arange(col - 1, 100, col): # | (s == 0)
remove_yticks(ax_external)
else:
set_ylabel_arrow(ax_external, xpos=xlabel_xpos_y_modelanddata(), ypos=0.87)
if s >= row * col - col:
set_xlabel_arrow(ax_external, ypos=ypos_x_modelanddata())
# ax.set_xlabel(F1_xlabel(), labelpad=20)
else:
remove_xticks(ax_external)
if len(cells) > 1:
a += 1
set_clim_same(ims, mats=mats, lim_type='up', nr_clim='perc', clims='', percnr=perc_model_full())#
#################################################
# Flowcharts
ax_ams, ax_external = plt_model_flowcharts(a_fr, ax_external, c, cs, grid, model, stimulus_length)
set_same_ylim(ax_ams, up='up')
axes = np.concatenate([ax_data, ax_model])
axes = [ax_ams[0], axes[1], axes[2], ax_ams[1], axes[3], axes[4], ] #ax_ams[2], axes[5], axes[6],
#axd1 = plt.subplot(grid[1, 1])
#axd2 = plt.subplot(grid[2, 1])
#ax_data.extend([,])
#axd1.show_spines('')
#axd2.show_spines('')
#embed()
#axes = [[ax_ams[0],ax_data[0],axes[2], axes[3]],[ax_ams[1],axd1,axes[4], axes[5]],[axd2,axd2, axes[6], axes[7]]]
fig.tag([ax_data], xoffs=-3, yoffs=1.6) # ax_ams[3],
fig.tag([axes[0:3]], xoffs=-3, yoffs=1.6) # ax_ams[3],
#fig.tag([[axes[4]]], xoffs=-3, yoffs=1.6, minor_index=0) # ax_ams[3],
fig.tag([axes[3:6]], xoffs=-3, yoffs=1.6) #, major_index = 1, minor_index = 2 ax_ams[3],
#fig.tag([[axes[7]]], xoffs=-3, yoffs=1.6, major_index=2,minor_index=0) # ax_ams[3],
#fig.tag([axes[8::]], xoffs=-3, yoffs=1.6, major_index=2, minor_index=2) # ax_ams[3],
#fig.tag([axes[7::]], xoffs=-3, yoffs=1.6) # ax_ams[3],
#fig.tag([ax_ams[0],ax_data[0],axes[2], axes[3]], xoffs=-3, yoffs=1.6)#ax_ams[3],
save_visualization(pdf=True)
if __name__ == '__main__':
model = resave_small_files("models_big_fit_d_right.csv", load_folder='calc_model_core')
cells = model.cell.unique()
params = {'cells': cells}
show = True
# if show == False:
# low CV: cells = ['2012-07-03-ak-invivo-1']
plot_style()
default_settings(lw=0.5, column=2, length=3.35) #8.5
redo = False
D_extraction_method = ['additiv_cv_adapt_factor_scaled']
# D_extraction_method = ['additiv_visual_d_4_scaled']
##########################
# hier printen wir die table Werte zum kopieren in den Text
path = 'print_table_suscept-model_params_suscept_table.csv'
if os.path.exists(path):
table = pd.read_csv(path)
table_printen(table)
path = 'print_table_all-model_params_suscept_table.csv'
if os.path.exists(path):
table = pd.read_csv()
print('model big')
table_printen(table)
#embed()
##########################
#embed()
#print('hi')
model_and_data2(eod_metrice=False, width=0.005, D_extraction_method=D_extraction_method) #r'$\frac{1}{mV^2S}$'