added plothist and plot_cuml_divisions_per_cell

plotting.py

+from copy import copy
+
 import matplotlib.cm
 import matplotlib.pylab as plt
 import numpy as np
@@ -290,3 +292,155 @@ def plot_replicate_array(
         plt.grid()
     if show:
         plt.show(block=False)
+
+
+def plothist(
+    t,
+    da,
+    num_fine=200,
+    vmax_factor=1,
+    norm="log",
+    edges=None,
+    title=None,
+    figsize=None,
+):
+    """
+    Plot the distribution of fluorescence as a function of time.
+
+    Time is on the x-axis; the support of the distribution is on the y-axis;
+    and shading represents the height of the distribution at each y value.
+
+    Parameters
+    ----------
+    t: 1D array
+        An array of time points.
+    da: 2D array
+        An array of the fluorescence data, with each row containing data from
+        a single cell.
+    num_fine: integer (optional)
+        The number of data points to include in a finer spacing of time points.
+        Linear interpolation is used to find the extra data values.
+    vmax_factor: float (optional)
+        Used to rescale the maximal data value, which is then passed to the vmax
+        argument in plt.pcolormesh to change the range of the colour bar.
+    norm: str (optional)
+        Passed to the norm argument in plt.pcolormesh for mapping the data
+        values to the colour map. Default is "log".
+    edges: list of arrays (optional)
+        Specifies the bins in time and in fluorescence values and is used to plot
+        different data sets on axes with the same range.
+    title: str (optional)
+        Title for the plot.
+    figsize: tuple (optional)
+        Sets the width and height of the figure.
+
+    Examples
+    --------
+    Load data:
+    >>> from wela.dataloader import dataloader
+    >>> from wela.plotting import plothist
+    >>> dl = dataloader()
+    >>> dl.load("1334_2023_03_28_pyrTo2Gal_01glc_00", use_tsv=True)
+    >>> dlc = dataloader()
+    >>> dlc.load("1005_2023_03_09_exp_00_pyrToGalInduction, use_tsv="True")
+
+    Get time-series data:
+    >>> t, d = dl.get_time_series("median_GFP")
+    >>> tc, dc = dlc.get_time_series("median_GFP")
+
+    Plot both data sets using axes with the same range:
+    >>> edges = plothist(t, dc, title="2% Gal", figsize=(4, 3))
+    >>> plothist(t, d, title="2% Gal and 0.1% Glu", edges=edges, figsize=(4, 3))
+    """
+    # interp to a new grid and replace NaN
+    t_fine = np.linspace(t.min(), t.max(), num_fine)
+    s_fine = np.empty((da.shape[0], num_fine), dtype=float)
+    for i in range(da.shape[0]):
+        s = da[i, :]
+        s_fine[i, :] = np.interp(t_fine, t[~np.isnan(s)], s[~np.isnan(s)])
+    t_fine = np.matlib.repmat(t_fine, da.shape[0], 1)
+    # make histogram
+    cmap = copy(plt.cm.viridis)
+    cmap.set_bad(cmap(0))
+    if edges is not None:
+        h, xedges, yedges = np.histogram2d(
+            t_fine.flatten(),
+            s_fine.flatten(),
+            bins=edges,
+        )
+    else:
+        h, xedges, yedges = np.histogram2d(
+            t_fine.flatten(),
+            s_fine.flatten(),
+            bins=[int(num_fine / 2), 50],
+        )
+    vmax = np.nanmax(da) / vmax_factor
+    # plot using pcolormesh
+    if figsize is not None:
+        plt.figure(figsize=figsize)
+    else:
+        plt.figure()
+    pcm = plt.pcolormesh(
+        xedges,
+        yedges,
+        h.T,
+        cmap=cmap,
+        rasterized=True,
+        vmax=vmax,
+        norm=norm,
+    )
+    plt.colorbar(pcm, label="number of cells", pad=0.02)
+    plt.xlabel("time (h)")
+    plt.ylabel("fluorescence")
+    plt.ylim(top=1200)
+    if title:
+        plt.title(title)
+    plt.show(block=False)
+    return [xedges, yedges]
+
+
+def plot_cuml_divisions_per_cell(t, buddings, nboots=30, col="b", label=None):
+    """
+    Plot the mean cumulative number of divisions per cell over time.
+
+    No figure is automatically generated. See example.
+
+    Parameters
+    ----------
+    t: 1D array
+        An array of time points.
+    buddings: 2D array
+        A binary array where a budding event is denoted by a one and each row is
+        a time series of budding events for one cell.
+    nboots: integer (optional)
+        The number of bootstraps to use to estimate errors.
+    col: str (optional)
+        Color of the line to plot.
+    label: str (optional)
+        Label for the legend.
+
+    Example
+    -------
+    >>> from wela.plotting import plot_cuml_divisions_per_cell
+    >>> plt.figure()
+    >>> plot_cuml_divisions_per_cell(t, b, label="Gal Glu")
+    >>> plt.show(block=False)
+    """
+
+    def find_cuml(b):
+        return (
+            np.array([np.nansum(b[:, :i]) for i in range(b.shape[1])])
+            / b.shape[0]
+        )
+
+    def sample(b):
+        return np.array(
+            [b[i, :] for i in np.random.randint(0, b.shape[0], b.shape[0])]
+        )
+
+    cuml = find_cuml(buddings)
+    err_cuml = 2 * np.std(
+        [find_cuml(sample(buddings)) for i in range(nboots)], 0
+    )
+    plt.plot(t, cuml, color=col, label=label)
+    plt.fill_between(t, cuml - err_cuml, cuml + err_cuml, color=col, alpha=0.2)