현재 Metropolis-Hastings 알고리즘 및 일부 수치 예제에 대한 개요를 기반으로하는 수학 학위 과정의 최종 프로젝트를 진행하고 있습니다. 지금까지 Gaussian으로 제안서 배포판을 사용하고 몇 가지 다른 배포판을 샘플링하여 좋은 결과를 얻었지만 다른 제안 배포판을 사용하여 한 걸음 더 나아가려고합니다.Metropolis-Hastings 알고리즘, Matlab의 가우시안 이외의 제안 분포 사용
지금까지 나는 다른 제안을 사용하는 것에 관해 온라인으로 제한된 리소스로이 코드를 사용했습니다 (Matlab을 사용하고 있습니다). 현실에서 나는 시도하는 방법을 너무 확신하지는 않습니다. 이것은 (특히 이것은 지금까지 유용한 데이터 출력을 제공하지 못함).
쉽게 접근 할 수있는 정보를 알고 있거나 전달하면 누군가가 손을 빌려주는 것이 환상적 일 것입니다. (나는 코딩 조언뿐만 아니라 수학도 알고 있습니다.)
그래서, 라플라스, 지금까지 내 코드의 제안 분포를 사용하여 가우시안에서 샘플링 할 :
n = 1000; %%%%number of iterations
x(1) = -3; %%%%Generate a starting point
%%%%Target distribution: Gaussian:
strg = '1/(sqrt(2*pi*(sig)))*exp(-0.5*((x - mu)/sqrt(sig)).^2)';
tnorm = inline(strg, 'x', 'mu', 'sig');
mu = 1; %%%%Gaussian Parameters (I will be estimating these from my markov chain x)
sig = 3;
%%%%Proposal distribution: Laplace:
strg = '(1/(2*b))*exp((-1)*abs(x - mu)/b)';
laplace = inline(strg, 'x', 'b', 'mu');
b = 2; %%%%Laplace parameter, I will be using my values for y and x(i-1) for mu
%%%%Generate markov chain by acceptance-rejection
for i = 2:n
%%%%Generate a candidate from the proposal distribution
y = laplace(randn(1), b, x(i-1));
%%%%Generate a uniform for comparison
u = rand(1);
alpha = min([1, (tnorm(y, mu, sig)*laplace(x(i-1), b, y))/(tnorm(x(i-1), mu, sig)*laplace(y, b, x(i-1)))]);
if u <= alpha
x(i) = y;
else
x(i) = x(i-1);
end
end
위에서 그것에 대해가는/완전히 잘못되면 누군가는 말해 줄 수있는 경우 틀린 길, 또는 단지 몇 가지 실수가 있습니다 (환상적인 for 루프가 완전히 잘못되었다는 것에 대한 제 세대에 대해서는 매우 조심 스럽습니다).
감사합니다, 톰 참고로