Parameter inputs from data

Turn point values, an exported draw matrix, or a weighted posterior sample into the parameter draw matrix run_psa expects.

Open In Colab

Not every analysis starts from ParameterSet distributions: some run at one base-case set of point values, some carry a draw matrix exported from another tool, and some carry a posterior sample with weights. This tutorial shows how to turn each of those three sources into a parameter draw matrix with single_draw (via ParameterSet.at_means), read_draws, and resample_posterior. It walks through examples/parameter_inputs.py step by step; bring your own outputs is the analogue on the outcome side. All three converge on the same rule: once the result is a draw matrix with an iteration index and numeric columns, run_psa runs it identically regardless of its source.

Specifying a model to run

Two interventions compete for a chronic disease, standard care and a new drug; cost and quality-adjusted life-years (QALYs) per iteration are simple functions of three parameters.

import numpy as np
import pandas as pd
from heormodel.models import Outcomes
from heormodel.params import (
    Beta, Gamma, ParameterSet, read_draws, resample_posterior, single_draw,
)
from heormodel.run import run_psa

def model(draws: pd.DataFrame) -> Outcomes:
    base_qaly = 8.0
    effect_drug = base_qaly + draws["u_gain"] * draws["p_response"] * 10
    costs = pd.DataFrame(
        {"Standard care": 40_000.0, "New drug": 40_000.0 + draws["c_drug"]},
        index=draws.index,
    )
    effects = pd.DataFrame(
        {"Standard care": base_qaly, "New drug": effect_drug}, index=draws.index,
    )
    return Outcomes.from_wide(costs, effects)

params = ParameterSet(
    {
        "p_response": Beta.from_mean_se(0.35, 0.05),
        "c_drug": Gamma.from_mean_se(12_000, 1_500),
        "u_gain": Beta.from_mean_se(0.12, 0.03),
    },
    correlation={("p_response", "u_gain"): 0.4},
)

Running a base case at one set of point values

ParameterSet.at_means returns the analytic means as a one-row draw matrix at iteration 0, giving the deterministic run that a probabilistic analysis is usually checked against. It is single_draw(params.means().to_dict()) underneath; call single_draw directly when the point values come from somewhere other than a ParameterSet, a published base case, for instance.

base = params.at_means()
run_psa(model, base).outcomes.summary().round(2)
cost qaly
intervention
Standard care 40000.0 8.00
New drug 52000.0 8.42

The summary holds the deterministic base case, evaluated once at the parameter means rather than over the sampled distribution the later sections use.

Reading a draw matrix from another tool

read_draws validates a CSV path or a DataFrame as a draw matrix. It uses an explicit iteration column when present and otherwise assigns a fresh index; a non-numeric column raises an error before it reaches the model, so a malformed export fails immediately rather than deep inside the analysis. The code below writes a sampled matrix to CSV to stand in for an external tool’s export, then reads it back as a real export would be.

external = params.sample(2_000, seed=1)
external.to_csv("external_draws.csv")
csv_draws = read_draws("external_draws.csv", iteration="iteration")
run_psa(model, csv_draws).outcomes.summary().round(2)
cost qaly
intervention
Standard care 40000.00 8.00
New drug 52002.43 8.42

The summary is the same probabilistic analysis a native ParameterSet sample would give; read_draws changes only where the draws come from, not how they are analyzed.

Resampling a weighted posterior

A calibration or a bootstrap may return parameter rows with a weight column instead of an already-uniform sample. resample_posterior draws whole rows with replacement, each row \(i\) with probability proportional to its weight,

\[\Pr(\text{row } i) = \frac{w_i}{\sum_j w_j},\]

so any correlation between parameters in the posterior survives, then drops the weight column, giving a plain draw matrix run_psa can use. Resampling to an n larger than the input adds no information; it only smooths Monte Carlo noise in downstream expectations.

grid = params.sample(500, seed=2)
grid["weight"] = np.exp(4.0 * grid["p_response"])  # favor higher response
posterior = resample_posterior(grid, n=2_000, seed=7)
print(f"grid mean p_response {grid['p_response'].mean():.3f}, "
      f"resample mean {posterior['p_response'].mean():.3f}")
run_psa(model, posterior).outcomes.summary().round(2)
grid mean p_response 0.348, resample mean 0.360
cost qaly
intervention
Standard care 40000.00 8.00
New drug 51898.82 8.44

The resampled mean of p_response is higher than the unweighted grid mean, because the weights favor higher values, and the reweighting raises the new drug’s expected QALYs.

The calibration workflow produces such a posterior with abc_calibrate and mixes it with literature draws through mix_draws. Next: deterministic sensitivity analysis sweeps parameters one at a time through the same analysis.

Back to top