Podsećanje

MCMC je kompjuterska metoda uzorkovanja koja omogućava da se karakterizuje raspodela bez poznavanja svih matematičkih svojstava raspodele, slučajnim uzorkovanjem vrednosti iz raspodele. Posebna snaga MCMC metoda je što se može koristiti za vađenje uzoraka iz raspodela, čak i kada je sve što je znamo o raspodeli kako izračunati gustinu za različite uzorke. Ime MCMC kombinuje dva važna pojma: Monte Carlo i Markovske lance.

Monte Carlo je način procenjivanja svojstava raspodele ispitivanjem slučajnih uzoraka iz raspodele. Na primer, umesto da tražimo srednju vrednost raspodele direktnim računanjem iz jednačina raspodele, Monte Carlo pristup bi bio da se izvuče veliki broj slučajnih uzoraka iz raspodele, a onda izračuna njihova srednja vrednost.

Prednost Monte Carlo pristupa je jasna: izračunavanje srednje vrednosti velikog uzorka brojeva je često mnogo lakše nego izračunavanje srednje vrednosti iz jednačina raspodele.

Ova prednost je najizraženija kada se lako vade slučajni uzorci iz raspodele, i kada je teško raditi sa jednačinama raspodele.

Deo koji se odnosi na Markovske lance u MCMC je ideja da se slučajni uzorci generišu posebnim nizovnim procesom. Svaki slučajan uzorak se koristi kao osnova da se generiše sledeći slučajan uzorak (odatle i lanac). Posebno svojstvo lanca je da, dok svaki novi uzorak zavisi od onog pre njega, novi uzorci ne zavise ni od jednog uzorka pre svog prethodnika (ovo je svojstvo Markova).

MCMC je naročito koristan u Bajesovskom pristupu jer je fokus na aposteriornoj raspodeli, koju je često teško odrediti analitičkim ispitivanjem. U ovim slučajevima, MCMC dozvoljava da se aproksimiraju karakteristike aposteriorne raspodele koje se ne mogu direktno izračunati (npr. slučajni uzorci iz aposteriorne raspodele, srednja vrednost aposteriorne raspodele, itd.).

Bajesovski pristup koristi se da se na osnovu informacija prikupljenih iz dobijenih podataka, unaprede prethodna verovanja o parametru (ili skupu parametara). Dobijamo aposteriornu raspodelu, tj. verovanja o parametru (ili skupu parametara) nakon što imamo nove informacije dobijene iz podataka. Formalno, Bajesovo pravilo se definiše na sledeći način: \[ p(\mu|D)\propto p(D|\mu)\cdot p(\mu)\] gde \(\mu\) označava parametar (ili skup parametara) od interesa, \(D\) označava podatke, \(p(\mu|D)\) označava aposteriornu verovatnoću od \(\mu\) za date podatke, \(p(D|\mu)\) označava verovatnoću dobijanja datih podataka ako je dat parametar \(\mu\), a \(p(\mu)\) predstavlja apriornu raspodelu od \(\mu\). (Simbol \(\propto\) predstavlja proporcionalnost.)

Dakle, način na koji se podaci koriste da unaprede naša prethodna verovanja o parametru (odnosno skupu parametara) je da se ispita verovatnoća dobijanja datih podataka za dati skup parametara od interesa. U idealnom slučaju, hteli bismo da procenimo ovu verovatnoću za svaku kombinaciju vrednosti parametara. Kada je dostupan analitički izraz, može se kombinovati sa apriornom raspodelom da bi se dobio analitički izraz za aposteriornu raspodelu. Međutim, u praksi često ne možemo dobiti lako taj analitički izraz. U Bajesovskom pristupu, ovaj problem se najčešće rešava putem MCMC-a: vadi se niz uzoraka iz aposteriorne, pa se računa njihova sredina, domen, itd.

Iako MCMC može zvučati komplikovano kada je apstraktno opisan, njegova praktična implementacija je uglavnom vrlo jednostavna.

Primer: In-class test

Pretpostavimo da nekog profesora interesuje srednja vrednost rezultata testa u populaciji studenata. Iako mu je srednja vrednost rezultata testa nepoznata, profesor zna da su rezultati testa normalno raspodeljeni sa standardnom devijacijom 15. Do sad, profesor je video rezultat testa za jednog studenta: 100. Takođe, profesor veruje da srednja vrednost ima \(\mathcal{N}(100,15^2)\) raspodelu.

On može iskoristiti MCMC za vađenje uzoraka iz ciljne raspodele, (u ovom slučaju aposteriorne), koja predstavlja verovatnoću svake moguće vrednosti populacijske sredine, za dati ovaj jedan podatak. Ovo je (preterano) pojednostavljen primer, jer postoji analitički izraz za aposteriornu raspodelu (\(\mathcal{N}(100,\frac{15^2}{2})\)), ali njegova svrha je da ilustruje MCMC.

Da bi vadio uzorke iz aposteriorne raspodele, MCMC počinje sa inicijalnom vrednošću (početnom pretpostavkom). Pretpostavimo da je inicijalna vrednost 110. MCMC se zatim koristi za pravljenje lanca uzoraka iz ove početne pretpostavke. Svaki novi uzorak se dobija u dva jednostavna koraka:

  1. Kreira se predlog za novi uzorak dodavanjem malog slučajnog šuma najskorijem uzorku.

  2. Ovaj novi predlog se ili prihvata kao novi uzorak, ili se odbacuje (u tom slučaju se zadržava stari uzorak).

Postoji mnogo načina za dodavanje slučajnog šuma pri kreiranju predloga, kao i različiti pristupi procesu prihvatanja i odbacivanja. Pominjali smo Metropolis algoritam, koji je jedan od najjednostavnih MCMC pristupa:

  1. Počnemo sa verodostojnom početnom vrednošću, 110 u ovom primeru.

  2. Generišemo novi predlog tako što uzimamo poslednji uzorak (110) i dodamo neki slučajan šum. Ovaj slučajni šum se generiše iz predložene raspodele, koja bi trebalo da bude simetrična i centrirana oko 0. U ovom primeru ćemo koristiti predloženu raspodelu koja je normalna sa parametrima 0 (srednja vrednost) i 5 (standardna devijacija). Ovo znači da je novi predlog 110 (poslednji uzorak) plus slučajan uzorak iz \(\mathcal{N}(0,5^2)\) raspodele. Pretpostavimo da to rezultira predlogom 108.

  3. Poredimo “visinu” aposteriorne raspodele novog predloga sa “visinom” za poslednji uzorak. Pošto je ciljna raspodela normalna sa očekivanjem 100 i disperzijom \(\frac{15^2}{2}\), to znači da poredimo \(\mathcal{N}(108|100,\frac{15^2}{2})\) protiv \(\mathcal{N}(110|100,\frac{15^2}{2})\). Ove dve verovatnoće govore nam koliko su predlog i najskoriji uzorak prihvatljivi za ciljnu raspodelu.

  4. Ako novi predlog ima veću aposteriornu verovatnoću od najskorijeg uzorka, onda prihvatamo novi predlog.

  5. Ako novi predlog ima manju aposteriornu verovatnoću od najskorijeg uzorka, na slučajan način se odlučuje da li se prihvata ili odbija novi predlog, sa verovatnoćom koja je jednaka odnosu aposteriornih verovatnoća za te dve vrednosti.

  6. Ako je predlog prihvaćen, on postaje sledeći uzorak u MCMC lancu, u suprotnom sledeći uzorak u MCMC lancu je samo kopija poslednjeg uzorka.

  7. Time se završava jedna MCMC iteracija. Sledeća iteracija se izvršava vraćanjem na korak 2.

  8. Stajemo kada dobijemo dovoljno uzoraka (npr. 500). Odluka o tome koliko je uzoraka dovoljno je posebno pitanje, o čemu će kasnije biti govora.

Ovaj jednostavni problem uzimanja uzoraka zahteva samo nekoliko redova koda u R-u. Rezultati pokretanja ovog uzorkovanja prikazani su u levoj koloni na sledećoj slici. Ovi se uzorci sada mogu koristiti u Monte Carlo svrhe. Na primer, srednja vrednost rezultata testa za populaciju studenata može se proceniti računanjem uzoračke srednje vrednosti od 500 dobijenih uzoraka.

samples = rep(0, 500) # kreiramo vektor koji ce da sadrzi uzorke
samples[1] = 110 # pocetna pretpostavka

# petlja ponavlja proces generisanja vrednosti predloga i odredjuje da li se vrednost predloga prihvata 
# ili se zadrzava trenutna vrednost
for(i in 2:500) {
  proposal = samples[i-1] + rnorm(1,0,5) #vrednost predloga
  
  r = dnorm(proposal,100,15/sqrt(2))/dnorm(samples[i-1],100,15/sqrt(2))
  
  if(r > runif(1))  
    samples[i] = proposal #prihvatamo predlog
  else (samples[i] = samples[i-1]) #odbacujemo predlog
}

Leva kolona: Lanac uzorkovanja koji počinje od dobre početne vrednosti Srednja kolona: Lanac uzorkovanja koji počinje od početne vrednosti u repovima prave raspodele. Desna kolona: Lanac uzorkovanja koji počinje od vrednosti daleko od prave raspodele.

Gornji red: Markovski lanac. Donji red: gustina uzorka. Analitička (prava) raspodela je označena isprekidanom linijom.

Gornji levi grafik na slici prikazuje 500 iteracija: ovo je Markovski lanac. Vrednosti uzorka su centrirane oko srednje vrednosti uzorka (\(100\)), ali takođe, ima i vrednosti koje su manje uobičajene. Donji levi grafik pokazuje gustinu vrednosti iz uzorka. Ponovo, vrednosti su centrirane oko sredine uzorka sa devijacijom koja je vrlo blizu populacijskoj disperziji od \(\frac{15^2}{2}=112.5\) (uzoračka disperzija za ovaj Markovski lanac je 114.71). Dakle, MCMC metoda je zadržala suštinu prave populacijske raspodele sa relativno malim brojem slučajnih uzoraka.

Ograničenja

Prisetimo se nekih ograničenja ovog algoritma, i pogledajmo kako ih možemo rešiti.

MCMC algoritam pruža moćan alat za vađenje uzoraka iz raspodele, kada je direktno uzorkovanje teško. Metoda će “raditi” (tj. raspodela uzoraka će zaista biti ciljna raspodela) dok god su zadovoljeni određeni uslovi.

Prvo, vrednosti verovatnoće izračunate u koracima \(4\) i \(5\), (da bismo prihvatili ili odbacili novi predlog), moraju precizno odražavati gustinu predloga i trenutne vrednosti u ciljnoj raspodeli.

Kada se MCMC primeni na Bajesov zaključak, to znači da izračunate vrednosti moraju biti aposteriorne verovatnoće, ili bar proporcionalne aposteriornim verovatnoćama (tj. izračunati odnos verovatnoća mora biti tačan).

Drugo, ako želimo da koristimo ovaj metod - Metropolis, raspodela predloga treba da bude simetrična (ili, ako se koristi asimetrična raspodela, onda koristimo složeniji Metropolis-Hastings algoritam).

Treće, pošto bi početna pretpostavka mogla biti veoma pogrešna, prvi deo Markovljevog lanca treba ignorisati - ne može se garantovati da će ovi rani uzorci biti izvučeni iz ciljne raspodele. Na ovo ćemo se vratiti kasnije.

Primer MCMC algoritma odozgo vadi uzorke iz normalne normalne \(\mathcal{N}(0,5^2)\) raspodele. Teorijski, bilo koja simetrična raspodela bi radila isto toliko dobro, ali u praksi izbor raspodele predloga može značajno da utiče na performanse samog uzorkovanja.

Ovo se može uočiti zamenjivanjem vrednosti devijacije za raspodelu predloga u gornjem primeru veoma velikom vrednošću, na primer 50. Tada bi mnogi predlozi bili daleko izvan ciljne raspodele (npr. predlozi sa negativnim rezultatima testova!) što bi dovelo do velike stope odbacivanja. S druge strane, za vrlo malu disperziju, kao što je 1, moralo bi se uzeti mnogo iteracija za konvergenciju od početne vrednosti do ciljne raspodele. Takođe, postoji rizik da se zaglavimo u oblastima lokalnih maksimuma: oblastima gde je verovatnoća veća za određenu vrednost nego za njene bliske susede, ali niža nego za neke udaljenije vrednosti.

Širina raspodele predloga se ponekad naziva parametar podešavanja (tuning parameter) ovog MCMC algoritma. Činjenica da praktični učinak uzorkovanja može da zavisi od vrednosti parametra podešavanja je jedno ograničenje standardnog Metropolis-Hastings algoritma za uzorkovanje, pa postoje mnoge proširene metode koje rešavaju ovaj problem. Na primer, auto-tuning algoritmi koji prilagođavaju širinu raspodele predloga prirodi podataka i raspodele.

Treći uslov, činjenica da početne uzorke treba ignorisati jer bi mogli biti veoma pogrešni, već smo spominjali, ovo predstavlja problem koji je poznat kao burn-in. Na primer, pretpostavimo da je početna pretpostavka bila takva da je jako malo verovatno da dolazi iz ciljne raspodele, kao što je na primer rezultat od \(250\) ili čak \(650\). Markovski lanci koji počinju ovim vrednostima, prikazani su u srednjoj i desnoj koloni prvog reda na gornjoj slici. Srednji grafik u prvom redu prikazuje Markovski lanac koji brzo ide ka pravoj aposteriornoj vrednosti. Nakon samo 80 iteracija, lanac se centrira oko prave populacijske sredine. Gornji desni grafik, gde je početna tačka još ekstremnija, pokazuje da je broj iteracija potreban da bi se došlo do prave populacijske sredine (oko 300) mnogo veći nego za bolje odabrane početne tačke. Ova dva primera jasno pokazuju da se za prvih nekoliko iteracija u bilo kom lancu ne može sa sigurnošću pretpostaviti da se izvlače iz ciljne raspodele. Na primer, uključivanje prvih 80 iteracija (u gornjem srednjem slučaju), ili prvih 300 iteracija (gornji desni) dovodi do pogrešnog zaključka o raspodeli populacije, koja je prikazana na slikama ispod.

Jedan od načina za ublažavanje ovog problema je bolji odabir početnih tačaka. Početne vrednosti koje su bliže modi aposteriorne raspodele će obezbediti brži burn-in i manje problema sa konvergencijom. U praksi može biti teško pronaći polazne tačke bliske modi aposteriorne raspodele, ali ocena maksimalne verodostojnosti (ili neka aproksimacija toga) može biti korisna u identifikaciji dobrih kandidata. Drugi pristup je korišćenje više lanaca: da se izvrši uzorkovanje mnogo puta sa različitim početnim vrednostima (npr. sa početnim vrednostima uzorkovanim iz apriorne raspodele). Razlike u raspodelama uzoraka iz različitih lanaca mogu ukazati na probleme sa burn-in-om i konvergencijom.

Drugi element rešavanja problema je uklanjanje ranih uzoraka: onih uzoraka iz nestacionarnih delova lanca. Kada ponovo pogledamo lance u gornjem redu na slici odozgo, možemo videti da je lanac u levom uglu došao u neku vrstu ravnoteže. Lanci u sredini i desno (gornji red) takođe konvergiraju, tj. ulaze u stanje ravnoteže, ali tek posle otprilike 80, odnosno 300 iteracija, redom. Prema tome, važno je napomenuti da svi uzorci pre konvergencije nisu uzorci iz ciljne raspodele i moraju biti odbačeni.

Odlučivanje o tome u kojoj tački lanac ulazi u ravnitežno stanje može biti teško, pa je često razlog zbunjenosti pri početku rada sa MCMC-om. Važan aspekt burn-in-a je razumevanje post-hoc prirode odluke - odluke o burn-in-u moraju biti donete nakon uzorkovanja, i nakon posmatranja lanca. Dobra ideja je da budemo oprezni: odbacivanje dodatnih uzoraka je bezbedno, jer je verovatnije da su preostali uzorci iz ravnotežnih delova lanca. Jedino ograničenje ovog opreznog pristupa je da mora imati dovoljno uzoraka nakon burn-in-a da bi se osigurala adekvatna aproksimacija raspodele.

MCMC primenjen na kognitivni model

Recimo da smo zainteresovani za parametre kognitivnih modela iz podataka o ponašanju. Kao što smo naveli u uvodu, MCMC metode obezbeđuju odličan način procenjivanja parametara u Bajesovskom pristupu.

Primer takvih kognitivnih modela su modeli zasnovani na detekciji signala (SDT). Modeli zasnovani na SDT su imali značajnu istoriju u kognitivnoj nauci, možda delimično zbog njihove intuitivne psihološke privlačnosti i računske jednostavnosti. Računska jednostavnost SDT-a čini ga dobrim kandidatom za procenu parametara preko MCMC-a.

Pretpostavimo da istraživač dobija podatke u obliku ispravno i pogrešno prepoznatih signala iz jednostavnog vizuelnog detekcionog eksperimenta. Ocena parametara SDT modela omogućava istraživaču da stekne uvid u način na koji ljudi donose odluke pod neizvešnošću.

Sproveden je eksperiment, i ispitanici su trebali da odlučuju o tome da li određeni obrazac prepoznaju kao signal (npr. sirena u bučnoj sredini) ili ne (samo buka).

Parametri SDT-a omogućavaju teorijsko razumevanje toga kako ljudi prepoznaju signale: osetljivost (ili \(d'\)) daje meru sposobnosti pojedinca da pravi razliku između signala i samo npr. buke; pristrasnost (ili \(C\)) meri pristrasnost pojedinca, na kom nivou su voljni da nešto prepoznaju kao signal.

Jedan od načina za procenjivanje SDT parametara iz podataka bi bio da se koristi Bajesov zaključak i ispita aposteriorna raspodela ovih parametara. Pošto SDT model ima dva parametra (\(d'\) i \(C\)), aposteriorna raspodela je dvodimenzionalna, tj. aposteriorna raspodela je definisana na skupu svih mogućih kombinacija vrednosti \(d'\) i \(C\). Željena raspodela je data jednačinom: \[p(d',C|D) \propto p(D|d',C)\cdot p(d',C)\]

Dakle, verovatnoća pogodaka i promašaja, za date SDT parametre, pomnožena sa prethodnim uverenjima za te SDT parametre. Sa ovom jednačinom, proces MCMC uzorkovanja iz aposteriorne raspodele za \(d'\) i \(C\) je relativno jednostavan, i zahteva samo malo promene u odnosu na primer testa odozgo.

Najvažnija promena je to što je sada lanac uzorkovanja dvodimenzionalan: svaki uzorak u Markovskom lancu sadrži dve vrednosti - jednu za \(d'\) i jednu za \(C\).

Ciljna raspodela je aposteriorna raspodela za parametre. Ovo omogućava istraživaču da odgovori na ključna pitanja, na primer, da li je \(d'\) pouzdano veće od nule, ili da li je \(C\) pouzdano različito od neke vrednosti. Da bi ciljna raspodela bila aposteriorna raspodela ovih parametara, odnos verovatnoća iz trećeg koraka mora se izračunati korišćenjem Bajesovske jednačine. Jednostavan primer takvog MCMC uzorkovanja za SDT model je sledeći:

# Metropolis uzorkovanje za procenu parametara SDT modela
# kod procenjuje zajednicku verovatnocu SDT parametara d' i C
# Novi predlozi za oba parametra se uzorkuju i procenjuju istovremeno

N.present = 100 # broj slucajeva kada je signal bio prisutan
N.absent = 100 # kada je signal bio odsutan
N.hits = 85 # tacni odgovori na prisutne signale
N.falsealarams = 12 # netacni odgovori kada nije bilo signala

posterior.density = function(parameters, h, fa, Np, Na) {
  # funkcija koja računa aposteriornu gustinu
  # "parameters" je dvodimenzionalni vektor, sa komponentama 
  # d'(osetljivost) i C (pristrasnost)
  # h i fa su brojaci pogodaka i promasaja
  # Np i Na su brojevi slucajeva sa signalom i bez signala
  
  # modelom predvidjena verovatnoca laznog alarma
  prob.fa = pnorm(-parameters["C"])
  # modelom predvidjena verovatnoca pogotka
  prob.h = pnorm(parameters["d'"]-parameters["C"])
  
  # log-likelihood opservacija fa laznih alarma
  loglike.fa = dbinom(x = fa, size = Na, prob = prob.fa, log = TRUE)
  # log-likelihood dobijanja h pravih
  loglike.h = dbinom(x = h, size = Np, prob = prob.h, log = TRUE)
  
  # logaritam apriorne verovatnoce parametara 
  # apriorna za oba parametra N(0,4)
  loglike.prior = dnorm(parameters, mean = 0, sd = 4, log = TRUE)
  
  # vracamo aposteriornu raspodelu
  exp(loglike.fa + loglike.h + sum(loglike.prior))
}

# broj uzoraka
nmc = 500

# vektor u koji smestamo uzorke
samples = array(dim = c(2, nmc), dimnames = list(c("C", "d'"), NULL))

# pocetna pretpostavka
samples[, 1] = c(0.5, 1)

# uzorkovanje
for(i in 2:nmc) {
  proposal = samples[, i-1] + rnorm(n = 2, mean = 0, sd = 0.1)
  
  new.likelihood = posterior.density(parameters = proposal, h = N.hits, fa = N.falsealarams, Np = N.present, Na = N.absent)
  
  old.likelihood = posterior.density(parameters = samples[, i-1], h = N.hits, fa = N.falsealarams, Np = N.present, Na = N.absent)
  
  likelihood.ratio = new.likelihood/old.likelihood
  
  if(runif(1) < likelihood.ratio) {
    samples[, i] = proposal
  }
  else {
    samples[, i] = samples[, i-1]
  }
}

Važna stvar kod SDT primera, na koju ranije nismo obratili pažnju, je da su parametri modela korelisani. Drugim rečima, relativna verovatnoća parametarskih vrednosti za \(d'\) će se razlikovati za različite vrednosti parametra \(C\). Korelisani modeli, u teoriji, ne deluju kao problem za MCMC, međutim, u praksi mogu uzorkovati velike poteškoće. Korelacije između parametara mogu dovesti do ekstremno spore konvergencije lanca uzorkovanja, nekad čak i do toga da nema konvergencije (bar u nekoj razumnoj količini vremena). Postoji više različitih pristupa uzorkovanju koji omogućavaju MCMC-u da se efikasno bavi takvim korelacijama.

Drugi tipovi uzorkovanja (pored osnovnog Metropolis-Hastings)

Metropolis-Hastings algoritam je veoma jednostavan i dovoljno moćan za mnoge probleme. Međutim, kada su parametri jako korelisani, bolje je koristiti kompleksniji pristup MCMC-u.

Gibbs-ovo uzorkovanje

Kad nam je data višedimenziona raspodela, kao što je primer SDT-a gore, često se koristi tip MCMC-a koji se naziva Gibbs-ovo uzorkovanje. Gibbs-ovo uzorkovanje razbija problem tako što izvlači uzorke za svaki parametar pojedinačno direktno iz uslovne raspodele tog parametra, tj. raspodele verovatnoće parametra, kada nam je drugi parametar poznat.

Ponekad se Gibbs-ovo uzorkovanje kombinuje sa Metropolis pristupom. Ovo ćemo ilustrovati koristeći SDT primer iz prethodnog dela. Ključ je u tome da se za višedimenzionu gustinu svaki parametar tretira odvojeno: predlog-prihvatanje/odbijanje koraci se preduzimaju parametar po parametar. Ovaj algoritam pokazuje kako se Metropolis i Gibbs mogu kombinvati za SDT primer:

  1. Biramo početne vrednosti za \(d'\) i \(C\), pretpostavimo da su ove vrednosti 1 i 0.5, redom.

  2. Generišemo novi predlog za \(d'\), analogno drugom koraku u Metropolis algortimu gore opisanom.

  3. Prihvatamo predlog ako je verovatnije da je došao iz populacione raspodele nego sadašnja vrednost \(d'\), za datu sadašnju vrednost \(C\). Dakle, za datu vrednost \(C = 0.5\), prihvatamo predlog, npr. \(d' = 1.2\), ako je to verovatnija vrednost od \(d'=1\) za ovu vrednost \(C\). Prihvatamo novi predlog sa verovatnoćom koja je jednaka odnosu verovatnoće predloženog \(d'\) (\(1.2\)) i sadašnjeg \(d'\) (\(1\)) za dato \(C = 0.5\). Pretpostavimo da je novi predlog (\(d' = 1.2\)) prihvaćen.

  4. Generišemo novi predlog za \(C\). Koristićemo istu raspodelu predloga kao za prvi parametar, u ovom slučaju \(\mathcal{N}(0,0.1)\) raspodelu. Pretpostavimo da je novi predlog za \(C\) \(0.6\).

  5. Prihvatamo novi predlog ako je verovatnije da je došao iz populacijske raspodele nego trenutna vrednost \(C\), za datu sadašnju vrednost \(d'\). Dakle, za datu vrednost \(d' = 1.2\), prihvatamo predlog \(C = 0.6\) ako je to verovatnija vrednost \(C\) nego \(0.5\) za fiksiranu vrednost \(d'\). Prihvatamo novu vrednost sa verovatnoćom koja je jednaka odnosu verovatnoće novog \(C\) \((0.6)\), i trenutnog \(C\) \((0.5)\), za datu vrednost \(d'\) od \(1.2\). Pretpostavimo u ovom slučaju da je predlog \(C = 0.6\) odbijen. Prema tome, ostaje \(C = 0.5\).

  6. Ovim se završava jedna iteracija uzorkovanja. Vraćamo se na korak 2 da bismo započeli novu iteraciju.

Sledi R kod za ovaj problem. Funkcija posterior.density je izostavljena, jer se definiše na isti način kao malopre. Ključna razlika između Metropolis uzorkovanja u prethodnom odeljku, i ovog Metropolis-a u sklopu Gibbs-a je da se predlog i procenjivanje odvijaju odvojeno za svaki parametar, a ne istovremeno za oba parametra. Petlja koja prolazi kroz broj parametara, for(j in rownames(samples)) dozvoljava da se parametaru \(d'\) daje nova predložena vrednost i da se odlučuje o prihvatanju, dok se parametar \(C\) drži na poslednjoj prihvaćenoj vrednosti i obrnuto.

# broj uzoraka
nmc = 1000

# broj parametara
n.pars = 2

# vektor u koji se smestaju uzorci
samples = array(dim = c(n.pars, nmc), dimnames = list(c("d'","C"), NULL))

# pocetne pretpostavke
samples[, 1] = c(1, 0.5)

# uzorkujemo
for(i in 2:nmc) {
  samples[, i] = samples[, i-1]
  
  for(j in rownames(samples)) {
    proposal = samples[, i]
    proposal[j] = proposal[j] + rnorm(1, 0, 0.1)
    
    new.likelihood = posterior.density(parameters = proposal, h = N.hits, fa = N.falsealarams, Np = N.present, Na = N.absent)
     old.likelihood = posterior.density(parameters = samples[, i], h = N.hits, fa = N.falsealarams, Np = N.present, Na = N.absent)
     likelihood.ratio = new.likelihood/old.likelihood
     
     if(runif(1) < likelihood.ratio) {
       samples[, i] = proposal
     }
  }
}

Dobijemo sledeće rezultate:

Leva kolona: Markovski lanac i i gustina uzorka za d'.

Srednja: Markovski lanac i gustina uzorka za C.

Desna kolona: Zajednički uzorci, jasno se vidi korelacija.

Leva i srednja kolona prikazuju promenljive \(d'\) i \(C\). Desna kolona pokazuje uzorke iz zajedničke aposteriorne raspodele, koja je dvodimenzionalna. Može se videti da su parametri korelisani. Ovakva korelacija je tipična kod parametara kognitivnih modela. Ovo može napraviti problem kod Metropolis-Hastings uzorkovanja, jer je korelisana ciljna raspodela vrlo slabo usklađena sa raspodelom predloga, koja ne uključuje nikakvu korelaciju među parametrima; uzorkovanje predloga iz nekorelisane zajedničke raspodele ignoriše činjenicu da se raspodela verovatnoće svakog parametra razlikuje u zavisnosti od vrednosti drugih parametara. Kombinacija Metropolis-a i Gibbs-a može ublažiti ovaj problem jer uklanja potrebu za razmatranjem višedimenzionih predloga, i umesto toga primenjuje korak prihvati/odbaci za svaki parametar zasebno.

Zaključak

Kroz jednostavne primere, upoznali smo se sa najjednostavnijim MCMC metodama uzorkovanja. Takođe, pomenuli smo i neke savete kako da se maksimalno iskoristi MCMC uzorkovanje (bez obzira na to koja se vrsta koristi), kao što je korišćenje više lanaca, procena burn-in-a i korišćenje parametara podešavanja. Opisali smo i neke scenarije u kojima je MCMC odličan alat za uzorkovanje iz zanimljivih raspodela. Primeri su se fokusirali na Bajesov zaključak, jer je MCMC moćan način da se izvede zaključak o takvim modelima i da se istraže aposteriorne raspodele parametara od interesa.

Dakle, razmotrili smo neke od detalja koji su svakako važni, ali postoji i mnogo drugih pitanja koja su ostala neistražena. Ovde smo želeli da uvidimo ideju MCMC-a kroz jednostavne metode. Drugi, napredniji, MCMC algoritmi kao što je Hamiltonian Monte Carlo zapravo rade veoma slično ovome, ali su samo mnogo pametniji u predlaganju sledećeg skoka.

Sve do sada opisano predstavlja samo uvod u izučavanju MCMC metoda uzorkovanja i njihovih primena. Postoji još mnogo MCMC metoda, a svaka metoda se razlikuje po svojoj složenosti i situacijama za koje je najprikladnija. Cilj nam je bio da se demistifikuje MCMC uzorkovanje i da se daju jednostavni primeri koji će ohrabriti da se usvoje MCMC metode u sopstvenim istraživanjima.