Home Artificial Intelligence Using Bayesian Networks to forecast ancillary service volume in hospitals

Using Bayesian Networks to forecast ancillary service volume in hospitals

1
Using Bayesian Networks to forecast ancillary service volume in hospitals

A Python example using diagnostic input variables

Photo from Unsplash, by EJ Strat

Since I’ve been working with healthcare data (almost 10 years now), forecasting future patient volume has been a tricky nut to crack. There are such a lot of dependencies to think about — patient requests and severity, administrative needs, exam room constraints, a provider just called out sick, a foul snow storm. Plus, unanticipated scenarios can have cascading impacts on scheduling and resource allocation that contradict even the perfect Excel projections.

These sorts of problems are really interesting to attempt to solve from an information perspective, one because they’re tough and you may chew on it for awhile, but in addition because even slight improvements can result in major wins (e.g., improve patient throughput, lower wait times, happier providers, lower costs).

The right way to solve it then? Well, Epic provides us with a number of data, including actual records of when patients arrived for his or her appointments. With historical outputs known, we’re mainly within the space of supervised learning, and Bayesian Networks (BNs) are good probabilistic graphical models.

While most decisions may be made on a single input (e.g., “should I bring a raincoat?”, if the input is “it’s raining”, then the choice is “yes”), BNs can easily handle more complex decision-making — ones involving multiple inputs, of various probability and dependencies. In this text, I’m going to “scratch pad” in python an excellent easy BN that may output a probability rating for a patient arriving in 2 months based on known probabilities for 3 aspects: symptoms, cancer stage, and treatment goal.

Understanding Bayesian Networks:

At its core, a Bayesian Network is a graphical representation of a joint probability distribution using a directed acyclic graph (DAG). Nodes within the DAG represent random variables, and directed edges denote causal relationships or conditional dependencies between these variables. As is true for all data science projects, spending a number of time with the stakeholder at first to properly map the workflows (e.g., variables) involved in decision-making is critical for high-quality predictions.

So, I’ll invent a scenario that we meet our Breast oncology partners they usually explain that three variables are critical for determining whether a patient will need an appointment in 2 months: their symptoms, cancer stage, and treatment goal. I’m making this up as I type, but let’s go together with it.

(In point of fact there shall be dozens of things that influence future patient volumes, a few of singular or multiple dependencies, others completely independent but still influencing).

I’ll say the workflow looks just like the above: Stage is dependent upon their symptom, but treatment type is independent of those and in addition influences the appointment occurring in 2 months.

Based on this, we might the fetch data for these variables from our data source (for us, Epic), which again, would contain known values for our rating node (Appointment_2months), labeled “yes” or “no”.

# install the packages
import pandas as pd # for data manipulation
import networkx as nx # for drawing graphs
import matplotlib.pyplot as plt # for drawing graphs

!pip install pybbn
# for creating Bayesian Belief Networks (BBN)
from pybbn.graph.dag import Bbn
from pybbn.graph.edge import Edge, EdgeType
from pybbn.graph.jointree import EvidenceBuilder
from pybbn.graph.node import BbnNode
from pybbn.graph.variable import Variable
from pybbn.pptc.inferencecontroller import InferenceController

# Create nodes by manually typing in probabilities
Symptom = BbnNode(Variable(0, 'Symptom', ['Non-Malignant', 'Malignant']), [0.30658, 0.69342])
Stage = BbnNode(Variable(1, 'Stage', ['Stage_III_IV', 'Stage_I_II']), [0.92827, 0.07173,
0.55760, 0.44240])
TreatmentTypeCat = BbnNode(Variable(2, 'TreatmentTypeCat', ['Adjuvant/Neoadjuvant', 'Treatment', 'Therapy']), [0.58660, 0.24040, 0.17300])
Appointment_2weeks = BbnNode(Variable(3, 'Appointment_2weeks', ['No', 'Yes']), [0.92314, 0.07686,
0.89072, 0.10928,
0.76008, 0.23992,
0.64250, 0.35750,
0.49168, 0.50832,
0.32182, 0.67818])

Above, let’s manually input some probability scores for levels in each variable (node). In practice, you’d use a crosstab to realize this.

For instance, for the symptom variable, I’ll get frequencies of their 2-levels, about 31% are non-malignant and 69% are malignant.

Photo by writer

Then, we consider the following variable, Stage, and crosstab that with Symptom to get these freqeuncies.

Photo by writer

And, so on and so forth, until all crosstabs between parent-child pairs are defined.

Now, most BNs include many parent-child relationships, so calculating probabilities can get tedious (and majorly error prone), so the function below can calculate the probability matrix for any child node corresponding with 0, 1 or 2 parents.

# This function helps to calculate probability distribution, which fits into BBN (note, can handle as much as 2 parents)
def probs(data, child, parent1=None, parent2=None):
if parent1==None:
# Calculate probabilities
prob=pd.crosstab(data[child], 'Empty', margins=False, normalize='columns').sort_index().to_numpy().reshape(-1).tolist()
elif parent1!=None:
# Check if child node has 1 parent or 2 parents
if parent2==None:
# Caclucate probabilities
prob=pd.crosstab(data[parent1],data[child], margins=False, normalize='index').sort_index().to_numpy().reshape(-1).tolist()
else:
# Caclucate probabilities
prob=pd.crosstab([data[parent1],data[parent2]],data[child], margins=False, normalize='index').sort_index().to_numpy().reshape(-1).tolist()
else: print("Error in Probability Frequency Calculations")
return prob

Then we create the actual BN nodes and the network itself:

# Create nodes through the use of our earlier function to mechanically calculate probabilities
Symptom = BbnNode(Variable(0, 'Symptom', ['Non-Malignant', 'Malignant']), probs(df, child='SymptomCat'))
Stage = BbnNode(Variable(1, 'Stage', ['Stage_I_II', 'Stage_III_IV']), probs(df, child='StagingCat', parent1='SymptomCat'))
TreatmentTypeCat = BbnNode(Variable(2, 'TreatmentTypeCat', ['Adjuvant/Neoadjuvant', 'Treatment', 'Therapy']), probs(df, child='TreatmentTypeCat'))
Appointment_2months = BbnNode(Variable(3, 'Appointment_2months', ['No', 'Yes']), probs(df, child='Appointment_2months', parent1='StagingCat', parent2='TreatmentTypeCat'))

# Create Network
bbn = Bbn()
.add_node(Symptom)
.add_node(Stage)
.add_node(TreatmentTypeCat)
.add_node(Appointment_2months)
.add_edge(Edge(Symptom, Stage, EdgeType.DIRECTED))
.add_edge(Edge(Stage, Appointment_2months, EdgeType.DIRECTED))
.add_edge(Edge(TreatmentTypeCat, Appointment_2months, EdgeType.DIRECTED))

# Convert the BBN to a join tree
join_tree = InferenceController.apply(bbn)

And we’re all set. Now let’s run some hypotheticals through our BN and evaluate the outputs.

Evaluating the BN outputs

First, let’s take a have a look at the probability of every node because it stands, without specifically declaring any conditions.

# Define a function for printing marginal probabilities
# Probabilities for every node
def print_probs():
for node in join_tree.get_bbn_nodes():
potential = join_tree.get_bbn_potential(node)
print("Node:", node)
print("Values:")
print(potential)
print('----------------')

# Use the above function to print marginal probabilities
print_probs()

Node: 1|Stage|Stage_I_II,Stage_III_IV
Values:
1=Stage_I_II|0.67124
1=Stage_III_IV|0.32876
----------------
Node: 0|Symptom|Non-Malignant,Malignant
Values:
0=Non-Malignant|0.69342
0=Malignant|0.30658
----------------
Node: 2|TreatmentTypeCat|Adjuvant/Neoadjuvant,Treatment,Therapy
Values:
2=Adjuvant/Neoadjuvant|0.58660
2=Treatment|0.17300
2=Therapy|0.24040
----------------
Node: 3|Appointment_2weeks|No,Yes
Values:
3=No|0.77655
3=Yes|0.22345
----------------

Meaning, all of the patients on this dataset have a 67% probability of being Stage_I_II, a 69% probability of being Non-Malignant, a 58% probability of requiring Adjuvant/Neoadjuvant treatment, and only 22% of them required an appointment 2 months from now.

We could easily get that from easy frequency tables and not using a BN.

But now, let’s ask a more conditional query: What’s the probability a patient would require care in 2 months on condition that they’ve Stage = Stage_I_II and have a TreatmentTypeCat = Therapy. Also, consider the undeniable fact that the provider knows nothing about their symptoms yet (perhaps they haven’t seen the patient yet).

We’ll run what we all know to be true within the nodes:

# So as to add evidence of events that happened so probability distribution may be recalculated
def evidence(ev, nod, cat, val):
ev = EvidenceBuilder()
.with_node(join_tree.get_bbn_node_by_name(nod))
.with_evidence(cat, val)
.construct()
join_tree.set_observation(ev)

# Add more evidence
evidence('ev1', 'Stage', 'Stage_I_II', 1.0)
evidence('ev2', 'TreatmentTypeCat', 'Therapy', 1.0)
# Print marginal probabilities
print_probs()

Which returns:

Node: 1|Stage|Stage_I_II,Stage_III_IV
Values:
1=Stage_I_II|1.00000
1=Stage_III_IV|0.00000
----------------
Node: 0|Symptom|Non-Malignant,Malignant
Values:
0=Non-Malignant|0.57602
0=Malignant|0.42398
----------------
Node: 2|TreatmentTypeCat|Adjuvant/Neoadjuvant,Treatment,Therapy
Values:
2=Adjuvant/Neoadjuvant|0.00000
2=Treatment|0.00000
2=Therapy|1.00000
----------------
Node: 3|Appointment_2months|No,Yes
Values:
3=No|0.89072
3=Yes|0.10928
----------------

That patient only has an 11% likelihood of arriving in 2 months.

A note in regards to the importance of quality input variables:

The success of a BN in providing a reliable future visit estimate depends heavily on an accurate mapping of workflows for patient care. Patients presenting similarly, in similar conditions, will typically require similar services. The permutation of those inputs, whose characteristics can span from the clinical to administrative, ultimately correspond to a somewhat deterministic path for service needs. However the more complicated or farther out the time projection, the upper the necessity for more specific, intricate BNs with high-quality inputs.

Here’s why:

  1. Accurate Representation: The structure of the Bayesian Network must reflect the actual relationships between variables. Poorly chosen variables or misunderstood dependencies can result in inaccurate predictions and insights.
  2. Effective Inference: Quality input variables enhance the model’s ability to perform probabilistic inference. When variables are accurately connected based on their conditional dependence, the network can provide more reliable insights.
  3. Reduced Complexity: Including irrelevant or redundant variables can unnecessarily complicate the model and increase computational requirements. Quality inputs streamline the network, making it more efficient.

Thanks for reading. Blissful to attach with anyone on LinkedIn! In case you are thinking about the intersection of information science and healthcare or if you may have interesting challenges to share, please leave a comment or DM.

Take a look at a few of my other articles:

Why Balancing Classes is Over-Hyped

Feature Engineering CPT Codes

7 Steps to Design a Basic Neural Network

1 COMMENT

LEAVE A REPLY

Please enter your comment!
Please enter your name here