egutierrez
47fac22230
- .claude/CLAUDE.md - .claude/commands/subagentes.md - .claude/rules/INDEX.md - .mcp.json - bash/functions/cybersecurity/analyze_dns.md - bash/functions/cybersecurity/audit_http_headers.md - bash/functions/cybersecurity/audit_ssh_config.md - bash/functions/cybersecurity/check_firewall.md - bash/functions/cybersecurity/detect_suspicious_users.md - bash/functions/cybersecurity/encrypt_file.md - ... Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
93 lines
3.8 KiB
Markdown
93 lines
3.8 KiB
Markdown
---
|
|
name: metropolis_hastings
|
|
kind: function
|
|
lang: cpp
|
|
domain: datascience
|
|
version: "1.0.0"
|
|
purity: pure
|
|
signature: "MHResult mh_run_1d(const std::function<double(double)>& log_pdf, double x0, double sigma, size_t n, double* out, Rng&); MHResult mh_run_nd(const std::function<double(const double*)>& log_pdf, const double* x0, const double* sigma, int d, size_t n, double* out, Rng&)"
|
|
description: "Metropolis-Hastings 1D y d-dimensional con proposal Gaussian symmetric. Target log-pdf inyectable via std::function (no necesita normalizarse). Devuelve cadena en out[] y accept rate. Pareja CPU del mc_metropolis_hastings_gpu."
|
|
tags: [mcmc, metropolis, hastings, sampling, bayesian, datascience, pendiente-usar]
|
|
uses_functions: ["rng_cpp_datascience"]
|
|
uses_types: []
|
|
returns: []
|
|
returns_optional: false
|
|
error_type: ""
|
|
imports: [cstddef, functional, cmath, vector]
|
|
tested: false
|
|
tests: []
|
|
test_file_path: ""
|
|
file_path: "cpp/functions/datascience/metropolis_hastings.cpp"
|
|
params:
|
|
- name: target_log_pdf
|
|
desc: "Functor con la log-densidad target (no normalizada). 1D recibe double; ND recibe const double* de tamano d."
|
|
- name: x0
|
|
desc: "Punto inicial (escalar 1D, array d-dim ND)."
|
|
- name: proposal_sigma
|
|
desc: "Stddev del proposal Gaussian. 1D escalar; ND array de d stddevs (uno por dimension)."
|
|
- name: d
|
|
desc: "(ND) dimension del espacio."
|
|
- name: n_samples
|
|
desc: "Total de samples a generar (incluyendo x0)."
|
|
- name: out_chain
|
|
desc: "Buffer destino. 1D: double[n]. ND: double[n*d] row-major (sample i, dim k = out[i*d+k])."
|
|
- name: r
|
|
desc: "Rng (de rng_cpp_datascience). State mutado durante el sampling."
|
|
output: "MHResult con n_samples, n_accepted y accept_rate. Tipico target: 0.20-0.40 accept rate (ajustar proposal_sigma)."
|
|
---
|
|
|
|
# metropolis_hastings
|
|
|
|
Random-walk Metropolis estandar. Replica el sampler que usan los 4 calculadores MCMC del set, en formato componible.
|
|
|
|
## Patron 1D (mcmc-bayes / mcmc-full Beta posterior)
|
|
|
|
```cpp
|
|
double k = 7.0, n_obs = 10.0; // observados
|
|
auto log_post = [k, n_obs](double theta) -> double {
|
|
if (theta <= 0.0 || theta >= 1.0) return -1e300;
|
|
// Beta(2, 2) prior + Binomial likelihood
|
|
return (k + 1.0) * std::log(theta) +
|
|
(n_obs - k + 1.0) * std::log(1.0 - theta);
|
|
};
|
|
|
|
fn::ds::Rng r;
|
|
fn::ds::rng_seed(r, 42);
|
|
|
|
std::vector<double> chain(10000);
|
|
auto res = fn::ds::mh_run_1d(log_post,
|
|
/*x0=*/0.5,
|
|
/*sigma=*/0.1,
|
|
chain.size(),
|
|
chain.data(),
|
|
r);
|
|
// Burn-in: descartar primeros 1000
|
|
// res.accept_rate ~ 0.3 si sigma esta bien ajustada
|
|
```
|
|
|
|
## Patron 2D (mcmc-visualizer)
|
|
|
|
```cpp
|
|
auto log_density = [](const double* xy) -> double {
|
|
double x = xy[0], y = xy[1];
|
|
// Bimodal target del visualizer
|
|
double d1 = (x-2)*(x-2) + (y-2)*(y-2);
|
|
double d2 = (x+2)*(x+2) + (y+2)*(y+2);
|
|
return std::log(std::exp(-0.5*d1) + std::exp(-0.5*d2));
|
|
};
|
|
|
|
double x0[2] = {0.0, 0.0};
|
|
double sigma[2] = {1.0, 1.0};
|
|
std::vector<double> chain(10000 * 2);
|
|
|
|
fn::ds::mh_run_nd(log_density, x0, sigma, 2,
|
|
10000, chain.data(), r);
|
|
```
|
|
|
|
## Notas
|
|
|
|
- **No hay tuning automatico** del proposal sigma. Si el accept rate es <0.1, baja sigma; si >0.6, sube. Para tunings serios usar adaptive MH (no incluido aqui).
|
|
- `target_log_pdf` se evalua 2 veces por sample en el peor caso (current cached, proposal nuevo). Si tu log-pdf es caro (10ms+ por eval), N=10^4 cuesta minutos — para esos casos preferir GPU MH (no caben dependencias caras en GLSL — tienes que portar la log-pdf).
|
|
- Compatible con `rhat_split` y `ess_basic`: corre M cadenas con seeds distintos, apila en `chains[m*n + i]`, llama los diagnosticos.
|
|
- 1D y ND comparten la maquinaria pero ND lleva el coste de los `std::vector<double>` por sample. Si rendimiento importa y conoces d en compile-time, fork-and-specialize.
|