How can I get the MAP estimate via MCMC? After inference finishes, should I take all of my samples, run them through the likelihood function and choose the one that has the highest likelihood or Pyro has some functions that I should call for that?
for MAP see here. the usual way to find MAP estimates is optimization, not MCMC.
But is there a way to get an approximate MAP estimate using the samples that MCMC has produced?
as you noted you could take the MCMC samples you have and compute k = argmax_i log p(latent_i) + log p(data|latent_i). then
your approximate MAP estimate would be latent_k.