{"spec_id":"probability-weibull","library":"bokeh","language":"python","code":"\"\"\" anyplot.ai\nprobability-weibull: Weibull Probability Plot for Reliability Analysis\nLibrary: bokeh 3.9.1 | Python 3.13.13\nQuality: 93/100 | Updated: 2026-06-07\n\"\"\"\n\nimport os\nimport sys\nimport time\nfrom pathlib import Path\n\n\n# bokeh.py is the script name — remove its directory from sys.path so that\n# `import bokeh` resolves to the installed package, not this file itself.\n_here = os.path.dirname(os.path.abspath(__file__))\nsys.path[:] = [p for p in sys.path if os.path.abspath(p) != _here]\n\nimport numpy as np\nfrom bokeh.io import output_file, save\nfrom bokeh.models import Band, ColumnDataSource, HoverTool, Label, Span\nfrom bokeh.plotting import figure\nfrom scipy import stats\nfrom selenium import webdriver\nfrom selenium.webdriver.chrome.options import Options\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 palette\nBRAND = \"#009E73\"  # position 1 — failures (primary data)\nLAVENDER = \"#C475FD\"  # position 2 — censored observations\nANYPLOT_AMBER = \"#DDCC77\"  # reference / caution line\n\n# Data - Turbine blade fatigue-life data (hours to failure)\nnp.random.seed(42)\nn_failures = 25\nn_censored = 5\nn_total = n_failures + n_censored\n\n# Generate Weibull-distributed failure times (shape=2.5, scale=5000)\ntrue_beta = 2.5\ntrue_eta = 5000\nfailure_times = np.sort(stats.weibull_min.rvs(true_beta, scale=true_eta, size=n_failures))\n\n# Censored observations (suspended tests, typically at higher times)\ncensored_times = np.sort(np.random.uniform(4000, 7000, size=n_censored))\n\n# Combine and compute median rank plotting positions for failures only\n# Using Bernard's approximation: (i - 0.3) / (n + 0.4)\nranks = np.arange(1, n_failures + 1)\nmedian_rank = (ranks - 0.3) / (n_total + 0.4)\n\n# Weibull linearized y-axis: ln(-ln(1 - F))\nweibull_y = np.log(-np.log(1 - median_rank))\nlog_failure_times = np.log(failure_times)\n\n# Fit line in linearized space: weibull_y = beta * ln(t) - beta * ln(eta)\nslope_fit, intercept_fit, r_value, _, _ = stats.linregress(log_failure_times, weibull_y)\nbeta_fit = slope_fit\neta_fit = np.exp(-intercept_fit / beta_fit)\nr_squared = r_value**2\n\n# Fitted line data\nlog_t_line = np.linspace(np.log(failure_times.min()) - 0.5, np.log(failure_times.max()) + 0.5, 200)\nweibull_y_line = beta_fit * log_t_line - beta_fit * np.log(eta_fit)\n\n# Reference probability levels for y-axis labels\nprob_levels = np.array([0.01, 0.05, 0.10, 0.20, 0.50, 0.632, 0.90, 0.99])\nweibull_y_levels = np.log(-np.log(1 - prob_levels))\n\n# B10 life calculation (10% failure probability)\nb10_weibull_y = np.log(-np.log(1 - 0.10))\nb10_life = eta_fit * (-np.log(1 - 0.10)) ** (1 / beta_fit)\n\n# 63.2% characteristic life reference\nchar_life_y = np.log(-np.log(1 - 0.632))\n\n# Compute censored plotting positions using adjusted ranks\nall_times = np.concatenate([failure_times, censored_times])\nall_censored = np.concatenate([np.zeros(n_failures), np.ones(n_censored)])\nsort_idx = np.argsort(all_times)\nall_times_sorted = all_times[sort_idx]\nall_censored_sorted = all_censored[sort_idx]\n\nreverse_ranks = np.arange(n_total, 0, -1)\nadjusted_rank = np.zeros(n_total)\nprev_adj = 0\nfor i in range(n_total):\n    if all_censored_sorted[i] == 0:\n        prev_adj += 1\n        adjusted_rank[i] = prev_adj\n    else:\n        increment = (n_total + 1 - prev_adj) / (reverse_ranks[i] + 1)\n        prev_adj += increment\n        adjusted_rank[i] = prev_adj\n\ncensored_mask = all_censored_sorted == 1\ncensored_adjusted_ranks = adjusted_rank[censored_mask]\ncensored_median_rank = (censored_adjusted_ranks - 0.3) / (n_total + 0.4)\ncensored_median_rank = np.clip(censored_median_rank, 0.001, 0.999)\ncensored_weibull_y = np.log(-np.log(1 - censored_median_rank))\n\n# Data sources\nfailure_source = ColumnDataSource(\n    data={\n        \"time\": failure_times,\n        \"log_time\": log_failure_times,\n        \"weibull_y\": weibull_y,\n        \"prob_pct\": [f\"{p * 100:.1f}%\" for p in median_rank],\n        \"time_fmt\": [f\"{t:.0f}\" for t in failure_times],\n    }\n)\n\ncensored_source = ColumnDataSource(\n    data={\n        \"time\": censored_times,\n        \"log_time\": np.log(censored_times),\n        \"weibull_y\": censored_weibull_y,\n        \"time_fmt\": [f\"{t:.0f}\" for t in censored_times],\n    }\n)\n\n# Confidence band (approximate ±2 SE)\nse_y = np.sqrt(\n    (1 - r_squared)\n    * np.var(weibull_y)\n    * (\n        1 / len(weibull_y)\n        + (log_t_line - np.mean(log_failure_times)) ** 2 / np.sum((log_failure_times - np.mean(log_failure_times)) ** 2)\n    )\n)\nband_source = ColumnDataSource(\n    data={\n        \"log_time\": log_t_line,\n        \"weibull_y\": weibull_y_line,\n        \"upper\": weibull_y_line + 2 * se_y,\n        \"lower\": weibull_y_line - 2 * se_y,\n    }\n)\nline_source = ColumnDataSource(data={\"log_time\": log_t_line, \"weibull_y\": weibull_y_line})\n\n# Title — descriptive prefix + spec-id · language · library · anyplot.ai\n# \"Turbine Blade Fatigue Life · probability-weibull · python · bokeh · anyplot.ai\" = 76 chars\n# fontsize = round(50 * 67 / 76) = 44pt\nTITLE = \"Turbine Blade Fatigue Life · probability-weibull · python · bokeh · anyplot.ai\"\n\n# Canvas: 3200×1800 (landscape 16:9) per bokeh.md hard rule\nW, H = 3200, 1800\np = figure(\n    width=W,\n    height=H,\n    title=TITLE,\n    x_axis_label=\"Time to Failure (hours)\",\n    y_axis_label=\"Cumulative Failure Probability\",\n    x_range=(log_t_line.min() - 0.2, log_t_line.max() + 0.2),\n    y_range=(weibull_y_levels[0] - 0.3, weibull_y_levels[-1] + 0.3),\n    toolbar_location=None,\n    min_border_bottom=160,\n    min_border_left=200,\n    min_border_top=110,\n    min_border_right=60,\n)\n\n# Theme-adaptive backgrounds\np.background_fill_color = PAGE_BG\np.border_fill_color = PAGE_BG\np.outline_line_color = INK_SOFT\n\n# 63.2% characteristic life reference line\nchar_life_span = Span(\n    location=char_life_y, dimension=\"width\", line_color=ANYPLOT_AMBER, line_width=4, line_dash=\"dashed\", line_alpha=0.8\n)\np.add_layout(char_life_span)\n\n# 95% confidence band\nband = Band(\n    base=\"log_time\",\n    lower=\"lower\",\n    upper=\"upper\",\n    source=band_source,\n    fill_color=BRAND,\n    fill_alpha=0.1,\n    line_color=BRAND,\n    line_alpha=0.2,\n    line_width=1,\n)\np.add_layout(band)\n\n# Fitted line\np.line(\"log_time\", \"weibull_y\", source=line_source, line_color=BRAND, line_width=4, legend_label=\"Weibull Fit\")\n\n# Failure data points (filled circles)\nfailure_renderer = p.scatter(\n    \"log_time\",\n    \"weibull_y\",\n    source=failure_source,\n    size=18,\n    color=BRAND,\n    alpha=0.9,\n    line_color=PAGE_BG,\n    line_width=2,\n    legend_label=\"Failures\",\n)\n\n# Censored data points (hollow triangles — distinct shape from failures)\ncensored_renderer = p.scatter(\n    \"log_time\",\n    \"weibull_y\",\n    source=censored_source,\n    size=20,\n    marker=\"triangle\",\n    color=PAGE_BG,\n    alpha=0.95,\n    line_color=LAVENDER,\n    line_width=3,\n    legend_label=\"Censored\",\n)\n\n# Hover tools\nhover_failure = HoverTool(\n    renderers=[failure_renderer],\n    tooltips=[(\"Time\", \"@time_fmt hours\"), (\"Cum. Probability\", \"@prob_pct\")],\n    mode=\"mouse\",\n)\np.add_tools(hover_failure)\n\nhover_censored = HoverTool(renderers=[censored_renderer], tooltips=[(\"Censored at\", \"@time_fmt hours\")], mode=\"mouse\")\np.add_tools(hover_censored)\n\n# Parameter annotation box — placed below the top tick so text flows downward fully visible\nparam_text = f\"β = {beta_fit:.2f}  (shape)\\nη = {eta_fit:.0f} h  (scale)\\nR² = {r_squared:.4f}\"\nparam_label = Label(\n    x=log_t_line.min() + 0.1,\n    y=weibull_y_levels[-1] - 0.6,\n    text=param_text,\n    text_font_size=\"30pt\",\n    text_color=INK,\n    text_font_style=\"bold\",\n    background_fill_color=ELEVATED_BG,\n    background_fill_alpha=0.95,\n    border_line_color=INK_SOFT,\n    border_line_alpha=0.4,\n    padding=14,\n)\np.add_layout(param_label)\n\n# 63.2% characteristic life label\nchar_label = Label(\n    x=log_t_line.max() - 0.5,\n    y=char_life_y + 0.15,\n    text=\"63.2% (Characteristic Life)\",\n    text_font_size=\"26pt\",\n    text_color=ANYPLOT_AMBER,\n    text_font_style=\"italic\",\n    text_align=\"right\",\n    background_fill_color=ELEVATED_BG,\n    background_fill_alpha=0.9,\n)\np.add_layout(char_label)\n\n# B10 life annotation (10% failure probability)\nb10_label = Label(\n    x=np.log(b10_life) + 0.1,\n    y=b10_weibull_y - 0.28,\n    text=f\"B10 = {b10_life:.0f} h\",\n    text_font_size=\"27pt\",\n    text_color=INK_SOFT,\n    text_font_style=\"italic\",\n)\np.add_layout(b10_label)\n\n# B10 reference marker on the fit line\np.scatter(\n    [np.log(b10_life)],\n    [b10_weibull_y],\n    size=18,\n    marker=\"diamond\",\n    color=INK_SOFT,\n    alpha=0.8,\n    line_color=PAGE_BG,\n    line_width=2,\n)\n\n# Custom y-axis tick labels showing probability percentages\np.yaxis.ticker = list(weibull_y_levels)\nprob_labels = {float(y): f\"{p * 100:.1f}%\" for y, p in zip(weibull_y_levels, prob_levels, strict=True)}\np.yaxis.major_label_overrides = prob_labels\n\n# Custom x-axis: show actual hours\nlog_tick_values = np.log([1000, 2000, 3000, 5000, 7000, 10000])\np.xaxis.ticker = list(log_tick_values)\ntime_labels = {float(lt): f\"{np.exp(lt):.0f}\" for lt in log_tick_values}\np.xaxis.major_label_overrides = time_labels\n\n# Title\np.title.text_font_size = \"44pt\"\np.title.text_color = INK\np.title.align = \"center\"\n\n# Axis labels and ticks\np.xaxis.axis_label_text_font_size = \"42pt\"\np.yaxis.axis_label_text_font_size = \"42pt\"\np.xaxis.axis_label_text_color = INK\np.yaxis.axis_label_text_color = INK\np.xaxis.axis_label_text_font_style = \"normal\"\np.yaxis.axis_label_text_font_style = \"normal\"\np.xaxis.major_label_text_font_size = \"34pt\"\np.yaxis.major_label_text_font_size = \"34pt\"\np.xaxis.major_label_text_color = INK_SOFT\np.yaxis.major_label_text_color = INK_SOFT\n\n# Legend\np.legend.label_text_font_size = \"34pt\"\np.legend.label_text_color = INK_SOFT\np.legend.location = \"bottom_right\"\np.legend.background_fill_color = ELEVATED_BG\np.legend.background_fill_alpha = 0.92\np.legend.border_line_color = INK_SOFT\np.legend.border_line_alpha = 0.4\np.legend.glyph_height = 40\np.legend.glyph_width = 40\np.legend.padding = 25\np.legend.spacing = 15\np.legend.margin = 20\n\n# Grid — minimal y-grid only for Weibull paper look\np.xgrid.grid_line_color = None\np.ygrid.grid_line_color = INK\np.ygrid.grid_line_alpha = 0.15\n\n# Axes chrome\np.xaxis.axis_line_color = INK_SOFT\np.yaxis.axis_line_color = INK_SOFT\np.xaxis.axis_line_width = 2\np.yaxis.axis_line_width = 2\np.xaxis.minor_tick_line_color = None\np.yaxis.minor_tick_line_color = None\np.xaxis.major_tick_line_color = INK_SOFT\np.yaxis.major_tick_line_color = INK_SOFT\np.xaxis.major_tick_out = 8\np.yaxis.major_tick_out = 8\n\np.outline_line_color = None\n\n# Save HTML (interactive artifact)\noutput_file(f\"plot-{THEME}.html\")\nsave(p)\n\n# Screenshot with headless Chrome (Selenium 4 / Selenium Manager)\n# Use CDP setDeviceMetricsOverride so the inner viewport is authoritative:\n# --window-size alone is eaten by Chrome in headless mode (gives 1661 instead of 1800).\nopts = Options()\nfor arg in (\n    \"--headless=new\",\n    \"--no-sandbox\",\n    \"--disable-dev-shm-usage\",\n    \"--disable-gpu\",\n    f\"--window-size={W},{H}\",\n    \"--hide-scrollbars\",\n):\n    opts.add_argument(arg)\ndriver = webdriver.Chrome(options=opts)\ndriver.execute_cdp_cmd(\n    \"Emulation.setDeviceMetricsOverride\", {\"width\": W, \"height\": H, \"deviceScaleFactor\": 1, \"mobile\": False}\n)\ndriver.get(f\"file://{Path(f'plot-{THEME}.html').resolve()}\")\ntime.sleep(3)\ndriver.save_screenshot(f\"plot-{THEME}.png\")\ndriver.quit()\n\n# Pin to exact canvas dims so the post-render gate passes\nfrom PIL import Image as _PILImage\n\n\n_img = _PILImage.open(f\"plot-{THEME}.png\").convert(\"RGB\")\nif _img.size != (W, H):\n    _norm = _PILImage.new(\"RGB\", (W, H), PAGE_BG)\n    _norm.paste(_img, ((W - _img.size[0]) // 2, (H - _img.size[1]) // 2))\n    _norm.save(f\"plot-{THEME}.png\")\n"}