<jupyter_start><jupyter_text>Subtract baseline/leak<jupyter_code>plt.figure(figsize=(10, 3))

plt.plot(gating[0][1, :, 0])
plt.axvspan(2300, 3300, alpha=0.4)
plt.axvspan(5500, 5900, alpha=0.4)
plt.xlim(0, 60000)
gating_tmp = []
for rec in gating:
    baseline_subtracted = subtract_baseline(rec, slice(2300, 3300), 0)
    leak_subtracted = subtract_leak(baseline_subtracted, slice(2300, 3300), slice(5500, 5900))
    gating_tmp.append(leak_subtracted)
gating = gating_tmp
del gating_tmp, baseline_subtracted, leak_subtracted<jupyter_output><empty_output><jupyter_text>Fit decay<jupyter_code>xlims = (28010, 35000)

plt.figure(figsize=(10, 6))

plt.subplot(211)
plt.plot(gating[0][0, :, -1])
plt.xlim(xlims)

plt.subplot(212)
plt.plot(gating[0][1, :, -1])
plt.axvspan(2300, 3300, alpha=0.4)
plt.axvspan(5500, 5900, alpha=0.4)
plt.xlim(xlims)
step_amplitude = 50  # 70mV (-70 to -20)

IA_inactivation = {
    'traces': [],
    'fitted_data': [],
    'fitted_curves': [], 
    'range_fitted': [],
    'inactivation_taus': [],
    'peak_latencies': [],
    'peak_amplitudes': [],
    'steady_state_amplitudes': [],
    'peak_exists': []
}

peak_fit_params = {
    'stimulus_start_time': 2801.6,
    'steady_state_time': 3801.6,
    'peak_slice': slice(28036, 28500),
    'decay_slice_end': 29000
}

for i in range(len(gating)):
    
    peak_inds = find_peaks(
        gating[i][0, peak_fit_params['peak_slice'], -1],
        distance=500,
        prominence=50.,
        width=100
    )[0]
    if len(peak_inds) > 0:
        if len(peak_inds) > 1:
            warnings.warn(
                '{} peaks detected for recording {}.'.format(len(peak_inds), i)
            )
        peak_exists = True
        peak_ind = np.argmax(gating[i][0, peak_fit_params['peak_slice'], -1])
    else:
        peak_exists = False
    IA_inactivation['peak_exists'].append(peak_exists)
    
    # Find peak latency
    if peak_exists:
        IA_inactivation['peak_latencies'].append(
            peak_ind * 0.1 
            + peak_fit_params['peak_slice'].start * 0.1 
            - peak_fit_params['stimulus_start_time']
        )
    else:
        IA_inactivation['peak_latencies'].append(np.nan)
    
    # Get steady-state amplitude
    steady_state_amplitude = gating[i][0, timeToIndex(peak_fit_params['steady_state_time'], 0.1)[0], -1]
    IA_inactivation['steady_state_amplitudes'].append(steady_state_amplitude)
    
    # Get peak amplitude.
    if peak_exists:
        peak_amplitude = gating[i][0, peak_ind + peak_fit_params['peak_slice'].start, -1]
        IA_inactivation['peak_amplitudes'].append(peak_amplitude)
    else:
        IA_inactivation['peak_amplitudes'].append(np.nan)
    
    # Fit decay
    if peak_exists:
        
        # Convert decay_slice to eighty twenty range
        decay_slice = slice(peak_ind + peak_fit_params['peak_slice'].start, peak_fit_params['decay_slice_end'])
        decay_slice = getIndicesByPercentile(gating[i][0, decay_slice, -1], [0.80, 0.20]) + decay_slice.start  # Get 80-20 interval
        decay_slice[1] = np.argmin(gating[i][0, decay_slice[0]:decay_slice[1], -1]) + decay_slice[0]  # Truncate to min.

        decay_slice += peak_ind
        decay_slice = slice(decay_slice[0], decay_slice[1])
        plt.plot(gating[i][0, decay_slice, -1])
        plt.show()
        IA_inactivation['range_fitted'].append([decay_slice.start * 0.1, decay_slice.stop * 0.1])
        

        t_vec = np.arange(0, gating[i].shape[1], 0.1)[:gating[i].shape[1]]
        fitted_tmp = gating[i][0, decay_slice, -1]
        IA_inactivation['fitted_data'].append(np.array(
            [fitted_tmp,
            t_vec[decay_slice]]
        ))

        IA_inactivation['traces'].append(np.array([gating[i][0, :, -1], gating[i][1, :, -1], t_vec]))

        p_tmp, fitted_tmp = fit_decay_curve(
            fitted_tmp,
            [fitted_tmp[0] - fitted_tmp[-1], fitted_tmp[-1], 20],
            dt = 0.1
        )

        IA_inactivation['inactivation_taus'].append(p_tmp[2])
        IA_inactivation['fitted_curves'].append(np.array(
            [fitted_tmp[0, :],
            np.linspace(decay_slice.start * 0.1, decay_slice.stop * 0.1, len(fitted_tmp[0, :]))]
        ))

        # Diagnostic plot of decay fit.
        plt.figure()
        plt.axhline(peak_amplitude)
        plt.axvline((peak_ind + peak_fit_params['peak_slice'].start) * 0.1)
        plt.plot(t_vec, gating[i][0, :, -1], 'k-', lw=0.5)
        #plt.plot(t_vec[peak_fit_params['decay_slice']], gating[i][0, peak_fit_params['decay_slice'], -1], 'r-')
        plt.plot(
            np.linspace(
                decay_slice.start * 0.1, decay_slice.stop * 0.1, len(fitted_tmp[0, :])
            ),
            fitted_tmp[0, :], '--', color = 'gray', lw=4
        )
        plt.xlim(2600, 4000)
        plt.show()
    else:
        IA_inactivation['range_fitted'].append([np.nan, np.nan])
        IA_inactivation['traces'].append(np.array([gating[i][0, :, -1], gating[i][1, :, -1], gating[i].time_supp]))
        IA_inactivation['inactivation_taus'].append(np.nan)
        IA_inactivation['fitted_curves'].append(None)
                                         


for key in IA_inactivation:
    if key != 'fitted_data':
        IA_inactivation[key] = np.array(IA_inactivation[key])
print('IA inactivation tau {:.1f} +/- {:.1f} (mean +/- SD)'.format(
    IA_inactivation['inactivation_taus'].mean(), IA_inactivation['inactivation_taus'].std()
))

#with open(PROCESSED_PATH + 'inactivation_fits.dat', 'wb') as f:
#    pickle.dump(IA_inactivation, f)
print(np.sum(IA_inactivation['peak_exists']) / len(IA_inactivation['peak_exists']))
IA_inactivation.keys()
for dataset in ['peak_amplitudes', 'peak_latencies', 'steady_state_amplitudes', 'inactivation_taus']:
    print('{:>25} {:>10.3} +/- {:>5.3}'.format(
        dataset, 
        np.nanmean(IA_inactivation[dataset]), 
        stats.sem(IA_inactivation[dataset], nan_policy='omit')
    ))<jupyter_output>peak_amplitudes        nan +/-   0.0
           peak_latencies        nan +/-   0.0
  steady_state_amplitudes    3.8e+02 +/-  80.4
        inactivation_taus        nan +/-   0.0<jupyter_text>Summary statistics for quantities in nS<jupyter_code>print('Peak amplitudes (nS) {:>20.3} +/- {:>5.3}'.format(
    np.nanmean(IA_inactivation['peak_amplitudes'] / step_amplitude),
    stats.sem(IA_inactivation['peak_amplitudes'] / step_amplitude, nan_policy='omit')
))
print('Steady state amplitudes (nS) {:>20.3} +/- {:>5.3}'.format(
    np.nanmean(IA_inactivation['steady_state_amplitudes'] / step_amplitude),
    stats.sem(IA_inactivation['steady_state_amplitudes'] / step_amplitude, nan_policy='omit')
))
IA_inactivation['steady_state_conductance'] = IA_inactivation['steady_state_amplitudes'] / step_amplitude
IA_inactivation['peak_conductance'] = IA_inactivation['peak_amplitudes'] / step_amplitude
statistics = [
    'peak_latencies', 
    'steady_state_amplitudes', 
    'steady_state_conductance', 
    'peak_conductance', 
    'inactivation_taus',
    'peak_exists'
]
stats_df = pd.DataFrame({key: IA_inactivation[key] for key in statistics})
stats_df
stats_df.mean()
stats_df.sem()
stats_df.to_csv(os.path.join('data', 'processed', 'mPFC', 'transient_current_parameters.csv'), index=False)<jupyter_output><empty_output><jupyter_text>Figures<jupyter_code>fit_example_no = 8
bg_tr_alpha = 0.8

tr_spec = gs.GridSpec(2, 1, height_ratios=[1, 0.2], hspace=0, top=0.97, bottom=0.05, right=0.97, left=0.1)

plt.figure(figsize=(1.5, 1))

wc_ax = plt.subplot(tr_spec[0, :])
plt.plot(
    (IA_inactivation['traces'][~IA_inactivation['peak_exists'], 2, :].T - peak_fit_params['stimulus_start_time']), 
    IA_inactivation['traces'][~IA_inactivation['peak_exists'], 0, :].T,
    '-', color='gray', lw=0.5, alpha=bg_tr_alpha
)
plt.plot(
    (IA_inactivation['traces'][IA_inactivation['peak_exists'], 2, :].T - peak_fit_params['stimulus_start_time']), 
    IA_inactivation['traces'][IA_inactivation['peak_exists'], 0, :].T,
    '-', color=colors['pyr'], lw=0.5, alpha=bg_tr_alpha
)
plt.xlim(-50, 1550)
plt.ylim(-250, 1600)
pltools.add_scalebar(y_units='pA', y_size=500, omit_x=True, anchor=(-0.05, 0.1), linewidth=sbarlw)
sns.despine(ax=plt.gca(), trim=True)

"""wc_ins = inset_axes(wc_ax, '60%', '50%', loc='upper right', borderpad=1)
plt.plot(
    (IA_inactivation['traces'][~IA_inactivation['peak_exists'], 2, :].T - peak_fit_params['stimulus_start_time']), 
    IA_inactivation['traces'][~IA_inactivation['peak_exists'], 0, :].T,
    '-', color='gray', lw=0.5, alpha=bg_tr_alpha
)
plt.plot(
    (IA_inactivation['traces'][IA_inactivation['peak_exists'], 2, :].T - peak_fit_params['stimulus_start_time']), 
    IA_inactivation['traces'][IA_inactivation['peak_exists'], 0, :].T,
    '-', color=colors['pyr'], lw=0.5, alpha=bg_tr_alpha
)
plt.xlim(-20, 150)
plt.ylim(-50, 1500)
pltools.add_scalebar(x_units='ms', x_size=50, omit_y=True, anchor=(0.95, 0.7), x_label_space=0.05, remove_frame=False, linewidth=2)
plt.xticks([])
plt.yticks([])
dashedBorder(wc_ins, lw=1.5)
mark_inset(wc_ax, wc_ins, 1, 4, ls='--', color='gray', lw=1.5)"""

plt.subplot(tr_spec[1, :])
plt.plot(
    (IA_inactivation['traces'][:, 2, :].T - peak_fit_params['stimulus_start_time']), 
    IA_inactivation['traces'][:, 1, :].T,
    '-', color=colors['input'], lw=0.5, alpha=bg_tr_alpha
)
plt.xlim(-50, 1550)
pltools.add_scalebar(x_units='ms', x_size=200, omit_y=True, anchor=(0.7, -0.05), x_label_space=0.05, linewidth=sbarlw)

plt.tight_layout()

if IMG_PATH is not None:
    plt.savefig(os.path.join(IMG_PATH, 'mPFC_kinetics_trace_only.png'))
    plt.savefig(os.path.join(IMG_PATH, 'mPFC_kinetics_trace_only.svg'))
    
plt.show()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1696687370600414, 'fraction_numerical': 0.04006211180124224, 'mean_word_length': 3.864551863041289, 'pattern_counts': {'":': 0, '<': 16, '<?xml version=': 0, '>': 25, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '4', 'dir': 'jupyter-structured-clean-dedup', 'id': '634492', 'license': 'MIT', 'n_tokens_mistral': 4146, 'n_tokens_neox': 3709, 'n_words': 518, 'path': 'figs/scripts/mpfc_currents.ipynb', 'repo_name': 'nauralcodinglab/raphegif', 'size': '9660'}