Skip to content
Snippets Groups Projects
Commit 954be0c2 authored by Guillaume Cerutti's avatar Guillaume Cerutti
Browse files

update visualization in notebook

parent d474e241
Branches
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
notebook/test_nuclei_quantification.ipynb
+ 10
− 11
View file @ 954be0c2
%% Cell type:markdown id:6704c6e8 tags: %% Cell type:markdown id:6704c6e8 tags:
# Nuclei Fie/Cas9 Quantification # Nuclei Fie/Cas9 Quantification
### Necessary imports ### Necessary imports
%% Cell type:code id:7002eb02 tags: %% Cell type:code id:7002eb02 tags:
``` python ``` python
import numpy as np import numpy as np
import scipy.ndimage as nd import scipy.ndimage as nd
import pandas as pd import pandas as pd
import matplotlib.pyplot as plt import matplotlib.pyplot as plt
import pyvista as pv import pyvista as pv
pv.set_jupyter_backend('ipyvtklink')
from timagetk import SpatialImage from timagetk import SpatialImage
from timagetk.io import imread from timagetk.io import imread
from timagetk.algorithms.peak_detection import detect_nuclei from timagetk.algorithms.peak_detection import detect_nuclei
from timagetk.algorithms.signal_quantification import quantify_nuclei_signal_intensity from timagetk.algorithms.signal_quantification import quantify_nuclei_signal_intensity
from visu_core.vtk.utils.image_tools import vtk_image_data_from_image from visu_core.vtk.utils.image_tools import vtk_image_data_from_image
``` ```
%% Cell type:markdown id:c196d5ba tags: %% Cell type:markdown id:c196d5ba tags:
### Data specification ### Data specification
%% Cell type:code id:ab076632 tags: %% Cell type:code id:ab076632 tags:
``` python ``` python
dirname = "../data/" dirname = "../data/"
filename = "230315_GK_86_9_T1_5" filename = "230315_GK_86_9_T1_5"
``` ```
%% Cell type:markdown id:a1023b95 tags: %% Cell type:markdown id:a1023b95 tags:
### Parameters ### Parameters
%% Cell type:code id:86777931 tags: %% Cell type:code id:86777931 tags:
``` python ``` python
detection_threshold = 100, detection_threshold = 100,
detection_radius_range = (0.8, 3) detection_radius_range = (0.8, 3)
quantification_radius = 2.5 quantification_radius = 2.5
``` ```
%% Cell type:markdown id:165d542f tags: %% Cell type:markdown id:165d542f tags:
### Image loading ### Image loading
%% Cell type:code id:8bfe671b tags: %% Cell type:code id:8bfe671b tags:
``` python ``` python
image_filename = f"{dirname}/{filename}.tif" image_filename = f"{dirname}/{filename}.tif"
img = imread(image_filename, channel_names=['FieGFP', 'Calco', 'Cas9RFP']) img = imread(image_filename, channel_names=['FieGFP', 'Calco', 'Cas9RFP'])
``` ```
%% Cell type:code id:2234e546 tags: %% Cell type:code id:2234e546 tags:
``` python ``` python
image_datas = {} image_datas = {}
image_data = pv.ImageData(dimensions=img.shape[::-1], spacing=img.voxelsize[::-1])
for channel_name in img.channel_names: for channel_name in img.channel_names:
image_data = vtk_image_data_from_image(img[channel_name].get_array(), img.voxelsize) image_data.point_data[channel_name] = img[channel_name].get_array().ravel()
image_data = pv.UniformGrid(image_data)
image_data.point_data.GetArray(0).SetName(channel_name)
image_datas[channel_name] = image_data
plotter = pv.Plotter(notebook=True) plotter = pv.Plotter(notebook=False)
plotter.set_background('w') plotter.set_background('w')
plotter.add_volume(image_datas['Cas9RFP'], cmap='Greys', clim=(0, 2000)) plotter.add_volume(image_data, scalars='Cas9RFP', cmap='Greys', clim=(0, 2000))
plotter.show() plotter.show()
``` ```
%% Cell type:code id:6195742e tags: %% Cell type:code id:6195742e tags:
``` python ``` python
figure = plt.figure(figsize=(20, 10)) figure = plt.figure(figsize=(20, 10))
extent = ( extent = (
-img.voxelsize[2]/2, -img.voxelsize[2]/2,
img.extent[2]+img.voxelsize[2]/2, img.extent[2]+img.voxelsize[2]/2,
img.extent[1]+img.voxelsize[1]/2, img.extent[1]+img.voxelsize[1]/2,
-img.voxelsize[1]/2, -img.voxelsize[1]/2,
) )
figure.add_subplot(1, 2 ,1) figure.add_subplot(1, 2 ,1)
figure.gca().imshow( figure.gca().imshow(
np.max([img['FieGFP'].get_array(), img['Cas9RFP'].get_array()], axis=0).max(axis=0), np.max([img['FieGFP'].get_array(), img['Cas9RFP'].get_array()], axis=0).max(axis=0),
extent=extent extent=extent
) )
figure.add_subplot(1, 2 ,2) figure.add_subplot(1, 2 ,2)
figure.gca().imshow( figure.gca().imshow(
np.sum([img['FieGFP'].get_array(), img['Cas9RFP'].get_array()], axis=0).max(axis=0), np.sum([img['FieGFP'].get_array(), img['Cas9RFP'].get_array()], axis=0).max(axis=0),
extent=extent extent=extent
) )
``` ```
%% Cell type:markdown id:7ab91cea tags: %% Cell type:markdown id:7ab91cea tags:
### Nuclei detection on fused channels ### Nuclei detection on fused channels
%% Cell type:code id:0d13fd28 tags: %% Cell type:code id:0d13fd28 tags:
``` python ``` python
nuclei_img = SpatialImage( nuclei_img = SpatialImage(
np.sum([img['FieGFP'].get_array(), img['Cas9RFP'].get_array()], axis=0), np.sum([img['FieGFP'].get_array(), img['Cas9RFP'].get_array()], axis=0),
voxelsize=img.voxelsize voxelsize=img.voxelsize
) )
nuclei_points = detect_nuclei( nuclei_points = detect_nuclei(
nuclei_img, nuclei_img,
threshold=detection_threshold, threshold=detection_threshold,
radius_range=detection_radius_range radius_range=detection_radius_range
) )
``` ```
%% Cell type:code id:66ee5adb tags: %% Cell type:code id:66ee5adb tags:
``` python ``` python
figure = plt.figure(figsize=(10, 10)) figure = plt.figure(figsize=(10, 10))
extent = ( extent = (
-img.voxelsize[2]/2, -img.voxelsize[2]/2,
img.extent[2]+img.voxelsize[2]/2, img.extent[2]+img.voxelsize[2]/2,
img.extent[1]+img.voxelsize[1]/2, img.extent[1]+img.voxelsize[1]/2,
-img.voxelsize[1]/2, -img.voxelsize[1]/2,
) )
figure.gca().imshow( figure.gca().imshow(
nuclei_img.get_array().max(axis=0), nuclei_img.get_array().max(axis=0),
extent=extent extent=extent
) )
figure.gca().scatter( figure.gca().scatter(
nuclei_points[:, 0], nuclei_points[:, 0],
nuclei_points[:, 1], nuclei_points[:, 1],
color='r' color='r'
) )
``` ```
%% Cell type:markdown id:6e171d20 tags: %% Cell type:markdown id:6e171d20 tags:
### Nuclei signal quantification ### Nuclei signal quantification
%% Cell type:code id:7c5f3bb1 tags: %% Cell type:code id:7c5f3bb1 tags:
``` python ``` python
nuclei_signals = { nuclei_signals = {
signal_name: quantify_nuclei_signal_intensity( signal_name: quantify_nuclei_signal_intensity(
img[signal_name], img[signal_name],
nuclei_points, nuclei_points,
nuclei_sigma=quantification_radius nuclei_sigma=quantification_radius
) )
for signal_name in ['FieGFP', 'Cas9RFP'] for signal_name in ['FieGFP', 'Cas9RFP']
} }
``` ```
%% Cell type:code id:a3d7b5c4 tags: %% Cell type:code id:a3d7b5c4 tags:
``` python ``` python
import matplotlib.pyplot as plt import matplotlib.pyplot as plt
figure = plt.figure(figsize=(20, 10)) figure = plt.figure(figsize=(20, 10))
extent = ( extent = (
-img.voxelsize[2]/2, -img.voxelsize[2]/2,
img.extent[2]+img.voxelsize[2]/2, img.extent[2]+img.voxelsize[2]/2,
img.extent[1]+img.voxelsize[1]/2, img.extent[1]+img.voxelsize[1]/2,
-img.voxelsize[1]/2, -img.voxelsize[1]/2,
) )
for i, signal_name in enumerate(['FieGFP', 'Cas9RFP']): for i, signal_name in enumerate(['FieGFP', 'Cas9RFP']):
figure.add_subplot(1, 2, i+1) figure.add_subplot(1, 2, i+1)
figure.gca().imshow( figure.gca().imshow(
img[signal_name].get_array().max(axis=0), img[signal_name].get_array().max(axis=0),
extent=extent, extent=extent,
) )
figure.gca().scatter( figure.gca().scatter(
nuclei_points[:, 0], nuclei_points[:, 0],
nuclei_points[:, 1], nuclei_points[:, 1],
c=nuclei_signals[signal_name], c=nuclei_signals[signal_name],
cmap='Reds' cmap='Reds'
) )
``` ```
%% Cell type:code id:6e7b756f tags: %% Cell type:code id:6e7b756f tags:
``` python ``` python
nuclei_mask = nuclei_points[:, 1]>200 nuclei_mask = nuclei_points[:, 1]>200
#nuclei_mask = np.ones_like(nuclei_points[:, 1]).astype(bool) #nuclei_mask = np.ones_like(nuclei_points[:, 1]).astype(bool)
figure = plt.figure(figsize=(10, 10)) figure = plt.figure(figsize=(10, 10))
figure.gca().scatter( figure.gca().scatter(
nuclei_signals['Cas9RFP'][nuclei_mask], nuclei_signals['Cas9RFP'][nuclei_mask],
nuclei_signals['FieGFP'][nuclei_mask], nuclei_signals['FieGFP'][nuclei_mask],
c=nuclei_points[:, 1][nuclei_mask], c=nuclei_points[:, 1][nuclei_mask],
cmap="Reds", cmap="Reds",
) )
``` ```
%% Cell type:markdown id:95a97827 tags: %% Cell type:markdown id:95a97827 tags:
### Data export ### Data export
%% Cell type:code id:0bea4d68 tags: %% Cell type:code id:0bea4d68 tags:
``` python ``` python
nuclei_data = { nuclei_data = {
'id': range(len(nuclei_points)), 'id': range(len(nuclei_points)),
} }
nuclei_data.update({ nuclei_data.update({
f'center_{dim}': nuclei_points[:, i] f'center_{dim}': nuclei_points[:, i]
for i, dim in enumerate('xyz') for i, dim in enumerate('xyz')
}) })
nuclei_data.update(nuclei_signals) nuclei_data.update(nuclei_signals)
nuclei_df = pd.DataFrame(nuclei_data) nuclei_df = pd.DataFrame(nuclei_data)
nuclei_df.head() nuclei_df.head()
``` ```
%% Cell type:code id:487153ec tags: %% Cell type:code id:487153ec tags:
``` python ``` python
nuclei_filename = f"{dirname}/{filename}_nuclei_data.csv" nuclei_filename = f"{dirname}/{filename}_nuclei_data.csv"
nuclei_df.to_csv(nuclei_filename, index=False) nuclei_df.to_csv(nuclei_filename, index=False)
``` ```
%% Cell type:markdown id:e7b650b8 tags: %% Cell type:markdown id:e7b650b8 tags:
### Individual nuclei slices ### Individual nuclei slices
%% Cell type:code id:95d74f9c tags: %% Cell type:code id:95d74f9c tags:
``` python ``` python
nucleus_window = 30 nucleus_window = 30
nucleus_extent = ( nucleus_extent = (
(-nucleus_window-0.5)*img.voxelsize[2], (-nucleus_window-0.5)*img.voxelsize[2],
(nucleus_window+0.5)*img.voxelsize[2], (nucleus_window+0.5)*img.voxelsize[2],
(nucleus_window+0.5)*img.voxelsize[1], (nucleus_window+0.5)*img.voxelsize[1],
(-nucleus_window-0.5)*img.voxelsize[1], (-nucleus_window-0.5)*img.voxelsize[1],
) )
for nucleus_id in range(len(nuclei_points)): for nucleus_id in range(len(nuclei_points)):
nucleus_coords = (nuclei_points[nucleus_id][::-1]/img.voxelsize).astype(int) nucleus_coords = (nuclei_points[nucleus_id][::-1]/img.voxelsize).astype(int)
figure = plt.figure(figsize=(20, 10)) figure = plt.figure(figsize=(20, 10))
for i, signal_name in enumerate(['FieGFP', 'Cas9RFP']): for i, signal_name in enumerate(['FieGFP', 'Cas9RFP']):
figure.add_subplot(1, 2, i+1) figure.add_subplot(1, 2, i+1)
nucleus_img = img[signal_name].get_array()[nucleus_coords[0], nucleus_img = img[signal_name].get_array()[
nucleus_coords[1]-40:nucleus_coords[1]+40, nucleus_coords[0],
nucleus_coords[2]-40:nucleus_coords[2]+40] nucleus_coords[1]-nucleus_window:nucleus_coords[1]+nucleus_window,
nucleus_coords[2]-nucleus_window:nucleus_coords[2]+nucleus_window
]
figure.gca().imshow( figure.gca().imshow(
nucleus_img, nucleus_img,
extent=nucleus_extent, extent=nucleus_extent,
) )
figure.gca().set_title(signal_name, size=18) figure.gca().set_title(signal_name, size=18)
figure.suptitle(f"Nucleus {nucleus_id}", size=24) figure.suptitle(f"Nucleus {nucleus_id}", size=24)
figure.tight_layout() figure.tight_layout()
``` ```
%% Cell type:code id:e32bc393 tags: %% Cell type:code id:e32bc393 tags:
``` python ``` python
``` ```
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment