<jupyter_start><jupyter_text>*CoNNear*: A convolutional neural-network model of human cochlear mechanics and filter tuning for real-time applicationsThis python notebook reproduces the evaluation results of the proposed CoNNear model. Please read through the prerequisites to install the correct packages.  Prerequisites- To run the trained CoNNear model, you should install numpy, scipy, keras and tensorflow, preferably through a Conda +3.6 environment. Open jupyter notebook to run this script.<jupyter_code>import os
#os.add_dll_directory(os.getcwd())

import numpy as np
from scipy import signal
import scipy.signal as sp_sig
import matplotlib.pyplot as plt
import keras
from keras.models import model_from_json
from keras.utils import CustomObjectScope
from keras.initializers import glorot_uniform

from tlmodel.get_tl_vbm_and_oae import tl_vbm_and_oae
from helper_ops import *

json_file = open("connear/Gmodel.json", "r")
loaded_model_json = json_file.read()
json_file.close()
with CustomObjectScope({'GlorotUniform': glorot_uniform()}):    
    connear = model_from_json(loaded_model_json)
with CustomObjectScope({'GlorotUniform': glorot_uniform()}):
     connear.load_weights("connear/Gmodel.h5")
connear.summary()<jupyter_output><empty_output><jupyter_text>Define parameters<jupyter_code># Define model specific variables
down_rate = 2
fs = 20e3
fs_tl = 100e3
factor_fs = int(fs_tl / fs)
p0 = 2e-5
right_context = 256
left_context = 256
# load CFs 
CF = np.loadtxt('tlmodel/cf.txt')
channels = CF.size<jupyter_output><empty_output><jupyter_text>Click responseShow the responses of the CoNNear model to a click stimulus.      **Notice that for all shown simulations the CoNNear model operates at 20kHz.**<jupyter_code># Define the click stimulus
dur = 128.0e-3 # for 2560 samples (2048 window length, 2x256 context)
stim = np.zeros((1, int(dur * fs)))
L = 70.0
samples = dur * fs
click_duration = 2 # 100 us click 
silence = 60 #samples in silence
samples = int(samples - right_context - left_context)

################ CoNNear ####################################
stim = np.zeros((1, int(dur * fs)))
stim[0, right_context + silence : right_context + silence + click_duration] = 2 * np.sqrt(2) * p0 * 10**(L/20)

stim = np.expand_dims(stim, axis=2)
connear_pred_click = connear.predict(stim.T, verbose=1)
bmm_click_connear = connear_pred_click[0,:,:].T * 1e-6<jupyter_output><empty_output><jupyter_text>Plotting the results.<jupyter_code>################ Plots ######################################
#Plot input stimulus
plt.plot(stim[0,256:-256]), plt.xlim(0,2000)
plt.show()

# Plot the CoNNear response
plt.imshow(bmm_click_connear, aspect='auto', cmap='jet')
plt.xlim(0,2000), plt.clim(-4e-7,5e-7)
plt.colorbar()
plt.show()<jupyter_output><empty_output><jupyter_text>Cochlear Excitation PatternsHere, we plot the simulated RMS levels of basilar memberane (BM) displacement across CF for tone stimuli presented at SPLs between 0 and 90 dB SPL.**You can change the `f_tone` variable to simulate tone stimuli of different frequencies, say 500Hz, 1kHz, 2kHz, etc..**<jupyter_code>f_tone = 1e3 # You can change this tone frequency to see how the excitation pattern changes with stimulus frequency
dur = 102.4e-3 # for 2048 samples
window_len = int(fs * dur)
L = np.arange(0., 91.0, 10.) # SPLs from 0 to 90dB
total_length = window_len + right_context + left_context #total length = 2560

################ CoNNear ####################################
t = np.arange(0., dur, 1./fs)
hanlength = int(10e-3 * fs) # 10ms length hanning window
stim_sin = np.sin(2 * np.pi * f_tone * t)
han = signal.windows.hann(hanlength)
stim_sin[:int(hanlength/2)] = stim_sin[:int(hanlength/2)] * han[:int(hanlength/2)]
stim_sin[-int(hanlength/2):] = stim_sin[-int(hanlength/2):] * han[int(hanlength/2):]
stim = np.zeros((len(L), total_length))
for j in range(len(L)):
    stim[j,right_context:window_len+right_context] = p0 * np.sqrt(2) * 10**(L[j]/20) * stim_sin
    
# prepare for feeding to the DNN
stim = np.expand_dims(stim, axis=2)

connear_pred_tone = connear.predict(stim, verbose=1)
bmm_tone_connear = connear_pred_tone

# Compute rms for each level
cochlear_pred_tone_rms = np.vstack([rms(bmm_tone_connear[i]) for i in range(len(L))])
################ Plots ######################################
# Plot the RMS for CoNNear
cftile=np.tile(CF, (len(L),1))
plt.semilogx((cftile.T), 20.*np.log10(cochlear_pred_tone_rms.T))
plt.xlim(0.25,8.), plt.grid(which='both'), 
plt.xticks(ticks=(0.25, 0.5, 1., 2., 4., 8.) , labels=(0.25, 0.5, 1., 2., 4., 8.))
plt.ylim(-80, 20)
plt.xlabel('CF (kHz)')
plt.ylabel('RMS of y_bm (dB)')
plt.title('CoNNear Predicted')
plt.show()<jupyter_output><empty_output><jupyter_text>QERB PlotsNext, the level dependent tuning properties of the cochlear filters (QERB) are shown. We chose click stimuli at three levels 0dB, 40dB and 70dB peSPL, computed the CoNNEar cochlear BM response and computed the QERB.<jupyter_code>L = [0., 40., 70.] # We will plot it for three SPLs
    
################ CoNNear ####################################
stim_con = np.zeros(((len(L),int(dur * fs)+right_context+left_context,1)))
QERB_connear = np.zeros(((len(L),channels)))

for i in range (len(L)):
    stim_con[i, right_context + silence : right_context + silence + click_duration, 0] = 2 * np.sqrt(2) * p0 * 10**(L[i]/20)
    
#Get CoNNear outputs
con_predicted = connear.predict(stim_con)

for i in range (len(L)):
    QERB_connear[i,:] = QERB_calculation(con_predicted[i, :, :].T, CF*1e3, fs)

################ Plots ######################################
# Plot QERB of CoNNear model
plt.semilogx(CF[0::5], (QERB_connear[0,0::5]),':gs', label='0dB')
plt.semilogx(CF[0::5], (QERB_connear[1,0::5]),'r', label='40dB')
plt.semilogx(CF[0::5], (QERB_connear[2,0::5]),':rs', label='70dB')
plt.xlim(0.25,8.), plt.grid(which='both'), 
plt.xticks(ticks=(0.25, 0.5, 1., 2., 4., 8.) , labels=(0.25, 0.5, 1., 2., 4., 8.))
plt.yticks(ticks=(5, 10, 15, 20) , labels=(5, 10, 15, 20))
plt.xlabel('CF (kHz)')
plt.ylim(2,20)
plt.ylabel('QERB')
plt.title('CoNNear Predicted')
plt.legend()
plt.show()<jupyter_output><empty_output><jupyter_text>Speech InputHere, a sentence from the Dutch speech matrix (unseen during training) will be input to both the TL and the CoNNear models. By adapting fragment_length parameters, various input lengths can be compared and visualised. **Notice that this part is computationally more expensive with a higher fragment_duration, both in terms of memory and time.**<jupyter_code>#load in speechfile
wavfile = 'dutch_sentence.wav'
signal_wav, fs_signal = wavfile_read(wavfile)
signalr = sp_sig.resample_poly(signal_wav, fs_tl, fs_signal)

L = np.array([70]) #sound-level of 70 dB SPL
stim_full = np.zeros((len(L), signalr.size))

for j in range(len(L)):
    stim_full[j, :] = p0 * 10**(L[j]/20) * signalr/rms(signalr)
    
fragment_length = 12345 #define fragment length (max 40000 for the included wav-file)
stim_length_init = factor_fs*(fragment_length+right_context+left_context)
stim_length = stim_length_init

#adapt fragment duration if no multiple of 16 (due to the CNN character of CoNNear, we need multiples of 16)
zero_pad = fragment_length%16
zeros = 0
if zero_pad != 0: 
    zeros = 16-zero_pad
    stim_length = factor_fs*(fragment_length+right_context+left_context+zeros)
        
################ CoNNear ####################################
stim = np.zeros((len(L), int(stim_length)))
stimrange = range(0, stim_length_init)
stim[:,stimrange] = stim_full[:,0:stim_length_init]
stim=sp_sig.resample_poly(stim, fs, fs_tl, axis=1)
stim=np.expand_dims(stim, axis=2)

tl_pred = connear.predict(stim)
tl_pred = tl_pred[0, :, :].T

################ Plots ######################################
fig, axarr = plt.subplots(2, sharex=True)
axarr[0].set_ylim(-0.35, 0.35)
axarr[0].plot(stim[0,(left_context):-(right_context),0])
axarr[0].set_title('Segment of Audio Input')

cax2 = axarr[1].imshow(tl_pred, cmap='bwr',aspect='auto', vmin=-0.5, vmax=0.5)  
axarr[1].set_title('Output of CoNNear')
axarr[1].set(ylabel='Center Frequency')

plt.show()<jupyter_output><empty_output><jupyter_text>DPOAE PlotsThe frequency response of the 12-kHz CF channel is evaluated as a proxy for the otoacoustic emissions recorded in the ear-canal. Frequency responses of model simulations are shown in response to two pure tones of $f_{1,2}$ of 2.0 and 2.4 kHz.The most pronounced distortion product in humans occurs at $2f_1 - f_2$ (1.6 kHz).**Notice that this part is computationally expensive (up to 20 mins), both in terms of memory and time.**<jupyter_code># Create the stimulus
L = [70.0]
f1 = 2000.
f2 = 1.2 * f1
L2 = 50.
L1 = 39 + 0.4 * L2 # scissors paradigm
print ("The tone frequencies are " + str(f1) + " and " + str(f2))
print ("with levels " + str(L1) + " and " + str(L2))
trailing_silence = 0.

# Here we will a pick stimulus longer than 2048 samples to get a better FFT
# We will prepare it with fs_tl sampling frequency
dur_sin_samples = np.lcm(int(f1), int(f2))
min_duration = 0.25 # in seconds
if dur_sin_samples < (min_duration * fs_tl):
    dur_sin = (((min_duration * fs_tl) // dur_sin_samples) + 1) * dur_sin_samples
else:
    dur_sin = dur_sin_samples
dur_sin = (dur_sin / fs_tl) 

t = np.arange(0., dur_sin, 1./fs_tl)

hanlength = int(10e-3 * fs_tl) # 10ms length hanning window
#f1
stim_sin1 = np.sin(2 * np.pi * f1 * t)
han = signal.windows.hann(hanlength)
stim_sin1 = p0 * np.sqrt(2) * 10**(L1/20) * stim_sin1
#f2
stim_sin2 = np.sin(2 * np.pi * f2 * t)
stim_sin2 = p0 * np.sqrt(2) * 10**(L2/20) * stim_sin2

stim_sin = stim_sin1 + stim_sin2

total_length =  int(trailing_silence * fs_tl) + len(stim_sin)
stim = np.zeros((1, int(total_length)))
stimrange = range(int(trailing_silence * fs_tl), int(trailing_silence * fs_tl) + len(stim_sin))
stim[0, stimrange] = stim_sin

################# CoNNear ####################################
# prepare for feeding to the CoNNear model
# first resample it to fs
shift_stim = 1
stim = stim[:, :]
stim = signal.decimate(stim, factor_fs, axis=1)
stim_1 = np.array(stim[0,:])
# window the signal into chunks of 2560 samples to be fed to the CoNNer model
stim = slice_1dsignal(stim_1, 2048, shift_stim, 256, left_context=256, right_context=256) 
connear_out_chunks = connear.predict(stim,verbose=1)
# undo the windowing to get back the full response
connear_out_full = undo_window (connear_out_chunks, 2048, shift_stim, ignore_first_set=0)
connear_out_full = connear_out_full[:,:stim_1.shape[0],:] * 1e-6

##############################################################
f_cf = 12000.
tone_index, tone_cf = min(enumerate(CF*1000), key=lambda x: abs( x [1]- f_cf))
print("CF nearest to " + str(f_cf) + " is " + str(CF[tone_index]*1000))
scale_val = (p0* np.sqrt(2))

################ Plots ######################################
# Plot the DPOAE CoNNear
connear_dpoae, nfft_connear = get_dpoae(connear_out_full, cf_location=tone_index)
    
freq_bins_connear = np.linspace(0, fs, num = nfft_connear)
plt.semilogx(freq_bins_connear[:int(nfft_connear/2)]/1000, 20 * np.log10(connear_dpoae/scale_val))
plt.title("DPOAEs - CoNNear")
plt.xlabel('Frequency [kHz]'), plt.ylabel('Magnitudes dbSPL'), plt.xlim((0.25, 8)), plt.grid()
plt.show()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1799252536038441, 'fraction_numerical': 0.040398647446164795, 'mean_word_length': 5.766405779650813, 'pattern_counts': {'":': 0, '<': 34, '<?xml version=': 0, '>': 33, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '8', 'dir': 'jupyter-structured-clean-dedup', 'id': '611275', 'license': 'RSA-MD', 'n_tokens_mistral': 4487, 'n_tokens_neox': 4078, 'n_words': 1045, 'path': 'connear_notebook_light.ipynb', 'repo_name': 'deepakbaby/CoNNear_cochlea', 'size': '11238'}