韦布尔分布的期望_weibull分布参数意义「建议收藏」

前言：

本章我们将利用算法工具Matlab，对常用的几个韦伯分布进行计算分析。
本章实例包括显示面板的数据分析、轮胎数据分析和纺线的数据分析。

文后附录，转分享韦伯分布的Python 和 C++的代码。

例一：一个纺线张力的韦伯分布模型和MT应用

1.1 图形的生成

这个例子，分析纺线的张力测量。

已知经验数据：参数shap β =2，参数scale η=0.5的案例
【Matlab的源码如下：】

beta =2;eta =0.5;
time = 0:0.1:5;
fx =(beta/eta)*(time/eta).^(beta-1).*exp(-(time/eta).^beta);
ft=1-exp(-(time./eta).^beta);
subplot(211);plot(time,fx,'linewidth',2.5);xlabel('Time');ylabel('PDF');
subplot(212);plot(time,ft,'linewidth',2.5);xlabel('Time');ylabel('CDF');

【如果改变参数：下一段代码，展示了连续变换参数的结果】

clc;clear;close all;
beta = [0.1 0.5 1 2]; eta = [1 2];time = 0:0.1:10;
for i = 1:length(beta)
    for j = 1:length(eta)
        fx =(beta(i)/eta(j))*(time/eta(j)).^(beta(i)-1).*exp(-(time/eta(j)).^beta(i));
        ft=1-exp(-(time./eta(j)).^beta(i));
        subplot(211);plot(time,fx,'linewidth',2.5);xlabel('Time');ylabel('PDF');
        subplot(212);plot(time,ft,'linewidth',2.5);xlabel('Time');ylabel('CDF');
    end
end

【案，大致解释一下，首先，前面给出的是韦伯分布的双参数值，分别是shap 和 scale，然后，时间段的步长为j0.1,那么就是51步，然后分别实现了PDF，CDF的公式，然后把他们画出来。】

1.2 weibull Matlab函数的使用

我们可以利用Matlab的函数，生成韦伯分布的随机数字。

1.2.1 利用wblrad函数生成韦伯分布数据

我们依旧用上面例子里面的参数，β =2，η=0.5的案例

a = wblrnd(0.5,2,[10000 1]);
time = linspace(0,1,10000);
plot(time,a,'.');
xlabel('time');
ylabel('weibull drandom number');

【如果我们变更一下】

a = wblrnd(2,2,[10000 1]);
time = linspace(0,1,10000);
plot(time,a,'.');
xlabel('time');
ylabel('weibull drandom number');

1.2.2 估算韦伯分布的参数

a = wblrnd(0.5,2,[10000 1]);
time = linspace(0,1,10000);
plot(time,a,'.');
xlabel('time');
ylabel('weibull drandom number');
parmhat = wblfit(a);

parmhat =
0.4989 2.0330
【案，验证，结果和我们之前设定结果一致】

1.2.3 绘制韦伯分布的图像

继续输入

wblplot(a)

图形接近线性，说明图像是韦伯分布。

1.2.4 获取韦伯分布的数学期望和均方差

[M V]= wblstat(0.5,2);

M =
0.4431
V =
0.0537

1.2.5 获取PDF

继续输入

fx = wblpdf(time,0.5,2);
plot(time,fx,'linewidth',2.5);xlabel('Time');ylabel('PDF');

1.2.6 获取CDF

继续输入

ft = wblcdf(time,0.5,2);
plot(time,ft,'linewidth',2.5);
xlabel('Time');
ylabel('CDF');

1.2.7 生成Histogram 数据：

histfit(a,100,'weibull');

1.3 weibull Matlab对象的使用

1.3.1 创建weibull分布对象

注册一个weibulll分布的对象，

pd = makedist('weibull','beta',2,'eta',0.5);
histfit(r,1000,'weibull');

然后，用这个对象可以做很多事情

pdf = pdf(pd,time);
plot(time,pdf,'linewidth',2.5);

例二：显示面板分析和源码

本节我们尝试用Matlab的统计工具和机器学习工具，对平面面板显示器的生产数据进行分析。

问题提出：已知平板显示器，经验数据为符合韦伯分布，且，缩放参数scale η=5000 ，形状参数shap β =0.5，求
1 面板寿命>20000hrs的概率
2 面板寿命<10000hrs的概率
3 10000<面板寿命<20000hrs的概率
4 面板报废的数学期望

问题转为求：
1 P(x > 20000)
2 P(x < 10000)
3 P(10000<x<20000)
4 E(x)

由第一章，我们知道故障的累计分布函数CDF公式为：
$\large\displaystyle F(t) = 1 – e^{-(\frac{t}{\eta })^{\beta }} (t>0)$
我们在Matlab中这么做，

crc;clear;
eta = 5000; beta = 0.5;
p2w_below = wblcdf(20000,eta ,beta);
p1w_below = wblcdf(10000,eta ,beta);
p2w_beyond = 1 - p2w_below ;
p1wto2w = wblcdf(20000,eta ,beta) - wblcdf(10000,eta ,beta);
[M V] = wblstat(eta ,beta);
time = 1:100:25000;
subplot(211); fx = wblpdf(time,eta,beta);
plot(time,fx,'linewidth',2.5);xlabel('Time Hour');ylabel('PDF');
subplot(212); ft = wblcdf(time,eta,beta);
plot(time,ft,'linewidth',2.5);xlabel('Time Hour');ylabel('CDF');

计算结果如下：

1 P(x > 20000) = 0.135335283236613 【p2w_beyond 】
2 P(x < 10000) = 0.756883265565786 【p1w_below 】
3 P(10000<x<20000) = 0.107781451197602 【p1wto2w 】
4 E(x）= 10000 【M】

例二：轮胎失效分析和源码

20000个轮胎，满足韦伯分布，历史数据有：缩放参数scale η=5000 ，形状参数shap β =3.5
求：
6000公里的时候，有多少轮胎报废？

crc;clear;
eta = 5000; beta = 3.5;
TotalTires = 20000;
cdf = wblcdf(6000,eta ,beta);


FailTireNumbert = TotalTires * cdf;

time = 1:10:10000;
subplot(211); fx = wblpdf(time,eta,beta);
plot(time,fx,'linewidth',2.5);xlabel('Kilometer');ylabel('PDF');
subplot(212); ft = wblcdf(time,eta,beta);
plot(time,ft,'linewidth',2.5);xlabel('Kilometer');ylabel('CDF');

FailTireNumbert = 1.698740114017983e+04
cdf = 0.849370057008992 【84.9%的轮胎可以保证到6000公里没问题】

用Matlab工具箱分析韦伯分布

Statistics and Machine Learning Toolbox
Probability Distribution Plotter
【直接上图吧，不解释了】

——————————————————————————————————————————————————————————

Matlab weibull 函数简介

wblrnd

生成韦伯分布的随机数值。Weibull random numbers

R = wblrnd(A,B)

参数 A, scale parameter
参数 B，shape parameter

[parmhat,parmci] = wblfit(data)

韦伯分布参数获取
parmhat = 缩放参数，
parmci = 形状参数

wblplot(x)

通过数据X,绘制韦伯概率分布图形.
图形如果是线性的话，说明是韦伯分布，否则不是。

[M,V] = wblstat(A,B)

输入：A, scale缩放参数； B，shape形状参数
返回：M，均值数学期望；V，均方差

Y = wblpdf(X,A,B)

输入：A, scale缩放参数； B，shape形状参数

p = wblcdf(x,a,b)

输入：X,样本周期空间，A, scale缩放参数； B，shape形状参数

WeibullDistribution

韦伯分布对象生成，这里他输入变量必须为，a，b

Matlab weibull 分布的函数列表：

wblcdf Weibull cumulative distribution function
wblpdf Weibull probability density function
wblinv Weibull inverse cumulative distribution function
wbllike Weibull negative log-likelihood
wblstat Weibull mean and variance
wblfit Weibull parameter estimates
wblrnd Weibull random numbers
wblplot Weibull probability plot

——————————————————————————————————————————————————————————

其他代码参考：

【这里COPY了其他开发者，用C++ 和 Python实现的源码，希望对分析有帮助】
韦伯分布python实现

import numpy as np
import matplotlib.pyplot as plt


# define the pdf of weibull distribution
def weib(x, scale, shape):
    return (shape / scale) * (x / scale)**(shape - 1) * np.exp(-(x / scale) ** shape)


scale = 50
shape = 1.5
x = np.arange(1, scale*2)
y = np.zeros(len(x))  # [0 for i in range(len(x))]
for i in range(len(x)):
    y[i] = weib(x[i], scale, shape)

scale = 50
shape = 2.5
y1 = np.zeros(len(x))  # [0 for i in range(len(x))]
for i in range(len(x)):
    y1[i] = weib(x[i], scale, shape)
scale = 50
shape = 4
y2 = np.zeros(len(x))  # [0 for i in range(len(x))]
for i in range(len(x)):
    y2[i] = weib(x[i], scale, shape)


scale = 30
shape = 2.5
y3 = np.zeros(len(x))  # [0 for i in range(len(x))]
for i in range(len(x)):
    y3[i] = weib(x[i], scale, shape)
scale = 70
shape = 2.5
y4 = np.zeros(len(x))  # [0 for i in range(len(x))]
for i in range(len(x)):
    y4[i] = weib(x[i], scale, shape)


plt.subplot(2, 1, 1)
plt.plot(x, y, 'r', label='scale=50, shape=1.5')
plt.plot(x, y1, 'b', label='scale=50, shape=2.5')
plt.plot(x, y2, 'g', label='scale=50, shape=4')
plt.legend()
plt.subplot(2, 1, 2)
plt.plot(x, y3, 'r', label='scale=30, shape=2.5')
plt.plot(x, y1, 'b', label='scale=50, shape=2.5')
plt.plot(x, y4, 'g', label='scale=70, shape=2.5')
plt.legend()
plt.show()
————————————————
版权声明：本文为CSDN博主「心态与做事习惯决定人生高度」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。
原文链接：https://blog.csdn.net/robert_chen1988/article/details/100785494

韦伯分布C++实现：

public class WeibullDistribution : ContinuousDistribution
    { 
   
        private class MMSystem
        { 
   
            private double mean;

            private double variance;

            public MMSystem(NumericalVariable variable)
            { 
   
                this.mean = variable.Mean;
                this.variance = variable.Variance;
            }

            private Vector Evaluate(Vector input, Vector output)
            { 
   
                double num;
                num = input[0];
                double num2;
                num2 = input[1];
                double num3;
                num3 = GammaFunctions.Gamma(1.0 + 1.0 / num2);
                output[0] = this.mean - num * num3;
                output[1] = this.variance - num * num * (GammaFunctions.Gamma(1.0 + 2.0 / num2) - num3 * num3);
                return output;
            }

            public void Solve(ref double a, ref double b)
            { 
   
                DoglegSystemSolver doglegSystemSolver;
                doglegSystemSolver = new DoglegSystemSolver();
                doglegSystemSolver.TargetFunction = Evaluate;
                doglegSystemSolver.InitialGuess = new GeneralVector(a, b);
                Vector vector;
                vector = doglegSystemSolver.Solve();
                a = vector[0];
                b = vector[1];
            }
        }

        private class MLSystem
        { 
   
            private NumericalVariable variable;

            private NumericalVariable log;

            public MLSystem(NumericalVariable variable)
            { 
   
                this.variable = variable;
                this.log = variable.Transforms.Log();
            }

            private Vector Evaluate(Vector input, Vector output)
            { 
   
                double num;
                num = input[0];
                double num2;
                num2 = input[1];
                int length;
                length = this.variable.Length;
                double num3;
                num3 = 0.0;
                double num4;
                num4 = 0.0;
                double num5;
                num5 = 0.0;
                double num6;
                num6 = Math.Log(num);
                for (int i = 0; i < length; i++)
                { 
   
                    double num7;
                    num7 = Math.Pow(this.variable[i] / num, num2);
                    double num8;
                    num8 = this.log[i] - num6;
                    num3 += num8;
                    num5 += num8 * num7;
                    num4 += num7;
                }
                output[0] = (double)length + num2 * num3 - num2 * num5;
                output[1] = (double)length - num4;
                return output;
            }

            public void Solve(ref double a, ref double b)
            { 
   
                DoglegSystemSolver doglegSystemSolver;
                doglegSystemSolver = new DoglegSystemSolver();
                doglegSystemSolver.TargetFunction = Evaluate;
                doglegSystemSolver.InitialGuess = new GeneralVector(a, b);
                doglegSystemSolver.Solve();
                a = doglegSystemSolver.CurrentPoint[0];
                b = doglegSystemSolver.CurrentPoint[1];
            }
        }

        private double scale;

        private double shape;

        private double location;

        public override double Mean => this.location + this.scale * GammaFunctions.Gamma(1.0 + 1.0 / this.shape);

        public override double Variance => this.scale * this.scale * (GammaFunctions.Gamma(1.0 + 2.0 / this.shape) - Math.Pow(GammaFunctions.Gamma(1.0 + 1.0 / this.shape), 2.0));

        public override double Skewness
        { 
   
            get
            { 
   
                double num;
                num = GammaFunctions.Gamma(1.0 + 1.0 / this.shape);
                double num2;
                num2 = GammaFunctions.Gamma(1.0 + 2.0 / this.shape);
                double num3;
                num3 = num2 - num * num;
                return (num * (2.0 * num * num - 3.0 * num2) + GammaFunctions.Gamma(1.0 + 3.0 / this.shape)) / (num3 * Math.Sqrt(num3));
            }
        }

        public override double Kurtosis
        { 
   
            get
            { 
   
                double num;
                num = GammaFunctions.Gamma(1.0 + 1.0 / this.shape);
                double num2;
                num2 = GammaFunctions.Gamma(1.0 + 2.0 / this.shape);
                double num3;
                num3 = num2 - num * num;
                return (6.0 * num * num * (2.0 * num2 - num * num) - 3.0 * num2 * num2 - 4.0 * num * GammaFunctions.Gamma(1.0 + 3.0 / this.shape) + GammaFunctions.Gamma(1.0 + 4.0 / this.shape)) / (num3 * num3);
            }
        }

        public double LocationParameter => this.location;

        public double ScaleParameter => this.scale;

        public double ShapeParameter => this.shape;

        public WeibullDistribution(double shape, double scale, double location)
        { 
   
            if (shape <= 0.0)
            { 
   
                ThrowException.ArgumentOutOfRange("shape");
            }
            if (scale <= 0.0)
            { 
   
                ThrowException.ArgumentOutOfRange("scale");
            }
            this.scale = scale;
            this.shape = shape;
            this.location = location;
        }

        public WeibullDistribution(double shape, double scale)
            : this(shape, scale, 0.0)
        { 
   
        }

        public WeibullDistribution(double shape)
            : this(shape, 1.0, 0.0)
        { 
   
        }

        public WeibullDistribution(NumericalVariable variable, EstimationMethod method)
        { 
   
            if (variable == null)
            { 
   
                ThrowException.ArgumentNull("variable");
            }
            double a;
            a = variable.Mean;
            double b;
            b = 1.0;
            switch (method)
            { 
   
            case EstimationMethod.Default:
            case EstimationMethod.MatchingMoments:
                new MMSystem(variable).Solve(ref a, ref b);
                break;
            case EstimationMethod.MaximumLikelihood:
                new MLSystem(variable).Solve(ref a, ref b);
                break;
            default:
                ThrowException.ArgumentOutOfRange("method");
                break;
            }
            this.scale = a;
            this.shape = b;
        }

        public WeibullDistribution(NumericalVariable variable)
            : this(variable, EstimationMethod.MatchingMoments)
        { 
   
        }

        public override double ProbabilityDensityFunction(double x)
        { 
   
            if (x < this.location)
            { 
   
                return 0.0;
            }
            if (x == this.location)
            { 
   
                if (this.shape == 1.0)
                { 
   
                    return 1.0 / this.scale;
                }
                if (this.shape < 1.0)
                { 
   
                    return double.PositiveInfinity;
                }
                return 0.0;
            }
            if (this.shape == 1.0)
            { 
   
                return Math.Exp((0.0 - (x - this.location)) / this.scale) / this.scale;
            }
            return this.shape / this.scale * Math.Exp(0.0 - Math.Pow((x - this.location) / this.scale, this.shape) + (this.shape - 1.0) * Math.Log((x - this.location) / this.scale));
        }

        public override double DistributionFunction(double x)
        { 
   
            if (x <= this.location)
            { 
   
                return 0.0;
            }
            return 1.0 - Math.Exp(0.0 - Math.Pow((x - this.location) / this.scale, this.shape));
        }

        public override double MomentFunction(int order, double x)
        { 
   
            if (order < 0)
            { 
   
                ThrowException.ArgumentOutOfRange("order");
            }
            if (order == 0)
            { 
   
                return this.DistributionFunction(x) - this.DistributionFunction(0.0);
            }
            if (x <= this.location)
            { 
   
                return 0.0;
            }
            if (this.shape == 1.0)
            { 
   
                return Math.Exp((0.0 - this.location) / this.scale) * ElementaryFunctions.Pow(0.0 - this.scale, order) * GammaFunctions.IncompleteGamma(1 + order, (0.0 - x) / this.scale);
            }
            if (this.location == 0.0)
            { 
   
                return ElementaryFunctions.Pow(this.scale, order) * (GammaFunctions.Gamma(1.0 + (double)order / this.shape) - GammaFunctions.IncompleteGamma(1.0 + (double)order / this.shape, Math.Pow(x / this.scale, this.shape)));
            }
            return base.MomentFunction(order, x);
        }

        public override double InverseDistributionFunction(double probability)
        { 
   
            if (probability < 0.0 || probability > 1.0)
            { 
   
                ThrowException.ArgumentOutOfRange("probability");
            }
            if (probability == 0.0)
            { 
   
                return this.location;
            }
            if (probability == 1.0)
            { 
   
                return double.PositiveInfinity;
            }
            return this.location + this.scale * Math.Pow(0.0 - Math.Log(1.0 - probability), 1.0 / this.shape);
        }

        public override double GetRandomVariate(System.Random random)
        { 
   
            if (random == null)
            { 
   
                ThrowException.ArgumentNull("random");
            }
            return this.location + this.scale * Math.Pow(0.0 - Math.Log(random.NextDouble()), 1.0 / this.shape);
        }
    }

上一章参考：

第1章
第2章
第3章

今天的文章韦布尔分布的期望_weibull分布参数意义「建议收藏」分享到此就结束了，感谢您的阅读。

版权声明：本文内容由互联网用户自发贡献，该文观点仅代表作者本人。本站仅提供信息存储空间服务，不拥有所有权，不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容，请发送邮件至举报，一经查实，本站将立刻删除。
如需转载请保留出处：https://bianchenghao.cn/89337.html

韦布尔分布的期望_weibull分布参数意义「建议收藏」

前言：

例一：一个纺线张力的韦伯分布模型和MT应用

1.1 图形的生成

1.2 weibull Matlab函数的使用

1.2.1 利用wblrad函数生成韦伯分布数据

1.2.2 估算韦伯分布的参数

1.2.3 绘制韦伯分布的图像

1.2.4 获取韦伯分布的数学期望和均方差

1.2.5 获取PDF

1.2.6 获取CDF

1.2.7 生成Histogram 数据：

1.3 weibull Matlab对象的使用

1.3.1 创建weibull分布对象

例二：显示面板分析和源码

例二：轮胎失效分析和源码

用Matlab工具箱分析韦伯分布

Matlab weibull 函数简介

wblrnd

R = wblrnd(A,B)

[parmhat,parmci] = wblfit(data)

wblplot(x)

[M,V] = wblstat(A,B)

Y = wblpdf(X,A,B)

p = wblcdf(x,a,b)

WeibullDistribution

Matlab weibull 分布的函数列表：

其他代码参考：

上一章参考：

相关推荐

发表回复