THE AUDITORY MODELING TOOLBOX
Go to function
BRUCE2018 - Auditory-nerve filterbank (improved synapse)
Program code:
function [output,info] = bruce2018(stim,fsstim, fc, varargin)
%BRUCE2018 Auditory-nerve filterbank (improved synapse)
% Usage: [output] = bruce2018(stim,fsstim, fc);
% [output,info] = bruce2018(stim, fsstim, fc, varargin);
%
%
% Input parameters:
% stim : Pressure waveform of stimulus (timeseries)
%
% fsstim : Sampling frequency of stimulus
%
% fc : Frequency vector containing the CFs.
% Use logspace(log10(flow), log10(fhigh),numCF) to
% replicate the results from Bruce et al. (2018).
%
% varargin : various flags and key-value pairs.
%
%
% Output parameters:
% output : a struct containing the various modelstage outputs
% The struct output contains:
%
% - cohcs: actually used COHCs in the simulations
% - 'cihcs' : actually used CIHCs in the simulations
% - sponts : actually used spontanous rates in the simulations
% - tabss : actually used absolute timings in the simulations
% - trels : actually used relative timings in the simulations
% - neurogram_ft : (fine-timing) neurogram calculated from the synapse output [time CFs]
% - t_ft : time axis of neurogram_ft (s)
% - neurogram_mr : average of spike count in each PSTH bin [time CF]
% - t_mr : time axis of neurogram_mr (s)
% - neurogram_Sout : synapse output [time CF]
% - t_Sout : time axis of neurogram_Sout (s)
% - fc : actually used CFs (Hz)
% - psth_ft : fine-timing PSTH [time CF]
% - meanrate : average spiking rate (spikes/s) [time CF]
% - varrate : variance of the spiking rate (spikes/s) [time CF]
%
%
% BRUCE2018(...) returns modeled responses of multiple AN fibers tuned to
% various characteristic frequencies characterstic.
%
% Please cite the references below if you use this model.
%
% This function takes the following optional key/value pairs:
%
% 'ag_fs',ag_fs Frequencies at which the audiogram should be
% evaluated. Used only in fitaudiogram.
%
% 'ag_db',ag_db Hearing loss (dB) at frequencies ag_fs.
% Used only in fitaudiogram.
%
% 'cohcs',cohcs OHC scaling factors: 1 denotes normal OHC function (default);
% 0 denotes complete OHC dysfunction. Can be a vector
% of the size of fc or a scalar. Not used in fitaudiogram.
%
% 'cihcs',cihcs IHC scaling factors: 1 denotes normal IHC function (default);
% 0 denotes complete IHC dysfunction. Can be a vector
% of the size of fc or a scalar. Not used in fitaudiogram.
%
% 'numL',nl number of nerve fibres with low SR. Used only in autoSR or specificSRautoTiming.
%
% 'numM',nm number of nerve fibres with medium SR. Used only in autoSR or specificSRautoTiming.
%
% 'numH',nh number of nerve fibres with high SR. Used only in autoSR or specificSRautoTiming.
%
% 'lossL',lls loss of low-SR fibres ranging from 0 (no fibres)
% to 1 (healthy, all fibres). Used only in autoSR or specificSRautoTiming.
%
% 'lossM',lms loss of medium-SR fibres ranging from 0 (no fibres)
% to 1 (healthy, all fibres). Used only in autoSR or specificSRautoTiming.
%
% 'lossH',lhs loss of high-SR fibres ranging from 0 (no fibres)
% to 1 (healthy, all fibres). Used only in autoSR or specificSRautoTiming.
%
% 'SRL',SRL SR of the low-SR fibres (spikes/s). Used only in specificSRautoTiming.
%
% 'SRM',SRM SR of the medium-SR fibres (spikes/s). Used only in specificSRautoTiming.
%
% 'SRH',SRH SR of the high-SR fibres (spikes/s). Used only in specificSRautoTiming.
%
% 'numsponts',n Overall numbers of fibers. Used only in specificSR.
%
% 'spont',spont SR (spikes/s). Can be scalar or size of fc. Used only in specificSR.
%
% 'tabs',tabs Absolute timings (s). Can be scalar or size of fc. Used only in specificSR.
%
% 'trel',trel Relative timings (s). Can be scalar or size of fc. Used only in specificSR.
%
% 'psthbinwidth_mr',psthbw mean-rate binwidth (s).
%
% 'windur_ft',winft fine-timing neurogram window length.
%
% 'windur_mr',winmr mean-rate neurogram window length.
%
% 'nrep',nrep Number of repetitions for the mean rate,
% rate variance & psth calculation. Default is 1.
%
% 'reptime',rt length of one repetition of the stimuli with pause.
% Default is 1.2 stimulus duration.
%
% 'fsmod',fsmod Model sampling rate. It is possible to run the model
% at a range of fsmod between 100 kHz and 500 kHz.
% Default value is 200 kHz for cats and 100 kHz for humans.
%
% BRUCE2018 accepts the following flags:
%
% 'fitaudiogram' Calculate the hearing-loss factors cihcs and
% cohcs from the frequencies ag_fs and threshold
% shifts ag_dbloss (dB) by using BRUCE2018_fitaudiogram.
% The default parameters reflect a healthy cochlea.
%
% 'no_fitaudiogram' Default. Use directly provided cihcs and cohcs
% either as scalars or size of fc.
% The default parameters reflect a healthy cochlea.
%
% 'human' Default. Use model parameters for humans.
%
% 'cat' Use model parameters for cats.
%
% 'fixedFGn' Default. Fractional Gaussian noise will be the same in every
% simulation.
%
% 'varFGn' Fractional Gaussian noise will be different in every
% simulation.
%
% 'approxPL' Default. Use approxiate implementation of the power-law
% functions.
%
% 'actualPL' Use actual implementation of the power-law functions.
%
% 'outputPerSynapse' Output the synapse output of each individual
% nerve fibre.
% This can considerably slow down the calculations.
%
% 'outputPerCF' Default. Output the average results over all synapses.
% This mode is faster than outputPerSynapse.
%
% 'autoSR' Generate the parameters for the AN
% population wil be generated by the function
% BRUCE2018_generateanpopulation based on
% the number of low, medium, and high SR fibres
% numL, numM, numH. Default.
%
% 'specificSRautoTiming' Generate the timing parameters for the AN
% population by the function
% BRUCE2018_generateanpopulation based on
% the number of low, medium, and high SR fibres
% numL, numM, numH. The spontanous rates are
% provided in SRL, SRM, and SRH.
%
% 'specificSR' Do not generate AN parameters. The overall number
% of fibres numsponts, spontanouse rate spont,
% absolute timing tabs, and relative timing trel
% must be provided.
%
%
% If 'outputPerSynapse' is specified, psth_ft, meanrate, and
% varrate have the dimensions [time CF Syn] (with Syn as the number
% of synapses) and output additionally contains:
%
% 'synout' the output of a synapse [time CF Syn]
%
% 'info' a struct containing the parameter settings applied
%
%
%
% See also: plot_bruce2018 demo_bruce2018_auditorynervemodel
% demo_bruce2018 bruce2018_synapse bruce2018_generateanpopulation
% bruce2018_innerhaircells bruce2018_fitaudiogram bruce2018_ffgn
% exp_bruce2018 demo_carney2015 carney2015_fitaudiogram carney2015_generateneurogram
% exp_osses2022 zilany2014
%
% References:
% I. C. Bruce, Y. Erfani, and M. S. R. Zilany. A phenomenological model
% of the synapse between the inner hair cell and auditory nerve:
% Implications of limited neurotransmitter release sites. Hearing
% Research, 360:40--54, 2018.
%
%
% Url: http://amtoolbox.org/amt-1.2.0/doc/models/bruce2018.php
% Copyright (C) 2009-2022 Piotr Majdak, Clara Hollomey, and the AMT team.
% This file is part of Auditory Modeling Toolbox (AMT) version 1.2.0
%
% This program is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% This program is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with this program. If not, see <http://www.gnu.org/licenses/>.
% #StatusDoc: Good
% #StatusCode: Perfect
% #Verification: Verified
% #Requirements: MATLAB MEX M-Signal
% #Author: Ian Bruce: basic code of the model
% #Author: Alejandro Osses (2020): original implementation
% #Author: Clara Hollomey (2021): adapted to the AMT 1.0
% #Author: Piotr Majdak (2021): adaptations to exp_osses2022; specificSRautoTiming added
%
% The other two modelstages 'bruce2018_innerhaircells' and
% 'bruce2018_synapse' are always active. '_innerhaircells' is called for
% each element in 'fcs' (each characteristic frequency), and '_synapse'
% is called for each nerve fiber. Per default and for execution speed,
% only results from outside of the nerve fiber loop are written to the
% struct 'output'. To retrieve results per fiber, set the flag
% 'ouputPerSynapse'. The actual parameters used in bruce2018 are
% output to the 'info' struct.
if nargin<3
error('%s: Too few input parameters.',upper(mfilename));
end;
% Define input flags and values
definput.import = {'bruce2018'}; % load defaults from arg_bruce2018
[flags,kv] = ltfatarghelper({},definput,varargin);
% Derive the number of CF fibres
numCF = length(fc);
fs=kv.fsmod; % sampling rate (Hz) of the model
tdres = 1/fs; % sampling interval (s), i.e., the reciprocal of fs
stim = resample(stim,fs,fsstim); % resample the stimulus to the model fs
T = length(stim)/fs; % actual stimulus duration in seconds
smw_ft = hamming(kv.windur_ft);
smw_mr = hamming(kv.windur_mr);
simdur = ceil(T*kv.reptime/kv.psthbinwidth_mr)*kv.psthbinwidth_mr; % time interval (s) between stimulus repetitions
% species is either "cat" (1) or "human" (2): "1" for cat and "2" for human
if flags.do_cat, species = 1; else species = 2; end
% noiseType is for fixed fGn (noise will be same in every simulation) or
% variable fGn: "0" for fixed fGn and "1" for variable fGn
if flags.do_varFGn, noiseType = 0; else noiseType = 1; end
% implnt is for "approxiate" or "actual" implementation of the power-law functions:
%"0" for approx. and "1" for actual implementation
if flags.do_actualPL, implnt = 1; else implnt = 0; end
if flags.do_fitaudiogram
% mixed loss
dbloss = interp1(kv.ag_fs,kv.ag_dbloss,fcs,'linear','extrap');
[cohcs,cihcs,OHC_loss, IHC_loss]=bruce2018_fitaudiogram(fc,dbloss,species);
info.OHC_loss = OHC_loss;
info.IHC_loss = IHC_loss;
else
if length(kv.cihcs) ~= numCF && length(kv.cihcs) ~= 1
error('cihc needs to either be a scalar or a vector of the same size as fc.')
end
if length(kv.cohcs) ~= numCF && length(kv.cohcs) ~= 1
error('cohc needs to either be a scalar or a vector of the same size as fc.')
end
cohcs = kv.cohcs.* ones(1, numCF);
cihcs = kv.cihcs.* ones(1, numCF);
info.OHC_loss = [];
info.IHC_loss = [];
end
% generate AN parameters?
if flags.do_autoSR,
% generate SR and timing for a population based on numL, numM, and numH
[sponts,tabss,trels] = bruce2018_generateanpopulation(numCF,[kv.numL kv.numM kv.numH]);
numsponts = round([kv.lossL kv.lossM kv.lossH].*[kv.numL kv.numM kv.numH]);
cntsponts = sum(numsponts);
end
if flags.do_specificSRautoTiming,
% use specific SRs but calculate the timing based on numL, numM, and numH
[~,tabss,trels] = bruce2018_generateanpopulation(numCF,[kv.numL kv.numM kv.numH]);
sponts.LS = kv.SRL*ones(1,kv.numL);
sponts.MS = kv.SRM*ones(1,kv.numM);
sponts.HS = kv.SRH*ones(1,kv.numH);
numsponts = round([kv.lossL kv.lossM kv.lossH].*[kv.numL kv.numM kv.numH]);
cntsponts = sum(numsponts);
end
if flags.do_specificSR,
% use predefined SR parameters for all fibers, used in e.g., exp_bruce2018('fig8b');
if ~isscalar(kv.numsponts), error('numsponts needs to a scalar when running bruce2018 in specificSR mode.'); end
cntsponts = kv.numsponts;
sponts = kv.spont;
tabss = kv.tabs;
trels = kv.trel;
end
% initializations
init_size = simdur*fs;
output_len=floor(simdur/tdres+0.5)*kv.nrep; % length of vihc, c1, c2, synout
psth_len=output_len/kv.nrep; % length of vr, mr, psth_ft
psthbins = round(kv.psthbinwidth_mr*fs); % number of psth_ft bins per psth bin
neurogram_ft = zeros(round(init_size),numCF);
neurogram_Sout = zeros(round(init_size)*kv.nrep,numCF);
neurogram_mr = zeros(round(init_size)/round(kv.psthbinwidth_mr*fs),numCF);
output.vihc=zeros(output_len,numCF);
output.C1=zeros(output_len,numCF);
output.C2=zeros(output_len,numCF);
if flags.do_outputPerSynapse
output.psth_ft=zeros(psth_len,numCF,cntsponts);
output.meanrate=zeros(psth_len,numCF,cntsponts);
output.varrate=zeros(psth_len,numCF,cntsponts);
output.psth=zeros(psth_len/psthbins,numCF,cntsponts);
end
if flags.do_outputPerCF
% Memory allocation:
output.psth_ft=zeros(psth_len,numCF);
output.meanrate=zeros(psth_len,numCF);
output.varrate=zeros(psth_len,numCF);
output.psth=zeros(psth_len/psthbins,numCF);
output.meanrate_LSR = zeros(size(output.meanrate));
output.meanrate_MSR = zeros(size(output.meanrate));
output.meanrate_HSR = zeros(size(output.meanrate));
output.psth_LSR = zeros(size(output.psth));
output.psth_MSR = zeros(size(output.psth));
output.psth_HSR = zeros(size(output.psth));
end
for ii = 1:numCF
%for each CF, collect the appropriate tabs, trels and sponts, and calculate
% the inner hair cell potential (for the whole stimulus)
% FC = fc(ii);
% cohc = cohcs(ii);
% cihc = cihcs(ii);
[vihc, C1, C2] = bruce2018_innerhaircells(stim,fc(ii),kv.nrep,tdres,simdur,cohcs(ii),cihcs(ii),species);
% [vihc, C1, C2] = zilany2014_innerhaircells(stim,FC,kv.nrep,tdres,simdur,cohc,cihc,species); equivalent to bruce2018
output.vihc(:,ii) = vihc;
output.C1(:,ii) = C1;
output.C2(:,ii) = C2;
if cntsponts>0
if flags.do_autoSR || flags.do_specificSRautoTiming, % use generated SR parameters
spont = [sponts.LS(ii,1:numsponts(1)) sponts.MS(ii,1:numsponts(2)) sponts.HS(ii,1:numsponts(3))];
tabs = [tabss.LS(ii,1:numsponts(1)) tabss.MS(ii,1:numsponts(2)) tabss.HS(ii,1:numsponts(3))];
trel = [trels.LS(ii,1:numsponts(1)) trels.MS(ii,1:numsponts(2)) trels.HS(ii,1:numsponts(3))];
end
if flags.do_specificSR % use the same SR parameters for all fibers
spont = sponts.*ones(1, cntsponts);
tabs = tabss.*ones(1, cntsponts);
trel = trels.*ones(1, cntsponts);
end
for jj = 1:cntsponts % calculate the synapse output for each fibre type
amt_disp(['CF = ' int2str(ii) '/' int2str(numCF) '; spont = ' int2str(jj) '/' int2str(cntsponts) '; SR = ' num2str(spont(jj)) ' spikes/s'],'volatile');
[psth_ft,mr,vr,synout] = bruce2018_synapse(vihc,fc(ii),kv.nrep,tdres,noiseType,implnt,spont(jj),tabs(jj),trel(jj));
neurogram_Sout(:,ii) = neurogram_Sout(:,ii)+synout;
cnt = sum(reshape(psth_ft,psthbins,length(psth_ft)/psthbins))';
neurogram_ft(:,ii) = neurogram_ft(:,ii)+filter(smw_ft,1,psth_ft);
neurogram_mr(:,ii) = neurogram_mr(:,ii)+filter(smw_mr,1,cnt);
if flags.do_outputPerSynapse
output.psth_ft(:, ii, jj) = psth_ft;
output.meanrate(:, ii, jj) = mr;
output.varrate(:, ii, jj) = vr;
output.synout(:, ii, jj) = synout;
output.psth(:,ii,jj)=cnt;
end
if flags.do_outputPerCF,
output.psth_ft(:,ii)=output.psth_ft(:,ii)+psth_ft;
output.meanrate(:,ii)=output.meanrate(:,ii)+mr;
output.varrate(:,ii)=output.varrate(:,ii)+vr;
output.psth(:,ii)=output.psth(:,ii)+cnt;
fiberType = find(jj <= cumsum([kv.numL kv.numM kv.numH]),1,'first');
switch fiberType
case 1
output.meanrate_LSR(:,ii) = output.meanrate_LSR(:,ii) + mr;
output.psth_LSR(:,ii) = output.psth_LSR(:,ii) + cnt;
case 2
output.meanrate_MSR(:,ii) = output.meanrate_MSR(:,ii) + mr;
output.psth_MSR(:,ii) = output.psth_MSR(:,ii) + cnt;
case 3
output.meanrate_HSR(:,ii) = output.meanrate_HSR(:,ii) + mr;
output.psth_HSR(:,ii) = output.psth_HSR(:,ii) + cnt;
end
end
end
amt_disp();
if flags.do_outputPerCF
output.psth_ft(:,ii)=output.psth_ft(:,ii)/cntsponts;
output.meanrate(:,ii)=output.meanrate(:,ii)/cntsponts;
output.varrate(:,ii)=output.varrate(:,ii)/cntsponts;
if kv.numL ~= 0
output.meanrate_LSR(:,ii)=output.meanrate_LSR(:,ii)/kv.numL;
end
if kv.numM ~= 0
output.meanrate_MSR(:,ii)=output.meanrate_MSR(:,ii)/kv.numM;
end
if kv.numH ~= 0
output.meanrate_HSR(:,ii)=output.meanrate_HSR(:,ii)/kv.numH;
end
end
end
end
% amt_disp();
neurogram_ft = neurogram_ft(1:kv.windur_ft/2:end,:); % 50% overlap in Hamming window
t_ft = 0:kv.windur_ft/2/fs:(size(neurogram_ft,1)-1)*kv.windur_ft/2/fs; % time vector for the fine-timing neurogram
neurogram_mr = neurogram_mr(1:kv.windur_mr/2:end,:); % 50% overlap in Hamming window
t_mr = 0:kv.windur_mr/2*kv.psthbinwidth_mr:(size(neurogram_mr,1)-1)*kv.windur_mr/2*kv.psthbinwidth_mr; % time vector for the mean-rate neurogram
t_Sout = 0:1/fs:(size(neurogram_Sout,1)-1)/fs; % time vector for the synapse output neurogram
%--------------------------------------------------------------------------
% write everything into an output struct
output.cohcs = cohcs;
output.cihcs = cihcs;
output.sponts = sponts;
output.tabss = tabss;
output.trels = trels;
output.t_ft = t_ft;
output.neurogram_ft = neurogram_ft;
output.neurogram_mr = neurogram_mr;
output.neurogram_Sout = neurogram_Sout;
output.t_mr = t_mr;
output.t_Sout = t_Sout;
output.fc = fc;
info.keyvals = kv;
info.flags = flags;
Build with Bootstrap