OSCR

Spatial cartography of human thymus enables the geopositioning of lineage transcription factors in rare mimetic thymic epithelial cells.

Code ↔ Paper

15 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 15 matches
  1. [1] § Method › Integration of spatial clusters across samples ↔ 2.DefineNiche/integrate_niches_functions.ipynb, lines 108–256 · score 0.98 · pyCombat, batch corrected, pseudobulk profiles, COL1A2, MS4A1, CD8A
  2. [2] § Method › Signalling pathways ↔ 3.CellCommunication/run_cellularCommunication.ipynb, lines 33–181 · score 0.86 · truncatedMean, cell communication, incoming, outgoing, probability, CellChat
  3. [3] § Method › Integration of spatial clusters across samples ↔ 1.ProcessSamples/scRNA_integrate.Rmd, lines 739–861 · score 0.85 · COL1A1, CD8A, FOXP3, IGHM, MKI67, TPSAB1
  4. [4] § Method › Spatial cell type annotation ↔ 2.DefineNiche/integrate_niches_functions.ipynb, lines 536–652 · score 0.84 · expressed genes identified, combined known cell, gene score, integrated fetal, 0–1, scRNA
  5. [5] § Method › Signalling pathways ↔ 1.ProcessSamples/scRNA_integrate.Rmd, lines 1228–1337 · score 0.81 · CellChat, communication probability, cell communication, incoming, outgoing, signals
  6. [6] § Method › Differential GO function ↔ 2.DefineNiche/integrate_niches_functions.ipynb, lines 776–797 · score 0.77 · GO_Biological_Process_2023, drug2cell, gene ontology
  7. [7] § Results › Thymus architecture and geopositioning of cell types ↔ 2.DefineNiche/integrate_niches_functions.ipynb, lines 108–256 · score 0.72 · PCA transformation, batch corrected, COL1A2, VWF, gene expression, CD34
  8. [8] § Method › sn/scRNA-seq integration with public data ↔ 1.ProcessSamples/scRNA_integrate.Rmd, lines 105–187 · score 0.69 · FindIntegrationAnchors, FindAllMarkers, UMAP, Seurat, sc, seq
  9. [9] § Results › Thymus architecture and geopositioning of cell types ↔ 2.DefineNiche/integrate_niches_functions.ipynb, lines 536–652 · score 0.69 · mTEC, known cell, spatial spot, expressed genes, scRNA, DP
  10. [10] § Results › High-resolution spatial transcriptomic atlas of human thymus ↔ 2.DefineNiche/integrate_niches_functions.ipynb, lines 654–720 · score 0.67 · COL1A2, MS4A1, CD8A, CD83, CD14, CD34
  11. [11] § Method › sc/snRNA-seq data processing ↔ 1.ProcessSamples/scRNA_integrate.Rmd, lines 33–100 · score 0.57 · SCTransform, Seurat, filtered, UMAP, RNA, PCA
  12. [12] § Method › Identification of mimeTFs in spatial multiome ↔ 4.MimeTFs/IdentifyMimeTFs.ipynb, lines 237–289 · score 0.55 · CITE spatial, mimeTFs, septa, multiome, protein, Stereo
  13. [13] § Method › Spatial transcriptomics data processing ↔ 3.CellCommunication/run_cellularCommunication.ipynb, lines 33–181 · score 0.54 · weighted, edge, mouse, diameter, networks, neighbours
  14. [14] § Results › Cell type mimicking function by mimeTFs expressing TECs ↔ 2.DefineNiche/integrate_niches_functions.ipynb, lines 776–797 · score 0.53 · Gene Ontology, biological processes, GO
  15. [15] § Results › Abundance and geopositioning of mimeTFs expressing cells ↔ 2.DefineNiche/Integrate_niche_celltype_samples.ipynb, lines 443–465 · score 0.50 · KRT1 high region, HC region, HCs, LQ, niche, spots

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 · 857 lines · 33 KB · CC-BY-4.0 · 7 matches

  1. # %% [markdown]
  2. # # Functions
  3. # ### 1. integrate spatial samples
  4. # ### 2. define niches
  5. # ### 3. determine cell type -> cell-type-specificity score and RCTD
  7. # %%
  8. import scanpy as sc
  9. import networkx as nx
  10. import pandas as pd
  11. import matplotlib.pyplot as plt
  12. import numpy as np
  13. import h5py
  14. from anndata._io.specs import read_elem
  15. import seaborn as sns
  16. from numpy import inf
  17. from scipy.interpolate import griddata
  18. from skimage import measure
  19. from shapely.geometry import MultiPoint, Polygon
  20. import alphashape
  21. import matplotlib.colors as mcol
  22. import decoupler as dc
  24. from combat.pycombat import pycombat
  25. import matplotlib.pyplot as plt
  26. import matplotlib.colors as mcolors
  28. from sklearn.cluster import KMeans
  29. from sklearn import datasets
  30. from sklearn.decomposition import PCA
  31. from sklearn.discriminant_analysis import LinearDiscriminantAnalysis
  32. import h5py
  33. from anndata._io.specs import read_elem
  34. import drug2cell as d2c
  35. import blitzgsea as blitz
  37. import warnings
  38. warnings.simplefilter(action='ignore', category=pd.errors.PerformanceWarning)
  39. warnings.filterwarnings("ignore", category=DeprecationWarning)
  40. warnings.simplefilter(action='ignore', category=FutureWarning)
  41. warnings.simplefilter(action='ignore', category=UserWarning)
  43. import matplotlib
  44. matplotlib.rcParams['pdf.fonttype'] = 42
  45. matplotlib.rcParams['ps.fonttype'] = 42
  46. plt.rcParams['pdf.fonttype'] = 42
  47. plt.rcParams['ps.fonttype'] = 42
  48. plt.style.use('seaborn-white')
  52. # %%
  53. # Load marker genes
  54. degenes_annot2 = pd.read_csv("/data/Combined_Analysis/DEG_scRNA-seq.csv", sep=",")
  55. custom_marker = pd.read_csv("/data/Combined_Analysis/Custom_markers.csv", sep=",")
  57. # %%
  58. ## Define colours scheme for the domain and niches
  59. domain_color = { 'C': '#4292c6',
  60. 'M':'#fc9272',
  61. 'S': '#41ab5d',
  62. 'LQ' : '#f0f0f0'
  63. }
  64. niche_color = {'C1': '#9ecae1','C2': '#4292c6','C3': '#08519c', 'M1':'#fc9272', 'M2':'#ef3b2c','M3':'#a50f15',
  65. 'S1': '#a1d99b','S2': '#74c476', 'S3': '#41ab5d', 'S4': '#006d2c', 'LQ' : '#f0f0f0', 'HC' : '#252525'
  66. }
  68. HC_niche_color = {'C1': '#9ecae1','C2': '#4292c6','C3': '#08519c', 'M1':'#fc9272', 'M2':'#ef3b2c','M3':'#a50f15',
  69. 'S1': '#a1d99b','S2': '#74c476', 'S3': '#41ab5d', 'S4': '#006d2c', 'LQ' : '#f0f0f0', 'KRT1+' : '#6a51a3','KRT1++' : '#3f007d'
  70. }
  72. # %%
  73. # Get pseudo-bulk profile for each leiden cluster within leiden cluster
  74. def findPseudobulk(adata, adata_raw, sample_name, save_folder, group='leiden'):
  75. """
  76. Calculate the pseudo bulk profile of the leiden clusters.
  78. Parameters:
  79. adata: spatial data with leiden cluster meta data
  80. adata_raw: raw counts spatial data
  81. sample_name: sample name
  82. save_folder: folder to save the file
  84. Returns:
  85. Write the pseudobulk profile to csv file
  86. """
  87. adata_raw.obs[group] = adata.obs.leiden
  88. adata_raw.obs.leiden = adata_raw.obs.leiden.astype('str')
  89. adata_raw.X = np.round(adata_raw.X)
  90. adata_raw.layers['raw_counts'] = adata_raw.X
  91. adata_raw.obs['orig.ident'] = sample_name
  92. pdata = dc.get_pseudobulk(
  93. adata_raw,
  94. sample_col='orig.ident',
  95. groups_col=group,
  96. layer='raw_counts',
  97. mode='sum',
  98. min_cells=0,
  99. min_counts=0
  100. )
  101. dc.plot_psbulk_samples(pdata, groupby=[group,'orig.ident'], figsize=(8, 4))
  102. pdata.T.to_df().to_csv(f"{save_folder}" + sample_name + group + '_sumgeneexp.csv')
  103. adata.obs.to_csv(f"{save_folder}" + sample_name + group + "_meta.csv")
  108. # %%
  109. def integrate_pca(df_exp, datasets, target_names, nPC=5, ncluster=20):
  110. """
  111. Integrate the sample cluster pseudobulk profiles using PCA
  112. - pyCombat batch correction
  113. - PCA transformation
  114. - Kmeans clustering
  115. - Plot marker gene expression
  117. Parameters:
  118. df_exp: pandas object with normalised
  119. datasets: list of data to correct for data
  120. target_names: sample names
  121. nPC: number of PCs to model
  122. ncluster: number of resulting groups
  124. Returns:
  125. Write the pseudobulk profile to csv file
  126. df_segment_avgSonly: Average profile of the kmeans group
  127. df_clusters: Contain sample-wise group information
  128. """
  129. batch = []
  130. for j in range(len(datasets)):
  131. batch.extend([j for _ in range(len(datasets[j].columns)-1)])
  132. # run pyComBat
  133. df_corrected = pycombat(df_exp,batch)
  134. df_corrected.to_csv(f"{save_folder}Psuedo_batchCorrected.csv")
  135. plt.boxplot(df_corrected)
  136. plt.show()
  137. X = df_corrected.transpose()
  138. y = np.array(batch)
  140. plt.style.use('seaborn-white')
  141. pca = PCA()
  142. pca.fit(X)
  143. per_var = np.round(pca.explained_variance_ratio_*100, decimals = 1)
  145. plt.figure(figsize = (10,6))
  146. plt.plot(range(1, len(per_var)+1), per_var.cumsum(), marker = "o", linestyle = "--")
  147. plt.grid()
  148. plt.ylabel("Percentage Cumulative of Explained Variance")
  149. plt.xlabel("Number of Components")
  150. plt.title("Explained Variance by Component")
  151. plt.show()
  153. pca = PCA(n_components = nPC)
  154. pca.fit(X)
  155. scores_pca = pca.transform(X)
  157. WCSS = []
  158. for i in range(1,30):
  159. kmeans_pca = KMeans(n_clusters = i, init = "k-means++", random_state = 42)
  160. kmeans_pca.fit(scores_pca)
  161. WCSS.append(kmeans_pca.inertia_)
  163. plt.figure(figsize = (10,6))
  164. plt.plot(range(1,30), WCSS, marker = "o", linestyle = "--")
  165. plt.grid()
  166. plt.title("Cluster using PCA Scores")
  167. plt.ylabel("WCSS")
  168. plt.xlabel("N Clusters")
  169. plt.show()
  171. kmeans_pca = KMeans(n_clusters = ncluster, init = "k-means++", random_state = 42)
  172. kmeans_pca.fit(scores_pca)
  174. df = pd.DataFrame(df_corrected.columns)
  176. # Concatening the original df with the components informations present in scores_pca
  177. df_clust_pca_kmeans = pd.concat([df, pd.DataFrame(scores_pca)], axis = 1)
  178. #print(df_clust_pca_kmeans)
  179. # Renaming the column label from each component
  180. df_clust_pca_kmeans.columns = ["sample_leiden", "comp1", "comp2", "comp3", "comp4", "comp5"]
  181. # Seting the cluster label to each observation, using the atribute .labels_
  182. df_clust_pca_kmeans["segment_kmeans_pca"] = kmeans_pca.labels_
  184. # Mapping each cluster segmentation and renaming their labels
  185. df_clust_pca_kmeans["segment"] = df_clust_pca_kmeans["segment_kmeans_pca"]
  187. df_clust_pca_kmeans.to_csv(f"{save_folder}Sample_psuedo_PCA.csv")
  188. sns.pairplot(df_clust_pca_kmeans[1:], hue='segment', palette='tab20')
  190. df_clust_pca_kmeans[['Sample', 'leiden_clusters']] = df_clust_pca_kmeans['sample_leiden'].str.split('_', n=1, expand=True)
  191. sns.pairplot(df_clust_pca_kmeans[1:], hue='Sample', palette='tab20')
  193. df_corrected_long = df_corrected.stack().reset_index().set_axis('Genes Sample_leiden avgExp'.split(), axis=1)
  194. df_corrected_long = pd.merge(df_corrected_long, df_clust_pca_kmeans, left_on='Sample_leiden', right_on='sample_leiden')
  195. df_corrected_long[['Sample', 'leiden_clusters']] = df_corrected_long['sample_leiden'].str.split('_', n=1, expand=True)
  196. df_corrected_long['Sample_segment'] = df_corrected_long['Sample'] + df_corrected_long['segment'].astype('str')
  198. df_segment_avg = df_corrected_long.groupby(['Sample_leiden', 'Genes']).mean('avgExp')
  199. df_segment_avg = df_segment_avg.reset_index()
  200. df_segment_avgM = df_segment_avg.pivot(index='Genes', columns='Sample_leiden', values='avgExp')
  202. sns.clustermap(df_segment_avgM[df_segment_avgM.axes[0].isin(['PRSS16',
  203. 'AIRE',
  204. 'COL1A2',
  205. 'HBB',
  206. 'CD4', 'CD8A',
  207. 'MS4A1', 'CD14', 'CD34'
  208. ])], annot=False, cmap="Reds", figsize=(50,3))
  211. df_segment_avg= df_corrected_long.groupby(['Sample_segment', 'Genes']).mean('avgExp')
  212. df_segment_avg = df_segment_avg.reset_index()
  213. df_segment_avgS = df_segment_avg.pivot(index='Genes', columns='Sample_segment', values='avgExp')
  214. df_segment_avgS.to_csv('df_segment_avgS.csv')
  216. sns.clustermap(df_segment_avgS[df_segment_avgS.axes[0].isin(['PRSS16', 'PSMB11','CCL25', 'LY75', 'TBATA',
  217. 'AIRE', 'FEZF2',
  218. 'KRT15',
  219. 'KRT1',
  220. 'EPCAM',
  221. 'LAMP3',
  222. 'ITGAX',
  223. 'CD74', 'CD83', 'MS4A1', 'CD34', 'VWF', 'COL1A2','CD4', 'CD8A', 'CD14','CDH5',
  224. 'PECAM1',
  225. 'HBB',
  226. 'LUM',
  227. 'PRX1'
  228. ])], annot=False, cmap="Reds", figsize=(15,8)).savefig(f"{save_folder}Sample_Segment.pdf")
  230. df_segment_avgSonly = df_corrected_long.groupby(['segment', 'Genes']).mean('avgExp')
  231. df_segment_avgSonly = df_segment_avgSonly.reset_index()
  232. df_segment_avgSonly = df_segment_avgSonly.pivot(index='Genes', columns='segment', values='avgExp')
  234. df_segment_avgSonly.to_csv(f"{save_folder}df_segment_avgSonlyF.csv")
  235. df_segment_avgSonly.sum().to_csv(f"{save_folder}segments_sumGenes.csv")
  236. df_clusters = df_clust_pca_kmeans
  237. df_clusters[['Sample', 'leiden_clusters']] = df_clusters['sample_leiden'].str.split('_', n=1, expand=True)
  239. sns.clustermap(df_segment_avgSonly[df_segment_avgSonly.axes[0].isin([
  240. 'PRSS16', 'PSMB11','CCL25', 'LY75', 'TBATA',
  241. 'AIRE', 'FEZF2',
  242. 'KRT15',
  243. 'KRT1',
  244. 'EPCAM',
  245. 'LAMP3',
  246. 'ITGAX',
  247. 'CD74', 'CD83', 'MS4A1', 'CD34', 'VWF', 'COL1A2','CD4', 'CD8A', 'CD14','CDH5',
  248. 'PECAM1',
  249. 'HBB',
  250. 'LUM',
  251. 'PRX1', 'VWF'
  252. ])], annot=False, cmap="viridis", z_score=True, figsize=(8,8)).savefig(f"{save_folder}Segments_markergenes_long.pdf")
  254. sns.clustermap(df_segment_avgSonly[df_segment_avgSonly.axes[0].isin(['AIRE',
  255. 'PRSS16', 'HBB'])], z_score=True, annot=False, cmap="Reds", figsize=(10,2)).savefig(f"{save_folder}Segments_markergenes.pdf")
  256. return df_segment_avgSonly, df_clusters
  259. # %%
  260. def addAnnotation(tmp_adata, Medulla, Cortex, Septa, LowQuality):
  261. """
  262. Add annotation based on the cluster/group numbers
  264. Parameters:
  265. tmp_adata: spatial data
  266. Medulla: list of medulla cluster numbers
  267. Cortex: list of cortex cluster numbers
  268. Septa: list of septa cluster numbers
  269. LowQuality: list of LQ cluster numbers
  271. Returns:
  272. spatial data with additional meta data
  273. """
  274. tmp_adata.obs["Domain"] = "NA"
  275. x = list(tmp_adata[tmp_adata.obs['segment'].isin(Medulla)].obs_names)
  276. tmp_adata.obs.loc[ x, 'Domain'] = 'M'
  277. x = list(tmp_adata[tmp_adata.obs['segment'].isin(Cortex)].obs_names)
  278. tmp_adata.obs.loc[ x, 'Domain'] = 'C'
  279. x = list(tmp_adata[tmp_adata.obs['segment'].isin(Septa)].obs_names)
  280. tmp_adata.obs.loc[ x, 'Domain'] = 'S'
  281. x = list(tmp_adata[tmp_adata.obs['segment'].isin(LowQuality)].obs_names)
  282. tmp_adata.obs.loc[ x, 'Domain'] = 'LQ'
  283. return tmp_adata
  286. # %%
  287. def compute_contour_boundary(points, alpha):
  288. """
  289. Compute the contour boundary (outline) of a given set of points.
  290. Compute alpha shape (concave hull)
  292. Parameters:
  293. points: x y points
  294. alpha: Large alpha values → More convex, simpler shape. Small alpha values → More concave, tighter fit around the points.
  296. Returns:
  297. Boundary coordinates
  298. """
  299. shape = alphashape.alphashape(points, alpha)
  300. if shape.is_empty:
  301. raise ValueError("Alpha shape failed; adjust alpha parameter.")
  302. if isinstance(shape, Polygon):
  303. boundary_coords = np.array(shape.exterior.coords)
  304. else:
  305. # If multiple polygons, combine their boundaries
  306. boundary_coords = np.vstack([np.array(poly.exterior.coords) for poly in shape.geoms])
  307. return boundary_coords
  309. def compute_contours(data, alpha):
  310. """
  311. Find the contour or boundary of the spatial groups
  313. Parameters:
  314. data: spatial data
  315. alpha: Large alpha values → More convex, simpler shape. Small alpha values → More concave, tighter fit around the points.
  317. Returns:
  318. The location of boundary position between spatial groups
  319. """
  320. df_all = pd.DataFrame(data)
  321. thres_len = 100
  322. category_mapping = {'C': 1, 'M': 2, 'S': 3}
  324. ###Outline
  325. outline_coords = compute_contour_boundary(df_all[['X', 'Y']].to_numpy(), alpha= alpha)
  327. boundary_data = []
  328. for point in outline_coords:
  329. boundary_data.append({
  330. 'Boundary': 'Outline',
  331. 'X': point[0],
  332. 'Y': point[1]
  333. })
  335. ### M
  336. df = df_all[df_all.Z.isin(['M'])]
  337. outline_coords = compute_contour_boundary(df[['X', 'Y']].to_numpy(), alpha=alpha)
  338. for point in outline_coords:
  339. boundary_data.append({
  340. 'Boundary': 'M',
  341. 'X': point[0],
  342. 'Y': point[1]
  343. })
  345. ### C-M
  346. df = df_all[df_all.Z.isin(['C','M'])]
  347. outline_coords = compute_contour_boundary(df[['X', 'Y']].to_numpy(), alpha=alpha)
  349. for point in outline_coords:
  350. boundary_data.append({
  351. 'Boundary': 'C-M',
  352. 'X': point[0],
  353. 'Y': point[1]
  354. })
  356. ### C-S
  357. df = df_all[df_all.Z.isin(['C','M', 'S'])]
  358. outline_coords = compute_contour_boundary(df[['X', 'Y']].to_numpy(), alpha= alpha)
  360. for point in outline_coords:
  361. boundary_data.append({
  362. 'Boundary': 'C-M-S',
  363. 'X': point[0],
  364. 'Y': point[1]
  365. })
  366. boundary_df = pd.DataFrame(boundary_data)
  367. return boundary_df
  371. def shortest_distance(point, coordinates):
  372. """
  373. Calculate the shortest distance from a given point to a list of coordinates.
  375. Parameters:
  376. point: one data x-y point
  377. coordinates: compare against these collection of points
  379. Returns:
  380. shortest distance value
  381. """
  382. point = np.array(point)
  383. coordinates = np.array(coordinates)
  384. # Compute squared distances (faster than directly computing Euclidean distances)
  385. squared_distances = np.sum((coordinates - point) ** 2, axis=1)
  386. # Find the index of the minimum distance
  387. min_index = np.argmin(squared_distances)
  388. # Compute the actual shortest distance
  389. shortest_dist= np.sqrt(squared_distances[min_index])
  390. # Get the closest point
  391. closest_point = tuple(coordinates[min_index])
  392. return shortest_dist
  395. def assign_subcategory(x, y, region, boundary_df): #contours):
  396. """
  397. Calculate distances and assign subcategories
  399. Parameters:
  400. region: C, M or S
  401. boundary_df: boundary positions
  403. Returns:
  404. Sub region assignment value
  405. """
  406. #distance_to_outline = shortest_distance((x,y), np.array(boundary_df[boundary_df.Boundary == 'Outline' ][['X', 'Y']].to_records(index=False)).tolist())
  407. distance_to_m = shortest_distance((x,y), np.array(boundary_df[boundary_df.Boundary == 'M' ][['X', 'Y']].to_records(index=False)).tolist())
  408. distance_to_cm = shortest_distance((x,y), np.array(boundary_df[boundary_df.Boundary == 'C-M' ][['X', 'Y']].to_records(index=False)).tolist())
  409. distance_to_cms = shortest_distance((x,y), np.array(boundary_df[boundary_df.Boundary == 'C-M-S' ][['X', 'Y']].to_records(index=False)).tolist())
  410. if region == 'C':
  411. if distance_to_m <= distance_to_cm:
  412. if distance_to_m <= distance_to_cms:
  413. return 'C1'
  414. elif distance_to_m > distance_to_cms:
  415. return 'C3'
  416. elif distance_to_m > distance_to_cm:
  417. if distance_to_m <= distance_to_cms:
  418. return 'C2'
  419. elif distance_to_m > distance_to_cms:
  420. return 'C3'
  421. elif region == 'M':
  422. if distance_to_m <= distance_to_cm:
  423. if distance_to_m <= 100:
  424. return 'M2'
  425. elif distance_to_m > 100:
  426. return 'M1'
  427. elif distance_to_m > distance_to_cm:
  428. return 'M3'
  429. elif region == 'S':
  430. if distance_to_cm <= distance_to_cms:
  431. if distance_to_m <= distance_to_cm:
  432. return 'S1'
  433. elif distance_to_m > distance_to_cm:
  434. return 'S2'
  435. else:
  436. return 'S3'
  439. def find_niche(tmp_adata, sample, s_value, alpha):
  440. """
  441. Main function to Find Niche based on boundaries and closeness to the boundaries
  443. Parameters:
  444. tmp_adata: spatial data
  445. s_value: spot size for spatial plotting function
  446. alpha: Large alpha values → More convex, simpler shape. Small alpha values → More concave, tighter fit around the points.
  447. sample: sample name
  449. Returns:
  450. spatial data with meta on Niche subregions based closeness to the boundaries
  451. """
  452. data = {
  453. 'X': tmp_adata.obs.x,
  454. 'Y': tmp_adata.obs.y,
  455. 'Z': tmp_adata.obs.Domain
  456. }
  457. df = pd.DataFrame(data)
  459. boundary_df = compute_contours(data, alpha)
  460. boundary_color = mcol.ListedColormap(["grey", "red", "blue","green"])
  461. plt.figure(figsize=(4, 3))
  462. plt.scatter(boundary_df['X'], boundary_df['Y'], c=pd.factorize(boundary_df['Boundary'])[0], s=0.05, cmap=boundary_color )
  463. plt.title("Boundaries")
  464. plt.xlabel("X")
  465. plt.ylabel("Y")
  466. plt.savefig(sample + "_boundary.pdf")
  467. print(boundary_df.Boundary.value_counts())
  468. df['Subcategory'] = df.apply(
  469. lambda row: assign_subcategory(row['X'], row['Y'], row['Z'], boundary_df),
  470. axis=1
  471. )
  473. df.to_csv(sample + '_subcategories.csv', index=False)
  474. boundary_df.to_csv(sample + '_boundaries.csv', index=False)
  475. tmp_adata.obs['Niche'] = df.Subcategory
  476. #sc.pl.embedding(tmp_adata, basis="spatial", color=['segment', 'Domain', 'Niche'], s=s_value, save= sample + "_domain_niche.pdf")
  477. tmp_adata.obs.to_csv(sample + '_meta_niche.csv', index=True)
  478. print(df.Subcategory.value_counts())
  480. # with plt.rc_context(): # Use this to set figure params like size and dpi
  481. # plt.rcParams['figure.dpi'] = 100
  482. # plt.rcParams["figure.figsize"] = (3,2)
  483. # sc.pl.dotplot(tmp_adata, gene_list, 'Niche', use_raw=False, cmap='viridis', show=True, dendrogram=False)
  484. return tmp_adata
  486. # %%
  489. # %%
  490. ## TissueTag
  491. import tissue_tag as tt
  492. def tissue_tag(adata, annotation_column, Number1, Number2):
  493. """
  494. Call the Tissue Tag python function to split the spatial Domains into concentric bins
  496. Parameters:
  497. adata: spatial data
  498. annotation_column: Domain C, M or S
  499. Number1: number of concentric bins - C
  500. Number2: number of concentric bins - M or S
  502. Returns:
  503. spatial data with meta on Tissue Tag based concentric niches
  504. """
  505. adata.obs['x'] = adata.obsm['spatial'][:, 0]
  506. adata.obs['y'] = adata.obsm['spatial'][:, 1]
  507. tt.dist2cluster_fast(
  508. df=adata.obs,
  509. annotation=annotation_column,
  510. KNN=10
  511. ) # calculate minimum mean distance of each spot to clusters
  512. structure = ['S','C','M']
  513. w = [0.2,0.8]
  514. df_anno = tt.calculate_axis_3p(adata.obs, anno=annotation_column, structure=structure, output_col='cma_3p', w=w)
  515. adata.obs = df_anno
  517. #cortex group
  518. labels_cortex = [f'C{i}' for i in range(1, Number1+1)]
  519. labels_cortex
  520. adata.obs['Cell_group'] = adata.obs[annotation_column].astype(str)
  521. adata.obs.loc[adata.obs[annotation_column] == 'C', 'Cell_group'] = pd.qcut(adata.obs.loc[adata.obs[annotation_column] == 'C', 'cma_3p'], q=Number1, labels=labels_cortex).astype(str)
  522. #Medulla group
  523. labels_medulla = [f'M{i}' for i in range(1, Number2+1)]
  524. labels_medulla
  525. adata.obs.loc[adata.obs[annotation_column] == 'M', 'Cell_group'] = pd.qcut(adata.obs.loc[adata.obs[annotation_column] == 'M', 'cma_3p'], q=Number2, labels=labels_medulla).astype(str)
  526. #Septa group
  527. labels_septa = [f'S{i}' for i in range(1, Number2+1)]
  528. labels_septa
  529. adata.obs.loc[adata.obs[annotation_column] == 'S', 'Cell_group'] = pd.qcut(adata.obs.loc[adata.obs[annotation_column] == 'S', 'cma_3p'], q=Number2, labels=labels_septa).astype(str)
  531. return adata
  536. # %%
  537. ## Cellspecificity score and RCTD cell type
  538. cts = [
  539. 'cTEC(P)','cTEC(Q)', 'cTEC(Q) SIRT1lo','cTEC DP',
  540. 'DN(Q)','DP(P) RUNX3hi', 'DP(P) FABP5hi', 'DP(Q) AQP3hi', 'DP(Q) FABP5lo', 'CD8aa','ab entry',
  541. 'mTEC I', 'mTEC II', 'mTEC III', 'mimeTEC', 'CD4/CD8/Treg', 'B pro(P)', 'B pro(Q)', 'B IGHMhi','B TRAF1hi', 'B TRAF5hi',
  542. 'DC1','DC2','aDC3','aDC1.2', 'pDC', 'Mono', 'Macro(Q)','Macro(Phago)', 'Mast',
  543. 'ETP', 'VSMC','Ery', 'Fib(P)', 'Fib(Q)', 'Endo']
  544. cts_rctd = ['ab.entry',
  545. 'aDC1.2', 'aDC3', 'B.IGHMhi', 'B.pro.P.', 'B.pro.Q.', 'B.TRAF1hi',
  546. 'B.TRAF5hi', 'CD4_CD8_Treg', 'CD8aa', 'cTEC.DP', 'cTEC.P.', 'cTEC.Q.',
  547. 'cTEC.Q..SIRT1lo', 'DC1', 'DC2', 'DN.Q.', 'DP.P..FABP5hi',
  548. 'DP.P..RUNX3hi', 'DP.Q..AQP3hi', 'DP.Q..FABP5lo', 'Endo', 'Ery', 'ETP',
  549. 'Fib.P.', 'Fib.Q.', 'Macro.Phago.', 'Macro.Q.', 'Mast', 'mimeTEC',
  550. 'Mono', 'mTEC.I', 'mTEC.II', 'mTEC.III', 'pDC', 'VSMC']
  552. def cell_specificity_score(adata, cts, sample):
  553. """
  554. cell-specificity-score calculates a gene score by combining known cell-type marker genes with the top ten differentially expressed genes
  555. identified from integrated fetal and pediatric scRNA-seq dataset
  557. Parameters:
  558. adata: spatial data
  559. cts: list of cell type names
  560. sample: sample name
  562. Returns:
  563. writes the cell type specificity score for all the spatial spots and cell type combination to a csv file
  564. """
  565. for c in cts:
  566. gl = degenes_annot2.gene[degenes_annot2.cluster == c][0:10]
  567. gl = list(set(gl).union(set(custom_marker.gene[custom_marker.celltype == c])))
  568. if len(gl) >0 :
  569. sc.tl.score_genes(adata,
  570. gl,
  571. ctrl_size=50,
  572. gene_pool=None, n_bins=25, score_name=c,
  573. random_state=0, copy=False, use_raw=None)
  574. adata.obs.loc[adata.obs[c] < 0.75, c] = 0
  575. adata.obs.loc[adata.obs[c] < 0.75, c] = 0
  576. adata.obs.loc[adata.obs[c] >= 0.75, c] = 1
  577. adata.obs.loc[adata.obs[c] >= 0.75, c] = 1
  578. adata.obs.to_csv(f'{save_folder}' + sample + "_CS_valuesMeta.csv")
  579. plot_CT_heatmaps(adata, cts, sample + 'CS')
  582. def plot_CT_heatmaps(adata, cts, sample):
  583. """
  584. Normalizes the cell type enrichment scores in sample based on Domains, Niche by boundary and Niche by tissuetag and visualises as heatmap
  586. Parameters:
  587. adata: spatial data
  588. cts: list of cell type names
  589. sample: sample name
  591. Returns:
  592. writes the cell type specificity score for all the spatial spots and cell type combination to a csv file
  593. """
  594. adata_meta = adata.obs
  595. meta_niche = adata_meta.groupby('Domain')[cts].sum()
  596. dfMax = meta_niche.max(axis=0) ### NOrmalize to 0-1 for each cell type
  597. meta_niche = meta_niche.divide(dfMax, axis=1)
  598. meta_niche[meta_niche ==np.nan] =0
  599. #print(meta_niche)
  600. t = meta_niche.sum(axis = 0)
  601. meta_niche.to_csv(f'{save_folder}' + sample + "_domain.csv")
  602. sub = meta_niche[meta_niche.index.isin (['C','M','S' ])]
  603. sns.clustermap(sub, row_cluster= False, col_cluster=False,
  604. cmap="viridis", figsize=(10,4)).savefig(f'{save_folder}' + sample + "_heatmap.pdf")
  606. adata_meta = adata.obs
  607. meta_niche = adata_meta.groupby('Niche')[cts].sum()
  608. dfMax = meta_niche.max(axis=0) ### NOrmalize to 0-1 for each cell type
  609. meta_niche = meta_niche.divide(dfMax, axis=1)
  610. meta_niche[meta_niche ==np.nan] =0
  611. #print(meta_niche)
  612. t = meta_niche.sum(axis = 0)
  613. meta_niche.to_csv(f'{save_folder}' + sample + "_niche.csv")
  614. # sub = meta_niche[meta_niche.index.isin (['C1','C2','C3','C4','M1','M2','M3','S1','S2','S3','S4' ])]
  615. sns.clustermap(sub, row_cluster= False, col_cluster=False,
  616. cmap="viridis", figsize=(10,4)).savefig(f'{save_folder}' + sample + "_niche_heatmap.pdf")
  618. adata_meta = adata.obs
  619. meta_niche = adata_meta.groupby('Niche_bysegment')[cts].sum() #Niche by boundary
  620. dfMax = meta_niche.max(axis=0) ### NOrmalize to 0-1 for each cell type
  621. meta_niche = meta_niche.divide(dfMax, axis=1)
  622. meta_niche[meta_niche ==np.nan] =0
  623. t = meta_niche.sum(axis = 0)
  624. meta_niche.to_csv(f'{save_folder}' + sample + "_Niche_bysegment.csv")
  625. sns.clustermap(meta_niche, row_cluster= False, col_cluster=False,
  626. cmap="viridis", figsize=(10,4)).savefig(f'{save_folder}' + sample + "_segment_heatmap.pdf")
  628. adata_meta = adata.obs
  629. meta_niche = adata_meta.groupby('Cell_group')[cts].sum() #Niche by tissuetag
  630. dfMax = meta_niche.max(axis=0) ### NOrmalize to 0-1 for each cell type
  631. meta_niche = meta_niche.divide(dfMax, axis=1)
  632. meta_niche[meta_niche ==np.nan] =0
  633. t = meta_niche.sum(axis = 0)
  634. meta_niche.to_csv(f'{save_folder}' + sample + "_Cell_group.csv")
  635. sns.clustermap(meta_niche, row_cluster= False, col_cluster=False,
  636. cmap="viridis", figsize=(10,4)).savefig(f'{save_folder}' + sample + "_Cell_group_heatmap.pdf")
  638. def rctd_score(adata, cts, sample, rctd):
  639. """
  640. Summarizes the RCTD cell type values
  642. Parameters:
  643. adata: spatial data
  644. cts: list of cell type names
  645. sample: sample name
  646. rctd: RCTD meta data
  648. Returns:
  649. calls cell type enrichment plots function for RCTD
  650. """
  651. adata.obs[cts] = rctd[cts].to_numpy()
  652. plot_CT_heatmaps(adata, cts, sample + 'RCTD')
  654. # %%
  655. def plot_spatialGeneExpression(adata, save_folder, sample, svalue):
  656. """
  657. Plotting functions for sptaial gene expression, clusters and niches.
  659. Parameters:
  660. adata: spatial data
  661. save_folder: folder to save plots
  662. sample: sample name
  663. svalue: spot size for spatial plotting function
  665. Returns:
  666. saves spatial plots in pdf and png for publication
  667. """
  668. with plt.rc_context(): # Use this to set figure params like size and dpi
  669. plt.rcParams['figure.dpi'] = 300
  670. plt.rcParams["figure.figsize"] = (3,3)
  672. sc.pl.embedding(adata, basis="spatial",
  673. color=['leiden_clusters'],
  674. s=svalue,
  675. show=False,
  676. color_map='Reds'
  677. )
  678. plt.savefig(f"{save_folder}" + sample + "leiden.png", bbox_inches="tight")
  679. plt.savefig(f"{save_folder}" + sample + "leiden.pdf", bbox_inches="tight")
  682. sc.pl.embedding(adata, basis="spatial",
  683. color=['Domain'],
  684. s=svalue, palette=domain_color,
  685. show=False,
  686. color_map='Reds'
  687. )
  688. plt.savefig(f"{save_folder}" + sample + "domain.png", bbox_inches="tight")
  689. plt.savefig(f"{save_folder}" + sample + "domain.pdf", bbox_inches="tight")
  691. sc.pl.embedding(adata, basis="spatial",
  692. color=[ 'Niche', 'Niche_bysegment', 'Cell_group'],
  693. s=svalue, palette=niche_color,
  694. show=False,
  695. color_map='Reds'
  696. )
  697. plt.savefig(f"{save_folder}" + sample + "groups.png", bbox_inches="tight")
  698. plt.savefig(f"{save_folder}" + sample + "groups.pdf", bbox_inches="tight")
  700. gene_list = ['CD4', 'CD8A', 'PRSS16', 'AIRE', 'CD34', 'MS4A1', 'CD83', 'CD14', 'COL1A2']
  701. nrow=1
  702. ncol=len(gene_list)
  704. fig,axs=plt.subplots(nrow,ncol,figsize=(18,2), dpi=100)
  705. # Plot expression for every marker on the corresponding Axes object
  706. col_idx=0
  707. row_idx=0
  708. for marker in gene_list:
  709. ax=axs[col_idx]
  710. sc.pl.embedding(adata, basis="spatial",
  711. color=marker,
  712. s=svalue, ax=ax,frameon=False,
  713. show=False,
  714. color_map='Reds')
  715. col_idx+=1
  717. # Alignment within the Figure
  718. fig.tight_layout()
  719. plt.savefig(f"{save_folder}" + sample + "markers.png", bbox_inches="tight")
  720. plt.savefig(f"{save_folder}" + sample + "marker.pdf", bbox_inches="tight")
  722. # %%
  723. def find_HC_regions(adata, save_folder, sample):
  724. """
  725. Annotating the spots with high KRT1 gene expression
  727. Parameters:
  728. adata: spatial data
  729. save_folder: folder to save plots
  730. sample: sample name
  732. Returns:
  733. saves meta data with KRT1 region information
  734. """
  735. cp_list = ['Sample', 'leiden_clusters', 'Domain', 'Niche', 'Niche_bysegment', 'Cell_group']
  736. adata.obs[cp_list] = meta[cp_list].to_numpy()
  737. adata.obs[cp_list] = adata.obs[cp_list].astype('category')
  739. adata.obs['Niche_bysegment_HC'] = adata.obs.Niche_bysegment.to_numpy()
  740. adata.obs['Niche_bysegment_HC'] = adata.obs.Niche_bysegment.to_numpy()
  741. adata.obs['Niche_bysegment_HC_M'] = adata.obs.Niche_bysegment.to_numpy()
  742. sub_adata = adata[adata[: , 'KRT1'].X > 0, :]
  743. sub_adata = sub_adata[sub_adata.obs.Domain.isin(['LQ', 'S']) ]
  744. sub_adataM = adata[adata[: , 'KRT1'].X > 0, :]
  745. sub_adataM = sub_adataM[sub_adataM.obs.Domain.isin(['LQ', 'S', 'M']) ]
  747. x = list(adata[adata.obs.index.isin(sub_adata.obs.index)].obs_names)
  748. adata.obs.loc[ x, 'Niche_bysegment_HC'] = 'KRT1+'
  749. x = list(adata[adata.obs.index.isin(sub_adataM.obs.index)].obs_names)
  750. adata.obs.loc[ x, 'Niche_bysegment_HC_M'] = 'KRT1+'
  752. sub_adata = adata[adata[: , 'KRT1'].X > 10, :]
  753. sub_adata = sub_adata[sub_adata.obs.Domain.isin(['LQ', 'S']) ]
  754. sub_adataM = adata[adata[: , 'KRT1'].X > 10, :]
  755. sub_adataM = sub_adataM[sub_adataM.obs.Domain.isin(['LQ', 'S', 'M']) ]
  757. x = list(adata[adata.obs.index.isin(sub_adata.obs.index)].obs_names)
  758. adata.obs.loc[ x, 'Niche_bysegment_HC'] = 'KRT1++'
  759. x = list(adata[adata.obs.index.isin(sub_adataM.obs.index)].obs_names)
  760. adata.obs.loc[ x, 'Niche_bysegment_HC_M'] = 'KRT1++'
  763. with plt.rc_context(): # Use this to set figure params like size and dpi
  764. plt.rcParams['figure.dpi'] = 300
  765. plt.rcParams["figure.figsize"] = (3,3)
  766. sc.pl.embedding(adata, basis="spatial", palette= HC_niche_color,
  767. color=['Niche_bysegment', 'Niche_bysegment_HC', 'Niche_bysegment_HC_M'], s=4)
  769. plt.savefig(f"{save_folder}" + sample + "groupsHC.png", bbox_inches="tight")
  770. plt.savefig(f"{save_folder}" + sample + "groupsHC.pdf", bbox_inches="tight")
  772. adata.obs.to_csv(f"{save_folder}" + sample + "_niche_HC.csv")
  773. print(adata.obs.Niche_bysegment_HC.value_counts())
  774. print(adata.obs.Niche_bysegment_HC_M.value_counts())
  776. # %%
  777. targets = blitz.enrichr.get_library("GO_Biological_Process_2023")
  778. def find_GO(adata, region):
  779. """
  780. Determine the differential gene ontology (GO) biological process terms y Domains or Niches
  782. Parameters:
  783. adata: spatial data
  784. save_folder: folder to save plots
  785. sample: sample name
  787. Returns:
  788. returns top 10 GO for each group as dataframe
  789. """
  790. d2c.score(adata, targets=targets, use_raw=True)
  791. sc.tl.rank_genes_groups(adata.uns['drug2cell'], groupby=region, method='wilcoxon')
  792. # Convert results to a pandas DataFrame
  793. df = pd.DataFrame({
  794. group: adata.uns['drug2cell'].uns['rank_genes_groups']['names'][group][:10]
  795. for group in adata.uns['drug2cell'].uns['rank_genes_groups']['names'].dtype.names
  796. })
  797. return (df)
  799. # %%
  800. def ThreeGene_plot(adata, gene, name):
  801. """
  802. Plot spatial three genes in RGB colours and also co-localization by the overlap of 2 colours
  804. Parameters:
  805. adata: spatial data
  806. gene: combination of genes
  807. name: file name
  809. Returns:
  810. spatial plots in png and pdf
  811. """
  812. combinations = gene
  813. for i in range(len(combinations)):
  814. cc = combinations[i]
  815. sc.tl.score_genes(adata, [cc[0]], ctrl_size=50,
  816. gene_pool=None, n_bins=25, score_name="1",
  817. random_state=0, copy=False, use_raw=None)
  818. sc.tl.score_genes(adata, [cc[1]], ctrl_size=50,
  819. gene_pool=None, n_bins=25, score_name="2",
  820. random_state=0, copy=False, use_raw=None)
  821. sc.tl.score_genes(adata, [cc[2]], ctrl_size=50,
  822. gene_pool=None, n_bins=25, score_name="3",
  823. random_state=0, copy=False, use_raw=None)
  825. adata.obs['1'] = adata.obs['1'] / max(adata.obs['1'])
  826. adata.obs['2'] = adata.obs['2'] / max(adata.obs['2'])
  827. adata.obs['3'] = adata.obs['3'] / max(adata.obs['3'])
  829. adata.obs.loc[adata.obs['1'] < 0.1, '1'] = 0
  830. adata.obs.loc[adata.obs['2'] < 0.1, '2'] = 0
  831. adata.obs.loc[adata.obs['3'] < 0.1, '3'] = 0
  833. adata.obs.loc[adata.obs['1'] >= 0.1, '1'] = 1
  834. adata.obs.loc[adata.obs['2'] >= 0.1, '2'] = 1
  835. adata.obs.loc[adata.obs['3'] >= 0.1, '3'] = 1
  837. adata.obs['Colour3'] = 'k'
  838. adata.obs['Colour3'][(adata.obs['1'] == 0) & (adata.obs['2'] == 0) & (adata.obs['3'] == 0)] = np.nan
  839. adata.obs['Colour3'][(adata.obs['1'] == 1) & (adata.obs['2'] == 0) & (adata.obs['3'] == 0)] = 'r'
  840. adata.obs['Colour3'][(adata.obs['1'] == 0) & (adata.obs['2'] == 1) & (adata.obs['3'] == 0)] = 'g'
  841. adata.obs['Colour3'][(adata.obs['1'] == 0) & (adata.obs['2'] == 0) & (adata.obs['3'] == 1)] = 'b'
  842. adata.obs['Colour3'][(adata.obs['1'] == 1) & (adata.obs['2'] == 1) & (adata.obs['3'] == 0)] = 'y'
  843. adata.obs['Colour3'][(adata.obs['1'] == 1) & (adata.obs['2'] == 0) & (adata.obs['3'] == 1)] = 'm'
  844. adata.obs['Colour3'][(adata.obs['1'] == 0) & (adata.obs['2'] == 1) & (adata.obs['3'] == 1)] = 'c'
  845. adata.obs['Colour3'][(adata.obs['1'] == 1) & (adata.obs['2'] == 1) & (adata.obs['3'] == 1)] = 'k'
  847. with plt.rc_context():
  848. plt.rcParams['figure.dpi'] = 4000
  849. plt.rcParams["figure.figsize"] = (2, 2)
  850. sc.pl.embedding(adata, basis="spatial",
  851. color="Colour3", s=0.35, show=False,
  852. palette=mcolors.BASE_COLORS, na_color='whitesmoke')
  853. plt.savefig(f"{name}_{cc}_50_4000dpi.png", bbox_inches="tight", format="png")
  854. plt.savefig(f"{name}_{cc}_50_4000dpi.pdf", bbox_inches="tight", format="pdf" )
  856. # %%

integrate_niches_functions.ipynb, under CC-BY-4.0 · at the source

Overview

Authors: Uma S. Kamaraj1, Ying Chen1, Junjie Lei1, Pradeep Gautam1, Pongsatorn Horcharoensuk1,2, Czaryna K. M. Clemente3, Katja G. Weinacht4, Nicholas R. J. Gascoigne5,6,7,8, Jinmiao Chen9,10, Ching Kit Chen3,11, Qingfeng Chen1,9, Qi-Jing Li1,9, Lai Guan Ng12, Yuin-Han Loh1,2,8,13
ORCID iDs: Uma S. Kamaraj, Katja G. Weinacht, Nicholas R. J. Gascoigne, Jinmiao Chen, Ching Kit Chen, Qingfeng Chen, Qi-Jing Li, Lai Guan Ng, Yuin-Han Loh
13 affiliations
  1. Institute of Molecular and Cell Biology (IMCB), Agency for Science, Technology and Research (A*STAR),Singapore, Republic of Singapore
  2. Department of Physiology, NUS Yong Loo Lin School of Medicine,Singapore, Singapore
  3. Department of Pediatrics, Yong Loo Lin School of Medicine, National University of Singapore,Singapore, Singapore
  4. Department of Pediatrics, Division of Stem Cell Transplantation and Regenerative Medicine, Stanford University,Stanford, CA USA
  5. Immunology Translational Research Programme, Yong Loo Lin School of Medicine, National University of Singapore,Singapore, Singapore
  6. Department of Microbiology and Immunology, Yong Loo Lin School of Medicine, National University of Singapore,Singapore, Singapore
  7. Cancer Translational Research Programme, Yong Loo Lin School of Medicine, National University of Singapore,Singapore, Singapore
  8. NUS Graduate School’s Integrative Sciences and Engineering Programme, National University of Singapore,Singapore, Singapore
  9. Singapore Immunology Network (SIgN), Agency for Science, Technology and Research (A*STAR),Singapore, Singapore
  10. Immunology Translational Research Program, Department of Microbiology and Immunology, Yong Loo Lin School of Medicine, National University of Singapore,Singapore, Singapore
  11. Division of Cardiology, Department of Paediatrics, Khoo Teck Puat-National University Children’s Medical Institute, National University Health System,Singapore, Singapore
  12. Shanghai Immune Therapy Institute, School of Medicine, Shanghai Jiao Tong University, Renji Hospital,Shanghai, China
  13. Department of Biological Sciences, National University of Singapore,Singapore, Singapore
Institutions: Institute of Molecular and Cell Biology (Singapore); National University of Singapore (Singapore); Stanford University (United States); Duke-NUS Medical School (Singapore); Singapore Immunology Network (Singapore); National University Health System (Singapore); Renji Hospital (China)
Journal: Nature communications, volume 17, issue 1, article 3721
Dates: received 25 June 2024; accepted 12 January 2026; published online 10 March 2026
Type: Research article · Language: English
License: CC BY-NC-ND
Identifiers: DOI 10.1038/s41467-026-68596-w · PMID 41807357 · PMCID PMC13102919 · OpenAlex W7134846924
Open access: gold, a free copy (OpenAlex)
Status: code verified
Categories: genetics / omics (modality), human (organism), cellular / molecular (subfield)
Methods: Smoothing, state filtering, decompositions, Machine learning, Preprocessing, fMRI & imaging
Keywords: Thymus, Systems analysis, Imaging the immune system, Transcriptomics
MeSH: Epithelial Cells*, Thymus Gland*, Transcription Factors*, B-Lymphocytes, Cell Lineage, Child, Child, Preschool, Dendritic Cells, Fetus, Gene Expression Profiling, Gene Expression Regulation, Developmental, Humans, Proteomics, Spatial Transcriptomics, Thymocytes (* major topic)
Topic: Myasthenia Gravis and Thymoma (Neurology, Medicine), according to OpenAlex
Funding: National Research Foundation (NRFI2018-02); National Research Foundation Singapore (H23J2a0097, NRFI2018-02); Ministry of Education, India; Medical Research Council; National Medical Research Council (OFIRG21nov-0088); Biomedical Research Council
Citations: cited by 2 papers (Europe PMC); 65 references in the paper

Abstract

The abstract is not reproduced here: the paper's license (CC BY-NC-ND) does not allow it. Read it in the paper, at the publisher or on Europe PMC.

Repositories

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

UmaSangumathi/mimeTFs

License: none: the authors keep all their rights
State: the link answers, verified on 30 September 2026
Evidence: files inventoried
Commit: 6d1e67b8430c94e17b7bbe2803ed279995eef5d7, 6 October 2025
Languages: Jupyter (11), R (1)
Size: 14 files, 12 scripts
Software Heritage: not archived
Found in: “Code availability”
Holds: README, environment (requirements.txt), 12 notebooks
Not found: license file, CITATION.cff, tests, continuous integration, documentation
Tools: Matplotlib (2 files), NumPy (2 files), Scanpy (2 files), pandas (1 file), PyTorch (1 file), Squidpy (1 file)
Availability: 1 check, the latest on 30 September 2026: the link answers
  • 30 September 2026: the link answers
2 files
At the source: github.com/UmaSangumathi/mimeTFs.git

Zenodo 17851510

License: CC-BY-4.0
State: the link answers, verified on 30 September 2026
Evidence: files inventoried
Size: 1 file
Software Heritage: not checked
Found in: “Code availability”
Not found: README, license file, CITATION.cff, environment file, tests, continuous integration, documentation
Tools: data.table (2 files), ggplot2 (2 files), Matplotlib (2 files), NumPy (2 files), pandas (2 files), Scanpy (2 files), Seurat (2 files), tidyverse (2 files), anndata (1 file), clusterProfiler (1 file), cowplot (1 file), h5py (1 file), NetworkX (1 file), PyTorch (1 file), reshape2 (1 file), scikit-image (1 file), scikit-learn (1 file), SciPy (1 file), seaborn (1 file)
Availability: 1 check, the latest on 30 September 2026: the link answers (HTTP 200)
  • 30 September 2026: the link answers (HTTP 200)
7 files
At the source:

Code availability statement

The paper has a code availability statement. Its license (CC BY-NC-ND) does not allow reproducing it here; in short, from what the harvester recognized in it:

Read it in the paper: doi.org/10.1038/s41467-026-68596-w.

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:

  • 2 repositories of the authors' code, each at its verified commit, with its license and how the link was found in the paper;
  • 9 scripts, each with its path and the digest of its content;
  • 15 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 paper has a data availability statement. Its license (CC BY-NC-ND) does not allow reproducing it here; in short, from what the harvester recognized in it:

Read it in the paper: doi.org/10.1038/s41467-026-68596-w.

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, 30 September 2026: the first record

Recorded: type, language, journal, volume, issue, pages, dates, 14 authors, 4 keywords, 15 MeSH terms, 6 funders, 64 references.

Cite

This paper

Kamaraj, U. S., Chen, Y., Lei, J., Gautam, P., Horcharoensuk, P., Clemente, C. K. M., Weinacht, K. G., Gascoigne, N. R. J., Chen, J., Chen, C. K., Chen, Q., Li, Q.-J., Ng, L. G., & Loh, Y.-H. (2026). Spatial cartography of human thymus enables the geopositioning of lineage transcription factors in rare mimetic thymic epithelial cells. Nature communications, 17(1), 3721. https://doi.org/10.1038/s41467-026-68596-w

BibTeX

@article{kamaraj2026spatial,
author = {Kamaraj, Uma S. and Chen, Ying and Lei, Junjie and Gautam, Pradeep and Horcharoensuk, Pongsatorn and Clemente, Czaryna K. M. and Weinacht, Katja G. and Gascoigne, Nicholas R. J. and Chen, Jinmiao and Chen, Ching Kit and Chen, Qingfeng and Li, Qi-Jing and Ng, Lai Guan and Loh, Yuin-Han},
title = {{Spatial cartography of human thymus enables the geopositioning of lineage transcription factors in rare mimetic thymic epithelial cells}},
journal = {Nature communications},
year = {2026},
month = mar,
volume = {17},
number = {1},
pages = {3721},
publisher = {Nature Publishing Group},
issn = {2041-1723},
doi = {10.1038/s41467-026-68596-w},
url = {https://doi.org/10.1038/s41467-026-68596-w},
pmid = {41807357},
pmcid = {PMC13102919}
}

RIS

TY - JOUR
AU - Kamaraj, Uma S.
AU - Chen, Ying
AU - Lei, Junjie
AU - Gautam, Pradeep
AU - Horcharoensuk, Pongsatorn
AU - Clemente, Czaryna K. M.
AU - Weinacht, Katja G.
AU - Gascoigne, Nicholas R. J.
AU - Chen, Jinmiao
AU - Chen, Ching Kit
AU - Chen, Qingfeng
AU - Li, Qi-Jing
AU - Ng, Lai Guan
AU - Loh, Yuin-Han
TI - Spatial cartography of human thymus enables the geopositioning of lineage transcription factors in rare mimetic thymic epithelial cells
T2 - Nature communications
J2 - Nat Commun
PY - 2026
DA - 2026/03/10
VL - 17
IS - 1
SP - 3721
SN - 2041-1723
PB - Nature Publishing Group
DO - 10.1038/s41467-026-68596-w
UR - https://doi.org/10.1038/s41467-026-68596-w
LA - en
ER -

CSL-JSON

{
"id": "10.1038/s41467-026-68596-w",
"type": "article-journal",
"title": "Spatial cartography of human thymus enables the geopositioning of lineage transcription factors in rare mimetic thymic epithelial cells",
"container-title": "Nature communications",
"author": [
{
"family": "Kamaraj",
"given": "Uma S."
},
{
"family": "Chen",
"given": "Ying"
},
{
"family": "Lei",
"given": "Junjie"
},
{
"family": "Gautam",
"given": "Pradeep"
},
{
"family": "Horcharoensuk",
"given": "Pongsatorn"
},
{
"family": "Clemente",
"given": "Czaryna K. M."
},
{
"family": "Weinacht",
"given": "Katja G."
},
{
"family": "Gascoigne",
"given": "Nicholas R. J."
},
{
"family": "Chen",
"given": "Jinmiao"
},
{
"family": "Chen",
"given": "Ching Kit"
},
{
"family": "Chen",
"given": "Qingfeng"
},
{
"family": "Li",
"given": "Qi-Jing"
},
{
"family": "Ng",
"given": "Lai Guan"
},
{
"family": "Loh",
"given": "Yuin-Han"
}
],
"container-title-short": "Nat Commun",
"volume": "17",
"issue": "1",
"page": "3721",
"DOI": "10.1038/s41467-026-68596-w",
"PMID": "41807357",
"PMCID": "PMC13102919",
"ISSN": "2041-1723",
"publisher": "Nature Publishing Group",
"URL": "https://doi.org/10.1038/s41467-026-68596-w",
"language": "en",
"issued": {
"date-parts": [
[
2026,
3,
10
]
]
}
}

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.1016/j.xcrm.2026.102766 [code]
A longitudinal single-cell and spatial multiomic atlas of pediatric high-grade glioma.
Journal: Cell reports. Medicine
In common: anndata, clusterProfiler, Scanpy, 15 other tools, genetics / omics, cellular / molecular, 5 references
[2] doi:10.1038/s42003-026-10034-0 [code]
Region- and cell type-specific changes in gene expression in the cerebellum after classical fear conditioning.
Journal: Communications biology
In common: Squidpy, anndata, clusterProfiler, 13 other tools, genetics / omics, cellular / molecular, 3 references
[3] doi:10.1186/s13059-026-04177-w [code]
Genomic sequence evolution underlying human neocortical interareal diversification.
Journal: Genome biology
In common: Squidpy, anndata, Scanpy, 11 other tools, genetics / omics, cellular / molecular, 5 references
[4] doi:10.1038/s41467-026-71803-3 [code]
Charting the transition from in vitro gliogenesis to the in vivo maturation of human glial progenitor cells transplanted into the hypomyelinated mouse brain.
Journal: Nature communications
In common: Squidpy, anndata, Scanpy, 13 other tools, genetics / omics, cellular / molecular, 2 references
[5] doi:10.1093/nar/gkag706 [code]
scDifformer: diffusion-based post-training for virtual cell modeling across large-scale single-cell data.
Journal: Nucleic acids research
In common: Squidpy, anndata, Scanpy, 12 other tools, 3 references
[6] doi:10.21203/rs.3.rs-9676637/v1 [code]
A Comprehensive Benchmarking of Spatial Deconvolution and Domain Detection Methods across Diverse Tissues and Spatial Transcriptomic Technologies
Journal: Research Square (preprint)
In common: Squidpy, anndata, Scanpy, 14 other tools, genetics / omics, 1 reference
[7] doi:10.64898/2026.03.30.714220 [code]
An integrated single cell and spatial omics atlas of human prenatal development
Journal: bioRxiv (preprint)
In common: Squidpy, anndata, Scanpy, 9 other tools, 5 references
[8] doi:10.1038/s41592-026-03194-8 [code]
Beyond benchmarking: an expert-guided consensus approach to spatially aware clustering.
Journal: Nature methods
In common: Squidpy, anndata, Scanpy, 13 other tools, genetics / omics, 2 references
[9] doi:10.1002/imt2.70163 [code]
Spatial multi-omics unveils sphingolipid metabolic reprogramming within the retinal pathological niche.
Journal: iMeta
In common: anndata, clusterProfiler, Scanpy, 13 other tools, genetics / omics, cellular / molecular, 2 references
[10] doi:10.1186/s13073-026-01704-z [code]
Gene expression profiling enables refined parcellation of cortical layers in the heterogeneous human cerebral cortex.
Journal: Genome medicine
In common: anndata, Scanpy, Seurat, 13 other tools, genetics / omics, cellular / molecular, 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.