"""Phase 4 — Project-defined high-charge flag analysis. The primary cohort is defined within each service as providers at or above the 75th percentile of submitted charge. IQR and z-score screens provide separate statistical comparisons. The outputs distinguish provider-service flag counts from payment-dollar measures and include explicit eligibility denominators for the cross-dataset provider overlap. """ from collections import OrderedDict import numpy as np import pandas as pd from .config import ( DATA_PROC, MIN_PROVIDERS_PER_SERVICE, IQR_MULTIPLIER, ZSCORE_MODERATE, ZSCORE_EXTREME, ) from .analysis_insights import save_json HIGH_CHARGE_METHOD = "project_high_charge_top_quartile" HIGH_CHARGE_QUANTILE = 0.75 CONCENTRATION_TARGET_SHARE = 0.40 FLAG_PREVIEW_LIMIT = 2000 COHORT_PREVIEW_LIMIT = 50 BREADTH_BIN_LABELS = ("1", "2-4", "5-9", "10+") DATASET_PERIODS = {"inpatient": "FY2023", "outpatient": "CY2023"} def _period_from_frames(*frames) -> str: """Read a period label when frames carry one; use the project period otherwise.""" period_columns = ("period", "data_period", "year", "data_year", "reporting_year") for frame in frames: if not isinstance(frame, pd.DataFrame): continue for column in period_columns: if column not in frame.columns: continue values = frame[column].dropna().astype(str).unique() if len(values) == 1: return values[0] return "2023" def _analysis_metadata(period: str = "dataset-specific") -> OrderedDict: """Return shared method metadata for every Phase 4 output.""" return OrderedDict( method=HIGH_CHARGE_METHOD, method_label="Project-defined high-charge cohort", inspired_by=( "OIG top-quartile screening approach; no external source data is " "required for this project method" ), flag_label="high-charge flag", flag_plural="high-charge flags", flag_definition=( "A high-charge flag is a provider-service row at or above the " "within-service 75th percentile of avg_submitted_charge." ), thresholds=OrderedDict( min_providers_per_service=int(MIN_PROVIDERS_PER_SERVICE), high_charge_quantile=HIGH_CHARGE_QUANTILE, iqr_multiplier=float(IQR_MULTIPLIER), zscore_moderate=float(ZSCORE_MODERATE), zscore_extreme=float(ZSCORE_EXTREME), concentration_target_share=CONCENTRATION_TARGET_SHARE, ), selection_metric="avg_submitted_charge", proxy_measure=OrderedDict( inpatient="avg_submitted_charge / avg_total_payment", outpatient="avg_submitted_charge / avg_allowed_amount", ), proxy_selection_note=( "The proxy reuses submitted charge as its numerator after the cohort " "is selected on submitted charge; proxy premiums are post-selection descriptions." ), payment_measure=OrderedDict( inpatient="avg_total_payment", outpatient="avg_allowed_amount", ), period=str(period), periods=OrderedDict(DATASET_PERIODS), period_note=( "Inpatient is FY2023 and outpatient is CY2023; cross-dataset " "results are descriptive rather than same-period or causal comparisons." ), ) def _stable_sort(df: pd.DataFrame, columns: list) -> pd.DataFrame: """Sort records deterministically without changing their values.""" if df.empty: return df.reset_index(drop=True) sort_columns = [column for column in columns if column in df.columns] if not sort_columns: return df.reset_index(drop=True) return df.sort_values(sort_columns, kind="mergesort", na_position="last").reset_index(drop=True) def _rounded(value, digits: int): if value is None or pd.isna(value): return None return float(round(float(value), digits)) def _mean(series: pd.Series, digits: int): if series is None or series.dropna().empty: return None return _rounded(series.mean(), digits) def _premium_pct(high_value, low_value): if high_value is None or low_value is None or pd.isna(high_value) or pd.isna(low_value): return None if low_value == 0: return None return _rounded(100 * (high_value / low_value - 1), 1) def _rate(numerator: int, denominator: int): if not denominator: return 0.0 return float(numerator / denominator) def _breadth_rate_profile( provider_ids: set, flagged_providers: set, opportunity_counts: dict, n_bins: int = 5, ) -> tuple[dict, dict]: """Estimate flag rates by provider opportunity-breadth bin.""" values = pd.Series( { provider: float(opportunity_counts.get(provider, 0)) for provider in sorted(provider_ids) }, dtype=float, ) if values.empty: return {}, {} rank = values.rank(method="first") bin_count = min(n_bins, len(values)) bins = pd.qcut(rank, q=bin_count, labels=False, duplicates="drop") assignments = {provider: int(bin_id) for provider, bin_id in bins.items()} rates = {} details = {} for bin_id in sorted(set(assignments.values())): members = [provider for provider, value in assignments.items() if value == bin_id] flagged_count = sum(provider in flagged_providers for provider in members) rates[bin_id] = _rate(flagged_count, len(members)) details[str(bin_id)] = OrderedDict( provider_count=int(len(members)), flagged_provider_count=int(flagged_count), flag_rate=float(rates[bin_id]), mean_eligible_service_opportunities=float( round(values.loc[members].mean(), 2) ), ) return assignments, details def _breadth_matched_expected_overlap( shared_providers: set, inp_flagged: set, out_flagged: set, inp_opportunities: dict, out_opportunities: dict, ) -> tuple[float, dict]: """Compute an independence null while matching opportunity breadth.""" inp_bins, inp_details = _breadth_rate_profile( shared_providers, inp_flagged, inp_opportunities ) out_bins, out_details = _breadth_rate_profile( shared_providers, out_flagged, out_opportunities ) expected = sum( inp_details[str(inp_bins[provider])]["flag_rate"] * out_details[str(out_bins[provider])]["flag_rate"] for provider in shared_providers ) return float(expected), OrderedDict( benchmark="independence conditional on shared provider opportunity-breadth quintiles", n_bins=max(len(inp_details), len(out_details)), inpatient_bins=inp_details, outpatient_bins=out_details, expected_flagged_provider_overlap=float(round(expected, 4)), ) def _flagged_rows(flags_df: pd.DataFrame, flag_column: str = "is_high_charge_flag") -> pd.DataFrame: """Return only rows carrying a true high-charge flag.""" if not isinstance(flags_df, pd.DataFrame) or flags_df.empty: return flags_df.copy() if isinstance(flags_df, pd.DataFrame) else pd.DataFrame() if flag_column in flags_df.columns: return flags_df[flags_df[flag_column].fillna(False)].copy() # This fallback keeps callers that already pass a flagged-only table useful. return flags_df.copy() def _records(df: pd.DataFrame, limit: int = None) -> list: """Convert a stable, optionally truncated frame to JSON records.""" if not isinstance(df, pd.DataFrame) or df.empty: return [] if limit is not None: df = df.head(limit) return df.to_dict(orient="records") # --------------------------------------------------------------------------- # Helpers # --------------------------------------------------------------------------- def _per_service(df: pd.DataFrame, service_col: str, billing_col: str, payment_col: str) -> list: """Return list of (service_code, sub_df) for services with ≥30 providers.""" out = [] for code, sub in df.groupby(service_col, as_index=False, sort=True): if len(sub) >= MIN_PROVIDERS_PER_SERVICE: out.append((code, sub)) return out # --------------------------------------------------------------------------- # 4a. Project-defined high-charge cohort — top quartile per service # --------------------------------------------------------------------------- def high_charge_cohort_method(df: pd.DataFrame, service_col: str, desc_col: str, billing_col: str, payment_col: str, volume_col: str, dataset_name: str) -> dict: """Flag high-charge provider-service rows and compare the two cohorts.""" flags = [] cohort_rows = [] for code, sub in _per_service(df, service_col, billing_col, payment_col): q3 = sub[billing_col].quantile(HIGH_CHARGE_QUANTILE) if pd.isna(q3): continue desc = sub[desc_col].iloc[0] sub = sub.copy() sub["is_high_charge_flag"] = sub[billing_col] >= q3 sub["outlier_method"] = HIGH_CHARGE_METHOD sub["service_code"] = code sub["service_desc"] = desc sub["eligible_service"] = True sub["eligible_provider_count"] = int(sub["provider_ccn"].nunique()) sub["high_charge_threshold"] = q3 proxy = ( sub["charge_to_payment_ratio"] if "charge_to_payment_ratio" in sub.columns else sub[billing_col].div(sub[payment_col]).replace([np.inf, -np.inf], np.nan) ) sub["charge_to_payment_proxy"] = proxy flag_columns = [ "provider_ccn", "provider_name", "state", "census_region", "urban_rural", "service_code", "service_desc", "is_high_charge_flag", "outlier_method", "eligible_service", "eligible_provider_count", "high_charge_threshold", billing_col, payment_col, volume_col, "charge_to_payment_proxy", ] flags.append(sub[flag_columns].rename(columns={ billing_col: "submitted_charge", payment_col: "payment", volume_col: "volume", })) # Cohort aggregate row hi = sub[sub["is_high_charge_flag"]] lo = sub[~sub["is_high_charge_flag"]] if len(hi) == 0 or len(lo) == 0: continue high_proxy = _mean(hi["charge_to_payment_proxy"], 3) low_proxy = _mean(lo["charge_to_payment_proxy"], 3) cohort_rows.append(OrderedDict( dataset=dataset_name, period=DATASET_PERIODS.get(dataset_name, _period_from_frames(df)), service_code=code, service_desc=desc, n_providers=int(len(sub)), n_high_charge_flags=int(len(hi)), n_high_outliers=int(len(hi)), high_mean_charge=_mean(hi[billing_col], 2), low_mean_charge=_mean(lo[billing_col], 2), charge_premium_pct=_premium_pct(hi[billing_col].mean(), lo[billing_col].mean()), high_mean_charge_to_payment_proxy=high_proxy, low_mean_charge_to_payment_proxy=low_proxy, charge_to_payment_proxy_premium_pct=_premium_pct(high_proxy, low_proxy), proxy_selection_note=( "Post-selection descriptive premium: cohort selected on submitted charge " "and proxy numerator is submitted charge." ), high_mean_volume=_mean(hi[volume_col], 1), low_mean_volume=_mean(lo[volume_col], 1), volume_ratio=_rounded(hi[volume_col].mean() / lo[volume_col].mean(), 3) if lo[volume_col].mean() > 0 else None, selection_metric=billing_col, payment_measure=payment_col, high_charge_quantile=HIGH_CHARGE_QUANTILE, high_charge_threshold=_rounded(q3, 2), )) flags_df = pd.concat(flags, ignore_index=True) if flags else pd.DataFrame() flags_df = _stable_sort(flags_df, ["service_code", "provider_ccn"]) cohort_df = _stable_sort(pd.DataFrame(cohort_rows), ["service_code", "service_desc"]) return { "metadata": _analysis_metadata( DATASET_PERIODS.get(dataset_name, _period_from_frames(df)) ), "flags": flags_df, "cohort_comparison": cohort_df, } # --------------------------------------------------------------------------- # 4b. Statistical methods — IQR + z-score (per service) # --------------------------------------------------------------------------- def iqr_method(df: pd.DataFrame, service_col: str, billing_col: str) -> pd.DataFrame: """Per service, flag providers outside [Q1 - 1.5*IQR, Q3 + 1.5*IQR] on billing.""" flags = [] for code, sub in _per_service(df, service_col, billing_col, billing_col): q1 = sub[billing_col].quantile(0.25) q3 = sub[billing_col].quantile(0.75) iqr = q3 - q1 lo_thr = q1 - IQR_MULTIPLIER * iqr hi_thr = q3 + IQR_MULTIPLIER * iqr sub = sub.copy() sub["is_iqr_outlier"] = (sub[billing_col] < lo_thr) | (sub[billing_col] > hi_thr) sub["iqr_direction"] = np.where(sub[billing_col] > hi_thr, "high", np.where(sub[billing_col] < lo_thr, "low", "none")) flags.append(sub[sub["is_iqr_outlier"]][[ "provider_ccn", "provider_name", "state", "census_region", "urban_rural", service_col, billing_col, "is_iqr_outlier", "iqr_direction"]].assign(method="iqr")) if not flags: return pd.DataFrame() out = pd.concat(flags, ignore_index=True) out = out.rename(columns={service_col: "service_code", billing_col: "submitted_charge"}) return _stable_sort(out, ["service_code", "provider_ccn"]) def zscore_method(df: pd.DataFrame, service_col: str, billing_col: str, threshold: float = None) -> pd.DataFrame: """Per service, flag providers |z| > threshold on billing.""" threshold = threshold or ZSCORE_MODERATE flags = [] for code, sub in _per_service(df, service_col, billing_col, billing_col): mu = sub[billing_col].mean() sigma = sub[billing_col].std(ddof=1) if sigma == 0 or pd.isna(sigma): continue sub = sub.copy() sub["z"] = (sub[billing_col] - mu) / sigma sub["is_z_outlier"] = sub["z"].abs() > threshold sub["z_direction"] = np.where(sub["z"] > threshold, "high", np.where(sub["z"] < -threshold, "low", "none")) flags.append(sub[sub["is_z_outlier"]][[ "provider_ccn", "provider_name", "state", "census_region", "urban_rural", service_col, billing_col, "z", "is_z_outlier", "z_direction"]].assign(method=f"zscore_{threshold}")) if not flags: return pd.DataFrame() out = pd.concat(flags, ignore_index=True) out = out.rename(columns={service_col: "service_code", billing_col: "submitted_charge"}) return _stable_sort(out, ["service_code", "provider_ccn"]) def method_agreement_matrix(inp_high_charge_flags, inp_iqr_flags, inp_z_flags, out_high_charge_flags, out_iqr_flags, out_z_flags) -> dict: """Return provider-level agreement counts for the three screening methods.""" def _providers(*flag_dfs): sets = [set(df["provider_ccn"].unique()) if isinstance(df, pd.DataFrame) and not df.empty else set() for df in flag_dfs] return sets results = {} for ds, high_charge, iqr, z in [ ("inpatient", inp_high_charge_flags, inp_iqr_flags, inp_z_flags), ("outpatient", out_high_charge_flags, out_iqr_flags, out_z_flags), ]: high_charge = _flagged_rows(high_charge) high_s, i_s, z_s = _providers(high_charge, iqr, z) three = len(high_s & i_s & z_s) any_one = len(high_s | i_s | z_s) results[ds] = { "high_charge_only": len(high_s - i_s - z_s), "iqr_only": len(i_s - high_s - z_s), "z_only": len(z_s - high_s - i_s), "high_charge_and_iqr": len((high_s & i_s) - z_s), "high_charge_and_z": len((high_s & z_s) - i_s), "iqr_and_z": len((i_s & z_s) - high_s), "all_three": three, "any_method": any_one, "high_charge_method": HIGH_CHARGE_METHOD, } return {"metadata": _analysis_metadata(), **results} # --------------------------------------------------------------------------- # 4c. Concentration analysis # --------------------------------------------------------------------------- def concentration_analysis(flags_df: pd.DataFrame, dataset_name: str, top_n: int = 16) -> dict: """Rank services by provider-service high-charge flag counts.""" flags = _flagged_rows(flags_df) if flags.empty: counts = pd.DataFrame(columns=[ "service_code", "service_desc", "provider_service_flag_count", "share_of_total_provider_service_flags", "cumulative_share", ]) else: pair_columns = [column for column in ("provider_ccn", "service_code") if column in flags.columns] flags = flags.drop_duplicates(pair_columns) if len(pair_columns) == 2 else flags counts = ( flags.groupby(["service_code", "service_desc"], dropna=False, as_index=False) .size() .rename(columns={"size": "provider_service_flag_count"}) ) counts = _stable_sort( counts, ["provider_service_flag_count", "service_code", "service_desc"], ) counts = counts.sort_values( ["provider_service_flag_count", "service_code", "service_desc"], ascending=[False, True, True], kind="mergesort", na_position="last", ).reset_index(drop=True) total = int(counts["provider_service_flag_count"].sum()) counts["share_of_total_provider_service_flags"] = ( counts["provider_service_flag_count"] / total if total else 0.0 ) counts["cumulative_share"] = counts["share_of_total_provider_service_flags"].cumsum() total = int(counts["provider_service_flag_count"].sum()) if not counts.empty else 0 top16_share = ( float(round(counts["share_of_total_provider_service_flags"].head(16).sum(), 4)) if total else 0.0 ) crossing_target = ( int((counts["cumulative_share"] < CONCENTRATION_TARGET_SHARE).sum()) + 1 if total else 0 ) crossing_target = min(crossing_target, len(counts)) if counts.size else 0 cumulative_top16 = ( float(round(counts["cumulative_share"].iloc[15] if len(counts) >= 16 else counts["cumulative_share"].iloc[-1], 4)) if total else 0.0 ) top_services = [] if not counts.empty: preview = counts.head(top_n).copy() preview["share_of_total"] = preview["share_of_total_provider_service_flags"] preview["outlier_count"] = preview["provider_service_flag_count"] top_services = preview.to_dict(orient="records") return { "dataset": dataset_name, "metadata": _analysis_metadata(), "flag_unit": "provider-service pair", "flag_label": "high-charge flag", "total_flagged_provider_service_pairs": total, "total_outlier_events": total, "n_services_with_high_charge_flags": int(len(counts)), "n_services_with_outliers": int(len(counts)), "top16_provider_service_flag_share": top16_share, "our_top16_share": top16_share, "smallest_n_accounting_for_target_share": int(crossing_target), "smallest_n_accounting_for_40_pct": int(crossing_target), "target_share": CONCENTRATION_TARGET_SHARE, "cumulative_share_at_top16": cumulative_top16, "payment_dollar_benchmark_comparable": False, "oig_payment_dollar_benchmark": { "top16_share": 0.41, "unit": "payment dollars", "comparable": False, "comparison_status": "not_comparable", "used_for_calculation": False, "reason": ( "The historical benchmark is a payment-dollar share, while " "this output counts provider-service high-charge flags." ), }, "top_services": top_services, } # --------------------------------------------------------------------------- # 4d. Cross-dataset high-charge overlap # --------------------------------------------------------------------------- def _eligible_opportunity_info(frame: pd.DataFrame, service_col: str, billing_col: str, desc_col: str) -> dict: """Build eligible-service and provider-service opportunity denominators.""" empty_pairs = pd.DataFrame(columns=["provider_ccn", "service_code"]) if not isinstance(frame, pd.DataFrame) or frame.empty: return { "available": False, "eligible_codes": set(), "eligible_provider_set": set(), "provider_eligible_counts": pd.Series(dtype=int), "eligible_pairs": empty_pairs, "by_service": [], } required = ["provider_ccn", service_col, billing_col] if any(column not in frame.columns for column in required): return { "available": False, "eligible_codes": set(), "eligible_provider_set": set(), "provider_eligible_counts": pd.Series(dtype=int), "eligible_pairs": empty_pairs, "by_service": [], } working = frame.dropna(subset=["provider_ccn", service_col, billing_col]).copy() working = working.drop_duplicates(["provider_ccn", service_col]) service_counts = ( working.groupby(service_col, sort=True)["provider_ccn"] .nunique() .rename("eligible_provider_count") .reset_index() ) eligible_codes = set( service_counts.loc[ service_counts["eligible_provider_count"] >= MIN_PROVIDERS_PER_SERVICE, service_col, ].tolist() ) eligible = working[working[service_col].isin(eligible_codes)].copy() eligible = eligible.rename(columns={service_col: "service_code"}) provider_counts = ( eligible.groupby("provider_ccn")["service_code"].nunique() if not eligible.empty else pd.Series(dtype=int) ) by_service = [] if not eligible.empty: desc_by_service = ( working.groupby(service_col, sort=True)[desc_col].first() if desc_col in working.columns else pd.Series(dtype=object) ) by_service_df = ( eligible.groupby("service_code", sort=True) .agg( eligible_provider_count=("provider_ccn", "nunique"), eligible_service_opportunities=("provider_ccn", "size"), ) .reset_index() ) for _, row in by_service_df.iterrows(): by_service.append(OrderedDict( service_code=row["service_code"], service_desc=desc_by_service.get(row["service_code"]), eligible_provider_count=int(row["eligible_provider_count"]), eligible_service_opportunities=int(row["eligible_service_opportunities"]), )) return { "available": True, "eligible_codes": eligible_codes, "eligible_provider_set": set(provider_counts.index), "provider_eligible_counts": provider_counts, "eligible_pairs": eligible[["provider_ccn", "service_code"]].drop_duplicates(), "by_service": by_service, } def _flag_opportunity_info(flags_df: pd.DataFrame) -> dict: """Use a flagged table as a clearly labeled denominator fallback.""" empty = { "available": False, "eligible_codes": set(), "eligible_provider_set": set(), "provider_eligible_counts": pd.Series(dtype=int), "eligible_pairs": pd.DataFrame(columns=["provider_ccn", "service_code"]), "by_service": [], } if not isinstance(flags_df, pd.DataFrame) or flags_df.empty: return empty required = {"provider_ccn", "service_code"} if not required.issubset(flags_df.columns): return empty pairs = flags_df.dropna(subset=["provider_ccn", "service_code"])[ ["provider_ccn", "service_code"] ].drop_duplicates() if pairs.empty: return empty desc_by_service = ( flags_df.groupby("service_code", sort=True)["service_desc"].first() if "service_desc" in flags_df.columns else pd.Series(dtype=object) ) by_service_df = ( pairs.groupby("service_code", sort=True) .agg( eligible_provider_count=("provider_ccn", "nunique"), eligible_service_opportunities=("provider_ccn", "size"), ) .reset_index() ) by_service = [OrderedDict( service_code=row["service_code"], service_desc=desc_by_service.get(row["service_code"]), eligible_provider_count=int(row["eligible_provider_count"]), eligible_service_opportunities=int(row["eligible_service_opportunities"]), ) for _, row in by_service_df.iterrows()] provider_counts = pairs.groupby("provider_ccn")["service_code"].nunique() return { "available": False, "eligible_codes": set(pairs["service_code"]), "eligible_provider_set": set(provider_counts.index), "provider_eligible_counts": provider_counts, "eligible_pairs": pairs, "by_service": by_service, } def _universe_stats(frame: pd.DataFrame, service_col: str, payment_col: str, eligible_codes: set, available: bool) -> dict: """Summarize full and cost-observed provider-service universes.""" if not available or not isinstance(frame, pd.DataFrame): return {"available": False} if frame.empty or "provider_ccn" not in frame.columns or service_col not in frame.columns: return {"available": True, "rows": 0, "providers": 0, "services": 0, "provider_service_pairs": 0, "eligible_service_opportunities": 0, "eligible_providers": 0, "eligible_services": 0} pairs = frame.dropna(subset=["provider_ccn", service_col]).drop_duplicates( ["provider_ccn", service_col] ) eligible = pairs[pairs[service_col].isin(eligible_codes)] cost_observed = ( frame[frame[payment_col].notna()] if payment_col in frame.columns else None ) cost_pairs = ( cost_observed.dropna(subset=["provider_ccn", service_col]).drop_duplicates( ["provider_ccn", service_col] ) if cost_observed is not None else None ) return OrderedDict( available=True, rows=int(len(frame)), providers=int(frame["provider_ccn"].nunique()), services=int(frame[service_col].nunique()), provider_service_pairs=int(len(pairs)), eligible_services=int( pairs.loc[pairs[service_col].isin(eligible_codes), service_col].nunique() ), eligible_service_opportunities=int(len(eligible)), eligible_providers=int(eligible["provider_ccn"].nunique()), cost_observed_rows=(int(len(cost_observed)) if cost_observed is not None else None), cost_observed_provider_service_pairs=( int(len(cost_pairs)) if cost_pairs is not None else None ), ) def _breadth_bin(count: int) -> str: if count <= 1: return "1" if count <= 4: return "2-4" if count <= 9: return "5-9" return "10+" def _breadth_bins(breadth: pd.Series) -> list: counts = breadth.value_counts() if isinstance(breadth, pd.Series) else pd.Series(dtype=int) total = int(len(breadth)) if isinstance(breadth, pd.Series) else 0 rows = [] for label in BREADTH_BIN_LABELS: provider_count = int(sum( int(frequency) for breadth_value, frequency in counts.items() if _breadth_bin(int(breadth_value)) == label )) rows.append(OrderedDict( breadth_bin=label, provider_count=provider_count, share_of_flagged_providers=_rate(provider_count, total), )) return rows def _method_flags(flags_df: pd.DataFrame, method: str) -> pd.DataFrame: if not isinstance(flags_df, pd.DataFrame) or flags_df.empty: return flags_df.copy() if isinstance(flags_df, pd.DataFrame) else pd.DataFrame() if "outlier_method" in flags_df.columns: return flags_df[flags_df["outlier_method"] == method].copy() return flags_df.copy() def cross_dataset_outliers(inp_flags: pd.DataFrame, out_flags: pd.DataFrame, method: str = HIGH_CHARGE_METHOD, inp_frame: pd.DataFrame = None, out_frame: pd.DataFrame = None, inp_full_frame: pd.DataFrame = None, out_full_frame: pd.DataFrame = None, metadata: dict = None, ) -> dict: """Measure provider overlap using high-charge flags and full denominators.""" inp_method_flags = _method_flags(inp_flags, method) out_method_flags = _method_flags(out_flags, method) inp_sub = _flagged_rows(inp_method_flags) out_sub = _flagged_rows(out_method_flags) inp_info = _eligible_opportunity_info( inp_frame, "drg_code", "avg_submitted_charge", "drg_desc" ) out_info = _eligible_opportunity_info( out_frame, "apc_code", "avg_submitted_charge", "apc_desc" ) # When raw frames are unavailable, retain a transparent flag-only fallback. inp_frame_available = inp_info["available"] out_frame_available = out_info["available"] if not inp_info["available"]: inp_info = _eligible_opportunity_info( inp_method_flags, "service_code", "submitted_charge", "service_desc" ) if not inp_info["available"]: inp_info = _flag_opportunity_info(inp_sub) if not out_info["available"]: out_info = _eligible_opportunity_info( out_method_flags, "service_code", "submitted_charge", "service_desc" ) if not out_info["available"]: out_info = _flag_opportunity_info(out_sub) inp_denominator_basis = ( "eligible provider-service rows in the input frame" if inp_frame_available else ( "complete cohort flag table" if inp_info["available"] else "flag table fallback; input frame unavailable" ) ) out_denominator_basis = ( "eligible provider-service rows in the input frame" if out_frame_available else ( "complete cohort flag table" if out_info["available"] else "flag table fallback; input frame unavailable" ) ) inp_p = set(inp_sub["provider_ccn"].dropna().unique()) if not inp_sub.empty else set() out_p = set(out_sub["provider_ccn"].dropna().unique()) if not out_sub.empty else set() common = inp_p & out_p inp_pairs = (inp_sub[["provider_ccn", "service_code"]].drop_duplicates() if not inp_sub.empty else pd.DataFrame(columns=["provider_ccn", "service_code"])) out_pairs = (out_sub[["provider_ccn", "service_code"]].drop_duplicates() if not out_sub.empty else pd.DataFrame(columns=["provider_ccn", "service_code"])) inp_breadth = (inp_pairs.groupby("provider_ccn")["service_code"].nunique() if not inp_pairs.empty else pd.Series(dtype=int)) out_breadth = (out_pairs.groupby("provider_ccn")["service_code"].nunique() if not out_pairs.empty else pd.Series(dtype=int)) inp_eligible_provider_count = len(inp_info["eligible_provider_set"]) out_eligible_provider_count = len(out_info["eligible_provider_set"]) inp_flagged_provider_rate = _rate(len(inp_p), inp_eligible_provider_count) out_flagged_provider_rate = _rate(len(out_p), out_eligible_provider_count) eligible_both = ( inp_info["eligible_provider_set"] & out_info["eligible_provider_set"] ) marginal_expected_overlap = ( len(eligible_both) * inp_flagged_provider_rate * out_flagged_provider_rate ) shared_inp_flagged = inp_p & eligible_both shared_out_flagged = out_p & eligible_both shared_inp_rate = _rate(len(shared_inp_flagged), len(eligible_both)) shared_out_rate = _rate(len(shared_out_flagged), len(eligible_both)) expected_overlap = len(eligible_both) * shared_inp_rate * shared_out_rate breadth_expected_overlap, breadth_detail = _breadth_matched_expected_overlap( eligible_both, inp_p, out_p, inp_info["provider_eligible_counts"], out_info["provider_eligible_counts"], ) inp_universe = OrderedDict( full=_universe_stats( inp_full_frame, "drg_code", "avg_total_payment", inp_info["eligible_codes"], inp_full_frame is not None, ), cost_observed=_universe_stats( inp_frame, "drg_code", "avg_total_payment", inp_info["eligible_codes"], inp_frame is not None, ), ) out_universe = OrderedDict( full=_universe_stats( out_full_frame, "apc_code", "avg_allowed_amount", out_info["eligible_codes"], out_full_frame is not None, ), cost_observed=_universe_stats( out_frame, "apc_code", "avg_allowed_amount", out_info["eligible_codes"], out_frame is not None, ), ) # Build provider metadata keyed by provider_ccn. meta_frames = [] for frame in (inp_sub, out_sub): if frame.empty or "provider_ccn" not in frame.columns: continue columns = [ column for column in ( "provider_name", "state", "census_region", "urban_rural", ) if column in frame.columns ] if columns: meta_frames.append(frame.drop_duplicates("provider_ccn").set_index("provider_ccn")[columns]) meta = ( pd.concat(meta_frames).groupby(level=0).first() if meta_frames else pd.DataFrame() ) inp_provider_eligible = inp_info["provider_eligible_counts"] out_provider_eligible = out_info["provider_eligible_counts"] rows = [] for ccn in sorted(common, key=str): inp_flag_count = int(inp_breadth.get(ccn, 0)) out_flag_count = int(out_breadth.get(ccn, 0)) inp_eligible_count = int(inp_provider_eligible.get(ccn, inp_flag_count)) out_eligible_count = int(out_provider_eligible.get(ccn, out_flag_count)) row = OrderedDict( provider_ccn=ccn, provider_name=meta.loc[ccn, "provider_name"] if ccn in meta.index and "provider_name" in meta.columns else None, state=meta.loc[ccn, "state"] if ccn in meta.index and "state" in meta.columns else None, census_region=meta.loc[ccn, "census_region"] if ccn in meta.index and "census_region" in meta.columns else None, urban_rural=meta.loc[ccn, "urban_rural"] if ccn in meta.index and "urban_rural" in meta.columns else None, inpatient_high_charge_services_count=inp_flag_count, outpatient_high_charge_services_count=out_flag_count, inpatient_eligible_services_count=inp_eligible_count, outpatient_eligible_services_count=out_eligible_count, inpatient_provider_flag_rate=_rate(inp_flag_count, inp_eligible_count), outpatient_provider_flag_rate=_rate(out_flag_count, out_eligible_count), inpatient_breadth_bin=_breadth_bin(inp_flag_count), outpatient_breadth_bin=_breadth_bin(out_flag_count), total_high_charge_services_count=inp_flag_count + out_flag_count, # Existing generic count fields remain useful to current consumers. inp_outlier_services_count=inp_flag_count, out_outlier_services_count=out_flag_count, ) rows.append(row) rows = sorted( rows, key=lambda row: (-row["total_high_charge_services_count"], str(row["provider_ccn"])), ) inp_summary = OrderedDict( flagged_provider_service_pairs=int(len(inp_pairs)), eligible_service_opportunities=int(len(inp_info["eligible_pairs"])), provider_service_flag_rate=_rate(len(inp_pairs), len(inp_info["eligible_pairs"])), provider_service_flag_rate_pct=float(round( 100 * _rate(len(inp_pairs), len(inp_info["eligible_pairs"])), 2 )), denominator_basis=inp_denominator_basis, eligible_service_count=int(len(inp_info["eligible_codes"])), flagged_provider_count=int(len(inp_p)), eligible_provider_count=int(inp_eligible_provider_count), provider_flag_rate=inp_flagged_provider_rate, provider_flag_rate_pct=float(round(100 * inp_flagged_provider_rate, 2)), ) out_summary = OrderedDict( flagged_provider_service_pairs=int(len(out_pairs)), eligible_service_opportunities=int(len(out_info["eligible_pairs"])), provider_service_flag_rate=_rate(len(out_pairs), len(out_info["eligible_pairs"])), provider_service_flag_rate_pct=float(round( 100 * _rate(len(out_pairs), len(out_info["eligible_pairs"])), 2 )), denominator_basis=out_denominator_basis, eligible_service_count=int(len(out_info["eligible_codes"])), flagged_provider_count=int(len(out_p)), eligible_provider_count=int(out_eligible_provider_count), provider_flag_rate=out_flagged_provider_rate, provider_flag_rate_pct=float(round(100 * out_flagged_provider_rate, 2)), ) expected_detail = OrderedDict( benchmark="conditional independence within the shared eligible-provider universe", eligible_provider_overlap=int(len(eligible_both)), inpatient_provider_flag_rate=shared_inp_rate, outpatient_provider_flag_rate=shared_out_rate, expected_flagged_provider_overlap=float(round(expected_overlap, 4)), marginal_expected_flagged_provider_overlap=float(round(marginal_expected_overlap, 4)), breadth_matched_expected_flagged_provider_overlap=float( round(breadth_expected_overlap, 4) ), observed_flagged_provider_overlap=int(len(common)), observed_to_expected_ratio=( float(round(len(common) / expected_overlap, 4)) if expected_overlap else None ), observed_to_breadth_matched_expected_ratio=( float(round(len(common) / breadth_expected_overlap, 4)) if breadth_expected_overlap else None ), assumption=( "Within the shared eligible-provider universe, the two flag indicators " "are treated as independent; service-opportunity breadth is not matched." ), breadth_matched_benchmark=breadth_detail, ) output_metadata = metadata or _analysis_metadata( _period_from_frames(inp_frame, out_frame, inp_full_frame, out_full_frame) ) return { "metadata": output_metadata, "summary": { "method": method, "high_charge_method": HIGH_CHARGE_METHOD, "inpatient_flagged_providers": int(len(inp_p)), "outpatient_flagged_providers": int(len(out_p)), "inpatient_eligible_providers": int(inp_eligible_provider_count), "outpatient_eligible_providers": int(out_eligible_provider_count), "inpatient_provider_flag_rate": inp_flagged_provider_rate, "outpatient_provider_flag_rate": out_flagged_provider_rate, "shared_eligible_provider_count": int(len(eligible_both)), "shared_inpatient_provider_flag_rate": shared_inp_rate, "shared_outpatient_provider_flag_rate": shared_out_rate, "inpatient_flagged_provider_service_pairs": int(len(inp_pairs)), "outpatient_flagged_provider_service_pairs": int(len(out_pairs)), "inpatient_eligible_service_opportunities": int(len(inp_info["eligible_pairs"])), "outpatient_eligible_service_opportunities": int(len(out_info["eligible_pairs"])), "inpatient_outlier_providers": int(len(inp_p)), "outpatient_outlier_providers": int(len(out_p)), "outliers_in_both": int(len(common)), "pct_of_inpatient": float(round(100 * len(common) / max(1, len(inp_p)), 2)), "pct_of_outpatient": float(round(100 * len(common) / max(1, len(out_p)), 2)), "expected_overlap_under_independence": float(round(expected_overlap, 4)), "marginal_expected_overlap_under_independence": float( round(marginal_expected_overlap, 4) ), "expected_overlap_under_breadth_matched_independence": float( round(breadth_expected_overlap, 4) ), "expected_overlap_under_independence_detail": expected_detail, }, "eligible_service_opportunity_denominators": { "inpatient": inp_info["by_service"], "outpatient": out_info["by_service"], }, "dataset_denominators": { "inpatient": inp_summary, "outpatient": out_summary, }, "breadth_bins": { "inpatient": _breadth_bins(inp_breadth), "outpatient": _breadth_bins(out_breadth), }, "universe_metadata": { "inpatient": inp_universe, "outpatient": out_universe, }, "top_cross_dataset_outliers": rows[:50], } # --------------------------------------------------------------------------- # 4e. Outlier profiling # --------------------------------------------------------------------------- def profile_outliers(inp_flags: pd.DataFrame, out_flags: pd.DataFrame) -> dict: """Geographic + breadth profile of outlier providers per dataset. For each dataset we report, by census_region and urban_rural, the count of distinct outlier providers + their mean number of services they're outliers in. """ def _prof(flags_df, ds_name): flags_df = _flagged_rows(flags_df) if flags_df.empty: return {"dataset": ds_name, "regions": [], "urban_rural": [], "breadth": {}} # Distinct providers per (geography) per_prov_breadth = flags_df.groupby(["provider_ccn", "census_region", "urban_rural"])["service_code"].nunique().reset_index(name="service_count") region_agg = (per_prov_breadth.groupby("census_region") .agg(n_outlier_providers=("provider_ccn", "nunique"), median_services_outlier_in=("service_count", "median")) .reset_index()) urban_agg = (per_prov_breadth.groupby("urban_rural") .agg(n_outlier_providers=("provider_ccn", "nunique"), median_services_outlier_in=("service_count", "median")) .reset_index()) breadth = per_prov_breadth["service_count"] return { "dataset": ds_name, "regions": region_agg.to_dict(orient="records"), "urban_rural": urban_agg.to_dict(orient="records"), "breadth": { "median_services_outlier_in": float(round(breadth.median(), 2)), "p90_services_outlier_in": float(round(breadth.quantile(0.9), 2)), "max_services_outlier_in": int(breadth.max()), "n_providers_outlier_in_one_service_only": int((breadth == 1).sum()), "n_providers_outlier_in_5plus_services": int((breadth >= 5).sum()), }, } return { "metadata": _analysis_metadata(), "inpatient": _prof(inp_flags, "inpatient"), "outpatient": _prof(out_flags, "outpatient"), } def _audit_frame(frame: pd.DataFrame, dataset: str) -> pd.DataFrame: """Add dataset identity and deterministic ordering to an audit table.""" if not isinstance(frame, pd.DataFrame): frame = pd.DataFrame() frame = frame.copy() if "dataset" in frame.columns: frame["dataset"] = dataset else: frame.insert(0, "dataset", dataset) return _stable_sort(frame, ["dataset", "service_code", "provider_ccn"]) def _write_audit_table(frame: pd.DataFrame, stem: str) -> str: """Write a complete audit table to processed data, with CSV fallback.""" parquet_path = DATA_PROC / f"{stem}.parquet" try: frame.to_parquet(parquet_path, index=False) return str(parquet_path) except Exception: csv_path = DATA_PROC / f"{stem}.csv" frame.to_csv(csv_path, index=False) return str(csv_path) def _preview_metadata(metadata: dict, full_counts: dict, preview_counts: dict, limit: int, sort_order: list) -> OrderedDict: """Attach explicit size and ordering information to a site preview.""" out = OrderedDict(metadata) out.update( preview=True, truncated=any( preview_counts.get(dataset, 0) < full_counts.get(dataset, 0) for dataset in full_counts ), preview_limit_per_dataset=int(limit), full_row_counts={key: int(value) for key, value in full_counts.items()}, preview_row_counts={key: int(value) for key, value in preview_counts.items()}, sort_order=sort_order, ) return out # --------------------------------------------------------------------------- # Orchestration # --------------------------------------------------------------------------- def run() -> dict: print("[Phase 4] Loading feature parquets ...") inp = pd.read_parquet(DATA_PROC / "inpatient_features.parquet") out = pd.read_parquet(DATA_PROC / "outpatient_features.parquet") try: out_full = pd.read_parquet(DATA_PROC / "outpatient_full.parquet") except FileNotFoundError: out_full = None print(f" inpatient: {len(inp):,} rows | outpatient: {len(out):,} rows") output_metadata = _analysis_metadata() # 4a — project-defined high-charge cohort (top quartile per service) print("[Phase 4] 4a — Project-defined high-charge cohort (top quartile per service) ...") inp_high_charge = high_charge_cohort_method( inp, "drg_code", "drg_desc", "avg_submitted_charge", "avg_total_payment", "total_discharges", "inpatient" ) out_high_charge = high_charge_cohort_method( out, "apc_code", "apc_desc", "avg_submitted_charge", "avg_allowed_amount", "apc_services", "outpatient" ) inp_high_charge_flags = _flagged_rows(inp_high_charge["flags"]) out_high_charge_flags = _flagged_rows(out_high_charge["flags"]) n_inp_flags = len(inp_high_charge_flags) n_out_flags = len(out_high_charge_flags) print(f" inpatient cohort: {n_inp_flags:,} provider-service high-charge flags") print(f" outpatient cohort: {n_out_flags:,} provider-service high-charge flags") # Cohort comparison headline inp_cohort = inp_high_charge["cohort_comparison"] out_cohort = out_high_charge["cohort_comparison"] median_inp_premium = inp_cohort["charge_premium_pct"].median() if not inp_cohort.empty else None median_out_premium = out_cohort["charge_premium_pct"].median() if not out_cohort.empty else None median_inp_proxy_prem = ( inp_cohort["charge_to_payment_proxy_premium_pct"].median() if not inp_cohort.empty else None ) median_out_proxy_prem = ( out_cohort["charge_to_payment_proxy_premium_pct"].median() if not out_cohort.empty else None ) print(f" inpatient high-charge charge premium (median): {median_inp_premium}%") print(f" outpatient high-charge charge premium (median): {median_out_premium}%") print(f" inpatient charge-to-payment proxy premium (median): {median_inp_proxy_prem}%") print(f" outpatient charge-to-payment proxy premium (median): {median_out_proxy_prem}%") # 4b — IQR + z-score methods print("[Phase 4] 4b — IQR + z-score (per service) ...") inp_iqr_flags = iqr_method(inp, "drg_code", "avg_submitted_charge") out_iqr_flags = iqr_method(out, "apc_code", "avg_submitted_charge") inp_z2_flags = zscore_method(inp, "drg_code", "avg_submitted_charge", ZSCORE_MODERATE) out_z2_flags = zscore_method(out, "apc_code", "avg_submitted_charge", ZSCORE_MODERATE) print(f" inpatient: IQR {len(inp_iqr_flags):,} | z>2 {len(inp_z2_flags):,}") print(f" outpatient: IQR {len(out_iqr_flags):,} | z>2 {len(out_z2_flags):,}") agreement = method_agreement_matrix( inp_high_charge_flags, inp_iqr_flags, inp_z2_flags, out_high_charge_flags, out_iqr_flags, out_z2_flags, ) agreement["metadata"] = output_metadata print(f" 3-way agreement (all three methods): " f"inpatient={agreement['inpatient']['all_three']} | " f"outpatient={agreement['outpatient']['all_three']}") print(f" total inpatient providers flagged by any method: {agreement['inpatient']['any_method']}") # 4c — Concentration analysis print("[Phase 4] 4c — Concentration analysis (provider-service high-charge flags) ...") inp_conc = concentration_analysis(inp_high_charge_flags, "inpatient") out_conc = concentration_analysis(out_high_charge_flags, "outpatient") inp_conc["metadata"] = output_metadata out_conc["metadata"] = output_metadata print(f" inpatient: top 16 provider-service flag share = " f"{inp_conc['top16_provider_service_flag_share'] * 100:.1f}% | " f"smallest N for 40% = {inp_conc['smallest_n_accounting_for_target_share']}") print(f" outpatient: top 16 provider-service flag share = " f"{out_conc['top16_provider_service_flag_share'] * 100:.1f}% | " f"smallest N for 40% = {out_conc['smallest_n_accounting_for_target_share']}") print(" historical payment-dollar benchmark is not comparable to these flag counts") # 4d — Cross-dataset high-charge overlap print("[Phase 4] 4d — Cross-dataset high-charge flag overlap ...") cross = cross_dataset_outliers( inp_high_charge_flags, out_high_charge_flags, method=HIGH_CHARGE_METHOD, inp_frame=inp, out_frame=out, inp_full_frame=inp, out_full_frame=out_full, metadata=output_metadata, ) print(f" inpatient flagged providers: {cross['summary']['inpatient_flagged_providers']:,}") print(f" outpatient flagged providers: {cross['summary']['outpatient_flagged_providers']:,}") print(f" providers flagged in BOTH: {cross['summary']['outliers_in_both']:,} " f"({cross['summary']['pct_of_inpatient']}% of inpatient / " f"{cross['summary']['pct_of_outpatient']}% of outpatient)") print(f" expected overlap under independence: " f"{cross['summary']['expected_overlap_under_independence']:.2f}") # 4e — Profiling print("[Phase 4] 4e — High-charge provider profiling (geography, breadth) ...") profile = profile_outliers(inp_high_charge_flags, out_high_charge_flags) print(f" inpatient high-charge breadth: median " f"{profile['inpatient']['breadth']['median_services_outlier_in']} services " f"(max {profile['inpatient']['breadth']['max_services_outlier_in']})") # Serialize print("[Phase 4] Serializing to site/data/ ...") cohort_sort = ["service_code", "service_desc"] cohort_preview_inp = _stable_sort(inp_cohort, cohort_sort) cohort_preview_out = _stable_sort(out_cohort, cohort_sort) cohort_full_counts = { "inpatient": len(cohort_preview_inp), "outpatient": len(cohort_preview_out), } cohort_preview_counts = { "inpatient": min(COHORT_PREVIEW_LIMIT, len(cohort_preview_inp)), "outpatient": min(COHORT_PREVIEW_LIMIT, len(cohort_preview_out)), } save_json({ "metadata": _preview_metadata( output_metadata, cohort_full_counts, cohort_preview_counts, COHORT_PREVIEW_LIMIT, cohort_sort, ), "inpatient": _records(cohort_preview_inp, COHORT_PREVIEW_LIMIT), "outpatient": _records(cohort_preview_out, COHORT_PREVIEW_LIMIT), "summary": { "inpatient_median_charge_premium_pct": _rounded(median_inp_premium, 2), "outpatient_median_charge_premium_pct": _rounded(median_out_premium, 2), "inpatient_median_charge_to_payment_proxy_premium_pct": _rounded( median_inp_proxy_prem, 2 ), "outpatient_median_charge_to_payment_proxy_premium_pct": _rounded( median_out_proxy_prem, 2 ), }, }, "outlier_cohort_comparison.json") save_json({**agreement, "metadata": output_metadata}, "outlier_method_agreement.json") save_json({ "metadata": output_metadata, "inpatient": inp_conc, "outpatient": out_conc, }, "outlier_concentration.json") save_json(cross, "cross_dataset_outliers.json") profile["metadata"] = output_metadata save_json(profile, "outlier_profile.json") flag_columns = [ "provider_ccn", "provider_name", "state", "census_region", "urban_rural", "service_code", "service_desc", "is_high_charge_flag", "outlier_method", "eligible_service", "eligible_provider_count", "high_charge_threshold", "submitted_charge", "payment", "volume", "charge_to_payment_proxy", ] inp_site_flags = _stable_sort(inp_high_charge_flags, ["service_code", "provider_ccn"]) out_site_flags = _stable_sort(out_high_charge_flags, ["service_code", "provider_ccn"]) inp_site_flags = inp_site_flags[[column for column in flag_columns if column in inp_site_flags.columns]] out_site_flags = out_site_flags[[column for column in flag_columns if column in out_site_flags.columns]] flag_full_counts = { "inpatient": len(inp_site_flags), "outpatient": len(out_site_flags), } flag_preview_counts = { "inpatient": min(FLAG_PREVIEW_LIMIT, len(inp_site_flags)), "outpatient": min(FLAG_PREVIEW_LIMIT, len(out_site_flags)), } save_json({ "metadata": _preview_metadata( output_metadata, flag_full_counts, flag_preview_counts, FLAG_PREVIEW_LIMIT, ["service_code", "provider_ccn"], ), "inpatient": _records(inp_site_flags, FLAG_PREVIEW_LIMIT), "outpatient": _records(out_site_flags, FLAG_PREVIEW_LIMIT), }, "outlier_flags.json") # Complete tables stay outside the website preview for auditability. high_charge_audit = pd.concat([ _audit_frame(inp_high_charge["flags"], "inpatient"), _audit_frame(out_high_charge["flags"], "outpatient"), ], ignore_index=True, sort=False) cohort_audit = pd.concat([ _audit_frame(inp_cohort, "inpatient"), _audit_frame(out_cohort, "outpatient"), ], ignore_index=True, sort=False) iqr_audit = pd.concat([ _audit_frame(inp_iqr_flags, "inpatient"), _audit_frame(out_iqr_flags, "outpatient"), ], ignore_index=True, sort=False) zscore_audit = pd.concat([ _audit_frame(inp_z2_flags, "inpatient"), _audit_frame(out_z2_flags, "outpatient"), ], ignore_index=True, sort=False) audit_artifacts = OrderedDict( outlier_flags_audit=_write_audit_table(high_charge_audit, "outlier_flags_audit"), outlier_cohort_comparison_audit=_write_audit_table( cohort_audit, "outlier_cohort_comparison_audit" ), outlier_iqr_flags_audit=_write_audit_table(iqr_audit, "outlier_iqr_flags_audit"), outlier_zscore_flags_audit=_write_audit_table( zscore_audit, "outlier_zscore_flags_audit" ), ) print("[Phase 4] DONE.") return { "cohort_comparison": {"inpatient": inp_cohort, "outpatient": out_cohort}, "agreement": agreement, "concentration": {"inpatient": inp_conc, "outpatient": out_conc}, "cross_dataset": cross, "profile": profile, "metadata": output_metadata, "audit_artifacts": audit_artifacts, } if __name__ == "__main__": # pragma: no cover run()