MarkovModelModule
StateTransitionGraph: A directed graph class for modeling continuous-time Markov chains (CTMCs).
This module defines the StateTransitionGraph
class, an extension of networkx.DiGraph
with specialized methods and attributes for analyzing, simulating, and visualizing
continuous-time Markov chains (CTMCs).
Key Features
- Symbolic and numerical parameter support using SymPy.
- Automatic generation and storage of the rate matrix (Q-matrix).
- Simulation of state transitions and computation of steady-state probabilities.
- Enhanced graph visualization with transition rates and labels.
- Support for dynamic import of constants via a custom module (default: "Constants").
Typical Use Case
- Define a set of model parameters (e.g., lambda, mu, n).
- Create a
StateTransitionGraph
instance with those parameters. - Add states and transitions with symbolic rates.
- Compute the rate matrix and steady-state probabilities.
- Run simulations and visualize the CTMC behavior.
Example
>>> parameters = {'lambda': lam, 'mu': 1.0, 'n': 5}
>>> G = StateTransitionGraph(parameters)
>>> # Define transitions
>>> for i in range(parameters['n']):
>>> G.addTransition(i, i + 1, G.sym("lambda"), tid=const.ARRIVAL)
>>> G.addTransition(i + 1, i, (i + 1) * G.sym("mu"), tid=const.DEPARTURE)
>>> # Set node properties
>>> G.setAllStateDefaultProperties()
>>> G.setStateColor(G.n, const.COLOR_NODE_BLOCKING)
>>> # Compute steady-state probabilities
>>> probs = G.calcSteadyStateProbs()
>>> # Plot bar chart of state probabilities
>>> import matplotlib.pyplot as plt
>>> plt.figure(1, clear=True)
>>> plt.bar(probs.keys(), probs.values())
>>> plt.xlabel('states')
>>> plt.ylabel('probability')
>>> # Layout for transition graph
>>> pos = {s: (s, 0) for s in G.states()}
>>> G.drawTransitionGraph(pos=pos, bended=True, label_pos=0.5, num=2)
>>> # Simulate and animate system dynamics
>>> states, times = G.simulateMarkovModel(startNode=0, numEvents=2000)
See Also
networkx.DiGraph : The base class that is extended. sympy.symbols : Used to define symbolic transition rates. matplotlib.pyplot : Used internally for graph visualization.
Author
Tobias Hossfeld tobias.hossfeld@uni-wuerzburg.de
Created
5 April 2025
1# -*- coding: utf-8 -*- 2""" 3StateTransitionGraph: A directed graph class for modeling continuous-time Markov chains (CTMCs). 4 5This module defines the `StateTransitionGraph` class, an extension of `networkx.DiGraph` 6with specialized methods and attributes for analyzing, simulating, and visualizing 7continuous-time Markov chains (CTMCs). 8 9Key Features 10------------ 11- Symbolic and numerical parameter support using SymPy. 12- Automatic generation and storage of the rate matrix (Q-matrix). 13- Simulation of state transitions and computation of steady-state probabilities. 14- Enhanced graph visualization with transition rates and labels. 15- Support for dynamic import of constants via a custom module (default: "Constants"). 16 17Typical Use Case 18---------------- 191. Define a set of model parameters (e.g., lambda, mu, n). 202. Create a `StateTransitionGraph` instance with those parameters. 213. Add states and transitions with symbolic rates. 224. Compute the rate matrix and steady-state probabilities. 235. Run simulations and visualize the CTMC behavior. 24 25 26Example 27------- 28>>> parameters = {'lambda': lam, 'mu': 1.0, 'n': 5} 29>>> G = StateTransitionGraph(parameters) 30 31>>> # Define transitions 32>>> for i in range(parameters['n']): 33>>> G.addTransition(i, i + 1, G.sym("lambda"), tid=const.ARRIVAL) 34>>> G.addTransition(i + 1, i, (i + 1) * G.sym("mu"), tid=const.DEPARTURE) 35 36>>> # Set node properties 37>>> G.setAllStateDefaultProperties() 38>>> G.setStateColor(G.n, const.COLOR_NODE_BLOCKING) 39 40>>> # Compute steady-state probabilities 41>>> probs = G.calcSteadyStateProbs() 42 43>>> # Plot bar chart of state probabilities 44>>> import matplotlib.pyplot as plt 45>>> plt.figure(1, clear=True) 46>>> plt.bar(probs.keys(), probs.values()) 47>>> plt.xlabel('states') 48>>> plt.ylabel('probability') 49 50>>> # Layout for transition graph 51>>> pos = {s: (s, 0) for s in G.states()} 52>>> G.drawTransitionGraph(pos=pos, bended=True, label_pos=0.5, num=2) 53 54>>> # Simulate and animate system dynamics 55>>> states, times = G.simulateMarkovModel(startNode=0, numEvents=2000) 56 57See Also 58-------- 59networkx.DiGraph : The base class that is extended. 60sympy.symbols : Used to define symbolic transition rates. 61matplotlib.pyplot : Used internally for graph visualization. 62 63 64Author 65------ 66Tobias Hossfeld <tobias.hossfeld@uni-wuerzburg.de> 67 68Created 69------- 705 April 2025 71 72""" 73 74 75import networkx as nx # For the magic 76import matplotlib.pyplot as plt # For plotting 77import numpy as np 78from scipy import linalg 79#import itertools 80import sympy as sp 81import warnings 82import random 83import importlib 84from scipy.linalg import expm 85from IPython.display import display, clear_output 86from openpyxl import load_workbook 87from openpyxl.styles import PatternFill 88import matplotlib.colors as mcolors 89import pandas as pd 90import os 91#import Constants as const 92#from collections import namedtuple 93#%% 94class StateTransitionGraph(nx.DiGraph): 95 """ 96 Initialize the StateTransitionGraph for continuous time Markov chains (CTMC). 97 Extends the networkx.Graph class to include additional methods and properties 98 for the analysis, visualization, and simulation of CTMCs. 99 The analysis allows deriving the state probabilities in the steady-state 100 and in the transient phase. 101 102 This constructor sets up symbolic and numeric parameters for the state transition graph, 103 imports a constants module dynamically, and initializes internal data structures for 104 rate matrix computation and steady-state probability analysis. 105 106 Parameters 107 ---------- 108 - **parameters** : dict 109 A dictionary mapping parameter names to their numerical values. 110 Example: {'lambda': 1.0, 'mu': 1.0}. 111 For each parameter key: 112 - A symbolic variable is created using sympy.symbols. Can be accessed as `G.sym("key")` or `G.sym_key` 113 - A variable is created which can be accessed via `G.key` (and the numerical value is assigned) 114 - For the key="lambda", `G.lam` and `G.sym_lambda` are created 115 116 The parameters can be accessed via property `parameters` 117 - **constants** : str, optional 118 The name of a Python module (as string) to import dynamically and assign to `self.const`. 119 This can be used to store global or user-defined constants. 120 Default is "Constants". 121 - `*args` : tuple 122 Additional positional arguments passed to the base `networkx.DiGraph` constructor. 123 - `**kwargs` : dict 124 Additional keyword arguments passed to the base `networkx.DiGraph` constructor. 125 126 Attributes 127 ---------- 128 - `_rate_matrix` : ndarray or None 129 The continuous-time rate matrix `Q` of the CTMC. Computed when calling `createRateMatrix()`. 130 - `_state_probabilities` : dict or None 131 Dictionary of steady-state probabilities, computed via `calcSteadyStateProbs()`. 132 - `_n2i` : dict or None 133 Mapping from node identifiers to matrix row/column indices used for rate matrix construction. 134 - `_sym` : dict 135 Mapping of parameter names (strings) to their symbolic `sympy` representations. 136 - `_subs` : dict 137 Mapping from symbolic parameters to their numerical values (used for substitutions). 138 - `_parameters` : dict 139 Copy of the original parameter dictionary passed at initialization. 140 - `_transitionProbabilitesComputed` : bool 141 Internal flag indicating whether transition probabilities for the embedded Markov chain 142 have been computed. 143 - `const` : module 144 Dynamically imported module (default `"Constants"`), used for storing user-defined constants 145 such as colors or transition IDs. 146 147 Inherits 148 -------- 149 networkx.Graph: The base graph class from NetworkX. 150 """ 151 152 def __init__(self, parameters, constants="Constants", *args, **kwargs): 153 """ 154 Initialize the state transition graph with symbolic parameters. 155 156 This constructor sets up the symbolic and numerical environment for modeling 157 a CTMC. Parameters are stored as both symbolic expressions and numerical values, 158 and convenient attribute access is provided for both forms. 159 160 Parameters 161 ---------- 162 - **parameters** : dict 163 A dictionary mapping parameter names to their numerical values. 164 Example: {'lambda': 1.0, 'mu': 1.0}. 165 For each parameter key: 166 - A symbolic variable is created using sympy.symbols. Can be accessed as `G.sym("key")` or `G.sym_key` 167 - A variable is created which can be accessed via `G.key` (and the numerical value is assigned) 168 - For the key="lambda", `G.lam` and `G.sym_lambda` are created 169 170 The parameters can be accessed via property `parameters` 171 - **constants** : str, optional 172 The name of a Python module (as string) to import dynamically and assign to `self.const`. 173 This can be used to store global or user-defined constants. 174 Default is "Constants". 175 - `*args` : tuple 176 Additional positional arguments passed to the base `networkx.DiGraph` constructor. 177 - `**kwargs` : dict 178 Additional keyword arguments passed to the base `networkx.DiGraph` constructor. 179 """ 180 super().__init__(*args, **kwargs) 181 self._rate_matrix = None # Rate Matrix Q of this transition graph 182 self._state_probabilities = None # steady state proabilities of this CTMC 183 self._n2i = None # internal index for generating the rate matrix Q: dictionary (node to index) 184 185 #parameters = {'lambda': 1.0, 'mu':1.0, 'n':5} 186 self._parameters = parameters 187 self._sym = {key: sp.symbols(key) for key in parameters} 188 self._subs = {} # contains symbolic parameter as key and provides numerical value for computation. 189 for k in parameters: 190 self._subs[ self._sym[k] ] = parameters[k] 191 192 # access to parameters and symbols using dot-notation, e.g. G.lam, G.sym_lambda, or G.sym("lambda") 193 for key, value in parameters.items(): 194 if key=='lambda': 195 setattr(self, 'lam', value) 196 else: 197 setattr(self, key, value) 198 199 setattr(self, "sym_"+key, self._sym[key]) 200 201 self._transitionProbabilitesComputed = False 202 203 self._const = importlib.import_module(constants) 204 205 206 207 def addState(self, state, color=None, label=None): 208 """ 209 Add a state (node) to the graph with optional color and label attributes. 210 211 If the node already exists, its attributes are updated with the provided values. 212 213 Parameters 214 ---------- 215 - **state** : hashable 216 The identifier for the state. Can be a string, int, or a tuple. If no label 217 is provided, the label will be auto-generated as: 218 - comma-separated string for tuples (e.g., `(1, 2)` -> `'1,2'`) 219 - stringified value for other types 220 221 - **color** : str, optional 222 The color assigned to the state (used for visualization). 223 If not specified, defaults to `self.const.COLOR_NODE_DEFAULT`. 224 225 - **label** : str, optional 226 A label for the node (used for display or annotation). If not provided, 227 it will be auto-generated from the state identifier. 228 229 Notes 230 ----- 231 - This method wraps `self.add_node()` from NetworkX. 232 - If the node already exists in the graph, this method does not raise an error; 233 it will update the node’s existing attributes. 234 """ 235 if label is None: 236 if isinstance(state, tuple): 237 label = ','.join(map(str, state)) # if state is a tuple 238 else: 239 label = str(state) 240 #self.nodes[state]["label"] = label 241 242 if color is None: 243 color = self._const.COLOR_NODE_DEFAULT 244 245 # no error if the node already exists.Attributes get updated if you pass them again. 246 self.add_node(state, color=color, label=label) 247 248 def setStateProperties(self, state, color=None, label=None, **kwargs): 249 """ 250 Update visual or metadata properties for an existing state (node) in the graph. 251 252 This method updates the node's attributes such as color, label, and any additional 253 custom key-value pairs. The node must already exist in the graph. 254 255 Parameters 256 ---------- 257 - **state** : hashable 258 Identifier of the state to update. Must be a node already present in the graph. 259 260 - **color** : str, optional 261 Color assigned to the node, typically used for visualization. If not provided, 262 the default node color (`self._const.COLOR_NODE_DEFAULT`) is used. 263 264 - **label** : str, optional 265 Label for the node. If not provided, it is auto-generated from the state ID 266 (especially for tuples). 267 268 - `**kwargs` : dict 269 Arbitrary additional key-value pairs to set as node attributes. These can include 270 any custom metadata. 271 272 Raises 273 ------ 274 KeyError 275 If the specified state is not present in the graph. 276 277 Notes 278 ----- 279 - This method ensures that required attributes (like `color` and `label`) are set for visualization. 280 - Existing node attributes are updated or extended with `kwargs`. 281 """ 282 if state not in self: 283 raise KeyError(f"State '{state}' does not exist in the graph.") 284 285 # mandatory for plotting 286 self.addState(state, color, label) 287 288 # # Add any additional attributes 289 self.nodes[state].update(kwargs) 290 291 def setAllStateDefaultProperties(self, **kwargs): 292 """ 293 Set default visual and metadata properties for all states (nodes) in the graph. 294 295 This method applies `setStateProperties()` to every node, ensuring each state 296 has at least the required default attributes (e.g., color and label). 297 Any additional keyword arguments are passed to `setStateProperties()`. 298 299 Parameters 300 ---------- 301 `**kwargs` : dict, optional 302 Optional keyword arguments to apply uniformly to all nodes, 303 such as color='gray', group='core', etc. 304 305 Notes 306 ----- 307 - This method is useful for initializing node attributes before plotting 308 or performing other graph-based computations. 309 - Nodes without existing attributes will receive defaults via `addState()`. 310 """ 311 for state in self.nodes(): 312 self.setStateProperties(state, **kwargs) 313 314 315 #def addTransition(self, origin_state, list_destination_state, sym_rate, color=plt.cm.tab10(5), tid = None): 316 def addTransition(self, origin_state, destination_state, sym_rate, tid = None, color=None, label=None): 317 """ 318 Add a directed transition (edge) between two states with a symbolic transition rate. 319 320 This method adds a directed edge from `origin_state` to `destination_state`, 321 stores both the symbolic and numerical rate, and sets optional display properties 322 such as color, label, and transition ID. If the edge already exists, a warning is issued. 323 324 Parameters 325 ---------- 326 - **origin_state** : hashable 327 The source node (state) of the transition. 328 329 - **destination_state** : hashable 330 The target node (state) of the transition. 331 332 - **sym_rate** : sympy.Expr 333 A symbolic expression representing the transition rate. 334 This will be numerically evaluated using internal substitutions (`self._subs`). 335 336 - **tid** : str or int, optional 337 A transition identifier (e.g. for grouping or styling). 338 If provided, the color will default to `self._const.COLOR_TRANSITIONS[tid]` if not set explicitly. 339 340 - **color** : str, optional 341 Color for the edge (used in visualization). Defaults to the transition-specific color 342 if `tid` is provided. 343 344 - **label** : str, optional 345 LaTeX-formatted string for labeling the edge. If not provided, a default label is 346 generated using `sympy.latex(sym_rate)`. 347 348 Notes 349 ----- 350 - The numeric rate is computed by substituting known parameter values into `sym_rate`. 351 - The following attributes are stored on the edge: 352 - `rate` (float): numeric value of the transition rate 353 - `sym_rate` (sympy.Expr): original symbolic expression 354 - `label` (str): display label 355 - `color` (str): edge color 356 - `tid` (str or int): transition ID 357 - If the edge already exists, a `UserWarning` is issued but the edge is still overwritten. 358 """ 359 if tid is not None: 360 if color is None: 361 color = self._const.COLOR_TRANSITIONS[ tid ] 362 363 rate = sym_rate.subs( self._subs ) 364 365 if label is None: string = "$"+sp.latex(sym_rate)+"$" 366 367 if self.has_edge(origin_state, destination_state): 368 theedge = self[origin_state][destination_state] 369 warningstr = f"Edge already exists: {origin_state} -> {destination_state} with rate {theedge['label']}" 370 warnings.warn(warningstr, category=UserWarning) 371 self.add_edge(origin_state, destination_state, rate=float(rate), label=string, 372 color=color, sym_rate=sym_rate, tid=tid) 373 374 def sym(self, key): 375 """ 376 Retrieve the symbolic representation of a model parameter. 377 378 Parameters 379 ---------- 380 **key** : str 381 The name of the parameter (e.g., 'lambda', 'mu', 'n'). 382 383 Returns 384 ------- 385 **sympy.Symbol** 386 The symbolic variable corresponding to the given key. 387 388 Raises 389 ------ 390 KeyError 391 If the key does not exist in the internal symbolic mapping (`_sym`). 392 393 Examples 394 -------- 395 >>> G.sym('lambda') 396 lambda 397 398 Notes 399 ----- 400 - The symbolic mapping is initialized during object construction based on the `parameters` dictionary. 401 - To access the symbol using dot-notation, use `.sym_key` attributes created during initialization. 402 """ 403 return self._sym[key] 404 405 406 @property 407 def parameters(self): 408 """ 409 Return the model parameters as a dictionary. 410 411 This includes the original numerical values provided during initialization. 412 413 Returns 414 ------- 415 **dict** 416 Dictionary of parameter names and their corresponding numeric values. 417 418 Examples 419 -------- 420 >>> G.parameters 421 {'lambda': 1.0, 'mu': 0.5, 'n': 3} 422 423 Notes 424 ----- 425 - This property does not return symbolic values. Use `self.sym(key)` for symbolic access. 426 - Parameter symbols can also be accessed using attributes like `sym_lambda` or the `sym()` method. 427 """ 428 return self._parameters 429 430 431 432 @property 433 def rate_matrix(self): 434 """ 435 Return the continuous-time rate matrix (Q matrix) of the state transition graph. 436 437 This matrix defines the transition rates between all pairs of states in the 438 continuous-time Markov chain (CTMC) represented by this graph. 439 440 Returns 441 ------- 442 **numpy.ndarray** or None 443 The rate matrix `Q`, where `Q[i, j]` is the transition rate from state `i` to `j`. 444 Returns `None` if the rate matrix has not been computed yet. 445 446 Notes 447 ----- 448 - The matrix is typically generated via a method like `computeRateMatrix()` or similar. 449 - Internal state indexing is handled by `_n2i` (node-to-index mapping). 450 - The rate matrix is stored internally as `_rate_matrix`. 451 """ 452 return self._rate_matrix 453 454 455 456 @property 457 def state_probabilities(self): 458 """ 459 Return the steady-state probabilities of the continuous-time Markov chain (CTMC). 460 461 These probabilities represent the long-run proportion of time the system spends 462 in each state, assuming the CTMC is irreducible and positive recurrent. 463 464 Returns 465 ------- 466 **dict** or None 467 A dictionary mapping each state to its steady-state probability. 468 Returns `None` if the probabilities have not yet been computed. 469 470 Notes 471 ----- 472 - The probabilities are stored internally after being computed via a dedicated method 473 `calcSteadyStateProbs()`. 474 - The result is stored in the `_state_probabilities` attribute. 475 """ 476 return self._state_probabilities 477 478 479 480 def printEdges(self, state): 481 """ 482 Print all outgoing edges from a given state in the graph. 483 484 For each edge from the specified `state`, this method prints the destination 485 state and the corresponding transition label (typically the symbolic rate). 486 487 Parameters 488 ---------- 489 **state** : hashable 490 The source node (state) from which outgoing edges will be printed. 491 492 Returns 493 ------- 494 None 495 This method prints output directly to the console and does not return a value. 496 497 Notes 498 ----- 499 - The transition label is retrieved from the edge attribute `'label'`. 500 - Assumes `state` exists in the graph; otherwise, a KeyError will be raised. 501 """ 502 print(f"from {state}") 503 node = self[state] 504 for e in node: 505 rate = node[e]['label'] 506 print(f" -> {e} with {rate}") 507 508 509 def createRateMatrix(self): 510 """ 511 Construct the continuous-time Markov chain (CTMC) rate matrix Q from the graph. 512 513 This method generates the generator matrix `Q` corresponding to the structure and 514 rates defined in the directed graph. Each node represents a system state, and each 515 directed edge must have a `'rate'` attribute representing the transition rate from 516 one state to another. 517 518 The matrix `Q` is constructed such that: 519 - `Q[i, j]` is the transition rate from state `i` to state `j` 520 - Diagonal entries satisfy `Q[i, i] = -sum(Q[i, j]) for j ≠ i`, so that each row sums to 0 521 522 The resulting matrix and the internal node-to-index mapping (`_n2i`) are stored 523 as attributes of the graph. 524 525 Returns 526 ------- 527 None 528 The method updates internal attributes: 529 - `_rate_matrix` : numpy.ndarray 530 The CTMC rate matrix Q 531 - `_n2i` : dict 532 Mapping from node identifiers to matrix row/column indices 533 534 Notes 535 ----- 536 - Assumes all edges in the graph have a `'rate'` attribute. 537 - The node order used in the matrix is consistent with the order from `list(self.nodes)`. 538 - Use `self.rate_matrix` to access the generated matrix after this method is called. 539 """ 540 n2i = {} 541 nodes = list(self.nodes) 542 for i,node in enumerate(nodes): 543 n2i[node] = i 544 Q = np.zeros((len(nodes),len(nodes))) # rate matrix 545 546 for edge in self.edges: 547 i0 = n2i[edge[0]] 548 i1 = n2i[edge[1]] 549 Q[i0,i1] = self[edge[0]][edge[1]]["rate"] 550 551 np.fill_diagonal(Q, -Q.sum(axis=1)) 552 self._rate_matrix = Q 553 self._n2i = n2i 554 #return Q, n2i 555 556 def calcSteadyStateProbs(self, verbose=False): 557 """ 558 Compute the steady-state probabilities of the continuous-time Markov chain (CTMC). 559 560 This method calculates the stationary distribution of the CTMC defined by the rate 561 matrix of the state transition graph. The computation solves the linear system 562 `XQ = 0` with the constraint that the probabilities X sum to 1, using a modified 563 version of the rate matrix `Q`. 564 565 Parameters 566 ---------- 567 `verbose` : bool, optional 568 If True, prints the rate matrix, the modified matrix, the right-hand side vector, 569 and the resulting steady-state probabilities. 570 571 Returns 572 ------- 573 **probs** : dict 574 A dictionary mapping each state (node in the graph) to its steady-state probability. 575 576 Notes 577 ----- 578 - The graph must already contain valid transition rates on its edges. 579 - If the rate matrix has not yet been created, `createRateMatrix()` is called automatically. 580 - The method modifies the last column of the rate matrix to impose the normalization constraint 581 `sum(X) = 1`. 582 - The inverse of the modified matrix is computed directly, so the graph should be small enough 583 to avoid numerical instability or performance issues. 584 - The resulting probabilities are stored internally in `_state_probabilities` and can be accessed 585 via the `state_probabilities` property. 586 """ 587 588 # compute transition matrix if not already done 589 if self._rate_matrix is None: self.createRateMatrix() 590 Q2 = self._rate_matrix.copy() 591 if verbose: 592 print("Rate matrix Q") 593 print(f'Q=\n{Q2}\n') 594 595 Q2[:, -1] = 1 596 if verbose: 597 print(f'Matrix is changed to\nQ2=\n{Q2}\n') 598 599 b = np.zeros(len(Q2)) 600 b[-1] = 1 601 if verbose: 602 print("Solve X = b @ Q2") 603 print(f'b=\n{b}\n') 604 605 # state probabilities 606 X = b @ linalg.inv(Q2) # compute the matrix inverse 607 608 # Generate a matrix with P(X,Y) 609 matrix = {} 610 for n in self.nodes: 611 matrix[n] = X[ self._n2i[n] ] 612 613 if verbose: 614 print("Steady-state probabilities X") 615 print(f'X=\n{matrix}\n') 616 self._state_probabilities = matrix 617 return matrix 618 619 def calcTransitionProbabilities(self): 620 """ 621 Compute and store transition probabilities for each node in the graph of the embedded discrete Markov chain. 622 623 This method processes each node in the graph, treating outgoing edges as 624 transitions in an embedded discrete-time Markov chain derived from a continuous-time 625 Markov chain (CTMC). Transition *rates* on outgoing edges are normalized into 626 *probabilities*, and relevant data is stored in the node attributes. 627 628 For each node, the following attributes are set: 629 - `"transitionProbs"` : np.ndarray 630 Array of transition probabilities to successor nodes. 631 - `"transitionNodes"` : list 632 List of successor node identifiers (ordered as in the probability array). 633 - `"sumOutgoingRates"` : float 634 Sum of all outgoing transition rates. Used to model the exponential waiting time 635 in continuous-time simulation (`T ~ Exp(sumOutgoingRates)`). 636 637 Notes 638 ----- 639 - The graph must be a directed graph where each edge has a `'rate'` attribute. 640 - Nodes with no outgoing edges are considered absorbing states. 641 - Useful for discrete-event simulation and Markov chain Monte Carlo (MCMC) modeling. 642 - Sets an internal flag `_transitionProbabilitesComputed = True` upon completion. 643 644 Raises 645 ------ 646 KeyError 647 If any outgoing edge lacks a `'rate'` attribute. 648 """ 649 for node in self: 650 successors = list(self.successors(node)) 651 rates = np.array([self[node][nbr]['rate'] for nbr in successors]) 652 probs = rates / rates.sum() # Normalize to probabilities 653 self.nodes[node]["transitionProbs"] = probs # probabilities 654 self.nodes[node]["transitionNodes"] = successors # list of nodes 655 self.nodes[node]["sumOutgoingRates"] = rates.sum() # time in state ~ EXP(sumOutgoingRates) 656 numel = len(successors) 657 if numel==0: 658 print(f"Absorbing state: {node}") 659 self._transitionProbabilitesComputed = True 660 661 def simulateMarkovModel(self, startNode=None, numEvents=10): 662 """ 663 Simulate a trajectory through the state transition graph as a continuous-time Markov chain (CTMC). 664 665 This method performs a discrete-event simulation of the CTMC by generating a random walk 666 through the graph, starting at `startNode` (or a default node if not specified). Transitions 667 are chosen based on the normalized transition probabilities computed from edge rates, and 668 dwell times are sampled from exponential distributions based on the total outgoing rate 669 from each state. 670 671 Parameters 672 ---------- 673 - **startNode** : hashable, optional 674 The node at which the simulation starts. If not specified, the first node in the graph 675 (based on insertion order) is used. 676 677 - **numEvents** : int, optional 678 The number of transition events (steps) to simulate. Default is 10. 679 680 Returns 681 ------- 682 **states** : list 683 A list of visited states, representing the sequence of nodes traversed in the simulation. 684 685 **times** : numpy.ndarray 686 An array of dwell times in each state, sampled from exponential distributions 687 with rate equal to the sum of outgoing transition rates from each state. 688 689 Notes 690 ----- 691 - This method automatically calls `calcTransitionProbabilities()` if transition probabilities 692 have not yet been computed. 693 - If an absorbing state (i.e., a node with no outgoing edges) is reached, the simulation stops early. 694 - The cumulative time can be obtained by `np.cumsum(times)`. 695 - State durations follow `Exp(sumOutgoingRates)` distributions per CTMC theory. 696 """ 697 if not self._transitionProbabilitesComputed: self.calcTransitionProbabilities() 698 699 # simulate random walk through Graph = Markov Chain 700 if startNode is None: 701 startNode = next(iter(self)) 702 703 704 states = [] # list of states 705 times = np.zeros(numEvents) # time per state 706 707 node = startNode 708 for i in range(numEvents): 709 outgoingRate = self.nodes[node]["sumOutgoingRates"] # time in state ~ EXP(sumOutgoingRates) 710 times[i] = random.expovariate(outgoingRate) # exponentially distributed time in state 711 states.append(node) 712 713 # get next node 714 probs = self.nodes[node]["transitionProbs"] # probabilities 715 successors = self.nodes[node]["transitionNodes"] # list of nodes 716 numel = len(successors) 717 if numel==0: 718 print(f"Absorbing state: {node}") 719 break 720 elif numel==1: 721 nextNode = successors[0] 722 else: 723 nextNode = successors[ np.random.choice(numel, p=probs) ] 724 # edge = G[node][nextNode] 725 node = nextNode 726 727 return states, times 728 729 def calcProbabilitiesSimulationRun(self, states, times): 730 """ 731 Estimate empirical state probabilities from a simulation run of the CTMC. 732 733 This method computes the fraction of total time spent in each state during a simulation. 734 It uses the sequence of visited states and corresponding sojourn times to estimate the 735 empirical (time-weighted) probability distribution over states. 736 737 Parameters 738 ---------- 739 - **states** : list of tuple 740 A list of state identifiers visited during the simulation. 741 - **times** : numpy.ndarray 742 An array of sojourn times corresponding to each visited state, same length as `states`. 743 744 Returns 745 ------- 746 **simProbs** : dict 747 A dictionary mapping each unique state to its estimated empirical probability, 748 based on total time spent in that state relative to the simulation duration. 749 750 Notes 751 ----- 752 - Internally, states are mapped to flat indices for efficient counting using NumPy. 753 - This method assumes that `states` contains only valid node identifiers from the graph. 754 - If `_n2i` (node-to-index mapping) has not been initialized, it will be computed here. 755 - This is especially useful for validating analytical steady-state probabilities 756 against sampled simulation results. Or in case of state-space explosions. 757 """ 758 if self._n2i is None: 759 n2i = {} 760 nodes = list(self.nodes) 761 for i,node in enumerate(nodes): 762 n2i[node] = i 763 self._n2i = n2i 764 765 766 # Map tuples to flat indices 767 data = np.array(states) 768 max_vals = data.max(axis=0) + 1 769 indices = np.ravel_multi_index(data.T, dims=max_vals) 770 771 # Weighted bincount 772 bincount = np.bincount(indices, weights=times) 773 774 # Decode back to tuples 775 tuple_indices = np.array(np.unravel_index(np.nonzero(bincount)[0], shape=max_vals)).T 776 timePerState = bincount[bincount > 0] 777 778 simProbs = {} 779 totalTime = times.sum() 780 for t, c in zip(tuple_indices, timePerState): 781 simProbs[tuple(t)] = c / totalTime 782 #print(f"Tuple {tuple(t)} → weighted count: {c}") 783 return simProbs 784 785 786 def _draw_networkx_edge_labels( 787 self, 788 pos, 789 edge_labels=None, label_pos=0.75, 790 font_size=10, font_color="k", font_family="sans-serif", font_weight="normal", 791 alpha=None, bbox=None, 792 horizontalalignment="center", verticalalignment="center", 793 ax=None, rotate=True, clip_on=True, rad=0): 794 """ 795 Draw edge labels for a directed graph with optional curved edge support. 796 797 This method extends NetworkX's default edge label drawing by supporting 798 curved (bended) edges through quadratic Bézier control points. It places 799 labels dynamically at a specified position along the curve and optionally 800 rotates them to follow the edge orientation. 801 802 Parameters 803 ---------- 804 pos : dict 805 Dictionary mapping nodes to positions (2D coordinates). Required. 806 807 edge_labels : dict, optional 808 Dictionary mapping edge tuples `(u, v)` to label strings. If None, 809 edge data is used and converted to strings automatically. 810 811 label_pos : float, optional 812 Position along the edge to place the label (0=head, 1=tail). 813 Default is 0.75 (closer to the target node). 814 815 font_size : int, optional 816 Font size of the edge labels. Default is 10. 817 818 font_color : str, optional 819 Font color for the edge labels. Default is 'k' (black). 820 821 font_family : str, optional 822 Font family for the edge labels. Default is 'sans-serif'. 823 824 font_weight : str, optional 825 Font weight for the edge labels (e.g., 'normal', 'bold'). Default is 'normal'. 826 827 alpha : float, optional 828 Opacity for the edge labels. If None, labels are fully opaque. 829 830 bbox : dict, optional 831 Matplotlib-style bbox dictionary to style the label background. 832 Default is a white rounded box. 833 834 horizontalalignment : str, optional 835 Horizontal alignment of the text. Default is 'center'. 836 837 verticalalignment : str, optional 838 Vertical alignment of the text. Default is 'center'. 839 840 ax : matplotlib.axes.Axes, optional 841 Axes on which to draw the labels. If None, uses the current axes. 842 843 rotate : bool, optional 844 Whether to rotate the label to match the edge direction. Default is True. 845 846 clip_on : bool, optional 847 Whether to clip labels at the plot boundary. Default is True. 848 849 rad : float, optional 850 Curvature of the edge (used to determine Bézier control point). 851 Positive values bend counterclockwise; negative values clockwise. Default is 0 (straight). 852 853 Returns 854 ------- 855 dict 856 Dictionary mapping edge tuples `(u, v)` to the corresponding matplotlib Text objects. 857 858 Notes 859 ----- 860 - If `rotate=True`, label text is automatically aligned with the direction of the edge. 861 - Curvature (`rad`) enables visualization of bidirectional transitions using bent edges. 862 - This method is intended for internal use and supports bended edge label placement to match custom edge rendering. 863 864 See Also 865 -------- 866 draw 867 draw_networkx 868 draw_networkx_nodes 869 draw_networkx_edges 870 draw_networkx_labels 871 """ 872 if ax is None: 873 ax = plt.gca() 874 if edge_labels is None: 875 labels = {(u, v): d for u, v, d in self.edges(data=True)} 876 else: 877 labels = edge_labels 878 text_items = {} 879 for (n1, n2), label in labels.items(): 880 (x1, y1) = pos[n1] 881 (x2, y2) = pos[n2] 882 (x, y) = ( 883 x1 * label_pos + x2 * (1.0 - label_pos), 884 y1 * label_pos + y2 * (1.0 - label_pos), 885 ) 886 pos_1 = ax.transData.transform(np.array(pos[n1])) 887 pos_2 = ax.transData.transform(np.array(pos[n2])) 888 #linear_mid = 0.5*pos_1 + 0.5*pos_2 889 linear_mid = (1-label_pos)*pos_1 + label_pos*pos_2 890 d_pos = pos_2 - pos_1 891 rotation_matrix = np.array([(0,1), (-1,0)]) 892 ctrl_1 = linear_mid + rad*rotation_matrix@d_pos 893 ctrl_mid_1 = 0.5*pos_1 + 0.5*ctrl_1 894 ctrl_mid_2 = 0.5*pos_2 + 0.5*ctrl_1 895 bezier_mid = 0.5*ctrl_mid_1 + 0.5*ctrl_mid_2 896 (x, y) = ax.transData.inverted().transform(bezier_mid) 897 898 if rotate: 899 # in degrees 900 angle = np.arctan2(y2 - y1, x2 - x1) / (2.0 * np.pi) * 360 901 # make label orientation "right-side-up" 902 if angle > 90: 903 angle -= 180 904 if angle < -90: 905 angle += 180 906 # transform data coordinate angle to screen coordinate angle 907 xy = np.array((x, y)) 908 trans_angle = ax.transData.transform_angles( 909 np.array((angle,)), xy.reshape((1, 2)) 910 )[0] 911 else: 912 trans_angle = 0.0 913 # use default box of white with white border 914 if bbox is None: 915 bbox = dict(boxstyle="round", ec=(1.0, 1.0, 1.0), fc=(1.0, 1.0, 1.0)) 916 if not isinstance(label, str): 917 label = str(label) # this makes "1" and 1 labeled the same 918 919 t = ax.text( 920 x, y, label, 921 size=font_size, color=font_color, family=font_family, weight=font_weight, 922 alpha=alpha, 923 horizontalalignment=horizontalalignment, verticalalignment=verticalalignment, 924 rotation=trans_angle, transform=ax.transData, 925 bbox=bbox, zorder=1,clip_on=clip_on) 926 text_items[(n1, n2)] = t 927 928 ax.tick_params(axis="both", which="both",bottom=False, left=False, labelbottom=False, labelleft=False) 929 930 return text_items 931 932 933 def drawTransitionGraph(self, pos=None, 934 bended=False, node_size=1000, num=1, rad=-0.2, 935 edge_labels=None, node_shapes='o', figsize=(14, 7), 936 clear=True, fontsize=10, 937 fontcolor='black', label_pos=0.75): 938 """ 939 Visualize the state transition graph with nodes and directed edges. 940 941 This method draws the graph using `matplotlib` and `networkx`, including labels, 942 node colors, and optional curved (bended) edges for better visualization of 943 bidirectional transitions. It supports layout customization and is useful for 944 understanding the structure of the Markov model. 945 946 Parameters 947 ---------- 948 - pos : dict, optional 949 A dictionary mapping nodes to positions. If None, a Kamada-Kawai layout is used. 950 951 - bended : bool, optional 952 If True, edges are drawn with curvature using arcs. Useful for distinguishing 953 bidirectional transitions. Default is False. 954 955 - node_size : int, optional 956 Size of the nodes in the plot. Default is 1000. 957 958 - num : int, optional 959 Figure number for matplotlib (useful for managing multiple figures). Default is 1. 960 961 - rad : float, optional 962 The curvature radius for bended edges. Only used if `bended=True`. Default is -0.2. 963 964 - edge_labels : dict, optional 965 Optional dictionary mapping edges `(u, v)` to labels. If None, the `'label'` 966 attribute from each edge is used. 967 968 - node_shapes : str, optional 969 Shape of the nodes (e.g., `'o'` for circle, `'s'` for square). Default is `'o'`. 970 971 - figsize : tuple, optional 972 Size of the matplotlib figure. Default is (14, 7). 973 974 - clear : bool, optional 975 If True, clears the figure before plotting. Default is True. 976 977 - fontsize : int, optional 978 Font size used for node labels. Default is 10. 979 980 - fontcolor : str, optional 981 Font color for node labels. Default is `'black'`. 982 983 - label_pos : float, optional 984 Position of edge labels along the edge (0=start, 1=end). Default is 0.75. 985 986 Returns 987 ------- 988 - fig : matplotlib.figure.Figure 989 The created matplotlib figure object. 990 991 - ax : matplotlib.axes.Axes 992 The axes object where the graph is drawn. 993 994 - pos : dict 995 The dictionary of node positions used for drawing. 996 997 Notes 998 ----- 999 - Node labels and colors must be set beforehand using `addState()` or `setStateProperties()`. 1000 - Edge attributes must include `'label'` and `'color'`. 1001 - This method modifies the current matplotlib figure and is intended for interactive or inline use. 1002 """ 1003 node_labels = {node: data["label"] for node, data in self.nodes(data=True)} 1004 node_colors = {node: data["color"] for node, data in self.nodes(data=True)} 1005 node_colors = [data["color"] for node, data in self.nodes(data=True)] 1006 #node_shapes = {(hw,sw): "o" if hw+sw<G.N else "s" for hw,sw in G.nodes()} 1007 1008 edge_cols = {(u,v): self[u][v]["color"] for u,v in self.edges} 1009 1010 if edge_labels is None: 1011 edge_labels = {(u,v): self[u][v]["label"] for u,v in self.edges} 1012 1013 if pos is None: 1014 pos = nx.kamada_kawai_layout(self) 1015 1016 plt.figure(figsize=figsize, num=num, clear=clear) 1017 nx.draw_networkx_nodes(self, pos, node_color=node_colors, node_shape = node_shapes, node_size=node_size) 1018 1019 nx.draw_networkx_labels(self, pos, labels=node_labels, font_size=fontsize, font_color=fontcolor) # Draw node labels 1020 if bended: 1021 nx.draw_networkx_edges(self, pos, width=1, edge_color=[edge_cols[edge] for edge in self.edges], 1022 node_size = node_size, 1023 arrows = True, arrowstyle = '-|>', 1024 connectionstyle=f"arc3,rad={rad}") 1025 self._draw_networkx_edge_labels(pos, ax=plt.gca(), edge_labels=edge_labels,rotate=False, rad = rad, label_pos=label_pos) 1026 else: 1027 nx.draw_networkx_edges(self, pos, width=1, edge_color=[edge_cols[edge] for edge in self.edges], 1028 node_size = node_size, 1029 #min_target_margin=17, arrowsize=15, 1030 arrows = True, arrowstyle = '-|>') 1031 nx.draw_networkx_edge_labels(self, pos, edge_labels, label_pos=label_pos) #, verticalalignment='center',) 1032 plt.axis('off'); 1033 1034 return plt.gcf(), plt.gca(), pos 1035 1036 def animateSimulation(self, expectedTimePerState=1, inNotebook=False, **kwargs): 1037 """ 1038 Animate a simulation trajectory through the state transition graph. 1039 The animation of the CTMC simulation run either in a Jupyter notebook or as a regular script. 1040 1041 This method visualizes the path of a simulated or predefined sequence of states by 1042 temporarily highlighting each visited state over time. The time spent in each state 1043 is scaled relative to the average sojourn time to produce a visually smooth animation. 1044 1045 Parameters 1046 ---------- 1047 - **expectedTimePerState** : float, optional 1048 Approximate duration (in seconds) for an average state visit in the animation. 1049 The actual pause duration for each state is scaled proportionally to its sojourn time. 1050 Default is 1 second. 1051 1052 - **inNotebook** : bool, optional 1053 If True, uses Jupyter-compatible animation (with display + clear_output). 1054 If False, uses standard matplotlib interactive animation. Default is True. 1055 1056 - `**kwargs` : dict 1057 Additional keyword arguments including: 1058 - states : list 1059 A list of visited states (nodes) to animate. 1060 - times : list or np.ndarray 1061 Sojourn times in each state (same length as `states`). 1062 - startNode : hashable 1063 Optional start node for automatic simulation if `states` and `times` are not provided. 1064 - numEvents : int 1065 Number of events (state transitions) to simulate. 1066 - Additional drawing-related kwargs are passed to `drawTransitionGraph()`. 1067 1068 Raises 1069 ------ 1070 ValueError 1071 If neither a `(states, times)` pair nor `(startNode, numEvents)` are provided. 1072 1073 Returns 1074 ------- 1075 None 1076 The animation is shown interactively using matplotlib but no data is returned. 1077 1078 Notes 1079 ----- 1080 - If `states` and `times` are not supplied, a simulation is run via `simulateMarkovModel()`. 1081 - The function uses matplotlib's interactive mode (`plt.ion()`) to animate transitions. 1082 - The average sojourn time across all states is used to normalize animation speed. 1083 - Each visited state is highlighted in red for its corresponding (scaled) dwell time. 1084 - Drawing arguments like layout, font size, or color can be customized via kwargs. 1085 """ 1086 # expTimePerState: in seconds 1087 if "states" in kwargs and "times" in kwargs: 1088 states = kwargs["states"] 1089 times = kwargs["times"] 1090 else: 1091 # Run default simulation if data not provided 1092 if "startNode" in kwargs and "numEvents" in kwargs: 1093 states, times = self.simulateMarkovModel(startNode=None, numEvents=10) 1094 else: 1095 raise ValueError("Must provide either ('states' and 'times') or ('startNode' and 'numEvents').") 1096 1097 #avg_OutGoingRate = np.mean([self.nodes[n].get("sumOutgoingRates", 0) for n in self.nodes]) 1098 avg_Time = np.mean([1.0/self.nodes[n].get("sumOutgoingRates", 0) for n in self.nodes]) 1099 1100 1101 # Remove simulation-related keys before drawing 1102 draw_kwargs = {k: v for k, v in kwargs.items() if k not in {"states", "times", "startNode", "numEvents"}} 1103 1104 fig, ax, pos = self.drawTransitionGraph(**draw_kwargs) 1105 plt.tight_layout() 1106 1107 if not inNotebook: 1108 plt.ion() 1109 1110 for node, time in zip(states,times): 1111 if inNotebook: 1112 clear_output(wait=True) 1113 artist = nx.draw_networkx_nodes(self, pos, nodelist=[node], ax=ax, 1114 node_color=self._const.COLOR_NODE_ANIMATION_HIGHLIGHT, node_size=1000) 1115 display(fig) 1116 else: 1117 artist = nx.draw_networkx_nodes(self, pos, nodelist=[node],ax=ax, 1118 node_color=self._const.COLOR_NODE_ANIMATION_HIGHLIGHT, node_size=1000) 1119 plt.draw() 1120 plt.pause(time/avg_Time*expectedTimePerState) 1121 1122 # Remove highlighted node 1123 artist.remove() 1124 1125 pass 1126 1127 1128 1129 def states(self, data=False): 1130 """ 1131 Return a view of the graph's states (nodes), optionally with attributes. 1132 1133 This method is functionally identical to `self.nodes()` in NetworkX, but is 1134 renamed to `states()` to reflect the semantics of a state transition graph or 1135 Markov model, where nodes represent system states. 1136 1137 Parameters 1138 ---------- 1139 **data** : bool, optional 1140 If True, returns a view of `(node, attribute_dict)` pairs. 1141 If False (default), returns a view of node identifiers only. 1142 1143 Returns 1144 ------- 1145 *networkx.classes.reportviews.NodeView* or NodeDataView 1146 A view over the graph’s nodes or `(node, data)` pairs, depending on the `data` flag. 1147 1148 Examples 1149 -------- 1150 >>> G.states() 1151 ['A', 'B', 'C'] 1152 1153 >>> list(G.states(data=True)) 1154 [('A', {'color': 'red', 'label': 'Start'}), ('B', {...}), ...] 1155 1156 Notes 1157 ----- 1158 - This is a convenience wrapper for semantic clarity in models where nodes represent states. 1159 - The view is dynamic: changes to the graph are reflected in the returned view. 1160 """ 1161 return self.nodes(data=data) 1162 1163 1164 def setStateColor(self, state, color): 1165 """ 1166 Set the color attribute of a specific state (node). 1167 1168 Parameters 1169 ---------- 1170 **state** : hashable 1171 The identifier of the state whose color is to be updated. 1172 1173 **color** : str 1174 The new color to assign to the state (used for visualization). 1175 1176 Returns 1177 ------- 1178 None 1179 """ 1180 self.nodes[state]["color"] = color 1181 1182 1183 def setStateLabel(self, state, label): 1184 """ 1185 Set the label attribute of a specific state (node). 1186 1187 Parameters 1188 ---------- 1189 **state** : hashable 1190 The identifier of the state whose label is to be updated. 1191 1192 **label** : str 1193 The new label to assign to the state (used for visualization or annotation). 1194 1195 Returns 1196 ------- 1197 None 1198 """ 1199 self.nodes[state]["label"] = label 1200 1201 def prob(self, state): 1202 """ 1203 Return the steady-state probability of a specific state. 1204 1205 Parameters 1206 ---------- 1207 **state** : hashable 1208 The identifier of the state whose probability is to be retrieved. 1209 1210 Returns 1211 ------- 1212 **float** 1213 The steady-state probability of the specified state. 1214 1215 Raises 1216 ------ 1217 KeyError 1218 If the state is not found in the steady-state probability dictionary. 1219 """ 1220 return self.state_probabilities[state] 1221 1222 1223 def probs(self): 1224 """ 1225 Return the steady-state probabilities for all states. 1226 1227 Returns 1228 ------- 1229 **dict** 1230 A dictionary mapping each state to its steady-state probability. 1231 1232 Notes 1233 ----- 1234 - The probabilities are computed via `calcSteadyStateProbs()` and stored internally. 1235 - Use this method to access the full steady-state distribution. 1236 """ 1237 return self.state_probabilities 1238 1239 def getEmptySystemProbs(self, state=None): 1240 """ 1241 Create an initial state distribution where all probability mass is in a single state. 1242 1243 By default, this method returns a distribution where the entire probability mass 1244 is placed on the first node in the graph. Alternatively, a specific starting state 1245 can be provided. 1246 1247 Parameters 1248 ---------- 1249 **state** : hashable, optional 1250 The state to initialize with probability 1. If None, the first node 1251 in the graph (based on insertion order) is used. 1252 1253 Returns 1254 ------- 1255 **dict** 1256 A dictionary representing the initial distribution over states. 1257 The selected state has probability 1, and all others have 0. 1258 1259 Examples 1260 -------- 1261 >>> X0 = G.getEmptySystem() 1262 >>> sum(X0.values()) # 1.0 1263 1264 Notes 1265 ----- 1266 - This method is useful for setting up deterministic initial conditions 1267 in transient simulations of a CTMC. 1268 - The return format is compatible with methods like `transientProbs()`. 1269 """ 1270 if state is None: 1271 first_node = next(iter(self.nodes)) 1272 X0 = {s: 1 if s == first_node else 0 for s in self} 1273 return X0 1274 1275 1276 1277 def transientProbs(self, initialDistribution, t): 1278 """ 1279 Compute the transient state probabilities at time `t` for the CTMC. 1280 1281 This method solves the forward equation of a continuous-time Markov chain (CTMC): 1282 X(t) = X0 · expm(Q·t), where `Q` is the generator (rate) matrix and `X0` is the 1283 initial state distribution. It returns a dictionary mapping each state to its 1284 probability at time `t`. 1285 1286 Parameters 1287 ---------- 1288 - **initialDistribution** : dict 1289 A dictionary mapping states (nodes) to their initial probabilities at time `t=0`. 1290 The keys must match the nodes in the graph, and the values should sum to 1. 1291 1292 - **t** : float 1293 The time at which the transient distribution is to be evaluated. 1294 1295 Returns 1296 ------- 1297 **dict** 1298 A dictionary mapping each state (node) to its probability at time `t`. 1299 1300 Notes 1301 ----- 1302 - The rate matrix `Q` is generated automatically if not yet computed. 1303 - The state order in the rate matrix corresponds to the internal `_n2i` mapping. 1304 - Internally, the dictionary `initialDistribution` is converted into a vector 1305 aligned with the state indexing. 1306 - The transient solution is computed using `scipy.linalg.expm` for matrix exponentiation. 1307 1308 Raises 1309 ------ 1310 ValueError 1311 If the input dictionary has missing states or probabilities that do not sum to 1. 1312 (This is not enforced but recommended for correctness.) 1313 """ 1314 if self._rate_matrix is None: self.createRateMatrix() 1315 Q = self._rate_matrix 1316 X0 = np.zeros(len(initialDistribution)) 1317 1318 for i,n in enumerate(self.nodes): 1319 X0[ self._n2i[n] ] = initialDistribution[n] 1320 1321 1322 X = X0 @ expm(Q*t) 1323 1324 matrix = {} 1325 for n in self.nodes: 1326 matrix[n] = X[ self._n2i[n] ] 1327 1328 return matrix 1329 1330 def symSolveMarkovModel(self): 1331 """ 1332 Symbolically solve the global balance equations of the CTMC. 1333 1334 This method constructs and solves the system of symbolic equations for the 1335 steady-state probabilities of each state in the CTMC. It uses symbolic transition 1336 rates stored on the edges to form balance equations for each node, along with 1337 the normalization constraint that all probabilities sum to 1. 1338 1339 Returns 1340 ------- 1341 - **solution** : dict 1342 A symbolic solution dictionary mapping SymPy symbols (e.g., x_{A}) to 1343 expressions in terms of model parameters. 1344 1345 - **num_solution** : dict 1346 A numerical dictionary mapping each state (node) to its steady-state probability, 1347 computed by substituting the numeric values from `self._subs` into the symbolic solution. 1348 1349 Notes 1350 ----- 1351 - For each node `s`, a symbolic variable `x_{label}` is created using the node's label attribute. 1352 - One balance equation is created per state: total inflow = total outflow. 1353 - An additional constraint `sum(x_i) = 1` ensures proper normalization. 1354 - The symbolic system is solved with `sympy.solve`, and the results are simplified. 1355 - Numeric values are computed by substituting known parameter values (`self._subs`) into the symbolic solution. 1356 1357 Examples 1358 -------- 1359 >>> sym_sol, num_sol = G.symSolveMarkovModel() 1360 >>> sym_sol # symbolic expressions for each state 1361 >>> num_sol # evaluated numerical values for each state 1362 """ 1363 X={} 1364 to_solve = [] 1365 for s in self.nodes(): 1366 X[s] = sp.Symbol(f'x_{{{self.nodes[s]["label"]}}}') 1367 to_solve.append(X[s]) 1368 #% 1369 eqs = [] 1370 for s in self.nodes(): 1371 succ = self.successors(s) # raus 1372 out_rate = 0 1373 for z in list(succ): 1374 #out_rate += X[s]*self.edges[s, z]["sym_rate"] 1375 out_rate += X[s]*self[s][z]["sym_rate"] 1376 1377 in_rate = 0 1378 for r in list(self.predecessors(s)): 1379 #in_rate += X[r]*self.edges[r,s]["sym_rate"] 1380 in_rate += X[r]*self[r][s]["sym_rate"] 1381 eqs.append( sp.Eq(out_rate, in_rate) ) 1382 #% 1383 Xsum = 0 1384 for s in X: 1385 Xsum += X[s] 1386 1387 eqs.append( sp.Eq(Xsum, 1) ) 1388 #% 1389 solution = sp.solve(eqs, to_solve) 1390 simplified_solution = sp.simplify(solution) 1391 #print(simplified_solution) 1392 num_solution = simplified_solution.subs(self._subs) 1393 #num_dict = {str(var): str(sol) for var, sol in simplified_solution.items()} 1394 num_dict = {s: num_solution[X[s]] for s in X} 1395 1396 return simplified_solution, num_dict 1397 1398 def exportTransitionsToExcel(self, excel_path = "output.xlsx", num_colors=9, open_excel=False, verbose=True): 1399 """ 1400 Export transition data from the state transition graph to an Excel file. 1401 1402 This method collects all transitions (edges) from the graph and writes them 1403 to an Excel spreadsheet, optionally applying background colors to distinguish 1404 outgoing transitions from different source states. This export is helpful for 1405 manually verifying the transition structure of the model. 1406 1407 The transition type is defined by the `tid` field of each edge, and mapped to 1408 a human-readable string via the `EXPORT_NAME_TRANSITIONS` dictionary in the 1409 constants module (e.g., `Constants.EXPORT_NAME_TRANSITIONS[tid]`). 1410 1411 See Also 1412 -------- 1413 StateTransitionGraph.addTransition : Method to add edges with symbolic rates and `tid`. 1414 StateTransitionGraph.__init__ : Initializes the graph and loads the constants module. 1415 1416 Parameters 1417 ---------- 1418 - **excel_path** : str, optional 1419 Path to the output Excel file. Default is "output.xlsx". 1420 1421 - num_colors : int, optional 1422 Number of unique background colors to cycle through for different states. 1423 Default is 9. 1424 1425 - open_excel : bool, optional 1426 If True, automatically opens the generated Excel file after export. 1427 Default is False. 1428 1429 - verbose : bool, optional 1430 If True, prints summary information about the export process. 1431 Default is True. 1432 1433 Returns 1434 ------- 1435 None 1436 The function writes data to disk and optionally opens Excel, 1437 but does not return a value. 1438 1439 Notes 1440 ----- 1441 - Each edge must contain a `tid` attribute for transition type and a `label` for the rate. 1442 - Background coloring groups all outgoing transitions from the same source state. 1443 - Requires the constants module to define `EXPORT_NAME_TRANSITIONS`. 1444 """ 1445 1446 rows = [] 1447 for node in self.nodes(): 1448 edges = list(self.edges(node)) 1449 if verbose: print(f"Edges connected to node '{node}':") 1450 1451 for i, (u, v) in enumerate(edges): 1452 latex = self[u][v]["label"] 1453 rate = sp.pretty(self[u][v]["sym_rate"]) 1454 tid = self[u][v]["tid"] 1455 1456 if tid: # 1457 if hasattr(self._const, 'EXPORT_NAME_TRANSITIONS'): 1458 tid_name = self._const.EXPORT_NAME_TRANSITIONS[tid] 1459 else: 1460 tid_name = '' 1461 1462 rows.append({"from": u, "type": tid_name, "to": v, "rate": rate, "latex": latex}) 1463 1464 df = pd.DataFrame(rows, columns=["from", "type", "to", "rate", "latex"]) 1465 df["checked"] = "" 1466 df.to_excel(excel_path, index=False) 1467 1468 # Load workbook with openpyxl 1469 wb = load_workbook(excel_path) 1470 ws = wb.active 1471 1472 # Define colors for different "from" values 1473 cmap = map(plt.cm.Pastel1, range(num_colors)) 1474 # Convert RGBA to hex 1475 color_list = [mcolors.to_hex(color, keep_alpha=False).upper().replace("#", "") for color in cmap] 1476 1477 # Assign color to each unique "from" value using modulo 1478 unique_from = (df["from"].unique()) # Sorted for consistent order 1479 from_color_map = { 1480 value: color_list[i % len(color_list)] for i, value in enumerate(unique_from) 1481 } 1482 1483 # Apply coloring row by row 1484 for i,row in enumerate(ws.iter_rows(min_row=2, max_row=ws.max_row)): # Skip header 1485 from_value = row[0].value # "from" is first column 1486 #fill_color = color_list[i % len(color_list)] 1487 if isinstance(from_value, int): 1488 fill_color = from_color_map[from_value] 1489 else: 1490 fill_color = from_color_map.get(eval(from_value)) 1491 1492 if fill_color: 1493 fill = PatternFill(start_color=fill_color, end_color=fill_color, fill_type="solid") 1494 for cell in row: 1495 cell.fill = fill 1496 1497 # Save styled Excel file 1498 wb.save(excel_path) 1499 if verbose: print(f"File {excel_path} created.") 1500 if open_excel: 1501 print(f"Opening the file {excel_path}.") 1502 os.startfile(excel_path) 1503
95class StateTransitionGraph(nx.DiGraph): 96 """ 97 Initialize the StateTransitionGraph for continuous time Markov chains (CTMC). 98 Extends the networkx.Graph class to include additional methods and properties 99 for the analysis, visualization, and simulation of CTMCs. 100 The analysis allows deriving the state probabilities in the steady-state 101 and in the transient phase. 102 103 This constructor sets up symbolic and numeric parameters for the state transition graph, 104 imports a constants module dynamically, and initializes internal data structures for 105 rate matrix computation and steady-state probability analysis. 106 107 Parameters 108 ---------- 109 - **parameters** : dict 110 A dictionary mapping parameter names to their numerical values. 111 Example: {'lambda': 1.0, 'mu': 1.0}. 112 For each parameter key: 113 - A symbolic variable is created using sympy.symbols. Can be accessed as `G.sym("key")` or `G.sym_key` 114 - A variable is created which can be accessed via `G.key` (and the numerical value is assigned) 115 - For the key="lambda", `G.lam` and `G.sym_lambda` are created 116 117 The parameters can be accessed via property `parameters` 118 - **constants** : str, optional 119 The name of a Python module (as string) to import dynamically and assign to `self.const`. 120 This can be used to store global or user-defined constants. 121 Default is "Constants". 122 - `*args` : tuple 123 Additional positional arguments passed to the base `networkx.DiGraph` constructor. 124 - `**kwargs` : dict 125 Additional keyword arguments passed to the base `networkx.DiGraph` constructor. 126 127 Attributes 128 ---------- 129 - `_rate_matrix` : ndarray or None 130 The continuous-time rate matrix `Q` of the CTMC. Computed when calling `createRateMatrix()`. 131 - `_state_probabilities` : dict or None 132 Dictionary of steady-state probabilities, computed via `calcSteadyStateProbs()`. 133 - `_n2i` : dict or None 134 Mapping from node identifiers to matrix row/column indices used for rate matrix construction. 135 - `_sym` : dict 136 Mapping of parameter names (strings) to their symbolic `sympy` representations. 137 - `_subs` : dict 138 Mapping from symbolic parameters to their numerical values (used for substitutions). 139 - `_parameters` : dict 140 Copy of the original parameter dictionary passed at initialization. 141 - `_transitionProbabilitesComputed` : bool 142 Internal flag indicating whether transition probabilities for the embedded Markov chain 143 have been computed. 144 - `const` : module 145 Dynamically imported module (default `"Constants"`), used for storing user-defined constants 146 such as colors or transition IDs. 147 148 Inherits 149 -------- 150 networkx.Graph: The base graph class from NetworkX. 151 """ 152 153 def __init__(self, parameters, constants="Constants", *args, **kwargs): 154 """ 155 Initialize the state transition graph with symbolic parameters. 156 157 This constructor sets up the symbolic and numerical environment for modeling 158 a CTMC. Parameters are stored as both symbolic expressions and numerical values, 159 and convenient attribute access is provided for both forms. 160 161 Parameters 162 ---------- 163 - **parameters** : dict 164 A dictionary mapping parameter names to their numerical values. 165 Example: {'lambda': 1.0, 'mu': 1.0}. 166 For each parameter key: 167 - A symbolic variable is created using sympy.symbols. Can be accessed as `G.sym("key")` or `G.sym_key` 168 - A variable is created which can be accessed via `G.key` (and the numerical value is assigned) 169 - For the key="lambda", `G.lam` and `G.sym_lambda` are created 170 171 The parameters can be accessed via property `parameters` 172 - **constants** : str, optional 173 The name of a Python module (as string) to import dynamically and assign to `self.const`. 174 This can be used to store global or user-defined constants. 175 Default is "Constants". 176 - `*args` : tuple 177 Additional positional arguments passed to the base `networkx.DiGraph` constructor. 178 - `**kwargs` : dict 179 Additional keyword arguments passed to the base `networkx.DiGraph` constructor. 180 """ 181 super().__init__(*args, **kwargs) 182 self._rate_matrix = None # Rate Matrix Q of this transition graph 183 self._state_probabilities = None # steady state proabilities of this CTMC 184 self._n2i = None # internal index for generating the rate matrix Q: dictionary (node to index) 185 186 #parameters = {'lambda': 1.0, 'mu':1.0, 'n':5} 187 self._parameters = parameters 188 self._sym = {key: sp.symbols(key) for key in parameters} 189 self._subs = {} # contains symbolic parameter as key and provides numerical value for computation. 190 for k in parameters: 191 self._subs[ self._sym[k] ] = parameters[k] 192 193 # access to parameters and symbols using dot-notation, e.g. G.lam, G.sym_lambda, or G.sym("lambda") 194 for key, value in parameters.items(): 195 if key=='lambda': 196 setattr(self, 'lam', value) 197 else: 198 setattr(self, key, value) 199 200 setattr(self, "sym_"+key, self._sym[key]) 201 202 self._transitionProbabilitesComputed = False 203 204 self._const = importlib.import_module(constants) 205 206 207 208 def addState(self, state, color=None, label=None): 209 """ 210 Add a state (node) to the graph with optional color and label attributes. 211 212 If the node already exists, its attributes are updated with the provided values. 213 214 Parameters 215 ---------- 216 - **state** : hashable 217 The identifier for the state. Can be a string, int, or a tuple. If no label 218 is provided, the label will be auto-generated as: 219 - comma-separated string for tuples (e.g., `(1, 2)` -> `'1,2'`) 220 - stringified value for other types 221 222 - **color** : str, optional 223 The color assigned to the state (used for visualization). 224 If not specified, defaults to `self.const.COLOR_NODE_DEFAULT`. 225 226 - **label** : str, optional 227 A label for the node (used for display or annotation). If not provided, 228 it will be auto-generated from the state identifier. 229 230 Notes 231 ----- 232 - This method wraps `self.add_node()` from NetworkX. 233 - If the node already exists in the graph, this method does not raise an error; 234 it will update the node’s existing attributes. 235 """ 236 if label is None: 237 if isinstance(state, tuple): 238 label = ','.join(map(str, state)) # if state is a tuple 239 else: 240 label = str(state) 241 #self.nodes[state]["label"] = label 242 243 if color is None: 244 color = self._const.COLOR_NODE_DEFAULT 245 246 # no error if the node already exists.Attributes get updated if you pass them again. 247 self.add_node(state, color=color, label=label) 248 249 def setStateProperties(self, state, color=None, label=None, **kwargs): 250 """ 251 Update visual or metadata properties for an existing state (node) in the graph. 252 253 This method updates the node's attributes such as color, label, and any additional 254 custom key-value pairs. The node must already exist in the graph. 255 256 Parameters 257 ---------- 258 - **state** : hashable 259 Identifier of the state to update. Must be a node already present in the graph. 260 261 - **color** : str, optional 262 Color assigned to the node, typically used for visualization. If not provided, 263 the default node color (`self._const.COLOR_NODE_DEFAULT`) is used. 264 265 - **label** : str, optional 266 Label for the node. If not provided, it is auto-generated from the state ID 267 (especially for tuples). 268 269 - `**kwargs` : dict 270 Arbitrary additional key-value pairs to set as node attributes. These can include 271 any custom metadata. 272 273 Raises 274 ------ 275 KeyError 276 If the specified state is not present in the graph. 277 278 Notes 279 ----- 280 - This method ensures that required attributes (like `color` and `label`) are set for visualization. 281 - Existing node attributes are updated or extended with `kwargs`. 282 """ 283 if state not in self: 284 raise KeyError(f"State '{state}' does not exist in the graph.") 285 286 # mandatory for plotting 287 self.addState(state, color, label) 288 289 # # Add any additional attributes 290 self.nodes[state].update(kwargs) 291 292 def setAllStateDefaultProperties(self, **kwargs): 293 """ 294 Set default visual and metadata properties for all states (nodes) in the graph. 295 296 This method applies `setStateProperties()` to every node, ensuring each state 297 has at least the required default attributes (e.g., color and label). 298 Any additional keyword arguments are passed to `setStateProperties()`. 299 300 Parameters 301 ---------- 302 `**kwargs` : dict, optional 303 Optional keyword arguments to apply uniformly to all nodes, 304 such as color='gray', group='core', etc. 305 306 Notes 307 ----- 308 - This method is useful for initializing node attributes before plotting 309 or performing other graph-based computations. 310 - Nodes without existing attributes will receive defaults via `addState()`. 311 """ 312 for state in self.nodes(): 313 self.setStateProperties(state, **kwargs) 314 315 316 #def addTransition(self, origin_state, list_destination_state, sym_rate, color=plt.cm.tab10(5), tid = None): 317 def addTransition(self, origin_state, destination_state, sym_rate, tid = None, color=None, label=None): 318 """ 319 Add a directed transition (edge) between two states with a symbolic transition rate. 320 321 This method adds a directed edge from `origin_state` to `destination_state`, 322 stores both the symbolic and numerical rate, and sets optional display properties 323 such as color, label, and transition ID. If the edge already exists, a warning is issued. 324 325 Parameters 326 ---------- 327 - **origin_state** : hashable 328 The source node (state) of the transition. 329 330 - **destination_state** : hashable 331 The target node (state) of the transition. 332 333 - **sym_rate** : sympy.Expr 334 A symbolic expression representing the transition rate. 335 This will be numerically evaluated using internal substitutions (`self._subs`). 336 337 - **tid** : str or int, optional 338 A transition identifier (e.g. for grouping or styling). 339 If provided, the color will default to `self._const.COLOR_TRANSITIONS[tid]` if not set explicitly. 340 341 - **color** : str, optional 342 Color for the edge (used in visualization). Defaults to the transition-specific color 343 if `tid` is provided. 344 345 - **label** : str, optional 346 LaTeX-formatted string for labeling the edge. If not provided, a default label is 347 generated using `sympy.latex(sym_rate)`. 348 349 Notes 350 ----- 351 - The numeric rate is computed by substituting known parameter values into `sym_rate`. 352 - The following attributes are stored on the edge: 353 - `rate` (float): numeric value of the transition rate 354 - `sym_rate` (sympy.Expr): original symbolic expression 355 - `label` (str): display label 356 - `color` (str): edge color 357 - `tid` (str or int): transition ID 358 - If the edge already exists, a `UserWarning` is issued but the edge is still overwritten. 359 """ 360 if tid is not None: 361 if color is None: 362 color = self._const.COLOR_TRANSITIONS[ tid ] 363 364 rate = sym_rate.subs( self._subs ) 365 366 if label is None: string = "$"+sp.latex(sym_rate)+"$" 367 368 if self.has_edge(origin_state, destination_state): 369 theedge = self[origin_state][destination_state] 370 warningstr = f"Edge already exists: {origin_state} -> {destination_state} with rate {theedge['label']}" 371 warnings.warn(warningstr, category=UserWarning) 372 self.add_edge(origin_state, destination_state, rate=float(rate), label=string, 373 color=color, sym_rate=sym_rate, tid=tid) 374 375 def sym(self, key): 376 """ 377 Retrieve the symbolic representation of a model parameter. 378 379 Parameters 380 ---------- 381 **key** : str 382 The name of the parameter (e.g., 'lambda', 'mu', 'n'). 383 384 Returns 385 ------- 386 **sympy.Symbol** 387 The symbolic variable corresponding to the given key. 388 389 Raises 390 ------ 391 KeyError 392 If the key does not exist in the internal symbolic mapping (`_sym`). 393 394 Examples 395 -------- 396 >>> G.sym('lambda') 397 lambda 398 399 Notes 400 ----- 401 - The symbolic mapping is initialized during object construction based on the `parameters` dictionary. 402 - To access the symbol using dot-notation, use `.sym_key` attributes created during initialization. 403 """ 404 return self._sym[key] 405 406 407 @property 408 def parameters(self): 409 """ 410 Return the model parameters as a dictionary. 411 412 This includes the original numerical values provided during initialization. 413 414 Returns 415 ------- 416 **dict** 417 Dictionary of parameter names and their corresponding numeric values. 418 419 Examples 420 -------- 421 >>> G.parameters 422 {'lambda': 1.0, 'mu': 0.5, 'n': 3} 423 424 Notes 425 ----- 426 - This property does not return symbolic values. Use `self.sym(key)` for symbolic access. 427 - Parameter symbols can also be accessed using attributes like `sym_lambda` or the `sym()` method. 428 """ 429 return self._parameters 430 431 432 433 @property 434 def rate_matrix(self): 435 """ 436 Return the continuous-time rate matrix (Q matrix) of the state transition graph. 437 438 This matrix defines the transition rates between all pairs of states in the 439 continuous-time Markov chain (CTMC) represented by this graph. 440 441 Returns 442 ------- 443 **numpy.ndarray** or None 444 The rate matrix `Q`, where `Q[i, j]` is the transition rate from state `i` to `j`. 445 Returns `None` if the rate matrix has not been computed yet. 446 447 Notes 448 ----- 449 - The matrix is typically generated via a method like `computeRateMatrix()` or similar. 450 - Internal state indexing is handled by `_n2i` (node-to-index mapping). 451 - The rate matrix is stored internally as `_rate_matrix`. 452 """ 453 return self._rate_matrix 454 455 456 457 @property 458 def state_probabilities(self): 459 """ 460 Return the steady-state probabilities of the continuous-time Markov chain (CTMC). 461 462 These probabilities represent the long-run proportion of time the system spends 463 in each state, assuming the CTMC is irreducible and positive recurrent. 464 465 Returns 466 ------- 467 **dict** or None 468 A dictionary mapping each state to its steady-state probability. 469 Returns `None` if the probabilities have not yet been computed. 470 471 Notes 472 ----- 473 - The probabilities are stored internally after being computed via a dedicated method 474 `calcSteadyStateProbs()`. 475 - The result is stored in the `_state_probabilities` attribute. 476 """ 477 return self._state_probabilities 478 479 480 481 def printEdges(self, state): 482 """ 483 Print all outgoing edges from a given state in the graph. 484 485 For each edge from the specified `state`, this method prints the destination 486 state and the corresponding transition label (typically the symbolic rate). 487 488 Parameters 489 ---------- 490 **state** : hashable 491 The source node (state) from which outgoing edges will be printed. 492 493 Returns 494 ------- 495 None 496 This method prints output directly to the console and does not return a value. 497 498 Notes 499 ----- 500 - The transition label is retrieved from the edge attribute `'label'`. 501 - Assumes `state` exists in the graph; otherwise, a KeyError will be raised. 502 """ 503 print(f"from {state}") 504 node = self[state] 505 for e in node: 506 rate = node[e]['label'] 507 print(f" -> {e} with {rate}") 508 509 510 def createRateMatrix(self): 511 """ 512 Construct the continuous-time Markov chain (CTMC) rate matrix Q from the graph. 513 514 This method generates the generator matrix `Q` corresponding to the structure and 515 rates defined in the directed graph. Each node represents a system state, and each 516 directed edge must have a `'rate'` attribute representing the transition rate from 517 one state to another. 518 519 The matrix `Q` is constructed such that: 520 - `Q[i, j]` is the transition rate from state `i` to state `j` 521 - Diagonal entries satisfy `Q[i, i] = -sum(Q[i, j]) for j ≠ i`, so that each row sums to 0 522 523 The resulting matrix and the internal node-to-index mapping (`_n2i`) are stored 524 as attributes of the graph. 525 526 Returns 527 ------- 528 None 529 The method updates internal attributes: 530 - `_rate_matrix` : numpy.ndarray 531 The CTMC rate matrix Q 532 - `_n2i` : dict 533 Mapping from node identifiers to matrix row/column indices 534 535 Notes 536 ----- 537 - Assumes all edges in the graph have a `'rate'` attribute. 538 - The node order used in the matrix is consistent with the order from `list(self.nodes)`. 539 - Use `self.rate_matrix` to access the generated matrix after this method is called. 540 """ 541 n2i = {} 542 nodes = list(self.nodes) 543 for i,node in enumerate(nodes): 544 n2i[node] = i 545 Q = np.zeros((len(nodes),len(nodes))) # rate matrix 546 547 for edge in self.edges: 548 i0 = n2i[edge[0]] 549 i1 = n2i[edge[1]] 550 Q[i0,i1] = self[edge[0]][edge[1]]["rate"] 551 552 np.fill_diagonal(Q, -Q.sum(axis=1)) 553 self._rate_matrix = Q 554 self._n2i = n2i 555 #return Q, n2i 556 557 def calcSteadyStateProbs(self, verbose=False): 558 """ 559 Compute the steady-state probabilities of the continuous-time Markov chain (CTMC). 560 561 This method calculates the stationary distribution of the CTMC defined by the rate 562 matrix of the state transition graph. The computation solves the linear system 563 `XQ = 0` with the constraint that the probabilities X sum to 1, using a modified 564 version of the rate matrix `Q`. 565 566 Parameters 567 ---------- 568 `verbose` : bool, optional 569 If True, prints the rate matrix, the modified matrix, the right-hand side vector, 570 and the resulting steady-state probabilities. 571 572 Returns 573 ------- 574 **probs** : dict 575 A dictionary mapping each state (node in the graph) to its steady-state probability. 576 577 Notes 578 ----- 579 - The graph must already contain valid transition rates on its edges. 580 - If the rate matrix has not yet been created, `createRateMatrix()` is called automatically. 581 - The method modifies the last column of the rate matrix to impose the normalization constraint 582 `sum(X) = 1`. 583 - The inverse of the modified matrix is computed directly, so the graph should be small enough 584 to avoid numerical instability or performance issues. 585 - The resulting probabilities are stored internally in `_state_probabilities` and can be accessed 586 via the `state_probabilities` property. 587 """ 588 589 # compute transition matrix if not already done 590 if self._rate_matrix is None: self.createRateMatrix() 591 Q2 = self._rate_matrix.copy() 592 if verbose: 593 print("Rate matrix Q") 594 print(f'Q=\n{Q2}\n') 595 596 Q2[:, -1] = 1 597 if verbose: 598 print(f'Matrix is changed to\nQ2=\n{Q2}\n') 599 600 b = np.zeros(len(Q2)) 601 b[-1] = 1 602 if verbose: 603 print("Solve X = b @ Q2") 604 print(f'b=\n{b}\n') 605 606 # state probabilities 607 X = b @ linalg.inv(Q2) # compute the matrix inverse 608 609 # Generate a matrix with P(X,Y) 610 matrix = {} 611 for n in self.nodes: 612 matrix[n] = X[ self._n2i[n] ] 613 614 if verbose: 615 print("Steady-state probabilities X") 616 print(f'X=\n{matrix}\n') 617 self._state_probabilities = matrix 618 return matrix 619 620 def calcTransitionProbabilities(self): 621 """ 622 Compute and store transition probabilities for each node in the graph of the embedded discrete Markov chain. 623 624 This method processes each node in the graph, treating outgoing edges as 625 transitions in an embedded discrete-time Markov chain derived from a continuous-time 626 Markov chain (CTMC). Transition *rates* on outgoing edges are normalized into 627 *probabilities*, and relevant data is stored in the node attributes. 628 629 For each node, the following attributes are set: 630 - `"transitionProbs"` : np.ndarray 631 Array of transition probabilities to successor nodes. 632 - `"transitionNodes"` : list 633 List of successor node identifiers (ordered as in the probability array). 634 - `"sumOutgoingRates"` : float 635 Sum of all outgoing transition rates. Used to model the exponential waiting time 636 in continuous-time simulation (`T ~ Exp(sumOutgoingRates)`). 637 638 Notes 639 ----- 640 - The graph must be a directed graph where each edge has a `'rate'` attribute. 641 - Nodes with no outgoing edges are considered absorbing states. 642 - Useful for discrete-event simulation and Markov chain Monte Carlo (MCMC) modeling. 643 - Sets an internal flag `_transitionProbabilitesComputed = True` upon completion. 644 645 Raises 646 ------ 647 KeyError 648 If any outgoing edge lacks a `'rate'` attribute. 649 """ 650 for node in self: 651 successors = list(self.successors(node)) 652 rates = np.array([self[node][nbr]['rate'] for nbr in successors]) 653 probs = rates / rates.sum() # Normalize to probabilities 654 self.nodes[node]["transitionProbs"] = probs # probabilities 655 self.nodes[node]["transitionNodes"] = successors # list of nodes 656 self.nodes[node]["sumOutgoingRates"] = rates.sum() # time in state ~ EXP(sumOutgoingRates) 657 numel = len(successors) 658 if numel==0: 659 print(f"Absorbing state: {node}") 660 self._transitionProbabilitesComputed = True 661 662 def simulateMarkovModel(self, startNode=None, numEvents=10): 663 """ 664 Simulate a trajectory through the state transition graph as a continuous-time Markov chain (CTMC). 665 666 This method performs a discrete-event simulation of the CTMC by generating a random walk 667 through the graph, starting at `startNode` (or a default node if not specified). Transitions 668 are chosen based on the normalized transition probabilities computed from edge rates, and 669 dwell times are sampled from exponential distributions based on the total outgoing rate 670 from each state. 671 672 Parameters 673 ---------- 674 - **startNode** : hashable, optional 675 The node at which the simulation starts. If not specified, the first node in the graph 676 (based on insertion order) is used. 677 678 - **numEvents** : int, optional 679 The number of transition events (steps) to simulate. Default is 10. 680 681 Returns 682 ------- 683 **states** : list 684 A list of visited states, representing the sequence of nodes traversed in the simulation. 685 686 **times** : numpy.ndarray 687 An array of dwell times in each state, sampled from exponential distributions 688 with rate equal to the sum of outgoing transition rates from each state. 689 690 Notes 691 ----- 692 - This method automatically calls `calcTransitionProbabilities()` if transition probabilities 693 have not yet been computed. 694 - If an absorbing state (i.e., a node with no outgoing edges) is reached, the simulation stops early. 695 - The cumulative time can be obtained by `np.cumsum(times)`. 696 - State durations follow `Exp(sumOutgoingRates)` distributions per CTMC theory. 697 """ 698 if not self._transitionProbabilitesComputed: self.calcTransitionProbabilities() 699 700 # simulate random walk through Graph = Markov Chain 701 if startNode is None: 702 startNode = next(iter(self)) 703 704 705 states = [] # list of states 706 times = np.zeros(numEvents) # time per state 707 708 node = startNode 709 for i in range(numEvents): 710 outgoingRate = self.nodes[node]["sumOutgoingRates"] # time in state ~ EXP(sumOutgoingRates) 711 times[i] = random.expovariate(outgoingRate) # exponentially distributed time in state 712 states.append(node) 713 714 # get next node 715 probs = self.nodes[node]["transitionProbs"] # probabilities 716 successors = self.nodes[node]["transitionNodes"] # list of nodes 717 numel = len(successors) 718 if numel==0: 719 print(f"Absorbing state: {node}") 720 break 721 elif numel==1: 722 nextNode = successors[0] 723 else: 724 nextNode = successors[ np.random.choice(numel, p=probs) ] 725 # edge = G[node][nextNode] 726 node = nextNode 727 728 return states, times 729 730 def calcProbabilitiesSimulationRun(self, states, times): 731 """ 732 Estimate empirical state probabilities from a simulation run of the CTMC. 733 734 This method computes the fraction of total time spent in each state during a simulation. 735 It uses the sequence of visited states and corresponding sojourn times to estimate the 736 empirical (time-weighted) probability distribution over states. 737 738 Parameters 739 ---------- 740 - **states** : list of tuple 741 A list of state identifiers visited during the simulation. 742 - **times** : numpy.ndarray 743 An array of sojourn times corresponding to each visited state, same length as `states`. 744 745 Returns 746 ------- 747 **simProbs** : dict 748 A dictionary mapping each unique state to its estimated empirical probability, 749 based on total time spent in that state relative to the simulation duration. 750 751 Notes 752 ----- 753 - Internally, states are mapped to flat indices for efficient counting using NumPy. 754 - This method assumes that `states` contains only valid node identifiers from the graph. 755 - If `_n2i` (node-to-index mapping) has not been initialized, it will be computed here. 756 - This is especially useful for validating analytical steady-state probabilities 757 against sampled simulation results. Or in case of state-space explosions. 758 """ 759 if self._n2i is None: 760 n2i = {} 761 nodes = list(self.nodes) 762 for i,node in enumerate(nodes): 763 n2i[node] = i 764 self._n2i = n2i 765 766 767 # Map tuples to flat indices 768 data = np.array(states) 769 max_vals = data.max(axis=0) + 1 770 indices = np.ravel_multi_index(data.T, dims=max_vals) 771 772 # Weighted bincount 773 bincount = np.bincount(indices, weights=times) 774 775 # Decode back to tuples 776 tuple_indices = np.array(np.unravel_index(np.nonzero(bincount)[0], shape=max_vals)).T 777 timePerState = bincount[bincount > 0] 778 779 simProbs = {} 780 totalTime = times.sum() 781 for t, c in zip(tuple_indices, timePerState): 782 simProbs[tuple(t)] = c / totalTime 783 #print(f"Tuple {tuple(t)} → weighted count: {c}") 784 return simProbs 785 786 787 def _draw_networkx_edge_labels( 788 self, 789 pos, 790 edge_labels=None, label_pos=0.75, 791 font_size=10, font_color="k", font_family="sans-serif", font_weight="normal", 792 alpha=None, bbox=None, 793 horizontalalignment="center", verticalalignment="center", 794 ax=None, rotate=True, clip_on=True, rad=0): 795 """ 796 Draw edge labels for a directed graph with optional curved edge support. 797 798 This method extends NetworkX's default edge label drawing by supporting 799 curved (bended) edges through quadratic Bézier control points. It places 800 labels dynamically at a specified position along the curve and optionally 801 rotates them to follow the edge orientation. 802 803 Parameters 804 ---------- 805 pos : dict 806 Dictionary mapping nodes to positions (2D coordinates). Required. 807 808 edge_labels : dict, optional 809 Dictionary mapping edge tuples `(u, v)` to label strings. If None, 810 edge data is used and converted to strings automatically. 811 812 label_pos : float, optional 813 Position along the edge to place the label (0=head, 1=tail). 814 Default is 0.75 (closer to the target node). 815 816 font_size : int, optional 817 Font size of the edge labels. Default is 10. 818 819 font_color : str, optional 820 Font color for the edge labels. Default is 'k' (black). 821 822 font_family : str, optional 823 Font family for the edge labels. Default is 'sans-serif'. 824 825 font_weight : str, optional 826 Font weight for the edge labels (e.g., 'normal', 'bold'). Default is 'normal'. 827 828 alpha : float, optional 829 Opacity for the edge labels. If None, labels are fully opaque. 830 831 bbox : dict, optional 832 Matplotlib-style bbox dictionary to style the label background. 833 Default is a white rounded box. 834 835 horizontalalignment : str, optional 836 Horizontal alignment of the text. Default is 'center'. 837 838 verticalalignment : str, optional 839 Vertical alignment of the text. Default is 'center'. 840 841 ax : matplotlib.axes.Axes, optional 842 Axes on which to draw the labels. If None, uses the current axes. 843 844 rotate : bool, optional 845 Whether to rotate the label to match the edge direction. Default is True. 846 847 clip_on : bool, optional 848 Whether to clip labels at the plot boundary. Default is True. 849 850 rad : float, optional 851 Curvature of the edge (used to determine Bézier control point). 852 Positive values bend counterclockwise; negative values clockwise. Default is 0 (straight). 853 854 Returns 855 ------- 856 dict 857 Dictionary mapping edge tuples `(u, v)` to the corresponding matplotlib Text objects. 858 859 Notes 860 ----- 861 - If `rotate=True`, label text is automatically aligned with the direction of the edge. 862 - Curvature (`rad`) enables visualization of bidirectional transitions using bent edges. 863 - This method is intended for internal use and supports bended edge label placement to match custom edge rendering. 864 865 See Also 866 -------- 867 draw 868 draw_networkx 869 draw_networkx_nodes 870 draw_networkx_edges 871 draw_networkx_labels 872 """ 873 if ax is None: 874 ax = plt.gca() 875 if edge_labels is None: 876 labels = {(u, v): d for u, v, d in self.edges(data=True)} 877 else: 878 labels = edge_labels 879 text_items = {} 880 for (n1, n2), label in labels.items(): 881 (x1, y1) = pos[n1] 882 (x2, y2) = pos[n2] 883 (x, y) = ( 884 x1 * label_pos + x2 * (1.0 - label_pos), 885 y1 * label_pos + y2 * (1.0 - label_pos), 886 ) 887 pos_1 = ax.transData.transform(np.array(pos[n1])) 888 pos_2 = ax.transData.transform(np.array(pos[n2])) 889 #linear_mid = 0.5*pos_1 + 0.5*pos_2 890 linear_mid = (1-label_pos)*pos_1 + label_pos*pos_2 891 d_pos = pos_2 - pos_1 892 rotation_matrix = np.array([(0,1), (-1,0)]) 893 ctrl_1 = linear_mid + rad*rotation_matrix@d_pos 894 ctrl_mid_1 = 0.5*pos_1 + 0.5*ctrl_1 895 ctrl_mid_2 = 0.5*pos_2 + 0.5*ctrl_1 896 bezier_mid = 0.5*ctrl_mid_1 + 0.5*ctrl_mid_2 897 (x, y) = ax.transData.inverted().transform(bezier_mid) 898 899 if rotate: 900 # in degrees 901 angle = np.arctan2(y2 - y1, x2 - x1) / (2.0 * np.pi) * 360 902 # make label orientation "right-side-up" 903 if angle > 90: 904 angle -= 180 905 if angle < -90: 906 angle += 180 907 # transform data coordinate angle to screen coordinate angle 908 xy = np.array((x, y)) 909 trans_angle = ax.transData.transform_angles( 910 np.array((angle,)), xy.reshape((1, 2)) 911 )[0] 912 else: 913 trans_angle = 0.0 914 # use default box of white with white border 915 if bbox is None: 916 bbox = dict(boxstyle="round", ec=(1.0, 1.0, 1.0), fc=(1.0, 1.0, 1.0)) 917 if not isinstance(label, str): 918 label = str(label) # this makes "1" and 1 labeled the same 919 920 t = ax.text( 921 x, y, label, 922 size=font_size, color=font_color, family=font_family, weight=font_weight, 923 alpha=alpha, 924 horizontalalignment=horizontalalignment, verticalalignment=verticalalignment, 925 rotation=trans_angle, transform=ax.transData, 926 bbox=bbox, zorder=1,clip_on=clip_on) 927 text_items[(n1, n2)] = t 928 929 ax.tick_params(axis="both", which="both",bottom=False, left=False, labelbottom=False, labelleft=False) 930 931 return text_items 932 933 934 def drawTransitionGraph(self, pos=None, 935 bended=False, node_size=1000, num=1, rad=-0.2, 936 edge_labels=None, node_shapes='o', figsize=(14, 7), 937 clear=True, fontsize=10, 938 fontcolor='black', label_pos=0.75): 939 """ 940 Visualize the state transition graph with nodes and directed edges. 941 942 This method draws the graph using `matplotlib` and `networkx`, including labels, 943 node colors, and optional curved (bended) edges for better visualization of 944 bidirectional transitions. It supports layout customization and is useful for 945 understanding the structure of the Markov model. 946 947 Parameters 948 ---------- 949 - pos : dict, optional 950 A dictionary mapping nodes to positions. If None, a Kamada-Kawai layout is used. 951 952 - bended : bool, optional 953 If True, edges are drawn with curvature using arcs. Useful for distinguishing 954 bidirectional transitions. Default is False. 955 956 - node_size : int, optional 957 Size of the nodes in the plot. Default is 1000. 958 959 - num : int, optional 960 Figure number for matplotlib (useful for managing multiple figures). Default is 1. 961 962 - rad : float, optional 963 The curvature radius for bended edges. Only used if `bended=True`. Default is -0.2. 964 965 - edge_labels : dict, optional 966 Optional dictionary mapping edges `(u, v)` to labels. If None, the `'label'` 967 attribute from each edge is used. 968 969 - node_shapes : str, optional 970 Shape of the nodes (e.g., `'o'` for circle, `'s'` for square). Default is `'o'`. 971 972 - figsize : tuple, optional 973 Size of the matplotlib figure. Default is (14, 7). 974 975 - clear : bool, optional 976 If True, clears the figure before plotting. Default is True. 977 978 - fontsize : int, optional 979 Font size used for node labels. Default is 10. 980 981 - fontcolor : str, optional 982 Font color for node labels. Default is `'black'`. 983 984 - label_pos : float, optional 985 Position of edge labels along the edge (0=start, 1=end). Default is 0.75. 986 987 Returns 988 ------- 989 - fig : matplotlib.figure.Figure 990 The created matplotlib figure object. 991 992 - ax : matplotlib.axes.Axes 993 The axes object where the graph is drawn. 994 995 - pos : dict 996 The dictionary of node positions used for drawing. 997 998 Notes 999 ----- 1000 - Node labels and colors must be set beforehand using `addState()` or `setStateProperties()`. 1001 - Edge attributes must include `'label'` and `'color'`. 1002 - This method modifies the current matplotlib figure and is intended for interactive or inline use. 1003 """ 1004 node_labels = {node: data["label"] for node, data in self.nodes(data=True)} 1005 node_colors = {node: data["color"] for node, data in self.nodes(data=True)} 1006 node_colors = [data["color"] for node, data in self.nodes(data=True)] 1007 #node_shapes = {(hw,sw): "o" if hw+sw<G.N else "s" for hw,sw in G.nodes()} 1008 1009 edge_cols = {(u,v): self[u][v]["color"] for u,v in self.edges} 1010 1011 if edge_labels is None: 1012 edge_labels = {(u,v): self[u][v]["label"] for u,v in self.edges} 1013 1014 if pos is None: 1015 pos = nx.kamada_kawai_layout(self) 1016 1017 plt.figure(figsize=figsize, num=num, clear=clear) 1018 nx.draw_networkx_nodes(self, pos, node_color=node_colors, node_shape = node_shapes, node_size=node_size) 1019 1020 nx.draw_networkx_labels(self, pos, labels=node_labels, font_size=fontsize, font_color=fontcolor) # Draw node labels 1021 if bended: 1022 nx.draw_networkx_edges(self, pos, width=1, edge_color=[edge_cols[edge] for edge in self.edges], 1023 node_size = node_size, 1024 arrows = True, arrowstyle = '-|>', 1025 connectionstyle=f"arc3,rad={rad}") 1026 self._draw_networkx_edge_labels(pos, ax=plt.gca(), edge_labels=edge_labels,rotate=False, rad = rad, label_pos=label_pos) 1027 else: 1028 nx.draw_networkx_edges(self, pos, width=1, edge_color=[edge_cols[edge] for edge in self.edges], 1029 node_size = node_size, 1030 #min_target_margin=17, arrowsize=15, 1031 arrows = True, arrowstyle = '-|>') 1032 nx.draw_networkx_edge_labels(self, pos, edge_labels, label_pos=label_pos) #, verticalalignment='center',) 1033 plt.axis('off'); 1034 1035 return plt.gcf(), plt.gca(), pos 1036 1037 def animateSimulation(self, expectedTimePerState=1, inNotebook=False, **kwargs): 1038 """ 1039 Animate a simulation trajectory through the state transition graph. 1040 The animation of the CTMC simulation run either in a Jupyter notebook or as a regular script. 1041 1042 This method visualizes the path of a simulated or predefined sequence of states by 1043 temporarily highlighting each visited state over time. The time spent in each state 1044 is scaled relative to the average sojourn time to produce a visually smooth animation. 1045 1046 Parameters 1047 ---------- 1048 - **expectedTimePerState** : float, optional 1049 Approximate duration (in seconds) for an average state visit in the animation. 1050 The actual pause duration for each state is scaled proportionally to its sojourn time. 1051 Default is 1 second. 1052 1053 - **inNotebook** : bool, optional 1054 If True, uses Jupyter-compatible animation (with display + clear_output). 1055 If False, uses standard matplotlib interactive animation. Default is True. 1056 1057 - `**kwargs` : dict 1058 Additional keyword arguments including: 1059 - states : list 1060 A list of visited states (nodes) to animate. 1061 - times : list or np.ndarray 1062 Sojourn times in each state (same length as `states`). 1063 - startNode : hashable 1064 Optional start node for automatic simulation if `states` and `times` are not provided. 1065 - numEvents : int 1066 Number of events (state transitions) to simulate. 1067 - Additional drawing-related kwargs are passed to `drawTransitionGraph()`. 1068 1069 Raises 1070 ------ 1071 ValueError 1072 If neither a `(states, times)` pair nor `(startNode, numEvents)` are provided. 1073 1074 Returns 1075 ------- 1076 None 1077 The animation is shown interactively using matplotlib but no data is returned. 1078 1079 Notes 1080 ----- 1081 - If `states` and `times` are not supplied, a simulation is run via `simulateMarkovModel()`. 1082 - The function uses matplotlib's interactive mode (`plt.ion()`) to animate transitions. 1083 - The average sojourn time across all states is used to normalize animation speed. 1084 - Each visited state is highlighted in red for its corresponding (scaled) dwell time. 1085 - Drawing arguments like layout, font size, or color can be customized via kwargs. 1086 """ 1087 # expTimePerState: in seconds 1088 if "states" in kwargs and "times" in kwargs: 1089 states = kwargs["states"] 1090 times = kwargs["times"] 1091 else: 1092 # Run default simulation if data not provided 1093 if "startNode" in kwargs and "numEvents" in kwargs: 1094 states, times = self.simulateMarkovModel(startNode=None, numEvents=10) 1095 else: 1096 raise ValueError("Must provide either ('states' and 'times') or ('startNode' and 'numEvents').") 1097 1098 #avg_OutGoingRate = np.mean([self.nodes[n].get("sumOutgoingRates", 0) for n in self.nodes]) 1099 avg_Time = np.mean([1.0/self.nodes[n].get("sumOutgoingRates", 0) for n in self.nodes]) 1100 1101 1102 # Remove simulation-related keys before drawing 1103 draw_kwargs = {k: v for k, v in kwargs.items() if k not in {"states", "times", "startNode", "numEvents"}} 1104 1105 fig, ax, pos = self.drawTransitionGraph(**draw_kwargs) 1106 plt.tight_layout() 1107 1108 if not inNotebook: 1109 plt.ion() 1110 1111 for node, time in zip(states,times): 1112 if inNotebook: 1113 clear_output(wait=True) 1114 artist = nx.draw_networkx_nodes(self, pos, nodelist=[node], ax=ax, 1115 node_color=self._const.COLOR_NODE_ANIMATION_HIGHLIGHT, node_size=1000) 1116 display(fig) 1117 else: 1118 artist = nx.draw_networkx_nodes(self, pos, nodelist=[node],ax=ax, 1119 node_color=self._const.COLOR_NODE_ANIMATION_HIGHLIGHT, node_size=1000) 1120 plt.draw() 1121 plt.pause(time/avg_Time*expectedTimePerState) 1122 1123 # Remove highlighted node 1124 artist.remove() 1125 1126 pass 1127 1128 1129 1130 def states(self, data=False): 1131 """ 1132 Return a view of the graph's states (nodes), optionally with attributes. 1133 1134 This method is functionally identical to `self.nodes()` in NetworkX, but is 1135 renamed to `states()` to reflect the semantics of a state transition graph or 1136 Markov model, where nodes represent system states. 1137 1138 Parameters 1139 ---------- 1140 **data** : bool, optional 1141 If True, returns a view of `(node, attribute_dict)` pairs. 1142 If False (default), returns a view of node identifiers only. 1143 1144 Returns 1145 ------- 1146 *networkx.classes.reportviews.NodeView* or NodeDataView 1147 A view over the graph’s nodes or `(node, data)` pairs, depending on the `data` flag. 1148 1149 Examples 1150 -------- 1151 >>> G.states() 1152 ['A', 'B', 'C'] 1153 1154 >>> list(G.states(data=True)) 1155 [('A', {'color': 'red', 'label': 'Start'}), ('B', {...}), ...] 1156 1157 Notes 1158 ----- 1159 - This is a convenience wrapper for semantic clarity in models where nodes represent states. 1160 - The view is dynamic: changes to the graph are reflected in the returned view. 1161 """ 1162 return self.nodes(data=data) 1163 1164 1165 def setStateColor(self, state, color): 1166 """ 1167 Set the color attribute of a specific state (node). 1168 1169 Parameters 1170 ---------- 1171 **state** : hashable 1172 The identifier of the state whose color is to be updated. 1173 1174 **color** : str 1175 The new color to assign to the state (used for visualization). 1176 1177 Returns 1178 ------- 1179 None 1180 """ 1181 self.nodes[state]["color"] = color 1182 1183 1184 def setStateLabel(self, state, label): 1185 """ 1186 Set the label attribute of a specific state (node). 1187 1188 Parameters 1189 ---------- 1190 **state** : hashable 1191 The identifier of the state whose label is to be updated. 1192 1193 **label** : str 1194 The new label to assign to the state (used for visualization or annotation). 1195 1196 Returns 1197 ------- 1198 None 1199 """ 1200 self.nodes[state]["label"] = label 1201 1202 def prob(self, state): 1203 """ 1204 Return the steady-state probability of a specific state. 1205 1206 Parameters 1207 ---------- 1208 **state** : hashable 1209 The identifier of the state whose probability is to be retrieved. 1210 1211 Returns 1212 ------- 1213 **float** 1214 The steady-state probability of the specified state. 1215 1216 Raises 1217 ------ 1218 KeyError 1219 If the state is not found in the steady-state probability dictionary. 1220 """ 1221 return self.state_probabilities[state] 1222 1223 1224 def probs(self): 1225 """ 1226 Return the steady-state probabilities for all states. 1227 1228 Returns 1229 ------- 1230 **dict** 1231 A dictionary mapping each state to its steady-state probability. 1232 1233 Notes 1234 ----- 1235 - The probabilities are computed via `calcSteadyStateProbs()` and stored internally. 1236 - Use this method to access the full steady-state distribution. 1237 """ 1238 return self.state_probabilities 1239 1240 def getEmptySystemProbs(self, state=None): 1241 """ 1242 Create an initial state distribution where all probability mass is in a single state. 1243 1244 By default, this method returns a distribution where the entire probability mass 1245 is placed on the first node in the graph. Alternatively, a specific starting state 1246 can be provided. 1247 1248 Parameters 1249 ---------- 1250 **state** : hashable, optional 1251 The state to initialize with probability 1. If None, the first node 1252 in the graph (based on insertion order) is used. 1253 1254 Returns 1255 ------- 1256 **dict** 1257 A dictionary representing the initial distribution over states. 1258 The selected state has probability 1, and all others have 0. 1259 1260 Examples 1261 -------- 1262 >>> X0 = G.getEmptySystem() 1263 >>> sum(X0.values()) # 1.0 1264 1265 Notes 1266 ----- 1267 - This method is useful for setting up deterministic initial conditions 1268 in transient simulations of a CTMC. 1269 - The return format is compatible with methods like `transientProbs()`. 1270 """ 1271 if state is None: 1272 first_node = next(iter(self.nodes)) 1273 X0 = {s: 1 if s == first_node else 0 for s in self} 1274 return X0 1275 1276 1277 1278 def transientProbs(self, initialDistribution, t): 1279 """ 1280 Compute the transient state probabilities at time `t` for the CTMC. 1281 1282 This method solves the forward equation of a continuous-time Markov chain (CTMC): 1283 X(t) = X0 · expm(Q·t), where `Q` is the generator (rate) matrix and `X0` is the 1284 initial state distribution. It returns a dictionary mapping each state to its 1285 probability at time `t`. 1286 1287 Parameters 1288 ---------- 1289 - **initialDistribution** : dict 1290 A dictionary mapping states (nodes) to their initial probabilities at time `t=0`. 1291 The keys must match the nodes in the graph, and the values should sum to 1. 1292 1293 - **t** : float 1294 The time at which the transient distribution is to be evaluated. 1295 1296 Returns 1297 ------- 1298 **dict** 1299 A dictionary mapping each state (node) to its probability at time `t`. 1300 1301 Notes 1302 ----- 1303 - The rate matrix `Q` is generated automatically if not yet computed. 1304 - The state order in the rate matrix corresponds to the internal `_n2i` mapping. 1305 - Internally, the dictionary `initialDistribution` is converted into a vector 1306 aligned with the state indexing. 1307 - The transient solution is computed using `scipy.linalg.expm` for matrix exponentiation. 1308 1309 Raises 1310 ------ 1311 ValueError 1312 If the input dictionary has missing states or probabilities that do not sum to 1. 1313 (This is not enforced but recommended for correctness.) 1314 """ 1315 if self._rate_matrix is None: self.createRateMatrix() 1316 Q = self._rate_matrix 1317 X0 = np.zeros(len(initialDistribution)) 1318 1319 for i,n in enumerate(self.nodes): 1320 X0[ self._n2i[n] ] = initialDistribution[n] 1321 1322 1323 X = X0 @ expm(Q*t) 1324 1325 matrix = {} 1326 for n in self.nodes: 1327 matrix[n] = X[ self._n2i[n] ] 1328 1329 return matrix 1330 1331 def symSolveMarkovModel(self): 1332 """ 1333 Symbolically solve the global balance equations of the CTMC. 1334 1335 This method constructs and solves the system of symbolic equations for the 1336 steady-state probabilities of each state in the CTMC. It uses symbolic transition 1337 rates stored on the edges to form balance equations for each node, along with 1338 the normalization constraint that all probabilities sum to 1. 1339 1340 Returns 1341 ------- 1342 - **solution** : dict 1343 A symbolic solution dictionary mapping SymPy symbols (e.g., x_{A}) to 1344 expressions in terms of model parameters. 1345 1346 - **num_solution** : dict 1347 A numerical dictionary mapping each state (node) to its steady-state probability, 1348 computed by substituting the numeric values from `self._subs` into the symbolic solution. 1349 1350 Notes 1351 ----- 1352 - For each node `s`, a symbolic variable `x_{label}` is created using the node's label attribute. 1353 - One balance equation is created per state: total inflow = total outflow. 1354 - An additional constraint `sum(x_i) = 1` ensures proper normalization. 1355 - The symbolic system is solved with `sympy.solve`, and the results are simplified. 1356 - Numeric values are computed by substituting known parameter values (`self._subs`) into the symbolic solution. 1357 1358 Examples 1359 -------- 1360 >>> sym_sol, num_sol = G.symSolveMarkovModel() 1361 >>> sym_sol # symbolic expressions for each state 1362 >>> num_sol # evaluated numerical values for each state 1363 """ 1364 X={} 1365 to_solve = [] 1366 for s in self.nodes(): 1367 X[s] = sp.Symbol(f'x_{{{self.nodes[s]["label"]}}}') 1368 to_solve.append(X[s]) 1369 #% 1370 eqs = [] 1371 for s in self.nodes(): 1372 succ = self.successors(s) # raus 1373 out_rate = 0 1374 for z in list(succ): 1375 #out_rate += X[s]*self.edges[s, z]["sym_rate"] 1376 out_rate += X[s]*self[s][z]["sym_rate"] 1377 1378 in_rate = 0 1379 for r in list(self.predecessors(s)): 1380 #in_rate += X[r]*self.edges[r,s]["sym_rate"] 1381 in_rate += X[r]*self[r][s]["sym_rate"] 1382 eqs.append( sp.Eq(out_rate, in_rate) ) 1383 #% 1384 Xsum = 0 1385 for s in X: 1386 Xsum += X[s] 1387 1388 eqs.append( sp.Eq(Xsum, 1) ) 1389 #% 1390 solution = sp.solve(eqs, to_solve) 1391 simplified_solution = sp.simplify(solution) 1392 #print(simplified_solution) 1393 num_solution = simplified_solution.subs(self._subs) 1394 #num_dict = {str(var): str(sol) for var, sol in simplified_solution.items()} 1395 num_dict = {s: num_solution[X[s]] for s in X} 1396 1397 return simplified_solution, num_dict 1398 1399 def exportTransitionsToExcel(self, excel_path = "output.xlsx", num_colors=9, open_excel=False, verbose=True): 1400 """ 1401 Export transition data from the state transition graph to an Excel file. 1402 1403 This method collects all transitions (edges) from the graph and writes them 1404 to an Excel spreadsheet, optionally applying background colors to distinguish 1405 outgoing transitions from different source states. This export is helpful for 1406 manually verifying the transition structure of the model. 1407 1408 The transition type is defined by the `tid` field of each edge, and mapped to 1409 a human-readable string via the `EXPORT_NAME_TRANSITIONS` dictionary in the 1410 constants module (e.g., `Constants.EXPORT_NAME_TRANSITIONS[tid]`). 1411 1412 See Also 1413 -------- 1414 StateTransitionGraph.addTransition : Method to add edges with symbolic rates and `tid`. 1415 StateTransitionGraph.__init__ : Initializes the graph and loads the constants module. 1416 1417 Parameters 1418 ---------- 1419 - **excel_path** : str, optional 1420 Path to the output Excel file. Default is "output.xlsx". 1421 1422 - num_colors : int, optional 1423 Number of unique background colors to cycle through for different states. 1424 Default is 9. 1425 1426 - open_excel : bool, optional 1427 If True, automatically opens the generated Excel file after export. 1428 Default is False. 1429 1430 - verbose : bool, optional 1431 If True, prints summary information about the export process. 1432 Default is True. 1433 1434 Returns 1435 ------- 1436 None 1437 The function writes data to disk and optionally opens Excel, 1438 but does not return a value. 1439 1440 Notes 1441 ----- 1442 - Each edge must contain a `tid` attribute for transition type and a `label` for the rate. 1443 - Background coloring groups all outgoing transitions from the same source state. 1444 - Requires the constants module to define `EXPORT_NAME_TRANSITIONS`. 1445 """ 1446 1447 rows = [] 1448 for node in self.nodes(): 1449 edges = list(self.edges(node)) 1450 if verbose: print(f"Edges connected to node '{node}':") 1451 1452 for i, (u, v) in enumerate(edges): 1453 latex = self[u][v]["label"] 1454 rate = sp.pretty(self[u][v]["sym_rate"]) 1455 tid = self[u][v]["tid"] 1456 1457 if tid: # 1458 if hasattr(self._const, 'EXPORT_NAME_TRANSITIONS'): 1459 tid_name = self._const.EXPORT_NAME_TRANSITIONS[tid] 1460 else: 1461 tid_name = '' 1462 1463 rows.append({"from": u, "type": tid_name, "to": v, "rate": rate, "latex": latex}) 1464 1465 df = pd.DataFrame(rows, columns=["from", "type", "to", "rate", "latex"]) 1466 df["checked"] = "" 1467 df.to_excel(excel_path, index=False) 1468 1469 # Load workbook with openpyxl 1470 wb = load_workbook(excel_path) 1471 ws = wb.active 1472 1473 # Define colors for different "from" values 1474 cmap = map(plt.cm.Pastel1, range(num_colors)) 1475 # Convert RGBA to hex 1476 color_list = [mcolors.to_hex(color, keep_alpha=False).upper().replace("#", "") for color in cmap] 1477 1478 # Assign color to each unique "from" value using modulo 1479 unique_from = (df["from"].unique()) # Sorted for consistent order 1480 from_color_map = { 1481 value: color_list[i % len(color_list)] for i, value in enumerate(unique_from) 1482 } 1483 1484 # Apply coloring row by row 1485 for i,row in enumerate(ws.iter_rows(min_row=2, max_row=ws.max_row)): # Skip header 1486 from_value = row[0].value # "from" is first column 1487 #fill_color = color_list[i % len(color_list)] 1488 if isinstance(from_value, int): 1489 fill_color = from_color_map[from_value] 1490 else: 1491 fill_color = from_color_map.get(eval(from_value)) 1492 1493 if fill_color: 1494 fill = PatternFill(start_color=fill_color, end_color=fill_color, fill_type="solid") 1495 for cell in row: 1496 cell.fill = fill 1497 1498 # Save styled Excel file 1499 wb.save(excel_path) 1500 if verbose: print(f"File {excel_path} created.") 1501 if open_excel: 1502 print(f"Opening the file {excel_path}.") 1503 os.startfile(excel_path)
Initialize the StateTransitionGraph for continuous time Markov chains (CTMC). Extends the networkx.Graph class to include additional methods and properties for the analysis, visualization, and simulation of CTMCs. The analysis allows deriving the state probabilities in the steady-state and in the transient phase.
This constructor sets up symbolic and numeric parameters for the state transition graph, imports a constants module dynamically, and initializes internal data structures for rate matrix computation and steady-state probability analysis.
Parameters
parameters : dict
A dictionary mapping parameter names to their numerical values.
Example: {'lambda': 1.0, 'mu': 1.0}. For each parameter key:- A symbolic variable is created using sympy.symbols. Can be accessed as
G.sym("key")
orG.sym_key
- A variable is created which can be accessed via
G.key
(and the numerical value is assigned) For the key="lambda",
G.lam
andG.sym_lambda
are createdThe parameters can be accessed via property
parameters
- A symbolic variable is created using sympy.symbols. Can be accessed as
- constants : str, optional
The name of a Python module (as string) to import dynamically and assign to
self.const
. This can be used to store global or user-defined constants. Default is "Constants". *args
: tuple Additional positional arguments passed to the basenetworkx.DiGraph
constructor.**kwargs
: dict Additional keyword arguments passed to the basenetworkx.DiGraph
constructor.
Attributes
_rate_matrix
: ndarray or None
The continuous-time rate matrixQ
of the CTMC. Computed when callingcreateRateMatrix()
._state_probabilities
: dict or None
Dictionary of steady-state probabilities, computed viacalcSteadyStateProbs()
._n2i
: dict or None
Mapping from node identifiers to matrix row/column indices used for rate matrix construction._sym
: dict
Mapping of parameter names (strings) to their symbolicsympy
representations._subs
: dict
Mapping from symbolic parameters to their numerical values (used for substitutions)._parameters
: dict
Copy of the original parameter dictionary passed at initialization._transitionProbabilitesComputed
: bool
Internal flag indicating whether transition probabilities for the embedded Markov chain have been computed.const
: module
Dynamically imported module (default"Constants"
), used for storing user-defined constants such as colors or transition IDs.
Inherits
networkx.Graph: The base graph class from NetworkX.
153 def __init__(self, parameters, constants="Constants", *args, **kwargs): 154 """ 155 Initialize the state transition graph with symbolic parameters. 156 157 This constructor sets up the symbolic and numerical environment for modeling 158 a CTMC. Parameters are stored as both symbolic expressions and numerical values, 159 and convenient attribute access is provided for both forms. 160 161 Parameters 162 ---------- 163 - **parameters** : dict 164 A dictionary mapping parameter names to their numerical values. 165 Example: {'lambda': 1.0, 'mu': 1.0}. 166 For each parameter key: 167 - A symbolic variable is created using sympy.symbols. Can be accessed as `G.sym("key")` or `G.sym_key` 168 - A variable is created which can be accessed via `G.key` (and the numerical value is assigned) 169 - For the key="lambda", `G.lam` and `G.sym_lambda` are created 170 171 The parameters can be accessed via property `parameters` 172 - **constants** : str, optional 173 The name of a Python module (as string) to import dynamically and assign to `self.const`. 174 This can be used to store global or user-defined constants. 175 Default is "Constants". 176 - `*args` : tuple 177 Additional positional arguments passed to the base `networkx.DiGraph` constructor. 178 - `**kwargs` : dict 179 Additional keyword arguments passed to the base `networkx.DiGraph` constructor. 180 """ 181 super().__init__(*args, **kwargs) 182 self._rate_matrix = None # Rate Matrix Q of this transition graph 183 self._state_probabilities = None # steady state proabilities of this CTMC 184 self._n2i = None # internal index for generating the rate matrix Q: dictionary (node to index) 185 186 #parameters = {'lambda': 1.0, 'mu':1.0, 'n':5} 187 self._parameters = parameters 188 self._sym = {key: sp.symbols(key) for key in parameters} 189 self._subs = {} # contains symbolic parameter as key and provides numerical value for computation. 190 for k in parameters: 191 self._subs[ self._sym[k] ] = parameters[k] 192 193 # access to parameters and symbols using dot-notation, e.g. G.lam, G.sym_lambda, or G.sym("lambda") 194 for key, value in parameters.items(): 195 if key=='lambda': 196 setattr(self, 'lam', value) 197 else: 198 setattr(self, key, value) 199 200 setattr(self, "sym_"+key, self._sym[key]) 201 202 self._transitionProbabilitesComputed = False 203 204 self._const = importlib.import_module(constants)
Initialize the state transition graph with symbolic parameters.
This constructor sets up the symbolic and numerical environment for modeling a CTMC. Parameters are stored as both symbolic expressions and numerical values, and convenient attribute access is provided for both forms.
Parameters
parameters : dict
A dictionary mapping parameter names to their numerical values.
Example: {'lambda': 1.0, 'mu': 1.0}. For each parameter key:- A symbolic variable is created using sympy.symbols. Can be accessed as
G.sym("key")
orG.sym_key
- A variable is created which can be accessed via
G.key
(and the numerical value is assigned) For the key="lambda",
G.lam
andG.sym_lambda
are createdThe parameters can be accessed via property
parameters
- A symbolic variable is created using sympy.symbols. Can be accessed as
- constants : str, optional
The name of a Python module (as string) to import dynamically and assign to
self.const
. This can be used to store global or user-defined constants. Default is "Constants". *args
: tuple Additional positional arguments passed to the basenetworkx.DiGraph
constructor.**kwargs
: dict Additional keyword arguments passed to the basenetworkx.DiGraph
constructor.
208 def addState(self, state, color=None, label=None): 209 """ 210 Add a state (node) to the graph with optional color and label attributes. 211 212 If the node already exists, its attributes are updated with the provided values. 213 214 Parameters 215 ---------- 216 - **state** : hashable 217 The identifier for the state. Can be a string, int, or a tuple. If no label 218 is provided, the label will be auto-generated as: 219 - comma-separated string for tuples (e.g., `(1, 2)` -> `'1,2'`) 220 - stringified value for other types 221 222 - **color** : str, optional 223 The color assigned to the state (used for visualization). 224 If not specified, defaults to `self.const.COLOR_NODE_DEFAULT`. 225 226 - **label** : str, optional 227 A label for the node (used for display or annotation). If not provided, 228 it will be auto-generated from the state identifier. 229 230 Notes 231 ----- 232 - This method wraps `self.add_node()` from NetworkX. 233 - If the node already exists in the graph, this method does not raise an error; 234 it will update the node’s existing attributes. 235 """ 236 if label is None: 237 if isinstance(state, tuple): 238 label = ','.join(map(str, state)) # if state is a tuple 239 else: 240 label = str(state) 241 #self.nodes[state]["label"] = label 242 243 if color is None: 244 color = self._const.COLOR_NODE_DEFAULT 245 246 # no error if the node already exists.Attributes get updated if you pass them again. 247 self.add_node(state, color=color, label=label)
Add a state (node) to the graph with optional color and label attributes.
If the node already exists, its attributes are updated with the provided values.
Parameters
state : hashable The identifier for the state. Can be a string, int, or a tuple. If no label is provided, the label will be auto-generated as:
- comma-separated string for tuples (e.g.,
(1, 2)
->'1,2'
) - stringified value for other types
- comma-separated string for tuples (e.g.,
color : str, optional The color assigned to the state (used for visualization). If not specified, defaults to
self.const.COLOR_NODE_DEFAULT
.label : str, optional A label for the node (used for display or annotation). If not provided, it will be auto-generated from the state identifier.
Notes
- This method wraps
self.add_node()
from NetworkX. - If the node already exists in the graph, this method does not raise an error; it will update the node’s existing attributes.
249 def setStateProperties(self, state, color=None, label=None, **kwargs): 250 """ 251 Update visual or metadata properties for an existing state (node) in the graph. 252 253 This method updates the node's attributes such as color, label, and any additional 254 custom key-value pairs. The node must already exist in the graph. 255 256 Parameters 257 ---------- 258 - **state** : hashable 259 Identifier of the state to update. Must be a node already present in the graph. 260 261 - **color** : str, optional 262 Color assigned to the node, typically used for visualization. If not provided, 263 the default node color (`self._const.COLOR_NODE_DEFAULT`) is used. 264 265 - **label** : str, optional 266 Label for the node. If not provided, it is auto-generated from the state ID 267 (especially for tuples). 268 269 - `**kwargs` : dict 270 Arbitrary additional key-value pairs to set as node attributes. These can include 271 any custom metadata. 272 273 Raises 274 ------ 275 KeyError 276 If the specified state is not present in the graph. 277 278 Notes 279 ----- 280 - This method ensures that required attributes (like `color` and `label`) are set for visualization. 281 - Existing node attributes are updated or extended with `kwargs`. 282 """ 283 if state not in self: 284 raise KeyError(f"State '{state}' does not exist in the graph.") 285 286 # mandatory for plotting 287 self.addState(state, color, label) 288 289 # # Add any additional attributes 290 self.nodes[state].update(kwargs)
Update visual or metadata properties for an existing state (node) in the graph.
This method updates the node's attributes such as color, label, and any additional custom key-value pairs. The node must already exist in the graph.
Parameters
state : hashable Identifier of the state to update. Must be a node already present in the graph.
color : str, optional Color assigned to the node, typically used for visualization. If not provided, the default node color (
self._const.COLOR_NODE_DEFAULT
) is used.label : str, optional Label for the node. If not provided, it is auto-generated from the state ID (especially for tuples).
**kwargs
: dict Arbitrary additional key-value pairs to set as node attributes. These can include any custom metadata.
Raises
KeyError If the specified state is not present in the graph.
Notes
- This method ensures that required attributes (like
color
andlabel
) are set for visualization. - Existing node attributes are updated or extended with
kwargs
.
292 def setAllStateDefaultProperties(self, **kwargs): 293 """ 294 Set default visual and metadata properties for all states (nodes) in the graph. 295 296 This method applies `setStateProperties()` to every node, ensuring each state 297 has at least the required default attributes (e.g., color and label). 298 Any additional keyword arguments are passed to `setStateProperties()`. 299 300 Parameters 301 ---------- 302 `**kwargs` : dict, optional 303 Optional keyword arguments to apply uniformly to all nodes, 304 such as color='gray', group='core', etc. 305 306 Notes 307 ----- 308 - This method is useful for initializing node attributes before plotting 309 or performing other graph-based computations. 310 - Nodes without existing attributes will receive defaults via `addState()`. 311 """ 312 for state in self.nodes(): 313 self.setStateProperties(state, **kwargs)
Set default visual and metadata properties for all states (nodes) in the graph.
This method applies setStateProperties()
to every node, ensuring each state
has at least the required default attributes (e.g., color and label).
Any additional keyword arguments are passed to setStateProperties()
.
Parameters
**kwargs
: dict, optional
Optional keyword arguments to apply uniformly to all nodes,
such as color='gray', group='core', etc.
Notes
- This method is useful for initializing node attributes before plotting or performing other graph-based computations.
- Nodes without existing attributes will receive defaults via
addState()
.
317 def addTransition(self, origin_state, destination_state, sym_rate, tid = None, color=None, label=None): 318 """ 319 Add a directed transition (edge) between two states with a symbolic transition rate. 320 321 This method adds a directed edge from `origin_state` to `destination_state`, 322 stores both the symbolic and numerical rate, and sets optional display properties 323 such as color, label, and transition ID. If the edge already exists, a warning is issued. 324 325 Parameters 326 ---------- 327 - **origin_state** : hashable 328 The source node (state) of the transition. 329 330 - **destination_state** : hashable 331 The target node (state) of the transition. 332 333 - **sym_rate** : sympy.Expr 334 A symbolic expression representing the transition rate. 335 This will be numerically evaluated using internal substitutions (`self._subs`). 336 337 - **tid** : str or int, optional 338 A transition identifier (e.g. for grouping or styling). 339 If provided, the color will default to `self._const.COLOR_TRANSITIONS[tid]` if not set explicitly. 340 341 - **color** : str, optional 342 Color for the edge (used in visualization). Defaults to the transition-specific color 343 if `tid` is provided. 344 345 - **label** : str, optional 346 LaTeX-formatted string for labeling the edge. If not provided, a default label is 347 generated using `sympy.latex(sym_rate)`. 348 349 Notes 350 ----- 351 - The numeric rate is computed by substituting known parameter values into `sym_rate`. 352 - The following attributes are stored on the edge: 353 - `rate` (float): numeric value of the transition rate 354 - `sym_rate` (sympy.Expr): original symbolic expression 355 - `label` (str): display label 356 - `color` (str): edge color 357 - `tid` (str or int): transition ID 358 - If the edge already exists, a `UserWarning` is issued but the edge is still overwritten. 359 """ 360 if tid is not None: 361 if color is None: 362 color = self._const.COLOR_TRANSITIONS[ tid ] 363 364 rate = sym_rate.subs( self._subs ) 365 366 if label is None: string = "$"+sp.latex(sym_rate)+"$" 367 368 if self.has_edge(origin_state, destination_state): 369 theedge = self[origin_state][destination_state] 370 warningstr = f"Edge already exists: {origin_state} -> {destination_state} with rate {theedge['label']}" 371 warnings.warn(warningstr, category=UserWarning) 372 self.add_edge(origin_state, destination_state, rate=float(rate), label=string, 373 color=color, sym_rate=sym_rate, tid=tid)
Add a directed transition (edge) between two states with a symbolic transition rate.
This method adds a directed edge from origin_state
to destination_state
,
stores both the symbolic and numerical rate, and sets optional display properties
such as color, label, and transition ID. If the edge already exists, a warning is issued.
Parameters
origin_state : hashable The source node (state) of the transition.
destination_state : hashable The target node (state) of the transition.
sym_rate : sympy.Expr A symbolic expression representing the transition rate. This will be numerically evaluated using internal substitutions (
self._subs
).tid : str or int, optional A transition identifier (e.g. for grouping or styling). If provided, the color will default to
self._const.COLOR_TRANSITIONS[tid]
if not set explicitly.color : str, optional Color for the edge (used in visualization). Defaults to the transition-specific color if
tid
is provided.label : str, optional LaTeX-formatted string for labeling the edge. If not provided, a default label is generated using
sympy.latex(sym_rate)
.
Notes
- The numeric rate is computed by substituting known parameter values into
sym_rate
. - The following attributes are stored on the edge:
rate
(float): numeric value of the transition ratesym_rate
(sympy.Expr): original symbolic expressionlabel
(str): display labelcolor
(str): edge colortid
(str or int): transition ID
- If the edge already exists, a
UserWarning
is issued but the edge is still overwritten.
375 def sym(self, key): 376 """ 377 Retrieve the symbolic representation of a model parameter. 378 379 Parameters 380 ---------- 381 **key** : str 382 The name of the parameter (e.g., 'lambda', 'mu', 'n'). 383 384 Returns 385 ------- 386 **sympy.Symbol** 387 The symbolic variable corresponding to the given key. 388 389 Raises 390 ------ 391 KeyError 392 If the key does not exist in the internal symbolic mapping (`_sym`). 393 394 Examples 395 -------- 396 >>> G.sym('lambda') 397 lambda 398 399 Notes 400 ----- 401 - The symbolic mapping is initialized during object construction based on the `parameters` dictionary. 402 - To access the symbol using dot-notation, use `.sym_key` attributes created during initialization. 403 """ 404 return self._sym[key]
Retrieve the symbolic representation of a model parameter.
Parameters
key : str The name of the parameter (e.g., 'lambda', 'mu', 'n').
Returns
sympy.Symbol The symbolic variable corresponding to the given key.
Raises
KeyError
If the key does not exist in the internal symbolic mapping (_sym
).
Examples
>>> G.sym('lambda')
lambda
Notes
- The symbolic mapping is initialized during object construction based on the
parameters
dictionary. - To access the symbol using dot-notation, use
.sym_key
attributes created during initialization.
407 @property 408 def parameters(self): 409 """ 410 Return the model parameters as a dictionary. 411 412 This includes the original numerical values provided during initialization. 413 414 Returns 415 ------- 416 **dict** 417 Dictionary of parameter names and their corresponding numeric values. 418 419 Examples 420 -------- 421 >>> G.parameters 422 {'lambda': 1.0, 'mu': 0.5, 'n': 3} 423 424 Notes 425 ----- 426 - This property does not return symbolic values. Use `self.sym(key)` for symbolic access. 427 - Parameter symbols can also be accessed using attributes like `sym_lambda` or the `sym()` method. 428 """ 429 return self._parameters
Return the model parameters as a dictionary.
This includes the original numerical values provided during initialization.
Returns
dict Dictionary of parameter names and their corresponding numeric values.
Examples
>>> G.parameters
{'lambda': 1.0, 'mu': 0.5, 'n': 3}
Notes
- This property does not return symbolic values. Use
self.sym(key)
for symbolic access. - Parameter symbols can also be accessed using attributes like
sym_lambda
or thesym()
method.
433 @property 434 def rate_matrix(self): 435 """ 436 Return the continuous-time rate matrix (Q matrix) of the state transition graph. 437 438 This matrix defines the transition rates between all pairs of states in the 439 continuous-time Markov chain (CTMC) represented by this graph. 440 441 Returns 442 ------- 443 **numpy.ndarray** or None 444 The rate matrix `Q`, where `Q[i, j]` is the transition rate from state `i` to `j`. 445 Returns `None` if the rate matrix has not been computed yet. 446 447 Notes 448 ----- 449 - The matrix is typically generated via a method like `computeRateMatrix()` or similar. 450 - Internal state indexing is handled by `_n2i` (node-to-index mapping). 451 - The rate matrix is stored internally as `_rate_matrix`. 452 """ 453 return self._rate_matrix
Return the continuous-time rate matrix (Q matrix) of the state transition graph.
This matrix defines the transition rates between all pairs of states in the continuous-time Markov chain (CTMC) represented by this graph.
Returns
numpy.ndarray or None
The rate matrix Q
, where Q[i, j]
is the transition rate from state i
to j
.
Returns None
if the rate matrix has not been computed yet.
Notes
- The matrix is typically generated via a method like
computeRateMatrix()
or similar. - Internal state indexing is handled by
_n2i
(node-to-index mapping). - The rate matrix is stored internally as
_rate_matrix
.
457 @property 458 def state_probabilities(self): 459 """ 460 Return the steady-state probabilities of the continuous-time Markov chain (CTMC). 461 462 These probabilities represent the long-run proportion of time the system spends 463 in each state, assuming the CTMC is irreducible and positive recurrent. 464 465 Returns 466 ------- 467 **dict** or None 468 A dictionary mapping each state to its steady-state probability. 469 Returns `None` if the probabilities have not yet been computed. 470 471 Notes 472 ----- 473 - The probabilities are stored internally after being computed via a dedicated method 474 `calcSteadyStateProbs()`. 475 - The result is stored in the `_state_probabilities` attribute. 476 """ 477 return self._state_probabilities
Return the steady-state probabilities of the continuous-time Markov chain (CTMC).
These probabilities represent the long-run proportion of time the system spends in each state, assuming the CTMC is irreducible and positive recurrent.
Returns
dict or None
A dictionary mapping each state to its steady-state probability.
Returns None
if the probabilities have not yet been computed.
Notes
- The probabilities are stored internally after being computed via a dedicated method
calcSteadyStateProbs()
. - The result is stored in the
_state_probabilities
attribute.
481 def printEdges(self, state): 482 """ 483 Print all outgoing edges from a given state in the graph. 484 485 For each edge from the specified `state`, this method prints the destination 486 state and the corresponding transition label (typically the symbolic rate). 487 488 Parameters 489 ---------- 490 **state** : hashable 491 The source node (state) from which outgoing edges will be printed. 492 493 Returns 494 ------- 495 None 496 This method prints output directly to the console and does not return a value. 497 498 Notes 499 ----- 500 - The transition label is retrieved from the edge attribute `'label'`. 501 - Assumes `state` exists in the graph; otherwise, a KeyError will be raised. 502 """ 503 print(f"from {state}") 504 node = self[state] 505 for e in node: 506 rate = node[e]['label'] 507 print(f" -> {e} with {rate}")
Print all outgoing edges from a given state in the graph.
For each edge from the specified state
, this method prints the destination
state and the corresponding transition label (typically the symbolic rate).
Parameters
state : hashable The source node (state) from which outgoing edges will be printed.
Returns
None This method prints output directly to the console and does not return a value.
Notes
- The transition label is retrieved from the edge attribute
'label'
. - Assumes
state
exists in the graph; otherwise, a KeyError will be raised.
510 def createRateMatrix(self): 511 """ 512 Construct the continuous-time Markov chain (CTMC) rate matrix Q from the graph. 513 514 This method generates the generator matrix `Q` corresponding to the structure and 515 rates defined in the directed graph. Each node represents a system state, and each 516 directed edge must have a `'rate'` attribute representing the transition rate from 517 one state to another. 518 519 The matrix `Q` is constructed such that: 520 - `Q[i, j]` is the transition rate from state `i` to state `j` 521 - Diagonal entries satisfy `Q[i, i] = -sum(Q[i, j]) for j ≠ i`, so that each row sums to 0 522 523 The resulting matrix and the internal node-to-index mapping (`_n2i`) are stored 524 as attributes of the graph. 525 526 Returns 527 ------- 528 None 529 The method updates internal attributes: 530 - `_rate_matrix` : numpy.ndarray 531 The CTMC rate matrix Q 532 - `_n2i` : dict 533 Mapping from node identifiers to matrix row/column indices 534 535 Notes 536 ----- 537 - Assumes all edges in the graph have a `'rate'` attribute. 538 - The node order used in the matrix is consistent with the order from `list(self.nodes)`. 539 - Use `self.rate_matrix` to access the generated matrix after this method is called. 540 """ 541 n2i = {} 542 nodes = list(self.nodes) 543 for i,node in enumerate(nodes): 544 n2i[node] = i 545 Q = np.zeros((len(nodes),len(nodes))) # rate matrix 546 547 for edge in self.edges: 548 i0 = n2i[edge[0]] 549 i1 = n2i[edge[1]] 550 Q[i0,i1] = self[edge[0]][edge[1]]["rate"] 551 552 np.fill_diagonal(Q, -Q.sum(axis=1)) 553 self._rate_matrix = Q 554 self._n2i = n2i 555 #return Q, n2i
Construct the continuous-time Markov chain (CTMC) rate matrix Q from the graph.
This method generates the generator matrix Q
corresponding to the structure and
rates defined in the directed graph. Each node represents a system state, and each
directed edge must have a 'rate'
attribute representing the transition rate from
one state to another.
The matrix Q
is constructed such that:
Q[i, j]
is the transition rate from statei
to statej
- Diagonal entries satisfy
Q[i, i] = -sum(Q[i, j]) for j ≠ i
, so that each row sums to 0
The resulting matrix and the internal node-to-index mapping (_n2i
) are stored
as attributes of the graph.
Returns
None
The method updates internal attributes:
- _rate_matrix
: numpy.ndarray
The CTMC rate matrix Q
- _n2i
: dict
Mapping from node identifiers to matrix row/column indices
Notes
- Assumes all edges in the graph have a
'rate'
attribute. - The node order used in the matrix is consistent with the order from
list(self.nodes)
. - Use
self.rate_matrix
to access the generated matrix after this method is called.
557 def calcSteadyStateProbs(self, verbose=False): 558 """ 559 Compute the steady-state probabilities of the continuous-time Markov chain (CTMC). 560 561 This method calculates the stationary distribution of the CTMC defined by the rate 562 matrix of the state transition graph. The computation solves the linear system 563 `XQ = 0` with the constraint that the probabilities X sum to 1, using a modified 564 version of the rate matrix `Q`. 565 566 Parameters 567 ---------- 568 `verbose` : bool, optional 569 If True, prints the rate matrix, the modified matrix, the right-hand side vector, 570 and the resulting steady-state probabilities. 571 572 Returns 573 ------- 574 **probs** : dict 575 A dictionary mapping each state (node in the graph) to its steady-state probability. 576 577 Notes 578 ----- 579 - The graph must already contain valid transition rates on its edges. 580 - If the rate matrix has not yet been created, `createRateMatrix()` is called automatically. 581 - The method modifies the last column of the rate matrix to impose the normalization constraint 582 `sum(X) = 1`. 583 - The inverse of the modified matrix is computed directly, so the graph should be small enough 584 to avoid numerical instability or performance issues. 585 - The resulting probabilities are stored internally in `_state_probabilities` and can be accessed 586 via the `state_probabilities` property. 587 """ 588 589 # compute transition matrix if not already done 590 if self._rate_matrix is None: self.createRateMatrix() 591 Q2 = self._rate_matrix.copy() 592 if verbose: 593 print("Rate matrix Q") 594 print(f'Q=\n{Q2}\n') 595 596 Q2[:, -1] = 1 597 if verbose: 598 print(f'Matrix is changed to\nQ2=\n{Q2}\n') 599 600 b = np.zeros(len(Q2)) 601 b[-1] = 1 602 if verbose: 603 print("Solve X = b @ Q2") 604 print(f'b=\n{b}\n') 605 606 # state probabilities 607 X = b @ linalg.inv(Q2) # compute the matrix inverse 608 609 # Generate a matrix with P(X,Y) 610 matrix = {} 611 for n in self.nodes: 612 matrix[n] = X[ self._n2i[n] ] 613 614 if verbose: 615 print("Steady-state probabilities X") 616 print(f'X=\n{matrix}\n') 617 self._state_probabilities = matrix 618 return matrix
Compute the steady-state probabilities of the continuous-time Markov chain (CTMC).
This method calculates the stationary distribution of the CTMC defined by the rate
matrix of the state transition graph. The computation solves the linear system
XQ = 0
with the constraint that the probabilities X sum to 1, using a modified
version of the rate matrix Q
.
Parameters
verbose
: bool, optional
If True, prints the rate matrix, the modified matrix, the right-hand side vector,
and the resulting steady-state probabilities.
Returns
probs : dict A dictionary mapping each state (node in the graph) to its steady-state probability.
Notes
- The graph must already contain valid transition rates on its edges.
- If the rate matrix has not yet been created,
createRateMatrix()
is called automatically. - The method modifies the last column of the rate matrix to impose the normalization constraint
sum(X) = 1
. - The inverse of the modified matrix is computed directly, so the graph should be small enough to avoid numerical instability or performance issues.
- The resulting probabilities are stored internally in
_state_probabilities
and can be accessed via thestate_probabilities
property.
620 def calcTransitionProbabilities(self): 621 """ 622 Compute and store transition probabilities for each node in the graph of the embedded discrete Markov chain. 623 624 This method processes each node in the graph, treating outgoing edges as 625 transitions in an embedded discrete-time Markov chain derived from a continuous-time 626 Markov chain (CTMC). Transition *rates* on outgoing edges are normalized into 627 *probabilities*, and relevant data is stored in the node attributes. 628 629 For each node, the following attributes are set: 630 - `"transitionProbs"` : np.ndarray 631 Array of transition probabilities to successor nodes. 632 - `"transitionNodes"` : list 633 List of successor node identifiers (ordered as in the probability array). 634 - `"sumOutgoingRates"` : float 635 Sum of all outgoing transition rates. Used to model the exponential waiting time 636 in continuous-time simulation (`T ~ Exp(sumOutgoingRates)`). 637 638 Notes 639 ----- 640 - The graph must be a directed graph where each edge has a `'rate'` attribute. 641 - Nodes with no outgoing edges are considered absorbing states. 642 - Useful for discrete-event simulation and Markov chain Monte Carlo (MCMC) modeling. 643 - Sets an internal flag `_transitionProbabilitesComputed = True` upon completion. 644 645 Raises 646 ------ 647 KeyError 648 If any outgoing edge lacks a `'rate'` attribute. 649 """ 650 for node in self: 651 successors = list(self.successors(node)) 652 rates = np.array([self[node][nbr]['rate'] for nbr in successors]) 653 probs = rates / rates.sum() # Normalize to probabilities 654 self.nodes[node]["transitionProbs"] = probs # probabilities 655 self.nodes[node]["transitionNodes"] = successors # list of nodes 656 self.nodes[node]["sumOutgoingRates"] = rates.sum() # time in state ~ EXP(sumOutgoingRates) 657 numel = len(successors) 658 if numel==0: 659 print(f"Absorbing state: {node}") 660 self._transitionProbabilitesComputed = True
Compute and store transition probabilities for each node in the graph of the embedded discrete Markov chain.
This method processes each node in the graph, treating outgoing edges as transitions in an embedded discrete-time Markov chain derived from a continuous-time Markov chain (CTMC). Transition rates on outgoing edges are normalized into probabilities, and relevant data is stored in the node attributes.
For each node, the following attributes are set:
"transitionProbs"
: np.ndarray Array of transition probabilities to successor nodes."transitionNodes"
: list List of successor node identifiers (ordered as in the probability array)."sumOutgoingRates"
: float Sum of all outgoing transition rates. Used to model the exponential waiting time in continuous-time simulation (T ~ Exp(sumOutgoingRates)
).
Notes
- The graph must be a directed graph where each edge has a
'rate'
attribute. - Nodes with no outgoing edges are considered absorbing states.
- Useful for discrete-event simulation and Markov chain Monte Carlo (MCMC) modeling.
- Sets an internal flag
_transitionProbabilitesComputed = True
upon completion.
Raises
KeyError
If any outgoing edge lacks a 'rate'
attribute.
662 def simulateMarkovModel(self, startNode=None, numEvents=10): 663 """ 664 Simulate a trajectory through the state transition graph as a continuous-time Markov chain (CTMC). 665 666 This method performs a discrete-event simulation of the CTMC by generating a random walk 667 through the graph, starting at `startNode` (or a default node if not specified). Transitions 668 are chosen based on the normalized transition probabilities computed from edge rates, and 669 dwell times are sampled from exponential distributions based on the total outgoing rate 670 from each state. 671 672 Parameters 673 ---------- 674 - **startNode** : hashable, optional 675 The node at which the simulation starts. If not specified, the first node in the graph 676 (based on insertion order) is used. 677 678 - **numEvents** : int, optional 679 The number of transition events (steps) to simulate. Default is 10. 680 681 Returns 682 ------- 683 **states** : list 684 A list of visited states, representing the sequence of nodes traversed in the simulation. 685 686 **times** : numpy.ndarray 687 An array of dwell times in each state, sampled from exponential distributions 688 with rate equal to the sum of outgoing transition rates from each state. 689 690 Notes 691 ----- 692 - This method automatically calls `calcTransitionProbabilities()` if transition probabilities 693 have not yet been computed. 694 - If an absorbing state (i.e., a node with no outgoing edges) is reached, the simulation stops early. 695 - The cumulative time can be obtained by `np.cumsum(times)`. 696 - State durations follow `Exp(sumOutgoingRates)` distributions per CTMC theory. 697 """ 698 if not self._transitionProbabilitesComputed: self.calcTransitionProbabilities() 699 700 # simulate random walk through Graph = Markov Chain 701 if startNode is None: 702 startNode = next(iter(self)) 703 704 705 states = [] # list of states 706 times = np.zeros(numEvents) # time per state 707 708 node = startNode 709 for i in range(numEvents): 710 outgoingRate = self.nodes[node]["sumOutgoingRates"] # time in state ~ EXP(sumOutgoingRates) 711 times[i] = random.expovariate(outgoingRate) # exponentially distributed time in state 712 states.append(node) 713 714 # get next node 715 probs = self.nodes[node]["transitionProbs"] # probabilities 716 successors = self.nodes[node]["transitionNodes"] # list of nodes 717 numel = len(successors) 718 if numel==0: 719 print(f"Absorbing state: {node}") 720 break 721 elif numel==1: 722 nextNode = successors[0] 723 else: 724 nextNode = successors[ np.random.choice(numel, p=probs) ] 725 # edge = G[node][nextNode] 726 node = nextNode 727 728 return states, times
Simulate a trajectory through the state transition graph as a continuous-time Markov chain (CTMC).
This method performs a discrete-event simulation of the CTMC by generating a random walk
through the graph, starting at startNode
(or a default node if not specified). Transitions
are chosen based on the normalized transition probabilities computed from edge rates, and
dwell times are sampled from exponential distributions based on the total outgoing rate
from each state.
Parameters
startNode : hashable, optional The node at which the simulation starts. If not specified, the first node in the graph (based on insertion order) is used.
numEvents : int, optional The number of transition events (steps) to simulate. Default is 10.
Returns
states : list A list of visited states, representing the sequence of nodes traversed in the simulation.
times : numpy.ndarray An array of dwell times in each state, sampled from exponential distributions with rate equal to the sum of outgoing transition rates from each state.
Notes
- This method automatically calls
calcTransitionProbabilities()
if transition probabilities have not yet been computed. - If an absorbing state (i.e., a node with no outgoing edges) is reached, the simulation stops early.
- The cumulative time can be obtained by
np.cumsum(times)
. - State durations follow
Exp(sumOutgoingRates)
distributions per CTMC theory.
730 def calcProbabilitiesSimulationRun(self, states, times): 731 """ 732 Estimate empirical state probabilities from a simulation run of the CTMC. 733 734 This method computes the fraction of total time spent in each state during a simulation. 735 It uses the sequence of visited states and corresponding sojourn times to estimate the 736 empirical (time-weighted) probability distribution over states. 737 738 Parameters 739 ---------- 740 - **states** : list of tuple 741 A list of state identifiers visited during the simulation. 742 - **times** : numpy.ndarray 743 An array of sojourn times corresponding to each visited state, same length as `states`. 744 745 Returns 746 ------- 747 **simProbs** : dict 748 A dictionary mapping each unique state to its estimated empirical probability, 749 based on total time spent in that state relative to the simulation duration. 750 751 Notes 752 ----- 753 - Internally, states are mapped to flat indices for efficient counting using NumPy. 754 - This method assumes that `states` contains only valid node identifiers from the graph. 755 - If `_n2i` (node-to-index mapping) has not been initialized, it will be computed here. 756 - This is especially useful for validating analytical steady-state probabilities 757 against sampled simulation results. Or in case of state-space explosions. 758 """ 759 if self._n2i is None: 760 n2i = {} 761 nodes = list(self.nodes) 762 for i,node in enumerate(nodes): 763 n2i[node] = i 764 self._n2i = n2i 765 766 767 # Map tuples to flat indices 768 data = np.array(states) 769 max_vals = data.max(axis=0) + 1 770 indices = np.ravel_multi_index(data.T, dims=max_vals) 771 772 # Weighted bincount 773 bincount = np.bincount(indices, weights=times) 774 775 # Decode back to tuples 776 tuple_indices = np.array(np.unravel_index(np.nonzero(bincount)[0], shape=max_vals)).T 777 timePerState = bincount[bincount > 0] 778 779 simProbs = {} 780 totalTime = times.sum() 781 for t, c in zip(tuple_indices, timePerState): 782 simProbs[tuple(t)] = c / totalTime 783 #print(f"Tuple {tuple(t)} → weighted count: {c}") 784 return simProbs
Estimate empirical state probabilities from a simulation run of the CTMC.
This method computes the fraction of total time spent in each state during a simulation. It uses the sequence of visited states and corresponding sojourn times to estimate the empirical (time-weighted) probability distribution over states.
Parameters
- states : list of tuple A list of state identifiers visited during the simulation.
- times : numpy.ndarray
An array of sojourn times corresponding to each visited state, same length as
states
.
Returns
simProbs : dict A dictionary mapping each unique state to its estimated empirical probability, based on total time spent in that state relative to the simulation duration.
Notes
- Internally, states are mapped to flat indices for efficient counting using NumPy.
- This method assumes that
states
contains only valid node identifiers from the graph. - If
_n2i
(node-to-index mapping) has not been initialized, it will be computed here. - This is especially useful for validating analytical steady-state probabilities against sampled simulation results. Or in case of state-space explosions.
934 def drawTransitionGraph(self, pos=None, 935 bended=False, node_size=1000, num=1, rad=-0.2, 936 edge_labels=None, node_shapes='o', figsize=(14, 7), 937 clear=True, fontsize=10, 938 fontcolor='black', label_pos=0.75): 939 """ 940 Visualize the state transition graph with nodes and directed edges. 941 942 This method draws the graph using `matplotlib` and `networkx`, including labels, 943 node colors, and optional curved (bended) edges for better visualization of 944 bidirectional transitions. It supports layout customization and is useful for 945 understanding the structure of the Markov model. 946 947 Parameters 948 ---------- 949 - pos : dict, optional 950 A dictionary mapping nodes to positions. If None, a Kamada-Kawai layout is used. 951 952 - bended : bool, optional 953 If True, edges are drawn with curvature using arcs. Useful for distinguishing 954 bidirectional transitions. Default is False. 955 956 - node_size : int, optional 957 Size of the nodes in the plot. Default is 1000. 958 959 - num : int, optional 960 Figure number for matplotlib (useful for managing multiple figures). Default is 1. 961 962 - rad : float, optional 963 The curvature radius for bended edges. Only used if `bended=True`. Default is -0.2. 964 965 - edge_labels : dict, optional 966 Optional dictionary mapping edges `(u, v)` to labels. If None, the `'label'` 967 attribute from each edge is used. 968 969 - node_shapes : str, optional 970 Shape of the nodes (e.g., `'o'` for circle, `'s'` for square). Default is `'o'`. 971 972 - figsize : tuple, optional 973 Size of the matplotlib figure. Default is (14, 7). 974 975 - clear : bool, optional 976 If True, clears the figure before plotting. Default is True. 977 978 - fontsize : int, optional 979 Font size used for node labels. Default is 10. 980 981 - fontcolor : str, optional 982 Font color for node labels. Default is `'black'`. 983 984 - label_pos : float, optional 985 Position of edge labels along the edge (0=start, 1=end). Default is 0.75. 986 987 Returns 988 ------- 989 - fig : matplotlib.figure.Figure 990 The created matplotlib figure object. 991 992 - ax : matplotlib.axes.Axes 993 The axes object where the graph is drawn. 994 995 - pos : dict 996 The dictionary of node positions used for drawing. 997 998 Notes 999 ----- 1000 - Node labels and colors must be set beforehand using `addState()` or `setStateProperties()`. 1001 - Edge attributes must include `'label'` and `'color'`. 1002 - This method modifies the current matplotlib figure and is intended for interactive or inline use. 1003 """ 1004 node_labels = {node: data["label"] for node, data in self.nodes(data=True)} 1005 node_colors = {node: data["color"] for node, data in self.nodes(data=True)} 1006 node_colors = [data["color"] for node, data in self.nodes(data=True)] 1007 #node_shapes = {(hw,sw): "o" if hw+sw<G.N else "s" for hw,sw in G.nodes()} 1008 1009 edge_cols = {(u,v): self[u][v]["color"] for u,v in self.edges} 1010 1011 if edge_labels is None: 1012 edge_labels = {(u,v): self[u][v]["label"] for u,v in self.edges} 1013 1014 if pos is None: 1015 pos = nx.kamada_kawai_layout(self) 1016 1017 plt.figure(figsize=figsize, num=num, clear=clear) 1018 nx.draw_networkx_nodes(self, pos, node_color=node_colors, node_shape = node_shapes, node_size=node_size) 1019 1020 nx.draw_networkx_labels(self, pos, labels=node_labels, font_size=fontsize, font_color=fontcolor) # Draw node labels 1021 if bended: 1022 nx.draw_networkx_edges(self, pos, width=1, edge_color=[edge_cols[edge] for edge in self.edges], 1023 node_size = node_size, 1024 arrows = True, arrowstyle = '-|>', 1025 connectionstyle=f"arc3,rad={rad}") 1026 self._draw_networkx_edge_labels(pos, ax=plt.gca(), edge_labels=edge_labels,rotate=False, rad = rad, label_pos=label_pos) 1027 else: 1028 nx.draw_networkx_edges(self, pos, width=1, edge_color=[edge_cols[edge] for edge in self.edges], 1029 node_size = node_size, 1030 #min_target_margin=17, arrowsize=15, 1031 arrows = True, arrowstyle = '-|>') 1032 nx.draw_networkx_edge_labels(self, pos, edge_labels, label_pos=label_pos) #, verticalalignment='center',) 1033 plt.axis('off'); 1034 1035 return plt.gcf(), plt.gca(), pos
Visualize the state transition graph with nodes and directed edges.
This method draws the graph using matplotlib
and networkx
, including labels,
node colors, and optional curved (bended) edges for better visualization of
bidirectional transitions. It supports layout customization and is useful for
understanding the structure of the Markov model.
Parameters
pos : dict, optional A dictionary mapping nodes to positions. If None, a Kamada-Kawai layout is used.
bended : bool, optional If True, edges are drawn with curvature using arcs. Useful for distinguishing bidirectional transitions. Default is False.
node_size : int, optional Size of the nodes in the plot. Default is 1000.
num : int, optional Figure number for matplotlib (useful for managing multiple figures). Default is 1.
rad : float, optional The curvature radius for bended edges. Only used if
bended=True
. Default is -0.2.edge_labels : dict, optional Optional dictionary mapping edges
(u, v)
to labels. If None, the'label'
attribute from each edge is used.node_shapes : str, optional Shape of the nodes (e.g.,
'o'
for circle,'s'
for square). Default is'o'
.figsize : tuple, optional Size of the matplotlib figure. Default is (14, 7).
clear : bool, optional If True, clears the figure before plotting. Default is True.
fontsize : int, optional Font size used for node labels. Default is 10.
fontcolor : str, optional Font color for node labels. Default is
'black'
.label_pos : float, optional Position of edge labels along the edge (0=start, 1=end). Default is 0.75.
Returns
fig : matplotlib.figure.Figure The created matplotlib figure object.
ax : matplotlib.axes.Axes The axes object where the graph is drawn.
pos : dict The dictionary of node positions used for drawing.
Notes
- Node labels and colors must be set beforehand using
addState()
orsetStateProperties()
. - Edge attributes must include
'label'
and'color'
. - This method modifies the current matplotlib figure and is intended for interactive or inline use.
1037 def animateSimulation(self, expectedTimePerState=1, inNotebook=False, **kwargs): 1038 """ 1039 Animate a simulation trajectory through the state transition graph. 1040 The animation of the CTMC simulation run either in a Jupyter notebook or as a regular script. 1041 1042 This method visualizes the path of a simulated or predefined sequence of states by 1043 temporarily highlighting each visited state over time. The time spent in each state 1044 is scaled relative to the average sojourn time to produce a visually smooth animation. 1045 1046 Parameters 1047 ---------- 1048 - **expectedTimePerState** : float, optional 1049 Approximate duration (in seconds) for an average state visit in the animation. 1050 The actual pause duration for each state is scaled proportionally to its sojourn time. 1051 Default is 1 second. 1052 1053 - **inNotebook** : bool, optional 1054 If True, uses Jupyter-compatible animation (with display + clear_output). 1055 If False, uses standard matplotlib interactive animation. Default is True. 1056 1057 - `**kwargs` : dict 1058 Additional keyword arguments including: 1059 - states : list 1060 A list of visited states (nodes) to animate. 1061 - times : list or np.ndarray 1062 Sojourn times in each state (same length as `states`). 1063 - startNode : hashable 1064 Optional start node for automatic simulation if `states` and `times` are not provided. 1065 - numEvents : int 1066 Number of events (state transitions) to simulate. 1067 - Additional drawing-related kwargs are passed to `drawTransitionGraph()`. 1068 1069 Raises 1070 ------ 1071 ValueError 1072 If neither a `(states, times)` pair nor `(startNode, numEvents)` are provided. 1073 1074 Returns 1075 ------- 1076 None 1077 The animation is shown interactively using matplotlib but no data is returned. 1078 1079 Notes 1080 ----- 1081 - If `states` and `times` are not supplied, a simulation is run via `simulateMarkovModel()`. 1082 - The function uses matplotlib's interactive mode (`plt.ion()`) to animate transitions. 1083 - The average sojourn time across all states is used to normalize animation speed. 1084 - Each visited state is highlighted in red for its corresponding (scaled) dwell time. 1085 - Drawing arguments like layout, font size, or color can be customized via kwargs. 1086 """ 1087 # expTimePerState: in seconds 1088 if "states" in kwargs and "times" in kwargs: 1089 states = kwargs["states"] 1090 times = kwargs["times"] 1091 else: 1092 # Run default simulation if data not provided 1093 if "startNode" in kwargs and "numEvents" in kwargs: 1094 states, times = self.simulateMarkovModel(startNode=None, numEvents=10) 1095 else: 1096 raise ValueError("Must provide either ('states' and 'times') or ('startNode' and 'numEvents').") 1097 1098 #avg_OutGoingRate = np.mean([self.nodes[n].get("sumOutgoingRates", 0) for n in self.nodes]) 1099 avg_Time = np.mean([1.0/self.nodes[n].get("sumOutgoingRates", 0) for n in self.nodes]) 1100 1101 1102 # Remove simulation-related keys before drawing 1103 draw_kwargs = {k: v for k, v in kwargs.items() if k not in {"states", "times", "startNode", "numEvents"}} 1104 1105 fig, ax, pos = self.drawTransitionGraph(**draw_kwargs) 1106 plt.tight_layout() 1107 1108 if not inNotebook: 1109 plt.ion() 1110 1111 for node, time in zip(states,times): 1112 if inNotebook: 1113 clear_output(wait=True) 1114 artist = nx.draw_networkx_nodes(self, pos, nodelist=[node], ax=ax, 1115 node_color=self._const.COLOR_NODE_ANIMATION_HIGHLIGHT, node_size=1000) 1116 display(fig) 1117 else: 1118 artist = nx.draw_networkx_nodes(self, pos, nodelist=[node],ax=ax, 1119 node_color=self._const.COLOR_NODE_ANIMATION_HIGHLIGHT, node_size=1000) 1120 plt.draw() 1121 plt.pause(time/avg_Time*expectedTimePerState) 1122 1123 # Remove highlighted node 1124 artist.remove() 1125 1126 pass
Animate a simulation trajectory through the state transition graph. The animation of the CTMC simulation run either in a Jupyter notebook or as a regular script.
This method visualizes the path of a simulated or predefined sequence of states by temporarily highlighting each visited state over time. The time spent in each state is scaled relative to the average sojourn time to produce a visually smooth animation.
Parameters
expectedTimePerState : float, optional Approximate duration (in seconds) for an average state visit in the animation. The actual pause duration for each state is scaled proportionally to its sojourn time. Default is 1 second.
inNotebook : bool, optional If True, uses Jupyter-compatible animation (with display + clear_output). If False, uses standard matplotlib interactive animation. Default is True.
**kwargs
: dict Additional keyword arguments including: - states : list A list of visited states (nodes) to animate. - times : list or np.ndarray Sojourn times in each state (same length asstates
). - startNode : hashable Optional start node for automatic simulation ifstates
andtimes
are not provided. - numEvents : int Number of events (state transitions) to simulate. - Additional drawing-related kwargs are passed todrawTransitionGraph()
.
Raises
ValueError
If neither a (states, times)
pair nor (startNode, numEvents)
are provided.
Returns
None The animation is shown interactively using matplotlib but no data is returned.
Notes
- If
states
andtimes
are not supplied, a simulation is run viasimulateMarkovModel()
. - The function uses matplotlib's interactive mode (
plt.ion()
) to animate transitions. - The average sojourn time across all states is used to normalize animation speed.
- Each visited state is highlighted in red for its corresponding (scaled) dwell time.
- Drawing arguments like layout, font size, or color can be customized via kwargs.
1130 def states(self, data=False): 1131 """ 1132 Return a view of the graph's states (nodes), optionally with attributes. 1133 1134 This method is functionally identical to `self.nodes()` in NetworkX, but is 1135 renamed to `states()` to reflect the semantics of a state transition graph or 1136 Markov model, where nodes represent system states. 1137 1138 Parameters 1139 ---------- 1140 **data** : bool, optional 1141 If True, returns a view of `(node, attribute_dict)` pairs. 1142 If False (default), returns a view of node identifiers only. 1143 1144 Returns 1145 ------- 1146 *networkx.classes.reportviews.NodeView* or NodeDataView 1147 A view over the graph’s nodes or `(node, data)` pairs, depending on the `data` flag. 1148 1149 Examples 1150 -------- 1151 >>> G.states() 1152 ['A', 'B', 'C'] 1153 1154 >>> list(G.states(data=True)) 1155 [('A', {'color': 'red', 'label': 'Start'}), ('B', {...}), ...] 1156 1157 Notes 1158 ----- 1159 - This is a convenience wrapper for semantic clarity in models where nodes represent states. 1160 - The view is dynamic: changes to the graph are reflected in the returned view. 1161 """ 1162 return self.nodes(data=data)
Return a view of the graph's states (nodes), optionally with attributes.
This method is functionally identical to self.nodes()
in NetworkX, but is
renamed to states()
to reflect the semantics of a state transition graph or
Markov model, where nodes represent system states.
Parameters
data : bool, optional
If True, returns a view of (node, attribute_dict)
pairs.
If False (default), returns a view of node identifiers only.
Returns
networkx.classes.reportviews.NodeView or NodeDataView
A view over the graph’s nodes or (node, data)
pairs, depending on the data
flag.
Examples
>>> G.states()
['A', 'B', 'C']
>>> list(G.states(data=True))
[('A', {'color': 'red', 'label': 'Start'}), ('B', {...}), ...]
Notes
- This is a convenience wrapper for semantic clarity in models where nodes represent states.
- The view is dynamic: changes to the graph are reflected in the returned view.
1165 def setStateColor(self, state, color): 1166 """ 1167 Set the color attribute of a specific state (node). 1168 1169 Parameters 1170 ---------- 1171 **state** : hashable 1172 The identifier of the state whose color is to be updated. 1173 1174 **color** : str 1175 The new color to assign to the state (used for visualization). 1176 1177 Returns 1178 ------- 1179 None 1180 """ 1181 self.nodes[state]["color"] = color
Set the color attribute of a specific state (node).
Parameters
state : hashable The identifier of the state whose color is to be updated.
color : str The new color to assign to the state (used for visualization).
Returns
None
1184 def setStateLabel(self, state, label): 1185 """ 1186 Set the label attribute of a specific state (node). 1187 1188 Parameters 1189 ---------- 1190 **state** : hashable 1191 The identifier of the state whose label is to be updated. 1192 1193 **label** : str 1194 The new label to assign to the state (used for visualization or annotation). 1195 1196 Returns 1197 ------- 1198 None 1199 """ 1200 self.nodes[state]["label"] = label
Set the label attribute of a specific state (node).
Parameters
state : hashable The identifier of the state whose label is to be updated.
label : str The new label to assign to the state (used for visualization or annotation).
Returns
None
1202 def prob(self, state): 1203 """ 1204 Return the steady-state probability of a specific state. 1205 1206 Parameters 1207 ---------- 1208 **state** : hashable 1209 The identifier of the state whose probability is to be retrieved. 1210 1211 Returns 1212 ------- 1213 **float** 1214 The steady-state probability of the specified state. 1215 1216 Raises 1217 ------ 1218 KeyError 1219 If the state is not found in the steady-state probability dictionary. 1220 """ 1221 return self.state_probabilities[state]
Return the steady-state probability of a specific state.
Parameters
state : hashable The identifier of the state whose probability is to be retrieved.
Returns
float The steady-state probability of the specified state.
Raises
KeyError If the state is not found in the steady-state probability dictionary.
1224 def probs(self): 1225 """ 1226 Return the steady-state probabilities for all states. 1227 1228 Returns 1229 ------- 1230 **dict** 1231 A dictionary mapping each state to its steady-state probability. 1232 1233 Notes 1234 ----- 1235 - The probabilities are computed via `calcSteadyStateProbs()` and stored internally. 1236 - Use this method to access the full steady-state distribution. 1237 """ 1238 return self.state_probabilities
Return the steady-state probabilities for all states.
Returns
dict A dictionary mapping each state to its steady-state probability.
Notes
- The probabilities are computed via
calcSteadyStateProbs()
and stored internally. - Use this method to access the full steady-state distribution.
1240 def getEmptySystemProbs(self, state=None): 1241 """ 1242 Create an initial state distribution where all probability mass is in a single state. 1243 1244 By default, this method returns a distribution where the entire probability mass 1245 is placed on the first node in the graph. Alternatively, a specific starting state 1246 can be provided. 1247 1248 Parameters 1249 ---------- 1250 **state** : hashable, optional 1251 The state to initialize with probability 1. If None, the first node 1252 in the graph (based on insertion order) is used. 1253 1254 Returns 1255 ------- 1256 **dict** 1257 A dictionary representing the initial distribution over states. 1258 The selected state has probability 1, and all others have 0. 1259 1260 Examples 1261 -------- 1262 >>> X0 = G.getEmptySystem() 1263 >>> sum(X0.values()) # 1.0 1264 1265 Notes 1266 ----- 1267 - This method is useful for setting up deterministic initial conditions 1268 in transient simulations of a CTMC. 1269 - The return format is compatible with methods like `transientProbs()`. 1270 """ 1271 if state is None: 1272 first_node = next(iter(self.nodes)) 1273 X0 = {s: 1 if s == first_node else 0 for s in self} 1274 return X0
Create an initial state distribution where all probability mass is in a single state.
By default, this method returns a distribution where the entire probability mass is placed on the first node in the graph. Alternatively, a specific starting state can be provided.
Parameters
state : hashable, optional The state to initialize with probability 1. If None, the first node in the graph (based on insertion order) is used.
Returns
dict A dictionary representing the initial distribution over states. The selected state has probability 1, and all others have 0.
Examples
>>> X0 = G.getEmptySystem()
>>> sum(X0.values()) # 1.0
Notes
- This method is useful for setting up deterministic initial conditions in transient simulations of a CTMC.
- The return format is compatible with methods like
transientProbs()
.
1278 def transientProbs(self, initialDistribution, t): 1279 """ 1280 Compute the transient state probabilities at time `t` for the CTMC. 1281 1282 This method solves the forward equation of a continuous-time Markov chain (CTMC): 1283 X(t) = X0 · expm(Q·t), where `Q` is the generator (rate) matrix and `X0` is the 1284 initial state distribution. It returns a dictionary mapping each state to its 1285 probability at time `t`. 1286 1287 Parameters 1288 ---------- 1289 - **initialDistribution** : dict 1290 A dictionary mapping states (nodes) to their initial probabilities at time `t=0`. 1291 The keys must match the nodes in the graph, and the values should sum to 1. 1292 1293 - **t** : float 1294 The time at which the transient distribution is to be evaluated. 1295 1296 Returns 1297 ------- 1298 **dict** 1299 A dictionary mapping each state (node) to its probability at time `t`. 1300 1301 Notes 1302 ----- 1303 - The rate matrix `Q` is generated automatically if not yet computed. 1304 - The state order in the rate matrix corresponds to the internal `_n2i` mapping. 1305 - Internally, the dictionary `initialDistribution` is converted into a vector 1306 aligned with the state indexing. 1307 - The transient solution is computed using `scipy.linalg.expm` for matrix exponentiation. 1308 1309 Raises 1310 ------ 1311 ValueError 1312 If the input dictionary has missing states or probabilities that do not sum to 1. 1313 (This is not enforced but recommended for correctness.) 1314 """ 1315 if self._rate_matrix is None: self.createRateMatrix() 1316 Q = self._rate_matrix 1317 X0 = np.zeros(len(initialDistribution)) 1318 1319 for i,n in enumerate(self.nodes): 1320 X0[ self._n2i[n] ] = initialDistribution[n] 1321 1322 1323 X = X0 @ expm(Q*t) 1324 1325 matrix = {} 1326 for n in self.nodes: 1327 matrix[n] = X[ self._n2i[n] ] 1328 1329 return matrix
Compute the transient state probabilities at time t
for the CTMC.
This method solves the forward equation of a continuous-time Markov chain (CTMC):
X(t) = X0 · expm(Q·t), where Q
is the generator (rate) matrix and X0
is the
initial state distribution. It returns a dictionary mapping each state to its
probability at time t
.
Parameters
initialDistribution : dict A dictionary mapping states (nodes) to their initial probabilities at time
t=0
. The keys must match the nodes in the graph, and the values should sum to 1.t : float The time at which the transient distribution is to be evaluated.
Returns
dict
A dictionary mapping each state (node) to its probability at time t
.
Notes
- The rate matrix
Q
is generated automatically if not yet computed. - The state order in the rate matrix corresponds to the internal
_n2i
mapping. - Internally, the dictionary
initialDistribution
is converted into a vector aligned with the state indexing. - The transient solution is computed using
scipy.linalg.expm
for matrix exponentiation.
Raises
ValueError If the input dictionary has missing states or probabilities that do not sum to 1. (This is not enforced but recommended for correctness.)
1331 def symSolveMarkovModel(self): 1332 """ 1333 Symbolically solve the global balance equations of the CTMC. 1334 1335 This method constructs and solves the system of symbolic equations for the 1336 steady-state probabilities of each state in the CTMC. It uses symbolic transition 1337 rates stored on the edges to form balance equations for each node, along with 1338 the normalization constraint that all probabilities sum to 1. 1339 1340 Returns 1341 ------- 1342 - **solution** : dict 1343 A symbolic solution dictionary mapping SymPy symbols (e.g., x_{A}) to 1344 expressions in terms of model parameters. 1345 1346 - **num_solution** : dict 1347 A numerical dictionary mapping each state (node) to its steady-state probability, 1348 computed by substituting the numeric values from `self._subs` into the symbolic solution. 1349 1350 Notes 1351 ----- 1352 - For each node `s`, a symbolic variable `x_{label}` is created using the node's label attribute. 1353 - One balance equation is created per state: total inflow = total outflow. 1354 - An additional constraint `sum(x_i) = 1` ensures proper normalization. 1355 - The symbolic system is solved with `sympy.solve`, and the results are simplified. 1356 - Numeric values are computed by substituting known parameter values (`self._subs`) into the symbolic solution. 1357 1358 Examples 1359 -------- 1360 >>> sym_sol, num_sol = G.symSolveMarkovModel() 1361 >>> sym_sol # symbolic expressions for each state 1362 >>> num_sol # evaluated numerical values for each state 1363 """ 1364 X={} 1365 to_solve = [] 1366 for s in self.nodes(): 1367 X[s] = sp.Symbol(f'x_{{{self.nodes[s]["label"]}}}') 1368 to_solve.append(X[s]) 1369 #% 1370 eqs = [] 1371 for s in self.nodes(): 1372 succ = self.successors(s) # raus 1373 out_rate = 0 1374 for z in list(succ): 1375 #out_rate += X[s]*self.edges[s, z]["sym_rate"] 1376 out_rate += X[s]*self[s][z]["sym_rate"] 1377 1378 in_rate = 0 1379 for r in list(self.predecessors(s)): 1380 #in_rate += X[r]*self.edges[r,s]["sym_rate"] 1381 in_rate += X[r]*self[r][s]["sym_rate"] 1382 eqs.append( sp.Eq(out_rate, in_rate) ) 1383 #% 1384 Xsum = 0 1385 for s in X: 1386 Xsum += X[s] 1387 1388 eqs.append( sp.Eq(Xsum, 1) ) 1389 #% 1390 solution = sp.solve(eqs, to_solve) 1391 simplified_solution = sp.simplify(solution) 1392 #print(simplified_solution) 1393 num_solution = simplified_solution.subs(self._subs) 1394 #num_dict = {str(var): str(sol) for var, sol in simplified_solution.items()} 1395 num_dict = {s: num_solution[X[s]] for s in X} 1396 1397 return simplified_solution, num_dict
Symbolically solve the global balance equations of the CTMC.
This method constructs and solves the system of symbolic equations for the steady-state probabilities of each state in the CTMC. It uses symbolic transition rates stored on the edges to form balance equations for each node, along with the normalization constraint that all probabilities sum to 1.
Returns
solution : dict A symbolic solution dictionary mapping SymPy symbols (e.g., x_{A}) to expressions in terms of model parameters.
num_solution : dict A numerical dictionary mapping each state (node) to its steady-state probability, computed by substituting the numeric values from
self._subs
into the symbolic solution.
Notes
- For each node
s
, a symbolic variablex_{label}
is created using the node's label attribute. - One balance equation is created per state: total inflow = total outflow.
- An additional constraint
sum(x_i) = 1
ensures proper normalization. - The symbolic system is solved with
sympy.solve
, and the results are simplified. - Numeric values are computed by substituting known parameter values (
self._subs
) into the symbolic solution.
Examples
>>> sym_sol, num_sol = G.symSolveMarkovModel()
>>> sym_sol # symbolic expressions for each state
>>> num_sol # evaluated numerical values for each state
1399 def exportTransitionsToExcel(self, excel_path = "output.xlsx", num_colors=9, open_excel=False, verbose=True): 1400 """ 1401 Export transition data from the state transition graph to an Excel file. 1402 1403 This method collects all transitions (edges) from the graph and writes them 1404 to an Excel spreadsheet, optionally applying background colors to distinguish 1405 outgoing transitions from different source states. This export is helpful for 1406 manually verifying the transition structure of the model. 1407 1408 The transition type is defined by the `tid` field of each edge, and mapped to 1409 a human-readable string via the `EXPORT_NAME_TRANSITIONS` dictionary in the 1410 constants module (e.g., `Constants.EXPORT_NAME_TRANSITIONS[tid]`). 1411 1412 See Also 1413 -------- 1414 StateTransitionGraph.addTransition : Method to add edges with symbolic rates and `tid`. 1415 StateTransitionGraph.__init__ : Initializes the graph and loads the constants module. 1416 1417 Parameters 1418 ---------- 1419 - **excel_path** : str, optional 1420 Path to the output Excel file. Default is "output.xlsx". 1421 1422 - num_colors : int, optional 1423 Number of unique background colors to cycle through for different states. 1424 Default is 9. 1425 1426 - open_excel : bool, optional 1427 If True, automatically opens the generated Excel file after export. 1428 Default is False. 1429 1430 - verbose : bool, optional 1431 If True, prints summary information about the export process. 1432 Default is True. 1433 1434 Returns 1435 ------- 1436 None 1437 The function writes data to disk and optionally opens Excel, 1438 but does not return a value. 1439 1440 Notes 1441 ----- 1442 - Each edge must contain a `tid` attribute for transition type and a `label` for the rate. 1443 - Background coloring groups all outgoing transitions from the same source state. 1444 - Requires the constants module to define `EXPORT_NAME_TRANSITIONS`. 1445 """ 1446 1447 rows = [] 1448 for node in self.nodes(): 1449 edges = list(self.edges(node)) 1450 if verbose: print(f"Edges connected to node '{node}':") 1451 1452 for i, (u, v) in enumerate(edges): 1453 latex = self[u][v]["label"] 1454 rate = sp.pretty(self[u][v]["sym_rate"]) 1455 tid = self[u][v]["tid"] 1456 1457 if tid: # 1458 if hasattr(self._const, 'EXPORT_NAME_TRANSITIONS'): 1459 tid_name = self._const.EXPORT_NAME_TRANSITIONS[tid] 1460 else: 1461 tid_name = '' 1462 1463 rows.append({"from": u, "type": tid_name, "to": v, "rate": rate, "latex": latex}) 1464 1465 df = pd.DataFrame(rows, columns=["from", "type", "to", "rate", "latex"]) 1466 df["checked"] = "" 1467 df.to_excel(excel_path, index=False) 1468 1469 # Load workbook with openpyxl 1470 wb = load_workbook(excel_path) 1471 ws = wb.active 1472 1473 # Define colors for different "from" values 1474 cmap = map(plt.cm.Pastel1, range(num_colors)) 1475 # Convert RGBA to hex 1476 color_list = [mcolors.to_hex(color, keep_alpha=False).upper().replace("#", "") for color in cmap] 1477 1478 # Assign color to each unique "from" value using modulo 1479 unique_from = (df["from"].unique()) # Sorted for consistent order 1480 from_color_map = { 1481 value: color_list[i % len(color_list)] for i, value in enumerate(unique_from) 1482 } 1483 1484 # Apply coloring row by row 1485 for i,row in enumerate(ws.iter_rows(min_row=2, max_row=ws.max_row)): # Skip header 1486 from_value = row[0].value # "from" is first column 1487 #fill_color = color_list[i % len(color_list)] 1488 if isinstance(from_value, int): 1489 fill_color = from_color_map[from_value] 1490 else: 1491 fill_color = from_color_map.get(eval(from_value)) 1492 1493 if fill_color: 1494 fill = PatternFill(start_color=fill_color, end_color=fill_color, fill_type="solid") 1495 for cell in row: 1496 cell.fill = fill 1497 1498 # Save styled Excel file 1499 wb.save(excel_path) 1500 if verbose: print(f"File {excel_path} created.") 1501 if open_excel: 1502 print(f"Opening the file {excel_path}.") 1503 os.startfile(excel_path)
Export transition data from the state transition graph to an Excel file.
This method collects all transitions (edges) from the graph and writes them to an Excel spreadsheet, optionally applying background colors to distinguish outgoing transitions from different source states. This export is helpful for manually verifying the transition structure of the model.
The transition type is defined by the tid
field of each edge, and mapped to
a human-readable string via the EXPORT_NAME_TRANSITIONS
dictionary in the
constants module (e.g., Constants.EXPORT_NAME_TRANSITIONS[tid]
).
See Also
StateTransitionGraph.addTransition : Method to add edges with symbolic rates and tid
.
StateTransitionGraph.__init__ : Initializes the graph and loads the constants module.
Parameters
excel_path : str, optional Path to the output Excel file. Default is "output.xlsx".
num_colors : int, optional Number of unique background colors to cycle through for different states. Default is 9.
open_excel : bool, optional If True, automatically opens the generated Excel file after export. Default is False.
verbose : bool, optional If True, prints summary information about the export process. Default is True.
Returns
None The function writes data to disk and optionally opens Excel, but does not return a value.
Notes
- Each edge must contain a
tid
attribute for transition type and alabel
for the rate. - Background coloring groups all outgoing transitions from the same source state.
- Requires the constants module to define
EXPORT_NAME_TRANSITIONS
.