MATLAB在可靠性数据分析中的应用
2019/07/30 00:00:21
1前言
可靠性(Reliability)数据是系统可靠性(Reliability)设计、研究、分析、评定和改进的基本依据，对它的分析则是一切可靠性(Reliability)工作的基础。通过有计划、有目的地收集和分析系统各组成单元的可靠性(Reliability)数据，可定量评定系统的可靠性(Reliability)水平，发现其可靠性(Reliability)的薄弱环节和产生原因，经过改进，使系统的性能与可靠性(Reliability)水平不断提高。所以，可靠性(Reliability)数据的分析在系统可靠性(Reliability)研究中具有十分重要的地位。
收集到的系统可靠性(Reliability)原始数据往往是杂乱无章的，因此不能直接用其对系统的可靠性(Reliability)进行分析，作出判断，必须运用概率论与数理统计知识对其进行一系列整理、分析后，才能透过大量现象掌握它的本质，使其数据显示出规律性。可靠性(Reliability)数据整理分析的步骤如图1所示。
对系统可靠性(Reliability)数据整理分析，需要处理大量的数据，用人工的方法完成这些工作需要花费很长的时间，而且容易出错，而MATLAB的统计工具箱是建立在MATLAB数值计算环境上的工具集，它的统计功能已经可以与SPSS，SAS等软件的统计功能相媲美，它可实现数据资料处理的自动化，保证数据分析工作及时、准确的完成。
2 运用MATLAB绘制可靠性(Reliability)数据分布频数直方图
把系统可靠性(Reliability)数据读入MATLAB的工作空间，便可绘制可靠性(Reliability)数据分布频数直方图。MATLAB统计工具箱提供的hist函数，是用于作直方图的Matlab命令 ，其命令格式如下：
HIST(data，k)；
其中：
Data——原始数据；
K——所分小区间数。
这个命令将区间[min(data),max(data)]分为k等份，并描出频数直方图，k的缺省值为10。如果需要事先给出小区间，则将区间的中点存放在向量nb中，然后使用命令hist(data，nb)。
有时并不需要直方图，只希望知道小区间的位置及数据落在各小区间的频数，可以用命令：
n返回k个小区间的频数，x返回小区间的中点。
3 参数估计
根据绘出的直方图的形状，假设可靠性(Reliability)数据服从某一分布，一般情况可靠性(Reliability)数据服从负指数分布。
指数分布可以通过命令expfit进行参数估计，该命令用极大似然法给出了常用的概率分布的参数的点估计和区域估计值，命令格式如下：
[muhat，muci]＝expfit（data）
式中：muhat——参数μ的估计值；
muci ——参数μ估计的置信区间。
正态分布可以通过命令normfit进行参数估计，该命令用极大似然法给出了常用的概率分布的参数的点估计和区域估计值，命令格式如下：
[muhat，sigmahat，muci，sigmaci]＝normfit（data）
式中：muhat ——参数μ的估计值；
muci ——参数μ估计的置信区间；
sigmahat ——参数σ的估计值；
sigmaci ——参数σ的置信区间。
4 假设检验
K-S检验法的基本思想是：将样本观测值的累积频率（Fn（x））与假设的理论频率分布（
Fx（x））相比较来建立统计量。它可直接处理原始数据，一般认为其功效比Chi-Square检验高，且在样本不大时也可用，命令格式如下：
式中：
P——原假设成立的概率；
KSSTAT——测试统计量的值；
CV ——是否接受假设的临界值。
说明：若H=0则不能拒绝原假设，H=1则可以拒绝原假设。
5 实例应用
下面用实例说明如何使用HIST命令及如何进行参数估计和假设检验。某综放面采煤机无故障运行时间、故障时间统计情况见表1。
据表1可以得出采煤机无故障工作时间data1和采煤机维修时间data2分别为：
data1=[110, 520, 645, 680, 330,75, 95, 480, 360, 575,1065, 170, 415, 15, 20,1275, 270, 90, 1500, 1923,715, 1523, 427, 730, 1120,390, 240, 0, 40, 220,673, 2397, 1032, 315]
data2=[30, 165, 30, 60, 60,10, 50, 50, 130, 15,20, 45, 35, 45, 30,120, 30, 105, 15, 147,15, 78, 25, 90, 30,170, 210, 30, 80, 60,79, 50,125]
这里我们自己划分小区间，把采煤机无故障工作时间data1的区间的中点存放在向量nb1中，把采煤机维修时间data1的区间的中点存放在向量nb2中，二者分别为：
nb1=[200 600 1000 1400 1800 2200 ]；
nb2=[15 45 75 105 135 165 195]；
然后分别使用命令hist(data1，nb1)和hist(data2，nb2)，得出采煤机无故障工作时间直方图如图2所示，采煤机维修时间直方图如图3所示。
据描出的频数直方图的形状，假设采煤机无故障工作时间和采煤机维修时间两者均服从负指数分布。
分别使用命令[muhat1，muci1]＝expfit(data1)和命令[muhat2，muci2]＝expfit(data2)进行参数估计可得：
1/λ＝muhat1＝601.029（min）
1/μ＝muhat2＝67.697（min）
式中：λ——故障率；
μ——修复率。
从而，
λ＝0.0017（1/min）；
μ＝0.0148（1/min）。
令x=data1'，下面的命令可对采煤机无故障工作时间的分布进行假设检验：
[H,p,ksstat,cv]=kstest(x,[x expcdf(x,601)],0.05)
计算结果如下：
H = 0
p =0.9965
ksstat =0.0680
cv =0.2275
对计算结果的说明：
H=0表示接受原假设，统计量ksstat小于临界值cv表示接受原假设。
再令x=data2'，下面的命令可对采煤机维修时间的分布进行假设检验：
[H,p,ksstat,cv]=kstest(x,[x expcdf(x,67)],0.05)
计算结果如下：
H =0
p =0.2166
ksstat =0.1791
cv =0.2308
对计算结果的说明：H=0表示接受原假设，统计量ksstat小于临界值cv表示接受原假设。
6 用MATLAB进行可靠性(Reliability)数据分析的意义
(1)充分利用计算机这一快速、有效的现代工具，减轻科研、管理人员的劳动强度，将信息处理、存储、综合分析等大量的繁重工作自动化；
(2)缩短资料处理周期，及时为生产、管理服务；
(3)综合利用系统信息，提管理工作的科学性和定量程度。
可靠性(Reliability)数据是系统可靠性(Reliability)设计、研究、分析、评定和改进的基本依据，对它的分析则是一切可靠性(Reliability)工作的基础。通过有计划、有目的地收集和分析系统各组成单元的可靠性(Reliability)数据，可定量评定系统的可靠性(Reliability)水平，发现其可靠性(Reliability)的薄弱环节和产生原因，经过改进，使系统的性能与可靠性(Reliability)水平不断提高。所以，可靠性(Reliability)数据的分析在系统可靠性(Reliability)研究中具有十分重要的地位。
收集到的系统可靠性(Reliability)原始数据往往是杂乱无章的，因此不能直接用其对系统的可靠性(Reliability)进行分析，作出判断，必须运用概率论与数理统计知识对其进行一系列整理、分析后，才能透过大量现象掌握它的本质，使其数据显示出规律性。可靠性(Reliability)数据整理分析的步骤如图1所示。
图1 系统可靠性(Reliability)数据整理分析逻辑图
对系统可靠性(Reliability)数据整理分析，需要处理大量的数据，用人工的方法完成这些工作需要花费很长的时间，而且容易出错，而MATLAB的统计工具箱是建立在MATLAB数值计算环境上的工具集，它的统计功能已经可以与SPSS，SAS等软件的统计功能相媲美，它可实现数据资料处理的自动化，保证数据分析工作及时、准确的完成。
2 运用MATLAB绘制可靠性(Reliability)数据分布频数直方图
把系统可靠性(Reliability)数据读入MATLAB的工作空间，便可绘制可靠性(Reliability)数据分布频数直方图。MATLAB统计工具箱提供的hist函数，是用于作直方图的Matlab命令 ，其命令格式如下：
HIST(data，k)；
其中：
Data——原始数据；
K——所分小区间数。
这个命令将区间[min(data),max(data)]分为k等份，并描出频数直方图，k的缺省值为10。如果需要事先给出小区间，则将区间的中点存放在向量nb中，然后使用命令hist(data，nb)。
有时并不需要直方图，只希望知道小区间的位置及数据落在各小区间的频数，可以用命令：
- [n,x]＝HIST（data）；
- [n,x]＝HIST（data，k）；
- [n,x]＝HIST（data，nb）；
3 参数估计
根据绘出的直方图的形状，假设可靠性(Reliability)数据服从某一分布，一般情况可靠性(Reliability)数据服从负指数分布。
指数分布可以通过命令expfit进行参数估计，该命令用极大似然法给出了常用的概率分布的参数的点估计和区域估计值，命令格式如下：
[muhat，muci]＝expfit（data）
式中：muhat——参数μ的估计值；
muci ——参数μ估计的置信区间。
正态分布可以通过命令normfit进行参数估计，该命令用极大似然法给出了常用的概率分布的参数的点估计和区域估计值，命令格式如下：
[muhat，sigmahat，muci，sigmaci]＝normfit（data）
式中：muhat ——参数μ的估计值；
muci ——参数μ估计的置信区间；
sigmahat ——参数σ的估计值；
sigmaci ——参数σ的置信区间。
4 假设检验
K-S检验法的基本思想是：将样本观测值的累积频率（Fn（x））与假设的理论频率分布（
Fx（x））相比较来建立统计量。它可直接处理原始数据，一般认为其功效比Chi-Square检验高，且在样本不大时也可用，命令格式如下：
- [H,P,KSSTAT,CV] = kstest(X,cdf,alpha)
P——原假设成立的概率；
KSSTAT——测试统计量的值；
CV ——是否接受假设的临界值。
说明：若H=0则不能拒绝原假设，H=1则可以拒绝原假设。
5 实例应用
下面用实例说明如何使用HIST命令及如何进行参数估计和假设检验。某综放面采煤机无故障运行时间、故障时间统计情况见表1。
表1 采煤机无故障运行时间、故障时间统计表
据表1可以得出采煤机无故障工作时间data1和采煤机维修时间data2分别为：
data1=[110, 520, 645, 680, 330,75, 95, 480, 360, 575,1065, 170, 415, 15, 20,1275, 270, 90, 1500, 1923,715, 1523, 427, 730, 1120,390, 240, 0, 40, 220,673, 2397, 1032, 315]
data2=[30, 165, 30, 60, 60,10, 50, 50, 130, 15,20, 45, 35, 45, 30,120, 30, 105, 15, 147,15, 78, 25, 90, 30,170, 210, 30, 80, 60,79, 50,125]
这里我们自己划分小区间，把采煤机无故障工作时间data1的区间的中点存放在向量nb1中，把采煤机维修时间data1的区间的中点存放在向量nb2中，二者分别为：
nb1=[200 600 1000 1400 1800 2200 ]；
nb2=[15 45 75 105 135 165 195]；
然后分别使用命令hist(data1，nb1)和hist(data2，nb2)，得出采煤机无故障工作时间直方图如图2所示，采煤机维修时间直方图如图3所示。
图2 采煤机无故障工作时间直方图
图3 采煤机维修时间直方图
据描出的频数直方图的形状，假设采煤机无故障工作时间和采煤机维修时间两者均服从负指数分布。
分别使用命令[muhat1，muci1]＝expfit(data1)和命令[muhat2，muci2]＝expfit(data2)进行参数估计可得：
1/λ＝muhat1＝601.029（min）
1/μ＝muhat2＝67.697（min）
式中：λ——故障率；
μ——修复率。
从而，
λ＝0.0017（1/min）；
μ＝0.0148（1/min）。
令x=data1'，下面的命令可对采煤机无故障工作时间的分布进行假设检验：
[H,p,ksstat,cv]=kstest(x,[x expcdf(x,601)],0.05)
计算结果如下：
H = 0
p =0.9965
ksstat =0.0680
cv =0.2275
对计算结果的说明：
H=0表示接受原假设，统计量ksstat小于临界值cv表示接受原假设。
再令x=data2'，下面的命令可对采煤机维修时间的分布进行假设检验：
[H,p,ksstat,cv]=kstest(x,[x expcdf(x,67)],0.05)
计算结果如下：
H =0
p =0.2166
ksstat =0.1791
cv =0.2308
对计算结果的说明：H=0表示接受原假设，统计量ksstat小于临界值cv表示接受原假设。
6 用MATLAB进行可靠性(Reliability)数据分析的意义
(1)充分利用计算机这一快速、有效的现代工具，减轻科研、管理人员的劳动强度，将信息处理、存储、综合分析等大量的繁重工作自动化；
(2)缩短资料处理周期，及时为生产、管理服务；
(3)综合利用系统信息，提管理工作的科学性和定量程度。