nhra_gt.sensitivity

Classes

Functions

generate_sensitivity_summary(morris_path, sobol_path, output_path)

Synthesizes Morris and Sobol results into a Markdown report.

Source code in src/nhra_gt/sensitivity.py

def generate_sensitivity_summary(morris_path: Path, sobol_path: Path, output_path: Path) -> None:
    """Synthesizes Morris and Sobol results into a Markdown report."""
    summary = "# Global Sensitivity Analysis Summary\n\n"
    summary += "This report summarizes the findings from the Morris screening and Sobol variance decomposition.\n\n"

    # Morris Section
    if morris_path.exists():
        df_m = pd.read_csv(morris_path, index_col=0)
        summary += "## 1. Morris Screening (Influence & Non-linearity)\n"
        summary += "The Morris method identifies parameters with the greatest overall influence (mu_star) and those with non-linear or interactive effects (sigma).\n\n"
        summary += df_m[["mu_star", "sigma"]].head(5).to_markdown()
        summary += "\n\n"

    # Sobol Section
    if sobol_path.exists():
        df_s = pd.read_csv(sobol_path)
        summary += "## 2. Sobol Analysis (Variance Decomposition)\n"
        summary += "The Sobol method quantifies the percentage of output variance attributable to each parameter (S1) and its total effect including interactions (ST).\n\n"
        summary += (
            df_s[["Parameter", "S1", "ST"]].sort_values("ST", ascending=False).head(5).to_markdown()
        )
        summary += "\n\n"

    summary += "## 3. Key Findings\n"
    # Logic to identify top driver
    if sobol_path.exists():
        top_param = df_s.sort_values("ST", ascending=False).iloc[0]["Parameter"]
        summary += f"- **Primary Driver:** The most influential parameter in the system is **{top_param}**.\n"

    summary += "- **Interactions:** High sigma values in Morris or gaps between ST and S1 in Sobol indicate strong parameter interactions.\n"

    output_path.write_text(summary)

get_parameter_lineage()

Returns a mapping of model parameters to their evidence sources in the context pack.

Source code in src/nhra_gt/sensitivity.py

def get_parameter_lineage() -> dict[str, str]:
    """Returns a mapping of model parameters to their evidence sources in the context pack."""
    return {
        "nominal_cth_share_target": "NHRA Section 127; Federal Financial Relations Agreement (2020-2025)",
        "nep_annual_growth": "IHACPA National Efficient Price Determination 2024-25 (Historical Indexation)",
        "bed_capacity_index": "AIHW Hospital Resources 2022-23 (Bed-to-Population Ratios)",
        "discharge_delay_base": "Medicare UCC Evaluation Report (2024); Aged Care/NDIS Interface Audit",
        "political_salience": "Model Assumption Log (05_assumptions_log.md) - Behavioural Weights",
        "audit_pressure": "NHRA Performance and Accountability Framework (Section 4.2)",
        "rurality_weight": "AIHW Hospital Activity Data (2024) - Regional/Remote peer group weights",
        "cost_shifting_intensity": "Productivity Commission (2023) - Report on Government Services (Health)",
    }

plot_sobol_indices(si, output_path)

Generates Sobol first-order and total sensitivity plots.

Source code in src/nhra_gt/sensitivity.py

def plot_sobol_indices(si: dict[str, Any], output_path: Path) -> None:
    """Generates Sobol first-order and total sensitivity plots."""
    import warnings

    from .visualization.base import PlotConfig, save_figure
    from .visualization.sensitivity import plot_sobol_indices as new_plot_sobol_indices

    warnings.warn(
        "plot_sobol_indices is deprecated, use nhra_gt.visualization.sensitivity instead",
        DeprecationWarning,
        stacklevel=2,
    )
    config = PlotConfig()

    # S1
    fig_s1 = new_plot_sobol_indices(si, config=config, total_order=False)
    save_figure(fig_s1, output_path.parent / (output_path.name + "_s1.png"), config)

    # ST
    fig_st = new_plot_sobol_indices(si, config=config, total_order=True)
    save_figure(fig_st, output_path.parent / (output_path.name + "_st.png"), config)

plot_sobol_heatmap(si, output_path)

Generates a heatmap of second-order interaction indices (S2).

Source code in src/nhra_gt/sensitivity.py

def plot_sobol_heatmap(si: dict[str, Any], output_path: Path) -> None:
    """Generates a heatmap of second-order interaction indices (S2)."""
    import warnings

    from .visualization.base import PlotConfig, save_figure
    from .visualization.sensitivity import plot_sobol_heatmap as new_plot_sobol_heatmap

    warnings.warn(
        "plot_sobol_heatmap is deprecated, use nhra_gt.visualization.sensitivity instead",
        DeprecationWarning,
        stacklevel=2,
    )
    config = PlotConfig()
    fig = new_plot_sobol_heatmap(si, config=config)
    if fig:
        save_figure(fig, output_path.with_suffix(".png"), config)

plot_morris_tornado(df, output_path)

Generates a Morris Tornado plot (mu_star ranking).

Source code in src/nhra_gt/sensitivity.py

def plot_morris_tornado(df: pd.DataFrame, output_path: Path) -> None:
    """Generates a Morris Tornado plot (mu_star ranking)."""
    import warnings

    from .visualization.base import PlotConfig, save_figure
    from .visualization.sensitivity import plot_morris_tornado as new_plot_morris_tornado

    warnings.warn(
        "plot_morris_tornado is deprecated, use nhra_gt.visualization.sensitivity instead",
        DeprecationWarning,
        stacklevel=2,
    )
    config = PlotConfig()
    fig = new_plot_morris_tornado(df, config=config)
    save_figure(fig, output_path.with_suffix(".png"), config)

get_salib_problem(param_names, bounds_override=None, default_variation=0.2)

Generates a SALib-compatible problem dictionary from the Params dataclass.

Parameters:

Name	Type	Description	Default
param_names	list[str]	List of parameter names to include in the GSA.	required
bounds_override	dict[str, list[float]] | None	Optional dictionary mapping parameter names to [min, max] bounds.	None
default_variation	float	If no override, bounds are set to [default * (1-var), default * (1+var)].	0.2

Returns:

Type	Description
dict[str, Any]	A dictionary with 'num_vars', 'names', and 'bounds'.

Source code in src/nhra_gt/sensitivity.py

def get_salib_problem(
    param_names: list[str],
    bounds_override: dict[str, list[float]] | None = None,
    default_variation: float = 0.20,
) -> dict[str, Any]:
    """Generates a SALib-compatible problem dictionary from the Params dataclass.

    Args:
        param_names: List of parameter names to include in the GSA.
        bounds_override: Optional dictionary mapping parameter names to [min, max] bounds.
        default_variation: If no override, bounds are set to [default * (1-var), default * (1+var)].

    Returns:
        A dictionary with 'num_vars', 'names', and 'bounds'.
    """
    bounds_override = bounds_override or {}
    defaults = Params().flatten()

    problem_names = []
    problem_bounds = []

    for name in param_names:
        if name not in defaults:
            raise ValueError(f"Parameter '{name}' not found in Params dataclass.")

        problem_names.append(name)

        if name in bounds_override:
            problem_bounds.append(bounds_override[name])
        else:
            val = float(defaults[name])
            # Handle boolean or binary-like flags if they exist (Params v8 is mostly floats)
            problem_bounds.append(
                [val * (1.0 - default_variation), val * (1.0 + default_variation)]
            )

    return {"num_vars": len(problem_names), "names": problem_names, "bounds": problem_bounds}

evaluate_parallel(model_func, param_values, n_procs=4)

Evaluates the model function in parallel across multiple processes.

Parameters:

Name	Type	Description	Default
model_func	Callable[[ndarray[Any, Any]], float]	Function that takes a parameter array and returns a scalar result.	required
param_values	ndarray[Any, Any]	2D array of shape (n_samples, n_vars).	required
n_procs	int	Number of worker processes to use.	4

Returns:

Type	Description
ndarray[Any, Any]	A 1D array of results.

Source code in src/nhra_gt/sensitivity.py

def evaluate_parallel(
    model_func: Callable[[np.ndarray[Any, Any]], float],
    param_values: np.ndarray[Any, Any],
    n_procs: int = 4,
) -> np.ndarray[Any, Any]:  # pragma: no cover
    """Evaluates the model function in parallel across multiple processes.

    Args:
        model_func: Function that takes a parameter array and returns a scalar result.
        param_values: 2D array of shape (n_samples, n_vars).
        n_procs: Number of worker processes to use.

    Returns:
        A 1D array of results.
    """
    with ProcessPoolExecutor(max_workers=n_procs) as executor:
        # map preserves order
        results = list(executor.map(model_func, param_values))

    return np.array(results)

run_morris_analysis(problem, model_func, n_trajectories=10, n_procs=4, seed=42)

Performs Morris analysis (Elementary Effects screening).

Returns:

Type	Description
DataFrame	A pandas DataFrame with mu_star and sigma indices.

Source code in src/nhra_gt/sensitivity.py

def run_morris_analysis(
    problem: dict[str, Any],
    model_func: Callable[[np.ndarray[Any, Any]], float],
    n_trajectories: int = 10,
    n_procs: int = 4,
    seed: int = 42,
) -> pd.DataFrame:  # pragma: no cover
    """Performs Morris analysis (Elementary Effects screening).

    Returns:
        A pandas DataFrame with mu_star and sigma indices.
    """
    param_values = morris_sampler.sample(problem, N=n_trajectories, seed=seed)

    # Run the model
    results = evaluate_parallel(model_func, param_values, n_procs=n_procs)

    # Perform analysis
    si = morris_analyzer.analyze(problem, param_values, results, conf_level=0.95, seed=seed)

    # Convert to DataFrame
    df = pd.DataFrame(
        {
            "mu": si["mu"],
            "mu_star": si["mu_star"],
            "sigma": si["sigma"],
            "mu_star_conf": si["mu_star_conf"],
        },
        index=problem["names"],
    )

    return df.sort_values("mu_star", ascending=False)

run_sobol_analysis(problem, model_func, n_samples=128, n_procs=4, seed=42)

Performs Sobol variance-based sensitivity analysis.

Parameters:

Name	Type	Description	Default
n_samples	int	The number of samples to generate (must be a power of 2).	128

Returns:

Type	Description
dict[str, Any]	A dictionary containing S1, ST, and S2 indices.

Source code in src/nhra_gt/sensitivity.py

def run_sobol_analysis(
    problem: dict[str, Any],
    model_func: Callable[[np.ndarray[Any, Any]], float],
    n_samples: int = 128,
    n_procs: int = 4,
    seed: int = 42,
) -> dict[str, Any]:  # pragma: no cover
    """Performs Sobol variance-based sensitivity analysis.

    Args:
        n_samples: The number of samples to generate (must be a power of 2).

    Returns:
        A dictionary containing S1, ST, and S2 indices.
    """
    param_values = sobol_sampler.sample(problem, N=n_samples, calc_second_order=True)

    # Run the model
    results = evaluate_parallel(model_func, param_values, n_procs=n_procs)

    # Perform analysis
    si = sobol_analyzer.analyze(
        problem, results, calc_second_order=True, conf_level=0.95, seed=seed
    )

    # SALib dict doesn't always have names, so we add them for our utilities
    if "names" not in si:
        si["names"] = problem["names"]

    return dict(si)

export_sensitivity_indices(si, output_path)

Exports S1 and ST indices to a CSV file.

Source code in src/nhra_gt/sensitivity.py

def export_sensitivity_indices(si: dict[str, Any], output_path: Path) -> None:
    """Exports S1 and ST indices to a CSV file."""
    df = pd.DataFrame(
        {
            "Parameter": si["names"],
            "S1": si["S1"],
            "S1_conf": si["S1_conf"],
            "ST": si["ST"],
            "ST_conf": si["ST_conf"],
        }
    )
    df.to_csv(output_path, index=False)

export_sobol_s2_to_csv(si, output_path)

Exports the S2 interaction matrix to a CSV file.

Source code in src/nhra_gt/sensitivity.py

def export_sobol_s2_to_csv(si: dict[str, Any], output_path: Path) -> None:
    """Exports the S2 interaction matrix to a CSV file."""
    if "S2" not in si or si["S2"] is None:
        return

    names = si["names"]
    s2 = si["S2"]
    # S2 is often a 2D array in modern SALib or a dict of pairs
    # If it's a 2D array, we can directly create a DF
    if isinstance(s2, np.ndarray) and s2.ndim == 2:
        df = pd.DataFrame(s2, index=names, columns=names)
        df.to_csv(output_path)
    else:
        # Some versions return it differently, handle basic case
        pd.DataFrame(str(s2), index=["S2_raw"], columns=["Value"]).to_csv(output_path)

run_psa(distributions, model_func, n_samples=1000, n_procs=4)

Performs Probabilistic Sensitivity Analysis (PSA).

Parameters:

Name	Type	Description	Default
distributions	dict[str, Callable[[int], ndarray[Any, Any]]]	Dict mapping param name to a sampler function (takes N, returns array).	required
model_func	Callable[[ndarray[Any, Any]], float]	Function taking param array (in order of dict keys) -> scalar result.	required
n_samples	int	Number of MC samples.	1000

Returns:

Type	Description
DataFrame	DataFrame with parameters and outcome.

Source code in src/nhra_gt/sensitivity.py

def run_psa(
    distributions: dict[str, Callable[[int], np.ndarray[Any, Any]]],
    model_func: Callable[[np.ndarray[Any, Any]], float],
    n_samples: int = 1000,
    n_procs: int = 4,
) -> pd.DataFrame:
    """Performs Probabilistic Sensitivity Analysis (PSA).

    Args:
        distributions: Dict mapping param name to a sampler function (takes N, returns array).
        model_func: Function taking param array (in order of dict keys) -> scalar result.
        n_samples: Number of MC samples.

    Returns:
        DataFrame with parameters and outcome.
    """
    param_names = list(distributions.keys())

    # Generate samples
    samples = {}
    for name, sampler in distributions.items():
        samples[name] = sampler(n_samples)

    # Create param array for model_func
    # shape: (n_samples, n_vars)
    param_matrix = np.column_stack([samples[name] for name in param_names])

    # Evaluate
    outcomes = evaluate_parallel(model_func, param_matrix, n_procs=n_procs)

    # Build result DF
    df = pd.DataFrame(samples)
    df["outcome"] = outcomes
    return df