-
Notifications
You must be signed in to change notification settings - Fork 0
/
Pkcorr3.m
290 lines (241 loc) · 6.65 KB
/
Pkcorr3.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
function output=Pkcorr(Sig,Pos,percent,skp)
% Imput:
% Sig= is an array of cell signals each column is one cell, rows are time points,
% Pos= 2-D array of x, y positions of each cell. order of cell positions needs to be
% the order of cells in columns of Sig. column 1 = X, Column 2= y.
% skp= is the number of frames that will be measured after each event to
% find neighboring evens, and plot the activity trigger average of the
% active cell and its neighbor. Note: The last skp frames of the data are
% disregarded.
Ncell=length(Sig(1,:));
Nfrm=length(Sig(:,1));
disp('begin Pkcorr3, creating intensity and peak graphs');
for x=1:Ncell % normalize signals to delta F/Fo, Fo is the minium value for that cell
m=min(Sig(:,x));
Sig(:,x)=(Sig(:,x)-m)/m;
end
rng=max(max(Sig)); %generate heatmap of cell activity, each row is activity of one cell
%HeatMap(Sig','colormap',jet(1000),'DisplayRange',rng,'Symmetric',false); %1000
%figure
%plot(mean(Sig,1))
%title('mean intensity per cell')
%xlabel('Cells')
%ylabel('Intensity (F/Fo)')
%figure
%plot(sum(Sig,1))
%title('summed intensity per cell')
%figure
%plot(sum(Sig,2))
%title('summed intensity over time')
highsig = zeros(Nfrm,Ncell);
lowsig = zeros(Nfrm,Ncell);
PK=zeros(Nfrm,Ncell);
highPK=zeros(Nfrm,Ncell);
lowPK=zeros(Nfrm,Ncell);
%Npt=40; %500
Npt = 15;
Aucor=zeros(Npt+1,Ncell);
lowcount = 1;
highcount = 1;
for x=1:Ncell % find peaks in each trace and setup locations to 1 in PK
[pks,locs]=findpeaks(Sig(:,x),'MinPeakProminence',percent); %1, set lower for more events
PK(locs,x)=1;
if sum(PK(:,x)) < 2;
lowsig(:,lowcount) = Sig(:,x);
lowPK(:,lowcount) = PK(:,x);
lowcount = lowcount + 1;
end
if sum(PK(:,x)) > 5;
highsig(:,highcount) = Sig(:,x);
highPK(:,highcount) = PK(:,x);
highcount = highcount + 1;
end
if sum(PK(:,x)) > 10 %if excessive events are detected in one cell, the cell is considered null.
PK(:,x)= 0;
end
Aucor(:,x)=autocorr(Sig(:,x),Npt);
end
figure
set(gcf,'position',[50,50,1100,420])
subplot(1,2,1);
plot(Sig)
hold on
plot(mean(Sig,2),'k','LineWidth',3)
title('Intensity over time')
xlabel('Frames')
ylabel('Intensity (F/Fo)')
subplot(1,2,2);
plot(mean(Sig,2))
title('mean intensity over time')
xlabel('Frames')
ylabel('Intensity (F/Fo)')
figure
set(gcf,'position',[50,50,1100,420])
subplot(1,2,1);
plot(lowsig)
hold on
plot(mean(lowsig,2),'k','LineWidth',3)
title('lowIntensity over time')
xlabel('Frames')
ylabel('Intensity (F/Fo)')
subplot(1,2,2);
plot(mean(lowsig,2))
title('mean lowintensity over time')
xlabel('Frames')
ylabel('Intensity (F/Fo)')
figure
set(gcf,'position',[50,50,1100,420])
subplot(1,2,1);
plot(highsig)
hold on
plot(mean(highsig,2),'k','LineWidth',3)
title('highIntensity over time')
xlabel('Frames')
ylabel('Intensity (F/Fo)')
subplot(1,2,2);
plot(mean(highsig,2))
title('mean highintensity over time')
xlabel('Frames')
ylabel('Intensity (F/Fo)')
figure
plot(mean(Sig,2))
hold on
plot(mean(highsig,2))
hold on
plot(mean(lowsig,2))
hold off
title('mean intensity over time comparison')
xlabel('Frames')
ylabel('Intensity (F/Fo)')
legend('Mean Intensity','Mean High Signaling Intensity', 'Mean Low Signaling Intensity')
figure
subplot(2,1,1);
PK1=sum(PK,1);
plot(PK1,'r','LineWidth',1)
title('sum Peaks per cell');
subplot(2,1,2);
PK1=sum(PK,2);
plot(PK1,'r','LineWidth',1)
title('sum Peaks per time');
detectedevents=[];
for x=1:Nfrm
if x > 1;
detectedevents(x) = PK1(x) + detectedevents(x-1);
else
detectedevents(x) = PK1(x);
end
end
figure
plot(detectedevents)
title('detected events over time');
%%
%{
figure
PK2=mean(PK);
plot(PK2,'g','LineWidth',3)
title('mean PK');
%}
%disp('for loop 2 done');
%{
plt=PK; %plot indentified events, each column is one cell
for x=1:Ncell
plt(:,x)=plt(:,x)+x;
end
figure
plot(plt)
title('detected events');
%}
%figure
%sumplt = sum(plt,2);
%sumplt = sumplt - min(sumplt);
%plot(sumplt)
%title('sum detected events over time');
%{
figure
plt2 = sum(plt)
plot(plt2)
title('sum detected events');
plt1=PK1; %plot indentified events, each column is one cell
for x=1:Ncell
plt1(:,x)=plt1(:,x)+x;
end
figure
Plt2=plt1;
plot(Plt2,'b','LineWidth',3)
title('Plt from sum PK1');
plt2=PK2; %plot indentified events, each column is one cell
for x=1:Ncell
plt2(:,x)=plt2(:,x)+x;
end
figure
Plt3=plt2;
plot(Plt3,'b','LineWidth',3)
title('Plt from mean PK2');
%}
%{
Dis=zeros(Ncell, Ncell);
NxtCeldis=85;
for r=1:Ncell % calculate distance between cells if dis< 85 pixel then flag as neighbor and set to 1
for c=1:Ncell
x=Pos(r,1)-Pos(c,1);
y=Pos(r,2)-Pos(c,2);
if c~=r
if sqrt((x*x)+(y*y))<NxtCeldis
Dis(r,c)=1;
end
end
end
end
%disp('for loop 3 done');
NBPK=0;
%skp=100; %skp is the number of frames to look at the cross correlations (i.e. how many frames to plot after an "event"
NBtrc=zeros(1,skp+1);
Autrc=zeros(1,skp+1);
for x=1:Ncell
locs=find(PK((1:Nfrm-skp),x));
NeBrs=find(Dis(:,x));
for l=1:length(locs)
if sum(sum(PK((locs(l):locs(l)+skp),NeBrs)))>0
NBPK=NBPK+1; %count if NBPK are active after a cell event
end
Autrc=cat(1,Autrc,Sig((locs(l):locs(l)+skp),x)'); %build array of autocorrelation activity of the active cell
for nb=1:length(NeBrs)
Sig((locs(l):locs(l)+skp),NeBrs(nb))';
NBtrc=cat(1,NBtrc,Sig((locs(l):locs(l)+skp),NeBrs(nb))'); %build array of autocorrelation activity of the active cell
end
end
end
disp('for loop 4 done');
%{
figure %plot event triggered average traces of neighbors
plot(NBtrc')
hold on
NBtrc=mean(NBtrc);
plot(NBtrc'*10,'r','LineWidth',3)
title('event triggered average trace of neighbors');
figure %plot event triggered average traces of that cell
plot(Autrc')
hold on
Autrc=mean(Autrc);
plot(Autrc'*10,'k','LineWidth',3)
title('event triggered average trace of that cell');
figure % plot together
plot(NBtrc','r','LineWidth',3)
hold on
plot(Autrc','k','LineWidth',3)
%}
%figure
%plot(Aucor)
%hold on
%Aucor=mean(Aucor');
%plot(Aucor','k','LineWidth',3)
disp('all figures done');
mnsig=mean(mean(Sig))*10;
totPK=sum(sum(PK((1:Nfrm-skp),:)));
NBPK;
prbNB=NBPK/sum(sum(PK((1:Nfrm-skp),:)));
rate=totPK/(Ncell*(Nfrm-skp));
%output=cat(1,NBtrc,Autrc);
%}
disp('all figures done');
disp ('Pkcorr3 done');