Repository: DavAug/ErlotinibGefitinib

Issue:
define a plot function with allows model exploration with parameter sliders

Discussion:

Patch:
diff --git a/pkpd/plots/__init__.py b/pkpd/plots/__init__.py
index 055e5f0..d251bd5 100644
--- a/pkpd/plots/__init__.py
+++ b/pkpd/plots/__init__.py
@@ -7,7 +7,8 @@
 
 from ._scatter import (  # noqa
     plot_measurements,
-    plot_measurements_and_predictions
+    plot_measurements_and_predictions,
+    plot_measurements_and_simulation
 )
 
 from ._box import (  # noqa
diff --git a/pkpd/plots/_scatter.py b/pkpd/plots/_scatter.py
index b7bc461..68852d1 100644
--- a/pkpd/plots/_scatter.py
+++ b/pkpd/plots/_scatter.py
@@ -33,13 +33,13 @@ def plot_measurements(data):
 
     if not isinstance(data, pd.DataFrame):
         raise ValueError(
-            'Input data <' + str(data) + '> has to be pandas.DataFrame.')
+            'Data has to be pandas.DataFrame.')
 
     keys = ['#ID', 'TIME in day', 'TUMOUR VOLUME in cm^3', 'BODY WEIGHT in g']
     for key in keys:
         if key not in data.keys():
             raise ValueError(
-                'Input data <' + str(data) + '> must have key <' + str(key) +
+                'Data must have key <' + str(key) +
                 '>.')
 
     # Get ids
@@ -211,13 +211,13 @@ def plot_measurements_and_predictions(data, model, parameters):
     # Check data has the correct type
     if not isinstance(data, pd.DataFrame):
         raise TypeError(
-            'Data <' + str(data) + '> has to be pandas.DataFrame.')
+            'Data has to be pandas.DataFrame.')
     # Check that data has the required keys
     keys = ['#ID', 'TIME in day', 'TUMOUR VOLUME in cm^3', 'BODY WEIGHT in g']
     for key in keys:
         if key not in data.keys():
             raise ValueError(
-                'Data <' + str(data) + '> must have key <' + str(key) +
+                'Data must have key <' + str(key) +
                 '>.')
     # Check that model has the correct type
     if not isinstance(model, pints.ForwardModel):
@@ -349,3 +349,316 @@ def plot_measurements_and_predictions(data, model, parameters):
     )
 
     return fig
+
+
+def plot_measurements_and_simulation(
+        data, model, default_parameters, min_parameters, max_parameters,
+        parameter_names=None, steps=50):
+    """
+    Returns a `plotly.graph_objects.Figure` containing the data and simulation
+    with interactive sliders.
+
+    This function assumes the following keys naming convention:
+        ids: '#ID'
+        time: 'TIME in day'
+        tumour volume: 'TUMOUR VOLUME in cm^3
+
+    The axis labels as well as the hoverinfo assume that time is measured in
+    day, volume is measured in cm^3.
+
+    Arguments:
+        data -- A pandas.DataFrame containing the measured time-series data of
+                the tumour volume and the mass.
+        model -- A `pints.ForwardModel`.
+        parameters -- An array-like object with the model parameters for each
+                      individual in the dataset.
+                      Shape: (n_individuals, n_parameters)
+    """
+    # Check data has the correct type
+    if not isinstance(data, pd.DataFrame):
+        raise TypeError(
+            'Data has to be pandas.DataFrame.')
+    # Check that data has the required keys
+    keys = ['#ID', 'TIME in day', 'TUMOUR VOLUME in cm^3', 'BODY WEIGHT in g']
+    for key in keys:
+        if key not in data.keys():
+            raise ValueError(
+                'Data must have key <' + str(key) +
+                '>.')
+    # Check that model has the correct type
+    if not isinstance(model, pints.ForwardModel):
+        raise TypeError(
+            'Model needs to be an instance of `pints.ForwardModel`.')
+    # Check that model has only one output dimension
+    if model.n_outputs() != 1:
+        raise ValueError(
+            'Model output dimension has to be 1.')
+    # Check that parameters have the correct dimensions
+    default_parameters = np.asarray(default_parameters)
+    if default_parameters.ndim != 1:
+        raise ValueError(
+            'Default parameters need to have dimension 1. `default_parameters`'
+            ' has dimension <' + str(default_parameters.ndim) + '>.')
+    min_parameters = np.asarray(min_parameters)
+    if min_parameters.ndim != 1:
+        raise ValueError(
+            'Minimum parameters need to have dimension 1. `min_parameters` has'
+            ' dimension <' + str(min_parameters.ndim) + '>.')
+    max_parameters = np.asarray(max_parameters)
+    if max_parameters.ndim != 1:
+        raise ValueError(
+            'Maximum parameters need to have dimension 1. `max_parameters` has'
+            ' dimension <' + str(max_parameters.ndim) + '>.')
+    # Check that parameters have the correct shape
+    if default_parameters.shape != (model.n_parameters(),):
+        raise ValueError(
+            'Default parameters do not have the correct shape. Expected '
+            'shape (n_parameters,) = ' + str((model.n_parameters(),)) + '.')
+    if min_parameters.shape != (model.n_parameters(),):
+        raise ValueError(
+            'Minimum parameters do not have the correct shape. Expected shape '
+            '(n_parameters,) = ' + str((model.n_parameters(),)) + '.')
+    if max_parameters.shape != (model.n_parameters(),):
+        raise ValueError(
+            'Maximum parameters do not have the correct shape. Expected shape '
+            '(n_parameters,) = ' + str((model.n_parameters(),)) + '.')
+
+    if parameter_names is None:
+        parameter_names = [
+            'Parameter %d' % n for n in range(model.n_parameters())]
+
+    # Check parameter names
+    if len(parameter_names) != model.n_parameters():
+        raise ValueError(
+            'Number of parameter names does not match number of parameters.')
+
+    # Define colorscheme
+    n_ids = len(data['#ID'].unique())
+    colors = plotly.colors.qualitative.Plotly[:n_ids]
+
+    # Create figure
+    fig = go.Figure()
+
+    # Plot simulation for each slider step
+    _add_slider_step_plots(
+        fig, data, model, default_parameters, min_parameters, max_parameters,
+        steps)
+
+    # Scatter plot of control growth data
+    _add_data(fig, data, colors)
+
+    # Set X, Y axis and figure size
+    fig.update_layout(
+        autosize=True,
+        xaxis_title=r'$\text{Time in day}$',
+        yaxis_title=r'$\text{Tumour volume in cm}^3$',
+        template="plotly_white")
+
+    # Create parameter sliders
+    sliders = _create_sliders(
+        parameter_names, min_parameters, max_parameters, steps, data)
+
+    # Add switch between linear and log y-scale
+    fig.update_layout(
+        sliders=sliders,
+        updatemenus=[
+            dict(
+                type="buttons",
+                direction="left",
+                buttons=list([
+                    dict(
+                        args=[{"yaxis.type": "linear"}],
+                        label="Linear y-scale",
+                        method="relayout"
+                    ),
+                    dict(
+                        args=[{"yaxis.type": "log"}],
+                        label="Log y-scale",
+                        method="relayout"
+                    )
+                ]),
+                pad={"r": 0, "t": -10},
+                showactive=True,
+                x=0.0,
+                xanchor="left",
+                y=1.15,
+                yanchor="top"
+            ),
+            dict(
+                type="buttons",
+                direction="down",
+                buttons=list([
+                    dict(
+                        args=[{
+                            "sliders[0].visible": True,
+                            "sliders[1].visible": False,
+                            "sliders[2].visible": False}],
+                        label="Initial volume",
+                        method="relayout"
+                    ),
+                    dict(
+                        args=[{
+                            "sliders[0].visible": False,
+                            "sliders[1].visible": True,
+                            "sliders[2].visible": False}],
+                        label="Exp. growth rate",
+                        method="relayout"
+                    ),
+                    dict(
+                        args=[{
+                            "sliders[0].visible": False,
+                            "sliders[1].visible": False,
+                            "sliders[2].visible": True}],
+                        label="Lin. growth rate",
+                        method="relayout"
+                    )
+                ]),
+                pad={"r": 0, "t": -10},
+                showactive=True,
+                x=1.07,
+                xanchor="left",
+                y=-0.1,
+                yanchor="top"
+            ),
+        ]
+    )
+
+    return fig
+
+
+def _add_data(fig, data, colors):
+    """
+    Adds a scatter plot of the data to the figure.
+    """
+    ids = data['#ID'].unique()
+    for index, id_m in enumerate(ids):
+        # Create mask for mouse
+        mask = data['#ID'] == id_m
+
+        # Get observed data for indiviudal
+        observed_times = data['TIME in day'][mask].to_numpy()
+        observed_data = data['TUMOUR VOLUME in cm^3'][mask]
+
+        # Plot data
+        fig.add_trace(go.Scatter(
+            x=observed_times,
+            y=observed_data,
+            legendgroup="ID: %d" % id_m,
+            name="ID: %d" % id_m,
+            showlegend=True,
+            hovertemplate=(
+                "<b>Measurement </b><br>" +
+                "ID: %d<br>" % id_m +
+                "Time: %{x:} day<br>" +
+                "Tumour volume: %{y:.02f} cm^3<br>" +
+                "Cancer type: LXF A677<br>" +
+                "<extra></extra>"),
+            mode="markers",
+            marker=dict(
+                symbol='circle',
+                opacity=0.7,
+                line=dict(color='black', width=1),
+                color=colors[index])
+        ))
+
+
+def _add_slider_step_plots(
+        fig, data, model, parameters, min_parameters, max_parameters, steps):
+    """
+    Add a plot for each slider step to the figure.
+    """
+    # Define time range based on data
+    start_experiment = data['TIME in day'].min()
+    end_experiment = data['TIME in day'].max()
+    simulated_times = np.linspace(
+        start=start_experiment, stop=end_experiment)
+
+    # Add plot for each slider step
+    for param_id in range(model.n_parameters()):
+        # Get default parameters
+        params = parameters.copy()
+
+        # Create line plot for each slider value
+        for value_id, value in enumerate(np.linspace(
+                start=min_parameters[param_id],
+                stop=max_parameters[param_id], num=steps)):
+            # Update parameters
+            params[param_id] = value
+
+            # Solve model
+            simulated_data = model.simulate(params, simulated_times)
+
+            # Show only first plot
+            visible = False
+            if param_id == 0 and value_id == 10:
+                visible = True
+
+            # Plot slider step
+            fig.add_trace(go.Scatter(
+                x=simulated_times,
+                y=simulated_data,
+                legendgroup="Model",
+                name="Model",
+                showlegend=True,
+                hovertemplate=(
+                    "<b>Simulation </b><br>" +
+                    "Time: %{x:.0f} day<br>" +
+                    "Tumour volume: %{y:.02f} cm^3<br>" +
+                    "Cancer type: LXF A677<br>" +
+                    "<br>" +
+                    "<b>Parameters </b><br>" +
+                    "Initial tumour volume: %.02f cm^3<br>" % params[0] +
+                    "Exp. growth rate: %.02f 1/day<br>" % params[1] +
+                    "Lin. growth rate: %.02f cm^3/day<br>" % params[2] +
+                    "<extra></extra>"),
+                visible=visible,
+                mode="lines",
+                line=dict(color='Black')
+            ))
+
+
+def _create_sliders(
+        parameter_names, min_parameters, max_parameters, steps, data):
+    """
+    Returns slider objects that can be used in plotly's `update_layout`.
+    """
+    # Get number of parameters and measured individuals
+    n_params = len(parameter_names)
+    n_ids = len(data['#ID'].unique())
+
+    # Create slider object for each parameter
+    sliders = []
+    for param_id in range(n_params):
+        # Compute parameter values
+        values = np.linspace(
+            start=min_parameters[param_id], stop=max_parameters[param_id],
+            num=steps)
+
+        # Display plot for slider position and data
+        slider_steps = []
+        for step in range(steps):
+            # All simulations False and data True
+            s = dict(
+                method="update",
+                args=[{
+                    "visible": [False] * (n_params * steps) +
+                    [True] * n_ids}],
+                label='%.02f' % values[step])
+
+            # Make correct simulation visible
+            s["args"][0]["visible"][param_id * steps + step] = True
+
+            # Safe slider steps
+            slider_steps.append(s)
+
+        # Save slider
+        sliders.append(
+            dict(
+                visible=True if param_id == 0 else False,
+                active=10,
+                currentvalue=dict(
+                    prefix=parameter_names[param_id] + ': '),
+                pad={"t": 50},
+                steps=slider_steps))
+
+    return sliders


Test patch:
diff --git a/pkpd/tests/test_plots_scatter.py b/pkpd/tests/test_plots_scatter.py
index 9e589c7..b9073e9 100644
--- a/pkpd/tests/test_plots_scatter.py
+++ b/pkpd/tests/test_plots_scatter.py
@@ -44,14 +44,16 @@ class TestPlotMeasurements(unittest.TestCase):
         data = np.ones(shape=(10, 4))
 
         self.assertRaisesRegex(
-            ValueError, 'Input data <', pkpd.plots.plot_measurements, data)
+            ValueError, 'Data has to be pandas.DataFrame.',
+            pkpd.plots.plot_measurements, data)
 
     def test_bad_column_keys(self):
 
         data = self.data.rename(columns={'TIME in day': 'SOMETHING ELSE'})
 
         self.assertRaisesRegex(
-            ValueError, 'Input data <', pkpd.plots.plot_measurements, data)
+            ValueError, 'Data must have key <',
+            pkpd.plots.plot_measurements, data)
 
     def test_create_figure(self):
 
@@ -90,7 +92,7 @@ class TestPlotMeasurementsAndPredictions(unittest.TestCase):
         data = np.ones(shape=(10, 4))
 
         self.assertRaisesRegex(
-            TypeError, 'Data <',
+            TypeError, 'Data has to be pandas.DataFrame.',
             pkpd.plots.plot_measurements_and_predictions, data, self.model,
             self.parameters)
 
@@ -99,7 +101,7 @@ class TestPlotMeasurementsAndPredictions(unittest.TestCase):
         data = self.data.rename(columns={'TIME in day': 'SOMETHING ELSE'})
 
         self.assertRaisesRegex(
-            ValueError, 'Data <',
+            ValueError, 'Data must have key <',
             pkpd.plots.plot_measurements_and_predictions, data, self.model,
             self.parameters)
 
@@ -159,5 +161,159 @@ class TestPlotMeasurementsAndPredictions(unittest.TestCase):
         self.assertIsInstance(fig, go.Figure)
 
 
+class TestPlotMeasurementsAndSimulation(unittest.TestCase):
+    """
+    Tests the `pkpd.plots.plot_mesurements_and_simulation` function.
+    """
+
+    @classmethod
+    def setUpClass(cls):
+        # Create test dataset
+        ids = [0, 0, 0, 1, 1, 1, 2, 2]
+        times = [0, 1, 2, 2, np.nan, 4, 1, 3]
+        volumes = [np.nan, 0.3, 0.2, 0.5, 0.1, 0.2, 0.234, 0]
+        masses = [1, 1, 1, 1, 1, 1, 1, 1]
+        cls.data = pd.DataFrame({
+            '#ID': ids,
+            'TIME in day': times,
+            'TUMOUR VOLUME in cm^3': volumes,
+            'BODY WEIGHT in g': masses})
+
+        # Create test model
+        path = pkpd.ModelLibrary().get_path('Tumour growth without treatment')
+        cls.model = pkpd.PharmacodynamicModel(path)
+
+        # Create test parameters
+        cls.default_params = [0.1, 0.1, 0.1]
+        cls.min_params = [0.01, 0.02, 0.03]
+        cls.max_params = [1, 2, 3]
+
+    def test_bad_data_type(self):
+
+        data = np.ones(shape=(10, 4))
+
+        self.assertRaisesRegex(
+            TypeError, 'Data has to be pandas.DataFrame.',
+            pkpd.plots.plot_measurements_and_simulation, data, self.model,
+            self.default_params, self.min_params, self.max_params)
+
+    def test_bad_column_keys(self):
+
+        data = self.data.rename(columns={'TIME in day': 'SOMETHING ELSE'})
+
+        self.assertRaisesRegex(
+            ValueError, 'Data must have key <',
+            pkpd.plots.plot_measurements_and_simulation, data, self.model,
+            self.default_params, self.min_params, self.max_params)
+
+    def test_bad_model_type(self):
+
+        model = 'bad model'
+
+        self.assertRaisesRegex(
+            TypeError, 'Model needs to be an instance of',
+            pkpd.plots.plot_measurements_and_simulation, self.data, model,
+            self.default_params, self.min_params, self.max_params)
+
+    def test_bad_model_n_output(self):
+
+        # Create model
+        path = pkpd.ModelLibrary().get_path('Tumour growth without treatment')
+        model = pkpd.PharmacodynamicModel(path)
+
+        # Change to multi-output model
+        model.set_outputs(['myokit.tumour_volume', 'myokit.tumour_volume'])
+
+        self.assertRaisesRegex(
+            ValueError, 'Model output dimension has to be 1.',
+            pkpd.plots.plot_measurements_and_simulation, self.data, model,
+            self.default_params, self.min_params, self.max_params)
+
+    def test_bad_default_parameters_dimension(self):
+
+        parameters = np.ones(shape=(3, 3))
+
+        self.assertRaisesRegex(
+            ValueError, 'Default parameters need to have dimension 1.',
+            pkpd.plots.plot_measurements_and_simulation, self.data, self.model,
+            parameters, self.min_params, self.max_params)
+
+    def test_bad_min_parameters_dimension(self):
+
+        parameters = np.ones(shape=(3, 3))
+
+        self.assertRaisesRegex(
+            ValueError, 'Minimum parameters need to have dimension 1.',
+            pkpd.plots.plot_measurements_and_simulation, self.data, self.model,
+            self.default_params, parameters, self.max_params)
+
+    def test_bad_max_parameters_dimension(self):
+
+        parameters = np.ones(shape=(3, 3))
+
+        self.assertRaisesRegex(
+            ValueError, 'Maximum parameters need to have dimension 1.',
+            pkpd.plots.plot_measurements_and_simulation, self.data, self.model,
+            self.default_params, self.min_params, parameters)
+
+    def test_bad_default_parameters_shape(self):
+
+        parameters = np.ones(shape=(10,))
+
+        self.assertRaisesRegex(
+            ValueError, 'Default parameters do not have the correct shape.',
+            pkpd.plots.plot_measurements_and_simulation, self.data, self.model,
+            parameters, self.min_params, self.max_params)
+
+    def test_bad_min_parameters_shape(self):
+
+        parameters = np.ones(shape=(10,))
+
+        self.assertRaisesRegex(
+            ValueError, 'Minimum parameters do not have the correct shape.',
+            pkpd.plots.plot_measurements_and_simulation, self.data, self.model,
+            self.default_params, parameters, self.max_params)
+
+    def test_bad_max_parameters_shape(self):
+
+        parameters = np.ones(shape=(10,))
+
+        self.assertRaisesRegex(
+            ValueError, 'Maximum parameters do not have the correct shape.',
+            pkpd.plots.plot_measurements_and_simulation, self.data, self.model,
+            self.default_params, self.min_params, parameters)
+
+    def test_parameter_names(self):
+
+        # Test bad parameter names
+        param_names = ['param 1', 'param 2']
+
+        self.assertRaisesRegex(
+            ValueError, 'Number of parameter names does not match number of',
+            pkpd.plots.plot_measurements_and_simulation, self.data, self.model,
+            self.default_params, self.min_params, self.max_params, param_names)
+
+        # Good parameter names
+        param_names = ['param 1', 'param 2', 'param 3']
+        fig = pkpd.plots.plot_measurements_and_simulation(
+            self.data, self.model, self.default_params, self.min_params,
+            self.max_params, param_names)
+
+        # Check all labels
+        sliders = fig.layout.sliders
+        self.assertEqual(len(sliders), 3)
+        self.assertEqual(sliders[0].currentvalue.prefix, param_names[0] + ': ')
+        self.assertEqual(sliders[1].currentvalue.prefix, param_names[1] + ': ')
+        self.assertEqual(sliders[2].currentvalue.prefix, param_names[2] + ': ')
+
+    def test_create_figure(self):
+
+        fig = pkpd.plots.plot_measurements_and_simulation(
+            self.data, self.model, self.default_params, self.min_params,
+            self.max_params)
+
+        self.assertIsInstance(fig, go.Figure)
+
+
 if __name__ == '__main__':
     unittest.main()