rmeddis@22: function showMAP () rmeddis@22: % defaults rmeddis@22: % options.showModelParameters=1; rmeddis@22: % options.showModelOutput=1; rmeddis@22: % options.printFiringRates=1; rmeddis@22: % options.showACF=1; rmeddis@22: % options.showEfferent=1; rmeddis@22: % options.surfProbability=0; rmeddis@22: % options.fileName=[]; rmeddis@22: rmeddis@22: dbstop if warning rmeddis@22: rmeddis@22: global dt ANdt saveAN_spikesOrProbability savedBFlist saveMAPparamsName... rmeddis@22: savedInputSignal TMoutput OMEoutput ARattenuation ... rmeddis@22: DRNLoutput IHC_cilia_output IHCrestingCiliaCond IHCrestingV... rmeddis@22: IHCoutput ANprobRateOutput ANoutput savePavailable tauCas ... rmeddis@22: CNoutput ICoutput ICmembraneOutput ICfiberTypeRates MOCattenuation rmeddis@22: global OMEParams DRNLParams IHC_cilia_RPParams IHCpreSynapseParams rmeddis@22: global AN_IHCsynapseParams MacGregorParams MacGregorMultiParams rmeddis@22: global showMapOptions rmeddis@22: rmeddis@22: rmeddis@22: restorePath=path; rmeddis@22: addpath ( ['..' filesep 'utilities'], ['..' filesep 'parameterStore']) rmeddis@22: rmeddis@22: if nargin<1 rmeddis@22: options.showModelParameters=1; rmeddis@22: options.showModelOutput=1; rmeddis@22: options.printFiringRates=1; rmeddis@22: options.showACF=0; rmeddis@22: options.showEfferent=1; rmeddis@22: options.surfProbability=0; rmeddis@22: options.fileName=[]; rmeddis@22: end rmeddis@22: rmeddis@22: if options.showModelParameters rmeddis@22: % Read parameters from MAPparams<***> file in 'parameterStore' folder rmeddis@22: % and print out all parameters rmeddis@22: cmd=['MAPparams' saveMAPparamsName ... rmeddis@22: '(-1, 1/dt, 1);']; rmeddis@22: eval(cmd); rmeddis@22: end rmeddis@22: rmeddis@22: if options.printFiringRates rmeddis@22: %% print summary firing rates rmeddis@22: fprintf('\n\n') rmeddis@22: disp('summary') rmeddis@22: disp(['AR: ' num2str(min(ARattenuation))]) rmeddis@22: disp(['MOC: ' num2str(min(min(MOCattenuation)))]) rmeddis@22: nANfiberTypes=length(tauCas); rmeddis@22: if strcmp(saveAN_spikesOrProbability, 'spikes') rmeddis@22: nANfibers=size(ANoutput,1); rmeddis@22: nHSRfibers=nANfibers/nANfiberTypes; rmeddis@22: duration=size(TMoutput,2)*dt; rmeddis@22: disp(['AN: ' num2str(sum(sum(ANoutput(end-nHSRfibers+1:end,:)))/... rmeddis@22: (nHSRfibers*duration))]) rmeddis@22: rmeddis@22: nCNneurons=size(CNoutput,1); rmeddis@22: nHSRCNneuronss=nCNneurons/nANfiberTypes; rmeddis@22: disp(['CN: ' num2str(sum(sum(CNoutput(end-nHSRCNneuronss+1:end,:)))... rmeddis@22: /(nHSRCNneuronss*duration))]) rmeddis@22: disp(['IC: ' num2str(sum(sum(ICoutput))/duration)]) rmeddis@22: % disp(['IC by type: ' num2str(mean(ICfiberTypeRates,2)')]) rmeddis@22: else rmeddis@22: disp(['AN: ' num2str(mean(mean(ANprobRateOutput)))]) rmeddis@22: [PSTH pointsPerBin]= UTIL_makePSTH(ANprobRateOutput, dt, 0.001); rmeddis@22: disp(['max max AN: ' num2str(max(max(... rmeddis@22: PSTH/pointsPerBin )))]) rmeddis@22: end rmeddis@22: end rmeddis@22: rmeddis@22: rmeddis@22: %% figure (99) summarises main model output rmeddis@22: if options.showModelOutput rmeddis@22: plotInstructions.figureNo=99; rmeddis@22: signalRMS=mean(savedInputSignal.^2)^0.5; rmeddis@22: signalRMSdb=20*log10(signalRMS/20e-6); rmeddis@22: rmeddis@22: % plot signal (1) rmeddis@22: plotInstructions.displaydt=dt; rmeddis@22: plotInstructions.numPlots=6; rmeddis@22: plotInstructions.subPlotNo=1; rmeddis@22: plotInstructions.title=... rmeddis@22: ['stimulus: ' num2str(signalRMSdb, '%4.0f') ' dB SPL']; rmeddis@22: r=size(savedInputSignal,1); rmeddis@22: if r==1, savedInputSignal=savedInputSignal'; end rmeddis@22: UTIL_plotMatrix(savedInputSignal', plotInstructions); rmeddis@22: rmeddis@22: % stapes (2) rmeddis@22: plotInstructions.subPlotNo=2; rmeddis@22: plotInstructions.title= ['stapes displacement']; rmeddis@22: UTIL_plotMatrix(OMEoutput, plotInstructions); rmeddis@22: rmeddis@22: % DRNL (3) rmeddis@22: plotInstructions.subPlotNo=3; rmeddis@22: plotInstructions.title= ['BM displacement']; rmeddis@22: plotInstructions.yValues= savedBFlist; rmeddis@22: UTIL_plotMatrix(DRNLoutput, plotInstructions); rmeddis@22: rmeddis@22: switch saveAN_spikesOrProbability rmeddis@22: case 'spikes' rmeddis@22: % AN (4) rmeddis@22: plotInstructions.displaydt=ANdt; rmeddis@22: plotInstructions.title='AN'; rmeddis@22: plotInstructions.subPlotNo=4; rmeddis@22: plotInstructions.yLabel='BF'; rmeddis@22: plotInstructions.yValues= savedBFlist; rmeddis@22: plotInstructions.rasterDotSize=1; rmeddis@22: plotInstructions.plotDivider=1; rmeddis@22: if sum(sum(ANoutput))<100 rmeddis@22: plotInstructions.rasterDotSize=3; rmeddis@22: end rmeddis@22: UTIL_plotMatrix(ANoutput, plotInstructions); rmeddis@22: rmeddis@22: % CN (5) rmeddis@22: plotInstructions.displaydt=ANdt; rmeddis@22: plotInstructions.subPlotNo=5; rmeddis@22: plotInstructions.title='CN spikes'; rmeddis@22: if sum(sum(CNoutput))<100 rmeddis@22: plotInstructions.rasterDotSize=3; rmeddis@22: end rmeddis@22: UTIL_plotMatrix(CNoutput, plotInstructions); rmeddis@22: rmeddis@22: % IC (6) rmeddis@22: plotInstructions.displaydt=ANdt; rmeddis@22: plotInstructions.subPlotNo=6; rmeddis@22: plotInstructions.title='IC'; rmeddis@22: if size(ICoutput,1)>3 rmeddis@22: if sum(sum(ICoutput))<100 rmeddis@22: plotInstructions.rasterDotSize=3; rmeddis@22: end rmeddis@22: UTIL_plotMatrix(ICoutput, plotInstructions); rmeddis@22: else rmeddis@22: plotInstructions.title='IC (HSR) membrane potential'; rmeddis@22: plotInstructions.displaydt=dt; rmeddis@22: plotInstructions.yLabel='V'; rmeddis@22: plotInstructions.zValuesRange= [-.1 0]; rmeddis@22: UTIL_plotMatrix(ICmembraneOutput, plotInstructions); rmeddis@22: end rmeddis@22: rmeddis@22: otherwise % probability (4-6) rmeddis@22: plotInstructions.displaydt=dt; rmeddis@22: plotInstructions.numPlots=2; rmeddis@22: plotInstructions.subPlotNo=2; rmeddis@22: plotInstructions.yLabel='BF'; rmeddis@22: if nANfiberTypes>1, rmeddis@22: plotInstructions.yLabel='LSR HSR'; rmeddis@22: plotInstructions.plotDivider=1; rmeddis@22: end rmeddis@22: plotInstructions.title='AN - spike probability'; rmeddis@22: UTIL_plotMatrix(ANprobRateOutput, plotInstructions); rmeddis@22: end rmeddis@22: end rmeddis@22: rmeddis@22: %% surface plot of probability rmeddis@22: if options.surfProbability rmeddis@22: figure(97), clf rmeddis@22: % select only HSR fibers at the bottom of the matrix rmeddis@22: ANprobRateOutput= ANprobRateOutput(end-length(savedBFlist)+1:end,:); rmeddis@22: [nY nX]=size(ANprobRateOutput); rmeddis@22: if nY>2 rmeddis@22: time=dt*(1:nX); rmeddis@22: surf(time, savedBFlist, ANprobRateOutput) rmeddis@22: shading interp rmeddis@22: set(gca, 'yScale','log') rmeddis@22: xlim([0 max(time)]), ylim([0 max(savedBFlist)]), zlim([0 1000]) rmeddis@22: xlabel('time (s)') rmeddis@22: ylabel('best frequency (Hz)') rmeddis@22: zlabel('spike rate') rmeddis@22: view([-20 60]) rmeddis@22: if isfield(options, 'fileName') rmeddis@22: title ([options.fileName ': ' num2str(signalRMSdb,'% 3.0f') ' dB']) rmeddis@22: else rmeddis@22: title ([ num2str(signalRMSdb,'% 3.0f') ' dB']) rmeddis@22: end rmeddis@22: rmeddis@22: end rmeddis@22: end rmeddis@22: rmeddis@22: rmeddis@22: %% plot efferent control values as dB rmeddis@22: if options.showEfferent rmeddis@22: plotInstructions=[]; rmeddis@22: plotInstructions.figureNo=98; rmeddis@22: figure(98), clf rmeddis@22: plotInstructions.displaydt=dt; rmeddis@22: plotInstructions.numPlots=2; rmeddis@22: plotInstructions.subPlotNo=1; rmeddis@22: plotInstructions.zValuesRange=[ -25 0]; rmeddis@22: plotInstructions.title= ['AR strength. Signal level= ' ... rmeddis@22: num2str(signalRMSdb,'%4.0f') ' dB SPL']; rmeddis@22: UTIL_plotMatrix(20*log10(ARattenuation), plotInstructions); rmeddis@22: rmeddis@22: plotInstructions.subPlotNo=2; rmeddis@22: plotInstructions.yValues= savedBFlist; rmeddis@22: plotInstructions.yLabel= 'BF'; rmeddis@22: plotInstructions.title= ['MOC strength']; rmeddis@22: plotInstructions.zValuesRange=[ -25 0]; rmeddis@22: subplot(2,1,2) rmeddis@22: % imagesc(MOCattenuation) rmeddis@22: UTIL_plotMatrix(20*log10(MOCattenuation), plotInstructions); rmeddis@22: colorbar rmeddis@22: end rmeddis@22: <<<<<<< HEAD rmeddis@22: rmeddis@22: if options.surfProbability rmeddis@22: %% surface plot of probability rmeddis@22: figure(97), clf rmeddis@22: % select only HSR fibers at the bottom of the matrix rmeddis@22: ANprobRateOutput= ANprobRateOutput(end-length(savedBFlist)+1:end,:); rmeddis@22: PSTHbinWidth=0.001; rmeddis@22: PSTH=UTIL_PSTHmakerb(ANprobRateOutput, ANdt, PSTHbinWidth); rmeddis@22: [nY nX]=size(PSTH); rmeddis@22: PSTH=PSTH(:,200:end); rmeddis@22: [nY nX]=size(PSTH); rmeddis@22: time=PSTHbinWidth*(1:nX); rmeddis@22: surf(time, savedBFlist, PSTH) rmeddis@22: shading interp rmeddis@22: set(gca, 'yScale','log') rmeddis@22: xlim([0 max(time)]) rmeddis@22: ylim([0 max(savedBFlist)]) rmeddis@22: zlim([0 1000]) rmeddis@22: xlabel('time (s)') rmeddis@22: ylabel('best frequency (Hz)') rmeddis@22: zlabel('spike rate') rmeddis@22: view([-20 60]) rmeddis@22: % view([0 90]) rmeddis@22: title ([options.fileName ': ' num2str(signalRMSdb,'% 3.0f') ' dB']) rmeddis@22: end rmeddis@22: rmeddis@22: rmeddis@22: %% ACF plot if required rmeddis@22: if options.showACF rmeddis@22: tic rmeddis@22: method.dt=dt; rmeddis@22: method.segmentNo=1; rmeddis@22: method.nonlinCF=savedBFlist; rmeddis@22: rmeddis@22: minPitch= 80; maxPitch= 4000; numPitches=100; % specify lags rmeddis@22: pitches=10.^ linspace(log10(minPitch), log10(maxPitch),numPitches); rmeddis@22: pitches=fliplr(pitches); rmeddis@22: filteredSACFParams.lags=1./pitches; % autocorrelation lags vector rmeddis@22: filteredSACFParams.acfTau= .003; % time constant of running ACF rmeddis@22: filteredSACFParams.lambda= 0.12; % slower filter to smooth ACF rmeddis@22: filteredSACFParams.lambda= 0.01; % slower filter to smooth ACF rmeddis@22: rmeddis@22: filteredSACFParams.plotACFs=0; % special plot (see code) rmeddis@22: filteredSACFParams.plotFilteredSACF=0; % 0 plots unfiltered ACFs rmeddis@22: filteredSACFParams.plotMoviePauses=.3; % special plot (see code) rmeddis@22: rmeddis@22: filteredSACFParams.usePressnitzer=0; % attenuates ACF at long lags rmeddis@22: filteredSACFParams.lagsProcedure= 'useAllLags'; rmeddis@22: % filteredSACFParams.lagsProcedure= 'useBernsteinLagWeights'; rmeddis@22: % filteredSACFParams.lagsProcedure= 'omitShortLags'; rmeddis@22: filteredSACFParams.criterionForOmittingLags=3; rmeddis@22: filteredSACFParams.plotACFsInterval=200; rmeddis@22: rmeddis@22: if filteredSACFParams.plotACFs rmeddis@22: % plot original waveform on ACF plot rmeddis@22: figure(13), clf rmeddis@22: subplot(4,1,1) rmeddis@22: ======= rmeddis@22: rmeddis@22: %% ACF plot if required rmeddis@22: if options.showACF rmeddis@22: tic rmeddis@22: method.dt=dt; rmeddis@22: method.segmentNo=1; rmeddis@22: method.nonlinCF=savedBFlist; rmeddis@22: rmeddis@22: minPitch= 80; maxPitch= 4000; numPitches=100; % specify lags rmeddis@22: pitches=10.^ linspace(log10(minPitch), log10(maxPitch),numPitches); rmeddis@22: pitches=fliplr(pitches); rmeddis@22: filteredSACFParams.lags=1./pitches; % autocorrelation lags vector rmeddis@22: filteredSACFParams.acfTau= .003; % time constant of running ACF rmeddis@22: filteredSACFParams.lambda= 0.12; % slower filter to smooth ACF rmeddis@22: filteredSACFParams.lambda= 0.01; % slower filter to smooth ACF rmeddis@22: rmeddis@22: filteredSACFParams.plotACFs=0; % special plot (see code) rmeddis@22: filteredSACFParams.plotFilteredSACF=0; % 0 plots unfiltered ACFs rmeddis@22: filteredSACFParams.plotMoviePauses=.3; % special plot (see code) rmeddis@22: rmeddis@22: filteredSACFParams.usePressnitzer=0; % attenuates ACF at long lags rmeddis@22: filteredSACFParams.lagsProcedure= 'useAllLags'; rmeddis@22: % filteredSACFParams.lagsProcedure= 'useBernsteinLagWeights'; rmeddis@22: % filteredSACFParams.lagsProcedure= 'omitShortLags'; rmeddis@22: filteredSACFParams.criterionForOmittingLags=3; rmeddis@22: filteredSACFParams.plotACFsInterval=200; rmeddis@22: rmeddis@22: if filteredSACFParams.plotACFs rmeddis@22: % plot original waveform on ACF plot rmeddis@22: figure(13), clf rmeddis@22: subplot(4,1,1) rmeddis@22: t=dt*(1:length(savedInputSignal)); rmeddis@22: plot(t,savedInputSignal) rmeddis@22: xlim([0 t(end)]) rmeddis@22: title(['stimulus: ' num2str(signalRMSdb, '%4.0f') ' dB SPL']); rmeddis@22: end rmeddis@22: rmeddis@22: % plot original waveform on summary/smoothed ACF plot rmeddis@22: figure(96), clf rmeddis@22: subplot(2,1,1) rmeddis@22: >>>>>>> 77bdf847c4be95767d3ea519d77878f9414a2686 rmeddis@22: t=dt*(1:length(savedInputSignal)); rmeddis@22: plot(t,savedInputSignal) rmeddis@22: xlim([0 t(end)]) rmeddis@22: title(['stimulus: ' num2str(signalRMSdb, '%4.0f') ' dB SPL']); rmeddis@22: rmeddis@22: rmeddis@22: % compute ACF rmeddis@22: switch saveAN_spikesOrProbability rmeddis@22: case 'probability' rmeddis@22: inputToACF=ANprobRateOutput.^0.5; rmeddis@22: otherwise rmeddis@22: inputToACF=ANoutput; rmeddis@22: end rmeddis@22: rmeddis@22: disp ('computing ACF...') rmeddis@22: [P, BFlist, sacf, boundaryValue] = ... rmeddis@22: filteredSACF(inputToACF, method, filteredSACFParams); rmeddis@22: disp(' ACF done.') rmeddis@22: rmeddis@22: % SACF rmeddis@22: subplot(2,1,2) rmeddis@22: imagesc(P) rmeddis@22: ylabel('periodicities (Hz)') rmeddis@22: xlabel('time (s)') rmeddis@22: title(['running smoothed (root) SACF. ' saveAN_spikesOrProbability ' input']) rmeddis@22: pt=[1 get(gca,'ytick')]; % force top xtick to show rmeddis@22: set(gca,'ytick',pt) rmeddis@22: set(gca,'ytickLabel', round(pitches(pt))) rmeddis@22: tt=get(gca,'xtick'); rmeddis@22: set(gca,'xtickLabel', round(100*t(tt))/100) rmeddis@22: end rmeddis@22: rmeddis@22: path(restorePath)