Skip to content

Commit 2a2e42f

Browse files
author
Corentin
committed
Notebook ATP Staining measurement explo
1 parent b393961 commit 2a2e42f

File tree

1 file changed

+282
-0
lines changed

1 file changed

+282
-0
lines changed

‎notebooks/atp_exploration.ipynb

Lines changed: 282 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,282 @@
1+
{
2+
"cells": [
3+
{
4+
"cell_type": "code",
5+
"execution_count": null,
6+
"metadata": {},
7+
"outputs": [],
8+
"source": [
9+
"import pandas as pd\n",
10+
"\n",
11+
"df = pd.read_csv('../data/muscle_atlas/muscle_atlas_2_7_filt_triple_full.csv')\n",
12+
"df.head()"
13+
]
14+
},
15+
{
16+
"cell_type": "code",
17+
"execution_count": null,
18+
"metadata": {},
19+
"outputs": [],
20+
"source": [
21+
"ATP_IMAGE_DF = df[df[\"Staining method\"] == \"ATP 9.4\"]\n",
22+
"ATP_IMAGE_NAME = ATP_IMAGE_DF[\"Number\"].tolist()\n",
23+
"IMAGE_INDEX = 189\n",
24+
"# show the first image in the list using the image name\n",
25+
"from IPython.display import Image\n",
26+
"Image(filename='../data/muscle_atlas/images/' + ATP_IMAGE_NAME[IMAGE_INDEX])"
27+
]
28+
},
29+
{
30+
"cell_type": "code",
31+
"execution_count": null,
32+
"metadata": {},
33+
"outputs": [],
34+
"source": [
35+
"from myoquant.src.common_func import is_gpu_availiable,load_cellpose,run_cellpose\n",
36+
"import matplotlib.pyplot as plt\n",
37+
"try:\n",
38+
" from imageio.v2 import imread\n",
39+
"except ImportError:\n",
40+
" from imageio import imread\n",
41+
"\n",
42+
"image_array = imread('../data/muscle_atlas/images/' + ATP_IMAGE_NAME[IMAGE_INDEX])\n",
43+
"model_cellpose = load_cellpose()\n",
44+
"mask_cellpose = run_cellpose(image_array, model_cellpose)\n",
45+
"plt.imshow(mask_cellpose)"
46+
]
47+
},
48+
{
49+
"cell_type": "code",
50+
"execution_count": null,
51+
"metadata": {},
52+
"outputs": [],
53+
"source": [
54+
"from skimage.measure import regionprops_table\n",
55+
"\n",
56+
"plt.imshow(mask_cellpose)\n",
57+
"props_cellpose = regionprops_table(\n",
58+
" mask_cellpose,\n",
59+
" properties=[\n",
60+
" \"label\",\n",
61+
" \"area\",\n",
62+
" \"centroid\",\n",
63+
" \"eccentricity\",\n",
64+
" \"bbox\",\n",
65+
" \"image\",\n",
66+
" \"perimeter\",\n",
67+
" \"feret_diameter_max\",\n",
68+
" ],\n",
69+
")\n",
70+
"df_cellpose = pd.DataFrame(props_cellpose)\n",
71+
"df_cellpose.head()"
72+
]
73+
},
74+
{
75+
"cell_type": "code",
76+
"execution_count": null,
77+
"metadata": {},
78+
"outputs": [],
79+
"source": [
80+
"import numpy as np\n",
81+
"all_cell_median_intensity = []\n",
82+
"for index in range(len(df_cellpose)):\n",
83+
" single_cell_img = image_array[\n",
84+
" df_cellpose.iloc[index, 5] : df_cellpose.iloc[index, 7],\n",
85+
" df_cellpose.iloc[index, 6] : df_cellpose.iloc[index, 8],\n",
86+
" ].copy()\n",
87+
"\n",
88+
" single_cell_mask = df_cellpose.iloc[index, 9].copy()\n",
89+
" single_cell_img[~single_cell_mask] = 0\n",
90+
" # Calculate median pixel intensity of the cell but ignore 0 values\n",
91+
" single_cell_median_intensity = np.median(single_cell_img[single_cell_img > 0])\n",
92+
" all_cell_median_intensity.append(single_cell_median_intensity)\n"
93+
]
94+
},
95+
{
96+
"cell_type": "code",
97+
"execution_count": null,
98+
"metadata": {},
99+
"outputs": [],
100+
"source": [
101+
"from scipy.stats import gaussian_kde\n",
102+
"\n",
103+
"# Build a \"density\" function based on the dataset\n",
104+
"# When you give a value from the X axis to this function, it returns the according value on the Y axis\n",
105+
"density = gaussian_kde(all_cell_median_intensity)\n",
106+
"density.covariance_factor = lambda : .25\n",
107+
"density._compute_covariance()\n",
108+
"\n",
109+
"# Create a vector of 256 values going from 0 to 256:\n",
110+
"xs = np.linspace(0, 255, 256)\n",
111+
"\n",
112+
"# Set the figure size\n",
113+
"plt.figure(figsize=(8, 4))\n",
114+
"\n",
115+
"# Make the chart\n",
116+
"# We're actually building a line chart where x values are set all along the axis and y value are\n",
117+
"# the corresponding values from the density function\n",
118+
"plt.plot(xs,density(xs))\n",
119+
"plt.show()"
120+
]
121+
},
122+
{
123+
"cell_type": "code",
124+
"execution_count": null,
125+
"metadata": {},
126+
"outputs": [],
127+
"source": [
128+
"# from sklearn.cluster import KMeans\n",
129+
"# import numpy as np\n",
130+
"\n",
131+
"# all_cell_median_intensity = np.array(all_cell_median_intensity)\n",
132+
"\n",
133+
"# # fit the k-means model to the data\n",
134+
"# kmeans = KMeans(n_clusters=2).fit(all_cell_median_intensity.reshape(-1, 1))\n",
135+
"\n",
136+
"# # get the threshold point between the two clusters\n",
137+
"# threshold = kmeans.cluster_centers_[0] if kmeans.cluster_centers_[0] < kmeans.cluster_centers_[1] else kmeans.cluster_centers_[1]\n"
138+
]
139+
},
140+
{
141+
"cell_type": "code",
142+
"execution_count": null,
143+
"metadata": {},
144+
"outputs": [],
145+
"source": [
146+
"print(len(density(xs)))\n",
147+
"print(len(all_cell_median_intensity))"
148+
]
149+
},
150+
{
151+
"cell_type": "code",
152+
"execution_count": null,
153+
"metadata": {},
154+
"outputs": [],
155+
"source": [
156+
"from sklearn.mixture import GaussianMixture\n",
157+
"\n",
158+
"# Fit the GMM\n",
159+
"\n",
160+
"gmm = GaussianMixture(n_components=2).fit(np.array(all_cell_median_intensity).reshape(-1, 1))\n",
161+
"\n",
162+
"# Find the x values of the two peaks\n",
163+
"peaks_x = gmm.means_.flatten()\n",
164+
"\n",
165+
"print('The x values of the two peaks are:', peaks_x)"
166+
]
167+
},
168+
{
169+
"cell_type": "code",
170+
"execution_count": null,
171+
"metadata": {},
172+
"outputs": [],
173+
"source": [
174+
"from scipy.stats import norm\n",
175+
"\n",
176+
"sorted_peaks = np.sort(peaks_x)\n",
177+
"# Find the minimum point between the two peaks\n",
178+
"min_index = np.argmin(density(xs)[(xs > sorted_peaks[0]) & (xs < sorted_peaks[1])])\n",
179+
"threshold = sorted_peaks[0]+xs[min_index]\n",
180+
"print(threshold)\n",
181+
"# Plot the data\n",
182+
"plt.plot(xs, density(xs), label='Density')\n",
183+
"plt.axvline(threshold, color='r', label='Threshold')\n",
184+
"plt.legend()\n",
185+
"plt.show()"
186+
]
187+
},
188+
{
189+
"cell_type": "code",
190+
"execution_count": null,
191+
"metadata": {},
192+
"outputs": [],
193+
"source": [
194+
"df_cellpose[\"cell_intensity\"] = all_cell_median_intensity\n",
195+
"df_cellpose[\"muscle_cell_type\"] = df_cellpose[\"cell_intensity\"].apply(\n",
196+
" lambda x: 1 if x < threshold else 2\n",
197+
")\n",
198+
"df_cellpose.head()"
199+
]
200+
},
201+
{
202+
"cell_type": "code",
203+
"execution_count": null,
204+
"metadata": {},
205+
"outputs": [],
206+
"source": [
207+
"df_cellpose[\"muscle_cell_type\"].value_counts(normalize=True)"
208+
]
209+
},
210+
{
211+
"cell_type": "code",
212+
"execution_count": null,
213+
"metadata": {},
214+
"outputs": [],
215+
"source": [
216+
"label_map = np.zeros(\n",
217+
" (image_array.shape[0], image_array.shape[1]), dtype=np.uint16\n",
218+
")\n",
219+
"# for index in track(range(len(df_cellpose)), description=\"Painting cells\"):\n",
220+
"for index in range(len(df_cellpose)):\n",
221+
" single_cell_mask = df_cellpose.iloc[index, 9].copy()\n",
222+
" if df_cellpose[\"muscle_cell_type\"][index] == 2:\n",
223+
" label_map[\n",
224+
" df_cellpose.iloc[index, 5] : df_cellpose.iloc[index, 7],\n",
225+
" df_cellpose.iloc[index, 6] : df_cellpose.iloc[index, 8],\n",
226+
" ][single_cell_mask] = 1\n",
227+
" elif df_cellpose[\"muscle_cell_type\"][index] == 1:\n",
228+
" label_map[\n",
229+
" df_cellpose.iloc[index, 5] : df_cellpose.iloc[index, 7],\n",
230+
" df_cellpose.iloc[index, 6] : df_cellpose.iloc[index, 8],\n",
231+
" ][single_cell_mask] = 2"
232+
]
233+
},
234+
{
235+
"cell_type": "code",
236+
"execution_count": null,
237+
"metadata": {},
238+
"outputs": [],
239+
"source": [
240+
"%config InlineBackend.figure_format = 'retina'\n",
241+
"from myoquant.src.common_func import label2rgb, blend_image_with_label\n",
242+
"labelRGB_map = label2rgb(image_array, label_map)\n",
243+
"overlay_img = blend_image_with_label(image_array, labelRGB_map)\n",
244+
"\n",
245+
"plt.figure(figsize=(10,20))\n",
246+
"\n",
247+
"f, axarr = plt.subplots(1,2) \n",
248+
"axarr[0].imshow(image_array)\n",
249+
"axarr[1].imshow(overlay_img)\n",
250+
"plt.tight_layout()\n",
251+
"plt.show()"
252+
]
253+
}
254+
],
255+
"metadata": {
256+
"kernelspec": {
257+
"display_name": ".venv",
258+
"language": "python",
259+
"name": "python3"
260+
},
261+
"language_info": {
262+
"codemirror_mode": {
263+
"name": "ipython",
264+
"version": 3
265+
},
266+
"file_extension": ".py",
267+
"mimetype": "text/x-python",
268+
"name": "python",
269+
"nbconvert_exporter": "python",
270+
"pygments_lexer": "ipython3",
271+
"version": "3.8.16"
272+
},
273+
"orig_nbformat": 4,
274+
"vscode": {
275+
"interpreter": {
276+
"hash": "da16f84656d11a3c096dd3524a83da95908ee8e4fba887e4173f286cf5f829c6"
277+
}
278+
}
279+
},
280+
"nbformat": 4,
281+
"nbformat_minor": 2
282+
}

0 commit comments

Comments
 (0)