| 54 | | Éste último test en principio sólo nos sirve para comparar las medianas pero tiene la ventaja de |
| 55 | | ser muy rápido y de que ya existe en TOL como la función [wiki:TolGuiaDelUsuarioSet#FunciónAlgLib.MannWhitneyUtest AlgLib.MannWhitneyUtest]. |
| 56 | | Podría usarse como filtro inicial a cualquier otro test de comparación de muestras, pues si las |
| 57 | | medianas son significativamente distintas es evidente que las distribuciones no pueden ser las |
| 58 | | mismas. Pero además de eso puede usar de forma recursiva binaria del siguiente modo: |
| | 55 | La principal conclusión tras la lectura de estos documentos es que la comparación no paramétrica de distribuciones es un problema tremendamente complejo desde el punto de vista estadísticoy todavía más en lo que se refiere a su implementación. De hecho sólo he encontrado un método libre implementado en C++, el [http://www.google.es/url?sa=t&rct=j&q=a%20c%2B%2B%20program%20for%20the%20cram%C2%B4er-von%20mises&source=web&cd=2&ved=0CCwQFjAB&url=http%3A%2F%2Fwww.jstatsoft.org%2Fv17%2Fi08%2Fpaper&ei=0B2lTpCkDZGWOvvY2K4N&usg=AFQjCNFaDAewPd7AcY7N3I8dmtQZAPMZIQ&sig2=ySS0WwmbqcButNNDMTDByQ Cramér-von Mises], que es realmente bastante robusto pero tremendamente lento para cadenas de apenas unos cientos de muestras. |
| | 56 | El [http://www.alglib.net/hypothesistesting/mannwhitneyu.php Mann-Whitney U-test]) sí que es rapidísimo pero en principio sólo nos sirve para comparar las medianas. |
| | 57 | Tiene la ventaja añadida de que ya existe en TOL como la función [wiki:TolGuiaDelUsuarioSet#FunciónAlgLib.MannWhitneyUtest AlgLib.MannWhitneyUtest]. |
| | 58 | Podría usarse como filtro inicial a cualquier otro test de comparación de muestras, pues si las medianas son significativamente distintas es evidente que las distribuciones no pueden ser las mismas. |
| | 59 | Pero además de eso puede usar de una forma sencilla, aunque un tanto heurística, para comprobar que las muestras se parecen no sólo en los valores centrales sino en todo el rango de valores posibles |
| 60 | | 1. Se calcula el test Mann-Whitney U-test sobre las muestras actuales |
| 61 | | 1. Si se rechaza el test se devuelve falso. |
| 62 | | 1. Si se sobrepasa el nivel máximo de recursividad se devuelve cierto |
| 63 | | 1. Si el tamaño de alguna de las muestras es menor que el mínimo requerido se devuelve cierto |
| 64 | | 1. En cualquier otro caso se divide cada muestra en dos partes, la primera con los valores |
| 65 | | superiores a la mediana y la otra con los inferiores. |
| 66 | | 1. Se llama recursivamente al paso 1 con las muestras superiores a la mediana. |
| | 61 | 1. Se crea el vector de muestra conjunta mediante la simple concatenación de las muestras a chequear. |
| | 62 | 1. Se divide el rango de valores de la muestra conjunta en un número dado no, más bien pequeño, de partes equi-frecuentes. |
| | 63 | 1. En cada intervalo se calcula el test Mann-Whitney U-test sobre las muestras actuales censuradas en dicho intervalo. |
| | 64 | 1. Si en cualquiera de los intervalos se supera el nivel de significación se rechaza la hipótesis nula. |
| | 65 | |
| | 66 | Habría que estudiar cómo se comporta en casos reales y simulados a ver si es posible parametrizarlo para que resulte eficaz en el ámbito de nuestros objetivos. |
| | 67 | |
| | 68 | == Mecanismo de generación de puntos sobre-dispersos == |
| | 69 | |
| | 70 | Lo primero que hace falta para poder simular varias cadenas es generar unos cuantos puntos de origen de la cadena lo más dispersos que sea posible dentro de la región factible. El problema es que para poder generar unos pocos puntos dispersos sin favorecer demasiado la zona cercana al punto factible de origen es necesario generar una muestra más o menos uniforme o sobre-disersa que tienda a recorrer todo la región factible y eso no es fácil hacerlo generando sólo unos pocos puntos, más bien hay que generar un paseo aleatorio más o menos largo. |
| | 71 | |
| | 72 | En https://www.tol-project.org/ticket/1306#comment:14 se presenta el algoritmo actualmente usado que funciona estupendamente con pocas dimensiones, pero que hay que mejorar cunado aumentan las variables porque se atasca un poco y hay que generar muchos miles de muestras para recorrer medianamente la región factible. Tampoco es que tarde mucho en generar porque sólo tiene que evaluar si el punto es factible, aquí no hay densidades ni Gibbs ni nada de eso, algún producto matriz x vector y cosas así bastante lineales. Comparado con lo que cuesta simular luego la distribución a posteriori no es demasiado, aunque lógicamente sobrecarga la habrá. |
| | 73 | |
| | 74 | Se trata de un algoritmo completamente genérico pero implica la existencia de una función {{{g(x):R^n->R^m}}} tal que la región factible quede definida por las inecuaciones {{{g(x)<=0}}}. Tengo algunas ideas sobre cómo se puede mejorar el rendimiento |
| | 75 | |
| | 76 | == Integración de métodos de evaluación de la factibilidad en BSR == |
| | 77 | |
| | 78 | Para poder aplicar el algoritmo del punto anterior, o cualquier variante futura del mismo, es necesario introducir en BSR un método que me devuelva g(x) para un punto arbitrario x. Como BSR maneja un modelo dividido en bloques de Gibbs, es cada uno de los bloques el que es capaz de evaluar su parte correspondiente de {{{g(x)}}}, a la que podemos llamar {{{g[i](x[i])}}} para el bloque i-ésimo. Luego el master de BSR se encarga de dividir la entrada x en los bloques {{{x[i]}}} y de concatenar los {{{g[i]}}} de salida. |
| | 79 | |
| | 80 | Esta parte está hecha aunque hace falta mucho chequeo para darlo por bueno, tanto en robustez como en eficacia, especialmente en los bloques no lineales definidos por el usuario que es la parte más complicada y a lo peor es necesario que el filtro incorpore algún método extra para mayor eficacia. |
| | 81 | |
| | 82 | Hay que tener en cuenta que la generación de múltiples cadenas es incompatible con algunas características de BSR |
| | 83 | 1. La continuación de cadenas es evidentemente imposible si queremos forzar distintos puntos de inicio[[BR]] |
| | 84 | {{{ |
| | 85 | #!java |
| | 86 | Real config::do.resume == 0; |
| | 87 | }}} |
| | 88 | 1. La simulación parcial no sería en principio incompatible al 100% pero la implementación de los métodos de restricción |
| | 89 | {{{g(x)}}} serían demasiado complicados de implementar pues dependerían de que partes de {{{x}}} están fijadas y cuáles |
| | 90 | son libres. Por otra parte no es demasiado problema pues la simulación parcial se usa para reestimar modelos estables o |
| | 91 | generar previsiones, por lo que se da por supuesto que ya no hay problemas de convergencia. |
| | 92 | |
| | 93 | == Generación de múltiples cadenas en BSR == |
| | 94 | |
| | 95 | Hasta ahora los pasos para la simulación con BSR eran |
| | 96 | 1. Crear la definición del modelo como un {{{Struct @BSR.ModelDef modelDef}}} |
| | 97 | 1. Parametrizar la configuración como una instancia de {{{BysMcmc::@Config config}}} |
| | 98 | 1. Crear el gestor del modelo BSR como una instancia de la clase {{{BysMcmc::@Cycler cycler}}} |
| | 99 | 1. Crear el gestor de la estimación BSR como una instancia de la clase {{{BysMcmc::@Estim estim}}} |
| | 100 | 1. Generar la cadena MCMC llamando al método {{{estim::Run}}} |
| | 101 | Si queremos generar múltiples cadenas, los 4 primeros pasos siguen siendo los mismos, salvo que en el paso 2 hay que guardar una copia de config para recuperar la configuración tras el proceso de mezcla |
| | 102 | 1. Generar unos pocos puntos iniciales muy alejados entre sí |
| | 103 | 1. Resetear el modelo llamando al método {{{cycler::Initialize}}} |
| | 104 | 1. Generar un paseo aleatorio sobre-disperso en la región factible, es decir, con más probabilidad cerca de la frontera que en el interior. |
| | 105 | 1. Extraer unos pocos puntos iniciales lo más alejados entre sí que sea posible partiendo del centro de masas. |
| | 106 | 1. Cambiar los parámetros de configuración de las estimaciones locales para que |
| | 107 | 1. no haga ''burn-in'' ni ''thinning'' estándar[[BR]] |
| | 108 | {{{ |
| | 109 | #!java |
| | 110 | Real config::mcmc.burnin := 0; |
| | 111 | Real config::mcmc.thinning := 0; |
| | 112 | }}} |
| | 113 | 1. no haga la diagnosis ni los informes de resultados usuales[[BR]] |
| | 114 | {{{ |
| | 115 | #!java |
| | 116 | Real config::do.report := False; |
| | 117 | Real config::do.eval := False; |
| | 118 | Real config::do.linear.effects := False; |
| | 119 | }}} |
| | 120 | 1. Para cada punto inicial {{{x0[k]}}} generar de forma independiente una cadena MCMC de la distribución a posteriori |
| | 121 | 1. Resetear el modelo llamando al método cycler::Initialize |
| | 122 | 1. Asignar el punto inicial con la sentencia {{{Real cycler::_.sampler::setStore(x0[k]);}}} |
| | 123 | 1. Crear el k-esimo gestor de la estimación BSR como una instancia local de la clase BysMcmc::@Estim local_estim |
| | 124 | 1. Generar la k-esima cadena local MCMC llamando al método {{{local_estim::Run}}} |
| | 125 | 1. Devolver la cadena llamando a {{{ cycler::loadFullMcmc(?) }}} |
| | 126 | |
| | 127 | == Generación en paralelo con BSR == |
