统计学_生存分析/Weibull Distribution韦布尔分布（python代码实现）

大鹏一日同风起，扶摇直上九万里。这篇文章主要讲述统计学_生存分析/Weibull Distribution韦布尔分布（python代码实现）相关的知识，希望能为你提供帮助。
威布尔分布（Weibull distribution），又称韦伯分布或韦布尔分布，是可靠性分析和寿命检验的理论基础。威布尔分布：在可靠性工程中被广泛应用，尤其适用于机电类产品的磨损累计失效的分布形式。由于它可以利用概率值很容易地推断出它的分布参数，被广泛应用于各种寿命试验的数据处理。随机变量分布之一。威布尔分布(Ⅲ型极值分布)记为W(k,a,b)。瑞典工程师威布尔从30年代开始研究轴承寿命，以后又研究结构强度和疲劳等问题。他采用了“链式”模型来解释结构强度和寿命问题。这个模型假设一个结构是由若干小元件(设为n个)串联而成，于是可以形象地将结构看成是由n个环构成的一条链条，其强度(或寿命)取决于最薄弱环的强度(或寿命)。单个链的强度(或寿命)为一随机变量，设各环强度(或寿命)相互独立，分布相同，则求链强度(或寿命)的概率分布就变成求极小值分布问题，由此给出威布尔分布函数。由于零件或结构的疲劳强度(或寿命)也应取决于其最弱环的强度(或寿命)，也应能用威布尔分布描述。根据1943年苏联格涅坚科的研究结果，不管随机变量的原始分布如何，它的极小值的渐近分布只能有三种，而威布尔分布就是第Ⅲ种极小值分布。由于威布尔分布是根据最弱环节模型或串联模型得到的，能充分反映材料缺陷和应力集中源对材料疲劳寿命的影响，而且具有递增的失效率，所以，将它作为材料或零件的寿命分布模型或给定寿命下的疲劳强度模型是合适的。威布尔分布有多种形式，包括一参数威布尔分布、二参数威布尔分布、三参数威布尔分布或混合威布尔分布。三参数的威布尔分布由形状、尺度（范围）和位置三个参数决定。其中形状参数是最重要的参数，决定分布密度曲线的基本形状，尺度参数起放大或缩小曲线的作用，但不影响分布的形状。通过改变形状参数可以表示不同阶段的失效情况；也可以作为许多其他分布的近似，如，可将形状参数设为合适的值以近似正态、对数正态、指数等分布。二参数的威布尔分布主要用于滚动轴承的寿命试验以及高应力水平下的材料疲劳试验，三参数的威布尔分布用于低应力水平的材料及某些零件的寿命试验，一般而言，它具有比对数正态分布更大的适用性。但是，威布尔分布参数的分析法估计较复杂，区间估计值过长，实践中常采用概率纸估计法，从而降低了参数的估计精度．这是威布尔分布目前存在的主要缺点，也限制了它的应用。
历史1. 1927年，Fréchet(1927)首先给出这一分布的定义。
2. 1933年，Rosin和Rammler在研究碎末的分布时，第一次应用了韦伯分布（Rosin, P.; Rammler, E. (1933), "The Laws Governing the Fineness of Powdered Coal", Journal of the Institute of Fuel 7: 29 - 36.）。
3. 1951年，瑞典工程师、数学家Waloddi Weibull（1887-1979）详细解释了这一分布，于是，该分布便以他的名字命名为Weibull Distribution。

文章图片

Weibull Distribution应用
1.生存分析
2.工业制造
研究生产过程和运输时间关系
3.极值理论
4.预测天气
【统计学_生存分析/Weibull Distribution韦布尔分布（python代码实现）】5.可靠性和失效分析
6.雷达系统
对接受到的杂波信号的依分布建模
7.拟合度
无线通信技术中，相对指数衰减频道模型，Weibull衰减模型对衰减频道建模有较好的拟合度
8.量化寿险模型的重复索赔
9.预测技术变革
10.风速
由于曲线形状与现实状况很匹配，被用来描述风速的分布
11.应用于保险业，病人治疗，信用卡诈骗,信用卡拖欠

文章图片

文章图片

测试脚本

文章图片

测试数据

文章图片

T is an array of durations, E is a either boolean or binary array representing whether the “death” was observed (alternatively an individual can be censored).

import lifelines from lifelines.datasets import load_waltons df = load_waltons() # returns a Pandas DataFrame T = df[\'T\'] E = df[\'E\'] from lifelines import KaplanMeierFitter kmf = KaplanMeierFitter() kmf.fit(T, event_observed=E) # more succiently, kmf.fit(T,E) kmf.survival_function_ \'\'\' Out[7]: KM_estimate timeline 0.01.000000 6.00.993865 7.00.987730 9.00.969210 13.00.950690 15.00.938344 17.00.932170 19.00.913650 22.00.888957 26.00.858090 29.00.827224 32.00.821051 33.00.802531 36.00.790184 38.00.777837 41.00.734624 43.00.728451 45.00.672891 47.00.666661 48.00.616817 51.00.598125 \'\'\' kmf.median_ \'\'\' Out[8]: 56.0 \'\'\' kmf.plot()

文章图片

import lifelines from lifelines.datasets import load_waltons from lifelines import KaplanMeierFitter df = load_waltons() # returns a Pandas DataFrame kmf = KaplanMeierFitter() T = df[\'T\'] E = df[\'E\'] groups = df[\'group\'] ix = (groups == \'miR-137\') kmf.fit(T[~ix], E[~ix], label=\'control\') ax = kmf.plot() kmf.fit(T[ix], E[ix], label=\'miR-137\') kmf.plot(ax=ax)

文章图片

import numpy as np import matplotlib.pyplot as plt from lifelines.plotting import plot_lifetimes from numpy.random import uniform, exponential N = 25 current_time = 10 actual_lifetimes = np.array([[exponential(12), exponential(2)][uniform()< 0.5] for i in range(N)]) observed_lifetimes = np.minimum(actual_lifetimes,current_time) observed= actual_lifetimes < current_time plt.xlim(0,25) plt.vlines(10,0,30,lw=2, linestyles="--") plt.xlabel(\'time\') plt.title(\'Births and deaths of our population, at $t=10$\') plot_lifetimes(observed_lifetimes, event_observed=observed) print "Observed lifetimes at time %d:\\n"%(current_time), observed_lifetimes

文章图片

import pandas as pd import lifelines from lifelines import KaplanMeierFitter import matplotlib.pyplot as plt data = https://www.songbingjia.com/android/lifelines.datasets.load_dd() kmf = KaplanMeierFitter() T = data["duration"] C = data["observed"] kmf.fit(T, event_observed=C ) plt.title(\'Survival function of political regimes\') kmf.survival_function_.plot() kmf.plot() kmf.median_

文章图片

import pandas as pd import lifelines from lifelines import KaplanMeierFitter import matplotlib.pyplot as plt data = https://www.songbingjia.com/android/lifelines.datasets.load_dd() kmf = KaplanMeierFitter() T = data["duration"] C = data["observed"] kmf.fit(T, event_observed=C ) plt.title(\'Survival function of political regimes\') kmf.survival_function_.plot() kmf.plot() ax = plt.subplot(111) dem = (data["democracy"] == "Democracy") kmf.fit(T[dem], event_observed=C[dem], label="Democratic Regimes") kmf.plot(ax=ax, ci_force_lines=True) kmf.fit(T[~dem], event_observed=C[~dem], label="Non-democratic Regimes") plt.ylim(0,1); plt.title("Lifespans of different global regimes") kmf.plot(ax=ax, ci_force_lines=True)

文章图片

测试代码

# -*- coding: utf-8 -*- # Import standard packages import numpy as np import matplotlib.pyplot as plt from numpy.random import uniform, exponential import os # additional packages from lifelines.plotting import plot_lifetimes import sys sys.path.append(os.path.join(\'..\', \'..\', \'Utilities\')) try: # Import formatting commands if directory "Utilities" is available from ISP_mystyle import setFontsexcept ImportError: # Ensure correct performance otherwise def setFonts(*options): return# Generate some dummy data np.set_printoptions(precision=2) N = 20 study_duration = 12 # Note: a constant dropout rate is equivalent to an exponential distribution! actual_subscriptiontimes = np.array([[exponential(18), exponential(3)][uniform()< 0.5] for i in range(N)]) observed_subscriptiontimes = np.minimum(actual_subscriptiontimes,study_duration) observed= actual_subscriptiontimes < study_duration # Show the data setFonts(18) plt.xlim(0,24) plt.vlines(12, 0, 30, lw=2, linestyles="--") plt.xlabel(\'time\') plt.title(\'Subscription Times, at $t=12$months\') plot_lifetimes(observed_subscriptiontimes, event_observed=observed) print("Observed subscription time at time %d:\\n"%(study_duration), observed_subscriptiontimes)

文章图片

文章图片

where k > 0 is the shape parameter and > 0 is the scale parameter of the distribution. (It is one of the rare cases where we use a shape parameter different from skewness and kurtosis.) Its complementary cumulative distribution function is a stretched exponential function. If the quantity x is a “time-to-failure,” theWeibull distribution gives a distribution for which the failure rate is proportional to a power of time. The shape parameter, k, is that power plus one, and so this parameter can be interpreted directly as follows: ? Avalueofk < 1 indicates that the failure rate decreases over time. This happens if there is significant “infant mortality,” or defective items failing early and the failure rate decreasing over time as the defective items are weeded out of the population. ? Avalueofk D 1 indicates that the failure rate is constant over time. This might suggest random external events are causing mortality, or failure. ? Avalueofk > 1 indicates that the failure rate increases with time. This happens if there is an “aging” process, or parts that are more likely to fail as time goes on. An example would be products with a built-in weakness that fail soon after the warranty expires. In the field of materials science, the shape parameter k of a distribution of strengths is known as the Weibull modulus.
其中 k > 0 是形状参数，> 0 是分布的尺度参数。（这是我们使用不同于偏度和峰度的形状参数的罕见情况之一。）它的互补累积分布函数是拉伸指数函数。如果数量 x 是“故障时间”，则威布尔分布给出了故障率与时间的幂成正比的分布。形状参数 k 是功率加 1，因此该参数可以直接解释如下： ? Avalueofk < 1 表示故障率随时间降低。如果存在显着的“婴儿死亡率”，或者有缺陷的产品过早失效，并且随着有缺陷的产品从总体中剔除，失效率会随着时间的推移而降低，就会发生这种情况。 ? Avalueofk D 1 表示故障率随时间保持恒定。这可能表明随机的外部事件正在导致死亡或失败。 ? Avalueofk > 1 表示故障率随时间增加。如果存在“老化”过程，或者随着时间的推移更可能出现故障的部件，就会发生这种情况。一个例子是具有内置弱点的产品，在保修期满后很快就会失效。在材料科学领域，强度分布的形状参数 k 被称为威布尔模量。

生存分析补充说明