0

我最近一直在尝试将一段 Matlab 代码转换为 Python 代码。

我已经做出了我需要的大部分更改，但是我遇到的问题是它说的那一行：

y(index(m)) = 1-x(index(m));

我得到错误：

“无法分配给函数调用”

但是，我不确定如何对其进行重组以消除此错误。

我环顾四周，人们提到“获取项目”和“设置项目”，但是我尝试使用它们，但我无法让它们工作（可能是因为我无法弄清楚结构）

这是完整的代码：

import numpy

N = 100;
B = N+1;
M = 5e4;
burnin = M;
Niter = 20;
p = ones(B,Niter+1)/B;
hit = zeros(B,1);

for j in range(1,Niter):
    x = double(rand(1,N)>0.5);
    bin_x = 1+sum(x);
    index = ceil(N*rand(1,M+burnin));
    acceptval = rand(1,M+burnin);
    for m in range(1,M+burnin):
        y = x;
        y(index(m)) = 1-x(index(m));
        bin_y = 1+sum(y);

        alpha = min(1, p(bin_x,j)/p(bin_y,j) );
        if acceptval(m)<alpha:
            x = y; bin_x = bin_y;
        end

        if m > burnin: hit(bin_x) = hit(bin_x)+1; end
    end

    pnew = p[:,j];
    for b in range(1,B-1):
        if (hit(b+1)*hit(b) == 0):
            pnew(b+1) = pnew(b)*(p(b+1,j)/p(b,j));
        else:
            g(b,j) = hit(b+1)*hit(b) / (hit(b+1)+hit(b));
            g_hat(b) = g(b,j)/sum(g(b,arange(1,j)));
            pnew(b+1) = pnew(b)*(p(b+1,j)/p(b,j))+((hit(b+1)/hit(b))^g_hat(b));
        end
    end
    p[:,j+1] = pnew/sum(pnew);
    hit[:] = 0;
end

提前致谢

4

1 回答 1

3

圆括号()表示一个函数。对于索引，您需要[]方括号-但这只是许多错误中的第一个……我目前正在逐行检查，但这需要一段时间。

这段代码至少可以运行......你需要弄清楚索引是否符合你的预期，因为 Python 数组是从零开始索引的，而 Matlab 数组从 1 开始。我试图在几个地方解决这个问题，但没有t 逐行遍历 - 那是调试。

一些关键的学习：

  • 没有end声明...只是停止缩进
  • 导入库时，需要引用它（numpy.zeros, not zeros）
  • 列表从零开始索引，而不是一
  • 索引是用 完成的[]，而不是()
  • 创建一个随机数数组是用 完成的[random.random() for r in xrange(N)]，而不是random(N)。
  • ...以及您在查看下面的代码时会发现的许多其他内容。

祝你好运！

import numpy
import random

N = int(100);
B = N+1;
M = 5e4;
burnin = M;
Niter = 20;
p = numpy.ones([B,Niter+1])/B;
hit = numpy.zeros([B,1]);
g = numpy.zeros([B, Niter]);
b_hat = numpy.zeros(B);

for j in range(1,Niter):
    x = [float(random.randint(0,1)>0.5) for r in xrange(N)];
    bin_x = 1+sum(x);
    index = [random.randint(0,N-1) for r in xrange(int(M+burnin))];
    #acceptval = rand(1,M+burnin);
    acceptval = [random.random() for r in xrange(int(M+burnin))];
    for m in range(1,int(M+burnin)):
        y = x;
        y[index[m]] = 1-x[index[m]];
        bin_y = 1+sum(y);

        alpha = min(1, p[bin_x,j]/p[bin_y,j] );
        if acceptval[m]<alpha:
            x = y; bin_x = bin_y;

        if m > burnin: 
            hit[bin_x] = hit[bin_x]+1;

    pnew = p[:,j];
    for b in range(1,B-1):
        if (hit[b+1]*hit[b] == 0):
            pnew[b+1] = pnew[b]*(p[b+1,j]/p[b,j]);
        else:
            g[b,j] = hit[b+1]*hit[b] / [hit[b+1]+hit[b]];
            g_hat[b] = g[b,j]/sum(g[b,numpy.arange(1,j)]);
            pnew[b+1] = pnew[b]*(p[b+1,j]/p[b,j])+((hit[b+1]/hit[b])^g_hat[b]);
    p[:,j+1] = pnew/sum(pnew);
    hit[:] = 0;
于 2013-03-13T15:38:41.327 回答