becksteinlab / kda Goto Github PK

Description

There is currently a test which verifies the correct number of partial and directional diagrams are generated for maximally-connected graphs of various sizes, but this can be expanded to include graphs of arbitrary complexity using Kirchhoff's Theorem. The paper that inspired this can be found here. This is important because we currently only test edge cases, but it would be better to include many of the graphs we use throughout our tests and analysis. Some of this is discussed in #22.

TODO

• Add a function to graph_utils.pydiagrams.py which calculates any cofactor of the Laplacian matrix for an input graph and outputs the number of spanning trees (partial diagrams) for that graph
• Add a test which leverages this function and compares expected values to the output values for many moderately-connected graphs

As mentioned in #34, the plotting function documentation could use some examples. This will have to be handled after #46 of course, but the examples would likely be useful for anyone planning on using the built-in plotting utilities.

KDA currently uses edge attributes to store/retrieve the edge weights and names. For example, say a graph has an edge connecting nodes 1 -> 2 with an edge weight of 10, the current approach will assign an attribute dictionary similar to the following:

attr_dict = {"name": "k12", "val": 10}

This works "fine", but at this point seems a little verbose/redundant when we have to use the edge indices (i.e. (0, 1) for this case) to retrieve this data in the first place. At that point we might as well just create the string using the indices (i.e. f"k{i+1}{j+1}").

Additionally, if we actually store the edge weights the way NetworkX intended, I believe we will have the ability to leverage more of their built-in utilities. That could prove useful if there is ever a spanning tree algorithm implemented there. I know there is already some work there, but last I checked there wasn't anything I could leverage. Most seemed to be min/max spanning tree-related things.

Lastly, I think we are currently iterating over the edges to collect their weights since there is no way to collect the values for a given attribute simultaneously. At least, of the top of my head, I believe that statement is correct. But the point is since we aren't storing that info as an attribute it should be easier to access and could potentially speedup some of the manual calculation functions, such as kda.calculations.calc_state_probs().

With the plotting code being changed substantially in gh-45, the documentation is out of date. There are several undocumented helper functions and several new parameters in the public functions that need to be detailed.

The idea is to use the new precedent set by kda.calculate_thermo_force() of incorporating the string-to-SymPy conversions into the calculation functions themselves instead of having separate functions for this. The original idea was to give the user the maximum amount of flexibility possible, but now it seems giving access to the string versions of functions doesn't really add any value. It would be simpler from a user's perspective to be given the SymPy functions that can then be substituted/simplified and converted to lambda functions.

As discussed in 6. of this comment, the single-diagram/cycle cases of draw_cycles() and draw_diagrams() probably belong elsewhere, potentially in a separate, public function. It seems a little awkward to input a single diagram/cycle into plural-named functions anyways, although a general plotting function might want to handle this case.

This will have to be resolved after #47, and likely after #48.

At the bottom of the kda.diagrams page, the _combine() function shows up in the list. This function should be hidden, and instead there should be an entry for kda.diagrams.enumerate_partial_diagrams().

Travis CI has changed their testing platform to only allow 1000 minutes of testing time per month for free users, which may become intractable once the testing suite is fully developed.

Can view a guide for the migration here.

Change License from LGPLv3 to GPLv3.

Motivation

Largely inspired by @ianmkenney's feedback, it would be useful to have a built-in thermodynamic consistency checker which verifies rates are internally consistent for each cycle in a kinetic diagram. If such a function can be written that is sufficiently low cost, it could be checked by a user at any time, and forced to run before any calculations are executed.

Proposed Solution(s)

We discussed a few solutions, each with their own difficulties:

1. Use sum of free energy differences around the cycle to check thermodynamic consistency. Since the sum of free energy differences around a kinetic cycle must be equal to zero, this would be one way to handle it. Something like this has already been implemented in the EmrE thermodynamic consistency module.
2. Another method would be to require a user to input units into their expressions and use SymPy units to find intrinsic rates based on their units. The external parameters should be easy to find since they should be the only rates expressed in units of molar concentration.

Potential Pitfalls

Solution 1 has a caveat -- I think the functions used for EmrE are operating on substitution dictionaries (for SymPy) which do not contain the external concentrations. A general approach would need to distinguish between internal/external parameters on its own so a user can implement complicated external parameters (e.g. membrane potential) without breaking the thermodynamic consistency checker. This might require us to make our API more explicit about internal/external parameters, which we make no distinction between at the moment.

Solution 2 may prove to be cumbersome if users have to input units for all expressions since the expressions will eventually have to be simplified by SymPy, which is already very costly.

Motivation

The KDA API could use a great deal of improvement. Right now things are a bit clunky and can lead easily to errors.

For example, most diagrams and calculations require the user to use many different functions that all require similar information (this is targeted more specifically in issue #58). There should be a way to minimize the burden on the user while providing the same output. Moving towards an object-oriented approach would likely help with this, where the user creates their "universe" (i.e. graph) and can easily retrieve info about the graph by just calling built-in methods.

Another example -- when generating the net cycle fluxes, you need to know the "order" (i.e. direction) for each cycle (see issue #57), which is input as a pair of nodes in the cycle. Then once the net cycle fluxes are generated, there is no real way of linking one cycle to the next, making it difficult to generate the operational fluxes for a given system. If we can find a way to input more data about the system (or create & store it in the diagram object) up front, we could make this process a lot easier.

API Ideas

Here is an example of a more approachable API where the graph is read in from a graph file (maybe something similar to the multibind graph.csv and states.csv files):

    import kda
    # retrieve graph data from file
    G = kda.read_graph("./test_graph.csv")
    # plot the directional diagrams
    G.directional.plot()
    # get the state probability expressions
    state_probs = G.probability.collect_expressions(norm=False)

For the net transition fluxes an API like the following would prove very useful:

    import kda
    # retrieve graph data from file
    G = kda.read_graph("./test_graph.csv")
    # get the expression for the one-way transition flux j_56
    j_56_expr = G.flux(i=5, j=6, net=False, expr=True)
    # get the expression for the net transition flux from 5->6
    J_56_expr = G.flux(i=5, j=6, net=True, expr=True)

As mentioned in 4. of this comment, there are several plotting functions that could use unit tests in plotting.py. At the moment the plotting test doesn't properly "test" anything, they primarily flush the code paths to make sure there aren't any errors raised. But adding some unit tests would be a good start.

It was found that for a given list of nodes, there can be a host of different orders that create unique cycles. This is illustrated quite plainly with the 10-state model of EmrE transport from the Henzler-Wildman Lab. This model is sufficiently complex (many possible pathways for most nodes) that when listing all possible cycles, the NetworkX simple cycles function is insufficient. For simple combinations of 4 nodes, it works fine, but for collections of 8 or 10 nodes it will only give 1 out of 6 or 1 out of 4 unique cycles, respectively. This may not be an issue for all use cases, but when calculating operational cycle fluxes this is an issue. There needs to be a method of finding every relevant cycle for a given operation so that these fluxes can be calculated and summed appropriately.

test_KDA.py keeps failing on Travis because it cannot import the packages properly. Make these package imports explicit (should be aware of package dependencies) and verify Travis is importing them correctly. Package import must be explicitly listed in both setup.py and Travis.yml.

kda.diagrams.generate_partial_diagrams() currently appends the partial diagram edge lists to a list instead of filling in an array because it was previously difficult to quantify the number of partial diagrams/spanning trees that could be generated from an input diagram. This is now pretty simple to calculate with the addition of kda.diagrams.enumerate_partial_diagrams(), given the rate matrix is provided.

To conserve the current call signature for generate_partial_diagrams(), a call to kda.graph_utils.retrieve_rate_matrix() should suffice for retrieving the rate matrix. With this it should be possible to preallocate an array for the partial diagram edges, which should result in a reasonable speedup and a small reduction in spacial complexity. The final array should have the shape P x (N-1) x M, where P is the number of partial diagrams (calculated from enumerate_partial_diagrams()), N is the number of nodes in the input diagram, and M is the length of the edge tuples (2 or 3).

The returned edge list is already converted to an array, so this should be a pretty straight-forward change.

When using cbt=True in kda.plotting.plot_diagrams() most of the nodes will be labeled correctly, but at least 1 will be mislabeled every time.

It would be nice to have a conftest.py for storing global pytest fixtures in test_kda.py, such as StateProbs3, StateProbs4, etc. This would just clean up the primary testing module a bit, and allow other modules access to them if need be.

Recent CI runs have been giving the following warning:

kda/tests/test_kda.py: 1340 warnings
  C:\Users\nikol\OneDrive\projects\kinetic_diagram_analysis\kda\kda\diagrams.py:342: DeprecationWarning: `product` is deprecated as of NumPy 1.25.0, and will be removed in NumPy 2.0. Please use `prod` instead.
    K_cof = _get_cofactor_matrix(K_laplace=K_laplace)

A simple change from np.product() to np.prod() should fix this.

Problem

kda.generate_flux_diagrams() currently operates by generating all possible flux diagrams for G, finding the flux diagrams that pertain to the input cycle, and then using those. The problem is this takes a long time because of all the permutations it has to sort through, and since every calculation is specific to a given cycle, there is no reason to generate all flux diagrams every time one runs this function.

Solution

• Remove most of the functionality from kda.generate_flux_diagrams() and create alternate function (kda.generate_all_flux_diagrams()) that simply takes in G as input and creates all the flux diagrams. This could be used mostly for confirmation/testing.
• Replace existing code with an alternate function that uses an input diagram G and an input cycle to construct only the flux diagrams for that particular cycle. This should run much quicker and be easier to read.

• Fix scaling issues with panels: automate scaling process
• Set up plot_diagrams() to allow for plotting a single plot with cbt=True

There are many testing recommendations in #3 that are covered in the verification script verification.py. Simply adding a test that executes this script will improve the rigor of the testing suite, as the verification script tests randomly generated diagrams.

Considerations

1. The verification script uses Multibind to handle the creation of the rate matrices for the randomly generated diagrams, so this will have to be installed to run the script.

The convention chosen for calculating pi differences and thermodynamic forces is CCW - CW, but the way they are explained in the documentation is not quite clear about why this is. This needs to be updated so it is clear what is going on and why.

Since there are already a host of examples worked out with solutions it would be prudent to throw these all into some sort of documentation where the code can be made available for others to tailor to their specific needs.

KDA is in dire need of a more user-friendly interface. Having functional code sorted into function-specific modules isn't always a bad thing, but at this point it is quite cumbersome to perform basic tasks, even for me. It would be nice if we could define a class to handle all of the function calls and data handling so a user can do something like the following, given a pre-defined rate matrix K:

import kda

KD = kda.KineticDiagram(K)

KD.calc_state_probs(store_diagrams=False, store_expressions=True)
KD.calc_net_cycle_fluxes(store_diagrams=False)

# print out results
print("State probabilities: ", KD.state_probs.values)
for cycle, net_cycle_flux in zip(KD.cycles, KD.net_cycle_fluxes):
    print(f"Net cycle flux for cycle {cycle}: {net_cycle_flux}")

# print out state probability expressions
for i, expression in enumerate(KD.state_probs.expressions):
    print(f"Expression for state {i+1}: {expression}")

Having a central object to store all of the various information associated with a kinetic diagram would substantially simplify the process to go from a simple rate matrix to a complex set of expressions or diagrams. I think this would make the user experience much more delightful, and prevent a lot of user errors since some of the functions are a bit tricky to use.

There should be basic documentation for the package, generated from the doc strings using sphinx.

The docs can also contain examples #6 .

At the moment the primary use case for generating the directional partial diagrams in NetworkX diagram for (instead of an array of edges) is to plot them using kda.plotting.draw_diagrams(), but kda.plotting.draw_diagrams() only draws many diagrams when a list of diagrams is input. This is an issue when running the following code:

dir_pars = kda.diagrams.generate_directional_partial_diagrams(G, return_edges=False)
kda.plotting.draw_diagrams(dir_pars)

because it raises the following error:

---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
<ipython-input-18-822cdbca70d2> in <module>
      1 # plotting.draw_diagrams(dir_pars, pos=pos, panel=True, panel_scale=1.5, cbt=True, rows=4, cols=8, font_size=font_size, path=PATH, label=dir_pars_label)
----> 2 plotting.draw_diagrams(dir_pars, pos=pos, panel=True, panel_scale=1.5, cbt=True, rows=4, cols=8, font_size=font_size)

c:\users\nikol\onedrive\phy495\kda\kda\plotting.py in draw_diagrams(diagrams, pos, panel, panel_scale, font_size, cbt, rows, cols, path, label)
    144         fig = plt.figure(figsize=(4, 4), tight_layout=True)
    145         ax = fig.add_subplot(111)
--> 146         node_list = list(G.nodes)
    147         node_size = 500
    148         nx.draw_networkx_nodes(

AttributeError: 'numpy.ndarray' object has no attribute 'nodes'

This is confusing because you expect to be able to feed the diagrams directly into the plotting function. There are basically 2 solutions for this:

1. change kda.diagrams.generate_directional_partial_diagrams(return_edges=False) to return a list of diagrams instead of an array. This should be done in a consistent way in comparison to kda.diagrams.generate_partial_diagrams() so they can both be used effectively.
2. change the way that kda.plotting.draw_diagrams() detects multiple diagrams (let it check for an array, for example).

When the pseudocode algorithms are finished for the partial, directional, and flux diagram generation functions it would be good to add references to these in the corresponding docstrings.

Create an initial release.

Put it on PyPi and let ReadTheDocs generate your initial docs #10

I strongly recommend using semantic versioning. Use a CHANGES or CHANGELOG file.

Pretty much exactly what the title says. It should be clear exactly what quantity is being calculated. Correcting kda.calc_cycle_flux() to kda.calc_net_cycle_flux() will make things more clear for future users.

Recent CI run failed on Python 3.9.

It appears there is an issue with test_kda.py::test_get_ordered_cycle_all_node_cycles (output pasted below) where a cycle is now considered ambiguous based on the original diagram.

FAILED kda/tests/test_kda.py::test_get_ordered_cycle_all_node_cycles - kda.exceptions.CycleError: Distinct ordered cycle could not be determined. Input diagram has multiple unique cycles that contain all nodes in the input cycle ([0, 4, 3, 2, 1]). To fix ambiguity, please input a cycle with the nodes in the correct orientation. Select one of the following possibilities: [[0, 1, 2, 3, 4], [0, 1, 2, 4, 3], [0, 1, 3, 2, 4], [0, 1, 3, 4, 2], [0, 1, 4, 2, 3], [0, 1, 4, 3, 2], [0, 2, 1, 3, 4], [0, 2, 1, 4, 3], [0, 2, 3, 1, 4], [0, 2, 4, 1, 3], [0, 3, 1, 2, 4], [0, 3, 2, 1, 4]]

I'm guessing this is due to a change in the NetworkX library. This should be fixed before addressing #70.

For the paper we should have

• 0.1 release at minimum (EDIT: I mean an official release, not just a 0.1.0 pre release)
• pypi package

In kda/tests/test_kda.py we have many tests that iterate over every combination of pairs of values, which take the majority of the test time. It might be worth using hypothesis instead to have more robust testing and potentially improve the test suite run time, since it currently takes ~150 s for tests that are not particularly rigorous.

When working on PR #25 I wrote a variant of generate_flux_diagrams() that used a while-loop to count how many cycles were in each diagram, and when it detected more than 1, it broke the loop. This seemed to be faster, but I didn't have a chance to test it thoroughly so I decided to go with a more straight-forward approach.

The code:

# count how many cycles are in the flux diagram
n_cycles = 0
contains_only_1_cycle = True
for cycle in nx.simple_cycles(diag):
    # don't count loops between pairs of nodes
    if len(cycle) > 2:
        n_cycles += 1
        # `nx.simple_cycles()` returns the forward and reverse cycles
        # so for a single-cycle diagram, there will be 2 cycles returned
        if n_cycles > 2:
            contains_only_1_cycle = False
            break
# if there is exactly 1 unique cycle in the
# generated diagram, it is valid
if contains_only_1_cycle:
    # add valid diagram to flux diagram list
    flux_diagrams.append(diag)

This works well enough, but the performance improvement depends on the order of the cycles returned by nx.simple_cycles(). So for some cases it could be substantially faster, but others could be the same. I would need to see solid evidence that it improves performance before switching to this, so for now I'm just going to dump it here.

Currently there are no instructions for how to install KDA.

Should just be

pip install kda

but we should probably ensure that everything works before posting.

Problem:
The issue stems from an oversight with using kda.get_ordered_cycle(). This function returns the NetworkX version of any list of cycle nodes so that a user can input a list of any arbitrary nodes that make up a cycle without having to conform to the strict order required for calculations. The oversight is that the NetworkX version is not always CCW or CW, but a mix of these two. So this directly affects the outputs of kda.calculate_pi_difference() and kda.calculate_thermo_force() because they have a functional form that requires the difference of the forward and reverse rate-products and the ratio of the forward and reverse rate-products, respectively.

Solution(s):

• A poor, but easy solution is to just have the user input the cycle in the correct order. This would immediately solve the problem but wouldn't allow a user to loop through all cycles as easily since they would have to make a list of all cycles manually
• Another solution would be to ask the user to input a list of two nodes that inform the direction of the CCW direction in that cycle. For a cycle [0, 1, 2, 3], the user would input [0, 1] or [1, 0], for example.
• The most difficult but also best solution would be to write a function that can determine which direction in a given cycle is the CCW direction. Then a user could input a random cycle, have it filtered by NetworkX, then oriented and output in the CCW format. To fix the aforementioned functions it would only require a single function call before the nodes are looped through.

It would be nice if the codecov behavior was more similar to some of the larger OSS projects, such as scipy.

Specific changes to make:

• Add a patch diff
• Turn the codecov PR bot comments off
• Turn annotations off (see here)
• Merge the coverage reports
• Double check that pytest is generating a report in the appropriate format so codecov can pick it up
• Double check that codecov isn't ever going to fail the CI by setting the minimum values to 1/0 for the project/patch coverage values.

The following lines in diagrams.enumerate_partial_diagrams() should be changed such that a specific error is raised:

    # get the cofactor matrix
    K_cof = _get_cofactor_matrix(K_laplace=K_laplace)
    # check that all values in the cofactor matrix are the same by
    # checking if minimum and maximum values are equal
    assert K_cof.min() == K_cof.max()

For this case it would be far better to raise an error with an error message and leave assertions for tests.

It was found that for large diagrams with many edges, the number of partial diagrams is very large (1e6 or more). With this number of diagrams the memory footprint of this code is large as well, and can easily consume Gb's or tens of Gb's of RAM. While the number of partial diagrams seems to scale as N!, there are some ways that might help reduce the memory footprint:

• Change functions to work in terms of lists of edge tuples instead of diagrams
• Save node labels as integers
• Find ways to reduce number of total graphs loaded in memory at a time
• Store values as float32 instead of float64

Note: no partial diagrams, directional partial diagrams, or flux diagrams contain edge data. They only contain node "istarget" data (True or False).

In T.L. Hill's book, "Free Energy Transduction and Biochemical Cycle Kinetics", he discusses simple relationships between transition fluxes, net cycle fluxes, and operational fluxes. Specifically, equations 3.3 and 3.4 show flux identities for a simple model. A test should be added that verifies these identities to further ensure consistency between various flux calculation types.

Test should be written to compare the output expressions, which can be compared with the normalization factor Sigma.

Improvements:

• Need to add check into flux diagram generation that makes sure each generated flux diagram has only 1 valid cycle in it
• Need to add more stringent sorting of generated flux diagrams: each relevant flux diagram must contain the entire cycle of interest, not just the nodes for that cycle (add edges check)

These issues have caused problems when evaluating diagrams with many more edges than nodes (nodes have multiple pathways away from them).

This issue was originally opened in GitLab and used to track testing ideas. I've consolidated them into a single list below, but future testing ideas will have their own issues.

Testing ideas:

• Add tests that sum probabilities and check if almost equal to 1
• Check if all generated flux diagrams contain only 1 cycle -- redundant
• For existing tests, add matrix method verification via kda.solve_matrix()
• Add tests that count numbers of partial and directional partial diagrams
• Add test for kda.calculate_pi_difference() to verify the fix for issue #15
• Add regression test that counts the number of flux diagrams generated for multiple models. This is already halfway implemented, but needs to be expanded to include more known cases
• Add additional unit tests for calculations.calc_thermo_force() (similar to test_thermo_force_4WL) where each test verifies the output of multiple models (4-state, 5-state).
• Add unit tests for calculations.calculate_pi_difference() where each test verifies the output of multiple models (4-state, 4-state with leakage, 5-state)
• Add assertion to test_flux_diagram_counts verifying each flux diagram contains only 1 unique cycle
• Add assertion to test_diagram_counts verifying the list of edges for each partial and directional diagram contains all nodes in the input diagram

Originally there were 3 types of diagrams: partial, directional partial, and flux diagrams. But the "partial" in "directional partial" was deemed unnecessary, and has since been removed from the function names (gh-41). It hasn't made its way out of all of the documentation though, and this should change to avoid confusion.

A simple "directional partial" -> "directional" change for each instance should be simple enough.

At the moment it is difficult to calculate operational fluxes because each cycle has to have its direction chosen "by hand". For example, the kda.calculations functions calc_thermo_force, calc_pi_difference, calc_net_cycle_flux, etc. require the user to input the "order" of the cycle, which is a deliberately ordered pair of nodes which is compared to the input cycle to determine which direction is the "forward" cycle direction. While this method allows for complete control over the cycle direction, it is very cumbersome and requires the user to manually record node pairs.

A better solution might be to create a function that operates on the input diagram G which goes through each cycle and labels them as "CCW" or "CW" (or some other naming convention) based on the node positions and a user-defined direction. That way a user would only be required to specify the "positive" direction once, and if they attempted to use a function that requires these to be labeled an error could be raised. The inputs to this function would likely be a dictionary containing the information of any ions/substrates that are being transported. For example, dict("H+"=(0, 1)) would mean a proton is transported in the positive direction through the edge (0, 1).

With this information, something like the following could be performed:

1. store all unique dict keys (substrate names) in the diagram object
2. iterate over all cycles
3. check if each cycle contains any of the edges (or their reverse):
   • if the cycle contains only 1 edge (either forward or reverse), label it according to the prescribed direction
   • if the cycle contains no forward/reverse edges, label None
   • if the cycle contains pairs of forward/reverse edges for the same species, label None
   • if the cycle contains a combination of edges from multiple species, add unique labels for each species (based on direction)

We might have to temporarily store all of the "positive" and "reverse" edge tuples for each cycle but hopefully that can be avoided.

This is probably a naïve solution, but storing this data along with the cycle data (maybe even in a new object) will help prevent user errors using the current method, and allow for easier cycle verification before calculations are performed. It might even warrant a special plotting function that draws each cycle and (where applicable) a CW/CCW arrow for each substrate it transports. That would make it very simple to double check things by eye.

• Break down overall cycle into several individual cycles, calculate their individual fluxes, recombine
• Construct analytic cycle flux function for both cycles (B and C)
• Input parameters from paper to verify paper results

As mentioned in 2.) of this comment, draw_cycles() can likely be further simplified by appropriate use of _plot_panel(). It would be nice to remove the redundant code there so there is only a single location to edit if further changes are needed.

As discussed in 3. of this comment, the structures of plotting.draw_diagrams() and plotting.draw_cycles() are very similar. It shouldn't be terribly difficult to combine these into a single plotting function that can handle both kinds of figures. This is probably best handled after the completion of #47.

As discussed in #33 (comment), the plotting.py module needs to be reorganized and cleaned up more.

TODO:

• Move the plotting code into functions targeted at specific plotting cases (single plot, multiple plots, panel, etc.)
• Update the defaults to use panel plots since these are generally the easiest to read
• Move remaining repeated code (like the cbt=True code) into functions
• Improve documentation overall, and add examples
• Add tests to test_plotting.py that more rigorously test the smaller functions

Problem

The primary limitation for graph complexity is the space complexity of functions like kda.diagrams.generate_partial_diagrams() and kda.diagrams.generate_directional_partial_diagrams(). This is inherently related to issue #2, where this issue is focused on one part of that issue.

This problem is 2 fold:

1. the algorithm for generating partial diagrams is inefficient
2. the functions store more information than they need to

While 2. does contribute to the issue, the main issue is 1., especially for simpler input graphs. At the moment, the partial diagrams are generated by finding every possible combination of N-1 edges then discarding subgraphs that don't fit the definition of a partial diagram. Complete input graphs contain N(N-1)/2 unique edges, and thus have N^(N-2) partial diagrams (after doing the appropriate combinatorics). These yield approximately N! possibilities (using Stirling's approximation), so generating all possibilities is probably appropriate. But for simpler input graphs the number of partial diagrams is much less numerous so generating every possibility is likely not necessary.

Solution(s)

Before I list the possible solutions, I have to mention that in graph theory a "partial diagram" is known as a "spanning tree", so the solutions below will be using the graph theory language.

There are many solutions to this problem, many of which have formulas for calculating the space complexity of the algorithm. At the moment it isn't obvious which method is the best method, so I'm just going to dump links below.

• This is the first useful article I found that gives a clear way to quantify the number of partial diagrams for a given graph using an adjacency matrix. They also discuss a host of different methods (9 I think) for generating the edge lists for all possible spanning trees.

I was going to list more, but the other articles/papers I have are all references in the above paper so for now I will leave it here.

The public plotting functions have parameters label and path that allow a user to specify a unique file name (without extension) and a save path to store their plot/plots. The way this is implemented is a bit awkward. It's probably sufficient to simply have a single parameter filename which, if specified, stores the plot/plots in the current working directory. It may not be the worst idea to give the option of specifying a separate path, but the current implementation requires both parameters, which is not ideal.

It would also be nice to save figures with different file types (e.g. .png, .pdf, .svg) based on user preference. For example, if the filename kinetic_diagram.svg is specified it will save it as this type, but if only kinetic_diagram is specified without an extension it defaults to .pdf. This will require some string parsing since not all file types would be supported, but it would be nice to have the option to choose the file type without too much hassle.

We are currently testing Python 3.7-3.9 but we should add Python 3.10 to the mix.