用 Matlab 绘制钙成像热图，并使用kmeans聚类

% 生成随机的神经元钙信号数据
data = readtable('suite2p_ROI_data_added_columns.csv','ReadVariableNames',true) % csv有header
cells = data(data.is_cell == 1,:)
cells_activity = cells(:, 2:1200)



% 使用k-means算法进行聚类
numClusters = 5;  % 聚类数量
data_t  = table2array(cells_activity);
[idx, centroids] = kmeans(data_t, numClusters,'Start', ones(5,1199));

class_name = string(1:numClusters);

[sortedResult, sortedIdx] = sort(idx);

sortedNeuronData = data_t(sortedIdx, :);


% 可视化聚类结果
fig=figure('Position',[100,100,800,800]);
axMain=axes('Parent',fig);
axMain.Position=[0.06,.1,.80,.80];
P=axMain.Position;
imagesc(axMain,sortedNeuronData);
axMain.YAxisLocation='left'; 
colormap(sky);
cb = colorbar('Location','eastoutside');
set(cb,"Position",cb.Position+[0.12,0,0,-0.6]);
xlabel('Frame');
ylabel('Neurons');
title('calcium imaging heatmap');

axBlockL=axes('Parent',fig);
axBlockL.Position=[0.86,P(2),P(3)/20,P(4)];
[X,Y]=SClusterBlock(sortedResult,'Orientation','left','Parent',axBlockL);
for i=1:length(X)
    text(axBlockL,X(i),Y(i),class_name(i),'FontSize',17,'HorizontalAlignment','center','FontName','Cambria')
end

% 绘制每个cluster 钙信号平均图
ColorList=...
    [0.5529    0.8275    0.7804
    1.00 0.07 0.65
    0.7451    0.7294    0.8549
    0.9843    0.5020    0.4471
    0.5020    0.6941    0.8275
    0.9922    0.7059    0.3843
    0.7020    0.8706    0.4118
    0.9882    0.8039    0.8980
    0.8510    0.8510    0.8510
    0.7373    0.5020    0.7412
    0.8000    0.9216    0.7725
    1.0000    0.9294    0.4353];

fig_signal = figure();
ax_signal = axes('Parent',fig_signal);
hold(ax_signal,"on");
for i_cluster = 1:numClusters
    % find indeces of traces where cluster label equals to Ki
    i_signal = data_t(labels == i_cluster,:);
    % calculate average activity trace for cluster Ki
    md  = mean(i_signal , 1);
    s = std(i_signal,1) 
    % ！注意，这里进行了标准化,方便看每类神经元的活动normalize so that each trace is between 0 and 1
    md= rescale(md)

  
    % plot
    xs = 1:length(md);

    %fill([xs fliplr(xs)], [md+i_cluster-1-s fliplr(md+i_cluster-1+s)],[0.2,0.2,0.2],"FaceAlpha",0.15,'EdgeColor','none');
    plot(xs,md+i_cluster-1,'Color',ColorList(i_cluster,:),'LineWidth',1.5)
end

ax_signal.YTick=[];
ax_signal.YLim=[0,numClusters];
title('Average activity of each cluster');
xlabel("Frame");
for i_cluster = 1:numClusters
    text(-46.5, i_cluster-0.5, ['c',num2str(i_cluster)]);
end