『求助』关于使用MATLAB实现可靠性的正态分布
浏览：143054
PS1：这个是脚本函数，无法直接使用MCC编译成为EXE文件，也就是无法脱离MATLAB独立运行。。。网上的解决方法有好几种：使用MATCOM，或者在代码前加FUNCTION，不是解决之道。。。反正我米有成功，有成功的留个言哈~~~^_^
PS2：对于TEST1.DAT，自己随便写两排数据就可以了。第一排是序号，第二排是电灯泡工作寿命（这里是这样就是了。。。）
以下是代码
<code>
clc;clf;clear %清楚原始数据干扰
a=0;
b=15; %初始变量
for k=1:1:15
i=[0.0001 0.001 0.01 0.05 0.1 0.2 0.4 0.5 0.6 0.8 0.9 0.95 0.99 0.999 0.9999];
y1=norminv(i(k));
x1=[a,b];
plot(x1,[y1,y1])
hold on
end
set(gca,'ytick',norminv(i(:))');
set(gca,'yticklabel',i);
xlabel('时间');
ylabel('概率');
title('正态坐标纸'); %绘制坐标图
load test1.dat %数据导入
m=size(test1);n=m(2); %求数据输入的数目
for k=1:1:n
l(k)=(k-0.3)/(n+0.4); %求失效概率中位秩
end
for i=1:1:n
plot(test1(2,i),norminv(l(i)),'o') %描点
hold on
end
x=test1(2,:); %赋值给x
y=norminv(l(:))'; %赋值给y
avgx=sum(x)/n; %x平均数
avgy=sum(y)/n; %y平均数
B=(sum(x.*y)-n*avgx*avgy)/(sum(x.*x)-n*avgx*avgx);A=avgy-B*avgx; %最小二乘法求线性直线
for p=0:0.01:15
plot(p,A+B*p)
end %绘制直线
rol=(sum(x.*y)-n*avgx*avgy)/(sqrt((sum(x.*x)-n*avgx*avgx)*(sum(y.*y)-n*avgy*avgy))); %求相关系数
rola=2.576/sqrt(n+1); %求标准相关系数
if rol>=rola %判断是否线性相关
disp('两者线性相关')
cu='期望等于 ';cr='标准差等于 ';cum='Kh';
cuu=(0.5-A)/B;crr=(-A+0.841)/B-(-A+0.5)/B;
disp([cu,num2str(cuu),cum])
disp([cr,num2str(crr)]) %参数估计
else
disp('两者线性不相关') %非线性相关，quit
end
a=0;
b=15; %初始变量
for k=1:1:15
i=[0.0001 0.001 0.01 0.05 0.1 0.2 0.4 0.5 0.6 0.8 0.9 0.95 0.99 0.999 0.9999];
y1=norminv(i(k));
x1=[a,b];
plot(x1,[y1,y1])
hold on
end
set(gca,'ytick',norminv(i(:))');
set(gca,'yticklabel',i);
xlabel('时间');
ylabel('概率');
title('正态坐标纸'); %绘制坐标图
load test1.dat %数据导入
m=size(test1);n=m(2); %求数据输入的数目
for k=1:1:n
l(k)=(k-0.3)/(n+0.4); %求失效概率中位秩
end
for i=1:1:n
plot(test1(2,i),norminv(l(i)),'o') %描点
hold on
end
x=test1(2,:); %赋值给x
y=norminv(l(:))'; %赋值给y
avgx=sum(x)/n; %x平均数
avgy=sum(y)/n; %y平均数
B=(sum(x.*y)-n*avgx*avgy)/(sum(x.*x)-n*avgx*avgx);A=avgy-B*avgx; %最小二乘法求线性直线
for p=0:0.01:15
plot(p,A+B*p)
end %绘制直线
rol=(sum(x.*y)-n*avgx*avgy)/(sqrt((sum(x.*x)-n*avgx*avgx)*(sum(y.*y)-n*avgy*avgy))); %求相关系数
rola=2.576/sqrt(n+1); %求标准相关系数
if rol>=rola %判断是否线性相关
disp('两者线性相关')
cu='期望等于 ';cr='标准差等于 ';cum='Kh';
cuu=(0.5-A)/B;crr=(-A+0.841)/B-(-A+0.5)/B;
disp([cu,num2str(cuu),cum])
disp([cr,num2str(crr)]) %参数估计
else
disp('两者线性不相关') %非线性相关，quit
end
</code>
