123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225 |
- import numpy as np
- import pandas as pd
- import mrestimator as mre
- import os
- import time
- import mr_utilities
- import math
- def analyse_timeresolved_mr(data_path, result_path, recordings,
- binnings='all', ksteps=(1,400),
- windowsize=20000, windowstep=500,
- fit_method="offset"):
- for recording in recordings:
- outputfile = '{}{}/pkl/activity_SU_4ms.pkl'.format(data_path, recording)
- if os.path.isfile(outputfile):
- data = pd.read_pickle(outputfile)
- else:
- print('Error: no datafile for recording {}'.format(recording))
- continue
- existing_binnings = data['binning'].unique().tolist()
- patient = int(data['patient'].unique()[0])
- if binnings == 'all':
- binnings = existing_binnings
- for binning in binnings:
- dt = data['dt'].unique()[0]
- activity_list = data[data['binning'] == binning]['activity'].tolist()
- if len(activity_list) == 0:
- continue
- if len(activity_list) > 1:
- print('Error: Multiple activity entries for binning.')
- data = data[~data['activity'].apply(tuple).duplicated()]
- activity_list = data[data['binning'] == binning][
- 'activity'].tolist()
- activity = activity_list[0]
- else:
- activity = activity_list[0]
- t_onset = 600 # in seconds, given 10min recordings
- # the recording is significantly less than 10min long,
- # it is unclear whether seizure onset is indeed at the end of
- # the recording, therefore better skip
- if len(activity) < 0.8 * t_onset * 10**3 / dt:
- continue
- save_path_recording = '{}{}/'.format(result_path, recording)
- timeresolved_mr(activity, windowsize, windowstep,
- ksteps[0], ksteps[1], fit_method,
- save_path_recording, patient, recording,
- binning, dt)
- print('break to cool cpu...')
- time.sleep(30)
- def timeresolved_mr(activity, windowsize, windowstep, kmin, kmax,
- fit_method, save_path_recording, patient, recording,
- binning, dt):
- """ performs MR estimation in sliding windows of provided activity and
- calls function to store results in pkl file
- """
- for t in range(0, len(activity) - windowsize, windowstep):
- windowed_data = activity[t:t + windowsize]
- prepared = mre.input_handler(windowed_data)
- rk = mre.coefficients(prepared, dt=dt, dtunit='ms', steps=(kmin, kmax))
- try:
- fit_res = mre.fit(rk, fitfunc=mr_utilities.string_to_func(
- fit_method))
- t0 = t * dt
- mt_savefile(save_path_recording, patient, recording, binning,
- dt, kmin, kmax, fit_method, windowsize, windowstep,
- t0=t0, activity=windowed_data, rk=rk, fitres=fit_res)
- except (RuntimeError, ValueError):
- continue
- def mt_savefile(result_path, patient, recording, binning, dt, kmin, kmax,
- fit_method, windowsize, windowstep, t0, activity, rk,
- fitres):
- """ save results of m(t) estimation performed by function
- timeresolved_mr to pkl file.
- """
- datafile = '{}mt_results_{}_{}_kmin{}_kmax{}_winsize{}_winstep{}_dt{}' \
- '.pkl'.format(result_path, binning, fit_method, kmin, kmax,
- windowsize, windowstep, dt)
- if math.isnan(fitres.tau):
- inconsistencies = ['not_converged']
- else:
- inconsistencies = mr_utilities.consistency_checks(activity, rk, fitres)
- try:
- fano = np.var(activity)/np.mean(activity)
- variance = np.var(activity)
- except ZeroDivisionError:
- fano = np.nan
- variance = np.nan
- save_res = dict()
- save_res['patient'] = patient
- save_res['recording'] = recording
- save_res['binning'] = binning
- save_res['kmin'] = kmin
- save_res['kmax'] = kmax
- save_res['windowsize'] = windowsize
- save_res['windowstep'] = windowstep
- save_res['t0'] = t0
- save_res['rk'] = [list(rk.coefficients)]
- save_res['rk_steps'] = [list(rk.steps)]
- save_res['dt'] = dt
- save_res['fit_method'] = fit_method
- save_res['fit_params'] = [list(fitres.popt)]
- save_res['m'] = fitres.mre
- save_res['tau'] = fitres.tau
- save_res['inconsistencies'] = inconsistencies
- save_res['fano'] = fano
- save_res['variance'] = variance
- save_res['offset'] = fitres.popt[2]
- if os.path.isfile(datafile):
- mr_results_patient = pd.read_pickle(datafile)
- mr_results_patient = mr_results_patient.append(save_res,
- ignore_index=True)
- else:
- mr_results_patient = pd.DataFrame.from_records([save_res])
- mr_results_patient.to_pickle(datafile)
- def mt_x_parts(x, outputfilename, data_path, result_path, recordings,
- binnings='all', ksteps=(1, 400), fit_method ="offset"):
- mx_allresults = None
- for recording in recordings:
- outputfile = '{}{}/{}'.format(data_path, recording, outputfilename)
- if os.path.isfile(outputfile):
- data = pd.read_pickle(outputfile)
- else:
- print(' File not found {}'.format(outputfile))
- continue
- existing_binnings = data['binning'].unique().tolist()
- patient = int(data['patient'].unique()[0])
- kmin = ksteps[0]
- kmax = ksteps[1]
- if binnings == 'all':
- binnings = existing_binnings
- for binning in binnings:
- dt = data['dt'].unique()[0]
- activity_list = data[data['binning'] == binning]['activity'].tolist()
- activity = activity_list[0]
- partlength = int(np.floor(len(activity)/x))
- tau = []
- m = []
- binning_focus = mr_utilities.binning_focus(binning, patient)
- for i in range(x):
- activity_part = activity[i*partlength:(i+1)*partlength]
- nspikes = np.sum(np.array(activity_part))
- n_nonzero = mr_utilities.count_nonzero(activity_part)
- prepared = mre.input_handler(activity_part)
- rk = mre.coefficients(prepared, dt=dt, dtunit='ms',
- steps=(kmin, kmax))
- try:
- fit_res = mre.fit(rk,
- fitfunc=mr_utilities.string_to_func(
- fit_method))
- if not math.isnan(fit_res.tau):
- tau.append(fit_res.tau)
- m.append(fit_res.mre)
- inconsistencies = mr_utilities.consistency_checks(
- activity_part, rk, fit_res)
- fano = np.var(activity_part)/np.mean(activity_part)
- variance = np.var(activity_part)
- mx_result = {'patient': patient,
- 'recording': recording,
- 'binning': binning,
- 'binning_focus': binning_focus,
- 'partno': i+1,
- 'm': fit_res.mre,
- 'tau': fit_res.tau,
- 'fano': fano,
- 'variance': variance,
- 'offset': fit_res.popt[2],
- 'rk': [list(rk.coefficients)],
- 'fit_params': [list(fit_res.popt)],
- 'nspikes': nspikes,
- 'n_nonzero': n_nonzero,
- 'inconsistencies': inconsistencies,
- 'kmin': kmin,
- 'kmax': kmax,
- 'fit_method': fit_method}
- if mx_allresults is None:
- mx_allresults = pd.DataFrame.from_records(
- [mx_result])
- else:
- mx_allresults = mx_allresults.append(
- mx_result, ignore_index=True)
- except (RuntimeError, ValueError):
- continue
- datafile_mx = '{}mt_{}_results_allrecords_kmin{}_kmax{}_{}.pkl'.format(
- result_path, x, ksteps[0], ksteps[1], fit_method)
- if mx_allresults is not None:
- mx_allresults.to_pickle(datafile_mx)
|