{
    "info": {
        "author": "Ally Hume",
        "author_email": "a.hume@epcc.ed.ac.uk",
        "bugtrack_url": null,
        "classifiers": [
            "Development Status :: 3 - Alpha",
            "License :: OSI Approved :: MIT License",
            "Operating System :: OS Independent",
            "Programming Language :: Python :: 3.5",
            "Programming Language :: Python :: 3.6",
            "Programming Language :: Python :: 3.7",
            "Programming Language :: Python :: 3.8",
            "Topic :: Scientific/Engineering"
        ],
        "description": "## Introduction\n\nBoolean Delay Equations (BDEs) can be used to model a variety of problems.  ```pybde``` allows\nto you write Boolean delay equations models in Python and simulate them.\n\nMore detailed documentation can be found at: https://github.com/EPCCed/pybde/wiki/pybde\n\nCode for the examples included here can be found at: https://github.com/EPCCed/pybde-examples\n\n## Install pybde\n\n```pybde``` requires Python 3.5 or above.  \n\nYou can install `pybde` using `pip`:\n\n```\npip install pybde\n```\n\n## Writing and simulating a model\n\nThe first model we will simulate has a single variable and a single delay.  The single equation \nfor the model is:\n\nx(t) = NOT x(t-\u00cf\u201e)\n\nwhere x is the model variable, t is time, and \u00cf\u201e is the time delay.  In our example \u00cf\u201e = 1.\n\nTo implement this model we must write a function that when given the state of the model's\nvariables at each of the time delays returns the state of the model's variables at the current time.\nThe argument to this function is a list of lists.  If the argument is called `z` then `z[i][j]` contains\nthe state of the _j_ th variable at the _i_ th delay (note that indexing starts from 0).\n\nSo in this case we can write our model in the following function:\n\n```\ndef my_model(z):\n   return [ not z[0][0] ]\n```\n\nTo simulate this model we must provide:\n* the history of the state variables prior to the start of the simulation,\n* the time delays, and\n* the end time for the simulation.\n\nOur model has only one variable and we will specify its history from t=0 until t=1.  We \ndefine this as a Boolean time series specifying:\n* a list of time points where the state changes, \n* a corresponding list of the new variable state at each of these time points, and\n* the final time point for the time series.\n\nThe code to do this is:\n\n```\nhistory = BooleanTimeSeries([0], [True], 1)\n```\n\nWe only have a single delay parameter in this model and its value is 1 so the delay_parameters list is:\n```\ndelay_parameters = [ 1 ]\n```\n\nOur simulation will run from the end of the defined history (t=1) and will end at time t=5:\n```\nend_time = 5\n```\n\nNote that the history must last at least a long as the maximum delay parameter.  In this case both are 1 seconds\nso this is valid.\n\nPutting this altogether gives:\n\n```\nfrom pybde import BDESolver, BooleanTimeSeries\n\n\ndef my_model(z):\n    return [not z[0][0]]\n\n\ndef main():\n    history = BooleanTimeSeries([0], [True], 1)\n    delay_parameters = [1]\n    end_time = 5\n\n    my_bde_solver = BDESolver(my_model, delay_parameters, [history])\n    my_bde_solver.solve(end_time)\n    my_bde_solver.show_result()\n\n\nif __name__ == \"__main__\":\n    main()\n```\n\nThis will display the following plot showing the state of the variable\nover the duration of the simulation.\n\n![One variable one delay output](https://github.com/EPCCed/pybde/wiki/images/v1.0/single_variable_output.png)\n\n## Multiple variables and delays\n\nIn this example our model will contain two variable and two delays.  The model\nequations are:\n\nx1(t) = x2(t-\u00cf\u201e1)\n\nx2(t) = NOT x1(t-\u00cf\u201e2)\n\nwhere x1 and x2 are the model variables, t is time, and \u00cf\u201e1 and \u00cf\u201e2 are the time delays.\nIn this example example \u00cf\u201e1 = 1 and \u00cf\u201e2 = 0.5.\n\nWe implement this model with the following function.  Note that the first index \nspecifies the delay and the second index specifies the variable.  Here we have\nexplicitly named index variables so the code looks more like the equations\nexpressed above.\n\n```\ndef my_two_variable_model(z):\n    x1 = 0\n    x2 = 1\n    tau1 = 0\n    tau2 = 1\n    return [z[tau1][x2], not z[tau2][x1]]\n```\n\nWe wish to start the simulation at t=2 with input history until this point as shown below:\n\n![Two variables, two delays history](https://github.com/EPCCed/pybde/wiki/images/v1.0/two_variables_history.png)\n\nSo we specify the history of variables x1 and x2 as:\n\n```\nx1_history = BooleanTimeSeries([0, 1.5], [True, False], 2)\nx2_history = BooleanTimeSeries([0, 1], [True, False], 2)\n```\n\nTo distinguish the variables when plotting results we can give them\nlabels and matlplotlib plotting styles:\n\n```\nx1_history.label = \"x1\"\nx1_history.style = \"-r\"\nx2_history.label = \"x2\"\nx2_history.style = \"-b\"\n```\n\nSo the full simulation is run with the following code:\n\n```\nfrom pybde import BDESolver, BooleanTimeSeries\n\n\ndef my_two_variable_model(z):\n    x1 = 0\n    x2 = 1\n    tau1 = 0\n    tau2 = 1\n    return [z[tau1][x2], not z[tau2][x1]]\n\n\ndef main():\n    x1_history = BooleanTimeSeries([0, 1.5], [True, False], 2)\n    x2_history = BooleanTimeSeries([0, 1], [True, False], 2)\n\n    x1_history.label = \"x1\"\n    x1_history.style = \"-r\"\n    x2_history.label = \"x2\"\n    x2_history.style = \"-b\"\n\n    delay_parameters = [1, 0.5]\n\n    end_time = 6\n\n    my_bde_solver = BDESolver(my_two_variable_model, delay_parameters,[x1_history, x2_history])\n    my_bde_solver.solve(end_time)\n    my_bde_solver.show_result()\n\n\nif __name__ == \"__main__\":\n    main()\n```\n\nThis will display the following plot showing the state of the variables\nover the duration of the simulation.\n\n![Two variables, two delays output](https://github.com/EPCCed/pybde/wiki/images/v1.0/two_variables_output.png)\n\n\n## Forcing inputs\n\nForcing inputs are input variables that must be specified for the whole duration of\nthe simulation.  These variables' state are not determined by model equations but\ncan be used within model equations to determine the state of other variables.\n\nAs an example of forced inputs consider the following model equations:\n\nx1(t) = 1 if t mod 1 >= 0.5, 0 otherwise   \n\nx2(t) = x1(t-\u00cf\u201e)                        \n\nwhere x1 and x2 are model state variables, t is the time and \u00cf\u201e is the delay.\nIn this case \u00cf\u201e is 0.3.\n\nHere we can model x1 as a forcing input as we can define the value of x1 for\nthe whole duration of the simulation. To specify x1 we must define the\nstarting state and switch points for the whole simulation.  \n\nFor a three second simulation this can be defined as:\n\n```\nx1_input = BooleanTimeSeries([0, 0.5, 1, 1.5, 2, 2.5, 3], [False], 3)\n```\n\nNote that in the code above we only specify the initial Boolean state rather than\nall seven Boolean states.  If the length of the state list is smaller than the length\nof the timepoint list then the state list is simply extended with alternative True\nand False values.  The above code is identical to:\n\n```\nx1_input = BooleanTimeSeries([0, 0.5, 1, 1.5, 2, 2.5, 3], [False, True, False, True, False, True, False], 3)\n```\n\n\nGiven that x1 is a forcing variable the only normal variable will be\nx2.  The following code initialises it to `True` for 0.5 time units:\n\n```\nx2_history = BooleanTimeSeries([0], [True], 0.5)\n```\n\nThe inputs to the simulation are shown in the following plot:\n\n![Forcing input before simulation](https://github.com/EPCCed/pybde/wiki/images/v1.0/forcing_input_before.png)\n\nWhen using forcing inputs the state of forcing inputs at the various\ntime delays is passed to the model function as a second \nargument.  The model function is therefore:\n\n```\ndef my_forced_input_model(z, forced_inputs):\n    tau = 0\n    x1 = 0\n    return [ forced_inputs[tau][x1] ]\n```\n\nTo run the simulation is very similar to before except the forcing inputs\nmust be passed when constructing the ```BDESolver``` object:\n\n```\nmy_bde_solver = BDESolver(my_forcing_input_model, \n                          delay_parameters, \n                          [x2_history], \n                          [x1_input])\n```\n\nThe whole code is:\n\n```\nfrom pybde import BDESolver, BooleanTimeSeries\n\n\ndef my_forcing_input_model(z, forcing_inputs):\n    tau = 0\n    x1 = 0\n    return [ forcing_inputs[tau][x1] ]\n\n\ndef main():\n    x2_history = BooleanTimeSeries([0], [True], 0.5)\n    x2_history.label = 'x2'\n    x2_history.style = '-r'\n\n    x1_input = BooleanTimeSeries([0, 0.5, 1, 1.5, 2, 2.5, 3], [False], 3)\n    x1_input.label = 'x1'\n    x1_input.style = '-b'\n\n    delay_parameters = [0.3]\n    end_time = 3\n\n    my_bde_solver = BDESolver(my_forcing_input_model, \n                              delay_parameters, \n                              [x2_history], \n                              [x1_input])\n    my_bde_solver.solve(end_time)\n\n    my_bde_solver.show_result()\n\n\nif __name__ == \"__main__\":\n    main()\n```\n\n\nRunning this simulation produces the following plot:\n\n![Forcing input after simulation](https://github.com/EPCCed/pybde/wiki/images/v1.0/forcing_input_after.png)\n\n## Plotting and printing result data\n\nThe `BDESolver` class provides basic methods to plot or print results. These\ncan be useful to quickly see the result of a simulation. For more detailed\nanalysis of the results see the `BooleanTimeSeries` convenience functions\nbelow.\n\n### `plot_result()`\n\n`plot_result` plots the variable state and any forcing inputs to a single\nplot. Line styles and labels are taken from the `BooleanTimeSeries` objects \npassed to the solver. The plot will not be displayed.  Use `matplotlib` functions\nto display or save the plot.  For example:\n\n```\nimport matplotlib.pyplot as plt\n\n...\n\nplt.figure(figsize=(5, 2))\nmy_bde_solver.plot_result()\nplt.savefig(\"result_plot.png\")\n```\n\n### `show_result()`\n\n`show_result` is similar to `plot_result` except the `show()` function is called\nto display the plot.\n\n### `print_result(file=sys.stdout)`\n\n`print_result` prints the state of the variables at each switch point producing\noutput such as:\n\n```\n    0.00 ->     1.00 : T T \n    1.00 ->     1.50 : T F \n    1.50 ->     2.00 : F F \n    2.00 ->     3.00 : F T \n    3.00 ->     3.50 : T T \n    3.50 ->     4.50 : T F \n    4.50 ->     5.00 : F F \n    5.00 ->     6.00 : F T \n    6.00 ->     6.00 : T T \n```\n\nBy default the method prints to standard output but alternative outputs can be\nspecified using the `file` argument.\n\n\n## Obtaining the result data\n\nThe `solve` method of `BDESolver` returns a list containing a `BooleanTimeSeries`\nobject for each of the variables.\n\nYou can obtain and process the results with code such as:\n\n```\nresult = my_bde_solver.solve(end_time)\nfor bts in result:\n    print(bts)\n```\n\nOr you can explicitly obtain the time series for each variable using code such as:\n\n```\n[x1_result, x2_result] = my_bde_solver.solve(end_time)\n```\n\n\n\n## `BooleanTimeSeries` convenience functions\n\nThe `BooleanTimeSeries` class includes various convenience functions that help\nprocessing and manipulating Boolean time series data.  These are documented\nhere.\n\n### BooleanTimeSeries(list_of_switch_point_times, list_of_variable_state, end_time, label=None, style=None)\n\nThe `BooleanTimeSeries` constructor takes a list of switch point times,\na list of the new variable state at each of these times and the end_time of the\ntime series. These values are represent the state of the Boolean time series.\n\nThe `list_of_switch_point_times` parameter may be a list of numeric values or\na numpy array of numeric values.\n\nThe `list_of_variable_state` may be a list of `bool` values or a numpy array\nof `bool` values. To save specifying a list of alternating `True` and `False`\nvalues it is possible to specify a list with just the first state value and\nthis will automatically be padded out with alternating Boolean values for\neach specified switch point.\n\nThe optional `label` parameter specifies a label to use when plotting the data.\nThe value also be accessed and set using the class's `label` attribute.\n\nThe optional `style` parameter specifies a style to use when plotting the data.\nThe value also be accessed and set using the class's `style` attribute.\n\n\n### plot(offset=0, scale=1)\n\nPlots the Boolean time series to a `matplotlib` plot. If present the plot\nlabel and line style are taken from the `label` and `style` attributes of this\n`BooleanTimeSeries` instance.\n\nThe plot will not be displayed. To show or save the plot use the appropriate\n`matplotlib` functionality.\n\nThe `offset` parameter can be used to specify an offset from 0 and 1 at which\nto plot the line.  This can be very useful if plotting multiple Boolean time\nseries on the same plot.\n\nThe `scale` parameter can be used to specify that the value to plot for `True`\nis a value other than 1.  This can be useful when plotting Boolean time series\nalongside experimental data.\n\n### show(offset=0, scale=1)\n\n`show` is similar to `plot` expect the matplotlib `show` method will be called\nto display the plot.\n\n### plot_many(list_of_time_series, offset=0.05)\n\nStatic method that plots multiple Boolean time series in a single plot. The offset\nparameter is used to specify the offset between plots in the y axis.\n\nExample of usage:\n\n```\nimport matplotlib.pyplot as plt\n\n...\n\nplt.figure(figsize=(5, 2))\nlist_of_boolean_time_series = my_bde_solver.solve(end_time)\nBooleanTimeSeries.plot_many(list_of_boolean_time_series, offset=0.1)\nplt.savefig(\"result_plot.png\")\n```\n\n\n### show_many(list_of_time_series, offset=0.05)\n\nStatic method that is similar to `plot_many` but calls the `matplotlib` `show`\nfunction to display the plot.\n\n### to_plot_data(offset=0, scale=1)\n\nThe `to_plot_data` method can use used to obtain the Boolean time series in a format suitable\nfor plotting as using various plotting libraries.  The method returns two\nlists: one for x (time) values and the other of y values.\n\nThis method is useful if you wish to take full control over how the results\nare plotted. \n\nThe `offset` parameter can be used to specify an offset from 0 and 1 at which\nto plot the line.  This can be very useful if plotting multiple Boolean time\nseries on the same plot.\n\nThe `scale` parameter can be used to specify that the value to plot for `True`\nis a value other than 1.  This can be useful when plotting Boolean time series\nalongside experimental data.\n\nExample of usage:\n\n```\nfrom pybde import BooleanTimeSeries\n\nbts = BooleanTimeSeries([0, 2, 6, 10], [True], 12)\n\nx, y = bts.to_plot_data()\n\nprint('x = {}'.format(x))\nprint('y = {}'.format(y))\n```\n\nOutputs:\n\n```\nx = [0, 2, 2, 6, 6, 10, 10, 12]\ny = [1, 1, 0, 0, 1, 1, 0, 0]\n```\n\n\n### absolute_threshold(t, y, threshold) \n\nThe static `absolute_threshold` method produces Boolean time series data from\nnumerical time series data. An absolute threshold value is specified above\nwhich the Boolean time series will be `True` and below which the Boolean time\nseries will be `False`.\n\nInput parameter `t` must be either a list of numeric values or a numpy array of\nnumeric values. Input parameter `y` must be either a list of `bool` values\nor a numpy array of `bool` values.\n\nLinear interpolation is used to determine the time at which the state\nchanges.\n\nFor example:\n\n```\nfrom pybde import BooleanTimeSeries\n\nt = [0,  1, 2, 3,  4]\ny = [0, 10, 8, 3, 12]\n\nbts = BooleanTimeSeries.absolute_threshold(t, y, 5)\n\nprint(bts)\n```\n\nproduces:\n\n```\nt=[0, 0.5, 2.6, 3.2222222222222223], y=[False, True, False, True], end=4\n```\n\n### relative_threshold(t, y, threshold)\n\nThe static `relative_threshold` method produces Boolean time series data from\nnumerical time series data. An threshold value is calculated specified above\nwhich the Boolean time series will be `True` and below which the Boolean time\nseries will be `False`.  The absolute threshold value used is calculated\nas `(max(y)-min(y))*threshold + min(y)`.  The specified threshold parameter\nshould be a number between 0 and 1.\n\nInput parameter `t` must be either a list of numeric values or a numpy array of\nnumeric values. Input parameter `y` must be either a list of `bool` values\nor a numpy array of `bool` values.\n\nLinear interpolation is used to determine the time at which the state\nchanges.\n\nFor example:\n\n```\nfrom pybde import BooleanTimeSeries\n\nt = [0,  1, 2, 3,  4]\ny = [4, 10, 8, 2, 12]\n\nbts = BooleanTimeSeries.relative_threshold(t, y, 0.5)\n\nprint(bts)\n```\n\nproduces:\n\n```\nt=[0, 0.5, 2.1666666666666665, 3.5], y=[False, True, False, True], end=4\n```\n\n### cut(new_start, new_end, keep_switch_on_end=False)\n\nThe `cut` method return a new `BooleanTimeSeries` which is a sub-series of the original\nseries. The returned series will run from the specified new start time to the specified\nnew end time. By default a state switch that occurs on the new end time will be omitted,\nthe `keep_switch_on_end` flag can be set to `True` to keep such state switches.\n\nFor example:\n\n```\n> from pybde import BooleanTimeSeries\n> bts = BooleanTimeSeries([0,1,2,3,4,5,6], [True], 7)\n> print(bts)\nt=[0, 1, 2, 3, 4, 5, 6], y=[True, False, True, False, True, False, True], end=7\n\n> print( bts.cut(0,6) )\nt=[0, 1, 2, 3, 4, 5], y=[True, False, True, False, True, False], end=6\n\n> print( bts.cut(0, 6, keep_switch_on_end=True) )\nt=[0, 1, 2, 3, 4, 5, 6], y=[True, False, True, False, True, False, True], end=6\n\n> print( bts.cut(1.5, 4.5) )\nt=[1.5, 2, 3, 4], y=[False, True, False, True], end=4.5\n```\n\n### hamming_distance(boolean_time_series)\n\nThe `hamming_distance` method compares the Boolean Time Series with another \nBoolean time series and returns the total duration for which they differ.\nTwo time series that are identical will have a Hamming distance of zero.\n\nFor example:\n\n```\n> from pybde import BooleanTimeSeries\n> bts = BooleanTimeSeries([0,1,2,3,4,5,6], [True], 7)\n> print(bts.hamming_distance(bts))\n0.0\n\n> bts2 = BooleanTimeSeries([0,1.5,2,3,4.3,5,6], [True], 7)\nprint(bts.hamming_distance(bts2))\n0.8\n```\n\n### merge(list_of_time_series)\n\nThe static `merge` method takes a list of BooleanTimeSeries objects and outputs\ntwo lists.  The first list is the switch point times and the second list is\na list of lists of the state variables at these time points.\n\nFor example:\n\n```\nfrom pybde import BooleanTimeSeries\n\nbts1 = BooleanTimeSeries([0, 1.0, 2.0], [True], 3)\nbts2 = BooleanTimeSeries([0, 1.5, 2.5], [True], 3)\n\nt, y = BooleanTimeSeries.merge([bts1, bts2])\n\nprint('t = {}'.format(t))\nprint('y = {}'.format(y))\n```\n\noutputs:\n\n```\nt = [0, 1.0, 1.5, 2.0, 2.5]\ny = [[True, True], [False, True], [False, False], [True, False], [True, True]]\n```\n\n### unmerge(list_of_switch_timepoints, list_of_lists_of_variable_states, end)\n\nThe static function `unmerge` is the opposite of `merge`. `unmerge` takes as input \na list a switch point times, a list of list of variable states at these\ntime points and the time series end time and returns a list of BooleanTimeSeries objects.\n\nFor example:\n\n```\nfrom pybde import BooleanTimeSeries\n\nt = [0, 1.0, 1.5, 2.0, 2.5]\ny = [[True, True], [False, True], [False, False], [True, False], [True, True]]\n\nfor bts in BooleanTimeSeries.unmerge(t, y, 3):\n    print(bts)\n```\n\noutputs\n\n```\nt=[0, 1.0, 2.0], y=[True, False, True], end=3\nt=[0, 1.5, 2.5], y=[True, False, True], end=3\n```\n\n## Do not include switch points at the end of variable's history\n\nWhen running a simulation the input history time series must not end on a switch point.\nThis is because when the simulation starts from the time point the model equations may\ncontradict the history state at this point. To avoid this simply remove the final switch\npoint from the history.  This can be easily achieved using the `cut` function which\nby default removes any switch point at the end of the time series duration.  For example:\n\n```\n> hist = BooleanTimeSeries([0,1,2], [True], 2)\n> print(hist)\nt=[0, 1, 2], y=[True, False, True], end=2\n> hist = hist.cut(0,hist.end)\n> print(hist)\nt=[0, 1], y=[True, False], end=2\n```\n\n## Logging\n\n`pybde` using Python's [logging library](https://docs.python.org/3/library/logging.html) to provide\nsome debug logging. For example, the following line can be used to turn own debug logging:\n\n```\nimport logging\n\nlogging.basicConfig(level=logging.DEBUG)\n```\n\n## Numerical accuracy\n\nThe implementation of `pydbe` has to compare possible switch times generated in \ndifferent ways to see if they are the same time.  For example, is t1+\u00cf\u201e2 the timepoint\nas t2+\u00cf\u201e1. To perform comparisons of floating point numbers `pydbe` uses [`math.isclose`](https://docs.python.org/3.5/library/math.html#math.isclose) function.  This function\ndefines the acceptable accuracy using the `rel_tol` and `abs_tol` arguments. To specify\nnon-default values for these arguments you can specify `rel_tol` and `abs_tol` arguments\nwhen constructing the `BDESolver` object.\n\nThe `BooleanTimeSeries` class also performs some floating point comparisons and adopts\nthe same approach a `BSESolver`.  To alter the default relative and absolute tolerances\nfor the `BooleanTimeSeries` class set the `rel_tol` and `abs_tol` static attributes of\nthe class.\n\n## Real world example\n\nIn this section show a real world example where a Boolean Delay Equation model is\nused model Circadian rhythms found in biology.  This example is available as an\niPython notebook (neurospora_example_notebook.ipybn) from https://github.com/EPCCed/pybde-examples.\n\nWe have experimentally obtained experssion levels for two genes.  The following\ncode obtains and plots these expression levels:\n\n```\nfrom pybde import BDESolver\nfrom pybde import BooleanTimeSeries\nfrom numpy import genfromtxt\nimport matplotlib.pyplot as plt\n\nexperiment_data = genfromtxt(\n    'https://raw.githubusercontent.com/EPCCed/pybde-examples/master/neurospora_data.csv', \n    delimiter=',')\n    \nplt.plot(experiment_data[:,2], experiment_data[:,0], 'b-', label=\"m\", )\nplt.plot(experiment_data[:,2], experiment_data[:,1], 'r-', label=\"ft\")\nplt.legend()\nplt.xlabel(\"time (hours)\")\nplt.ylabel(\"Expression levels (AU)\")\nplt.title(\"Experiment data\")\nplt.show()\n```\n\n![Experiment data](https://github.com/EPCCed/pybde/wiki/images/v1.0/neurospora_1.png)\n\nTo use this data in a Boolean Delay Equation we need to convert it into Boolean data. We use the `BooleanTimeSeries` class to store Boolean time series data. We can create Boolean time series data from this experiment data by applying relative or absolute thresholding.\n\nRelative thresholding thresholds the relative to the range of the values in the data. A threshold of 0.5 will correspond to a threshold value midway between the minimum and maximum value.\n\nHere we choose a relative thresholdhold values of 0.3.\n\n```\n# Turn experiment data into Boolean time series\nm_bts = \\\n    BooleanTimeSeries.relative_threshold(experiment_data[:,2], experiment_data[:,0], 0.3)\nft_bts = \\\n    BooleanTimeSeries.relative_threshold(experiment_data[:,2], experiment_data[:,1], 0.3)\n\n# Add labels and plot styles to the switch points objects\nft_bts.label = \"ft\"\nm_bts.label = \"m\"\nft_bts.style = \"r-\"\nm_bts.style = \"b-\"\n\n# Plot the experiment data as Boolean time series\nBooleanTimeSeries.plot_many([m_bts, ft_bts])\nplt.xlabel(\"Time (hours)\")\nplt.show()\n```\n\n![Boolean version of experiment data](https://github.com/EPCCed/pybde/wiki/images/v1.0/neurospora_2.png)\n\nWe can plot this over the original data to see how it matches. Note that we rescale the Boolean data purely for display purposes they approximately match the same scale as the experiment data. Remember that the Boolean data only has states True (1) and False (0).\n\n```\n# Now plot the experiment data Boolean time series with the original data\n\nplt.plot(experiment_data[:,2], experiment_data[:,0], 'b-', label=\"m\", )\nm_bts.plot(scale=150)\nplt.xlabel(\"time (hours)\")\nplt.ylabel(\"Expression levels (AU)\")\nplt.title(\"m\")\nplt.show()\n\n\nplt.plot(experiment_data[:,2], experiment_data[:,1], 'r-', label=\"ft\")\nft_bts.plot(scale=300)\nplt.xlabel(\"time (hours)\")\nplt.ylabel(\"Expression levels (AU)\")\nplt.title(\"ft\")\nplt.show()\n```\n\n![Experiment data on Boolean version](https://github.com/EPCCed/pybde/wiki/images/v1.0/neurospora_3.png)\n\n![Experiment data on Boolean version](https://github.com/EPCCed/pybde/wiki/images/v1.0/neurospora_4.png)\n\nWe wish to use the first 24 hours' worth of data as history for our simulation. So we can cut this data to extract the first 24 hours:\n\n```\n# We wish to use the first 24 hours as our history\n\nhist_ft = ft_bts.cut(0, 24)\nhist_m = m_bts.cut(0, 24)\n\n# Plot the history\nBooleanTimeSeries.plot_many([hist_m, hist_ft])\nplt.xlabel(\"Time (hours)\")\nplt.legend()\nplt.show()\n```\n\n![Experiment history data](https://github.com/EPCCed/pybde/wiki/images/v1.0/neurospora_5.png)\n\nThe Circadian model uses light as an input so we need to prepare the light input:\n\n```\n# The model uses a forcing input for light\n\nlight_t = [0] + list(range(6,120,12))\nlight_y = []\nfor t in light_t:\n    light_y.append(6 <= t % 24 < 18)\nlight_bts = BooleanTimeSeries(light_t, light_y, 118)\nlight_bts.label = \"light\"\nlight_bts.style = \"-g\"\n\n# Plot light\nlight_bts.plot()\nplt.show()\n```\n\n![Light data](https://github.com/EPCCed/pybde/wiki/images/v1.0/neurospora_6.png)\n\nSo the inputs to our simulation are can be plotted on one graph.\n\n```\n# Plot light and the histories\nBooleanTimeSeries.plot_many([hist_m, hist_ft, light_bts])\nplt.show()\n```\n\n![Light data](https://github.com/EPCCed/pybde/wiki/images/v1.0/neurospora_7.png)\n\nNow we have to define our Boolean Delay Equation model. The model has two simulated states (ft and m) and one input state (light). The model has three delays.\n\n```\ndef neurospora_eqns(z, forced_inputs):\n    m = 0\n    ft = 1\n    light = 0\n    tau1 = 0\n    tau2 = 1\n    tau3 = 2\n\n    return [ (not z[tau2][ft]) or forced_inputs[tau3][light], z[tau1][m] ]\n```\n\nFor given values of the delays (determined by parameter fitting algorithms such as [pynmmso](https://github.com/EPCCed/pynmmso/wiki/pynmmso)) we can run the simulation for 118 hours:\n\n```\n# Run the simulator\n\ntau1 = 5.0752\ntau2 = 6.0211\ntau3 = 14.5586\ndelays = [tau1, tau2, tau3]\n\nsolver = BDESolver(neurospora_eqns, delays, [hist_m, hist_ft], [light_bts])\n[m_output, ft_output] = solver.solve(118)\n```\n\nWe can now plot the output of the simulation:\n\n```\nplt.title(\"Simulation output\")\nBooleanTimeSeries.plot_many([m_output , ft_output, light_bts])\nplt.xlabel(\"Time (hours)\")\nplt.legend()\nplt.show()\n```\n\n![Simulation result](https://github.com/EPCCed/pybde/wiki/images/v1.0/neurospora_8.png)\n\nPlotting these outputs over the original experiment data:\n\n```\nplt.plot(experiment_data[:,2], experiment_data[:,0], 'b-', label=\"m\", )\nm_output.plot(scale=150)\nplt.xlabel(\"time (hours)\")\nplt.ylabel(\"Expression levels (AU)\")\nplt.title(\"m\")\nplt.show()\n\n\nplt.plot(experiment_data[:,2], experiment_data[:,1], 'r-', label=\"ft\")\nft_output.plot(scale=300)\nplt.xlabel(\"time (hours)\")\nplt.ylabel(\"Expression levels (AU)\")\nplt.title(\"ft\")\nplt.show()\n```\n\n![Results and experiment data](https://github.com/EPCCed/pybde/wiki/images/v1.0/neurospora_9.png)\n\n![Results and experiment data](https://github.com/EPCCed/pybde/wiki/images/v1.0/neurospora_10.png)\n\nPlotting the simulated data alongside the thresholded data gives:\n\n```\n# Plot simulated data alongside thresholded data\n\nm_output.label = \"m sim\"\nm_output.style = 'c-'\nBooleanTimeSeries.plot_many([m_bts, m_output])\nplt.xlabel(\"time (hours)\")\nplt.title(\"m\")\nplt.legend()\nplt.show()\n\nft_output.label = \"ft sim\"\nft_output.style = 'c-'\nBooleanTimeSeries.plot_many([ft_bts, ft_output])\nplt.xlabel(\"time (hours)\")\nplt.title(\"ft\")\nplt.legend()\nplt.show()\n```\n\n![Results and experiment data](https://github.com/EPCCed/pybde/wiki/images/v1.0/neurospora_11.png)\n\n![Results and experiment data](https://github.com/EPCCed/pybde/wiki/images/v1.0/neurospora_12.png)\n\n\nWe can calculate the Hamming distance which gives a measure of what duration of time two Boolean time series have differing signals. For 96 hours of simulated time the Hamming distance measures are low:\n\n```\n# Calculate Hamming distances\n\nprint(\"Hamming of m : {:.4f}\".format(m_bts.hamming_distance(m_output)))\nprint(\"Hamming of ft : {:.4f}\".format(ft_bts.hamming_distance(ft_output)))\n\nprint(\"Hamming of m as %age of simulated time: {:.2f}%\".format(\n    m_bts.hamming_distance(m_output)/0.96))\nprint(\"Hamming of ft as %age of simulated time : {:.2f}%\".format(\n    ft_bts.hamming_distance(ft_output)/0.96))\n```\n\nWhich outputs:\n\n```\nHamming of m : 9.3353\nHamming of ft : 14.7796\nHamming of m as %age of simulated time: 9.72%\nHamming of ft as %age of simulated time : 15.40%\n```\n\n## Using pybde with pynmmso\n\nIn the above real world example the delay values were explicitly given. In practice they often be obtained\nas the results of a parameter optimisation exercise.  Here we show how to use [`pynmmso`](https://github.com/EPCCed/pynmmso/wiki/) can be used to determine the values of the delays that result in the best fit to the experiment data.  \n\nThis code follows on from the above real world example.  The details of pynmmso are explained in the \n[`pynmmso` documentation](https://github.com/EPCCed/pynmmso/wiki). It is important to know that \n`pynmmso` seeks to maximise a fitness function and returns multiple maxima at peaks of the fitness\nlandscape rather then simply the global maximum.\n\nDefine our problem class that provides the fitness function and the parameter bounds:\n\n```\nclass MyProblem:\n    def __init__(self, equations, history, inputs, experiment_data):\n        self.equations = equations\n        self.history = history\n        self.inputs = inputs\n        self.experiment_data = experiment_data\n        \n    def fitness(self, delays):\n        # Simulate with the given delay parameters\n        solver = BDESolver(self.equations, delays, self.history, self.inputs)\n        [m_output, ft_output] = solver.solve(118)\n        \n        # Compare the result with the Boolean version of the experiment data\n        cost = m_output.hamming_distance(self.experiment_data[0]) + ft_output.hamming_distance(self.experiment_data[1])\n        \n        # Because pynmmso looks for maxima we negate the cost so it looks for the smallest\n        # cost\n        return -cost\n        \n    def get_bounds(self):\n        # For each of the three delay values explore the range 1..20\n        return [1, 1, 1],[20,20,20]\n```\n\nNow use `pynmmso` to find the maxima:\n\n```\nnumber_of_fitness_evaluations = 50000\n\nnmmso = Nmmso(MyProblem(neurospora_eqns, [hist_m, hist_ft], [light_bts], [m_bts, ft_bts]))\nmy_result = nmmso.run(number_of_fitness_evaluations)\nfor mode_result in my_result:\n    print(\"Mode at {} has value {}\".format(mode_result.location, mode_result.value))\n```\n\nThis will produce output similar to the following although the actual values returned\nmay be slightly different:\n\n```\nMode at [ 9.99999999 19.99999998 15.17078072] has value -48.860801508507805\nMode at [19.99999996  7.27463383 15.46914558] has value -73.18395433102842\nMode at [ 5.64100973  7.90184814 13.70187647] has value -19.722341133743335\nMode at [19.95614228  3.5720228   1.03870254] has value -90.91490635269514\nMode at [20.         11.40358844 16.59000251] has value -74.00207006980949\nMode at [19.99915781 15.64418847 15.86933091] has value -75.92303065045431\nMode at [19.99998653 19.44037368  1.07194677] has value -95.11612551719436\nMode at [19.99999884  9.72558647  1.00212142] has value -86.6759957429595\nMode at [20.         11.40358844 18.00923951] has value -74.00207006980948\nMode at [19.87677477  3.40047693 18.00923951] has value -76.36160166053793\nMode at [19.87677477  3.40047693 16.94265386] has value -75.85875891998987\nMode at [ 4.45554284 19.27618803 16.02055511] has value -50.81186043628735\nMode at [1.00000066 4.99999932 1.        ] has value -83.91461625103807\nMode at [ 5.3953705   4.27897001 15.24473766] has value -18.939534918789544\n```\n\nThe two solutions with a value of -18.9 and -19.7 look primising and have a lower combined Hamming distance score than the parameters originally used in the example above.\n\n## Acknowledgements\n\nThis work was supported by the Engineering and Physical Sciences Research Council (grant number [EP/N018125/1](https://gow.epsrc.ukri.org/NGBOViewGrant.aspx?GrantRef=EP/N018125/1))",
        "description_content_type": "text/markdown",
        "docs_url": null,
        "download_url": "",
        "downloads": {
            "last_day": -1,
            "last_month": -1,
            "last_week": -1
        },
        "home_page": "https://github.com/EPCCed/pybde/wiki/pybde",
        "keywords": "BDE,boolean,delay,equations,solver,simulator",
        "license": "MIT",
        "maintainer": "",
        "maintainer_email": "",
        "name": "pybde",
        "package_url": "https://pypi.org/project/pybde/",
        "platform": "",
        "project_url": "https://pypi.org/project/pybde/",
        "project_urls": {
            "Homepage": "https://github.com/EPCCed/pybde/wiki/pybde"
        },
        "release_url": "https://pypi.org/project/pybde/0.4.3/",
        "requires_dist": null,
        "requires_python": ">=3.5",
        "summary": "Boolean Delay Equation simulator",
        "version": "0.4.3"
    },
    "last_serial": 5505763,
    "releases": {
        "0.1": [
            {
                "comment_text": "",
                "digests": {
                    "md5": "c3be7dd30fc7a9966e6fc672a14c2eea",
                    "sha256": "7d33e64e5007b24cfb8ac99e061eda33a4d09ec4e29a43511c8fb69aa743546e"
                },
                "downloads": -1,
                "filename": "pybde-0.1.tar.gz",
                "has_sig": false,
                "md5_digest": "c3be7dd30fc7a9966e6fc672a14c2eea",
                "packagetype": "sdist",
                "python_version": "source",
                "requires_python": ">=3.5",
                "size": 14569,
                "upload_time": "2019-06-05T14:11:22",
                "url": "https://files.pythonhosted.org/packages/ee/00/f2b955438cb53b586046b8f318574528d3ac87f8a435b0aeed52d80eb315/pybde-0.1.tar.gz"
            }
        ],
        "0.2": [
            {
                "comment_text": "",
                "digests": {
                    "md5": "f7043098cab570e5a5bce6a66afaf2a8",
                    "sha256": "9e33c7f7bbf72b2d05e3de8eefc22da2255466346ce08796cdf5e4382c72fffa"
                },
                "downloads": -1,
                "filename": "pybde-0.2.tar.gz",
                "has_sig": false,
                "md5_digest": "f7043098cab570e5a5bce6a66afaf2a8",
                "packagetype": "sdist",
                "python_version": "source",
                "requires_python": ">=3.5",
                "size": 14550,
                "upload_time": "2019-06-05T14:13:57",
                "url": "https://files.pythonhosted.org/packages/9d/3c/eae39c3c344da3c3e7770adb8686077097f149f45f46c65305d12e52dc0f/pybde-0.2.tar.gz"
            }
        ],
        "0.3": [
            {
                "comment_text": "",
                "digests": {
                    "md5": "4c92a628378cb27e4fd6e820f4fd8053",
                    "sha256": "1d369664fba1fb73699a184ff17e82861d47e8cd5164f36d650e6a5f05111bfa"
                },
                "downloads": -1,
                "filename": "pybde-0.3.tar.gz",
                "has_sig": false,
                "md5_digest": "4c92a628378cb27e4fd6e820f4fd8053",
                "packagetype": "sdist",
                "python_version": "source",
                "requires_python": ">=3.5",
                "size": 14545,
                "upload_time": "2019-06-05T14:19:20",
                "url": "https://files.pythonhosted.org/packages/26/38/bb5681b7cf5ebfa326d1052a1d91fa8295178ffe75c2838d190d115513ce/pybde-0.3.tar.gz"
            }
        ],
        "0.4": [
            {
                "comment_text": "",
                "digests": {
                    "md5": "2e0364ca3708fa8dd3a7b4a55e85e414",
                    "sha256": "dc1172c7c21d8e446e28ec3ca0343fc0d84c0e61ee2f542581dd7597b74dcfae"
                },
                "downloads": -1,
                "filename": "pybde-0.4.tar.gz",
                "has_sig": false,
                "md5_digest": "2e0364ca3708fa8dd3a7b4a55e85e414",
                "packagetype": "sdist",
                "python_version": "source",
                "requires_python": ">=3.5",
                "size": 21500,
                "upload_time": "2019-06-28T15:47:00",
                "url": "https://files.pythonhosted.org/packages/3b/87/bb419c9de6acf922ee6155a4df685a5ee83e3873fbf979e1aca32a4f9ff4/pybde-0.4.tar.gz"
            }
        ],
        "0.4.1": [
            {
                "comment_text": "",
                "digests": {
                    "md5": "96a092f8131bf6de867d445fe93fd959",
                    "sha256": "9f9d405886564f7bd88e753ee5219e66e9a15f069b4b5659903abd0f600d7d81"
                },
                "downloads": -1,
                "filename": "pybde-0.4.1.tar.gz",
                "has_sig": false,
                "md5_digest": "96a092f8131bf6de867d445fe93fd959",
                "packagetype": "sdist",
                "python_version": "source",
                "requires_python": ">=3.5",
                "size": 21469,
                "upload_time": "2019-06-28T15:53:00",
                "url": "https://files.pythonhosted.org/packages/1f/dd/30e8a45a8d8704441569f0e9d43f23864300e90a0cc0691787eed43a940d/pybde-0.4.1.tar.gz"
            }
        ],
        "0.4.2": [
            {
                "comment_text": "",
                "digests": {
                    "md5": "9b1cfa621f5339fdd6da2264d9acf642",
                    "sha256": "80e1a36c902c154007f8b651352c0ad636816e3c70b1f24fec2306dde5b0c1f5"
                },
                "downloads": -1,
                "filename": "pybde-0.4.2.tar.gz",
                "has_sig": false,
                "md5_digest": "9b1cfa621f5339fdd6da2264d9acf642",
                "packagetype": "sdist",
                "python_version": "source",
                "requires_python": ">=3.5",
                "size": 29880,
                "upload_time": "2019-07-02T15:49:55",
                "url": "https://files.pythonhosted.org/packages/96/af/1cb7ae7820aedb8692016511791e8e0749bb1d6efb8d684683088bf4484b/pybde-0.4.2.tar.gz"
            }
        ],
        "0.4.3": [
            {
                "comment_text": "",
                "digests": {
                    "md5": "6a49b59f0bd45298079389bc23ab378e",
                    "sha256": "f265bec773db68ee7d01e2cb4ac0b6c83bedf7062b655e346bd498a82abd64c8"
                },
                "downloads": -1,
                "filename": "pybde-0.4.3.tar.gz",
                "has_sig": false,
                "md5_digest": "6a49b59f0bd45298079389bc23ab378e",
                "packagetype": "sdist",
                "python_version": "source",
                "requires_python": ">=3.5",
                "size": 35780,
                "upload_time": "2019-07-09T09:21:15",
                "url": "https://files.pythonhosted.org/packages/ae/c7/c40f4a9d057d345b5164fd64f35de7799abec5f880cf184cfe0e9bc85e1e/pybde-0.4.3.tar.gz"
            }
        ]
    },
    "urls": [
        {
            "comment_text": "",
            "digests": {
                "md5": "6a49b59f0bd45298079389bc23ab378e",
                "sha256": "f265bec773db68ee7d01e2cb4ac0b6c83bedf7062b655e346bd498a82abd64c8"
            },
            "downloads": -1,
            "filename": "pybde-0.4.3.tar.gz",
            "has_sig": false,
            "md5_digest": "6a49b59f0bd45298079389bc23ab378e",
            "packagetype": "sdist",
            "python_version": "source",
            "requires_python": ">=3.5",
            "size": 35780,
            "upload_time": "2019-07-09T09:21:15",
            "url": "https://files.pythonhosted.org/packages/ae/c7/c40f4a9d057d345b5164fd64f35de7799abec5f880cf184cfe0e9bc85e1e/pybde-0.4.3.tar.gz"
        }
    ]
}