OSCR

A modular multi-color fluorescence microscope for simultaneous tracking of cellular activity and behavior.

Code ↔ Paper

6 matches between paragraphs of the paper and lines of its authors' code, computed by the harvester (lexical-v1). Click a colored paragraph or line to see its counterpart.

The 6 matches
  1. [1] § Methods › Statistical analysis of calcium imaging ↔ TRN_MixedLinearModel.ipynb, lines 78–91 · score 0.64 · mixed linear model, post, calcium, worms, reversals, signal
  2. [2] § Methods › Gait identification of tardigrades ↔ tardigrade_cluster.ipynb, lines 299–323 · score 0.62 · pairwise distance, UMAP, embedding, clustering, tardigrades
  3. [3] § Methods › Statistical analysis of prey signal ↔ GFP_preycontact.ipynb, lines 243–299 · score 0.59 · biting events, bite onset, baseline, ratio, GFP, signal
  4. [4] § Methods › Pose estimation of tardigrades ↔ tardigrade_cluster.ipynb, lines 128–152 · score 0.59 · stance duration, duty factor, swing, tardigrades
  5. [5] § Results › Brightfield tracking of tardigrades ↔ tardigrade_cluster.ipynb, lines 336–365 · score 0.50 · silhouette score, embedding, cycles, clustering, leg, tardigrades
  6. [6] § Results › Ratiometric imaging of single neurons in moving animals ↔ Dros_fig2.ipynb, lines 14–137 · score 0.50 · motion correction, ratiometric imaging, R0, tracked

Paper

Loaded from Europe PMC by your browser, not stored by OSCR: doi.org · Europe PMC

The paper is loaded when this pane is shown.

The authors' code

Jupyter notebook · 602 lines · 24 KB · no license · 3 matches

  1. # %%
  2. import pandas as pd
  3. import numpy as np
  4. from scipy.signal import savgol_filter
  5. import os
  6. import yaml
  7. import matplotlib.pyplot as plt
  8. from tslearn import metrics, barycenters
  9. import matplotlib.pyplot as plt
  12. from tslearn.preprocessing import TimeSeriesScalerMeanVariance, \
  13. TimeSeriesResampler
  15. from scipy.signal import find_peaks
  16. import tardigrade_functions as tg
  18. from sklearn.cluster import AgglomerativeClustering
  19. from sklearn.metrics import silhouette_score
  20. from matplotlib.gridspec import GridSpec
  21. from scipy import stats
  22. from umap import UMAP
  26. # %%
  27. # settings for dataset
  28. date = 20251209
  29. home = os.path.expanduser("~")
  30. expID = 'LB016'
  31. inpath = os.path.join(home, f'data/{expID}/out_{date}')
  32. pose_path = os.path.join(home, f'data/{expID}/poseestimation')
  33. pose_file = 'DLC_Resnet50_TardiLocoOct23shuffle1_snapshot_200.csv'
  35. # %%
  36. # load config
  37. config_path = "config.yaml"
  38. config = yaml.safe_load(open(config_path, "r"))
  40. fps = config['video']['fps']
  42. bodyparts = config['analysis']['bodyparts']
  43. limbs = config['analysis']['limbs']
  44. limbs_ord = config['analysis']['limbs_ord']
  45. pairs = config['analysis']['pairs']
  46. bp_pairs = config['analysis']['bp_pairs']
  48. bp_color_dict = config['color']['bp_color']
  49. roi_defined = config['analysis']['roi_defined']
  51. idx = pd.IndexSlice
  53. # %%
  54. Trkangle_all = pd.DataFrame([])
  55. CMSvelo_all = pd.DataFrame([])
  56. CLbp_all = pd.DataFrame([])
  57. swing_all = pd.DataFrame([])
  58. turns_all = pd.DataFrame([])
  59. straight_all = pd.DataFrame([])
  60. lastmaxidx = 0
  61. for f in os.listdir(inpath):
  62. if expID in f and not '.' in f:
  64. # load angle of tardigrade track
  65. fangle = pd.read_csv(os.path.join(inpath,f,f'{f}_trackangle.csv'), index_col=0).fillna(0)
  66. _ = savgol_filter(fangle, fps, 2, axis=0)
  67. fangle_smooth = pd.DataFrame(_, columns=fangle.columns, index=fangle.index)
  68. fangle_smooth.index = pd.MultiIndex.from_product([[f],fangle_smooth.index])
  69. Trkangle_all = pd.concat([Trkangle_all,fangle_smooth])
  71. # load cms coordinates
  72. cms = pd.read_csv(os.path.join(inpath,f,f'{f}_stage.csv'),index_col=0).loc[:,['Xcms','Ycms']]
  73. cmsvelo = np.sqrt(cms.iloc[:,0].diff(5)**2 + cms.iloc[:,1].diff(5)**2).interpolate(limit_direction='both')/5*fps
  74. cmsvelo.index = pd.MultiIndex.from_product([[f],cmsvelo.index])
  75. CMSvelo_all = pd.concat([CMSvelo_all,cmsvelo.to_frame()])
  77. # load leg position on centerline
  78. f_clbp = pd.read_csv(os.path.join(inpath,f,f'{f}_CLbp.csv'), index_col=0, header=[0,1,2])
  79. f_clbp_norm = (f_clbp - f_clbp.mean(axis=0)) / f_clbp.std(axis=0)
  80. f_clbp_x = f_clbp_norm.loc[:,idx[:,:,'x']].droplevel(level=[0,2],axis=1).fillna(0)
  81. CLbp_all = pd.concat([CLbp_all,f_clbp_x])
  83. # load bool for swings
  84. fswing = pd.read_csv(os.path.join(inpath,f,f'{f}_swing.csv'), index_col=0)
  85. fswing.index = pd.MultiIndex.from_product([[f],fswing.index])
  86. swing_all = pd.concat([swing_all,fswing])
  87. # load bool for turns
  88. fturns = pd.read_csv(os.path.join(inpath,f,f'{f}_turns.csv'), index_col=0)+lastmaxidx
  89. fturns.index = pd.MultiIndex.from_product([[f],fturns.index])
  90. turns_all = pd.concat([turns_all,fturns])
  91. # load bool for straight walks
  92. fstraight = pd.read_csv(os.path.join(inpath,f,f'{f}_straight.csv'), index_col=0)
  93. fstraight.index = pd.MultiIndex.from_product([[f],fstraight.index])
  94. st_offset = pd.DataFrame([False]*(len(f_clbp)-len(fstraight)), columns=fstraight.columns)
  95. st_offset.index = pd.MultiIndex.from_product([[f],st_offset.index])
  96. straight_all = pd.concat([straight_all, fstraight, st_offset])
  98. lastmaxidx += len(f_clbp)
  100. # %%
  101. lastmaxidx = 0
  102. rec_cycinfo = {}
  103. cyc_bouts = pd.Series([])
  104. cyc_rec = pd.Series([])
  105. # split recordings into bins based on swings
  106. for i,rec in enumerate(swing_all.index.get_level_values(0).unique()):
  107. # get swing onsets
  108. swing_dict = tg.get_bouts(swing_all.loc[rec][limbs])
  109. swing_on = {k:[s[0] for s in t] for k,t in swing_dict.items()}
  110. swing_dur = {k:[s[1] for s in t] for k,t in swing_dict.items()}
  111. on_df = pd.DataFrame().from_dict(swing_on, orient='index').T
  112. # use the first new swing onset of leg to determine bin end
  113. cyc = on_df.min(axis=1).astype(int)
  114. cyc = pd.concat([pd.Series([0]), cyc, pd.Series(len(swing_all.loc[rec]))])
  115. cyc_ = cyc + lastmaxidx
  116. # keep bin duration and bin index for each rec
  117. cyc_info = cyc_.diff().rename('dur').dropna().reset_index(drop=True).to_frame()
  118. cyc_info['idx'] = cyc_info.index + len(cyc_bouts) #use +len(cyc_bouts) to have one index for all recs
  119. rec_cycinfo[rec] = cyc_info
  120. # keep rec index for bins
  121. cyc_bouts = pd.concat([cyc_bouts,cyc_[:-1]])
  122. cyc_rec = pd.concat([cyc_rec,pd.Series(np.full_like(cyc_[:-1], i))])
  123. lastmaxidx += len(swing_all.loc[rec])
  125. cyc_bouts = cyc_bouts.reset_index(drop=True)
  126. len(cyc_bouts), cyc_bouts.diff().mean(), cyc_bouts.diff().std()
  128. # %%
  129. # duty factor
  130. stancedur_all = pd.DataFrame([])
  131. phasedur_all = pd.DataFrame([])
  132. phasevelo_all = pd.DataFrame([])
  133. for i,rec in enumerate(swing_all.index.get_level_values(0).unique()):
  134. stance_dict = tg.get_bouts(~swing_all.loc[rec][limbs])
  135. stance_on = {k:[s[0] for s in t] for k,t in stance_dict.items()}
  136. stance_dur = {k:[s[1] for s in t] for k,t in stance_dict.items()}
  137. on_df = pd.DataFrame().from_dict(stance_on, orient='index').T
  139. stancedur = pd.DataFrame().from_dict(stance_dur, orient='index').T
  140. stancedur.index = pd.MultiIndex.from_product([[rec],stancedur.index])
  141. stancedur_all = pd.concat([stancedur_all,stancedur])
  143. phasedur = pd.DataFrame([np.diff(on_df[bp], append=len(swing_all.loc[rec]), axis=0) for bp in limbs], index=limbs).T
  144. phasedur.index = pd.MultiIndex.from_product([[rec],phasedur.index])
  145. phasedur_all = pd.concat([phasedur_all,phasedur])
  147. bp_grouper = pd.DataFrame([on_df[bp].dropna().set_axis(on_df[bp].dropna()).reindex(np.arange(len(swing_all.loc[rec])), method='ffill') for bp in limbs], index=limbs).T
  148. phasevelo = pd.concat([CMSvelo_all.loc[rec].groupby(bp_grouper[bp].values).mean().reset_index(drop=True) for bp in limbs], axis=1)
  149. phasevelo.columns = limbs
  150. phasevelo.index = pd.MultiIndex.from_product([[rec],phasevelo.index])
  151. phasevelo_all = pd.concat([phasevelo_all,phasevelo])
  152. duty_factor = stancedur_all / phasedur_all
  154. # %%
  155. # set base outpath
  156. outpath = os.path.join(inpath, 'cluster', 'out_cluster_050226')
  157. if not os.path.exists(outpath):
  158. os.makedirs(outpath)
  160. # %%
  161. np.sum(duty_factor.count(axis=0).values.reshape(4,2),axis=1), duty_factor.count(axis=0)
  163. # %%
  164. duty_factor.mean(axis=0)
  166. # %%
  167. #pair_color = ['r', 'g', 'b', 'y']
  168. fig, axs = plt.subplots(1,5, figsize=(10,3), width_ratios=[1,1,1,1,.05])
  169. for i,p in enumerate(['leg_1', 'leg_2', 'leg_3', 'hindleg']):
  170. pair = [bp for bp in duty_factor.columns if p in bp]
  172. or_nan = (phasevelo_all[pair].isna().any(axis=1) | duty_factor[pair].isna().any(axis=1))
  173. phasevelo_pair = phasevelo_all[pair][~or_nan]#.values.flatten()
  174. phasevelo_pair = pd.concat([phasevelo_pair[c] for c in pair])
  175. duty_factor_pair = duty_factor[pair][~or_nan]#.values.flatten()
  176. duty_factor_pair = pd.concat([duty_factor_pair[c] for c in pair])
  178. # test individually
  179. rs, pvs = [],[]
  180. for rec in phasevelo_pair.index.get_level_values(0).unique():
  181. r,pv = stats.spearmanr(phasevelo_pair.loc[rec], duty_factor_pair.loc[rec])
  182. rs.append(r)
  183. pvs.append(pv)
  185. # combine after fisher, using fisher-z transformation for rho and fisher method for p values
  186. fisher_z = np.arctanh(rs)
  187. fisher_z_mean = np.mean(fisher_z)
  188. r_mean = np.tanh(fisher_z_mean)
  189. chi2, p_comb = stats.combine_pvalues(pvs)
  190. print(p, r_mean, p_comb, chi2)
  192. Z, xedges, yedges = np.histogram2d(phasevelo_pair, duty_factor_pair, bins=[50,30], range=[[0,150],[0,1]], density=True)
  193. _ = axs[i].pcolormesh(xedges, yedges, Z.T, cmap='magma', vmin=0, vmax=.21)
  194. axs[i].set_title(p)
  195. axs[i].text(1,.01,f'r: {np.round(r_mean, 2)}; p: {np.format_float_scientific(p_comb, precision=2)}', c='w', )
  197. plt.colorbar(_, cax=axs[-1])
  198. plt.savefig(os.path.join(outpath,f'duty_factor.pdf'), bbox_inches='tight')
  199. plt.show()
  201. # %%
  202. # example bins
  203. fig, axs = plt.subplots(1,10, figsize=(20,2), sharex=True)
  204. s = 100
  205. for i,j in enumerate(range(s,s+10)):
  206. axs[i].imshow(swing_all.iloc[cyc_bouts[j]:cyc_bouts[j+1]][limbs_ord].T, aspect='auto', interpolation='None')
  207. plt.savefig(os.path.join(outpath,'examplecycle.png'), bbox_inches='tight')
  209. # %% [markdown]
  210. # ### cross-similarity
  212. # %%
  213. # wrangle bins for each data type into lists
  214. X_train = [CLbp_all[limbs].iloc[cyc_bouts[j]:cyc_bouts[j+1]].values for j in range(len(cyc_bouts)-1)] #limbs / frontlegs
  215. swing_train = [swing_all[limbs].iloc[cyc_bouts[j]:cyc_bouts[j+1]].values for j in range(len(cyc_bouts)-1)] #limbs / frontlegs
  216. straight_train = [straight_all.iloc[cyc_bouts[j]:cyc_bouts[j+1]].values for j in range(len(cyc_bouts)-1)]
  217. angle_train = [Trkangle_all.iloc[cyc_bouts[j]:cyc_bouts[j+1]].values for j in range(len(cyc_bouts)-1)]
  218. cmsvelo_train = [CMSvelo_all.iloc[cyc_bouts[j]:cyc_bouts[j+1]].values for j in range(len(cyc_bouts)-1)]
  219. cyc_idx_turns = np.concatenate([np.where((cyc_bouts[:-1]<t).values&(cyc_bouts[1:]>t).values)[0] for t in turns_all.values.flatten()])
  221. cyc_trainidx = np.arange(len(cyc_bouts)-1)#[np.array(legs_in_swing) >= 8]
  223. dur_train = [len(x) for x in X_train]
  224. len(dur_train), np.mean(dur_train), np.std(dur_train)
  226. # %%
  227. # resample into same length
  228. ts_resampler = TimeSeriesResampler(sz=15)
  229. X_train = ts_resampler.fit_transform(X_train)
  230. X_train = TimeSeriesScalerMeanVariance().fit_transform(X_train)
  232. # transform bools for swing and straight accordingly
  233. swing_train_new = ts_resampler.transform(swing_train) #TODO check if necessary when below step
  234. swing_train_new = np.round(swing_train_new).astype(int)
  235. straight_train_new = ts_resampler.transform(straight_train) #TODO check if necessary when below step
  236. straight_train_new = np.round(straight_train_new).astype(int)
  239. # %%
  240. # get mean angle and velocity for each bin
  241. angle_train_new = np.vstack([np.nanmean(v) for v in angle_train]).flatten()
  242. cmsvelo_train_new = np.vstack([np.nanmean(v) for v in cmsvelo_train]).flatten()
  243. cmsvelo_train_new.shape
  245. # %%
  246. # compute cross similarity matrix for different sets of bodyparts, using dynamical time warping
  247. dtw_list = []
  248. dtw_range = []
  249. for r in [(0,8), (0,6)]:
  250. x_train = X_train[:,:,r[0]:r[1]]
  251. dtw_dist = metrics.cdist_dtw((x_train), (x_train))
  252. dtw_list.append(dtw_dist)
  253. dtw_range.append(r)
  255. # %%
  256. # plot cross-similarity in comparable way between sets
  257. triu_idx = np.triu_indices(dtw_list[0].shape[0])
  258. dtw_list_max = np.max([np.max(d) for d in dtw_list])
  259. dtw_list_diff_max = np.max([np.max(abs(dtw_list[0]-d)) for d in dtw_list])+1
  261. for d,r in zip(dtw_list,dtw_range):
  262. plt.imshow(d, vmin=0, vmax=dtw_list_max)
  263. plt.colorbar()
  264. plt.savefig(os.path.join(outpath,f'dtw_crosssimilaritydistance_{r[0]}-{r[1]}.pdf'), bbox_inches='tight')
  265. plt.show()
  267. plt.imshow(dtw_list[0]-d,vmin=0, vmax=dtw_list_diff_max)
  268. plt.colorbar()
  269. plt.savefig(os.path.join(outpath,f'dtw_crosssimilaritydistance_{r[0]}-{r[1]}_difffull.pdf'), bbox_inches='tight')
  270. plt.show()
  272. # %%
  273. # calculate the agreement between pairwise cross-similarity matrices using spearman r
  274. agree_mat = np.zeros((len(dtw_list),len(dtw_list)))
  275. for i in range(len(dtw_list)):
  276. for j in range(len(dtw_list)):
  277. dtw_i = dtw_list[i]
  278. dtw_j = dtw_list[j]
  279. rho, p = stats.spearmanr(dtw_i[triu_idx], dtw_j[triu_idx])
  280. agree_mat[i,j] = rho
  281. # plot agreement
  282. plt.imshow(agree_mat, vmin=0, vmax=1)
  283. plt.colorbar()
  284. plt.savefig(os.path.join(outpath,f'dtw_spearmanr.pdf'), bbox_inches='tight')
  285. np.savetxt(os.path.join(outpath,f'dtw_spearmanr.csv'), agree_mat, delimiter=',')
  287. # %% [markdown]
  288. # ### clustering
  290. # %%
  291. # define x_train for one set of bodyparts, redo for others
  292. dtw_dist = dtw_list[0] # [1]
  293. dtw_name = dtw_range[0] # [1]
  294. x_train = X_train[:,:,dtw_name[0]:dtw_name[1]]
  295. outpath_leg = os.path.join(outpath,f'legs{dtw_name[0]}-{dtw_name[1]}')
  296. if not os.path.exists(outpath_leg):
  297. os.makedirs(outpath_leg)
  299. # %%
  300. # embedd bins into lower dim space using umap and similarity matrix as pairwise distance
  302. rng = np.random.default_rng(seed=42) #set random seeds to evaluate randomness in UMAP
  304. # UMAP paramaters
  305. neigh = 30
  306. neg_rate = 30
  307. min_dist = 0
  308. disc_dist = 12
  309. rep_str = 1
  310. # define UMAP
  311. emb = UMAP(metric='precomputed', init='random', n_neighbors=neigh, min_dist=min_dist,
  312. negative_sample_rate=neg_rate,
  313. disconnection_distance=disc_dist,
  314. repulsion_strength=rep_str,
  315. )
  317. # perform embedding multiple times
  318. Xemb_l = []
  319. while len(Xemb_l) < 6:
  320. emb.random_state = rng.integers(1000)
  321. print(emb.random_state)
  322. Xemb = emb.fit_transform(dtw_dist)
  323. Xemb_l.append(Xemb)
  325. # %%
  326. # plot UMAP embeddings
  327. fig, axs = plt.subplots(2, len(Xemb_l)//2, figsize=(10,5))
  328. axs = axs.T.flatten()
  329. for i in range(len(Xemb_l)):
  330. sc = axs[i].scatter(*Xemb_l[i].T, alpha=1, s=1, c=cmsvelo_train_new)
  331. axs[i].axis('equal')
  332. plt.colorbar(sc)
  333. plt.savefig(os.path.join(outpath_leg,f'UMAP_neigh{neigh}_mindist{min_dist}_negrate{neg_rate}_discdist{disc_dist}_repstr{rep_str}_legs{dtw_name[0]}-{dtw_name[1]}.pdf'), bbox_inches='tight')
  334. plt.show()
  336. # %%
  337. # perform clustering on embeddings
  338. # cluster settings
  339. nclust = 4
  340. metr = 'euclidean'
  341. link = 'ward'
  343. # cluster fr all embeddings
  344. cluster_l = []
  345. silscore_l = []
  346. for i in range(len(Xemb_l)):
  347. cluster = AgglomerativeClustering(n_clusters=nclust, metric=metr, linkage=link).fit(Xemb_l[i])
  348. cluster_l.append(cluster)
  349. silscore_l.append(silhouette_score(Xemb_l[i], cluster.labels_))
  352. # plot clustering
  353. plt.rcParams["axes.prop_cycle"] = plt.cycler("color", plt.cm.viridis(np.linspace(0,1,len(np.unique(cluster.labels_)))))
  354. fig, axs = plt.subplots(4, len(Xemb_l)//2, figsize=(10,5), height_ratios=(1,.1,1,.1))
  355. axs = axs.T.flatten()
  356. for i,j in enumerate(range(0,len(Xemb_l)*2,2)):
  357. for c in np.unique(cluster_l[i].labels_):
  358. axs[j].scatter(*Xemb_l[i][cluster_l[i].labels_==c].T, alpha=1,s=2, label=c)
  359. axs[j].axis('equal')
  360. axs[j+1].text(0,0, f'silscore: {np.round(float(silscore_l[i]),4)}')
  361. axs[j+1].axis('off')
  362. plt.legend()
  363. plt.savefig(os.path.join(outpath_leg,f'AgglomClust_ncl{nclust}_metric{metr}_link{link}.pdf'), bbox_inches='tight')
  364. plt.show()
  365. plt.rcdefaults()
  368. # %%
  369. # calculate agreement matrix between clusters
  370. agree_mat = np.zeros((len(Xemb_l),len(Xemb_l)))
  371. for i in range(len(Xemb_l)):
  372. for j in range(len(Xemb_l)):
  373. # get labels for two clusters
  374. labels_i = cluster_l[i].labels_.copy()
  375. labels_j = cluster_l[j].labels_.copy()
  377. # calcualte overlap between clusters to map clusters with same identity
  378. maparr = np.zeros((len(np.unique(labels_i)),len(np.unique(labels_i))))
  379. for c in np.unique(labels_i):
  380. n = np.zeros(len(np.unique(labels_i)))
  381. c_, n_ = np.unique(labels_j[labels_i==c], return_counts=True) # how many of the rows contain the same labels in j as in i
  382. n[c_] = n_
  383. maparr[c] = n/np.sum(n)
  385. # maps clusters with highest overlap, create dict of corresponding clusters
  386. mapper = {k:np.nan for k in np.unique(labels_i)}
  387. for c in np.argsort(np.max(maparr, axis=1))[::-1]: # sorts cluster with highest overlap to map first
  388. m_cand = np.argsort(maparr[c])[::-1] # sorts index of clusters in j, according to overlap with clusters in i
  389. for m in m_cand:
  390. if m in mapper.values(): # checks if cluster already assigned
  391. continue
  392. mapper[c] = m # assigns cluster
  393. break
  395. if not all([o in np.unique(list(mapper.values())) for o in np.unique(labels_i)]):
  396. raise ValueError(f"Not all values are assigned properly, current mapping dictionary is {mapper}")
  397. # rename labels in i
  398. for c in np.unique(labels_i):
  399. labels_i[labels_i==c] = -mapper[c] #first to negative to not merge different clusters
  400. labels_i = abs(labels_i) #revert negative with abs
  402. # after clusters are mapped, calculate agreement
  403. agree_mat[i,j] = np.sum(labels_i == labels_j)/len(labels_i)
  405. # plot agreement
  406. plt.imshow(agree_mat, vmin=0, vmax=1)
  407. plt.colorbar()
  408. plt.savefig(os.path.join(outpath_leg,f'agreement.pdf'), bbox_inches='tight')
  409. np.savetxt(os.path.join(outpath_leg,f'agreement.csv'), agree_mat, delimiter=',')
  410. print(np.mean(agree_mat), np.std(agree_mat))
  412. # %%
  413. # decide for best embedding
  414. best_embedd = 0
  415. Xemb = Xemb_l[best_embedd]
  416. cluster = cluster_l[best_embedd]
  418. # %%
  419. np.savetxt(os.path.join(outpath_leg,f'Xemb.csv'), Xemb, delimiter=',')
  420. np.savetxt(os.path.join(outpath_leg,f'clusterlabels.csv'), cluster.labels_, delimiter=',')
  422. # %%
  423. # get barycenter for different clusters
  424. barycent_clust = {}
  425. for c in np.unique(cluster.labels_):
  426. x_c = X_train[cluster.labels_==c]
  427. barycent_clust[c] = {}
  428. for bp in range(X_train.shape[2]):
  429. barycent_clust[c][bp] = barycenters.softdtw_barycenter(x_c[:,:,bp])
  431. # %%
  432. # plot hildebrandt style heatmap of clusters, after each bin has been time warped to barycenter
  433. fig, axs = plt.subplots(len(np.unique(cluster.labels_)),2, figsize=(6,len(np.unique(cluster.labels_))*2.2), sharey=True)
  434. gridspec = axs.flatten()[0].get_subplotspec().get_gridspec()
  435. subfigs = [fig.add_subfigure(gs) for gs in gridspec]
  437. limbs_ord_num = [6,4,2,0,1,3,5,7]
  438. for c in np.unique(cluster.labels_):
  439. x_train_c = X_train[cluster.labels_==c]
  440. bary_c = np.hstack(list(barycent_clust[c].values()))
  441. x_shift = []
  442. d_shift = []
  443. print(len(x_train_c))
  444. # shift each trace in comparision to barycenter
  445. for x_ in x_train_c:
  446. path,d = metrics.dtw_path(bary_c, x_)
  447. path = np.array(path)
  448. u,ui = np.unique(path[:,0], return_index=True)
  449. path = path[ui]
  450. x_shift.append(x_[path[:,1]])
  451. d_shift.append(d)
  452. x_shift = np.stack(x_shift)
  453. x_shift_mean = np.mean(x_shift, axis=0)
  455. # plot mean
  456. subfigs[0].suptitle(f'mean')
  457. axs[c,0].set_title(f'Cluster {c}')
  458. im = axs[c,0].imshow(x_shift_mean.T, vmin=-1, vmax=1) #[limbs_ord_num]
  459. axs[c,0].set_yticks(range(x_shift_mean.shape[1]))
  460. axs[c,0].set_yticklabels(limbs) #limbs_ord
  461. plt.colorbar(im)
  463. # plot barycenter
  464. subfigs[1].suptitle(f'barycenter')
  465. axs[c,1].set_title(f'Cluster {c}')
  466. im = axs[c,1].imshow(bary_c.T, vmin=-1, vmax=1) #[limbs_ord_num]
  467. plt.colorbar(im)
  468. plt.savefig(os.path.join(outpath_leg,f'hildeheatmap_c{c}.pdf'), bbox_inches='tight')
  469. plt.show()
  471. # %%
  472. # plot the movement wave of each leg for the different clusters
  473. limb_color_list = list(pd.Series(bp_color_dict)[limbs])
  474. for c in np.unique(cluster.labels_):
  476. rndm = np.random.choice(np.arange(len(cluster.labels_))[cluster.labels_==c], min(500, len(cluster.labels_[cluster.labels_==c])))
  477. bary_c = np.hstack(list(barycent_clust[c].values()))
  478. x_rndm = []
  479. d_rndm = []
  480. for r in rndm:
  481. path,d = metrics.dtw_path(bary_c[:,dtw_name[0]:dtw_name[1]], x_train[r])
  482. path = np.array(path)
  483. u,ui = np.unique(path[:,0], return_index=True)
  484. path = path[ui]
  485. x_rndm.append(X_train[r,path[:,1]])
  486. d_rndm.append(d)
  488. plt.rcParams["axes.prop_cycle"] = plt.cycler("color", plt.cm.viridis(np.linspace(0,1,len(rndm))))
  489. fig = plt.figure(layout="constrained", figsize=(4,10))
  490. gs = GridSpec(5, 2, figure=fig)
  491. axs = []
  492. for i in range(4):
  493. axs.append(fig.add_subplot(gs[i, 0]))
  494. axs.append(fig.add_subplot(gs[i, 1]))
  495. axs.append(fig.add_subplot(gs[i+1, :2]))
  496. axs = np.array(axs)
  498. plt.suptitle(f'Cluster {c}')
  499. r_i = np.argsort(d_rndm)[::-1]
  500. for i in range(X_train.shape[2]):
  501. for j,r in enumerate(r_i):
  502. x_r = x_rndm[r][:,i]
  503. axs[i].plot(x_r, alpha=.1, zorder=j)
  504. axs[i].plot(barycent_clust[c][i], ls=':', c=limb_color_list[i], zorder=len(rndm)+1, lw=3)
  505. axs[-1].plot(barycent_clust[c][i], ls=':', c=limb_color_list[i], lw=5)
  507. plt.savefig(os.path.join(outpath_leg,f'gaitwave_c{c}.pdf'), bbox_inches='tight')
  508. plt.show()
  509. plt.rcdefaults()
  511. # %% [markdown]
  512. # ### Statistics
  514. # %%
  515. # create a new dataframe o extract statistics
  516. turns_train_new = np.zeros(len(cluster.labels_))
  517. turns_train_new[cyc_idx_turns] = 1
  518. label_tracktype = pd.DataFrame([cluster.labels_, cmsvelo_train_new, angle_train_new, np.round(np.mean(straight_train_new, axis=1)).astype(int).flatten(), turns_train_new], index=['label','velocity','angle','straight','turn']).T
  519. label_tracktype_rec = label_tracktype.groupby([cyc_rec.reset_index(drop=True),label_tracktype['label']]).mean() # group by biological replicant and label
  521. # %%
  522. # ANOVA
  523. velosamples = []
  524. anglesamples = []
  525. for name, gr in label_tracktype_rec.groupby('label'):
  526. if name < 0:
  527. continue
  528. velosamples.append(gr['velocity'].dropna())
  529. anglesamples.append(gr['angle'].dropna())
  530. stats.f_oneway(*velosamples), stats.f_oneway(*anglesamples),
  532. # %%
  533. # tukey_hsd
  534. velores = stats.tukey_hsd(*velosamples)
  535. velo_sig_comb=np.unique(np.array(np.where(velores.pvalue < .05)).T, axis=0)
  536. velo_sig_combuniq=np.unique((np.sort(velo_sig_comb, axis=1)), axis=0)
  537. velo_pv_combuniq=velores.pvalue[velo_sig_combuniq[:,0],velo_sig_combuniq[:,1]]
  539. angleres = stats.tukey_hsd(*anglesamples)
  540. angle_sig_comb=np.unique(np.array(np.where(angleres.pvalue < .05)).T, axis=0)
  541. angle_sig_combuniq=np.unique((np.sort(angle_sig_comb, axis=1)), axis=0)
  542. angle_pv_combuniq=angleres.pvalue[angle_sig_combuniq[:,0],angle_sig_combuniq[:,1]]
  543. velores.pvalue, angleres.pvalue
  545. # %%
  546. def sigbracket_loc(x0, x1, y, h):
  547. return (x0,y), (x0,y+h), (x1,y+h), (x1,y)
  549. # %%
  550. # Plot boxplot for cluster (animal mean) for velocity and angle
  551. fig, axs = plt.subplots(2, figsize=(7,7))
  552. for name, gr in label_tracktype_rec.groupby('label'):
  553. if name < 0:
  554. continue
  555. axs[0].boxplot(gr['velocity'].dropna(), positions=[name], showfliers=False)
  556. axs[0].scatter(np.full_like(gr['velocity'].dropna(),name-.3), gr['velocity'].dropna(), c='k')
  557. axs[1].boxplot(gr['angle'].dropna(), positions=[name], showfliers=False)
  558. axs[1].scatter(np.full_like(gr['angle'].dropna(),name-.3), gr['angle'].dropna(), c='k')
  559. for pv in range(len(velo_pv_combuniq)):
  560. p0,p1,p2,p3 = sigbracket_loc(velo_sig_combuniq[pv][0],velo_sig_combuniq[pv][1], 90+(10*velo_sig_combuniq[pv][0]),5)
  561. axs[0].plot(*np.array([p0,p1,p2,p3]).T)
  562. for pv in range(len(angle_sig_combuniq)):
  563. p0,p1,p2,p3 = sigbracket_loc(angle_sig_combuniq[pv][0],angle_sig_combuniq[pv][1], .15+(.1*angle_sig_combuniq[pv][0]),.03)
  564. axs[1].plot(*np.array([p0,p1,p2,p3]).T)
  565. axs[0].set_ylim(0,100)
  566. axs[1].set_ylim(0,0.4)
  567. plt.savefig(os.path.join(outpath_leg,f'clusterdist_velo-angle.pdf'), bbox_inches='tight')
  568. plt.show()
  570. # %%
  571. # plot usage and normalized usage of cluster during turn
  572. counts = pd.DataFrame([], index=np.unique(cluster.labels_))
  573. for name, gr in label_tracktype.groupby('turn'):
  574. count, bin_edge = np.histogram(gr['label'], bins=len(np.unique(cluster.labels_)), range=(0,max(np.unique(cluster.labels_))+1))
  575. counts.loc[bin_edge[:-1].astype(int),name] = count
  576. if name == 0:
  577. continue
  579. plt.hist(gr['label'], histtype='bar', bins=len(np.unique(cluster.labels_)), range=(0,max(np.unique(cluster.labels_))+1))
  580. plt.savefig(os.path.join(outpath_leg,f'clusteruse_turn.pdf'), bbox_inches='tight')
  581. plt.show()
  583. plt.bar(counts.loc[0:].index, counts.loc[0:,1]/counts.loc[0:,0], .9)
  584. plt.savefig(os.path.join(outpath_leg,f'clusteruse_turn_norm.pdf'), bbox_inches='tight')
  586. # %%
  587. # save cluster labels for each recording
  588. all_label_org = pd.DataFrame([])
  589. for f in os.listdir(inpath):
  590. if expID in f and '.' not in f:
  591. rec_label = np.full(len(rec_cycinfo[f]), np.nan)
  592. rec_idx = [t for t in rec_cycinfo[f]['idx'] if t <= len(cluster.labels_)-1]
  593. rec_label[:len(rec_idx)] = cluster.labels_[rec_idx]
  594. rec_label_org = np.concatenate([np.repeat(rec_label[i], rec_cycinfo[f]['dur'][i]) for i in range(len(rec_label))]).astype(int)
  595. print(f, np.unique(rec_label_org, return_counts=True))
  596. np.savetxt(os.path.join(outpath_leg,f'{f}_gait.csv'), rec_label_org, delimiter=',', fmt='%i')
  597. rec_label_org = pd.DataFrame(rec_label_org)
  598. rec_label_org.index = pd.MultiIndex.from_product([[f],rec_label_org.index])
  599. all_label_org = pd.concat([all_label_org, rec_label_org])
  601. # %%

tardigrade_cluster.ipynb at commit 6841dc7, no license · at the source

Overview

Authors: Euphrasie Ramahefarivo1,2, Leonard Böger1,2,3, Takkasila Saichol1, Behzad Shomali1, Luis Alvarez1, Monika Scholz1
ORCID iDs: Euphrasie Ramahefarivo, Leonard Böger, Takkasila Saichol, Behzad Shomali, Monika Scholz
  1. Max Planck Research Group Neural Information Flow, Max Planck Institute for Neurobiology of Behavior—Caesar, Bonn, Germany
  2. International Max Planck Research School for Brain and Behavior, Bonn, Germany
  3. Max Planck Research Group Genetics of Behavior, Max Planck Institute for Neurobiology of Behavior—Caesar, Bonn, Germany
Institutions: Max Planck Society (Germany); Max Planck Institute for Neurobiology of Behavior – caesar (Germany)
Journal: Nature communications, volume 17, issue 1, article 4412
Dates: received 2 May 2025; accepted 23 April 2026; published online 19 May 2026
Type: Research article · Language: English
License: CC BY
Identifiers: DOI 10.1038/s41467-026-72710-3 · PMID 42156370 · PMCID PMC13187449 · OpenAlex W4409044563
Open access: gold, a free copy (OpenAlex)
Preprint: doi:10.1101/2025.04.01.646316
Status: code verified
Categories: histology / microscopy (modality), optical imaging (calcium, voltage, 2-photon) (modality), drosophila (organism), C. elegans (organism)
Methods: Statistics, Smoothing, state filtering, decompositions, Machine learning, Evoked potentials, Connectivity, fMRI & imaging, Single-unit activity, calcium imaging
Keywords: Wide-field fluorescence microscopy, Ca2+ imaging, Behavioural methods
MeSH: Microscopy, Fluorescence*, Animals, Caenorhabditis elegans, Calcium, Drosophila melanogaster, Larva, Muscles, Neurons (* major topic)
Topic: Advanced Fluorescence Microscopy Techniques (Biophysics, Biochemistry, Genetics and Molecular Biology), according to OpenAlex
Funding: EC | Horizon 2020 Framework Programme (EU Framework Programme for Research and Innovation H2020) (101098722); NSF | Directorate for Mathematical &amp; Physical Sciences | Division of Physics (PHY) (PHY-1748958)
Citations: cited by 1 paper (Europe PMC); 73 references in the paper

Abstract

We present a modular epifluorescence tracking microscope which enables ratiometric imaging of muscles, neurons, and other structures in moving animals. The microscope is assembled entirely from commercial parts within 3 h, making the system broadly accessible. Leveraging the improved brightness and bleaching characteristics of recent genetically encoded indicators and fluorophores, the simple microscope is even suitable for calcium imaging of neurons in behaving animals, as we demonstrate in C. elegans. We also show how muscle dynamics in D. melanogaster larvae can be analyzed and how dual color fluorescence tracking elucidates inter-species interactions by visualizing both predatory nematodes and their prey. Finally, we showcase a configuration for brightfield imaging by tracking tardigrade gait as an example of utility for non-labeled species. The affordability of the hardware and ease of use of the accompanying software make this a suitable tool for education in addition to its use in research.

Reproduced under the paper's license (CC BY), from the paper cited above.

Repositories

Its files are read in the Code ↔ Paper reader above, with 6 matches between paragraphs and lines of code.

scholz-lab.github.io/glowtracker

License: none: the authors keep all their rights
State: the link answers, verified on 28 September 2026
Evidence: the link answers
Software Heritage: not checked
Found in: the text, “Hardware: Modular design enables multiple imagin”
Not found: README, license file, CITATION.cff, environment file, tests, continuous integration, documentation
Availability: 1 check, the latest on 28 September 2026: the link answers (HTTP 200)
  • 28 September 2026: the link answers (HTTP 200)
At the source: scholz-lab.github.io/GlowTracker/

scholz-lab/Macroscope-Paper

License: none: the authors keep all their rights
State: the link answers, verified on 28 September 2026
Evidence: files inventoried
Commit: 6841dc791e28c7d1bcfdcc2ab4940f0b6899e933, 23 March 2026
Languages: Jupyter (10), Python (5)
Size: 18 files, 15 scripts
Software Heritage: not archived
Found in: “Data availability”
Holds: README, 10 notebooks
Not found: license file, CITATION.cff, environment file, tests, continuous integration, documentation
Tools: NumPy (15 files), pandas (14 files), Matplotlib (13 files), SciPy (10 files), scikit-image (7 files), scikit-learn (4 files), Pillow (1 file), statsmodels (1 file), UMAP (1 file)
Availability: 1 check, the latest on 28 September 2026: the link answers
  • 28 September 2026: the link answers
16 files
At the source: github.com/scholz-lab/Macroscope-Paper

scholz-lab/GlowTracker

License: GPL-3.0
State: the link answers, verified on 28 September 2026
Evidence: files inventoried
Commit: 5585048feccaee94ba5a2fd4e7f2c9059f0126a8, 15 July 2026
Languages: Python (29)
Size: 72 files, 29 scripts
Software Heritage: not archived
Found in: “Code availability”
Holds: README, license file, environment (pyproject.toml), continuous integration
Not found: CITATION.cff, tests, documentation
Tools: NumPy (17 files), Matplotlib (10 files), OpenCV (8 files), scikit-image (5 files), pandas (3 files), SciPy (2 files), Pillow (1 file), tifffile (1 file)
Availability: 1 check, the latest on 28 September 2026: the link answers
  • 28 September 2026: the link answers
31 files
At the source: github.com/scholz-lab/GlowTracker

Code availability

The Glowtracker software is available on GitHub https://github.com/scholz-lab/GlowTracker and installable via the Python package manager pip as “glowtracker”. Further documentation and instructions are available at https://scholz-lab.github.io/GlowTracker. At the time of publication, the version of GlowTracker was glowtracker 0.11.0.

Reproduced under the paper's license (CC BY), from the paper cited above.

Tracing map

Proposed by the machine: these links were found in the paper and verified at the source, without human review. The map will receive a Zenodo DOI once one of the paper's authors has validated it with their ORCID.

What the map holds:

  • 3 repositories of the authors' code, each at its verified commit, with its license and how the link was found in the paper;
  • 44 scripts, each with its path and the digest of its content;
  • 6 matches between paragraphs of the paper and lines of the code (method lexical-v1);
  • neither the text of the paper nor the code itself.

Its JSON (tracing-map.json) is deposited on Zenodo with its DOI once the map is validated.

Data

Datasets cited

Data Availability Statement

The data underlying this manuscript are shown in the figures and supplementary tables. Extended behavioral and imaging data are available at https://osf.io/rpfny/. Code for generating the figures is deposited under https://github.com/scholz-lab/Macroscope-Paper.

The Glowtracker software is available on GitHub https://github.com/scholz-lab/GlowTracker and installable via the Python package manager pip as “glowtracker”. Further documentation and instructions are available at https://scholz-lab.github.io/GlowTracker. At the time of publication, the version of GlowTracker was glowtracker 0.11.0.

Reproduced under the paper's license (CC BY), from the paper cited above.

Versions

The history of this record: each version stored by the harvester or made by a correction of its authors or of the maintainers of its code, and what changed in its facts. The texts of the paper (its abstract, its availability statements) are not part of it; versions that changed only those are not listed.

Version 1, 28 September 2026: the first record

Recorded: type, language, journal, volume, issue, pages, dates, 6 authors, 3 keywords, 8 MeSH terms, 2 funders, 71 references.

Cite

This paper

Ramahefarivo, E., Böger, L., Saichol, T., Shomali, B., Alvarez, L., & Scholz, M. (2026). A modular multi-color fluorescence microscope for simultaneous tracking of cellular activity and behavior. Nature communications, 17(1), 4412. https://doi.org/10.1038/s41467-026-72710-3

BibTeX

@article{ramahefarivo2026modular,
author = {Ramahefarivo, Euphrasie and Böger, Leonard and Saichol, Takkasila and Shomali, Behzad and Alvarez, Luis and Scholz, Monika},
title = {{A modular multi-color fluorescence microscope for simultaneous tracking of cellular activity and behavior}},
journal = {Nature communications},
year = {2026},
month = may,
volume = {17},
number = {1},
pages = {4412},
publisher = {Nature Publishing Group},
issn = {2041-1723},
doi = {10.1038/s41467-026-72710-3},
url = {https://doi.org/10.1038/s41467-026-72710-3},
pmid = {42156370},
pmcid = {PMC13187449}
}

RIS

TY - JOUR
AU - Ramahefarivo, Euphrasie
AU - Böger, Leonard
AU - Saichol, Takkasila
AU - Shomali, Behzad
AU - Alvarez, Luis
AU - Scholz, Monika
TI - A modular multi-color fluorescence microscope for simultaneous tracking of cellular activity and behavior
T2 - Nature communications
J2 - Nat Commun
PY - 2026
DA - 2026/05/19
VL - 17
IS - 1
SP - 4412
SN - 2041-1723
PB - Nature Publishing Group
DO - 10.1038/s41467-026-72710-3
UR - https://doi.org/10.1038/s41467-026-72710-3
LA - en
ER -

CSL-JSON

{
"id": "10.1038/s41467-026-72710-3",
"type": "article-journal",
"title": "A modular multi-color fluorescence microscope for simultaneous tracking of cellular activity and behavior",
"container-title": "Nature communications",
"author": [
{
"family": "Ramahefarivo",
"given": "Euphrasie"
},
{
"family": "Böger",
"given": "Leonard"
},
{
"family": "Saichol",
"given": "Takkasila"
},
{
"family": "Shomali",
"given": "Behzad"
},
{
"family": "Alvarez",
"given": "Luis"
},
{
"family": "Scholz",
"given": "Monika"
}
],
"container-title-short": "Nat Commun",
"volume": "17",
"issue": "1",
"page": "4412",
"DOI": "10.1038/s41467-026-72710-3",
"PMID": "42156370",
"PMCID": "PMC13187449",
"ISSN": "2041-1723",
"publisher": "Nature Publishing Group",
"URL": "https://doi.org/10.1038/s41467-026-72710-3",
"language": "en",
"issued": {
"date-parts": [
[
2026,
5,
19
]
]
}
}

The tracing map gets a citation of its own once an author has validated it and it has a DOI.

Similar papers

The papers with a page that share the most with this one: the tools found in their code, their categories, datasets, cited references and authors, the rarest counting most.

[1] doi:10.1038/s41467-026-72709-w [code]
An epifluorescence microscope design for naturalistic behavior and cellular activity in freely moving Caenorhabditis elegans.
Journal: Nature communications
In common: tifffile, OpenCV, scikit-image, 5 other tools, C. elegans, histology / microscopy, 12 references
[2] doi:10.1038/s41593-026-02257-5 [code]
Neural sequences underlying directed turning in Caenorhabditis elegans.
Journal: Nature neuroscience
In common: tifffile, Pillow, pandas, 3 other tools, C. elegans, 7 references
[3] doi:10.1038/s41467-026-72437-1 [code]
High-speed whole-brain imaging in Drosophila.
Journal: Nature communications
In common: OpenCV, scikit-image, Pillow, 6 other tools, drosophila, 5 references
[4] doi:10.1371/journal.pcbi.1013441 [code]
Large vision model framework for automated C. elegans analysis: From static morphometry to dynamic neural activity.
Journal: PLoS computational biology
In common: tifffile, OpenCV, scikit-image, 5 other tools, C. elegans, optical imaging (calcium, voltage, 2-photon), 3 references
[5] doi:10.1111/ejn.70582 [code]
Multifiber Array-Based Photometry System for Multiregional Functional Mapping in the Mouse Brain.
Journal: The European journal of neuroscience
In common: tifffile, OpenCV, scikit-image, 5 other tools, optical imaging (calcium, voltage, 2-photon), 3 references
[6] doi:10.1016/j.isci.2026.116206 [code]
Gut distension evokes rapid neural dynamics in vagal and hindbrain populations of larval zebrafish.
Journal: iScience
In common: tifffile, OpenCV, scikit-image, 7 other tools, optical imaging (calcium, voltage, 2-photon), 1 reference
[7] doi: [code]
Real-time closed-loop feedback system for mouse mesoscale cortical signal and movement control
Journal: eLife
In common: tifffile, OpenCV, scikit-image, 7 other tools, optical imaging (calcium, voltage, 2-photon), 1 reference
[8] doi:10.1038/s41467-026-74823-1 [code]
Cerebellar activity is triggered by reach endpoint during learning of a complex locomotor task.
Journal: Nature communications
In common: tifffile, UMAP, OpenCV, 7 other tools, 1 reference
[9] doi:10.7554/elife.109717 [code]
Retrosplenial cortex enables context-dependent goal-directed sensorimotor transformation.
Journal: eLife
In common: tifffile, OpenCV, scikit-image, 7 other tools, 1 reference
[10] doi:10.1038/s41467-026-75352-7 [code]
Mechanosensory encoding of surface mechanics optimizes locomotion.
Journal: Nature communications
In common: tifffile, OpenCV, scikit-image, 4 other tools, C. elegans, 2 references

Contribute

The authors of this paper can claim it, correct its record and validate its tracing map, and the maintainers of its code (its owner, or a public member of its organization) correct what it says of their repository; anyone signed in can ask for its removal. Every request goes to OSCR's own machine, which answers it; your account page follows them.

Sign in with ORCID to claim this paper as one of its authors, correct its record or validate its tracing map: when the paper's metadata lists your ORCID iD, you are recognized at once. Maintainers of its code: sign in with GitHub, then claim the repository on your account page.

Request its removal

To ask OSCR to remove this record, the copies of its authors' scripts or its tracing map, use the removal request page: signed in, you say who you are, what to remove and why, then review and confirm the request. Published rules decide every request (how).

Discussion, reproductions, activity

Discussion: questions and error reports about this paper and its code, from signed-in readers and its authors. It opens with sign-in.

Reproductions: reports from readers who ran the authors' code: what they reproduced, with which environment, commit and data. It opens with sign-in.

Activity: what happens around this paper: new versions of its record, its map's validation, discussions and reproductions. It opens with sign-in.