function [pred_output]=BC_pred(STA_output) %[pred_output]=BC_pred(STA_output) %this function shows the prinicpal calculation for the the prediction to %the white noise stimulus with the linear-nonlinear model (LN-model)for %continous voltage signals of bipolar cells (BC). %Inputs: % STA_output = output of the BC_STA function; containing the % STA and stimulus signal % %Outpus: pred_output = structure containing the output of the % prediction % % % Note: the function is built only for cells recorded without frozen frames. % But additional comments are made for cells with frozen frames. % code by HS, last modified March 2021 %------------------------------------------------------------------------------- %% 1. recompute STA and nonlinearity without the stimulus signal that we want to predict %%Important we want to predict a stimulus sequence that was not used to %%build the model (filter and nonlinearity). Therefore we have to redo the %%filter(STA) and nonlinearity without the stimulus-sequence that we want %to predict the response(see Methods of Schreyer&Gollisch,2021). %In the Manuscript we predict 200 sequences of 300 frames long that where %randomly chosen (see Method section). Every time we rebuilt the STA and %NL. The overall prediction performance is then the average R^2 over all %200 sequences. Here, for simplicity I only show the prediction to one %sequence. Note however, to get a good estimate of how well the LN-model can predict %the response of the cell it is important to run it over multiple %sequences. %If we recorded Frozen frames, we predict the response to the frozen %frames from the model (filter and nonlinearity) we build from the running %frames. %% 1.1 STA %to keep the code easy: I here show the principal idea, by predicting the %response to the last 300frames of the white noise stimulus. %(in the manuscript we randomly choose 200 sequences and redo this part 200 times) %take out the stimulus sequence you want to predict stimPred=STA_output.stimulus2D_zoom(:,end-299:end); %the last 300frames %voltage trace to predict voltPred=STA_output.membranPotAVG(end-299:end); %the response to the last 300 frames %build the STA without the last 300 frames [STA_norm]=buildSTA ... (STA_output.stimulus2D_zoom(:,1:end-300),... STA_output.membranPotAVG(1:end-300),STA_output.STA_fnbr); STA3D_norm=reshape(STA_norm,size(STA_output.rangeY_zoom,2),size(STA_output.rangeX_zoom,2), ... size(STA_norm,2));% Y-zoom checkers (monitor axis), X-zoom checkers (monitor axis) and Z (STA_fnbr) %% 1.2 Nonlinearity %build the NL from the data without the last 300 frames generator_signal=filter2(STA_norm,STA_output.stimulus2D_zoom(:,1:end-300),'valid'); %see also Point 7 %in the function BC_STA.m and BC_NL of how you do it when you have a best %spatial and temporal component and reduce noice by cutting 3-sigma %contour around the receptive field region. %get the corresponding membrane potential responses y_axisNL_unsorted=STA_output.membranPotAVG(STA_output.STA_fnbr:end-300); %the first value is at the length of the filter %sort the x axis [NL_xAxis,indx]=sort(generator_signal); %sort the corresponding y axis NL_yAxis=y_axisNL_unsorted(indx); %do a binning with percentile to have same amount of data in each bin %make 40 bins nbrbins=40; bins= doBin(NL_xAxis,NL_yAxis,nbrbins); %------------------------------------------------------------------------------- %% 2. compute the prediction %% 2.1 convolution %convolve the new stimulus sequence (not used for prediciton) with the %normalized STA stim2pred=filter2(STA_norm,stimPred,'valid'); %% 2.2 compute the prediction through linear interpolation NL_x=bins.NL_xbin; NL_y=bins.NL_ybin; %example of how to make the prediction through linear interpolation predMemP_y=nan(1,size(stim2pred,2)); for kk=1:size(stim2pred,2) nbrTemp=stim2pred(kk); %first x-axis value to predict [xIndxS,xIndxL]= getIndex(NL_x,nbrTemp,1); %get those values that are on left and rigth xValueS=NL_x(xIndxS); xValueL=NL_x(xIndxL); %do the linear interpolation [polyCoeff] = polyfit([xValueS;xValueL],[NL_y(xIndxS);... NL_y(xIndxL)],1); predMemP_y(kk)=polyCoeff(1)*nbrTemp+polyCoeff(2); end %------------------------------------------------------------------------------- %% 3. evaluate the prediction performance %pearson correlation corrcoeff_temp=corrcoef(voltPred(STA_output.STA_fnbr:end),predMemP_y); %the first value is at the length of the filter pearCorrSqrt=corrcoeff_temp(2)^2; %------------------------------------------------------------------------------- %% 4. plot the output figure; %STA subplot(2,2,1) [~,indx_Zoom]=max(abs(STA3D_norm(:))); %max intensity value of STA [~,~,posZ]=ind2sub(size(STA3D_norm),indx_Zoom); %get the frame position imagesc(STA3D_norm(:,:,posZ)) colormap('gray'); %how the sitmulus was shown axis equal; axis tight; title('STA best frame') %NL subplot(2,2,2) plot(NL_xAxis,NL_yAxis,'.','color','k') %plot non-binned signal hold on; plot(bins.NL_xbin,bins.NL_ybin,'o','markerfacecolor','r','markeredgecolor','r') %plot binned signal title('NL') axis tight xlabel('generator signal') ylabel('voltage (mV)') %prediction subplot(2,2,3:4) plot(voltPred(STA_output.STA_fnbr:end),'k') %original membrane potential hold on; plot(predMemP_y,'r')%predicted membrane potential title(['pCorr:' mat2str(round(pearCorrSqrt,2)) ' Filename: ' STA_output.filename]) xlabel('stimulus frames') ylabel('voltage (ms)') legend('original response','predicted response') %------------------------------------------------------------------------------- %% 5. add variables into the output pred_output.pearCorrSqrt=pearCorrSqrt; end function [IndS,IndL]= getIndex(x,valConv,Points) IndS=find(x=valConv,Points); if isempty(IndS)|| length(IndS)=valConv,Points*2); IndS=IndL(1); IndL=IndL(2:end); else IndL=find(x>=valConv,Points*2-length(IndS)); end elseif isempty(IndL)|| length(IndL)