余弦波的贝叶斯拟合比预期花费更长的时间 [英] Bayesian fit of cosine wave taking longer than expected

问题描述

在最近的作业中，我被要求使用Metropolis算法对一组数据 a 和 b 进行贝叶斯拟合.给出 a 和 b 之间的关系:

In a recent homework, I was asked to perform a Bayesian fit over a set of data a and b using a Metropolis algorithm. The relationship between a and b is given:

e(t)= e_0 * cos(w * t)

e(t) = e_0*cos(w*t)

w = 2 * pi

w = 2 * pi

Metropolis算法是(在其他拟合中也可以正常工作):

The Metropolis algorithm is (it works fine with other fit):

def metropolis(logP, args, v0, Nsteps, stepSize):

    vCur = v0
    logPcur = logP(vCur, *args)
    v = []
    Nattempts = 0

    for i in range(Nsteps):
        while(True):
            #Propose step:
            vNext = vCur + stepSize*np.random.randn(*vCur.shape)
            logPnext = logP(vNext, *args)
            Nattempts += 1

            #Accept/reject step:
            Pratio = (1. if logPnext>logPcur else np.exp(logPnext-logPcur))

            if np.random.rand() < Pratio: 
                vCur = vNext
                logPcur = logPnext
                v.append(vCur)
                break

    acceptRatio = Nsteps*(1./Nattempts)
    return np.array(v), acceptRatio

我尝试使贝叶斯拟合余弦波并使用上面的Metropolis算法:

I have tried to Bayesian fit the cosine wave and used the Metropolis algorithm above:

e_0 = -0.00155 

def strain_t(e_0,t):
    return e_0*np.cos(2*np.pi*t)

data = pd.read_csv('stressStrain.csv')
t = np.array(data['t'])
e = strain_t(e_0,t)

def logfitstrain_t(params,t,e):
    e_0 = params[0]
    sigmaR = params[1]
    strainModel = strain_t(e_0,t)
    return np.sum(-0.5*((e-strainModel)/sigmaR)**2 - np.log(sigmaR))


params0 = np.array([-0.00155,np.std(t)]) 
params, accRatio = metropolis(logfitstrain_t, (t,e), params0, 1000, 0.042)
print('Acceptance ratio:', accRatio)

e0 = np.mean(params[0])
print('e0=',e0)

e_t = e0*np.cos(2*np.pi*t)
sns.jointplot(t, e_t, kind='hex',color='purple')

.csv中的数据看起来像

The data in .csv looks like

在我运行run之后没有任何错误消息显示，但是python给我输出输出需要花费很多时间. 我在这里做什么错了?

There isn't any error message showing after I hit run, but it takes forever for python to give me an output. What did I do wrong here?

推荐答案

为什么永远收拾"

您的算法被设计为可以运行直到它接受给定数量的提案(示例中为1000).因此，如果运行了很长时间，您可能会拒绝一堆建议.当步长太大时，可能会发生这种情况，导致新建议最终出现在似然空间的遥远，低概率区域中. 尝试减小步长.这可能还需要您增加样本数量，以确保充分探究后方空间.

Why it might "take forever"

Your algorithm is designed to run until it accepts a given number of proposals (1000 in the example). Thus, if it's running for a long time, you're likely rejecting a bunch of proposals. This can happen when the step size is too large, leading new proposals to end up in distant, low probability regions of the likelihood space. Try reducing your step size. This may require you to also increase the number of samples to ensure the posterior space becomes adequately explored.

因为您只将接受的提案附加到链v，所以您实际上并没有实现Metropolis算法，而是获得了一组有偏见的样本，这些样本往往会过度代表不大可能的区域.后空间.每当新提案被拒绝时，真正的Metropolis实施都会重新附加上一个提案.您仍然可以强制执行最少数量的已接受提案，但实际上每次都必须添加一些内容.

Because you only append accepted proposals to the chain v, you haven't actually implemented the Metropolis algorithm, and instead obtain a biased set of samples that will tend to overrepresent less likely regions of the posterior space. A true Metropolis implementation re-appends the previous proposal whenever the new proposal is rejected. You can still enforce a minimum number of accepted proposals, but you really must append something every time.

这篇关于余弦波的贝叶斯拟合比预期花费更长的时间的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！