Time series analysis

Time series analysis#

Download this Notebook Open in Colab Open in Colab 
# /// script
# requires-python = ">=3.10"
# dependencies = [
#     "matplotlib",
#     "numpy",
#     "scikit-image",
#     "scipy",
#     "tifffile",
#     "imagecodecs",
#     "pandas",
#     "seaborn",
#     "bobiac_tools @ git+https://github.com/bobiac/bobiac-tools.git"
# ]
# ///

from pathlib import Path

import matplotlib.pyplot as plt
import pandas as pd
import seaborn as sns
import skimage
import tifffile
from bobiac_tools import overlay_labels

Overview#

In this notebook, we analyse time resolved data.

We have data of cells where the location of our gene of interes oscillates between the nucleus and the cytoplasm. Our question is What is the frequency of the oscillations?

We have a total of 16 time points with a temporal spacing of 30 min (🤔). The background of the images was measured to be 4700.

You can find the data here.

Step 3#

  1. Create a new code cell in the notebook.

  2. Download the dataset folder that contains the time-lapse file F01_1615_tcyx.tif.

  3. For the F01_1615_tcyx.tif file, display the time points of the last channel as a grid of images (or use ndv if you have it in the dependencies).

  4. In the same dataset folder, there is a F01_1615_tcyx_cell_labels.tif file (labeled mask of the cells) and a F01_1615_tcyx_nuclei_labels.tif file (labeled mask of the nuclei). Visualize these masks as an overlay

image_path = Path("../../_static/images/quant/03_time_series/images/F01_1615_tcyx.tif")

time_stack = tifffile.imread(image_path)
print(f"Time stack shape: {time_stack.shape}")  # (16, 3, 1040, 1392)

t_max = time_stack.shape[0]
ch = 2  # channel index for the oscillating gene of interest
fig = plt.figure(figsize=(10, 8), layout="constrained")
for t in range(t_max):
    image = time_stack[t, ch, :, :]
    plt.subplot(4, 4, t + 1)
    plt.imshow(image, cmap="gray")
    plt.title(f"t={t}", fontsize=20)
    plt.axis("off")
plt.show()

These are the first eight time points.

mask_cells = tifffile.imread(
    "../../_static/images/quant/03_time_series/masks/F01_1615_tcyx_cell_labels.tif"
)
mask_nuclei = tifffile.imread(
    "../../_static/images/quant/03_time_series/masks/F01_1615_tcyx_nuclei_labels.tif"
)

overlay_labels(label_mask=[mask_cells, mask_nuclei])

This are the two label masks for the nuclei and the cells.

Step 4#

Measure nuclear intensity for the first time point

In this step, you will analyse the first image and quantify the intensity of each nucleus.

  1. Load the first time point of the last channel (index 2) in the F01_1615_tcyx.tif file.

  2. Remove objects touching the border of the image.

  3. Measure the intensity of each nucleus with skimage.measure.regionprops_table.

  4. Subtract the background intensity from the measurements and inspect the resulting values.

1. Measure nuclear intensity for the first time point#

First, we measure the intensities of all the nuclei in the image.

Load the first time point.

ch = 2  # channel index for the oscillating gene of interest
image_t0 = time_stack[0, ch, :, :]
overlay_labels(image=image_t0)

Remove boundary objects.

mask_nuclei = skimage.segmentation.clear_border(mask_nuclei)
overlay_labels(image=image_t0, label_mask=mask_nuclei)

Measure intensities with regionprops_table ans save measurements as DataFrame.

properties = ["area", "intensity_mean", "label"]
props = skimage.measure.regionprops_table(
    mask_nuclei, intensity_image=image_t0, properties=properties
)

df_nuc = pd.DataFrame(props)
print(df_nuc)

Subtract background intensity from measurements. Background was measured to 460.

background_intensity = 460
df_nuc["intensity_bg_cor"] = df_nuc["intensity_mean"] - background_intensity
print(df_nuc["intensity_bg_cor"].head())
print(
    "\nAre there any negative values after background subtraction?",
    (df_nuc["intensity_bg_cor"] < 0).any(),
)

Step 5#

Perform quality control and analyse all time points

Now you will turn the workflow into a batch analysis over all 16 images.

  1. Check the measurements for artifacts, for example by plotting the distribution of background-corrected intensities.

  2. Overlay the measured values on the image and inspect whether any nuclei look suspicious.

  3. Loop through all time points, measure the intensity of each nucleus, and add a time column.

  4. Combine all measurements into a single DataFrame.

Quality control on background-corrected intensities#

sns.histplot(
    data=df_nuc,
    x="intensity_bg_cor",
    bins=20,
)

Seems unsuspicious.

overlay_labels(
    image=image_t0,
    label_mask=mask_nuclei,
    df=df_nuc,
    id_col="label",
    measurement_col="intensity_bg_cor",
)

Looks fine.

Batch processing#

Add time column.

t = 0  # index of the first time point in the stack
time = t * 0.5
print(time)

df_nuc["time"] = time

list_df = []

for t in range(t_max):
    image = time_stack[t, ch, :, :]
    props = skimage.measure.regionprops_table(
        mask_nuclei, intensity_image=image, properties=properties
    )
    sdf = pd.DataFrame(props)
    sdf["intensity_bg_cor"] = sdf["intensity_mean"] - background_intensity
    sdf["time"] = t * 0.5
    list_df.append(sdf)

df = pd.concat(list_df)
print(df)

Step 6#

Visualise the oscillations and estimate the frequency

In this final step, you will answer the main question: what is the frequency of the oscillations?

  1. Use sns.lineplot to visualise the intensity of each nucleus over time.

  2. Add a second plot showing the average signal across cells.

  3. Estimate the oscillation frequency from the pattern you observe.

  4. Write a short conclusion in a markdown cell.

Plotting#

df.groupby(["label", "time"])["intensity_bg_cor"].mean().reset_index()

sns.lineplot(data=df, x="time", y="intensity_bg_cor", errorbar="sd")

sns.lineplot(data=df, x="time", y="intensity_bg_cor", hue="label")