Example output from stroke outcome model

In this notebook we provide an example of the output from the stroke outcome model assuming IVT is delivered at 90 mins and MT is delivered at 120 mins after stroke onset.

The model provides a sample distribution of mRS scores for 1,000 patients.

Load packages and data file

import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
from outcome_utilities.clinical_outcome import Clinical_outcome

import warnings
warnings.filterwarnings("ignore")

# Load mRS distributions
mrs_dists = pd.read_csv(
    './outcome_utilities/mrs_dist_probs_cumsum.csv', index_col='Stroke type')

View the loaded mRS distributions

For each stroke type (by row) the the imported table shows the cumulative proportion of patients with each mRS score (0-6)

mrs_dists

	0	1	2	3	4	5	6
Stroke type
pre_stroke_nlvo	0.582881	0.745419	0.848859	0.951082	0.993055	1.000000	1.0
pre_stroke_nlvo_ivt_deaths	0.576469	0.737219	0.839522	0.940620	0.982131	0.989000	1.0
pre_stroke_lvo	0.417894	0.560853	0.679283	0.843494	0.957269	1.000000	1.0
pre_stroke_lvo_ivt_deaths	0.403644	0.541728	0.656119	0.814731	0.924626	0.965900	1.0
pre_stroke_lvo_mt_deaths	0.402850	0.540662	0.654829	0.813128	0.922807	0.964000	1.0
no_treatment_nlvo	0.197144	0.460000	0.580032	0.707768	0.855677	0.917702	1.0
no_effect_nlvo_ivt_deaths	0.197271	0.460000	0.577583	0.702252	0.845244	0.904454	1.0
t0_treatment_nlvo_ivt	0.429808	0.630000	0.738212	0.848427	0.929188	0.956300	1.0
no_treatment_lvo	0.050000	0.129000	0.265000	0.429000	0.676000	0.811000	1.0
no_effect_lvo_ivt_deaths	0.047898	0.123576	0.253858	0.410962	0.647576	0.776900	1.0
no_effect_lvo_mt_deaths	0.047781	0.123274	0.253237	0.409957	0.645993	0.775000	1.0
t0_treatment_lvo_ivt	0.112916	0.200000	0.327377	0.484757	0.698212	0.811443	1.0
t0_treatment_lvo_mt	0.314082	0.436315	0.554431	0.712335	0.853604	0.916750	1.0

Set up outcome model and get output

# Set up outcome model
outcome_model = Clinical_outcome(mrs_dists)

# Get outputs
time_to_ivt = 90
time_to_mt = 120
outcomes = outcome_model.calculate_outcomes(
    time_to_ivt, time_to_mt, patients=10000, random_spacing=False)

Show raw model output

The model output is a dictionary of results.

outcomes

{'lvo_untreated_probs': array([0.05 , 0.079, 0.136, 0.164, 0.247, 0.135, 0.189]),
 'nlvo_untreated_probs': array([0.1972, 0.2628, 0.12  , 0.1278, 0.1478, 0.0621, 0.0823]),
 'lvo_ivt_probs': array([0.0926, 0.0865, 0.1298, 0.1581, 0.2195, 0.1171, 0.1964]),
 'lvo_mt_probs': array([0.2085, 0.1272, 0.1377, 0.1697, 0.1703, 0.0783, 0.1083]),
 'nlvo_ivt_probs': array([0.366 , 0.2248, 0.1128, 0.1164, 0.0941, 0.033 , 0.0529]),
 'lvo_untreated_mean_utility': 0.33261,
 'nlvo_untreated_mean_utility': 0.5993989999999999,
 'lvo_ivt_mean_utility': 0.37059999999999993,
 'lvo_mt_mean_utility': 0.528597,
 'nlvo_ivt_mean_utility': 0.7128859999999999,
 'lvo_ivt_added_utility': 0.03798999999999991,
 'lvo_mt_added_utility': 0.19598699999999997,
 'nlvo_ivt_added_utility': 0.113487,
 'lvo_untreated_cum_probs': array([0.05 , 0.129, 0.265, 0.429, 0.676, 0.811, 1.   ]),
 'nlvo_untreated_cum_probs': array([0.1972, 0.46  , 0.58  , 0.7078, 0.8556, 0.9177, 1.    ]),
 'lvo_ivt_cum_probs': array([0.0926, 0.1791, 0.3089, 0.467 , 0.6865, 0.8036, 1.    ]),
 'lvo_mt_cum_probs': array([0.2085, 0.3357, 0.4734, 0.6431, 0.8134, 0.8917, 1.    ]),
 'nlvo_ivt_cum_probs': array([0.366 , 0.5908, 0.7036, 0.82  , 0.9141, 0.9471, 1.    ]),
 'lvo_ivt_shift': array([ 0.0426,  0.0075, -0.0062, -0.0059, -0.0275, -0.0179,  0.0074]),
 'lvo_mt_shift': array([ 0.1585,  0.0482,  0.0017,  0.0057, -0.0767, -0.0567, -0.0807]),
 'nlvo_ivt_shift': array([ 0.1688, -0.038 , -0.0072, -0.0114, -0.0537, -0.0291, -0.0294]),
 'lvo_untreated_mean_mRS': 3.64,
 'nlvo_untreated_mean_mRS': 2.2817,
 'lvo_ivt_mean_mRS': 3.4623,
 'lvo_mt_mean_mRS': 2.6342,
 'nlvo_ivt_mean_mRS': 1.6584,
 'lvo_ivt_mean_shift': -0.1777,
 'lvo_mt_mean_shift': -1.0058,
 'nlvo_ivt_mean_shift': -0.6233,
 'lvo_ivt_improved': 0.1851,
 'lvo_mt_improved': 0.8088,
 'nlvo_ivt_improved': 0.6125}

Plot mRS distributions

mRS distributions

fig = plt.figure(figsize=(10,6))

# nLVO
x = np.arange(7)
width = 0.4
ax1 = fig.add_subplot(121)
y = outcomes['nlvo_untreated_probs']
ax1.bar(x - width/2, y, width = width, label='Untreated', color='k')
y = outcomes['nlvo_ivt_probs']
ax1.bar(x + width/2, y, width = width, label='IVT', color='y')
title = f'nLVO\nTime to IVT {time_to_ivt} mins.'
ax1.set_title(title)
ax1.set_xlabel('mRS')
ax1.set_ylabel('Probability')
ax1.grid()
ax1.legend()

# LVO
width = 0.25
x = np.arange(7)
ax2 = fig.add_subplot(122)
y = outcomes['lvo_untreated_probs']
ax2.bar(x - width, y, width = width, label='Untreated', color='k')
y = outcomes['lvo_ivt_probs']
ax2.bar(x, y, width = width, label='IVT', color='y')
y = outcomes['lvo_mt_probs']
ax2.bar(x + width, y, width = width, label='MT', color='r')
title = f'LVO\nTime to IVT {time_to_ivt} mins; Time to MT {time_to_mt} mins.'
ax2.set_title(title)
ax2.set_xlabel('mRS')
ax2.set_ylabel('Probability')
ax2.grid()
ax2.legend()

plt.tight_layout(pad=2)
plt.savefig('./images/demo_mrs_dists.jpg', dpi=300)
plt.show()

Cumulative mRS distributions

fig = plt.figure(figsize=(10,6))

# nLVO
x = np.arange(7)
width = 0.4
ax1 = fig.add_subplot(121)
y = np.cumsum(outcomes['nlvo_untreated_probs'])
ax1.bar(x - width/2, y, width = width, label='Untreated', color='k')
y = np.cumsum(outcomes['nlvo_ivt_probs'])
ax1.bar(x + width/2, y, width = width, label='IVT', color='y')
title = f'nLVO\nTime to IVT {time_to_ivt} mins.'
ax1.set_title(title)
ax1.set_xlabel('mRS')
ax1.set_ylabel('Cumulative probability')
ax1.grid()
ax1.legend()

# LVO
width = 0.25
x = np.arange(7)
ax2 = fig.add_subplot(122)
y = np.cumsum(outcomes['lvo_untreated_probs'])
ax2.bar(x - width, y, width = width, label='Untreated', color='k')
y = np.cumsum(outcomes['lvo_ivt_probs'])
ax2.bar(x, y, width = width, label='IVT', color='y')
y = np.cumsum(outcomes['lvo_mt_probs'])
ax2.bar(x + width, y, width = width, label='MT', color='r')
title = f'LVO\nTime to IVT {time_to_ivt} mins; Time to MT {time_to_mt} mins.'
ax2.set_title(title)
ax2.set_xlabel('mRS')
ax2.set_ylabel('Cumulative probability')
ax2.grid()
ax2.legend()

plt.tight_layout(pad=2)
plt.savefig('./images/demo_cum_mrs_dists.jpg', dpi=300)
plt.show()

Plot changes in mRS proportions with treatment

fig = plt.figure(figsize=(10,6))

# nLVO
x = np.arange(7)
width = 0.8
ax1 = fig.add_subplot(121)
y = outcomes['nlvo_ivt_shift']
ax1.bar(x, y, width = width, label='IVT', color='y')
title = f'nLVO\nTime to IVT {time_to_ivt} mins.'
ax1.set_title(title)
ax1.set_xlabel('mRS')
ax1.set_ylabel('Change in probability')
ax1.grid()
ax1.legend()

# LVO
width = 0.4
x = np.arange(7)
ax2 = fig.add_subplot(122)
y = outcomes['lvo_ivt_shift']
ax2.bar(x - width/2, y, width = width, label='IVT', color='y')
y = outcomes['lvo_mt_shift']
ax2.bar(x + width/2, y, width = width, label='MT', color='r')
title = f'LVO\nTime to IVT {time_to_ivt} mins; Time to MT {time_to_mt} mins.'
ax2.set_title(title)
ax2.set_xlabel('mRS')
ax2.set_ylabel('Probability')
ax2.grid()
ax2.legend()

plt.tight_layout(pad=2)
plt.savefig('./images/demo_mrs_shifts.jpg', dpi=300)
plt.show()

Other stats

Mean mRS

print('Mean mRS')
print('--------')
print('LVO untreated:', outcomes['lvo_untreated_mean_mRS'])
print('LVO IVT:', outcomes['lvo_ivt_mean_mRS'])
print('LVO MT:', outcomes['lvo_mt_mean_mRS'])
print('nLVO untreated:', outcomes['nlvo_untreated_mean_mRS'])
print('nLVO IVT:', outcomes['nlvo_ivt_mean_mRS'])

Mean mRS
--------
LVO untreated: 3.64
LVO IVT: 3.4623
LVO MT: 2.6342
nLVO untreated: 2.2817
nLVO IVT: 1.6584

Mean shift in mRS

print('Mean mRS shift')
print('--------------')
print('LVO IVT:', outcomes['lvo_ivt_mean_shift'])
print('LVO MT:', outcomes['lvo_mt_mean_shift'])
print('nLVO IVT:', outcomes['nlvo_ivt_mean_shift'])

Mean mRS shift
--------------
LVO IVT: -0.1777
LVO MT: -1.0058
nLVO IVT: -0.6233

The proportion of patients with improved mRS

Assuming all patients move up the mRS.

print('Proportion improved')
print('-------------------')
print('LVO IVT:', outcomes['lvo_ivt_improved'])
print('LVO MT:', outcomes['lvo_mt_improved'])
print('nLVO IVT:', outcomes['nlvo_ivt_improved'])

Proportion improved
-------------------
LVO IVT: 0.1851
LVO MT: 0.8088
nLVO IVT: 0.6125

Utility-weighted mRS outcomes

In addition to mRS, we may calculate utility-weighted mRS (UW-mRS).

UW-mRS incorporates both treatment effect and patient perceived quality of life as a single outcome measure for stroke trials.

UW-mRS scores are based on a pooled analysis of 2,000+ patients. From Wang X, Moullaali TJ, Li Q, Berge E, Robinson TG, Lindley R, et al. Utility-Weighted Modified Rankin Scale Scores for the Assessment of Stroke Outcome. Stroke. 2020 Aug 1;51(8):2411-7.

mRS Score	0	1	2	3	4	5	6
Utility	0.97	0.88	0.74	0.55	0.20	-0.19	0.00

x = outcomes['lvo_untreated_mean_utility']
print(f'LVO untreated UW-mRS: {x:0.3f}')

x1 = outcomes['lvo_ivt_mean_utility']
x2 = outcomes['lvo_ivt_added_utility']
print(f'LVO IVT UW-mRS: {x1:0.3f} (added UW-mRS: {x2:0.3f})')

x1 = outcomes['lvo_mt_mean_utility']
x2 = outcomes['lvo_mt_added_utility']
print(f'LVO MT UW-mRS: {x1:0.3f} (added UW-mRS: {x2:0.3f})')

x = outcomes['nlvo_untreated_mean_utility']
print(f'nLVO untreated UW-mRS: {x:0.3f}')

x1 = outcomes['nlvo_ivt_mean_utility']
x2 = outcomes['nlvo_ivt_added_utility']
print(f'nLVO IVT UW-mRS: {x1:0.3f} (added UW-mRS: {x2:0.3f})')

LVO untreated UW-mRS: 0.333
LVO IVT UW-mRS: 0.371 (added UW-mRS: 0.038)
LVO MT UW-mRS: 0.529 (added UW-mRS: 0.196)
nLVO untreated UW-mRS: 0.599
nLVO IVT UW-mRS: 0.713 (added UW-mRS: 0.113)

An example showing how untreated and treated mRS are compared at a patient level

In the example below we look at the treatment effect of 5 LVO patients treated with MT.

After calculating the treated mRS distribution at the specified treatment time, we can sample random patients by sampling from a uniform 0-1 distribution and using that same sampled value for each patient compare their location on untreated and treated distributions.

For illustration we use more evenly spaced patient values (rather than random), and we can see:

• Patient #1 (P=0.1): mRS untreated = 1, mRS treated = 0

• Patient #2 (P=0.3): mRS untreated = 3, mRS treated = 1

• Patient #3 (P=0.5): mRS untreated = 4, mRS treated = 3

• Patient #4 (P=0.7): mRS untreated = 5, mRS treated = 4

• Patient #5 (P=0.9): mRS untreated = 6, mRS treated = 6

This model is likely a simplification of actual effects, but should capture the average effect of treatment well, and provide a good guide to the proportion of patients who will move at least one mRS unit with treatment.

def draw_horizontal_bar(dist,label=''):
    """
    Draw a stacked horizontal bar chart of the values in 'dist'.
    
    dist  - list or np.array. The probability distribution 
            (non-cumulative).
    label - string. The name printed next to these stacked bars.
    """
    colour_list = plt.rcParams['axes.prop_cycle'].by_key()['color']
    # The first bar will start at this point on the x-axis:
    left = 0
    for i in range(len(dist)):
        # Draw a bar starting from 'left', the end of the previous bar,
        # with a width equal to the probability of this mRS:
        plt.barh(label, width=dist[i], left=left, height=0.5, 
                 label=f'{i}', edgecolor='k', color=colour_list[i%6])
        # Update 'left' with the width of the current bar so that the 
        # next bar drawn will start in the correct place.    
        left += dist[i]

def draw_connections(dist_t0, dist_tne, top_tne=0.25, bottom_t0=0.75):
    """
    Draw lines connecting the mRS bins in the top and bottom rows.
    
    dist_t0, dist_tne - lists or arrays. Probability distributions.
    top_tne, bottom_t0 - floats. y-coordinates just inside the bars. 
    """
    left_t0   = 0.0
    left_tne  = 0.0
    for i, d_t0 in enumerate(dist_t0):
        left_t0  += dist_t0[i]
        left_tne += dist_tne[i]
        plt.plot([left_t0,left_tne],[bottom_t0,top_tne],color='k')

bar_1 = outcomes['lvo_untreated_probs']
bar_2 = outcomes['lvo_mt_probs']

# Draw no effect distribution
draw_horizontal_bar(bar_2, 'Treated (MT)')

# Add legend now to prevent doubling all the labels:
plt.legend(loc='center',ncol=7, title='mRS', 
           bbox_to_anchor=[0.5,0.0,0.0,-0.5])   # Legend below axis.

# Draww t=0 distribution
draw_horizontal_bar(bar_1, 'Untreated')

# Darw connecting lines
draw_connections(bar_1, bar_2)

plt.vlines(0.1, -0.25, 1.25, colors='y', linestyles='dashed')
plt.vlines(0.3, -0.25, 1.25, colors='y', linestyles='dashed')
plt.vlines(0.5, -0.25, 1.25, colors='y', linestyles='dashed')
plt.vlines(0.7, -0.25, 1.25, colors='y', linestyles='dashed')
plt.vlines(0.9, -0.25, 1.25, colors='y', linestyles='dashed')

# Add general content
plt.xlabel('Probability')
plt.title('Treatment effect in 5 LVO patients (dashed lines)')
plt.xlim(0,1)
plt.savefig(f'./images/treatment_shift.jpg', dpi=300, bbox_inches='tight', 
    pad_inches=0.2)
plt.show()