2014-10-30 172 views
1

最初的问题是模拟一个24/7全天候使用的灯泡,并且通常会持续25天。一盒灯泡包含12.盒子可能持续一年以上的概率是多少?使用Matlab查找高斯分布的概率

我不得不使用MATLAB来建模一个基于指数变量的高斯曲线。 下面的代码生成一个平均值= 300和std = sqrt(12)* 25的高斯模型。 我不得不使用如此多的不同变量并加上它们的原因是因为我应该证明中心极限定理。高斯曲线表示一盒灯泡持续一天的概率,其中300是盒子将持续的平均天数。

我在使用生成的高斯函数和查找天> 365的概率时遇到了问题。声明1-normcdf(365,300,sqrt(12)* 25)是试图找出概率的预期值,我得到了.2265。 有关如何根据生成的高斯I来查找天> 365的概率的任何提示将不胜感激。

谢谢!

clear all 
samp_num=10000000; 
param=1/25; 
a=-log(rand(1,samp_num))/param; 
b=-log(rand(1,samp_num))/param; 
c=-log(rand(1,samp_num))/param; 
d=-log(rand(1,samp_num))/param; 
e=-log(rand(1,samp_num))/param; 
f=-log(rand(1,samp_num))/param; 
g=-log(rand(1,samp_num))/param; 
h=-log(rand(1,samp_num))/param; 
i=-log(rand(1,samp_num))/param; 
j=-log(rand(1,samp_num))/param; 
k=-log(rand(1,samp_num))/param; 
l=-log(rand(1,samp_num))/param; 
x=a+b+c+d+e+f+g+h+i+j+k+l; 


mean_x=mean(x); 
std_x=std(x); 
bin_sizex=.01*10/param; 
binsx=[0:bin_sizex:800]; 
u=hist(x,binsx); 
u1=u/samp_num; 

1-normcdf(365,300, sqrt(12)*25) 
bar(binsx,u1) 
legend(['mean=',num2str(mean_x),'std=',num2str(std_x)]); 
+0

你的高斯代表什么?这是一个灯泡将持续多少小时?剩余在一个盒子里或其他东西的灯泡数量?请更新您的信息。 – kkuilla 2014-10-30 10:22:05

+0

我相信高斯表示一盒灯泡持续某些天数的概率。 300通常会持续多长时间,平均值。 – azumakazuma 2014-10-30 10:31:11

+0

为什么你的'std = sqrt(12)* 25'?你确定答案.2265是错误的吗? – kkuilla 2014-10-30 11:14:04

回答

1

[f, y]=ecdf(x)x为数据创建一个empirical cdf。然后你可以找到第一次穿过365的概率来得到你的答案。

+0

我试着用[f,y] = ecdf(binsx),它给出了概率与天数的关系图,但是当我看365时,概率大约是4579,并且因为我想> 365, .4579 = .5421,这真的没有道理。如果我看高斯曲线,> 365显然远小于0.5。我不知道我是在做还是在想这个错误。 – azumakazuma 2014-10-30 13:13:17

+0

我得到0.7911 => p = 0.2089。你想找到y> 365的第一个点,而不是f的元素365。 – nivag 2014-10-30 14:18:02

+0

你的ecdf()的输入是什么? – azumakazuma 2014-10-30 14:24:18

0

生成N重复的x,其中N应该是数千或数万。然后p-hat = count(x > 365)/N,并且有一个标准错误sqrt[p-hat * (1 - p-hat)/N]。重复次数越多,估计的误差幅度越小。

当我在JMP上做了这个N = 10,000时,我以[0.2039, 0.2199]作为95%可信区间,结果是一盒灯泡持续时间超过一年。与您的价值0.2265的差异以及10,000个结果的直方图表明,实际分配仍然有些偏差。换句话说,对12个指数的总和使用CLT近似可以给出稍微偏离的答案。