MATLAB library developed to produce the results for the paper
and for the previous version on
The code has been tested successfully with older versions of MATLAB (at least R2016b (9.1)) and with one of the latest version (R2020a (9.8)) on Windows.
Good starting point is running the figure-generating or movie-generating scripts, as they reproduce the main results and simulations. In each figure class the subfigure plots are generated by a separate function where applicable.
In general, there are two major ways for running simulations - through a lineage tree, where cell divisions are taking place, or directly on a fixed-size grid. For the first scenario, the default way of setting up would be the following:
neigh_class = @neighbourhood.hop_x_grid; % set the neighbourhood class
neigh_range = 2; % communication range
model_class = @models.u_v_dn_model; % model/circuit
% the cell_division object that combines the above settings and starts from a single cell, running for 7 cell divisions
cd = cell_division(1,7,'neigh_class',neigh_class,'model_class',model_class,'neigh_range',neigh_range);
cd.grid_simulation_cell_division(); % simulation
fig = cd.merge_lintree_ratio_plot(); % plot the lineage tree with a ratio plot aboveThe simulations require choosing a 'neighbourhood', which defines how the cells on a grid interact with each other, hence how the intercellular communication is defined. More specifically, for each cell, it defines from which neighbors on the grid the cell is receiving signal, typically defined by the communicating range. In the example above communication is established between any cells on distance less or equal than 2, i.e. at most two hops away. The simulations also require choosing a 'model', i.e. a causality circuit that describes how the genes interact with each other, how the signaling input affects the circuit dynamics and in turn, how the signal secretion is regulated. In the example above the default model of u-v mutual inhibition is chosen.
The state changes and the evolution of spatial patterns on the grids during the cell divisions can also be animated:
cd.animation();Alternatively, summary information can be plotted, containing the stable spatial configurations at the end of each cell cycle, the ratio plot, and the fraction of u+ cells and the u+ cluster radius at the end of each cell cycle:
cd.make_results_plot();To run a simulation directly on a grid with fixed size (in this case 4x8 grid), the following code can be used:
neigh = neighbourhood.hop_2_grid(4,8);
model = models.u_v_dn_model(neigh);
ms = model_simulation(model);
ms.stochastic = true; % switch between deterministic and stochastic simulation
ms.integrate();The results can be visualized showing the time traces of each cell in its respective location (for grids up to 8x8 size), together with the marginal projections on the u-v state space:
ms.plot_grid();The parameters of the model can also be changed to set the system in the appropriate regime, for example:
model.update_parameters(struct('alfa1', 2.52));Likewise, the distribution from which initial conditions are sampled can be defined, for example:
ics_mean = [1.0,2.0,1.5]';
ics_std = 0.15;
model.set_initial_cond(ics_mean, ics_std);