SFig7B.m 10 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318
  1. close all;
  2. clear;
  3. % plot spiekcount and spontaneous for paired data set
  4. load('sp_context_pair.mat')
  5. dur_echo=1.3311;
  6. dur_distress=1.9266;
  7. c_=sp_context(:,2)==0; % unit classification
  8. aw_=sp_context(:,3)==0; % awake data
  9. particular1=logical(aw_.*c_);
  10. an_units=sp_context(particular1,1)+1; % anesthetized units
  11. particular2=ismember(sp_context(:,1),unique(an_units));
  12. spont_=sp_context(:,7)==0; % spontaneous firing rate
  13. cont_=sp_context(:,7)==1; % context firing rate intact
  14. % submatrices
  15. aw_spont=sp_context(logical(particular1.*spont_),:);
  16. an_spont=sp_context(logical(particular2.*spont_),:);
  17. aw_cont=sp_context(logical(particular1.*cont_),:);
  18. an_cont=sp_context(logical(particular2.*cont_),:);
  19. %% calculate average firing rate from spike counts across trials
  20. % awake data - spontaneous & context
  21. units=unique(aw_spont(:,1));
  22. aw_spontR=nan(numel(units),1);
  23. aw_contR=nan(numel(units),2);
  24. for u=1:numel(units)
  25. u_=aw_spont(:,1)==units(u);
  26. u_2=aw_cont(:,1)==units(u);
  27. rate=mean(aw_spont(u_,8));
  28. rate2=mean(aw_cont(u_2,5)); % nav
  29. rate3=mean(aw_cont(u_2,6)); % com
  30. aw_spontR(u,1)=rate;
  31. aw_contR(u,1)=rate2;
  32. aw_contR(u,2)=rate3;
  33. end
  34. aw_spontR=aw_spontR./((3.5-0.2));
  35. aw_contR(:,1)=aw_contR(:,1);
  36. aw_contR(:,2)=aw_contR(:,2);
  37. % anesthetized data - spontaneous & context
  38. units=unique(an_spont(:,1));
  39. an_spontR=nan(numel(units),1);
  40. an_contR=nan(numel(units),2);
  41. for u=1:numel(units)
  42. u_=an_spont(:,1)==units(u);
  43. u_2=an_cont(:,1)==units(u);
  44. rate=mean(an_spont(u_,8));
  45. rate2=mean(an_cont(u_2,5)); % nav
  46. rate3=mean(an_cont(u_2,6)); % com
  47. an_spontR(u,1)=rate;
  48. an_contR(u,1)=rate2;
  49. an_contR(u,2)=rate3;
  50. end
  51. an_spontR=an_spontR./((3.5-0.2));
  52. an_contR(:,1)=an_contR(:,1);
  53. an_contR(:,2)=an_contR(:,2);
  54. y=[aw_spontR,an_spontR,aw_contR(:,1),an_contR(:,1),aw_contR(:,2),an_contR(:,2)];
  55. % firing rate in response to probes
  56. load('pair_data.mat');
  57. c_=all_di(:,2)==0;
  58. aw_=all_di(:,3)==0;
  59. part=logical(c_.*aw_);
  60. an_=[0; part];
  61. an_(end)=[];
  62. probe_aw=all_di(part,18:19).*20;
  63. probe_an=all_di(logical(an_),18:19).*20;
  64. % compare paired samples
  65. pv1=signrank(aw_spontR,an_spontR);
  66. pv2=signrank(aw_contR(:,1),an_contR(:,1));
  67. pv3=signrank(aw_contR(:,2),an_contR(:,2));
  68. pv4=signrank(aw_contR(:,1),aw_contR(:,2));
  69. pv5=signrank(an_contR(:,1),an_contR(:,2));
  70. pv6=signrank(probe_aw(:,1),probe_an(:,1));
  71. pv7=signrank(probe_aw(:,2),probe_an(:,2));
  72. pv8=signrank(probe_aw(:,1),probe_aw(:,2));
  73. pv9=signrank(probe_an(:,1),probe_an(:,2));
  74. %% plot spontaneous and firing during context in same plot
  75. red =[0.8500,0.3250, 0.0980];
  76. blue=[0,0.4470,0.7410];
  77. red_light=[1,0.5,0.5].*red;
  78. blue_light=[1,0.5,0.5].*blue;
  79. lines_color=[0.7 0.7 0.7]; % light gray
  80. bar_colors=[[0.3 0.3 0.3];[0.5 0.5 0.5];blue;red;blue_light;red_light];
  81. x = 1:6;
  82. figure(1); set(gcf,'Position',[400 200 250 170])
  83. ax = axes();
  84. hold(ax);
  85. for i=1:6
  86. h=boxchart(x(i)*ones(size(y(:,i))), y(:,i), 'BoxFaceColor', bar_colors(i,:),'Notch','off');
  87. h.MarkerStyle='.';
  88. h.MarkerColor=bar_colors(i,:);
  89. end
  90. % add lines per unit
  91. hold on
  92. for i=1:numel(units)
  93. l=plot([1,2], [aw_spontR(i), an_spontR(i)],'color',lines_color);
  94. uistack(l,'bottom')
  95. l=plot([3,4], [aw_contR(i,1), an_contR(i,1)],'color',lines_color);
  96. uistack(l,'bottom')
  97. l=plot([5,6], [aw_contR(i,2), an_contR(i,2)],'color',lines_color);
  98. uistack(l,'bottom')
  99. end
  100. hold off
  101. xlim([0.5 6.5])
  102. ylabel('# spikes')
  103. xticklabels({'awake','anesthetized','awake','anesthetized','awake','anesthetized'});
  104. set(gca, 'box', 'off')
  105. set(gca, 'Color','none')
  106. set(gca,'linewidth',1);set(gca,'fontsize',8);
  107. xt = get(gca, 'XTick');
  108. line1=180;
  109. line2=200;
  110. line3=220;
  111. sl=0.75;
  112. % plot lines for comparison
  113. hold on
  114. plot(xt(1:2), [1 1]*line1, '-k','LineWidth',sl)
  115. plot([1;1]*xt(1),[line1,line1-3],'-k', 'LineWidth',sl);
  116. plot([1;1]*xt(2),[line1,line1-3],'-k', 'LineWidth',sl);
  117. plot(xt(3:4), [1 1]*line1, '-k','LineWidth',sl)
  118. plot([1;1]*xt(3),[line1,line1-3],'-k', 'LineWidth',sl);
  119. plot([1;1]*xt(4),[line1,line1-3],'-k', 'LineWidth',sl);
  120. plot(xt(5:6), [1 1]*line1, '-k','LineWidth',sl)
  121. plot([1;1]*xt(5),[line1,line1-3],'-k', 'LineWidth',sl);
  122. plot([1;1]*xt(6),[line1,line1-3],'-k', 'LineWidth',sl);
  123. % plot([xt(3) xt(5)], [1 1]*line2, '-k','LineWidth',sl)
  124. % plot([1;1]*xt(3),[line2,line2-3],'-k', 'LineWidth',sl);
  125. % plot([1;1]*xt(5),[line2,line2-3],'-k', 'LineWidth',sl);
  126. %
  127. % plot([xt(4) xt(6)], [1 1]*line3, '-k','LineWidth',sl)
  128. % plot([1;1]*xt(4),[line3,line3-3],'-k', 'LineWidth',sl);
  129. % plot([1;1]*xt(6),[line3,line3-3],'-k', 'LineWidth',sl);
  130. ylim([-5 180])
  131. % plot level of significance
  132. if pv1 <= 0.05 && pv1 > 0.01
  133. text(mean(xt(1:2))-0.05, line1+5, '*','HorizontalAlignment','center')
  134. elseif pv1 <=0.01 && pv1 > 0.001
  135. text(mean(xt(1:2))-0.05, line1+5, '**','HorizontalAlignment','center')
  136. elseif pv1 <= 0.001
  137. text(mean(xt(1:2))-0.05, line1+5, '***','HorizontalAlignment','center')
  138. else
  139. text(mean(xt(1:2))-0.05, line1+12, 'ns','FontSize',7,'HorizontalAlignment','center')
  140. end
  141. if pv2 <= 0.05 && pv2 > 0.01
  142. text(mean(xt(3:4))-0.05, line1+5, '*','HorizontalAlignment','center')
  143. elseif pv2 <=0.01 && pv2 > 0.001
  144. text(mean(xt(3:4))-0.05, line1+5, '**','HorizontalAlignment','center')
  145. elseif pv2 <= 0.001
  146. text(mean(xt(3:4))-0.05, line1+5, '***','HorizontalAlignment','center')
  147. else
  148. text(mean(xt(3:4))-0.05, line1+12, 'ns','FontSize',7,'HorizontalAlignment','center')
  149. end
  150. if pv3 <= 0.05 && pv3 > 0.01
  151. text(mean(xt(5:6))-0.05, line1+5, '*','HorizontalAlignment','center')
  152. elseif pv3 <=0.01 && pv3 > 0.001
  153. text(mean(xt(5:6))-0.05, line1+5, '**','HorizontalAlignment','center')
  154. elseif pv3 <= 0.001
  155. text(mean(xt(5:6))-0.05, line1+5, '***','HorizontalAlignment','center')
  156. else
  157. text(mean(xt(5:6))-0.05, line1+12, 'ns','FontSize',7,'HorizontalAlignment','center')
  158. end
  159. % if pv4 <= 0.05 && pv4 > 0.01
  160. % text(mean(xt(3:5))-0.05, line2+2, '*','HorizontalAlignment','center')
  161. % elseif pv4 <=0.01 && pv4 > 0.001
  162. % text(mean(xt(3:5))-0.05, line2+2, '**','HorizontalAlignment','center')
  163. % elseif pv4 <= 0.001
  164. % text(mean(xt(3:5))-0.05, line2+2, '***','HorizontalAlignment','center')
  165. % else
  166. % text(mean(xt(3:5))-0.05, line2+8, 'ns','FontSize',7,'HorizontalAlignment','center')
  167. % end
  168. %
  169. % if pv5 <= 0.05 && pv5 > 0.01
  170. % text(mean(xt(4:6))-0.05, line3+2, '*','HorizontalAlignment','center')
  171. % elseif pv5 <=0.01 && pv5 > 0.001
  172. % text(mean(xt(4:6))-0.05, line3+2, '**','HorizontalAlignment','center')
  173. % elseif pv5 <= 0.001
  174. % text(mean(xt(4:6))-0.05, line3+2, '***','HorizontalAlignment','center')
  175. % else
  176. % text(mean(xt(4:6))-0.05, line3+8, 'ns','FontSize',7,'HorizontalAlignment','center')
  177. % end
  178. exportgraphics(gcf,'E:\Users\User\Desktop\delay paper\anesthetized\figures_paper\sameunit_inj_rate.pdf',...
  179. 'Resolution',300','ContentType','vector','BackgroundColor','none')
  180. %% plotting probes
  181. colors=[blue_light;blue;red_light;red];
  182. y2=nan(58,4);
  183. y2(1:numel(probe_aw(:,1)),1)=probe_aw(:,1);
  184. y2(1:numel(probe_an(:,1)),2)=probe_an(:,1);
  185. y2(1:numel(probe_aw(:,2)),3)=probe_aw(:,2);
  186. y2(1:numel(probe_an(:,2)),4)=probe_an(:,2);
  187. x = 1:4;
  188. figure(2); set(gcf,'Position',[700 200 170 170])
  189. ax = axes();
  190. hold(ax);
  191. for i=1:4
  192. h=boxchart(x(i)*ones(size(y2(:,i))), y2(:,i), 'BoxFaceColor', colors(i,:),'Notch','off');
  193. h.MarkerStyle='.';
  194. h.MarkerColor=colors(i,:);
  195. end
  196. % add lines per unit
  197. hold on
  198. for i=1:numel(units)
  199. l=plot([1,2], [probe_aw(i,1), probe_an(i,1)],'color',lines_color);
  200. uistack(l,'bottom')
  201. l=plot([3,4], [probe_aw(i,2), probe_an(i,2)],'color',lines_color);
  202. uistack(l,'bottom')
  203. end
  204. hold off
  205. xlim([0.5 4.5])
  206. ylabel('firing rate ')
  207. xticks(1:4)
  208. xticklabels({'awake','anesthetized','awake','anesthetized'});
  209. set(gca, 'box', 'off')
  210. set(gca, 'Color','none')
  211. set(gca,'linewidth',1);set(gca,'fontsize',8);
  212. xt = get(gca, 'XTick');
  213. line1=220;
  214. % line2=340;
  215. % line3=390;
  216. sl=0.75;
  217. % plot lines for comparison
  218. hold on
  219. plot(xt(1:2), [1 1]*line1, '-k','LineWidth',sl)
  220. plot([1;1]*xt(1),[line1,line1-10],'-k', 'LineWidth',sl);
  221. plot([1;1]*xt(2),[line1,line1-10],'-k', 'LineWidth',sl);
  222. plot(xt(3:4), [1 1]*line1, '-k','LineWidth',sl)
  223. plot([1;1]*xt(3),[line1,line1-10],'-k', 'LineWidth',sl);
  224. plot([1;1]*xt(4),[line1,line1-10],'-k', 'LineWidth',sl);
  225. % plot([xt(1) xt(3)], [1 1]*line2, '-k','LineWidth',sl)
  226. % plot([1;1]*xt(1),[line2,line2-10],'-k', 'LineWidth',sl);
  227. % plot([1;1]*xt(3),[line2,line2-10],'-k', 'LineWidth',sl);
  228. %
  229. % plot([xt(2) xt(4)], [1 1]*line3, '-k','LineWidth',sl)
  230. % plot([1;1]*xt(2),[line3,line3-10],'-k', 'LineWidth',sl);
  231. % plot([1;1]*xt(4),[line3,line3-10],'-k', 'LineWidth',sl);
  232. ylim([-5 240])
  233. % plot level of significance
  234. if pv6 <= 0.05 && pv6 > 0.01
  235. text(mean(xt(1:2))-0.05, line1+10, '*','HorizontalAlignment','center')
  236. elseif pv6 <=0.01 && pv6 > 0.001
  237. text(mean(xt(1:2))-0.05, line1+10, '**','HorizontalAlignment','center')
  238. elseif pv6 <= 0.001
  239. text(mean(xt(1:2))-0.05, line1+10, '***','HorizontalAlignment','center')
  240. else
  241. text(mean(xt(1:2))-0.05, line1+18, 'ns','FontSize',8,'HorizontalAlignment','center')
  242. end
  243. if pv7 <= 0.05 && pv7 > 0.01
  244. text(mean(xt(3:4))-0.05, line1+10, '*','HorizontalAlignment','center')
  245. elseif pv7 <=0.01 && pv7 > 0.001
  246. text(mean(xt(3:4))-0.05, line1+10, '**','HorizontalAlignment','center')
  247. elseif pv7 <= 0.001
  248. text(mean(xt(3:4))-0.05, line1+10, '***','HorizontalAlignment','center')
  249. else
  250. text(mean(xt(3:4))-0.05, line1+18, 'ns','FontSize',8,'HorizontalAlignment','center')
  251. end
  252. % if pv10 <= 0.05 && pv10 > 0.01
  253. % text(mean(xt(1:3))-0.05, line2+20, '*','HorizontalAlignment','center')
  254. % elseif pv10 <=0.01 && pv10 > 0.001
  255. % text(mean(xt(1:3))-0.05, line2+20, '**','HorizontalAlignment','center')
  256. % elseif pv10 <= 0.001
  257. % text(mean(xt(1:3))-0.05, line2+20, '***','HorizontalAlignment','center')
  258. % else
  259. % text(mean(xt(1:3))-0.05, line2+28, 'ns','FontSize',8,'HorizontalAlignment','center')
  260. % end
  261. %
  262. % if pv11 <= 0.05 && pv11 > 0.01
  263. % text(mean(xt(2:4))-0.05, line3+20, '*','HorizontalAlignment','center')
  264. % elseif pv11 <=0.01 && pv11 > 0.001
  265. % text(mean(xt(2:4))-0.05, line3+20, '**','HorizontalAlignment','center')
  266. % elseif pv11 <= 0.001
  267. % text(mean(xt(2:4))-0.05, line3+20, '***','HorizontalAlignment','center')
  268. % else
  269. % text(mean(xt(2:4))-0.05, line3+28, 'ns','FontSize',8,'HorizontalAlignment','center')
  270. % end
  271. exportgraphics(gcf,'E:\Users\User\Desktop\delay paper\anesthetized\figures_paper\sameunit_rate_probe.pdf',...
  272. 'Resolution',300','ContentType','vector','BackgroundColor','none')