Specify the subplot and hold on at the starting of for loop and hold off at the end of for loop.If size(mean_fr,1) is assumed to be 10,consider subplot(2,5,c);
fit_thresh = .25 ; %set a threshold for the fit of the linear regression, just randomly picked this number, could be more stringent
%could change the fit_thresh
for c = 1:size(mean_fr,1)
subplot(2,5,c);
hold on;
%Then do a simple linear regression, look for a good fit and a
%positive slope
x = 1:numel(smoothdata(bin_spikes(c,:),'movmean',5));x=x'; %X in this case is just time (so 1 : number of bins)
X = [ones(length(x),1) x];
y = smoothdata(bin_spikes(c,:),'movmean',5); y=y';
b = X\y ; % can use "mldividde" instead of "\"
yfit = X *b ; %the line that is fit
%assess fit using an R^2 measure
r = 1 - sum((y - yfit).^2)/sum((y - mean(y)).^2) ; %goes from 0 to 1, with 1 being perfect fit
if r > fit_thresh
if b > 0 %looking for positive slope, could eventually put a threshold here
thisfig = figure();
plot(smoothdata(bin_spikes(c,:),'movmean',5)) %look for positive correlation, on the smoothened, binned data
title(['Channel ' num2str(recording_info.channel_numbers(c)) ' Area ' recording_info.area{c} ])
fig = gca;
fig.XTickLabel = num2cell(fig.XTick *10);
fig.XLabel.String = 'Time (ms) after sample offset';
fig.YLabel.String = 'FR (spikes/s)';
plot (x, yfit, 'r-.');
end
I hope this works fine.
