# Modeling patient flow in ER with Markov chains

Playing with AI and name generating routines got me started with markov chains (https://en.wikipedia.org/wiki/Markov_chain) and I decided to do a quick experiment about modeling patient flow in an emergency room. As markov chains can be arbitrarily long, or even infinite, I’m modeling them with generators. With them it’s possible to traverse chain one element at a time until chain is exhausted.

First step was to define a factory function that takes a configuration for chain as an input and returns a generator as an output. I’m not particularly concerned on how easy to use configuration format is, as the configuration can always be generated by another function. In fact, that’s probably a smart idea as soon as one starts doing non-trivial chains. But here we’re doing it by hand and concentrating on simple model of ER flow (I got the data from http://www.slideshare.net/akolker/shs-asq-2010-conference-presentation-hospital-system-patient-flow).

Code is a bit gnarly, as current yield doesn’t play nicely with let and if in current version of Hy. Routine first selects starting element from a list of possible starting elements and yields that. After that subsequent items are selected from a list defined by previous element. In both cases the configuration states what are the odds that a certain element will be selected.

```(import random)

(defn chain-factory [start-elements elements]
"create factory function that can create markov chain instances"
(fn []
"create generator for chain"
(defn select-next-element [elements-list]
"select element"
(let [[high (max (list-comp upper [, (element lower upper) elements-list]))]
[value (.randint random 0 high)]
[matches (list-comp element [, (element lower upper) elements-list]
(> lower value upper))]]
(if matches
(first (.choice random matches))
(first (.choice random elements-list)))))

(setv current-element (select-next-element start-elements))
(setv running true)
(yield current-element)
(while running
(setv next-elements (get elements current-element))
(if next-elements
(setv current-element (select-next-element next-elements))
(setv running false))
(when (not running) (break))
(yield current-element))))
```

Now that we can define our chains, lets model the flow of patients in our example ER. The chain has two possible start elements: Ambulance and Walk In. Both are as likely to occur. From there the chain proceeds through several links, ending at Ambulance diversion, Left (patient decided to leave without receiving treatment) or Discharge. Discharge can be reached after either Emergency room or Medical unit.

```(setv patient-events
(chain-factory :start-elements [(, "Ambulance" 1 50)
(, "Walk in" 51 100)]
:elements {"Ambulance" [(, "Waiting room" 1 100)]
"Walk in" [(, "Waiting room" 1 100)]
"Waiting room" [(, "Emergency room" 1 100)]
"Emergency room" [(, "Ambulance diversion" 1 5)
(, "Left" 6 10)
(, "Discharge" 11 15)
(, "Operating room" 16 60)
(, "Intensive care unit" 61 90)
(, "Medical unit" 91 100)]
"Ambulance diversion" []
"Left" []
"Discharge" []
"Operating room" [(, "Intensive care unit" 1 30)
(, "Medical unit" 31 100)]
"Intensive care unit" [(, "Medical unit" 1 100)]
"Medical unit" [(, "Intensive care unit" 1 5)
(, "Discharge" 6 100)]}))
```

Lets generate some sample chains by creating instance of generator and forcing it to a list:

```=> (list (patient-events))
['Ambulance', 'Waiting room', 'Emergency room', 'Discharge']

=> (list (patient-events))
['Ambulance', 'Waiting room', 'Emergency room', 'Intensive care unit',
'Medical unit', 'Intensive care unit', 'Medical unit', 'Intensive care unit',
'Medical unit', 'Discharge']

=> (list (patient-events))
['Walk in', 'Waiting room', 'Emergency room', 'Ambulance diversion']

=> (list (patient-events))
['Walk in', 'Waiting room', 'Emergency room', 'Left']
```

We can even ask a chain of at least 20 elements to be generated (this is still pretty fast, but for example chain of 50 elements might take a long time to be generated):

```=> (setv history [])
=> (while (< (len (list history)) 20)
...  (setv history (list (patient-events))))
=> history
['Ambulance', 'Waiting room', 'Emergency room', 'Medical unit', 'Intensive care unit',
'Medical unit', 'Intensive care unit', 'Medical unit', 'Intensive care unit',
'Medical unit', 'Intensive care unit', 'Medical unit', 'Intensive care unit',
'Medical unit', 'Intensive care unit', 'Medical unit', 'Intensive care unit',
'Medical unit', 'Intensive care unit', 'Medical unit', 'Intensive care unit',
'Medical unit', 'Discharge']
```

Same principle can be used to create hidden markov models (https://en.wikipedia.org/wiki/Hidden_Markov_model). Instead of having output of process directly in state of markov chain, we use functions that produce the output. For sake of the brevity following code assumes that each output is equally likely to happen than any other. In the example we’re modeling how the actual emergency room of the previous example functions with different amounts of patients. With low flow of the patients, each patient can be treated in a timely manner: nobody will leave without receiving treatment and no ambulance diversions are necessary. However, when the amount of patients is high, some patients may decide to leave without receiving treatment and some may be diverted to different hospital with an ambulance. Model begings with low amount of patients and there’s 25% change of it switching to high amount of patients. When amount of patients is high, there’s 50% chance that it will return back to low.

```(defn low-flow []
(.choice random ["Discharge" "Operating room" "Intensive care unit" "Medical unit"]))

(defn high-flow []
(.choice random ["Ambulance diversion" "Left" "Discharge" "Operating room"
"Intensive care unit" "Medical unit"]))

(setv er-events
(chain-factory :start-elements [(, low-flow 1 100)]
:elements {low-flow [(, low-flow 1 75)
(, high-flow 76 100)]
high-flow [(, low-flow 1 50)
(, high-flow 51 100)]}))

(setv events (genexpr (event) [event (er-events)]))
```

Since er-events returns functions, I’m creating generator events that will automatically call the resulting function, in order to easily get strings. One has to remember that this is an infinite sequence, so forcing it to a list is a bad idea (you’ll run out of memory or out of time). Lets try it out by sampling 4 items at a time:

```=> (list (take 4 events))
['Intensive care unit', 'Medical unit', 'Ambulance diversion', 'Operating room']

=> (list (take 5 events))
['Intensive care unit', 'Discharge', 'Discharge', 'Medical unit']

=> (list (take 5 events))
['Ambulance diversion', 'Intensive care unit', 'Ambulance diversion', 'Discharge']

=> (list (take 5 events))
['Medical unit', 'Discharge', 'Ambulance diversion', 'Ambulance diversion']

=> (list (take 5 events))
['Medical unit', 'Medical unit', 'Intensive care unit', 'Intensive care unit']
```

There’s equal amount of discharges and ambulance diversions. I was expecting discharges to be more frequent than ambulance diversions, but with such a small amount of samples it’s hard to say if the system is functioning correctly. So lets run model for 1 000 000 steps and tally up amount of ambulance diversions and discharges:

```=> (reduce (fn [acc item]
...          (cond [(= item "Ambulance diversion")
...                 {:ambulance-diversion (inc (:ambulance-diversion acc))
...                  :discharge (:discharge acc)}]
...                [(= item "Discharge")
...                 {:ambulance-diversion (:ambulance-diversion acc)
...                  :discharge (inc (:discharge acc))}]
...                [true acc]))
...        (take 1000000 events)
...        {:ambulance-diversion 0
...         :discharge 0})
{'﷐ambulance-diversion': 83220, '﷐discharge': 207921}
```

With 207 921 against 83 220, discharge is over 100% more likely to occur than ambulance diversion. That sounds about the correct ball park to me.

As mentioned in the beginning of this post, creating configuration by hand only works for trivial cases. More interesting results could be obtained, if the percentages could be calculated automatically from data retrieved from actual hospital information system. Then they could be updated periodically and would better reflect current situation. Even more interesting (albeit probably not fitting this case particularly well as process is hopefully well defined) would be to define process automatically by analyzing data in hospital information system. This is the method that is often used to generate random names. One just needs a big enough sample of names that are broken into parts and analyzed to create configuration for markov chain generator.