ForSys: non-invasive stress inference from time-lapse microscopy

Augusto Borges 1, 2, Jerónimo R. Miranda-Rodríguez 1, 3, Alberto Sebastián Ceccarelli 4, Guilherme Ventura5, Jakub Sedzinski 5, Hernán López-Schier 1,2,6 & Osvaldo Chara 7, 8, 9

1. Unit Sensory Biology, Helmholtz Zentrum München, Munich, Germany

2. Graduate School of Quantitative Biosciences (QBM), Munich, Germany

3. Instituto de Neurobiología, Universidad Nacional Autónoma de México (UNAM), Boulevard Juriquilla 3001, Juriquilla, México

4. Systems Biology Group (SysBio), Institute of Physics of Liquids and Biological Systems (IFLySIB), National Scientific and Technical Research Council (CONICET), University of La Plata, La Plata, Argentina

5. The Novo Nordisk Foundation Center for Stem Cell Medicine (reNEW), University of Copenhagen, Blegdamsvej 3B, 2200, Copenhagen, Denmark

6. Division of Science, New York University Abu Dhabi, Saadiyat Island, United Arab Emirates

7. School of Biosciences, University of Nottingham, Sutton Bonington Campus, Nottingham, LE12 5RD, UK

8. Instituto de Tecnología, Universidad Argentina de la Empresa, Buenos Aires, Argentina

9. Corresponding author: osvaldo.chara@nottingham.ac.uk

This notebook will reproduce the contents of dynamic_in_vivo.ipynb to show how to incorporate a connections JSON file into the pipeline for the first two frames

Generate an inference in a time series of microscopy frame

import sys
sys.path.append('..')
import forsys as fs
import os
import matplotlib.pyplot as plt
import numpy as np

DATA_FOLDER = os.path.join("data", "in_vivo")
RESULTS_FOLDER = os.path.join("results")

max_time = 2

# This is only necessary if you wish to create outputs
# if not os.path.exists(RESULTS_FOLDER):
#     os.makedirs(RESULTS_FOLDER)

Create all the frames by iterating in time.

The real_time is used to get the \(\Delta t\) in the velocity calculation

frames = {}
for time in range(max_time):
    real_time = 1 * time

    current_tif_file = os.path.join(DATA_FOLDER, f"t_{time}", f"t_{time}.tif")
    skeleton = fs.skeleton.Skeleton(current_tif_file, mirror_y=False)
    vertices, edges, cells = skeleton.create_lattice()
    vertices, edges, cells, _ = fs.virtual_edges.generate_mesh(vertices, edges, cells, ne=5)
    frames[time] = fs.frames.Frame(time,
                                   vertices,
                                   edges,
                                   cells,
                                   time=real_time)

initial_guess = fs.auxiliar.load_initial_guess(os.path.join("data", "in_vivo", "connections.json"), 0, max_time)

forsys = fs.ForSys(frames, initial_guess=initial_guess)

Build and solve the system of equations for the force and pressure

You could solve only the time of interest or all the times with a loop. The b_matrix parameter turns the Dynamic instance on

for time in range(max_time):
    forsys.build_force_matrix(when=time)
    forsys.solve_stress(when=time, allow_negatives=False, b_matrix="velocity")
    forsys.build_pressure_matrix(when=time)
    forsys.solve_pressure(when=time, method="lagrange_pressure")

all_velocities = forsys.get_system_velocity_per_frame()
ave_velocity_to_normalize = np.mean(all_velocities)
for time in range(max_time):
    forsys.build_force_matrix(when=time)
    forsys.solve_stress(when=time,
                        b_matrix="velocity",
                        velocity_normalization= 0.1 / ave_velocity_to_normalize,
                        adimensional_velocity=False,
                        method="nnls",
                        use_std=False,
                        allow_negatives=False)
    forsys.build_pressure_matrix(when=time)
    forsys.solve_pressure(when=time, method="lagrange_pressure")

Create system’s plots

fig, axs = plt.subplots(max_time, 1, figsize=(8, 3 * max_time))

for time, current_ax in enumerate(axs.reshape(-1)):
    _, ax = fs.plot.plot_inference(forsys.frames[time],
                                    normalized="max",
                                    colorbar=False,
                                    pressure=True,
                                    ax=current_ax)
    current_ax = ax
    current_ax.set_title(f"Time: {time}")
    current_ax.axis("off")


plt.tight_layout()
plt.show()