 Hola a todos, bienvenidos al tutorial de Modelado Valleciano, en R con R Estanaro. Mi nombre es Fernando, Fernando Antonio C.P. Derrera. Yo soy estadístico, estudio en la Universidad de Warwick en Inglaterra. Mi carrera la estudié en actuaría en el ITAM en México y bueno, pues tengo algo de experiencia en investigación de mercados y en consultoría estadística en México, en dos empresas, una llamada numérica, la otra Think Data MX y ese es un poco mi perfil, también para que puedan ustedes conocerlo y pues cualquier duda que pudieran tener también contactarme y con toda confianza. Vamos a empezar. Básicamente el tutorial de hoy y que va a ser un poquito más teórico tiene tres grandes componentes. Primero una es una introducción al paradigma valleciano. Sé que algunos probablemente conozcan ya de estadística valleciana que hayan utilizado, probablemente ya de estadística valleciana en la práctica y dentro de R, pero no queremos asumir que ya lo conocen. Entonces intentaré un poco desde el inicio presentarles el paradigma valleciano. Vamos a hablar también un poquito de modelos lineales generalizados desde esa perspectiva valleciana y en particular cómo ajustar modelos lineales generalizados en R con el paquete de restanar y un poco también de cómo funciona el paquete y qué hay detrás de él, en términos de los métodos, porque evidentemente es importante o al menos desde mi punto de vista es siempre muy importante conocer los métodos estadísticos detrás de los paquetes que nos resuelven la vida o nos ayudan a resolver nuestros problemas del día a día. En cualquier momento, por favor, si tienen alguna pregunta, hay algún comentario. Siéntanse con toda la confianza de preguntarlo mediante el chat o haciendo la pregunta. Creo que somos un tutorial lo suficientemente pequeño para que se pueda dar un poco esa conversación. Entonces, con toda confianza me pueden preguntar o detener. Trataré de estar al pendiente del chat por cualquier cosa. Entonces, vamos a empezar y en primer lugar, pues, hablar un poco de este paradigma valleciano y a qué nos vamos a referir durante toda la sesión del día de hoy. Quiero empezar con esta pregunta que puede ser quizás muy abstracta o quizás como que muy trivial, al mismo tiempo, pero que es importante. Porque probabilidad es un término que tiene un significado matemático que quizás ese precisamente no es trivial o sencillo, pero también es cotidiano el uso del término probabilidad. En ese sentido, lo utilizamos muchas veces en nuestro día a día y valdría la pena preguntarnos cuál es la relación entre estos dos formas de pensar la probabilidad desde el punto de vista matemático y en nuestro día a día. En matemáticas, sin ser 100% formales, podemos decir que es una función que es no negativa, que suma uno y dependiendo de los escenarios podemos hablar de varias cosas relacionadas con esta función. Masa de probabilidad, función de densidad, de distribución, pero quiero que pensemos o que empecemos pensando en esa última vez que utilizamos el término probabilidad en el día a día, en nuestras relaciones cotidianas. Por ejemplo, podríamos preguntarnos cuál es la probabilidad de que esté lloviendo afuera de mi ventana. Esa es una pregunta quizás que no se podemos hacer todos, ¿no? Va llover hoy, no va llover hoy, qué tan probable es que si voy a hacer al mercado y me tengo que llevar el paraguas o no, porque va llover. Y no necesariamente estamos pensando en matemáticas, sino simplemente como que en nuestro lenguaje. También nos podemos preguntar la probabilidad de que un equipo gane, ¿no? Yo soy aficionado del necaxa, entonces pues quizás quiera yo o espere yo que mi equipo gane muchas veces, ¿no? Pero en general pues es usual hablar en términos de probabilidad cuando nos referimos a eventos deportivos, ¿no? Probabilidad de victoria, de empate, no es muy poco probable que tu equipo gane o es muy probable. Y no sé si algunos de ustedes tengan algunos otros ejemplos que con toda confianza si nos quieren compartir se les ocurren cómo utilizan en el día a día ustedes el término probabilidad. Si no, no se sientan obligados, pero pueden utilizar el chat o pueden comentarlo. Yo pienso por ejemplo que repito un poco en mi background de quizás en el trabajo, podríamos estar pensando en un producto, ¿no? Queremos lanzar al mercado y alguien puede decir o puede preguntarse sobre qué tan probable es que pegue en el mercado, ¿no? Que a los consumidores, que a nuestros clientes les guste que las personas lo compren, se les haga un producto atractivo o no, usamos muchas veces ese término. Y aquí es donde viene un poco el corazón, digamos, de la probabilidad vallesiana, perdón, la interpretación vallesiana de la probabilidad frente a las más tradicionales que conocemos en nuestros cursos de inferencia estadística o de probabilidad incluso, ¿no? Que son las interpretaciones clásicas o frecuentistas, así llamadas clásicas o frecuentistas, que se enfocan en cierto sentido ver la probabilidad como una medida de variabilidad y no tanto, perdón, y exclusivamente como variabilidad. Por ejemplo, si tiramos un dado, hemos escuchado, ¿no? La probabilidad de un sexto y la lógica desde un punto de vista clásico es que es un caso favorable, por ejemplo, de que salga seis, ¿no? Estoy hablando de la probabilidad de que salga seis, es un sexto, porque es una de las seis caras del dado. Y pensamos en que cada que lancemos un dado estamos hablando de un resultado variable y por eso hablamos de probabilidad. La interpretación frecuentista simplemente es que la probabilidad de un evento es un límite de frecuencias relativas, en el sentido de que si repetimos muchas veces un experimento, esos resultados van a ir variando y podemos calcular en el límite cuando hagamos muchos, muchos, muchos experimentos la proporción, la frecuencia relativa de las veces que ocurrió el evento frente a las que no ocurrió el evento. Sin embargo, hoy en cuando hablé de yo de interpretación vallesiana, los invito a pensar en la probabilidad como una medida de incertidumbre y en el sentido muy, muy, muy similar o prácticamente en el mismo sentido en el que lo usamos en ese día a día. Cuando usamos el término probabilidad en todas estas tres ejemplos que quería y en aquellos que ustedes puedan pensar, yo considero que estamos hablando, estamos utilizando el término de probabilidad como medida de incertidumbre. Y hay mucha literatura desde el punto de vista vallesiano que se enfoca en este concepto. La probabilidad es una medida de incertidumbre. Si yo estoy muy seguro de algo, voy a decir que ese evento tiene probabilidad prácticamente uno. Y usamos, estoy 80% seguro, estoy 70% seguro, 60% seguro para referirnos cada vez a eventos un punto más inciertos. Y este concepto de medida de incertidumbre como probabilidad, lo que nos dice es que entonces podemos asignarle probabilidad y en particular distribuciones de probabilidad a todas aquellas cosas o cantidades que desconozcamos en una situación da. Y ese es como el corazón de la estadística vallesiano. Ahora en la práctica, cuando queremos hacer inferencia, terminamos haciéndolo o un recurso importante es lo que Gutiérrez Peña llamaría la única receta o la inferencia, la receta de la inferencia vallesiano que es encontrar la distribución condicional de todas aquellas cantidades de interés cuyo valor desconocemos dado el valor conocido de las variables observadas. Y aquí voy a tratar de ir por partes digamos de esta receta y es de nuevo estamos pensando en hablar de probabilidad para aquellas cantidades que no conozcamos. Estas son nuestras cantidades de interés cuyo valor desconocemos. Pero cuando hacemos inferencia, cuando hacemos estadística o ciencia de datos, pues queremos datos y queremos entonces hacer estas estimaciones o inferencias condicionando en los valores de datos que obtuvimos. Queremos relacionar esos posibles cantidades desconocidas con datos que sí hayamos observado, por lo que podemos hablar en inferencia de un propósito de inversión y el propósito de un análisis estadístico sería digamos que recuperar las causas pensando en inferencia paramétrica estaríamos hablando de los parámetros que generan cierto conjunto de datos, un mecanismo generador de datos y a esos parámetros los podríamos identificar con un concepto de causas que son parámetros desconocidos, pero por eso queremos deducir de algo o inferir algo más bien sobre esos parámetros dadas unas observaciones, dados unos efectos y entonces tenemos que la inferencia queremos invertir la probabilidad, queremos basar la inferencia en la distribución de teta, es decir los parámetros que no los conocemos, condicional en x y esta distribución de teta condicional en x se conoce como la distribución posterior. Ahora bien, ¿cómo hacemos esta inversión en la práctica? Normalmente utilizamos lo que se conoce como el teorema de valles y este teorema de valles tiene cuatro componentes básicos, evidentemente de este valles es que viene el apodo o el sobrenombre de vallesiano, pero recordemos entonces que el concepto vallesiano es más bien referido a la interpretación de la probabilidad y el teorema simplemente es una herramienta que utilizamos mucho en la práctica. Los cuatro componentes básicos son la distribución posterior, que es nuestro objetivo y aquí la denoto como p de teta dado x en color azul y esa es recordando la receta lo que queremos, queremos la distribución de aquello que desconocemos el parámetro dados datos observados. Tenemos una distribución que le vamos a llamar inicial o a priori que llamo aquí p de teta y esa es la manera o la distribución en la que reflejamos la incertidumbre que tengamos sobre el parámetro, pero antes de observar los datos y por eso se llama inicial o a priori, porque es a priori de haber observado los datos. La verosimilitud que es la función de verosimilitud que conocemos de inferencia estadística y que podemos pensarla como la distribución condicional de los datos dado el parámetro y aquí es donde viene esta conexión de inversión, pero la vemos como una función del parámetro y con los datos fijos. Entonces si recuerdan o conocen de estimación de máxima verosimilitud, recordarán que fijamos o la vemos como una función del parámetro y buscamos maximizar el parámetro que, perdón, encontrar el parámetro que maximiza esa función para los datos observados en un contexto frecuentista y finalmente la distribución marginal de los datos y esta distribución es un poco como un promedio ponderado de este denominador de los numeradores, es decir, del producto entre la verosimilitud y la inicial. Es un promedio ponderado de manera informal entre estos dos componentes que tan probable es el dato que observamos dado un valor fijo del parámetro, pero ponderado por la probabilidad de que efectivamente ese y esimo valor posible del parámetro sea el verdadero valor. Como no lo conocemos, tenemos que promediar por la probabilidad de que un valor específico del parámetro sea el correcto. Entonces con esto quisiera, digamos, aterrizarlo un poco más en un primer ejemplo, poco podríamos un ejemplo de juguete, pero que en la vida real es un poco más complicado, pero que nos ayudará, espero, a entender un poquito mejor estos estos conceptos. Supongamos, por ejemplo, que queremos saber si un tratamiento médico es efectivo o no contra alguna enfermedad en general y sabemos que probablemente queremos tener n voluntarios en un ensayo clínico y vamos a estimar la efectividad que no conocemos, no, theta, no sabemos qué porcentaje de pacientes van a responder satisfactoriamente a ese tratamiento y eso es en lo que estamos interesados en estimar ese porcentaje. ¿Cómo vamos a representarlo? Si una persona responde al tratamiento, podemos usar este dibujito del corazón verde o numéricamente como un uno, pensar como un éxito, si no usamos el cero y aquí más adelante voy a utilizar este color más naran aranjado del corazón. ¿Cómo modelaríamos esta situación? Pues el modelo más sencillo es un modelo Bernoulli, por ejemplo, para cada paciente con una probabilidad de éxito común. Ese es, digamos, el pensar que cada aplicación del tratamiento es un polado, entre comillas, y la probabilidad de éxito es este parámetro theta, es esta efectividad del tratamiento y que no conocemos. En este sentido, x igual a 1 significaría respuesta satisfactoria, x igual a 0, que desafortunadamente no sirvió el tratamiento en esa persona. Ahora, esa es la función de verosimilitud. Si conocemos esto, esto es, aplica para la función de verosimilitud, pero ¿cuál sería la distribución inicial en este contexto vallesiano? Pues, pensemos primero que queremos modelar una probabilidad o un porcentaje de efectividad que está entre 0 y 1. Si no tenemos mayor información, antes de observar los datos, digamos que podríamos pensar en que cualquier posible porcentaje sea igualmente probable. Y una forma usual de reflejar esto es mediante el uso de una distribución uniforme. En este caso, la distribución uniforme o la densidad uniforme en el 01 asigna la misma densidad a todos los posibles valores de theta, con esta función de densidad constante igual a 1. Y lo veríamos así. Esto se conoce como distribución no informativa o frecuentemente van a encontrar el término de distribución no informativa. Y la idea es que busca reflejar ignorancia en el parámetro. No sabemos nada del parámetro. Entonces, más o menos has sentido que cualquier posible valor tenga la misma, entre comillas, probabilidad inicial de ser el verdadero valor. Digo, entre comillas, probabilidad, porque aquí estamos hablando de densidades. Pero en realidad hay que preguntarnos si esto en realidad está reflejando falta de información total. Por ejemplo, aquí resalto dos dos regiones, una de el punto 4 al punto 7 y otra del punto 7 al 1, en colores verde y naranja. Y resulta que estas dos regiones le estamos asignando la misma probabilidad de que el verdadero valor del parámetro esté entre punto 4 y punto 7. Es la misma probabilidad que esté entre punto 7 y 100 por ciento y uno. Eso hasta cierto punto es información que estamos diciendo. ¿Por qué? Porque, por ejemplo, podría llegar una persona que es experta en medicina, que conoce tanto del tratamiento como de la enfermedad, que sabe cómo se desarrolló este medicamento. Y entonces que tiene una noción de la posible efectividad de la evolución de los en la evolución de los pacientes. Y pues nos puede decir que quizás es muy poco probable que el tratamiento sea demasiado efectivo, no 80 por ciento o más efectivo y nos dice que pues es más probable que esté alrededor del del 45 por ciento. Si ese es el conocimiento, entonces como que no querríamos utilizar una distribución uniforme, porque esa es información. Y si estamos pensando en la probabilidad como una medida de incertidumbre, pues queremos incorporar la información que tenemos. Si no estaríamos utilizando, no estaríamos, digamos, reflejando la incertidumbre correctamente. Una distribución que se puede utilizar es la distribución beta. Y en particular, por ejemplo, podemos utilizar esta distribución beta con dos parámetros, alfa igual a 2.5, beta igual a 3. Y la distribución beta es una distribución muy flexible entre 0 y 1. No, no me quiero detener demasiado en esa distribución en particular. Pero aquí resalto, por ejemplo, el corte en la línea roja divide, digamos, la región menor al 80 por ciento de efectividad que esta persona experta nos nos dijo que no creía que el verdadero, la verdadera efectividad del tratamiento fuera a superar esa efectividad. Y pongo en la línea punteada como referencia esa distribución uniforme para que se vea este contraste entre cómo estamos asignándole mayor probabilidad, mayor densidad y mayor masa a las regiones centrales alrededor del prácticamente el 45 por ciento que este experto nos estaba mencionando. Esta persona nos comentó. Ahora, esa es la distribución inicial. Y en particular, es algo que es una distribución que se conoce, así como la uniforme, hablábamos de no informativa. Este es un ejemplo de algo que se llama distribución conjugada, porque la distribución posterior cuando aplicamos el sistema de valles nos va a devolver otra distribución beta pertenece a la misma familia que la inicial. Entonces tenemos un inicial beta de esta forma y la verosimilitud Bernoulli son las dos expresiones que ya había puesto en las diapositivas anteriores. No me voy a meter, digamos, o detener mucho en en la marginal, pero la posterior, si ustedes arrastraran el lápiz y quisieran hacer todos los cálculos, podrían ver que al hacer o aplicar el trim de valles, esta sería la forma de la distribución posterior. Y resalto en verde y en naranja un poco algunos de los términos para que veamos de dónde vienen. Entonces, por ejemplo, el primer factor tenemos teta elevado a la X y si se fijan esa teta a la X viene de la verosimilitud en verde y teta a la alfa menos uno viene de la inicial y lo mismo para el uno menos uno menos teta. Entonces, esta distribución conjugada es muy útil, es muy pragmática, porque si se fijan la X o el uno menos X nos dice si hubo un éxito o hubo un fracaso y entonces esta distribución conjugada es una forma rápida de ir actualizando nuestras creencias, porque si hubo un éxito le vamos a sumar el uno de la X al parámetro alfa de la distribución beta inicial y si hubo un fracaso le vamos a sumar el uno menos cero al parámetro beta de la inicial. Por ejemplo, pensemos que el primer paciente tiene un éxito. Entonces, empezamos con nuestra distribución inicial en azul oscuro, es la beta 2.53, observamos un corazoncito verde y actualizamos nuestra distribución posterior. Le sumamos un éxito al primer parámetro de la beta y cero fracasos al segundo parámetro de la beta y vemos cómo la distribución se mueve y si se fijan se desplaza un poco a la derecha, es decir, aumenta la probabilidad posterior de efectividades mayores. Luego tenemos un segundo paciente y podemos repetir el ciclo, porque ahora la distribución posterior es el conocimiento que tenemos cuando aplicamos un segundo tratamiento, esa es nuestra distribución inicial y resulta que vuelve a ser efectivo el tratamiento y volvemos a actualizar la distribución posterior y vemos que se vuelve a desplazar a la derecha de la curva verde y entonces empezamos a dar más probabilidad a valores más altos de la efectividad. Pero el tercer paciente no responde al tratamiento y entonces actualizamos las creencias aumentándole uno digamos al contador a este segundo parámetro. Así podríamos seguir con el cuarto paciente, el quinto paciente y supongamos por ejemplo que hasta el décimo paciente este es el escenario que tenemos y si se fijan conforme vamos aumentando el número de pacientes las distribuciones empiezan a hacerse un poquito más pícudas y eso quiere decir que estamos cada vez más seguros o estamos estimando con mayor certeza el verdadero valor del parámetro teta de la efectividad. Podríamos seguir el ensayo clínico, llegar al paciente 20 y de nuevo están cada vez más estrechas las distribuciones y si se fijan ya el cada nuevo dato ya no se mueve tanto la distribución conforme más datos tenemos pues estamos ya más seguros del conocimiento que hemos ido adquiriendo. En total por ejemplo podríamos pensar que tuvimos 30 pacientes y que los resultados fueron los siguientes no hay primero los dos éxitos del inicio luego tuvimos tres pacientes que no respondieron la sexta persona sí respondió la séptima no así sucesivamente en total tuvimos 20 personas que respondieron al tratamiento y 10 que no respondieron al tratamiento. Al final tendríamos la distribución posterior beta 22.5 y 13 y notemos que aquí podríamos haber hecho si lo tenemos los 30 datos en un solo jalón podríamos hacer esta actualización en un solo paso pero también la podemos ir haciendo paso a paso. Esta es una de las ventajas o de las cosas atractivas del paradigma vallesiano en el sentido de que podemos ir actualizando nuestras creencias conforme van vamos recopilando información. Ahora esto perdón fue fue una introducción más menos general y un ejemplo específico de la inferencia de Diana. Alguien tiene alguna duda quisiera hacer una pregunta. Yo mismo. Sí, claro. Pongo la cámara. La probabilidad del denominador la que dice es que es la de los datos. No me queda del todo claro lo que represento. Ok. Estoy regresando un poquito esta marginal verdad. Exacto. Ok. Digamos que pensemos que el en este ejemplo médico la X es si un paciente específico respondió al tratamiento o no. Entonces la pdx nos referimos a la probabilidad de que un paciente dado si tomáramos a alguien al azar. Voy a permitir por ejemplo que este aplicáramos. Si la primera persona que llega a la clínica le aplicamos el tratamiento e queremos ver la probabilidad de que responda esa persona al tratamiento. Y en este modelo estamos pensando que esa probabilidad depende de un parámetro que le llamamos efectividad del tratamiento que le llamamos theta. Ese es un número fijo, pero que desconocemos. Si supiéramos el verdadero valor de theta utilizaríamos la función verde, la verosimilitud. Porque ya tendríamos la probabilidad de respuesta, satisfactoria, el tratamiento X igual a uno, pero dado un valor específico de theta. El problema es que no conocemos el verdadero valor de theta. Entonces esta distribución roja, digamos, lo que hace es ponderar o promediar todos los posibles valores de theta. Pensemos que hubieran sólo K posibles valores, 3 que fuera 0.5 o 1. Entonces nunca funciona el tratamiento, no sirve de nada el tratamiento. 0.5 de efectividad es un volado el tratamiento o es un tratamiento extremadamente efectivo y siempre cura a la enfermedad. Pues si fuera 0 tendríamos pdx dado theta igual a theta 1 igual a 0. Si fuera 0.5 sería pdx dado theta igual a 0.5. Y si fuera 1 igual, pero no sabemos cuál de las tres es. Vamos a ponderar cada una de esas tres probabilidades por la probabilidad que le estamos asignando a theta. Si nuestro experto médico nos dice, aquí me voy a mover de nuevo perdonar la distribución inicial. Nuestro experto médico nos dice que es poco probable que rebase el 80 por ciento. Eso quiere decir que a las verosimilitudes con theta mayor a 0.8 las va a disminuir en esa suma y va a considerar mayor. Le va a dar mayor peso a las verosimilitudes en el centro de esta distribución. No sé perdón si sí más o menos. Al final la cosa también es que la verosimilitud no se no suma uno cuando la integras y supongo que el denominador esta te ayuda. Exacto. Efectivamente cuando obviamente cuando haces la la un solo producto en el de arriba eso no necesariamente te integra uno. Entonces lo que queremos hacer es normalizar para que la posterior efectivamente sea una una probabilidad y cumpla estos réxitos matemáticos. Pero digamos que el. De hecho se conoce muchas veces como constante de normalización esta distribución marginal y le llamamos constante porque no depende de theta. No una vez que hacemos la suma o que hacemos la integral estamos eliminando la dependencia de theta. Entonces si efectivamente es una forma digamos de hacer que todo integre o sume a uno como bien dices. Perfecto gracias. A ti por la pregunta y espero que también los demás la hayan encontrado de provecho. Alguien más tiene tiene algún comentario duda. Hasta ahora. Si yo yo quisiera preguntar acerca del pues el ejemplo que diste que lo presentaste era con la pinomial. Cierto. O es ser o uno y algo en el medio. Y que pasa con con distribuciones como normal o otros tipos. Hay alguna suposición. Si por ejemplo más adelante podemos en en otro contexto ya de modelos generalizados hablar de otro tipo de distribuciones. Pero pero la idea sería sería más o menos la la misma. Por ejemplo tendríamos que. Plantear otro otro este. Otro otro experimento en el que el número sea continuo aquí como nuestra observación es un éxito o un fracaso. Naturalmente pensamos como como dices en una binomial o en una Bernoulli. Pero pero también podríamos hablar de distribuciones normales. Y habrán distribuciones conjugadas para cada uno de los casos. De hecho la distribución normal es una distribución conjugada de una verosimilitud normal. Entonces básicamente si observamos datos que vienen de una distribución normal cada vez vamos a estar actualizando el promedio. No la la me el promedio es el primer parámetro de la de la distribución normal. Y vamos a ir disminuyendo nuestra incertidumbre como como se ve aquí en esta en esta gráfica. No al inicio teníamos una distribución muy ancha o relativamente no muy no pero relativamente ancha. Y cuando tuvimos 30 puntos la actualizamos y la la estrechamos. Este fenómeno se podría dar también para una distribución normal. Donde los la efectividad por ejemplo ya fuera un número real. Este digamos ya no no es este con este ejemplo en particular pero pero sí sí que podría hacerse este digamos el equivalente. Gracias. De nada. Bueno entonces siguiendo un poco en esta línea y vamos a hablar de modelos lineales generalizados y de regresión que aquí ya por ejemplo es algo más parecido a este este pregunta que hace a Samuel de una distribución normal por ejemplo. En en regresión lineal estamos pensando en describir una variable de interés ya a partir de otras variables explicativas X y normalmente pueden observarse distintos valores de Y para un mismo valor de X. Entonces no no estamos planteando un modelo de igualdad de que es igual una función determinística de de X. Pero si pensamos o en empíricamente podemos ver algunas asociación por ejemplo que cuando aumenta X esperamos que aumente Y también no. Y tiende a aumentar y hay muchas formas digamos de expresar este tipo de asociaciones. Es una forma de ver las regresiones es que están pensar que estamos modelando el valor esperado y aquí uso digamos con la negocia y ventaja el término esperar de la esperanza matemática. El valor esperado de nuestra variable de interés mediante alguna función de las variables explicativas y la función entonces perdón y aquí vamos a estar hablando entonces de incertidumbre de nuevo cuenta porque no tenemos esta función determinística no de que si supiéramos X hay un solo valor posible de Y no tenemos incertidumbre y esa incertidumbre desde la perspectiva vallesiana la reflejamos como una distribución de probabilidad condicional de Y pero dado lo que conocemos y que conocemos pues los datos y supongamos que conociéramos la relación específica que vamos a llamar o indexar con unos parámetros teta. Si hablamos de esta interpretación de modelar el valor esperado entonces tendríamos que la forma en la que X afecta a Y o está asociada a Y no no no hablo aquí de de causalidad no este sino simplemente de asociación. Es mediante el valor esperado entonces tenemos un modelo de probabilidad v de Y dado X y teta pero lo que sí sabemos es que la esperanza de Y dada X y teta es una función es así determinística de un subconjunto de parámetros teta D y de los datos de las variables explicativas X y la forma más sencilla es de esta función determinística es la una función lineal en los coeficientes y esta es la forma tradicional digamos de la regresión en la que el valor esperado de Y dado X es beta 0 más beta 1 por X 1 más beta 2 por X 2 así hasta beta de menos 1 por X de menos 1 no y en total tenemos determinó de coeficientes incluyendo el intercepto beta 0 podríamos tener otros parámetros que aquí le llamé teta R que no determinan la esperanza pero que también influyen en en G por ejemplo la regresión lineal tradicional una de las formas más usuales de expresarla es G y cada para cada y es igual a beta 0 más beta 1 y al final le sumamos un error normal así es como muchas veces se presenta o se explica la distribución la regresión lineal tradicional pues esa distribución la normal nos permite también reexpresar este mismo modelo de la de la siguiente manera no G y dado X y beta y sigma es una normal con una media muy y una varianza sigma cuadrada y esa media muy es de la misma forma que en la de positiva anterior no la ed y dado X y los parámetros entonces tenemos que beta juegue el papel de teta D y sigma cuadrada juegue el papel de teta R no hay un subconjunto de parámetros que determinan el valor esperado de la variable de interés y hay otros parámetros que podemos conocer o desconocer no y recordarán que hay veces que hacemos regresiones asumiendo un valor de teta perdón de sigma cuadrada y hay situaciones en las que también desconocemos el valor de sigma cuadrada ahora en el contexto vallesiano además tendríamos que asignar la distribución inicial de los parámetros desconocidos sean estos todos beta y sigma cuadrada o simplemente beta esos también serían el siguiente paso que tendríamos que hacer desde el punto de vista vallesiano pero los modelos lineales generalizados van un poquito más allá de la regresión lineal sencilla o múltiple y permiten modelar otro tipo de datos o situaciones más complejas por ejemplo datos de conteos de proporciones porque recordemos que la distribución la regresión lineal normal es para valores continuos no ya tiene que ser un valor continuo y hay veces que nuestra ejemplo o nuestro modelo exige otro tipo de de modelos estadísticos que no sean la distribución normal pero todos estos modelos generalizados siguen la misma noción de modelar el valor esperado de la variable de interés como una función de las variables explicativas en particular son tres elementos y perdón se me juntan quería tenerlos un poco más separados los los bullets para que no se enzimar el texto pero ya va a ser un miembro de algo que se conoce como familia exponencial no me quiero detener digamos en cuál es la forma específica de la familia exponencial pero esto incluye a la distribución normal a la poason a la Bernoulli a la beta y incluso otras distribuciones bajo algunas condiciones como la binomial la multinomial la binomial negativa entonces la familia de modelos generalizados es efectivamente amplia pero conceptualmente la idea es vamos a elegir una familia de distribuciones luego vamos a tener un predictor lineal esta función más sencilla que podíamos pensar de x beta 0 más beta 1 x 1 más beta de menos 1 x de menos 1 y que es la misma expresión digamos que en el modelo de lineal normal pero ahora vamos a tener una función liga que vaya que va a vincular el valor esperado de la variable de interés con el predictor lineal y esta función liga puede tomar muchas formas en el caso de la regresión normal la función liga es la identidad entonces no le hacemos nada pero por ejemplo si quisiéramos usar una regresión poason la función liga más utilizada es el logaritmo y entonces se habla de regresión log poason en ese caso tendríamos que la g está asociada al logaritmo o la regresión logística no sé si algunos de ustedes han escuchado o han tenido contacto con la regresión logística es un modelo lineal generalizado en el que la g es bernoulli tenemos un predictor lineal y ese predictor lineal está relacionado con yé mediante la función liga logística y de ahí viene el nombre de regresión logística y decimos que el valor esperado de yé la por la probabilidad o la proporción es el logit no de un predictor lineal que depende de las variables explicativas ahora el paquete r están arm que es lo que lo que nos permite pues lo que nos va a permitir es aplicar o modelar de manera más o menos automatizada y con vamos a decir valores de default sensatos o sensibles este tipo de de modelos una regresión lineal normal pero también otro tipo de regresiones como la poason la bernoulli e incluso aunque hoy no hablaremos de eso pero regresión beta y como como cuál es la función digamos clave del paquete r están arm bueno pensemos por ejemplo en en este ejemplo este este gráfico tenemos una variable que es log disk estos datos por cierto vienen de del conjunto de datos mt cars que seguramente ya conocen algunos de ustedes porque es muy utilizado no en el contexto de r y tenemos la variable de disk y la variable de mpg nunca recuerdo que que es disk creo que mpg es miles per gallon no de los de los carros pero pero bueno digamos que tenemos esta esta relación y este podría ser un ejemplo de de una regresión que quisiéramos hacer un ejemplo de juguete pero que nos ayuda digamos a ver el código específico de cómo se se ajustan estos modelos en r el enfoque clásico tradicional es utilizar la función glem del paquete stats creo pero es de r base y es muy utilizada esta esta función y es muy sencilla entre comillas en el sentido de que simplemente le pasamos una fórmula y un conjunto de datos esos son los dos parámetros de la función glem la fórmula está pensada para que trate de reflejar este predictor lineal del que hablaba o la relación entre y x entonces utilizamos el latilde y del lado izquierdo tenemos la respuesta entonces mpg va a ser nuestra variable de interés la y y del lado derecho de la tilde las variables explicativas en este caso solamente tenemos una variable explicativa entonces decimos que y tilde x y en este caso mpg tilde log disk y son los nombres de las variables y le tenemos que pasar un segundo parámetro que son los datos donde provienen estas estas variables nuestras observaciones en este caso tendríamos una matriz data frame o un tibel dependiendo de lo que utilicen ustedes en dentro de r pero los nombres los de las columnas tendrían que ser estos nombres de la fórmula que vayamos a utilizar y por ejemplo si ustedes correrían correrían esto en en en el estudio este tendría que ser el resultado bueno menos lo que tenemos donde tenemos estimadores de máxima velocidad similitud este es un enfoque frecuentista clásico y es muy utilizado no en mer pero si quisiéramos hacer lo mismo desde el punto de vista valles ya no podemos utilizar el paquete r sanar y este paquete creo que está teniendo mucho éxito puede ser muy útil para muchas personas porque uno de sus aciertos es que busca prácticamente ser un paso directo de esta función glem o de este tipo de funciones de modelos lineales generalizados que se han utilizado en el repórmulo tiempo la función en r sanar perdón es la siguiente y si se fijan lo único que cambia es además del del paquete es que le agregamos un tan guión bajo antes del glem pero la sintaxis es la misma hasta ahorita no la fórmula y los datos eso no cambio sólo utilizamos una función distinta que se llama tan guión bajo glem ahora que nos dice están glem una vez que termina no de correr nos imprime un poco de información distinta primero nos dice la familia y ese si se fijan es el primer elemento que les mencionaba de los modelos lineales generalizados dice que la familia es gaussian si estamos hablando de una regresión lineal normal más aún entre corchetes o nos dice que hay identity esa es la función liga que es la identidad decir no no estamos haciendo otra cosa que no sea una regresión lineal simple nos da la fórmula cuántas observaciones y cuántos predictores hubo aquí predictores es de tenemos 2 porque tenemos un intercepto y la variable si se fijan el intercepto lo pone de default no tuvimos que ponerlo en la fórmula y ese es el default en la función de glem y por eso lo hereda es restanar en no tenemos que especificar forzosamente el intercepto que también se podría poner como uno más en la fórmula pero no es necesario y si se fijan los resultados son prácticamente los mismos e los coeficientes de glem nos decía que el intercepto era 69.2 más o menos y de la variable explicativa log disk menos 9.3 si redondeamos que son el 69 y el menos 9.3 los numericamente tenemos algo más o menos equivalente pero también nos da los parámetros auxiliares los que yo llamaba teta r en este caso sigma no estoy asumiendo que que se conoce sigma entonces también lo estima y nos da un estimado de 2.6 y la diferencia es que aquí está dándonos y se fijan la mediana porque porque tenemos este gráfico y esta es creo una de las ventajas o de los motivos por los que queríamos utilizar el paradigma vallesiano y queríamos hacer regresiones desde el punto de vista vallesiano y es que tenemos o reflejamos la incertidumbre sobre el parámetro en los enfoques tradicionales de máxima verosimilitud solo queremos dar un estimador puntual no el estimador de máxima verosimilitud y eso es bueno y eso es muy útil pero el añadido del paradigma vallesiano es que también contamos con lo que decía al inicio que era la receta con la distribución completa es decir reflejamos la incertidumbre que tenemos sobre ese valor estimado en su máxima expresión entonces aquí por ejemplo estoy graficando las dos distribuciones marginales posteriores es decir la distribución del betas 0 el intercepto y la distribución de beta 1 dados los datos que se observaron y en el centro está la línea vertical que refleja esta mediana para el intercepto 69 para beta 1 menos 9 punto 3 y sombreado no sé si se alcanzan a ver en en azul está el intervalo central que sería un poco el equivalente a la noción de intervalo de confianza pero aquí también hay digamos una distinción que a las personas que utilizan estadística vallesiana ponen mucho énfasis en ello y es que estos se conocen como intervalos de probabilidad y que efectivamente corresponden a la noción que tenemos muchas veces de manera intuitiva de lo que queremos no sé si les ha pasado que cuando quieren hacer un análisis estadístico de datos como que pensamos cuando decimos bueno y dame un intervalo de estimación voy a utilizar el término como para no usar ni intervalo de confianza ni intervalo de probabilidad pero dame un intervalo como que intuitivamente pensamos en cuál es la probabilidad o queremos un intervalo que tenga mucha probabilidad de que el verdadero valor esté dentro de ese intervalo y esa interpretación directa de la probabilidad de que el parámetro esté entre a y b sólo se puede dar digamos de manera coherente dentro de este paradigma vallesiano porque aquí sí tenemos toda esta distribución y entonces podemos dar este tipo de enunciados probabilísticos del parámetro la probabilidad de y estos por cierto son al 95 por ciento entonces quiere decir que hay un 95 por ciento de probabilidad que el verdadero valor del intercepto esté en esta región sombreada y la mediana posterior es decir ese es un estimador puntual es ese menos 69 podríamos haber utilizado otro estimador puntual el valor esperado o la moda por ejemplo pero por default el restanar nos da o nos imprime directamente la mediana y nos da de hecho este math sd que es min absolute deviation perdón median absolute deviation standard deviation y es prácticamente digamos un estimación del error digamos que tenemos en este o de la incertidumbre que tenemos sobre ese estimador puntual de la mediana hasta aquí perdón más o menos todo claro no sé si algunos de ustedes ya habían utilizado este esta sintaxis de fórmula o el o el la función glem o si es nuevo para todos ustedes veo ahí este pulgares arriba algunos ok ahora pasé mucho tiempo o bueno pasé tiempo no sé si mucho poco pero al hablando de la distribución inicial al inicio y aquí en la función de restanar parece que dónde quedó la distribución inicial no qué está pasando bueno cómo podemos determinar esa distribución en una regresión y si puedo hacer una pregunta si si claro yo acabo de hacer lo mismo que que presentaste en el slide en mi en mi sesión y aquí y lo que ver es que hay además de los dos estimados intercepto y el log disc me parece también un sigma aquí es qué significa el sigma si este de aquí pero no sé si en su ustedes que están viendo mi pantalla también se ve cuando con el cursor este resalto estoy tratando de resaltar justo ese sigma que mencionas en mi pantalla si si se ve en donde lo estoy resaltando ese sigma es no voy a regresar es esta sigma de de la regresión lineal digamos la el valor de la varianza de los errores por ejemplo si es si están ustedes acostumbrados a pensar en la regresión lineal como este y igual a la suba este predictor lineal y un error el la varianza de ese error ese es el valor de sigma y le llaman parámetro auxiliar porque no pertenece a ese predictor lineal que determina el valor esperado en este concepto de modelos lineales generalizados pero es este esta sigma de estas de estas formulas ok gracias y perdón el mi por tuñol no pero todas las preguntas que tengan este con mucho gusto y si puedo responderlas encantado para para eso también es es el tutorial y ok entonces realizamos un poquito a la distribución inicial este si se fijan y samuel había hecho haber hecho la pregunta ya estamos hablando de distribuciones normales podemos pensar en en hacer un poco lo mismo que en el ejemplo médico y es a ver qué sabemos sobre los parámetros sobre beta sobre sigma antes de observar los datos no supongamos que no conocemos esta gráfica y sólo sólo sabemos que vamos a hacer una regresión lineal de dos variables no pues sabemos que los coeficientes son números reales y que sigma no puede ser negativa no porque es una varianza pero esto si quisiéramos hacer lo mismo que hicimos en el ejemplo médico y asignar una distribución uniforme como esa noción de no informativa nos va a dar un pequeño problema relacionado con lo que yo sé creo fue quien quien mencionaba ya no va a sumar a uno que es una distribución impropia no podemos asignar una distribución uniforme en todos los reales o no podemos asignar una distribución uniforme en todos los números no negativos pongo problema entre comillas porque existen digamos justificaciones teóricas y casos en los que se podrían utilizar y mucha gente lo utiliza pero no es la única forma de asignar distribuciones no informativas que pueden encontrar en la práctica y en particular muchas distribuciones inusual y iniciales usadas son asignar una distribución normal a los coeficientes por ejemplo o distribuciones te de estudent si queremos reflejar mayor incertidumbre normalmente porque la te de estudent permite digamos mayores errores tiene colas más pesadas para los coeficientes y distribuciones digamos en los números no negativos son por ejemplo la gama o el caso particular de una exponencial y muchas veces esas son las distribuciones que se utilizan para la regresión normal de hecho son también conjugadas es decir existe fórmulas como vimos en el caso de la beta donde nada más había que irle sumando a los parámetros los resultados de éxitos o fracasos para el caso de la regresión lineal normal existen fórmulas si las distribuciones iniciales son normal es gama normal para los coeficientes y gama para la sigma y hay una fórmula pero digamos que la noción de utilizar luego una distribución no uniforme perdón si no uniforme pero que se parezca lo más a la uniforme para representar esta ignorancia y se ha utilizado muchas veces la idea de buscar distribuciones planas y en inglés van a encontrarlas como flat non informative flat distributions por ejemplo si quisiéramos asignar de distribuciones normales la forma en la que la gente trata de aplanar la distribución normal es aumentando el valor de la desviación estándar de esa normal entonces por ejemplo aquí les grafico 4 de esas distribuciones normales y cada una con una desviación estándar mayor una normal 0 sigma 2.5 luego una 0.5 una 0.10 y una 0.50 y conforme aumentamos el parámetro de la desviación estándar de la normal vemos que vamos de una distribución muy estrecha a una muy plana sin embargo esto es en realidad un poco un espejismo de la escala no si se fijan aquí está y graficando de menos 100 a 100 que pasa por ejemplo si graficamos no sé de menos 2500 a 2500 pues las distribuciones de desviaciones estándar pequeñas se siguen viendo estrechas de hecho se ven mucho más estrechas pero la que parecía la curva roja o la naranja que parecían más planas bueno sotó la roja en este primer gráfico empezamos a ver esta misma forma de campana en realidad seguimos teniendo una distribución normal que sigue teniendo la forma tradicional de la campana pero más preocupante entre comillas es que si tienen como algunas implicaciones extrañas recordemos la distribución inicial nos da el conocimiento que tenemos sobre el parámetro antes de observar los datos supongamos que alguien nos dice voy a utilizar una beta normal 0 sigma 100 porque sigma 100 pues es una desviación estándar muy muy amplia entonces estoy haciendo una distribución normal súper súper variable súper plana en esta noción porque quiero hacerla muy plana pues en realidad por ejemplo si yo tengo este modelo resulta que la probabilidad de que beta está entre menos 25 y 25 es de alrededor de punto 97 eso equivale que a priori estamos diciendo que hay más de 80 por ciento de probabilidad de que el efecto entiendo efecto aquí como el coeficiente beta de esa variable a la que le estamos asignando esa distribución inicial es mayor a 25 digamos si si la probabilidad de que beta sea en valor absoluto menor a 25 es menor a 20 por ciento quiere decir que estoy asignando más de 80 por ciento de probabilidad para para efectos de magnitud mayor a 25 y muchas veces sobre todo en procesos de ciencia de datos o de estadística donde hacemos estandarización de variables variables entradas vemos todo esto de estas pasos que muchas veces se siguen no en los pipelines de análisis de datos encontrarse con un efecto de magnitud 25 ya es muy extremo a veces obviamente dependerá de las unidades pero por eso hablo después de estandarización muchas veces si inmediatamente ponemos distribuciones con valores de sigma enormes iniciales queriendo hacer esto de una distribución plana tenemos que ser conscientes de que estamos diciendo que a priori creemos que hay un efecto enorme es decir la varianza a priori muy grande no significa menos certeza sino más bien que creemos probables efectos de mayor magnitud entonces qué es porque les comentaba que r estanar busca tener defaults sensatos o digamos generalizables lo más seguros posibles para las distintas aplicaciones que se pueden utilizar cuáles son las distribuciones iniciales que utilizo el restanar en particular este es el mismo ejemplo y si utilizo la función priori de r estanar nos imprime esta información y nos dice que hay priors es decir las distribuciones a priori iniciales para el intercepto para los coeficientes intercepto coeficientes y para el bar para metro auxiliar de sigma de los interceptos nos dice que son después de que los predictores sean centrados entonces esto es algo frecuente también que puede puedes asignar distribuciones iniciales una vez que ya se centraron los los predictores nos dice que la distribución especificada es una normal con parámetro de locación 20 y parámetro de escala o desviación estándar 2.5 y luego nos habla de una distribución ajustar y lo que hace de restanar es ajusta esta distribución que si se fijan era la más estrecha que yo había puesto acá de la de 2.5 si una distribución que refleja que estamos pensando en que hay efectos muy pequeños decir que casi no hay efectos pero la escala la ajusta por la magnitud o el orden de magnitud de los datos específicos que le que le estamos pasando entonces toma las desviaciones estándares o los órdenes de magnitud de los datos específicos que le pasemos a la función y reescala esta distribución de default de 2.5 al valor específico que se requiera no lo multiplica por en particular la desviación estándar de y entre la desviación estándar de x y eso lo hace tras bambalinas y ese es como el default del parámetro y la locación toma la media empírica si ustedes le presentan esto a un vallesiano muy purista muy purista hay algunos que no les gusta esto porque estás en teoría una distribución inicial tendría que ser antes de observar los datos y aquí como que estás utilizando los datos para fijar la distribución inicial pero digamos que aquí la necesidad es dar unos valores como que de default que funcionen y lo que está buscando la distribución inicial aquí es solo identificar la escala los órdenes de magnitud precisamente para que no estemos creyendo que tenemos una distribución súper plana cuando en realidad no es plana sino que estamos asignando probabilidades poco poco sensatas y los coeficientes siempre lo centra en en cero de default es decir a priori no asume que hay un efecto positivo o negativo en la asociación de las variables y da el mismo fenómeno de tomar como base una distribución que se conoce como mínimo informativa escalando por la la escala de los datos específicos que le estemos pasando ahora estos son los defaults que utilizó el restanar pero sería digamos poco si no pudiéramos nosotros especificar nuestras distribuciones iniciales digamos que sería una limitante del paquete pero afortunadamente nosotros tenemos una manera relativamente directa que no es completamente flexible pero tiene cierto grado de flexibilidad de especificar nuestras propias distribuciones iniciales en particular los tres tipos de distribuciones iniciales que imprimíamos en la diapositiva anterior podemos especificar una distribución inicial para el intercepto para los coeficientes o para signa y lo hacemos con parámetros adicionales recordemos que esta primer línea del código están glem la fórmula y los datos era el la llamada que habíamos hecho originalmente y ahora podríamos agregarle tres parámetros adicionales prior intercept para el intercepto prior sin apellido para los coeficientes y prior aux para estos parámetros auxiliares que en el caso de la regresión lineal es sigma y si se fijan la sintaxis es la misma que nos imprimía acá normal location y scale entonces podemos decirle sabes que el intercepto quiero que sea nada más para ejemplificar sea una normal con locación cero y escala uno y le podemos dar este parámetro extra de lógico si queremos que haga esta este escalamiento automático o no entonces digamos que la primer parte location y scale es esta especificada esta primera línea donde imprimíamos specified y si queremos que se salte de la specified al ajuste le tendríamos que dar autoscale true sino falls o bien nosotros podemos digamos que fijar los parámetros de locación y de escala que consideramos en saptos para nuestro problema en específico la distribución de sigma les comentaba normalmente se utilizan gamas o exponenciales en este caso podemos hacer un exponencial con un parámetro de de tasa 5 igual con autoscale true of false dependiendo de si queremos que haga el ajuste por la escala de los datos o no si ya nosotros directamente le estamos dando la escala que consideramos podríamos hacer esta llamada y obtener el siguiente gráfico de las distribuciones posteriores y si se fijan cambia bastante de del gráfico que teníamos antes no creo que ya son varias dispositivas atrás pero si cambia entonces aquí me voy a ir digamos a al otro extremo de digamos del abogado del diablo y es que una de las críticas que se le hace normalmente a los procedimientos vallesianos es que pueden llegar a ser sensibles a la decisión que tomemos sobre las distribuciones iniciales sobre todo cuando no hay muchos datos recordemos que eran 32 observaciones aquí en este ejemplo de juguete entonces una distribución inicial puede afectar distinta puede cambiar las inferencias pero desde el punto de vista vallesiano esto no es un problema porque la probabilidad es una medida de incertidumbre y distintas personas podrían tener distintas conocimiento recordemos el caso inicial donde alguien que no sabe nada puede decir no sé nada pero una persona experta en el tema tiene ya una opinión distinta formada y puede dar una distribución distinta informativa o mínimo informativa no tan uniforme esto no se trata de moverle a los datos ni que no existe una respuesta correcta sino simplemente digamos que los vallesianos o muchas personas dentro del paradigma vallesiano buscan transparentar que las los y las personas que hacemos estadística tenemos ciertos grados de libertad donde decidimos cosas de manera subjetiva por qué decidimos un modelo normal por qué decidimos un modelo lineal generalizado por qué estamos haciendo regresión y no algún método de machine learning en la práctica digamos existen muchas decisiones que hasta cierto punto podrían ser arbitrarias y es importante reconocer que existe esa arbitraridad y lo mejor que podamos hacer es ser transparentes y decir cuáles fueron esas decisiones que tomamos y eso lo permite permite que en los modelos que presentamos también sean criticables cuestionables en caso de que alguien considere que no es un buen modelo no se trata de manipular los datos sino de tener buenos modelos y ser transparente en las decisiones de modelado que estemos tomando y la distribución inicial es una decisión de modelado como pueden haber otras como mencionaba el modelo específico no utilizar regresión utilizar bosques aleatorios utilizar regresiones regularizadas como lazo reach etcétera entonces la lección digamos más que lección el mensaje yo quisiera darles es no hay que tenerle miedo a cuando utilicemos el paquete nosotros proporcionarle las distribuciones informativas o mínimas informativas que queramos sobre todo por ejemplo en el caso del ejemplo médico decía que íbamos teniendo información secuencial puede ser que tomemos como distribuciones iniciales los resultados de un modelo anterior que hayamos hecho de un ejercicio anterior si estamos trabajando con datos anuales quizás tenemos cada año que presentarle en nuestras empresas un modelo un dashboard no y cada año tenemos que hacer estos procedimientos pues podríamos utilizar como distribuciones iniciales de el año de cada año las estimaciones del año pasado por ejemplo y una buena práctica es aumentar esta escala un poco para que tampoco estemos este casándonos con el resultado del año anterior no sino que sí si permitamos a priori variabilidad pero digamos que no necesariamente tenemos que siempre poner coeficiente cero y escalas muy amplias hay que hay que ser digamos transparentes en lo que estamos haciendo pero sin miedo a utilizar la información que sí tenemos no nosotros estamos realizando un análisis podemos tener esa información hay contextos distintos donde por ejemplo por ética quizás en un ejemplo médico si sea requerido buscar una distribución no informativa desde el principio porque por ética tenemos que no asumir digamos que yo soy un experto y entonces confío demasiado en mi opinión o por ejemplo en contextos políticos hay veces que autoridades electorales pueden estar utilizando modelos y tienen que ser árbitros y entonces a priori no podrían estar quieren estimar la proporción de votos de los candidatos y puede ser que a priori aunque sepamos que van a haber candidaturas que son más probables de ganar no y personas que que no van a recibir tantos votos a priori un árbitro obtendrá que poner distribuciones uniformes o equiprobables entonces lo importante es ser honestos con el concepto contexto perdón y transparentes ahora perdón hasta aquí algún otra otra duda ahora esta fue una regreso un ejemplo de una regresión normal y niña pero mencionaba que los modelos lineales generalizados incluyen otro tipo de regresiones por ejemplo la regresión logística y tanto el contexto base de r con la función glem como su equivalente vallesiano en están glem permite utilizar modelos de regresión logística entonces la regresión logística recordemos nuestra variable de respuesta son bernulís o binomiales ceros y unos éxitos y fracasos y nuestra función liga es el logit entonces la sintaxis es la siguiente tenemos los mismos dos parámetros anteriores una fórmula y unos datos aquí por ejemplo ya voy a utilizar dos variables explicativas además del intercepto ya no solo una para que también veamos el uso de este más en la fórmula entonces tenemos j switch de lado izquierdo de la tilde y del lado derecho nuestras variables x dist 100 y arsenic este es un ejemplo también muy utilizado creo que originalmente es de gelman y sus colaboradores los datos de pozos en vangladesh y la idea es ver si las personas se cambian de de pozo de agua debido a la distancia en metros o en kilómetros que tienen a distintos pozos disponibles y las concentraciones de arsenico que están reportadas en los pozos entonces la idea es que uno buscaría tener acceso a agua más limpia y que no se contaminara pero como somos los seres humanos buscando conveniencia pues hay veces que preferimos el pozo más cercano no sobre todo si estamos hablando de tener que ir y sacar agua y cargar no las cubetas regresar pues hay veces que ni modo no que haya un poco más de contaminantes prefieren las personas un pozo más cercano y en ese contexto es el modelo no de los de los datos la única diferencia digamos desde el punto de vista del código es que ahora agregamos un parámetro de familia estamos cambiando el modelo ya no es un modelo normal ahora es un modelo binomial entonces el primer componente modelo lineal generalizado es lo que estamos cambiando ya la familia ya es binomial y la función liga ya no es la identidad sino el logic aquí por ejemplo existen otro tipo de regresiones con modelo binomial por ejemplo la probit en ese caso la liga sería probit en lugar de logit no una regresión probit en lugar de logística por si algunos de ustedes conocen o las regresiones probit y logística es digamos que la misma sintaxis necesitamos nada más ajustar de logit a probit el link dentro de la familia y de nuevo la belleza de r estanar es que permite tener la misma sintaxis que la función frecuentista de r y es exactamente lo mismo excepto que con estan de digamos de prem de prefijo y si corramos si corremos esto podríamos también tomar los resultados y las distribuciones posteriores serían las siguientes no para el intercepto para esta distancia en kilómetros o en 100 metros no no no recuerdo 100 metros y para la concentración de arsénico no entonces digamos que podemos hablar ya de digamos las consecuencias del modelo vemos que a mayores concentraciones de arsénico aumentan la probabilidad de que se cambien las personas de pozo si se dan cuenta de que está contaminado el pozo al que van normalmente pues tienden a cambiar de pozo pero entre más lejos esté el otro pozo pues disminuye la probabilidad de que se cambien no si está muy lejos la alternativa no me no me quiero no me quiero cambiar una pregunta claro y que quiere decir cuando el intercepto se sobrepone al decir al cero muy buena muy buena pregunta justo ahí habría digamos que dentro de la regresión logística siempre es como truculento a veces el interpretar los coeficientes el intercepto en este caso recordemos que está en el predictor lineal entonces un intercepto de cero supongamos que tenemos una distancia de cero es decir estamos plantados en un pozo y la concentración de arsénico es cero entonces es un pozo completamente limpio y la alternativa está también ahí no no es no es algo muy realista no porque tendríamos como que dos pozos exactamente en el mismo lugar en este contexto uno completamente limpio el suelo y ahí el el pozo alternativo pero tendríamos que interpretar ese coeficiente de esa manera para que la contribución al predictor lineal sea solo la del intercepto si queremos hacer cero la variable dist cien y cero la variable arsénico para que el producto sea cero y sólo nos quedemos con el intercepto y al transformar el cero en la escala logística con la función inversa del logit inverso a la probabilidad resulta que la p de la binomial o la teta de la binomial es punto cinco entonces cero en la escala logit es una probabilidad de punto cinco que quizás en este mundo repito el problema aquí es que luego los interceptos en relaciones logísticas no tienen una interpretación muy intuitiva pero quizás aquí tiene algo de sentido en escaso no es pues si yo tengo un pozo limpio al lado del otro pues pues quizás es un volado no sea punto cinco de que me cambio no porque son prácticamente intercambiables los los dos pozos no sé si eso perdón responde tu tu pregunta es un poco complicado pero creo que sí pero sí pero voy voy a digamos aprovecharla para mencionar algo no específico en este caso del intercepto pero supongamos que hubiera sido uno de los de los otros coeficientes que estuviera cerca del cero esa es también una de las ventajas de esta distribución posterior y es que no sé si les ha pasado que las pruebas de hipótesis por ejemplo luego son complicadas y son en mi pueblo decimos rebobujadas no una prueba de hipótesis y los valores pe como que son conceptos muy útiles pero que sí son un poco técnicos y a veces la de nuevo como que el impulso que tenemos es de decir oye cuál es la probabilidad de que sea cero o de que sea muy cercano a cero el el el efecto con esta distribución posterior si podemos tener esta interpretación directa entonces si un coeficiente su distribución está alrededor del cero si tenemos un intervalo de probabilidad que contiene al cero podríamos decir digamos que ese es un la prueba que la probabilidad no de que el verdadero valor esté dentro del intervalo es 95 por ciento entonces más bien cuando tenemos estos dos ejemplos dist 100 y arsenic con distribuciones completamente separadas del cero podríamos hablar del equivalente del concepto de que son variables significativas no no se utiliza tanto en el concepto vallesiano el término significativo porque no hablamos de un valor pe específico ni un nivel de significancia pero intuitivamente si si es esta si que podemos decir esta cosa de pues la probabilidad de que sea menor a cero por ejemplo el coeficiente de dist 100 es prácticamente uno no todas nuestras nuestra distribución está a la izquierda del cero no en los negativos obviamente esta es una representación de la distribución la distribución es sobre todos los números entonces técnicamente no es 100 por ciento pero es prácticamente 100 por ciento la probabilidad de que el coeficiente sea negativo o de que sea la variable significativa pero pues aquí por la hora creo que llamamos no sé si prefieran hacer una pausa de cinco minutitos o prefieren que me que continúe no sé cómo estamos por qué podemos seguir pero me adaptó si si seguimos aquí he hablado de las distribuciones posteriores o les hemos mostrado estas distribuciones pero no de dónde salieron no si nos regresamos pues nada más hemos hablado del del objeto pero no explorado el pero de las sintaxis de de están glem pero no hemos dicho cuál es el objeto que regresa esta función de están glem para hablar de eso voy un poquito a que está detrás de restanar y lo que está detrás es un método computacional que se conoce como marcof chain montecar cadenas de marcof o montecar lo vía cadenas de marcof justamente el estimador puntual de el valor esperado de una de un parámetro es una meta común en la estadística como que querríamos o es muy natural utilizar el valor esperado de algo como nuestro estimador pero en modelos a veces muy sencillos la regresión logística es uno de los casos es decir es un modelo digamos muy sencillos en el sentido de que es de los primeros modelos estadísticos que podríamos pensar los valores esperados son integrales y hay veces que no podemos calcular esas integrales por ejemplo también esa constante de normalización que que yucep mencionaba de la marginal en el teorema de valles ese denominador que es un integral en modelos como la regresión logística y más complejos no la podemos encontrar analíticamente podrán estar ahí arrastrando el lápiz todo lo que quieran y no tiene una solución analítica no qué podemos hacer pues podemos utilizar simulación y las computadoras eso es lo que nos permiten no necesitamos o no es 100% necesario hacer todos los cálculos y han explotado no los métodos estadísticos gracias al poder de cómputo en particular montecarlo y esta noción de cadenas de marco la idea de montecarlo es que vamos a tener simulaciones de variables independientes en la computadora y las cadenas de marco o montecarlo vía cadenas de marco esas simulaciones van a tener una correlación van a estar encadenadas y por eso llaman cadenas de de marco la idea que está detrás de estos métodos es la ley de grandes números para el caso de montecarlo y se llama teoremas ergódicos para mc mc y es muy intuitiva no si nosotros simulamos un montón de eventos podemos estimar con mayor precisión lo que queremos en este caso estoy graficando por ejemplo cuatro simulaciones independientes más bien corridas independientes de volados y aquí es muy útil la interpretación frecuentista de la probabilidad si tenemos mil volados en el eje y voy calculando las frecuencias relativas cuántas han caído en méxico decimos águila sol no sé caro cruz no en algunos otros países creo es es como se dice no el volado si vemos la frecuencia relativa cuando lanzamos solo un volado y cuán que proporciona caído cuando lanzamos 2 3 4 10 así hasta mil podemos ver cómo va acercándose a esta probabilidad teórica de punto 5 y no importa si hacemos una cadena o una secuencia u otras no aquí hice 4 corridas independientes de de mil y cada una de esas va convergiendo no no son exactamente iguales y obviamente esta convergencia es mejor entre más simulemos pero la idea es que está detrás de esto es si yo tengo una muestra de realizaciones independientes o correlacionadas suficientemente grande puedo estimar con mucha precisión lo que yo quiero en particular el objetivo va a ser estimar estos valores esperados o asociar a una distribución teórica un histograma o puntos en realidad ustedes imaginarán que estas distribuciones que yo he estado graficando pues en realidad son una colección finita de puntos no es no la curva específica sino y se nota que un poquito con las formas pues en la computadora tenemos 10 mil numeritos y estamos graficando frecuencias o alturas como si fuera un histograma o un estimador de densidad y estamos asociando esa muestra finita a una distribución pero lo vamos a hacer de una manera teóricamente válida que es Montecarlo o mc mc y en particular restan o restan arm que está detrás de yo el lenguaje de programación están utiliza una técnica en particular que se llama jamiltonia montecarlo una familia de métodos que se llama jamiltonia montecarlo y que me gusta pensar o a la gente lo explica normalmente como un sistema de física imaginaria y es una cuestión que a mí me parece fenomenal y me encanta hablar de eso porque es mind blowing dirían la idea es que pensemos en esa distribución o esa densidad que que tenemos la distribución posterior recordamos ese es nuestro nuestro objetivo como una colina y la idea es que querríamos tener más muestra más realizaciones en el centro cerca de la moda de la media de la mediana bueno de la mediana pero de la moda y de la y de la del valor esperado más muestra y menos muestra en las colas de las distribuciones equivalentemente podríamos verla como un tazón si le aplicamos el menos logaritmo estamos volteando la distribución una solo una transformación y podemos pensar que queremos más muestra en el fondo del tazón que está equi equi parado con la moda de la distribución y menos muestra hacia los extremos hacia los bordes de ese tazón que es nuestra distribución densidad volteada ahora imagínense esto es completamente imaginario que ese tazón fuera real y que fuera un sistema físico y tuviéramos una partícula o una pelotita que rodara en ese tazón siguiendo las leyes de la física pues si se juegan de la primaria o la algunos de ustedes quizás sabrán mucho más física que que yo yo no soy físico el el la pelotita estaría intercambiando energía cinética y potencial y si no hubiera fricción estaría subiendo y bajando subiendo y bajando de acuerdo a las leyes de la podríamos pensar en ver este tazón pero también graficar y el tazón perdón estoy graficando en el eje x una variable teta que ese es nuestro parámetro que estábamos interesados la la distribución posterior estamos pensando que es el tazón pero también el momentum de esa partícula la inercia que tiene esa pelotita que está subiendo y bajando pues resulta que si graficáramos conjuntamente la posición en el eje x de la pelotita junto con el momentum los físicos le llaman espacio de fases se vuelven unas elipses o unas trayectorias de energía constante porque siempre está cambiando energía cinética por energía potencial no no quiero perdón aquí entrar como en estos detalles de física porque repito es completamente imaginario ya no se preocupen pues por tanto por el tema de la de la física acá pero la idea básica es que podemos encontrar o graficar estas trayectorias o estas órbitas asociadas a una partícula que está bailando o rodando en el tazón que es nuestra distribución y vamos a pensar en un proceso de 12 etapas empezamos en una posición dada ponemos la partícula la pelotita en algún lado de nuestro tazón y a esa posición le vamos a llamar teta 0 y vamos a pensar en un momento maleatorio eso quiere decir que le vamos a aplicar una fuerza aleatoria poquita fuerza o mucha fuerza y una dirección a la pelotita la empujamos y eso determina el momento y entonces colocamos a la partícula en una órbita en particular esa órbita define una trayectoria que las leyes de la física determinan esa órbita está completamente determinada por las leyes de la física entonces vamos a seguir esa trayectoria por un tiempo tau y cuando pasa el tiempo tau vemos en dónde está la pelotita en qué posición entonces con la flecha aquí naranja se fue moviendo la pelotita y registramos donde cayó dado perdón donde cayó teta y vamos a hacer este proceso una y otra y otra vez empezamos aquí a la derecha no sé si se ve mi cursor muestreamos empujando la pelotita en una dirección aleatoria con una fuerza aleatoria seguimos a la pelotita un tiempo y vemos en dónde está en el momento en que la detenemos y vemos en dónde está le aplicamos otra fuerza completamente aleatoria y eso hace que salte a otra de estas órbitas que también están determinadas la trayectoria se mueve por otro tiempo tau vemos donde cayó y registramos esa posición y vamos repitiendo ese proceso y procedimiento estas son 10 iteraciones de ese procedimiento pues lo mágico es que esta idea funciona para lo que queremos y es que si nosotros hacemos muchas de estas iteraciones y vemos en dónde estaban cayendo las pelotitas después de los tiempos tau que íbamos observando resulta que podemos reconstruir el tazón original vamos a tener más muestra en el hacia el centro o hacia el fondo del tazón menos muestra hacia los bordes del tazón es más difícil que una pelotita suba mucho no y más bien como siempre está bajando y rodando pasando por el centro pues es más probable que cuando la detenemos esté cerca del fondo del tazón y resulta que existe un procedimiento para que es para estas iteraciones nos den una descripción de la distribución original que estamos en las que estamos interesados y esto nos sirve porque podemos pensar en nuestros modelos estadísticos como completamente imaginario pero como sistemas físicos y la constante de normalización justamente que que yo sé mencionaba al inicio en el teorema de valles no nos que siempre es la que nos meten problemas porque es un integral complicada resulta que como no depende de teta no afecta la dinámica de esta partícula inventada y si recuerdan a veces decían que la energía potencial es proporcional a la altura y estamos pensando en la altura como el producto de la verosimilitud y la distribución inicial volteada por el logaritmo para que tenga esta forma de tazón y entonces efectivamente hablemos de una altura pues resulta que esta analogía funciona y podemos pedirle a la computadora que simule empujar pelotitas en estas curvas que son productos de las verosimilitudes y las distribuciones iniciales que les quedemos y siga las trayectorias de la física y nos devuelva entonces realizaciones del modelo posterior que estamos interesados de teta dado los los datos y o x esto es digamos conceptualmente la idea de jambiltonia nonte carlo obviamente hay mucha teoría y justificaciones de por qué funciona cómo funciona la implementación específica que se utiliza y muchos trucos de programación detrás del paquete pero es importante como mencionarlo para que tengamos la idea de que lo que estamos pidiendo le hay restan arma es que simule este sistema físico y nos regrese una secuencia de posiciones estas iteraciones forman la cadena de marcos por eso llama mc mc cadenas de marco marcos chain monte carlo porque nos va a ir generando estos rayitas naranjas que puse en el eje x que están correlacionadas porque dependen la una de las otras de las otras no están aquí vinculadas entonces son una cadena pero al final esa cadena si juntamos todas las observaciones y la pensamos como una muestra nos delinea las distribuciones posteriores y lo bueno es que tenemos una muestra de esa distribución y nos va a permitir propagar a equiparar la muestra con la distribución y hacer estimaciones como si tuviéramos una muestra como si hubiéramos hecho no un ejercicio físico de los ejemplos de sus primeros cursos de probabilidad de por ejemplo tengo una urna con muchas pelotitas y saco una al azar no este sería un poco el el equivalente no aquí nuestra urna tiene todas estas rayitas que fueron las posiciones que generó el la computadora y si nos preguntamos por ejemplo la probabilidad de que beta sea mayor a cero pues podemos ver de todas las betas que generó nuestro mecanismo en la computadora cuantas fueron mayores a cero y cuantas no y con eso podemos dar un estimado de esa probabilidad y por eso podemos hablar directamente de la probabilidad de que beta sea menor a cero la probabilidad o que beta esté entre 5 y menos 5 con probabilidad 80 por ciento y por eso es muy útil digamos o muy atractiva hoy día la estadística vallesian ahorita pero este sí ya es si se fijan la última diapositiva entonces aquí creo que si es momento de hacer una breve pausa porque todavía tenemos si no me equivocó 40 minutos más o menos también para perdón servirme un vaso de agua y me voy a pasar después de esta pausa de 5 minutos a r studio para mostrarles digamos un objeto de restanar y cómo se va es como se ve digamos cuando le damos clic a esa a esa función porque aquí si es diferente digamos a la tradicional de glem y el motivo por el que todo el output de restanar es diferente es esto no entonces si quería mencionar un poco esta noción antes de ver en r studio el output que nos genera el restanar pero antes de hacer una pausa es alguien tiene alguna duda sé perdón que esto es como digamos como que muy volado y teórico pero pero si alguien tiene una duda con mucho gusto podemos podemos platicar de ello bueno solo tal vez que aquí es una dimensión no supongo que con la regresión es lo mismo pero multidimensional y la posición de la pelotita es un vector de tantas componentes como tantos parámetros exactamente es es correcto la la aquí y de hecho la interpretación como la pelotita es sobre todo muy intuitiva en pocas dimensiones no obviamente si tenemos muchos parámetros ya como que imaginarte un tazón de muchas dimensiones es es como pues poco poco realista no o fuera de nuestras capacidades y en realidad pasan cosas extravagantes de hecho la muestra ya no se concentra o ya no se debería concentrar solo en el fondo del tazón sino en realidad se forman como serillos o o como unas geometrías muy raras pero pero digamos que en pocas dimensiones esta sería la interpretación y si tienes toda la razón obviamente en muchas dimensiones lo que tenemos es un vector de teta no tenemos un cada una de estas posiciones se convierte en un vector de beta 0 beta 1 beta 2 así todas las realizaciones de tantos parámetros como como tengas perfecto y una cosa es lo objetivo al final de esto es sacar las muestras de donde hay más densidad en la distribución posterior no exacto es el cometido del método es no perder tiempo generando muestras de no donde la probabilidad es baja exactamente porque si justo lo que no quieres es perder tu tiempo en regiones donde no es probable que esté el parámetro quieres tenerlas en donde si es probable pero usaste el término correcto de alguna manera densidad o sea ahí es donde es un poquito como técnico el asunto porque el objetivo es o más bien lo que busca tener mayor eficiencia estos métodos es en poder estimar valores esperados que son estas integrales entonces combina un poquito la la probabilidad con el valor de la de la función que estés que estés buscando entonces como que digamos es si busca no perder el tiempo el método en lugares o en regiones de poca contribución al valor esperado pero pero sí lo que me parece un poquito raro es el logaritmo pero bueno supongo que tiene que ver con la transformación que sea menos el logaritmo si bueno es en realidad es un artificio digamos un poquito de para darles de intuición de un tazón la volteas pero también es por motivos numéricos no para que la computadora funcione y porque si tenemos solo probabilidades sin el menos logaritmo podemos tener underflow este y se nos van a cero todos los números entonces la estabilidad numérica es pero digamos que tiene el truco de que se convierte en algo intuitivo o más intuitivo yo no sé qué tan intuitivo les pareció verdad no es correcto pero bueno digamos que que tiene esas dos roles no el principal es que le da estabilidad numérica al algoritmo pero también tiene la ventaja de que se le puede dar una interpretación perfecto gracias ok entonces voy a hacer una una pausa creo que tengo que post-recording y perdón en cuatro minutos al cuarto para la hora regreso en factos mencionábamos que teníamos aquí los datos ya modificados y este es el el comando de glem que tradicionalmente utilizaríamos para hacer regresión frecuentista en r no y el resultado es lo que lo que habíamos visto en en xaringan y obviamente podríamos también explorar el la estructura de de todo lo que nos regresa el objeto de de la regresión y tenemos coeficientes residuales valores ajustados efectos las medidas no de r todos los digamos los resultados que nos da la pero creo que sí pero el glem pero qué pasa cuando queremos usar rstan arma no es decir la función que tiene el prefijo de stan le damos la fórmula los datos la familia y voy a utilizar aquí los defaults de distribuciones iniciales que auto ajustan a la escala de los datos de rstan si corremos empezamos a ver que nos imprime un montón de cosas en computadora entonces entenderán todo esto no no lo podía poner en xaringan y qué es lo que nos va imprimiendo aquí lo está imprimiendo esto de manera secuencial no no estoy utilizando cómpluto paralelo pero si se fijan nos habla de cadenas cadena 1 cadena 2 cadena 3 y cadena 4 y estas serían las cadenas que podrían paralelizarse pues podríamos correr distintas cadenas en paralelo ya que nos referimos con cadenas pues a una secuencia de estas iteraciones del algoritmo de jamiltonia monte estos primeros mensajes para hacer el el algoritmo hay que hacer unas evaluaciones de gradientes y nos habla de unas transiciones que toman algunos segundos y entonces nos piden que ajustemos nuestras expectativas no no le dan mucho caso digamos al valor específico de tiempo que les digan aquí sino más bien relativamente piensen en que si ven un número grande o relativamente grande quiere decir que el algoritmo va a tardar en terminar y obviamente querríamos poco tiempo y y algoritmos eficientes y nos va imprimiendo cómo van las iteraciones hay una etapa que se llama de calentamiento eso digamos que escapa un poco los objetivos de hoy pero la idea es que va calentando el algoritmo y ya eventualmente empieza a muestrear correctamente o esperamos que correctamente de la distribución posterior en las que estamos interesados al final es para estos ejemplos sencillos es bastante eficiente de los algoritmos de mc mc tardó menos de 2 segundos en esta cadena completamente tanto el calentamiento como ya el muestreo y si se fijan hubo 2 mil interacciones en total de las cuales mil fueron de calentamiento entonces tenemos mil resultados mil posiciones al final y efectivamente como decía yuseb en realidad tenemos mil pares ordenados no pues un vector porque te perdón tres porque tenemos el intercepto y cuatro con las y me andaba haciendo esas son es todo un vector no el intercepto beta 1 beta 2 más sigma y hace lo mismo cuatro veces no por default los defaults de restanar son hacer cuatro cadenas cada una de 2 mil con mil iteraciones de calentamiento y mil de muestreo solo vamos a quedarnos con las de muestreo las observaciones de calentamiento no nos sirven digamos para los gráficos de las densidades porque queremos cuatro cadenas y no sólo una digamos para tener un poco mayor deficiencia tanto desde el punto de vista computacional como les mencionaba si lo podemos hacer en paralelo podemos ir generando más muestra al mismo tiempo y no tenemos que hacer una sola cadena para contar cuatro mil si queremos una muestra final de cuatro mil pues no tenemos que esperar hacer las cuatro mil sino que en paralelo hacemos mil mil mil y entonces tenemos más en menos tiempo pero también por validez como estamos haciendo cadenas de de marcov no queremos que una sola cadena domine digamos las inferencias que vayamos a hacer con esa muestra específica y si empezamos con distintos valores las cadenas y tenemos varias cadenas nos estamos protegiendo un poquito de de ese posible problema digamos de depender solo de una semilla en particular no del valor del algoritmo con el que empezamos y tiene digamos también su componente teórico matemático o analítico de porque queremos varias varias cadenas si ustedes tienen una computadora con muchos cores y pueden hacer muchas cadenas lo pueden lo pueden aprovechar tampoco es completamente necesario hacer todo en el extremo no mil cadenas paralelas se puede encontrar un balance y normalmente por eso el default es cuatro seis siete ocho cadenas dependiendo de qué tan todos recursos computacionales tengamos y al mismo tiempo también podemos aumentar el número de iteraciones si queremos inferencias cada vez más precisas el final recordemos que pues la la estimación que podemos hacer tanto de la densidad como de las preguntas que nos hagamos por ejemplo de la probabilidad de que el coeficiente sea menor a cero o mayor a cero que tan buenas sean nuestras respuestas con esta muestra dependerá de cuánta muestra tengamos entre más muestra tengamos más iteraciones generemos mejores van a ser las las estimaciones pero eso tiene el costo computacional entonces de nuevo como todo en muy frecuentemente en ciencia de datos y en estadística pues tenemos que hacer estos trade-offs entonces este es como el el output el verboz no que que le llaman de están glem ahora podemos obtener lo que aquí le vamos a llamar esta posterior como una matriz y en realidad también como es para ver qué es lo que qué es lo que tenemos y si se fijan pues tenemos la los los coeficientes no no perdón que quise perdón huelz modifier quiero posterior que a la correcta tenemos los cuervos ahora sí ya los tres los tres coeficientes de nuestro predictor lineal el intercepto el de dist 100 y el del arsénico y son un vector cada una de las iteraciones y en total tenemos 4.000 entradas aquí apiló las cuatro cadenas una vez que corremos en paralelo ya podemos juntar las muestras y el código para generar las densidades viene de esta de este otro paquete de biasplot la función mc mc áreas en el que le pasamos esa muestra esa matriz de de simulaciones de la distribución posterior podemos elegir qué parámetros mostrar entonces por ejemplo si no quisiéramos mostrar el intercepto y sólo queremos las distribuciones posteriores de dist 100 y de arsénico se lo pasamos al parámetro de paris y la probabilidad de nuestro de nuestro intervalo de probabilidad sombrado entonces en el ejemplo de la presentación era con 95 aquí se ve el el gráfico podríamos cambiarle por ejemplo a un intervalo del 80 por ciento muchas veces es es necesario cuando hay modelos más complejos tenemos mayor incertidumbre a veces no podemos dar intervalos de o si podemos pero lo no queremos porque no es lo más útil no lo más operacionable dar intervalos al 95 por ciento necesitamos dar intervalos más más cortos con el riesgo que eso que eso implica estamos bajando la probabilidad de que es la cobertura del intervalo pero digamos como para que vean también el parámetro de esta de esta función que es muy útil y ahora les perdón hasta aquí alguna duda o comentarios sobre este output de las cadenas ahora aquí tengo otros datos ficticios y se ven un poco así es un gráfico muy muy elaborado ni los datos son muy elaborados pero se ven se ven de esta manera y tenemos una variable zeta que la estamos tratando como como factor entonces como grupos una variable x1 y una y una respuesta la respuesta son conteos son son discretas no si se fijan están en el cero en el uno dos etcétera y entonces podríamos pensar en un modelo poason pero así como hablamos de modelos lineales generalizados también es común utilizar modelos jerárquicos o multinivel también se conocen como modelos de efectos mixtos aunque muchas de las personas involucradas con estan y y desde el punto de vista vayasiano no no les gusta tanto el término o no utilizan mucho el término efectos fijos y efectos mixtos y efectos aleatorios que son los los conceptos que quizás algunos de ustedes hayan escuchado antes y más bien prefieren los términos modelos jerárquicos o multinivel y aquí la idea es que permites que tienen unas variables de grupo y que los coeficientes sean distintos o varían a través de los diferentes grupos obviamente podemos interpretar los modelos de distintas maneras y hay literatura al respecto y se pueden ver los modelos de distintas formas como interacciones por eso efectos fijos y aleatorios a veces pero quería mencionarlo aquí porque es digamos una función distinta en caso de que quisiéramos utilizar modelos multinivel y la fórmula o la sintaxis de la fórmula cambia también un poquito esta función si se fijan ahora dice glem en lugar de glem y esto viene de que también existe el el equivalente para para frecundista y lo que estamos haciendo es la primera parte de la fórmula y la sintaxis es la misma nada más cambia la función si queremos modelos multinivel o defectos mixtos vamos a usar glem en lugar de glem pero la primer parte de la fórmula es igual pero le estoy agregando un término de entre paréntesis y con esta barra y esta sintaxis de incorporar entre paréntesis una variable condicional la barra podemos interpretarla como condicional en otra refleja que queremos que el coeficiente de x1 varíe o sea distinto para cada grupo de las de las etas entonces esta es la la sintaxis de de un modelo muy básico de multinivel obviamente esto se podría complicar y ustedes después podrían podrían viendo la documentación y haciendo más más ejemplos más realistas y que se adapten a sus necesidades podrían agregarle más términos a la fórmula pero digamos que para que podamos ver un ejemplo mínimo es este este sería de nueva cuenta le pasamos los datos le pasamos la familia ahora en lugar de la regresión logística o en lugar de una regresión lineal normal como sabemos que son conteos un modelo usual para lineal generalizado para conteos es el modelo lo poason o poason entonces el miembro de la familia exponencial ahora va a ser poason ya no binomial o normal y en la función liga el link va a ser el logarito lo si se fijan la sintaxis es la misma y lo que estamos cambiando aquí entonces es el modelo específico la forma paramétrica y la función liga y estoy agregando a propósito un un parámetro adicional que se llama adapt delta y aquí quizás no es la mejor práctica pero lo estoy poniendo en punto 4 para echarlo a perder a propósito porque quiero que vean lo que puede pasar y pasa muy seguido en la práctica de cuando algo falla entonces voy a correr y empieza igual no genera las cuatro cadenas por 2000 pero ahora al terminar nos aparecen muchos errores o bueno muchos mensajes en rosa de guarniz el más importante y en el por el que lo quería mostrar es el primero que dice que hubo transiciones divergentes después del calentamiento 1490 y nos da una liga en la que nos explican también un poquito qué hacer o qué qué sucede con este con este warning es un warning pero desde el punto de vista teórico es un error entonces por eso no no lo quería omitir y aunque aquí a propósito puse que fallara o busqué que aparecieran las transiciones divergentes la realidad es que pues en la práctica nos va a pasar muchas veces no a propósito no puede ser que que nos genere este error y normalmente es un error no del algoritmo en sí sino nuestro en el sentido de que el modelo que estamos proponiendo no es el mejor para los datos entonces por eso es importante verlo porque si en el pro en el pipeline que tengamos si nos encontramos con modelos ajustados por están arm que nos generan este mensaje de transiciones divergentes es muy importante que nos detengamos a pensar en el modelo que estamos proponiendo y es un grito que nos da están de que algo está no del todo bien entonces algo que podría parecer como muy escandaloso y malo de porque el repito pasa muy frecuentemente no no no sé si muy frecuentemente pero no es extraño que que suceda y podríamos decir no pues entonces eso no me sirve ese algoritmo la verdad es que es una ventaja porque muchas veces en los procedimientos que hagamos podemos irnos como en silencio y no estar seguros de si algo está fallando en los algoritmos o si el modelo que estamos proponiendo ajusta correctamente a los datos y el hecho de que la computadora se queje tan estruendozamente cuando algo falla al final es una como bendición oculta en el sentido de que nos avisa a tiempo de que tenemos que ver el modelo y de que no dejemos pasar algo que quizás no no es este no es lo mejor no pero no siempre es un problema del modelo por eso fui capaz aquí de hacerlo a propósito muchas veces se corrige con este parámetro de adept delta y que está relacionado con la forma en la que la computadora simula las trayectorias de la partícula en estas órbitas que ya como es computadora no puede simular exactamente con toda la precisión del mundo las trayectorias de las ecuaciones de la física tiene que hacerlas de estimar las trayectorias de manera numérica y hay un tamaño de paso o de adaptación en el que va generando esas trayectorias entre más alto sea el parámetro adept delta el tamaño de paso es más chiquito entonces eso le permite al algoritmo seguir de mejor manera las trayectorias y dar buenos resultados cuando es relativamente alto el tamaño de paso es decir relativamente bajo este parámetro adapt delta el algoritmo puede no seguir correctamente la trayectoria se sale no da un paso demasiado grande y se va diverge de la trayectoria y por eso se llaman divergentes y entonces falla no el algoritmo el default no es punto 4 el default es mucho más alto es decir un tamaño de paso más chiquito precisamente para asegurar buenos resultados o intentar asegurar buenos resultados yo a propósito lo puse en un número punto 4 que no es el correcto pero digamos que puede ser que cuando intenten ustedes hacer algún modelo les hace con el valor de default de la adaptación y aquí podemos ver si se fijan la segunda el segundo mensaje es que nos dice que veamos pairs plot pueden utilizar el persplot directo o utilizar base plot tiene esta otra función que se llama msmc pairs que es un poco más estética y les permite ver algunas otras cosas necesita que le pasemos los estos parámetros nudes params para para poder a esta no cargué si perdón quizás no había quedado y le podemos pasar también los los parámetros a perdón no imprimí el modelo de divergencias para que veamos cuáles son los los resultados aquí los nombres de los parámetros son un poco más complicados por este tema de del modelo multinivel en el que tenemos parámetros que cambian por grupo de de zeta y la forma en la que restanar los llama es nos dice b de de grupo y este es el intercepto del nivel uno de zeta esto es digamos el costo que pagamos por por un una función digamos muy general y muy flexible es que pues luego los nombres de los parámetros son son poco complicados pero bueno si si tenemos la estructura de nuestros nuestras simulaciones posteriores aquí vienen los los parámetros y ustedes digamos que podrían ponerle nombres más más interpretables o modificarlos ya de acuerdo a las necesidades no aquí por efectos ilustrativos pues nos vamos con con el nombre feo no este automático y vemos este gráfico de de pares en el que podemos ver las distribuciones posteriores como histogramas aquí del x1 y de marginales y las distribuciones de variadas porque como tenemos como distribución posterior una distribución multivariada el vector para cada iteración tenemos un vector podemos ver coordenada coordenada o podemos ver combinaciones y nos marcan en crucecitas rojas los lugares donde hubo esas transiciones divergentes esto es un diagnóstico muy útil porque cuando no vemos un patrón claro de dónde están localizadas las divergencias eso nos nos da un poco más de confianza en que si aumentamos el parámetro de delta a delta vamos a resolver el problema en cambio si vemos alguna región específica donde están concentrándose no los puntitos rojos nos está apuntando el modelo que ahí hay un problema con el modelo y que probablemente aunque le aumentemos el parámetro de delta y queramos hacer pasos muy chiquitos se va a seguir equivocando el algoritmo y eso nos ayuda repito digamos como que es un modelo quejumbroso pero esa queja del propio modelo es algo muy útil porque preferible que se queje y nos señale errores de nuestro modelo a que se nos vaya digamos por debajo del agua y que es un riesgo que corremos todos los días con muchos modelos entonces este sería el como el diagnóstico que podríamos hacer aquí si se fijan el patrón como que no es muy claro y esto tiene sentido porque el motivo por el que está fallando es que a propósito utilicé un parámetro un mal parámetro de paso qué pasa si aumento el adapt delta entonces si a la delta aquí vamos a ver si pasamos de punto 4 a punto 6 se redujeron las distintas transiciones divergentes a 46 pero siguen habiendo entonces podemos buscar un modelo más con punto 9 también si se fijan conforme aumento el tamaño de este parámetro disminuyó el tamaño de paso y tarda un poquito más entonces ese también es un tema a tener en consideración y tres no podríamos aumentar un poco más creo que el default de están gmr es punto 95 si no me equivoco punto 94 entonces digamos que si corremos el normal voy a correr nada más así ahorita sin perdón los daron un poco más estos por cierto también fueron datos generados un poco a propósito que sean no tan buenos para también hablar de la diferencia en la en la eficiencia no se desacuerdan los primeros modelos que estábamos haciendo tardaban un segundo dos segundos y aquí cuando ponemos datos o modelos que no son perfectos tardan tardan más no si hay que tener en cuenta que el costo de poder tener esta expresión total de la distribución posterior y de poder hacer todos este tipo de enunciados sobre los parámetros interpretándolos como como como objetos probabilísticos es en mc mc en particular pues un costo del algoritmo y de y que tarda un poco más o a veces mucho más y el reto es encontrar buenas parametrizaciones y buenos modelos para hacerlos eficientes el original también tenía pocas pocas divergencias podríamos seguir aumentando el parámetro hasta eliminarlas pero dado que nos quedan 15 minutitos y quiero dejar al final un poquito de tiempo para para preguntas les quiero enseñar también esta otra función de trace que es muy utilizada y que se llama el trace plot o o en español lo conocen muchas personas como como gráfico de de orugas vamos a ver si puedo hacer más grande la no sé qué tan que también se vea para para verlo a ustedes pero tenemos para cada parámetro unos gráficos de líneas que son las cadenas que fuimos generando y en el eje x está el número de iteración de la cero a mil recordemos que sueron mil simulaciones después del calentamiento y en el eje x ya perdón el valor de esa simulación y entonces se generan estas estos gráficos de cadenas y el objetivo para saber que esto es un que como gráfico de diagnóstico lo que se busca es que las cadenas se mezclan es decir que están explorando las mismas regiones si las cadenas se mezclan quiere decir que el algoritmo está explorando las mismas regiones y eso es una señal positiva no no completamente definitiva pero si es bueno y y que no tenga estos patrones claros sino que queremos que se vean como ruido y que esté bailando digamos en una banda que determina esa esa región en la que estamos interesados cuando vemos un patrón muy claro quiere decir que algo algo está mal entonces este es evidentemente un problema en el gráfico de diagnóstico porque vemos que pues esta cadena está como que muy clara de de identificar y se va para para otro lado aquí vemos este patrón como descalerita y lo que queremos es más más ruido por por por cuestiones de tiempo pues no no podríamos como bueno no en tiempo de encontrar otros o presentarles otros ejemplos y a un dar un poco más a detalle de estos diagnósticos pero en las referencias que les daré darles muy buenas referencias la verdad es que la comunidad que desarrolla el paquete de restán de restan arm y en general están son son muy muy buenos y generan muchos recursos sin los cuales pues obviamente no hubiera podido hacer esta esta presentación y que les recomiendo totalmente para que sigan el camino de explorar este paquete de explorar el modelado vallesiano y pues con eso pues les les agradezco mucho la participación espero que haya sido pues de su agrado de utilidad y si tienen alguna pregunta o comentario con muchísimo gusto para para cerrar estos últimos cinco seis minutitos que tenemos Fernando que te parece de los demás paquetes como arm y otros que hacen también vallation análisis por ejemplo otros otros paquetes como v vrms también está basado en están no no no escuché bien cual cual paquete mencionaste no sé si nos lo puedes poner en el en el chat si yo creo arm arm simplemente arm por escribir aquí no sé de hecho si tengo el paquete de rm si creo que no tengo el paquete de de arm entonces la respuesta es no no lo conozco perdón no sé si es una dependencia pero pero pero sí creo que la respuesta más correcta es no no conozco el paquete perdón ok gracias porque porque tiene de hecho también las gráficas que iguales a estas que que la presenta de ahorita en cientos debe ser otra implementación u otro de están este tipo de diagnósticos están muy basados en mc mc otro tipo de mc mc es un poco más general que que jamil tonia en montecarlo que es el algoritmo o la familia de algoritmos que utiliza están hay otros paquetes como jags y como winbox que utilizan otros métodos de mc mc pero yo creo que que lo más state of the art hoy por hoy es es están entonces también por eso digamos decidí presentarles este este paquete y y en las referencias espero que también puedan puedan explorarlo perdón no sé si porque sé que algunos tienen que tienen que irse pero gracias natalia y yonatan es por los comentarios y sí sí con mucho gusto les les pasaré el el código y las slides para para que tengan las tengan ustedes perdón samuel no sé si respondí tu tu pregunta sé que dije que no lo conozco entonces no gracias a pero no más una última pregunta y cómo cómo lo haces para reportar los los los outputs los resultados de una análisis en de un análisis paisiano en comparación a una normal frecuentista o o sí como una hora la creo que una de las cosas es que que se presentan estas estas distribuciones como el ejemplo de las de las necesarias de desplot cambié creo y entonces este tipo de gráficos se presentan y también se pueden dar los intervalos no en lugar de de graficar la distribución completa una forma usual de presentar los resultados es dar los intervalos que serían el equivalente al intervalo de confianza pero pero desde el punto de vista vallesiano son son intervalos de probabilidad y eso permite también que a la hora de explicar tus resultados podamos entre comillas hablar directamente y con este lenguaje coloquial de correctamente o sea no sin riesgo de equivocarnos desde el punto de vista técnico si podemos decir hay un 80 por ciento de probabilidad o un 90 y 5 por ciento de probabilidad de que el parámetro esté entre tal como lo lo generamos por ejemplo también podemos encontrar igual dentro de la función de valles plot por ejemplo hay hay hay otra de intervalos no de la del área sino este intervalos rápido y te puede generar como que los los intervalos entonces existen muchos muchos diagnósticos y como todo pues la verdad es que siempre hay mucha libertad digamos entonces pero el el el la diferencia clave creo es que estamos contando aquí con la distribución posterior completa entonces eso nos permite incluso hacer por ejemplo predicciones más más sencillas porque propagamos la incertidumbre es decir si queremos predecir una nueva observación de podemos utilizar cada una de las simulaciones posteriores para hacer una nueva predicción y tenemos entonces también una distribución completa para la predicción de una nueva observación por ejemplo entonces ese es un poco la diferencia espero a ver si va más o menos claro ok y como y qué le pasa con los priors cómo lo haces para incluir en un modelo por ejemplo por ejemplo tienes tienes un estudio a priori y luego tienes como una idea más o menos de cómo se distribuye tu tus datos y luego haces un estudio más grande después y quieres hacer el primer estudio como un prior como un una información a priori cómo lo incluís en en en un modelo por ejemplo en arstan por ejemplo si si haces una prueba piloto digamos no si si eso quiere decir y entonces podrías utilizar por ejemplo para los resultados de la prueba piloto usas los defaults si quieres de de restanar entonces no especificas iniciales ajustas el modelo a los datos de la prueba piloto y una vez que tienes la estimación de esa prueba piloto es decir tienes el valor esperado posterior y puedes tener una un intervalo de confianza del 99 por ciento incluso del 100 por ciento o sea el mínimo y el máximo de las simulaciones posteriores puedes utilizar esos dos datos para fijar los parámetros de una distribución inicial que utilices en el en para los para el estudio o su secuente para el estudio completo perdón creo que justo se nos se nos sacado ya un poquito el tiempo pero podemos creo que crear el canal de deslac y si podemos seguir platicando un ratito si tienen algunas otras dudas por el canal de deslac les parecen muchas gracias y y pueden incluso hacer ahí una pues un chat ahora le pongo el nombre estudio en bajo el restan arm perfecto muchas gracias a ti y a todos por asistir a justo natalia y a los mil mil gracias a todos voy a parar la grabación y gracias a todos por haber asistido