rmeddis@29: function vectorStrength=testAN(targetFrequency,BFlist, levels, ... rmeddis@29: paramsName,paramChanges) rmeddis@29: % testIHC used either for IHC I/O function ... rmeddis@29: % or receptive field (doReceptiveFields=1) rmeddis@29: rmeddis@29: global IHC_VResp_VivoParams IHC_cilia_RPParams IHCpreSynapseParams rmeddis@29: global AN_IHCsynapseParams rmeddis@29: rmeddis@29: global ANoutput ANdt CNoutput ICoutput ICmembraneOutput tauCas rmeddis@29: global ARattenuation MOCattenuation rmeddis@29: rmeddis@29: dbstop if error rmeddis@29: restorePath=path; rmeddis@29: rmeddis@29: addpath (['..' filesep 'MAP'], ['..' filesep 'utilities'], ... rmeddis@29: ['..' filesep 'parameterStore'], ['..' filesep 'wavFileStore'],... rmeddis@29: ['..' filesep 'testPrograms']) rmeddis@29: rmeddis@29: if nargin<5 rmeddis@29: paramChanges=[]; rmeddis@29: end rmeddis@29: rmeddis@29: if nargin<4 rmeddis@29: paramsName='Normal'; rmeddis@29: end rmeddis@29: rmeddis@29: if nargin<3 rmeddis@29: levels=-10:10:80; rmeddis@29: % levels=80:10:90; rmeddis@29: end rmeddis@29: rmeddis@29: nLevels=length(levels); rmeddis@29: rmeddis@29: toneDuration=.2; rmeddis@29: rampDuration=0.002; rmeddis@29: silenceDuration=.02; rmeddis@29: localPSTHbinwidth=0.001; rmeddis@29: rmeddis@29: % Use only the first frequency in the GUI targetFrequency box to defineBF rmeddis@29: % targetFrequency=stimulusParameters.targetFrequency(1); rmeddis@29: % BFlist=targetFrequency; rmeddis@29: rmeddis@29: AN_HSRonset=zeros(nLevels,1); rmeddis@29: AN_HSRsaturated=zeros(nLevels,1); rmeddis@29: AN_LSRonset=zeros(nLevels,1); rmeddis@29: AN_LSRsaturated=zeros(nLevels,1); rmeddis@29: CNLSRrate=zeros(nLevels,1); rmeddis@29: CNHSRsaturated=zeros(nLevels,1); rmeddis@29: ICHSRsaturated=zeros(nLevels,1); rmeddis@29: ICLSRsaturated=zeros(nLevels,1); rmeddis@29: vectorStrength=zeros(nLevels,1); rmeddis@29: rmeddis@29: AR=zeros(nLevels,1); rmeddis@29: MOC=zeros(nLevels,1); rmeddis@29: rmeddis@29: % ANoutput=zeros(200,200); rmeddis@29: rmeddis@29: figure(15), clf rmeddis@29: set(gcf,'position',[980 356 401 321]) rmeddis@29: figure(5), clf rmeddis@29: set(gcf,'position', [980 34 400 295]) rmeddis@29: drawnow rmeddis@29: rmeddis@29: %% guarantee that the sample rate is at least 10 times the frequency rmeddis@29: sampleRate=50000; rmeddis@29: while sampleRate< 10* targetFrequency rmeddis@29: sampleRate=sampleRate+10000; rmeddis@29: end rmeddis@29: rmeddis@29: %% adjust sample rate so that the pure tone stimulus has an integer rmeddis@29: %% numver of epochs in a period rmeddis@29: dt=1/sampleRate; rmeddis@29: period=1/targetFrequency; rmeddis@29: ANspeedUpFactor=5; %anticipating MAP (needs clearing up) rmeddis@29: ANdt=dt*ANspeedUpFactor; rmeddis@29: ANdt=period/round(period/ANdt); rmeddis@29: dt=ANdt/ANspeedUpFactor; rmeddis@29: sampleRate=1/dt; rmeddis@29: epochsPerPeriod=sampleRate*period; rmeddis@29: rmeddis@29: %% main computational loop (vary level) rmeddis@29: levelNo=0; rmeddis@29: for leveldB=levels rmeddis@29: levelNo=levelNo+1; rmeddis@29: rmeddis@29: fprintf('%4.0f\t', leveldB) rmeddis@29: amp=28e-6*10^(leveldB/20); rmeddis@29: rmeddis@29: time=dt:dt:toneDuration; rmeddis@29: rampTime=dt:dt:rampDuration; rmeddis@29: ramp=[0.5*(1+cos(2*pi*rampTime/(2*rampDuration)+pi)) ... rmeddis@29: ones(1,length(time)-length(rampTime))]; rmeddis@29: ramp=ramp.*fliplr(ramp); rmeddis@29: rmeddis@29: silence=zeros(1,round(silenceDuration/dt)); rmeddis@29: rmeddis@29: % create signal (leveldB/ targetFrequency) rmeddis@29: inputSignal=amp*sin(2*pi*targetFrequency'*time); rmeddis@29: inputSignal= ramp.*inputSignal; rmeddis@29: inputSignal=[silence inputSignal]; rmeddis@29: rmeddis@29: %% run the model rmeddis@29: AN_spikesOrProbability='spikes'; rmeddis@29: showPlotsAndDetails=0; rmeddis@29: rmeddis@29: rmeddis@29: MAP1_14(inputSignal, 1/dt, BFlist, ... rmeddis@29: paramsName, AN_spikesOrProbability, paramChanges); rmeddis@29: rmeddis@29: nTaus=length(tauCas); rmeddis@29: rmeddis@29: %LSR (same as HSR if no LSR fibers present) rmeddis@29: [nANFibers nTimePoints]=size(ANoutput); rmeddis@29: rmeddis@29: numLSRfibers=nANFibers/nTaus; rmeddis@29: numHSRfibers=numLSRfibers; rmeddis@29: rmeddis@29: LSRspikes=ANoutput(1:numLSRfibers,:); rmeddis@29: PSTH=UTIL_PSTHmaker(LSRspikes, ANdt, localPSTHbinwidth); rmeddis@29: PSTHLSR=mean(PSTH,1)/localPSTHbinwidth; % across fibers rates rmeddis@29: PSTHtime=localPSTHbinwidth:localPSTHbinwidth:... rmeddis@29: localPSTHbinwidth*length(PSTH); rmeddis@29: AN_LSRonset(levelNo)= max(PSTHLSR); % peak in 5 ms window rmeddis@29: AN_LSRsaturated(levelNo)= mean(PSTHLSR(round(length(PSTH)/2):end)); rmeddis@29: rmeddis@29: % HSR rmeddis@29: HSRspikes= ANoutput(end- numHSRfibers+1:end, :); rmeddis@29: PSTH=UTIL_PSTHmaker(HSRspikes, ANdt, localPSTHbinwidth); rmeddis@29: PSTH=mean(PSTH,1)/localPSTHbinwidth; % sum across fibers (HSR only) rmeddis@29: AN_HSRonset(levelNo)= max(PSTH); rmeddis@29: AN_HSRsaturated(levelNo)= mean(PSTH(round(length(PSTH)/2): end)); rmeddis@29: rmeddis@29: figure(5), subplot(2,2,2) rmeddis@29: hold off, bar(PSTHtime,PSTH, 'b') rmeddis@29: hold on, bar(PSTHtime,PSTHLSR,'r') rmeddis@29: ylim([0 1000]) rmeddis@29: xlim([0 length(PSTH)*localPSTHbinwidth]) rmeddis@29: set(gcf,'name',[num2str(BFlist), ' Hz: ' num2str(leveldB) ' dB']); rmeddis@29: rmeddis@29: % AN - CV rmeddis@29: % CV is computed 5 times. Use the middle one (3) as most typical rmeddis@29: cvANHSR= UTIL_CV(HSRspikes, ANdt); rmeddis@29: rmeddis@29: % AN - vector strength rmeddis@29: PSTH=sum(HSRspikes); rmeddis@29: [PH, binTimes]=UTIL_periodHistogram... rmeddis@29: (PSTH, ANdt, targetFrequency); rmeddis@29: VS=UTIL_vectorStrength(PH); rmeddis@29: vectorStrength(levelNo)=VS; rmeddis@29: disp(['sat rate= ' num2str(AN_HSRsaturated(levelNo)) ... rmeddis@29: '; phase-locking VS = ' num2str(VS)]) rmeddis@29: title(['AN HSR: CV=' num2str(cvANHSR(3),'%5.2f') ... rmeddis@29: 'VS=' num2str(VS,'%5.2f')]) rmeddis@29: rmeddis@29: % CN - first-order neurons rmeddis@29: rmeddis@29: % CN LSR rmeddis@29: [nCNneurons c]=size(CNoutput); rmeddis@29: nLSRneurons=round(nCNneurons/nTaus); rmeddis@29: CNLSRspikes=CNoutput(1:nLSRneurons,:); rmeddis@29: PSTH=UTIL_PSTHmaker(CNLSRspikes, ANdt, localPSTHbinwidth); rmeddis@29: PSTH=sum(PSTH)/nLSRneurons; rmeddis@29: CNLSRrate(levelNo)=mean(PSTH(round(length(PSTH)/2):end))/localPSTHbinwidth; rmeddis@29: rmeddis@29: %CN HSR rmeddis@29: MacGregorMultiHSRspikes=... rmeddis@29: CNoutput(end-nLSRneurons:end,:); rmeddis@29: PSTH=UTIL_PSTHmaker(MacGregorMultiHSRspikes, ANdt, localPSTHbinwidth); rmeddis@29: PSTH=sum(PSTH)/nLSRneurons; rmeddis@29: PSTH=mean(PSTH,1)/localPSTHbinwidth; % sum across fibers (HSR only) rmeddis@29: rmeddis@29: CNHSRsaturated(levelNo)=mean(PSTH(length(PSTH)/2:end)); rmeddis@29: rmeddis@29: figure(5), subplot(2,2,3) rmeddis@29: bar(PSTHtime,PSTH) rmeddis@29: ylim([0 1000]) rmeddis@29: xlim([0 length(PSTH)*localPSTHbinwidth]) rmeddis@29: cvMMHSR= UTIL_CV(MacGregorMultiHSRspikes, ANdt); rmeddis@29: title(['CN CV= ' num2str(cvMMHSR(3),'%5.2f')]) rmeddis@29: rmeddis@29: % IC LSR rmeddis@29: [nICneurons c]=size(ICoutput); rmeddis@29: nLSRneurons=round(nICneurons/nTaus); rmeddis@29: ICLSRspikes=ICoutput(1:nLSRneurons,:); rmeddis@29: PSTH=UTIL_PSTHmaker(ICLSRspikes, ANdt, localPSTHbinwidth); rmeddis@29: ICLSRsaturated(levelNo)=mean(PSTH(round(length(PSTH)/2):end))/localPSTHbinwidth; rmeddis@29: rmeddis@29: %IC HSR rmeddis@29: MacGregorMultiHSRspikes=... rmeddis@29: ICoutput(end-nLSRneurons:end,:); rmeddis@29: PSTH=UTIL_PSTHmaker(MacGregorMultiHSRspikes, ANdt, localPSTHbinwidth); rmeddis@29: PSTH=sum(PSTH)/nLSRneurons; rmeddis@29: PSTH=mean(PSTH,1)/localPSTHbinwidth; % sum across fibers (HSR only) rmeddis@29: rmeddis@29: ICHSRsaturated(levelNo)=mean(PSTH(length(PSTH)/2:end)); rmeddis@29: rmeddis@29: AR(levelNo)=min(ARattenuation); rmeddis@29: MOC(levelNo)=min(MOCattenuation(length(MOCattenuation)/2:end)); rmeddis@29: rmeddis@29: time=dt:dt:dt*size(ICmembraneOutput,2); rmeddis@29: figure(5), subplot(2,2,4) rmeddis@29: plot(time,ICmembraneOutput(2, 1:end),'k') rmeddis@29: ylim([-0.07 0]) rmeddis@29: xlim([0 max(time)]) rmeddis@29: title(['IC ' num2str(leveldB,'%4.0f') 'dB']) rmeddis@29: drawnow rmeddis@29: rmeddis@29: figure(5), subplot(2,2,1) rmeddis@29: plot(20*log10(MOC), 'k'), rmeddis@29: title(' MOC'), ylabel('dB attenuation') rmeddis@29: ylim([-30 0]) rmeddis@29: rmeddis@29: rmeddis@29: end % level rmeddis@29: figure(5), subplot(2,2,1) rmeddis@29: plot(levels,20*log10(MOC), 'k'), rmeddis@29: title(' MOC'), ylabel('dB attenuation') rmeddis@29: ylim([-30 0]) rmeddis@29: xlim([0 max(levels)]) rmeddis@29: rmeddis@29: fprintf('\n') rmeddis@29: toneDuration=2; rmeddis@29: rampDuration=0.004; rmeddis@29: silenceDuration=.02; rmeddis@29: nRepeats=200; % no. of AN fibers rmeddis@29: fprintf('toneDuration %6.3f\n', toneDuration) rmeddis@29: fprintf(' %6.0f AN fibers (repeats)\n', nRepeats) rmeddis@29: fprintf('levels') rmeddis@29: fprintf('%6.2f\t', levels) rmeddis@29: fprintf('\n') rmeddis@29: rmeddis@29: rmeddis@29: % ---------------------------------------------------- display parameters rmeddis@29: rmeddis@29: figure(15), clf rmeddis@29: nRows=2; nCols=2; rmeddis@29: rmeddis@29: % AN rate - level ONSET functions rmeddis@29: subplot(nRows,nCols,1) rmeddis@29: plot(levels,AN_LSRonset,'ro'), hold on rmeddis@29: plot(levels,AN_HSRonset,'ko'), hold off rmeddis@29: ylim([0 1000]), xlim([min(levels) max(levels)]) rmeddis@29: ttl=['tauCa= ' num2str(IHCpreSynapseParams.tauCa)]; rmeddis@29: title( ttl) rmeddis@29: xlabel('level dB SPL'), ylabel('peak rate (sp/s)'), grid on rmeddis@29: text(0, 800, 'AN onset', 'fontsize', 16) rmeddis@29: rmeddis@29: % AN rate - level ADAPTED function rmeddis@29: subplot(nRows,nCols,2) rmeddis@29: plot(levels,AN_LSRsaturated, 'ro'), hold on rmeddis@29: plot(levels,AN_HSRsaturated, 'ko'), hold off rmeddis@29: ylim([0 400]) rmeddis@29: set(gca,'ytick',0:50:300) rmeddis@29: xlim([min(levels) max(levels)]) rmeddis@29: set(gca,'xtick',[levels(1):20:levels(end)]) rmeddis@29: % grid on rmeddis@29: ttl=[ 'spont=' num2str(mean(AN_HSRsaturated(1,:)),'%4.0f')... rmeddis@29: ' sat=' num2str(mean(AN_HSRsaturated(end,1)),'%4.0f')]; rmeddis@29: title( ttl) rmeddis@29: xlabel('level dB SPL'), ylabel ('adapted rate (sp/s)') rmeddis@29: text(0, 340, 'AN adapted', 'fontsize', 16), grid on rmeddis@29: rmeddis@29: % CN rate - level ADAPTED function rmeddis@29: subplot(nRows,nCols,3) rmeddis@29: plot(levels,CNLSRrate, 'ro'), hold on rmeddis@29: plot(levels,CNHSRsaturated, 'ko'), hold off rmeddis@29: ylim([0 400]) rmeddis@29: set(gca,'ytick',0:50:300) rmeddis@29: xlim([min(levels) max(levels)]) rmeddis@29: set(gca,'xtick',[levels(1):20:levels(end)]) rmeddis@29: % grid on rmeddis@29: ttl=[ 'spont=' num2str(mean(CNHSRsaturated(1,:)),'%4.0f') ' sat=' ... rmeddis@29: num2str(mean(CNHSRsaturated(end,1)),'%4.0f')]; rmeddis@29: title( ttl) rmeddis@29: xlabel('level dB SPL'), ylabel ('adapted rate (sp/s)') rmeddis@29: text(0, 350, 'CN', 'fontsize', 16), grid on rmeddis@29: rmeddis@29: % IC rate - level ADAPTED function rmeddis@29: subplot(nRows,nCols,4) rmeddis@29: plot(levels,ICLSRsaturated, 'ro'), hold on rmeddis@29: plot(levels,ICHSRsaturated, 'ko'), hold off rmeddis@29: ylim([0 400]) rmeddis@29: set(gca,'ytick',0:50:300) rmeddis@29: xlim([min(levels) max(levels)]) rmeddis@29: set(gca,'xtick',[levels(1):20:levels(end)]), grid on rmeddis@29: rmeddis@29: rmeddis@29: ttl=['spont=' num2str(mean(ICHSRsaturated(1,:)),'%4.0f') ... rmeddis@29: ' sat=' num2str(mean(ICHSRsaturated(end,1)),'%4.0f')]; rmeddis@29: title( ttl) rmeddis@29: xlabel('level dB SPL'), ylabel ('adapted rate (sp/s)') rmeddis@29: text(0, 350, 'IC', 'fontsize', 16) rmeddis@29: set(gcf,'name',' AN CN IC rate/level') rmeddis@29: rmeddis@29: peakVectorStrength=max(vectorStrength); rmeddis@29: rmeddis@29: fprintf('\n') rmeddis@29: disp('levels vectorStrength') rmeddis@29: fprintf('%3.0f \t %6.4f \n', [levels; vectorStrength']) rmeddis@29: fprintf('\n') rmeddis@29: fprintf('Phase locking, max vector strength=\t %6.4f\n\n',... rmeddis@29: max(vectorStrength)) rmeddis@29: rmeddis@29: allData=[ levels' AN_HSRonset AN_HSRsaturated... rmeddis@29: AN_LSRonset AN_LSRsaturated ... rmeddis@29: CNHSRsaturated CNLSRrate... rmeddis@29: ICHSRsaturated ICLSRsaturated]; rmeddis@29: fprintf('\n levels \tANHSR Onset \tANHSR adapted\tANLSR Onset \tANLSR adapted\tCNHSR\tCNLSR\tICHSR \tICLSR \n'); rmeddis@29: UTIL_printTabTable(round(allData)) rmeddis@29: fprintf('VS (phase locking)= \t%6.4f\n\n',... rmeddis@29: max(vectorStrength)) rmeddis@29: rmeddis@29: UTIL_showStruct(IHC_cilia_RPParams, 'IHC_cilia_RPParams') rmeddis@29: UTIL_showStruct(IHCpreSynapseParams, 'IHCpreSynapseParams') rmeddis@29: UTIL_showStruct(AN_IHCsynapseParams, 'AN_IHCsynapseParams') rmeddis@29: rmeddis@29: fprintf('\n') rmeddis@29: disp('levels vectorStrength') rmeddis@29: fprintf('%3.0f \t %6.4f \n', [levels; vectorStrength']) rmeddis@29: fprintf('\n') rmeddis@29: fprintf('Phase locking, max vector strength= \t%6.4f\n\n',... rmeddis@29: max(vectorStrength)) rmeddis@29: rmeddis@29: allData=[ levels' AN_HSRonset AN_HSRsaturated... rmeddis@29: AN_LSRonset AN_LSRsaturated ... rmeddis@29: CNHSRsaturated CNLSRrate... rmeddis@29: ICHSRsaturated ICLSRsaturated]; rmeddis@29: fprintf('\n levels \tANHSR Onset \tANHSR adapted\tANLSR Onset \tANLSR adapted\tCNHSR\tCNLSR\tICHSR \tICLSR \n'); rmeddis@29: UTIL_printTabTable(round(allData)) rmeddis@29: fprintf('VS (phase locking)= \t%6.4f\n\n',... rmeddis@29: max(vectorStrength)) rmeddis@29: rmeddis@29: path(restorePath)