我想知道是否有人能帮我把下面的代码从MatLab翻译成Python。该方程用于确定截断正态分布的99%置信区间。在
function sigma = var_truncNormal( a, b, mu, sigma, data )
x1 = (a-mu)/sigma * normpdf( a, mu, sigma );
x2 = (b-mu)/sigma * normpdf( b, mu, sigma );
cx = normcdf( b, mu, sigma) - normcdf( a, mu, sigma );
yhat = var( data(data>(mu-3000)&data<(mu+3000)) );
sigma2 = yhat/((1+(x1-x2)/cx - ((x1-x2)/cx)^2));
sigma = sqrt( sigma2 );
return;
function ci99 = GetCI99( data )
mu = median( data );
sigma = std( data );
fprintf( 1, 'initial sigma = %.1f\n', sigma );
sigma = var_truncNormal( mu-3000, mu+3000, mu, sigma, data );
fprintf( 1, 'updated sigma = %.1f\n', sigma );
sigma = var_truncNormal( mu-3000, mu+3000, mu, sigma, data );
fprintf( 1, 'updated sigma = %.1f\n', sigma );
ci99 = 2*mu-norminv( 0.01, mu, sigma );
figure( 'visible', 'off' );
hist( data, 5000:200:20000 );
axis( [5000 35000 0 550] );
hold;
[n2, xx] = ksdensity( data, 'npoints', 100 );
plot( xx, n2*length(data)*200, 'r' );
hdl = plot( xx, normpdf( xx, mu, sigma )*length(data)*200, 'k' );
set( hdl, 'linewidth', 2 );
line( [ci99 ci99], [0 550] );
print( '-dpdf', 'testFigure' );
close;
return;
我会很感激你的帮助。在
我想你会很高兴从matlab到python/numpy。你只需要一行一行地过。在
例如函数的第一行:
x1=(a-mu)/西格玛*标准差(a,mu,sigma)
python中的normpd是:1/(sigma*sqrt(2*pi))exp(-1(a-mu)**2/2*sigma**2)。在
所以我们可以在python中定义一个小函数:
Numpy/Scipy有一套很深的统计包,您应该始终检查预先存在的函数,以满足您的需要。在本例中,http://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.norm.html#scipy.stats.normpdf(x,loc=0,scale=1)很接近,但可能不够,因为它被定义为:
^{pr2}$正确的公式是:
你忘了(2*sigma**2)中的括号
顺便说一句:规范.pdf(x,loc,scale)结果相同
相关问题 更多 >
编程相关推荐