mhsample
Выборка Гастингса столицы
Синтаксис
smpl = mhsample(start,nsamples,'pdf',pdf,'proppdf',proppdf,
'proprnd',proprnd)
smpl = mhsample(...,'symmetric',sym)
smpl = mhsample(...,'burnin',K)
smpl = mhsample(...,'thin',m)
smpl = mhsample(...,'nchain',n)
[smpl,accept] = mhsample(...)
Описание
smpl = mhsample(start,nsamples,'pdf',pdf,'proppdf',proppdf,
'proprnd',proprnd) чертит случайные выборки nsamples от целевого стационарного распределения pdf с помощью алгоритма Гастингса Столицы.
start является вектором - строкой, содержащим значение запуска Цепи Маркова, nsamples является целым числом, задающим количество выборок, которые будут сгенерированы, и pdf, proppdf, и proprnd является указателями на функцию, созданными с помощью @. proppdf задает плотность распределения предложения, и proprnd задает генератор случайных чисел для распределения предложения. pdf и proprnd берут один аргумент в качестве входа с тем же типом и размером как start. proppdf берет два аргумента в качестве входных параметров с тем же типом и размером как start.
smpl является вектор-столбцом или матрицей, содержащей выборки. Если логарифмическая функция плотности предпочтена, 'pdf' и 'proppdf' могут быть заменены 'logpdf' и 'logproppdf'. Функции плотности, используемые в алгоритме Гастингса Столицы, не обязательно нормированы.
q распределения предложения (x, y) дает плотность вероятности для выбора x как следующий вопрос, когда y является текущей точкой. Это иногда пишется как q (x |y).
Если proppdf или logproppdf удовлетворяют q (x, y) = q (y, x), то есть, распределение предложения симметрично, mhsample реализует Случайный Обход выборка Гастингса Столицы. Если proppdf или logproppdf удовлетворяют q (x, y) = q (x), то есть, распределение предложения независимо от текущих значений, mhsample реализует Независимую выборку Гастингса Столицы.
smpl = mhsample(...,'symmetric',sym) чертит случайные выборки nsamples от целевого стационарного распределения pdf с помощью алгоритма Гастингса Столицы. sym является логическим значением, которое указывает, симметрично ли распределение предложения. Значение по умолчанию является ложным, который соответствует асимметричному распределению предложения. Если sym верен, например, распределение предложения симметрично, proppdf и logproppdf являются дополнительными.
smpl = mhsample(...,'burnin',K) генерирует Цепь Маркова со значениями между отправной точкой и k th точка, не использованная в сгенерированной последовательности. Значения вне k th точка сохранены. k является неотрицательным целым числом со значением по умолчанию 0.
smpl = mhsample(...,'thin',m) генерирует Цепь Маркова с m-1 из значений m, не использованных в сгенерированной последовательности. m является положительным целым числом со значением по умолчанию 1.
smpl = mhsample(...,'nchain',n) генерирует Цепи Маркова n с помощью алгоритма Гастингса Столицы. n является положительным целым числом со значением по умолчанию 1. smpl является матрицей, содержащей выборки. Последняя размерность содержит индексы для отдельных цепочек.
[smpl,accept] = mhsample(...) также возвращает accept, пропускную способность предложенного распределения. accept является скаляром, если одна цепочка сгенерирована и является вектором, если кратные цепи сгенерированы.
Примеры
Оцените моменты Используя независимую выборку Гастингса столицы
Используйте Независимую выборку Гастингса Столицы, чтобы оценить момент второго порядка Гамма распределения.
rng default; % For reproducibility alpha = 2.43; beta = 1; pdf = @(x)gampdf(x,alpha,beta); % Target distribution proppdf = @(x,y)gampdf(x,floor(alpha),floor(alpha)/alpha); proprnd = @(x)sum(... exprnd(floor(alpha)/alpha,floor(alpha),1)); nsamples = 5000; smpl = mhsample(1,nsamples,'pdf',pdf,'proprnd',proprnd,... 'proppdf',proppdf);
Постройте график результатов.
xxhat = cumsum(smpl.^2)./(1:nsamples)'; figure; plot(1:nsamples,xxhat)
Случайный обход выборка Гастингса столицы
Используйте Случайный Обход выборка Гастингса Столицы, чтобы сгенерировать выборочные данные от стандартного нормального распределения.
rng default % For reproducibility delta = .5; pdf = @(x) normpdf(x); proppdf = @(x,y) unifpdf(y-x,-delta,delta); proprnd = @(x) x + rand*2*delta - delta; nsamples = 15000; x = mhsample(1,nsamples,'pdf',pdf,'proprnd',proprnd,'symmetric',1);
Постройте выборочные данные.
figure; h = histfit(x,50); h(1).FaceColor = [.8 .8 1];
