{"spec_id":"probability-weibull","library":"plotnine","language":"python","code":"\"\"\" anyplot.ai\nprobability-weibull: Weibull Probability Plot for Reliability Analysis\nLibrary: plotnine 0.15.5 | Python 3.13.13\nQuality: 91/100 | Updated: 2026-06-07\n\"\"\"\n\nimport os\nimport sys\n\nimport numpy as np\nimport pandas as pd\n\n\n# Work around naming conflict between plotnine.py script and plotnine package\nscript_dir = os.path.dirname(os.path.abspath(__file__))\nif script_dir in sys.path:\n    sys.path.remove(script_dir)\nif \"\" in sys.path:\n    sys.path.remove(\"\")\nif \".\" in sys.path:\n    sys.path.remove(\".\")\n\nfrom plotnine import (\n    aes,\n    annotate,\n    element_blank,\n    element_line,\n    element_rect,\n    element_text,\n    geom_line,\n    geom_point,\n    geom_ribbon,\n    geom_segment,\n    ggplot,\n    guide_legend,\n    guides,\n    labs,\n    scale_color_manual,\n    scale_fill_manual,\n    scale_shape_manual,\n    scale_x_log10,\n    scale_y_continuous,\n    theme,\n    theme_minimal,\n)\nfrom scipy import stats\n\n\n# Theme-adaptive chrome tokens (Imprint palette)\nTHEME = os.getenv(\"ANYPLOT_THEME\", \"light\")\nPAGE_BG = \"#FAF8F1\" if THEME == \"light\" else \"#1A1A17\"\nELEVATED_BG = \"#FFFDF6\" if THEME == \"light\" else \"#242420\"\nINK = \"#1A1A17\" if THEME == \"light\" else \"#F0EFE8\"\nINK_SOFT = \"#4A4A44\" if THEME == \"light\" else \"#B8B7B0\"\nINK_MUTED = \"#6B6A63\" if THEME == \"light\" else \"#A8A79F\"\n\n# Imprint categorical palette — canonical hybrid-v3 order\nIMPRINT_PALETTE = [\"#009E73\", \"#C475FD\", \"#4467A3\", \"#BD8233\", \"#AE3030\", \"#2ABCCD\", \"#954477\", \"#99B314\"]\nFIT_COLOR = IMPRINT_PALETTE[2]  # blue — Weibull regression line & band\nCHAR_COLOR = IMPRINT_PALETTE[4]  # matte red — characteristic life focal point\n\n# Data — turbine blade fatigue-life: 25 failures + 7 censored\nnp.random.seed(42)\nn_failures = 25\nn_censored = 7\n\nfailure_times = np.sort(stats.weibull_min.rvs(c=2.8, scale=5000, size=n_failures))\ncensored_times = np.sort(stats.uniform.rvs(loc=1000, scale=5000, size=n_censored))\n\nall_times = np.concatenate([failure_times, censored_times])\nis_failure = np.array([True] * n_failures + [False] * n_censored)\n\nsort_idx = np.argsort(all_times)\nall_times = all_times[sort_idx]\nis_failure = is_failure[sort_idx]\n\nfailure_order = np.cumsum(is_failure)\nmedian_rank = np.where(is_failure, (failure_order - 0.3) / (n_failures + 0.4), np.nan)\nweibull_y = np.where(is_failure, np.log(-np.log(1 - median_rank)), np.nan)\n\ndf = pd.DataFrame({\"time\": all_times, \"weibull_y\": weibull_y, \"status\": np.where(is_failure, \"Failure\", \"Censored\")})\n\n# Weibull regression on ln(time) vs ln(-ln(1-F))\nfailures_df = df[df[\"status\"] == \"Failure\"].dropna()\nlog_times = np.log(failures_df[\"time\"].values)\nslope, intercept, r_value, _, _ = stats.linregress(log_times, failures_df[\"weibull_y\"].values)\nbeta = slope\neta = np.exp(-intercept / beta)\n\n# Fitted line and 95% prediction band (in raw time units for scale_x_log10)\nlog_time_range = np.linspace(np.log(all_times.min() * 0.7), np.log(all_times.max() * 1.3), 200)\nfitted_y = beta * log_time_range + intercept\n\nn_fit = len(failures_df)\nx_mean = log_times.mean()\nx_ss = ((log_times - x_mean) ** 2).sum()\nresiduals = failures_df[\"weibull_y\"].values - (beta * log_times + intercept)\nmse = (residuals**2).sum() / (n_fit - 2)\nse_pred = np.sqrt(mse * (1 + 1 / n_fit + (log_time_range - x_mean) ** 2 / x_ss))\nt_crit = stats.t.ppf(0.975, n_fit - 2)\n\nfit_df = pd.DataFrame(\n    {\n        \"time\": np.exp(log_time_range),\n        \"weibull_y\": fitted_y,\n        \"ymin\": fitted_y - t_crit * se_pred,\n        \"ymax\": fitted_y + t_crit * se_pred,\n    }\n)\n\n# Y-axis: Weibull probability scale (linearised CDF)\nprob_levels = np.array([0.01, 0.05, 0.10, 0.20, 0.40, 0.632, 0.80, 0.90, 0.99])\nweibull_ticks = np.log(-np.log(1 - prob_levels))\nprob_labels = [f\"{p * 100:.1f}%\".replace(\".0%\", \"%\") for p in prob_levels]\n\n# X-axis breaks (raw time, scale_x_log10 handles the transform)\nx_tick_values = [1000, 2000, 3000, 5000, 7000, 10000]\nx_labels = [f\"{v:,}\" for v in x_tick_values]\n\n# Reference values\nref_y = np.log(-np.log(1 - 0.632))  # ≈ 0: characteristic life on Weibull y-scale\n\n# Scatter dataframes (raw time units)\nfailures_plot = failures_df[[\"time\", \"weibull_y\"]].assign(Status=\"Failure\")\ncensored_plot = df[df[\"status\"] == \"Censored\"][[\"time\"]].assign(weibull_y=weibull_ticks[0] + 0.18, Status=\"Censored\")\nscatter_df = pd.concat([failures_plot, censored_plot], ignore_index=True)\n\nhighlight_df = pd.DataFrame({\"time\": [eta], \"weibull_y\": [ref_y]})\n\n# Annotation positions in data space (time units for x; Weibull-y for y)\n# Shift further right to avoid crowding with high-probability scatter points\nannot_time = eta * 1.55\nannot_y = ref_y + 0.42\n\nplot = (\n    ggplot()\n    # 95% prediction band\n    + geom_ribbon(fit_df, aes(x=\"time\", ymin=\"ymin\", ymax=\"ymax\"), fill=FIT_COLOR, alpha=0.12)\n    # Weibull regression line\n    + geom_line(fit_df, aes(x=\"time\", y=\"weibull_y\"), color=FIT_COLOR, size=0.9, alpha=0.85)\n    # 63.2% horizontal reference (dashed guide)\n    + geom_segment(\n        aes(x=x_tick_values[0] * 0.75, xend=eta, y=ref_y, yend=ref_y),\n        linetype=\"dashed\",\n        color=INK_MUTED,\n        size=0.35,\n        alpha=0.65,\n    )\n    # Vertical drop guide to x-axis\n    + geom_segment(\n        aes(x=eta, xend=eta, y=weibull_ticks[0] - 0.1, yend=ref_y),\n        linetype=\"dashed\",\n        color=INK_MUTED,\n        size=0.35,\n        alpha=0.65,\n    )\n    # Failure and censored data points\n    + geom_point(\n        scatter_df, aes(x=\"time\", y=\"weibull_y\", color=\"Status\", shape=\"Status\"), size=2.8, alpha=0.82, stroke=0.3\n    )\n    # Characteristic life diamond marker\n    + geom_point(\n        highlight_df,\n        aes(x=\"time\", y=\"weibull_y\"),\n        color=CHAR_COLOR,\n        fill=CHAR_COLOR,\n        size=4.5,\n        shape=\"D\",\n        alpha=0.95,\n        show_legend=False,\n    )\n    # Eta annotation label\n    + annotate(\n        \"text\",\n        x=annot_time,\n        y=annot_y,\n        label=f\"η = {eta:,.0f} hrs (63.2%)\",\n        size=3.5,\n        ha=\"left\",\n        color=CHAR_COLOR,\n        fontweight=\"bold\",\n        fontstyle=\"italic\",\n    )\n    # Thin connector from annotation to diamond\n    + geom_segment(\n        aes(x=annot_time * 0.98, xend=eta * 1.02, y=annot_y - 0.10, yend=ref_y + 0.05),\n        color=CHAR_COLOR,\n        size=0.28,\n        alpha=0.50,\n    )\n    # Weibull parameters summary (upper-left)\n    + annotate(\n        \"text\",\n        x=x_tick_values[0] * 1.06,\n        y=weibull_ticks[-2] + 0.22,\n        label=f\"β = {beta:.2f}  ·  η = {eta:,.0f} hrs  ·  R² = {r_value**2:.3f}\",\n        size=3.5,\n        ha=\"left\",\n        va=\"top\",\n        color=FIT_COLOR,\n        fontweight=\"bold\",\n    )\n    # Scales — Imprint palette: failures=green (first), censored=lavender (second)\n    + scale_color_manual(\n        values={\"Failure\": IMPRINT_PALETTE[0], \"Censored\": IMPRINT_PALETTE[1]},\n        name=\"Observation\",\n        limits=[\"Failure\", \"Censored\"],\n    )\n    + scale_shape_manual(values={\"Failure\": \"o\", \"Censored\": \"^\"}, name=\"Observation\", limits=[\"Failure\", \"Censored\"])\n    + scale_fill_manual(\n        values={\"Failure\": IMPRINT_PALETTE[0], \"Censored\": IMPRINT_PALETTE[1]},\n        name=\"Observation\",\n        limits=[\"Failure\", \"Censored\"],\n    )\n    # scale_x_log10: idiomatic plotnine/ggplot2 log axis (no manual ln transform needed)\n    + scale_x_log10(breaks=x_tick_values, labels=x_labels)\n    + scale_y_continuous(breaks=weibull_ticks.tolist(), labels=prob_labels)\n    + guides(color=guide_legend(override_aes={\"size\": 3.5, \"alpha\": 1}))\n    + labs(\n        x=\"Time to Failure (hours)\",\n        y=\"Cumulative Failure Probability\",\n        title=\"probability-weibull · plotnine · anyplot.ai\",\n        caption=\"Turbine blade fatigue-life — 25 failures, 7 censored observations\",\n    )\n    + theme_minimal(base_family=\"sans-serif\")\n    + theme(\n        figure_size=(8, 4.5),\n        plot_title=element_text(size=12, weight=\"bold\", color=INK, margin={\"b\": 4}),\n        plot_caption=element_text(size=7.5, color=INK_MUTED, fontstyle=\"italic\", ha=\"center\"),\n        axis_title_x=element_text(size=10, color=INK, margin={\"t\": 8}),\n        axis_title_y=element_text(size=10, color=INK, margin={\"r\": 8}),\n        axis_text_x=element_text(size=8, color=INK_SOFT),\n        axis_text_y=element_text(size=8, color=INK_SOFT),\n        panel_grid_major_x=element_line(color=INK, size=0.2, alpha=0.12),\n        panel_grid_major_y=element_line(color=INK, size=0.2, alpha=0.12),\n        panel_grid_minor=element_blank(),\n        panel_background=element_rect(fill=PAGE_BG, color=\"none\"),\n        plot_background=element_rect(fill=PAGE_BG, color=\"none\"),\n        legend_title=element_text(size=9, weight=\"bold\", color=INK),\n        legend_text=element_text(size=8, color=INK_SOFT),\n        legend_position=\"right\",\n        legend_background=element_rect(fill=ELEVATED_BG, color=INK_SOFT, size=0.3),\n        legend_key=element_rect(fill=\"none\", color=\"none\"),\n        axis_ticks=element_line(color=INK_SOFT, size=0.2),\n        plot_margin=0.04,\n    )\n)\n\nplot.save(f\"plot-{THEME}.png\", dpi=400, width=8, height=4.5, units=\"in\", verbose=False)\n"}