Actividad evaluativa número 1.

Modelos lineales generalizados

Presentado por:
Julio César Riascos
Diego Fernando Muñoz
María Alejandra Marmol
Michael Steven Morales
Pedro Pablo García Lasso

Presentado a:
PhD. Sergio Alexander Gómez Noguera

UNIVERSIDAD DE NARIÑO
Facultad de Ciencias Naturales y Matemáticas
Maestría en Estadística Aplicada
San Juan de Pasto 2024

1
 1. Punto 1
En el archivo [Link] son descritas características de una muestra aleatoria de 38 vinos
de la marca "Pinot Noir". El objetivo del estudio es relacionar la variable qualidade (calidad)
del vino con las variables explicativas (i) claridade, (ii) aroma, (iii) corpo(cuerpo), (iv) sabor,
(v) aromac (aroma de la barrica de roble) (vi) regiao (1: región 1, 2: región 2, 3: región 3). Una
vez haya cargado los datos use el comando regiao<- factor(regiao) para transformar la variable
región en factor. Haga un análisis descriptivo construyendo por ejemplo boxplots y diagramas
de dispersión de cada variable explicativa contra la variable respuesta qualidade (calidad).
Calcule también las correlaciones lineales entre las variables. Proponga un modelo normal
lineal con todas las variables explicativas. Use el AIC para seleccionar un submodelo.
Interprete los coeficientes estimados.
Análisis descriptivo
De manera inicial, se procedió con un análisis gráfico de la distribución de la variable
cualidad (Qualidade), variable respuesta del modelo que se pretende adelantar. Se observa una
distribución con un buen grado de simetría, por lo cual la postulación del modelo lineal normal
pareciera ser una opción adecuada:
Figura 1.
Distribución de la variable Qualidade

2
 Al abordar el análisis de la variable “qualidade” bajo el diagnóstico propuesto por Cullen
& Frey, se tiene que su distribución bajo bootstraping se aproxima a una beta, una Normal y una
Gamma, lo que posibilita el desarrollo de una aproximación alrededor de la naturaleza descrita
por la variable en cuestión, que permite el desarrollo de los modelos a proponer.
Figura 2.
Collen and Frey

En la base Wine existe la variable de interés y se encuentran seis variables relacionadas

con características de los vinos: cinco de ellas numéricas y una categórica. En el siguiente
gráfico, se resume la inspección visual de la relación entre las variables que serán tomadas como
posibles predictores, y la variable respuesta de calidad:

3
Figura 3.
Relación grafica de variables y calidad del vino

De manera preliminar, se observa que no existe relación entre la claridad y la calidad del
vino, así como entre el aroma de la barrica y calidad. Por otro lado, parece existir una relación
positiva entre el aroma y la calidad, el cuerpo y la calidad, el sabor y la calidad, y finalmente la
región de procedencia y la calidad. Si bien esta exploración visual indica algunas tendencias en
los datos, para poder establecer la magnitud de la relación, se procedió a adelantar análisis
estadísticos de la asociación entre variables: correlaciones entre las variables numéricas, y
prueba de hipótesis para la variable categórica de región y la variable calidad. Los resultados de
los análisis de correlación se muestran en las gráficas subsiguientes:

4
Figura 4.
Correlaciones

Dado que no todas las variables siguen una distribución normal (el caso de la variable
claridad), se procedió a calcular tanto las correlaciones de Pearson como las de Spearman (para
el caso del contraste Qualidade - Claridade), en las cuales se identifican correlaciones de
magnitud considerable entre aroma y calidad (r=0.71), cuerpo y calidad (r=0.55), y finalmente
sabor y calidad (0.79).
En términos de la asociación de la variable respuesta “Claridade” (numérica) y región
(categórica), se adelantó un análisis ANOVA. Con un valor de probabilidad inferior al 5%,
(p=6.59e-08) se encontró evidencia para rechazar 𝐻0 , por lo cual se identifica que existen
diferencias en términos de la calidad del producto dada la región de procedencia. Para verificar la
pertinencia de la prueba de hipótesis elegida, se comprobaron los supuestos de normalidad de los
residuos (con un valor de probabilidad correspondiente al 13.2%) y de homogeneidad de
varianzas (90.04%), comprobando la adecuación de la prueba ANOVA.
Al identificar las diferencias en las comparaciones dos a dos, a través de la prueba post
hoc de Tukey, se identifica la existencia de diferencias significativas entre la calidad de los vinos
5
procedentes de la región 1 y 2 (p=0.02), y entre la calidad de los vinos procedentes de la región 1
y 3 (p<0.0001), así como entre la calidad de los vinos procedentes de las regiones 2 y tres
(p<0.0001). Es decir, cada región imprime un sello diferencial en la calidad del vino, siendo los
mejores aquellos procedentes de la región tres, seguidos de la región 1, y finalmente los de
menor calidad son aquellos de la región 2, como se observa en el siguiente gráfico de valores
post hoc:

Figura 5.
Valores post hoc

Modelo Lineal Propuesto

Considerando que, dadas las características de las variables, inicialmente, se postula un
modelo lineal normal:
Distribución variable de respuesta
𝑌𝑖𝑗 ~𝑖𝑑 𝑁(𝜇𝑖𝑗 ; 𝜎 2 )
Componente sistemático
𝜇𝑖𝑗 = 𝛼 + 𝛽1 𝑥1𝑖 + 𝛽2 𝑥2𝑖 + 𝛽3 𝑥3𝑖 + 𝛽4 𝑥4𝑖 + 𝛽5 𝑥5𝑖 + 𝛾𝑗 𝑥𝑗
Donde: 𝑌𝑖𝑗 = 𝐶𝑎𝑙𝑖𝑑𝑎𝑑 𝑜𝑏𝑠𝑒𝑟𝑣𝑎𝑑𝑎 𝑑𝑒𝑙 𝑖 − 𝑒𝑠𝑖𝑚𝑜 𝑣𝑖𝑛𝑜 𝑒𝑛 𝑙𝑎 𝑗 − 𝑒𝑠𝑖𝑚𝑎 𝑟𝑒𝑔𝑖ó𝑛
𝜇𝑖𝑗 = 𝐶𝑎𝑙𝑖𝑑𝑎𝑑 𝑚𝑒𝑑𝑖𝑎 𝑑𝑒𝑙 𝑖 − 𝑒𝑠𝑖𝑚𝑜 𝑣𝑖𝑛𝑜 𝑒𝑛 𝑙𝑎 𝑗 − 𝑒𝑠𝑖𝑚𝑎 𝑟𝑒𝑔𝑖ó𝑛
𝑥1𝑖 = 𝐶𝑙𝑎𝑟𝑖𝑑𝑎𝑑 del i-ésimo vino
6
 𝑥2𝑖 = 𝐴𝑟𝑜𝑚𝑎 𝑑𝑒𝑙 𝑖 − é𝑠𝑖𝑚𝑜 𝑣𝑖𝑛𝑜
𝑥3𝑖 = 𝐶𝑢𝑒𝑟𝑝𝑜 𝑑𝑒𝑙 𝑖 − é𝑠𝑖𝑚𝑜 𝑣𝑖𝑛𝑜
𝑥4𝑖 = 𝑆𝑎𝑏𝑜𝑟 del i-ésimo vino
𝑥5𝑖 = 𝐴𝑟𝑜𝑚𝑎 𝑑𝑒 𝑙𝑎 𝑏𝑎𝑟𝑟𝑖𝑐𝑎 𝑑𝑒 𝑟𝑜𝑏𝑙𝑒 𝑑𝑒𝑙 𝑖 − é𝑠𝑖𝑚𝑜 𝑣𝑖𝑛𝑜
𝑥𝑗 = 𝑅𝑒𝑔𝑖𝑎𝑜; j=1,2,3.
𝑗 = 1 (Regiao 1)

𝑗 = 2 (Regiao 2)

𝑗 = 3 (Regiao 3)

𝛼 = 𝑖𝑛𝑡𝑒𝑟𝑐𝑒𝑝𝑡𝑜
𝛽1 = Mide el cambio que la claridad ejerce sobre la calidad media del vino manteniendo todo lo demás constante
𝛽2 = Mide el cambio que el aroma ejerce sobre la calidad media del vino manteniendo todo lo demás constante
𝛽3 = Mide el cambio que el cuerpo ejerce sobre la calidad media del vino manteniendo todo lo demás constante
𝛽4 = Mide el cambio que el sabor ejerce sobre la calidad media del vino manteniendo todo lo demás constante
𝛽5
= Mide el cambio que el aroma de la barrica de roble ejerce sobre la calidad media del vino manteniendo todo lo demás constante
𝛾𝑗 = Mide el cambio que la región ejerce sobre la calidad media del vino manteniendo todo lo demás constante

Una vez postulado tal modelo, y previo al proceso de posible interpretación, se procedió
con la identificación del mejor ajuste del modelo, a través del criterio AIC. El mejor modelo
seleccionado fue aquel que contenía las variables aroma, sabor y región como predictores, ya que
fue el modelo con el menor AIC registrado. La formulación de este nuevo modelo lineal normal
obedece a la siguiente postulación:
Distribución variable de respuesta
𝑌𝑖𝑗 ~𝑖𝑑 𝑁(𝜇𝑖𝑗 , 𝜎 2 )
Componente sistemático
𝜇𝑖𝑗 = 𝛼 + 𝛽1 𝑥1𝑖 + 𝛽2 𝑥2𝑖 + 𝛾𝑗 𝑥𝑗
Donde: 𝑌𝑖𝑗 = 𝐶𝑎𝑙𝑖𝑑𝑎𝑑 𝑜𝑏𝑠𝑒𝑟𝑣𝑎𝑑𝑎 𝑑𝑒𝑙 𝑖 − é𝑠𝑖𝑚𝑜 𝑣𝑖𝑛𝑜 𝑒𝑛 𝑙𝑎 𝑗 − 𝑒𝑠𝑖𝑚𝑎 𝑟𝑒𝑔𝑖ó𝑛
𝜇𝑖𝑗 = 𝐶𝑎𝑙𝑖𝑑𝑎𝑑 𝑚𝑒𝑑𝑖𝑎 𝑑𝑒𝑙 𝑖 − é𝑠𝑖𝑚𝑜 𝑣𝑖𝑛𝑜 𝑒𝑛 𝑙𝑎 𝑗 − 𝑒𝑠𝑖𝑚𝑎 𝑟𝑒𝑔𝑖ó𝑛
𝑥1𝑖 = 𝑎𝑟𝑜𝑚𝑎
𝑥2𝑖 = 𝑠𝑎𝑏𝑜𝑟
𝑥𝑗 = 𝑅𝑒𝑔𝑖𝑎𝑜; j=1,2,3.
𝑗 = 1 (Regiao 1), como casilla de referencia.

𝑗 = 2 (Regiao 2)
7
 𝑗 = 3 (Regiao 3)
𝛼 = 𝑖𝑛𝑡𝑒𝑟𝑐𝑒𝑝𝑡𝑜
𝛽1 = Mide el cambio que el aroma ejerce sobre la calidad media del vino manteniendo todo lo demás constante
𝛽2 = Mide el cambio que el sabor ejerce sobre la calidad media del vino manteniendo todo lo demás constante
𝛾𝑗 = Mide el cambio que la región ejerce sobre la calidad media del vino manteniendo todo lo demás constante

Considerando tal modelo, se procedió a adelantar análisis diagnóstico respecto al ajuste

del modelo mediante la identificación del envelope para la distribución normal, así como la
identificación de posibles datos influyentes. El resultado gráfico del ajuste envelope y de las
distancias se muestra a continuación:
Figura 6.
Envelope, Distancia de Cook y residuales

El ajuste del envelope muestra adecuación, y se identifica que las observaciones más
distantes son la 12 y la 20. Se observa un patrón en los residuos estandarizados, similar a un
cono, por lo cual se procede a adelantar el diagnóstico de residuos a través de la función de

8
modelo aditivos generalizados de localización escala y forma (GAMLSS). En tal proceso, se
observaron los siguientes resultados:

Figura 7.
Plot Gamlss

El modelo presenta un patrón de cono identificable en el gráfico de los cuantiles

residuales versus los valores ajustados, lo cual parece sugerir la presencia de heterocedasticidad
residual. Por tal razón, el equipo de trabajo consideró pertinente la postulación de un modelo
Gamma, que permita modelar de una manera más adecuada la calidad del vino.
Modelo Gamma
Considerando que, dadas las características de las variables y teniendo en cuenta que el
modelo lineal tiene ciertas dificultades frente a la distribución de sus residuos y patrones
heterocedásticos, se puede pensar en la pertinencia de postular un modelo tipo Gamma, a
continuación, se presenta la formulación del modelo correspondiente:
Distribución variable de respuesta
𝑌𝑖𝑗 ~𝑖𝑑 𝐺(𝜇𝑖𝑗 , Φ)
9
 Componente sistemático
log (𝜇𝑖𝑗 ) = 𝛼 + 𝛽1 𝑥1𝑖 + 𝛽2 𝑥2𝑖 + 𝛽3 𝑥3𝑖 + 𝛽4 𝑥4𝑖 + 𝛽5 𝑥5𝑖 + 𝛾𝑗 𝑥𝑗
Donde: 𝑌𝑖𝑗 = 𝐶𝑎𝑙𝑖𝑑𝑎𝑑 𝑜𝑏𝑠𝑒𝑟𝑣𝑎𝑑𝑎 𝑑𝑒𝑙 𝑖 − é𝑠𝑖𝑚𝑜 𝑣𝑖𝑛𝑜 𝑒𝑛 𝑙𝑎 𝑗 − é𝑠𝑖𝑚𝑎 𝑟𝑒𝑔𝑖ó𝑛
Log (𝜇𝑖𝑗 ) = 𝐶𝑎𝑙𝑖𝑑𝑎𝑑 𝑚𝑒𝑑𝑖𝑎 𝑑𝑒𝑙 𝑖 − 𝑒𝑠𝑖𝑚𝑜 𝑣𝑖𝑛𝑜 𝑒𝑛 𝑙𝑎 𝑗 − 𝑒𝑠𝑖𝑚𝑎 𝑟𝑒𝑔𝑖ó𝑛
Donde: 𝑌𝑖𝑗 = 𝐶𝑎𝑙𝑖𝑑𝑎𝑑 𝑜𝑏𝑠𝑒𝑟𝑣𝑎𝑑𝑎 𝑑𝑒𝑙 𝑖 − 𝑒𝑠𝑖𝑚𝑜 𝑣𝑖𝑛𝑜 𝑒𝑛 𝑙𝑎 𝑗 − 𝑒𝑠𝑖𝑚𝑎 𝑟𝑒𝑔𝑖ó𝑛
𝜇𝑖𝑗 = 𝐶𝑎𝑙𝑖𝑑𝑎𝑑 𝑚𝑒𝑑𝑖𝑎 𝑑𝑒𝑙 𝑖 − 𝑒𝑠𝑖𝑚𝑜 𝑣𝑖𝑛𝑜 𝑒𝑛 𝑙𝑎 𝑗 − 𝑒𝑠𝑖𝑚𝑎 𝑟𝑒𝑔𝑖ó𝑛
𝑥1𝑖 = 𝐶𝑙𝑎𝑟𝑖𝑑𝑎𝑑 del i-ésimo vino
𝑥2𝑖 = 𝐴𝑟𝑜𝑚𝑎 𝑑𝑒𝑙 𝑖 − é𝑠𝑖𝑚𝑜 𝑣𝑖𝑛𝑜
𝑥3𝑖 = 𝐶𝑢𝑒𝑟𝑝𝑜 𝑑𝑒𝑙 𝑖 − é𝑠𝑖𝑚𝑜 𝑣𝑖𝑛𝑜
𝑥4𝑖 = 𝑆𝑎𝑏𝑜𝑟 del i-ésimo vino
𝑥5𝑖 = 𝐴𝑟𝑜𝑚𝑎 𝑑𝑒 𝑙𝑎 𝑏𝑎𝑟𝑟𝑖𝑐𝑎 𝑑𝑒 𝑟𝑜𝑏𝑙𝑒 𝑑𝑒𝑙 𝑖 − é𝑠𝑖𝑚𝑜 𝑣𝑖𝑛𝑜
𝑥𝑗 = 𝑅𝑒𝑔𝑖𝑎𝑜; j=1,2,3.
𝑗 = 1 (Regiao 1)

𝑗 = 2 (Regiao 2)

𝑗 = 3 (Regiao 3)
𝛼 = 𝑖𝑛𝑡𝑒𝑟𝑐𝑒𝑝𝑡𝑜
𝛽1 = Mide el cambio que la claridad ejerce sobre la calidad media del vino manteniendo todo lo demás constante
𝛽2 = Mide el cambio que el aroma ejerce sobre la calidad media del vino manteniendo todo lo demás constante
𝛽3 = Mide el cambio que el cuerpo ejerce sobre la calidad media del vino manteniendo todo lo demás constante
𝛽4 = Mide el cambio que el sabor ejerce sobre la calidad media del vino manteniendo todo lo demás constante
𝛽5
= Mide el cambio que el aroma de la barrica de roble ejerce sobre la calidad media del vino manteniendo todo lo demás constante
𝛾𝑗 = Mide el cambio que la región ejerce sobre la calidad media del vino manteniendo todo lo demás constante

Se propuso inicialmente un modelo gamma que incluyó todos los predictores disponibles,
el cual fue evaluado mediante el AIC, para determinar el ajuste más óptimo. Los resultados
indicaron que, el modelo con el menor AIC y, por ende, el más adecuado para los datos, fue el
que incluyó únicamente los predictores Sabor y Región. En consecuencia, se plantea el siguiente
modelo basado en esta estructura específica:
Distribución variable de respuesta
𝑌𝑖𝑗 ~𝑖𝑑 𝐺(𝜇𝑖𝑗 , Φ)
Componente sistemático
log (𝜇𝑖𝑗 ) = 𝛼 + 𝛽1 𝑥1𝑖 + 𝛾𝑗 𝑥𝑗
10
 Donde: 𝑌𝑖𝑗 = 𝐶𝑎𝑙𝑖𝑑𝑎𝑑 𝑜𝑏𝑠𝑒𝑟𝑣𝑎𝑑𝑎 𝑑𝑒𝑙 𝑖 − é𝑠𝑖𝑚𝑜 𝑣𝑖𝑛𝑜 𝑒𝑛 𝑙𝑎 𝑗 − é𝑠𝑖𝑚𝑎 𝑟𝑒𝑔𝑖ó𝑛
Log (𝜇𝑖 𝑖) = 𝐶𝑎𝑙𝑖𝑑𝑎𝑑 𝑚𝑒𝑑𝑖𝑎 𝑑𝑒𝑙 𝑑𝑒𝑙 𝑖 − é𝑠𝑖𝑚𝑜 𝑣𝑖𝑛𝑜 𝑒𝑛 𝑙𝑎 𝑗 − é𝑠𝑖𝑚𝑎 𝑟𝑒𝑔𝑖ó𝑛
𝑥1𝑖 = 𝑠𝑎𝑏𝑜𝑟
𝑥𝑗 = 𝑅𝑒𝑔𝑖𝑎𝑜; j=1,2,3.
𝑗 = 1 (Regiao 1)

𝑗 = 2 (Regiao 2)

𝑗 = 3 (Regiao 3)
𝛼 = 𝑖𝑛𝑡𝑒𝑟𝑐𝑒𝑝𝑡𝑜
𝛽1 = Mide el cambio que el sabor ejerce sobre la calidad media del vino manteniendo todo lo demás constante
𝛾𝑗 = Mide el cambio sobre la calidad media del vinoque ejercido por la región manteniendo todo lo demás constante

Previo a adelantar procesos de interpretación, se procedió a realizar el diagnóstico a

través de envelope, identificación de puntos influyentes y análisis de residuos. Los primeros dos
resultados gráficos se muestran a continuación:
Figura 8.
Envelope, Distancia de Cook y residuales

11
 El envelope muestra un comportamiento adecuado, con un solo punto que se encuentra
parte fuera y parte dentro de la banda de confianza. Así mismo, se identifican dos puntos
posiblemente influyentes, a saber: 12 y 20. Los datos identificados no generan cambios
inferenciales al eliminarlos del modelo, por lo cual se conservan. Por otro lado, al calcular la
función de desvío, se identifica que con un valor p de 0.2904, no se rechaza 𝐻0 , lo cual da cuenta
de que el modelo gamma posee un buen ajuste a los datos. Al adelantar el diagnóstico de los
residuos a través de gamlss, se identificó el siguiente comportamiento:

Figura 9.
Plot Gamlss

Nuevamente, se identifica un patrón de cono en los residuos, lo cual muestra que siguen
existiendo problemas de heterocedasticidad aun cuando el modelo se modifica. Dado que, a
través de los modelos aditivos de localización, escala y forma se puede modelar más de un
parámetro, se procedió a ajustar la fórmula de sigma. Al efectuar dicha modificación se
identificó la siguiente estructura de los residuos:

12
Figura 11.
Plot Galmss

Así mismo, el wormplot del modelo ajustado a través de Gamlss muestra un ajuste
adecuado:
Figura 10.
Wormplot

13
 Considerando que, con los ajustes realizados, el comportamiento del modelo muestra
mejora, se procedió a realizar las interpretaciones correspondientes. En términos de
interpretación, se considera la siguiente salida:

Se identifica que, en cuanto al valor promedio de calidad, respecto al sabor, un aumento

unitario en tal variable se asocia con un incremento del 8.48% en el valor promedio de la calidad,
teniendo las demás variables constantes, dado que 𝑒 0.08148 ≈1.0848. En cuanto a la región, el
hecho de que un vino provenga de la región 2 se asocia con una disminución del 11.67% en el
valor medio de la calidad respecto a los vinos provenientes de la región 1, manteniendo las
demás variables constantes, puesto que 𝑒 −0.12409≈0.8833. Finalmente, los vinos provenientes de

14
la región tres se asocian con un incremento del 8.2% en el valor promedio de calidad respecto a
la región 1, y manteniendo todo lo demás constante, considerando que 𝑒 0.07886 ≈1.082.
Puesto que se pretendió modelar el parámetro de dispersión de los datos, a continuación,
se presenta la postulación del modelo respecto a sigma, donde:
log (𝜎𝑖𝑗 ) = 𝛼 + 𝜏1 𝑥1𝑖 + 𝜏𝑗 𝑥𝑗
Log (𝜇𝑖𝑗 ) = 𝐿𝑜𝑔𝑎𝑟𝑖𝑡𝑚𝑜 𝑑𝑒 𝑙𝑎 𝑑𝑖𝑠𝑝𝑒𝑟𝑠𝑖ó𝑛 𝑑𝑒 𝑙𝑎 𝑐𝑎𝑙𝑖𝑑𝑎𝑑 𝑑𝑒𝑙 𝑖 − 𝑒𝑠𝑖𝑚𝑜 𝑣𝑖𝑛𝑜 𝑑𝑒 𝑙𝑎 𝑗 −
𝑒𝑠𝑖𝑚𝑎 𝑟𝑒𝑔𝑖ó𝑛
𝑥1𝑖 = 𝑠𝑎𝑏𝑜𝑟
𝑥𝑗 = 𝑅𝑒𝑔𝑖𝑎𝑜; j=1,2,3.
𝑗 = 1 (Regiao 1)

𝑗 = 2 (Regiao 2)

𝑗 = 3 (Regiao 3)
𝛼 = 𝑖𝑛𝑡𝑒𝑟𝑐𝑒𝑝𝑡𝑜
𝜏1
= 𝑚𝑖𝑑𝑒 𝑒𝑙 𝑐𝑎𝑚𝑏𝑖𝑜 𝑒𝑛 la dispersión de la calidad del vino por cada unidad adicional de sabor, manteniendo todo lo demás constante

𝜏𝑗
= mide el cambio en la dispersión de la calidad del vino entre cada región, manteniendo las demás variables constantes

Bajo esta consideración, se identificó que solamente sabor presenta una relación
significativa con la variación de Calidad. Los datos muestran que un aumento unitario en sabor
se asocia con una disminución de un 40.77% en la dispersión identificada en la variable.

15
 2. Punto 2
El archivo [Link] describe datos de un estudio prospectivo con 100 personas de al
menos 65 años en buenas condiciones físicas. El objetivo del estudio es tratar de relacionar el
número medio de caídas en un periodo de seis meses con algunas variables explicativas. Los
datos se describen en el siguiente orden: caidas (número de caídas en el período), intervencion
(=0 solo educación, =1 educación y ejercicios físicos), genero (=0 femenino, =1 masculino),
equilibrio (puntuación) y fuerza (puntuación). Para las variables equilibrio y fuerza, a mayor
valor, mayor equilibrio y fuerza del individuo, respectivamente. Inicialmente hacer un análisis
descriptivo de los datos. Ajustar un modelo de Poisson log-lineal que incluye efectos principales
e interacciones hasta el primer orden. Seleccione un submodelo utilizando los criterios AIC.
Interpretar los resultados y realizar un análisis de diagnóstico (Use glm y gamlss). ¿Si el modelo
Poisson seleccionado no es adecuado que modelo usted recomendaría? Realice un análisis
completo de dicho modelo
Análisis descriptivo
En primera instancia, se realizó un análisis gráfico de la distribución de la variable caídas,
variable respuesta del modelo que se desea formular. Se observa una distribución con mayor
frecuencia en el sector izquierdo del gráfico y por ende una cola a derecha. Teniendo en cuenta la
naturaleza de la variable respuesta (cuantitativa discreta) y la distribución de la misma, se sugiere
que un modelo adecuado puede ser un modelo Poisson.

16
Figura 12.
Distribución de la variable caídas

La base Geriatra está constituida por cuatro variables explicativas. Dos de ellas de tipo
cualitativa nominal (intervención y genero) y, las otras dos de tipo cuantitativo (equilibrio y
fuerza). En el siguiente gráfico, se resume la exploración gráfica de la relación entre las variables
que serán tomadas como posibles predictores, y la variable respuesta de número medio de caídas.
La gráfica no sugiere una relación clara entre fuerza y equilibrio de los pacientes con el
número de caídas, sin embargo, se observa que aquellos pacientes con menos fuerza y equilibrio
pueden experimentar un mayor número de caídas. En cuanto a la variable intervención se
observa que aquellos pacientes que solo recibieron educación tienen una mediana superior

17
respecto a aquellos pacientes que recibieron educación y ejercicios, finalmente la diferencia de
numero de caídas clasificada por género es leve.
Figura 13.
Caídas vs las otras variables

Modelo Lineal Propuesto

Considerando que, dadas las características de las variables, se puede pensar en la
pertinencia de postular un modelo Poisson log-lineal, a continuación, se presenta la formulación
del modelo correspondiente:
Distribución variable de respuesta
𝑌𝑖𝑗𝑘 ~𝑖𝑑 𝑃𝑜𝑖𝑠𝑠𝑜𝑛(𝜇𝑖𝑗𝑘 )
Componente sistemático
log (𝜇𝑖𝑗𝑘 ) = 𝛼 + 𝛿𝑗 + 𝛾𝑘 + 𝛽1 + 𝛽2 + 𝜓1 + 𝜓2 + 𝜓3 + 𝜓4 + 𝜓5 + 𝜓6
Donde:

18
 Donde:
𝑌𝑖𝑗𝑘 = 𝑁ú𝑚𝑒𝑟𝑜 𝑑𝑒 𝑐𝑎𝑖𝑑𝑎𝑠 𝑜𝑏𝑠𝑒𝑟𝑣𝑎𝑑𝑜 𝑑𝑒𝑙 𝑖 − é𝑠𝑖𝑚𝑜 𝑝𝑎𝑐𝑖𝑒𝑛𝑡𝑒 𝑐𝑜𝑛 𝑙𝑎 𝑗 −
é𝑠𝑖𝑚𝑎 𝑖𝑛𝑡𝑒𝑟𝑣𝑒𝑛𝑐𝑖ó𝑛 𝑑𝑒𝑙 𝑘 − é𝑠𝑖𝑚𝑜 𝑔𝑒𝑛𝑒𝑟𝑜
𝜇𝑖𝑗 = 𝑁ú𝑚𝑒𝑟𝑜 𝑑𝑒 𝑐𝑎í𝑑𝑎𝑠 𝑚𝑒𝑑𝑖𝑎 𝑑𝑒𝑙 𝑖 − é𝑠𝑖𝑚𝑜 𝑝𝑎𝑐𝑖𝑒𝑛𝑡𝑒 𝑐𝑜𝑛 𝑙𝑎 𝑗 −
é𝑠𝑖𝑚𝑎 𝑖𝑛𝑡𝑒𝑟𝑣𝑒𝑛𝑐𝑖ó𝑛 𝑑𝑒𝑙 𝑘 − é𝑠𝑖𝑚𝑜 𝑔𝑒𝑛𝑒𝑟𝑜.
𝑥𝑗 = intervenciones

𝑗 = 0 (solo educación) casilla de referencia

𝑗 = 1(Educación y ejercicios)
𝑥𝑘 = género
𝑘 = 0 (femenino) casilla de referencia

𝑘 = 1(masculino)
𝑥1𝑖 = 𝑓𝑢𝑒𝑟𝑧𝑎
𝑥1𝑗 = 𝑒𝑞𝑢𝑖𝑙𝑖𝑏𝑟𝑖𝑜
𝛿𝑗 = Mide el cambio que ejerce el tipo de intervención sobre el número de caídas

manteniendo todo lo demás constante

𝛾𝑘 = Mide el cambio que el género tiene sobre el número de caídas manteniendo todo lo

demás constante

𝛽1 = Mide el cambio que la fuerza tiene sobre el número de caídas manteniendo todo lo

demás constante

𝛽2 = Mide el cambio que el equilibrio tiene sobre el número de caídas manteniendo todo lo

demás constante

𝜓1 = Mide el cambio que ejerce la interacción entre intervención y equilibrio sobre el

número de caídas manteniendo todo lo demás constante

𝜓2 = Mide el cambio que ejerce la interacción entre intervención y fuerza sobre el número

de caídas manteniendo todo lo demás constante

19
 𝜓3 = Mide el cambio que ejerce la interacción entre intervención y género sobre el número

de caídas manteniendo todo lo demás constante

𝜓4 = Mide el cambio que ejerce la interacción entre equilibrio y fuerza sobre el número de

caídas manteniendo todo lo demás constante

𝜓5 = Mide el cambio que ejerce la interacción entre equilibrio y genero sobre el número de

caídas manteniendo todo lo demás constante

𝜓6 = Mide el cambio que ejerce la interacción entre fuerza y género sobre el número de

caídas manteniendo todo lo demás constante

𝛼 = 𝑖𝑛𝑡𝑒𝑟𝑐𝑒𝑝𝑡𝑜
A través del AIC, se identificó que el mejor modelo seleccionado fue aquel que contenía
las variables intervención, equilibrio y fuerza como predictores. Se resalta que, todos los
modelos en los cuales se incluyó interacciones reportaron AICs superiores.
Adicionalmente se formularon seis modelos en los cuales se incluyen todas las posibles
variables explicativas y cada una de las posibles interacciones, lo cual arrojó un resultado
semejante al anteriormente mencionado. Por lo tanto, la formulación de este nuevo modelo
obedece, entonces, a la siguiente postulación:
Distribución variable de respuesta
𝑌𝑖𝑗 ~𝑖𝑑 𝑃𝑜𝑖𝑠𝑠𝑜𝑛(𝜇𝑖𝑗 )
Componente sistemático
log (𝜇𝑖𝑗 ) = 𝛼 + 𝛿𝑗 + 𝛽1 + 𝛽2
𝑌𝑖𝑗 = 𝑁ú𝑚𝑒𝑟𝑜 𝑑𝑒 𝑐𝑎𝑖𝑑𝑎𝑠 𝑜𝑏𝑠𝑒𝑟𝑣𝑎𝑑𝑜 𝑑𝑒𝑙 𝑖 − é𝑠𝑖𝑚𝑜 𝑝𝑎𝑐𝑖𝑒𝑛𝑡𝑒 𝑐𝑜𝑛 𝑙𝑎 𝑗
− é𝑠𝑖𝑚𝑎 𝑖𝑛𝑡𝑒𝑟𝑣𝑒𝑛𝑐𝑖ó𝑛
𝜇𝑖𝑗 = 𝑁ú𝑚𝑒𝑟𝑜 𝑑𝑒 𝑐𝑎í𝑑𝑎𝑠 𝑚𝑒𝑑𝑖𝑎 𝑑𝑒𝑙 𝑖 − é𝑠𝑖𝑚𝑜 𝑝𝑎𝑐𝑖𝑒𝑛𝑡𝑒 𝑐𝑜𝑛 𝑙𝑎 𝑗 −
é𝑠𝑖𝑚𝑎 𝑖𝑛𝑡𝑒𝑟𝑣𝑒𝑛𝑐𝑖ó𝑛 .
𝑥𝑗 = intervenciones

𝑗 = 0 (solo educación) casilla de referencia

𝑗 = 1(Educación y ejercicios)
20
 𝑥1𝑖 = 𝑓𝑢𝑒𝑟𝑧𝑎
𝑥1𝑗 = 𝑒𝑞𝑢𝑖𝑙𝑖𝑏𝑟𝑖𝑜
𝛿𝑗 = Mide el cambio que ejerce el tipo de intervención sobre el número de caídas

manteniendo todo lo demás constante

𝛽1 = Mide el cambio que la fuerza tiene sobre el número de caídas manteniendo todo lo

demás constante

𝛽2 = Mide el cambio que el equilibrio tiene sobre el número de caídas manteniendo todo lo

demás constante

𝛼 = 𝑖𝑛𝑡𝑒𝑟𝑐𝑒𝑝𝑡𝑜
Una vez seleccionado el modelo, se realizó el análisis diagnóstico respecto al ajuste del
mismo mediante la identificación del envelope para la distribución Poisson, así como la
identificación de posibles datos influyentes. El resultado gráfico del ajuste envelope y de las
distancias se muestra a continuación:

21
Figura 14.
Envelope, distancia de Cook y Residuos, Modelo Poisson

El ajuste del envelope muestra que el modelo no tiene el mejor comportamiento esto se
evidencia por medio de puntos por fuera de las bandas. Por lo tanto, el modelo es susceptible de
mejorar empleando otra distribución en su formulación, y se identifica que las observaciones más
distantes son la 37, 42 y la 52. No se observan patrones en los residuos estandarizados.
Adicionalmente, se obtuvo el diagnóstico de residuos a través de la función de modelo aditivos
generalizados de localización escala y forma (GAMLSS) en R. En tal proceso, se observaron los
siguientes resultados:

22
Figura 15.
Plot Gamlss

El modelo evidencia que la distribución de los residuos cuantílicos no es completamente

simétrica, lo cual confirma que, su comportamiento no es óptimo. A continuación, se presenta el
Worm plot

23
 Teniendo en cuenta que, el modelo se formuló empleando una distribución de Poisson se
obtienen residuos aleatorizados por tal motivo se obtuvieron cuatro worm plot. En las gráficas
mencionadas se observa la presencia de puntos al interior de las zonas restringidas,
específicamente en el área posicionada en el primer cuadrante.
Posteriormente se retiraron los puntos influyentes previamente identificados y se
determinó si dicha extracción generó cambios inferenciales en el modelo ya formulado. A
continuación, se presenta las salidas obtenidas:

24
2.2.1 Modelo seleccionado

2.2.2 Modelo seleccionado sin puntos influyentes

Una vez retirados los puntos influyentes se observa que la variable fuerza no es
estadísticamente significativa al 5%, no hay cambio en cuanto a futuras interpretaciones.

25
 Interpretación
Manteniendo todo lo demás constante, la presencia de la intervención en educación y
ejercicio produce una disminución del 65.96% en el número promedio de caídas, con respecto a
la intervención basada solo en educación ya que: 𝑒 −1.077770 = 0.3403. Manteniendo las demás
variables constantes, por cada unidad adicional en que se aumenta el equilibrio se produce un
incremento 0.95% en las caídas(𝑒 0.009471 = 1.0095159), si bien lo anterior tiene un sentido
estadístico lo que se esperaría teóricamente es que a mayor equilibrio se reduzca el número de
caídas. Manteniendo todo lo demás constante por cada unidad que se incremente la fuerza el
número de caídas aumenta en un 0.90% (0.008979 = 1.0090), aunque esta conclusión tiene
sentido estadístico, se esperaría que el número de caídas disminuyera a mayor fuerza.

26
 3. Punto 3
(Valor 2,0) El archivo [Link] describe las siguientes variables observadas en una
muestra de 186 descargas de pesca en Bahía de Todos los Santos (costa noreste de Brasil), en el
período de enero de 2013, referentes a la captura de la raya blanca por el método artesanal. Las
características en estudio fueron (i) período (período de pesca, seco o lluvioso), (ii) lugar (lugar
de pesca, área1, área2, área3 y área4), (iii) marea (marea, cuadrado o sicigia), (iv) viento
(velocidad del viento en m /s), (v)tmax (temperatura máxima, en C), (vi) tmin (temperatura
mínima en C), (vii) ins (insolación en horas) y (viii) cpue (captura por unidad de esfuerzo, en
kg). Las variables (iii) a (vii) se observaron en el sitio de pesca. El objetivo principal del estudio
es relacionar la Cpue media con las demás variables explicativas. Realizar inicialmente un
análisis descriptivo de los datos. Proponga un modelo gama con enlace logarítmico con todas las
variables explicativas y por medio del criterio de AKAIKE seleccione un submodelo adecuado.
Evalúe si es necesario introducir interacciones de primer orden a un nivel de significancia de 10
%. Realice un análisis completo, que debe inclur análisis de residuos y análisis de sensibilidad.
Análisis Exploratorio
Inicialmente se realizó un análisis exploratorio univariable y luego se relacionan las
variables independientes con la variable dependiente (Cpue).
Figura 16.
Distribución de periodos

27
 La figura muestra la distribución de períodos de pesca en términos de frecuencia, la barra
roja representa el período lluvioso, con una frecuencia que alcanza aproximadamente el 50.5%,
mientras que la barra azul representa el período seco, con una frecuencia cercana al 49.5%. Esto
indica una distribución casi uniforme entre los días lluviosos y los días secos.
La siguiente figura ilustra la distribución de observaciones en cuatro áreas distintas. En el
área 1, se observa una frecuencia del 29.0%, mientras que el área 2 presenta una proporción
superior, alcanzando el 45.7% de las observaciones. Por otro lado, el área 3 muestra una
frecuencia del 15.1%, y el área 4 registra la menor frecuencia con un 10.2%. Este análisis
permite identificar el área 2 como la de mayor actividad en comparación con las demás, lo que
podría indicar condiciones más favorables para las actividades realizadas en dicho lugar.

Figura 17.
Distribución de lugar de pesca

La siguiente figura, visualiza el histograma de la distribución de la velocidad del viento,

medida en metros por segundo. Se observa que la mayoría de las velocidades registradas se
concentran en el rango de 1 a 3 metros por segundo. Específicamente, la mayor frecuencia se
encuentra entre 2 y 3 metros por segundo, destacando una alta incidencia de vientos de
intensidad moderada. A su vez, las velocidades cercanas a 0 y mayores a 4 metros por segundo
muestran frecuencias significativamente menores, lo que indica que son menos comunes. Este
28
patrón sugiere que el área estudiada generalmente experimenta vientos de baja a moderada
intensidad, con eventos de alta velocidad siendo poco frecuentes.
Figura 18.
Histograma de velocidad del Viento

Figura 19.
Distribución por marea

La primera barra representa la marea tipo cuadrado, que abarca el 55.4% de las
observaciones. En contraste, la marea tipo sicigia es representada por la segunda barra con un
44.6% de las observaciones. Esta distribución sugiere una ligera predominancia de las mareas
cuadradas sobre las sicigias en el área de estudio.

29
Figura 20.
Histograma, Temperatura Máxima

Las temperaturas oscilan principalmente entre 25 y 32.5 grados Celsius. La densidad más
alta se observa alrededor de los 27.5 grados Celsius, lo que sugiere que esta es la temperatura
máxima más frecuente en el área estudiada. La curva de densidad, marcada con una línea azul,
indica un comportamiento bimodal con dos picos pronunciados, uno cerca de los 27.5 grados y
otro alrededor de los 30 grados. Esto puede indicar la presencia de dos condiciones climáticas
distintas dentro del mismo marco temporal o geográfico, con temperaturas que tienden a
agruparse alrededor de estos dos valores.

30
Figura 21.
Histograma Temperatura Mínima

Se visualiza que las temperaturas se agrupan principalmente en torno a dos rangos: cerca
de 21 grados Celsius y alrededor de 22.5 grados Celsius. La densidad alcanza su punto más alto
en estos dos valores, lo que sugiere que estas temperaturas mínimas son las más comunes en el
área estudiada. La presencia de un patrón bimodal puede indicar variaciones en las condiciones
ambientales nocturnas durante el periodo analizado.
Figura 22.
Histograma de insolación en horas

31
 Se observa una distribución variada en la frecuencia de diferentes niveles de intensidad
solar. La intensidad cerca de 7.5 horas representa el valor más frecuente, seguido por los niveles
cerca de 10 y 2.5 horas, respectivamente. La presencia de una frecuencia más baja alrededor de 5
indica una menor ocurrencia de esta intensidad específica.
Figura 23.
Histograma de CPUE

El histograma presentado muestra la distribución de la Captura por Unidad de Esfuerzo

(CPUE), expresada en kilogramos, se observa una densidad más alta en valores bajos de CPUE,
con un pico significativo cerca de 0 kg, que disminuye gradualmente a medida que aumenta la
CPUE. El área sombreada en rosa sobre la curva de densidad indica que la mayoría de las
capturas se concentran en valores bajos, con una caída notable en la frecuencia a medida que la
CPUE aumenta hacia 50 kg y más allá. En este punto, la distribución de CPUE indica que su
distribución posiblemente sea Gamma.
En las siguientes figuras se relaciona las variables independientes su relación con CPUE.
La sucesiva figura, exhibe la distribución de la Captura por Unidad de Esfuerzo (CPUE),
diferenciada por períodos de clima lluvioso y seco, utilizando colores rojo y azul
respectivamente para indicar cada período. Se observa que, durante el período lluvioso,
representado en rojo, la densidad de CPUE es más alta en valores bajos y decrece rápidamente a

32
medida que aumenta la CPUE. Por otro lado, en el período seco, indicado en azul, la densidad de
CPUE también empieza alta en valores bajos, pero mantiene una leve caída comparada con el
período lluvioso.
Esta representación permite deducir que las capturas en períodos lluviosos tienden a ser
más consistentes en valores bajos de esfuerzo, mientras que en períodos secos hay una mayor
dispersión en la eficiencia del esfuerzo de captura.
Figura 24.
Histograma de CPUE vs periodo

33
Figura 25.
Histograma de CPUE por Lugar

En el área 1, la densidad es más alta en valores muy bajos de CPUE, lo que sugiere que
las capturas son generalmente bajas en esta área. El segundo lugar presenta una distribución más
variada, con un pico significativo en un rango intermedio de CPUE y menores frecuencias a
medida que la CPUE aumenta. En el tercer lugar, la distribución comienza con un pico en
valores bajos, seguido de una disminución y un pequeño aumento en la densidad para valores
intermedios de CPUE. Finalmente, en el área 4, muestra una densidad más equilibrada con un
pico principal en un rango de CPUE intermedio y valores más bajos tanto en los extremos
inferiores como superiores.

34
Figura 26.
Relación entre Velocidad del viento, Marea y CPUE

Observando la distribución de los puntos, se aprecia que la mayoría de las capturas

registradas, tanto en marea cuadrado como en marea sicigia, se concentran en velocidades del
viento de 1 a 3 metros por segundo. Sin embargo, no se detecta una correlación directa entre la
velocidad del viento y la CPUE, dado que los valores de CPUE varían ampliamente dentro de
este rango de velocidad.
Además, en velocidades del viento superiores a 2,5 metros por segundo, la cantidad de
datos es menor y las capturas tienden a ser bajas, lo que podría indicar un impacto negativo de
vientos fuertes sobre la eficacia de las actividades pesqueras. La representación gráfica no
muestra diferencias marcadas entre las mareas cuadrado y sicigia en relación con la velocidad del
viento, sugiriendo que el tipo de marea no influye significativamente en el CPUE para las
diferentes velocidades de viento analizadas.

35
Figura 27.
Relación temperatura máxima, período y CPUE

La mayoría de las observaciones de CPUE se distribuyen en un rango de temperatura de

27.5°C a 32.5°C. Se detecta que los valores más altos de CPUE suelen presentarse entre las
temperaturas de 30.0°C a 32.5°C y están predominantemente asociados con el período seco. Esto
sugiere que condiciones más cálidas y secas podrían ser más propicias para las capturas.
La intensidad solar también muestra un patrón interesante: puntos que representan
mayores períodos de insolación tienden a agruparse en temperaturas más altas y están
correlacionados con valores más elevados de CPUE, especialmente durante el período seco.

36
Figura 28.
Temperatura máxima, período, lugar, intensidad solar y C PUE

En cada panel, se observa una variabilidad considerable en la CPUE a lo largo de un

rango de temperaturas máximas entre 25.0°C y 32.5°C. Aunque no se evidencia una correlación
lineal clara entre la temperatura máxima y la CPUE, es notorio que, en los períodos secos,
representados por puntos color “cian”, se tienden a registrar valores más altos de CPUE en
temperaturas más elevadas, especialmente en los lugares 2 y 4. En contraste, durante el período
lluvioso, los valores de CPUE tienden a ser más bajos y menos dispersos a lo largo de las
temperaturas.
La intensidad solar también parece jugar un papel importante, especialmente en los
lugares 2 y 4, donde los puntos con mayor tamaño, que indican una mayor insolación, coinciden
frecuentemente con temperaturas más altas y mayores valores de CPUE. Esta tendencia podría
sugerir que una mayor exposición solar, combinada con temperaturas altas durante períodos
secos, facilita una mayor productividad pesquera
37
Figura 29.
Distribución de CPUE por marea y periodo

En ambas mareas, los períodos secos (en color cian) exhiben una mediana de CPUE
superior a los períodos lluviosos (en color rosa), lo que sugiere que las condiciones secas podrían
ser más favorables para la pesca o que la actividad de los peces aumenta en estas condiciones.
Además, en ambas mareas, los rangos de variabilidad de CPUE, indicados por la longitud de las
cajas y los bigotes, son comparables, pero los valores extremos, representados por los puntos
fuera de los bigotes, indican que existen algunas capturas excepcionalmente altas, especialmente
durante el período seco.

38
3.1.1 Análisis descriptivo numérico
Tabla 1.
Exploración numérica
Tmax Tmax Tmax Tmin Tmin Tmin Ins Ins Ins Cpue Cpue Cpue
Periodo Lugar Marea
Media Min Max Media Min Max Media Min Max Media Min Max
Chivoso área 1 cuadrado 26.6 23.2 30.3 20.8 19.1 23 3.96 0 8.5 7.99 1.19 37.5
Chivoso área 2 sicigia 27 23.8 29.3 20.6 19 22.4 5.75 0.2 8.3 11.4 0.714 27
Chivoso área 2 cuadrado 27.5 23.2 31.2 20.8 18.8 23.5 5.05 0 9.1 20.1 4.69 55
Chivoso área 2 sicigia 27 24.2 29.5 20.6 19.6 22.4 6.44 0.1 8.9 17.2 1.43 78
Chivoso área 3 cuadrado 29.2 28.2 31.2 21.4 20.2 23.5 5.67 2.7 8.8 11.6 1.25 25.2
Chivoso área 3 sicigia 26.8 25.5 28.3 20.4 19 21.5 6.78 5.1 7.6 15.8 5.12 36.1
Chivoso área 4 cuadrado 23.8 23.8 23.8 19.9 19.9 19.9 0 0 0 25 25 25
Chivoso área 4 sicigia 28.6 27.9 29.5 21 19.6 21.5 7.45 7.1 7.7 42.9 9.83 66.7
Seco área 1 cuadrado 30.3 27 32.3 23.1 21.4 24.8 4.43 0.2 9.5 11.8 4.06 38
Seco área 2 sicigia 31.4 29.2 33.7 23.4 22.8 24.3 7.57 2.6 9.5 8.09 0.596 35
Seco área 2 cuadrado 30 27 32.3 22.9 21.9 24.8 5.99 0.2 9.3 18.9 4.26 80
Seco área 2 sicigia 31.4 27.1 33.7 23.1 22 24.3 7.86 1.6 9.5 34.5 3.08 167
Seco área 3 cuadrado 31 30 32.3 23.2 22.5 24.8 6.72 2.3 9.3 13.7 2.92 26.2
Seco área 3 sicigia 30.3 27.1 33.7 22.8 22.3 23.7 5.28 1.6 8.9 23.2 4.33 66.7
Seco área 4 cuadrado 30.2 27.7 31.3 22.9 22.6 23.3 5.33 0.5 9.3 15.3 1.5 32.2
Seco área 4 sicigia 30.6 30.4 30.8 22.5 22.3 22.9 7.1 3.5 9.3 18.3 5.13 26.2

La eficiencia más alta se observa en el área 2 durante la marea de sicigia en el período seco, con un CPUE máximo de 167 kg.
Este valor sugiere condiciones óptimas para la pesca en esta área y período, posiblemente debido a la confluencia de condiciones
climáticas favorables. En contraste, el área 4 durante la marea cuadrada en el período Chuvoso registra un CPUE uniforme de 25 kg.
Además, el área 4 en la marea de sicigia durante el mismo período muestra un incremento en la eficiencia con un CPUE de 66.7 kg.
39
 Al observar los datos, se nota que la marea de sicigia generalmente produce un CPUE
más alto en comparación con la marea cuadrada. Por ejemplo, el área 3 durante la marea de
sicigia en el período Chuvoso alcanza un CPUE de 36.1 kg, superior al registrado durante la
marea cuadrada (25.2 kg), lo que apunta a la sicigia como un factor que potencialmente mejora la
captura de los peces. La variabilidad en la temperatura máxima y mínima no muestra una
correlación directa en los valores de CPUE, aunque la temperatura podría influir en la
distribución de las especies objetivo. Por ejemplo, las temperaturas más altas en el área 3 durante
el período seco no necesariamente coinciden con un aumento proporcional en CPUE, lo que
indica que otras variables podrían tener un rol más determinante.
Propuesta del modelo
Para modelar la relación entre el esfuerzo de pesca (CPUE) y varios factores ambientales
y temporales, el modelo se ajusta asumiendo una distribución Gamma y utilizando una función
de enlace logarítmico. En primera instancia, este modelo incorpora todas las variables
independientes, cuya expresión es la siguiente:

Distribución variable de respuesta

𝑌𝑖𝑗𝑘𝑙 ~𝐺𝑎𝑚𝑚𝑎(𝜇𝑖𝑗𝑘𝑙 ; Φ)
Componente sistemático

Log (𝜇𝑖𝑗𝑘𝑙 ) = 𝛼 + 𝛽𝑗 + 𝛾𝑘 + θ𝑙 + 𝛿1 + 𝛿2 + 𝛿3 + 𝛿4

Donde: 𝑌𝑖𝑗𝑘𝑙 = 𝐶𝑎𝑝𝑡𝑢𝑟𝑎 𝑝𝑜𝑟 𝑢𝑛𝑖𝑑𝑎𝑑 𝑑𝑒 𝑒𝑠𝑓𝑢𝑒𝑟𝑧𝑜 𝑜𝑏𝑠𝑒𝑟𝑣𝑎𝑑𝑜 𝑒𝑛 𝐾𝑔 𝑝𝑎𝑟𝑎 𝑙𝑎 𝑖 −

é𝑠𝑖𝑚𝑎 𝑒𝑚𝑏𝑎𝑟𝑐𝑎𝑐𝑖ó𝑛, 𝑑𝑒𝑙 𝑗 − é𝑠𝑖𝑚𝑜 𝑝𝑒𝑟𝑖𝑜𝑑𝑜, 𝑑𝑒 𝑙𝑎 𝑘 − é𝑠𝑖𝑚𝑎 𝑚𝑎𝑟𝑒𝑎, 𝑒𝑛 𝑒𝑙 𝑙 − é𝑠𝑖𝑚𝑜 𝑙𝑢𝑔𝑎𝑟
𝜇𝑖𝑗𝑘𝑙 = 𝐶𝑎𝑝𝑡𝑢𝑟𝑎 𝑚𝑒𝑑𝑖𝑎 𝑑𝑒𝑙 𝑒𝑠𝑓𝑢𝑒𝑟𝑧𝑜 𝑜𝑏𝑠𝑒𝑟𝑣𝑎𝑑𝑜 𝑒𝑛 𝐾𝑔 𝑝𝑎𝑟𝑎 𝑙𝑎 𝑖
− é𝑠𝑖𝑚𝑎 𝑒𝑚𝑏𝑎𝑟𝑐𝑎𝑐𝑖ó𝑛, 𝑑𝑒𝑙 𝑗 − é𝑠𝑖𝑚𝑜 𝑝𝑒𝑟𝑖𝑜𝑑𝑜, 𝑑𝑒 𝑙𝑎 𝑘
− é𝑠𝑖𝑚𝑎 𝑚𝑎𝑟𝑒𝑎, 𝑒𝑛 𝑒𝑙 𝑙 − é𝑠𝑖𝑚𝑜 𝑙𝑢𝑔𝑎𝑟
𝑥𝑗 = 𝑃𝑒𝑟𝑖𝑜𝑑𝑜
𝑥𝑗 = 𝑃𝑒𝑟𝑖𝑜𝑑𝑜; j=1,2
𝑗 = 1; 𝐶ℎ𝑢𝑣𝑜𝑠𝑜
𝑗 = 2; 𝑆𝑒𝑐𝑜

𝑥𝑘 = 𝑀𝑎𝑟𝑒𝑎
𝑥𝑘 = 𝑀𝑎𝑟𝑒𝑎; k=1,2
𝑘 = 1; 𝐶𝑢𝑎𝑑𝑟𝑎𝑑𝑎
𝑘 = 2; 𝑆𝑖𝑠𝑖𝑔𝑖𝑎

40
 𝑥𝑙 = 𝐿𝑢𝑔𝑎𝑟
𝑥𝑙 = 𝐿𝑢𝑔𝑎𝑟; l=1,2,3,4
𝑙 = 1; Á𝑟𝑒𝑎 1
𝑙 = 2; Á𝑟𝑒𝑎 2
𝑙 = 3; Á𝑟𝑒𝑎 3
𝑙 = 4; Á𝑟𝑒𝑎 4

𝑥1𝑖 = 𝑉𝑖𝑒𝑛𝑡𝑜 𝑑𝑒 𝑙𝑎 𝑖 − é𝑠𝑖𝑚𝑎 𝑒𝑚𝑏𝑎𝑟𝑐𝑎𝑐𝑖ó𝑛

𝑥2𝑖 = 𝑇𝑒𝑚𝑝𝑒𝑟𝑎𝑡𝑢𝑟𝑎 𝑚á𝑥𝑖𝑚𝑎 𝑑𝑒 𝑙𝑎 𝑖 − é𝑠𝑖𝑚𝑎 𝑒𝑚𝑏𝑎𝑟𝑐𝑎𝑐𝑖ó𝑛
𝑥3𝑖 = 𝑇𝑒𝑚𝑝𝑒𝑟𝑎𝑡𝑢𝑟𝑎 𝑚í𝑛𝑖𝑚𝑎 𝑑𝑒 𝑙𝑎 𝑖 − é𝑠𝑖𝑚𝑎 𝑒𝑚𝑏𝑎𝑟𝑐𝑎𝑐𝑖ó𝑛
𝑥4𝑖 = 𝐼𝑛𝑠𝑜𝑙𝑎𝑐𝑖ó𝑛 de la i-ésima embarcación
𝛼 = 𝑖𝑛𝑡𝑒𝑟𝑐𝑒𝑝𝑡𝑜

𝛽𝑗 = Mide el cambio que el periodo ejerce sobre la captura por unidad de esfuerzo manteniendo todo lo demás constante

𝛾𝑘 = Mide el cambio que la marea ejerce sobre la captura por unidad de esfuerzo manteniendo todo lo demás constante

θ𝑙 = Mide el cambio que el lugar ejerce sobre la captura por unidad de esfuerzo manteniendo todo lo demás constante

𝛿1 = Mide el cambio que el viento ejerce sobre la captura por unidad de esfuerzo manteniendo todo lo demás constante
𝛿2
= Mide el cambio que la temperatura máxima ejerce sobre la captura por unidad de esfuerzo manteniendo todo lo demás constante
𝛿3
= Mide el cambio que la temperatura mínima ejerce sobre la captura por unidad de esfuerzo manteniendo todo lo demás constante

𝛿4
= Mide el cambio que la insolación ejerce sobre la captura por unidad de esfuerzo manteniendo todo lo demás constante

41
 Modelo 1
Primero, se propone el modelo lineal generalizado utilizando la función glm() de R. cuyo
resultado es el siguiente:

Se observa que los estimadores de las variables no son significativos, a excepción de la

variable lugar 2 y 4, esto implica que diferentes lugares pueden tener un impacto positivo
significativo en el CPUE. El modelo tiene un AIC de 1416.5.

42
 Modelo 2
Alternativamente, se realizó un modelo aditivo generalizado para posición, escala y
forma mediante la función gamlss (), cuyo resultado es el siguiente:

El intercepto y los coeficientes para las variables como periodo seco, lugar, viento,
marea, temperaturas máxima y mínima, y la variable insolación, se mantienen sin cambios
significativos en términos de su relación con la variable respuesta, CPUE. La variable lugar 2, 3
y 4 sigue siendo significativa, indicando un efecto positivo en el CPUE.
Además, en este modelo se incluye un parámetro adicional, sigma, que controla la
variabilidad de los datos. El intercepto para sigma es significativo, sugiriendo que el modelo
captura una variación en los datos (-0.23977 con un valor p significativo) indica una fuerte
influencia en la dispersión de los datos. El número total de observaciones es de 186, con un total
de 11 grados de libertad para el ajuste del modelo, dejando 175 grados de libertad residuales.

43
 El AIC es de 1415.542, relativamente inferior al AIC del modelo lineal generalizado.

Figura 30.
Plot Gamlss

Los residuos parecen estar distribuidos de manera relativamente uniforme a lo largo de

los rangos de valores ajustados, sin patrones claros que sugieran problemas de
heteroscedasticidad. No se observan tendencias claras o patrones sistemáticos en los residuos
respecto al índice de las observaciones, lo que sugiere que el modelo no deja de captar alguna
estructura sistemática en los datos. Las marcas rojas en la base indican valores extremos, pero
parecen ser pocos. Los puntos siguen de cerca la línea roja, lo que indica que los residuos tienen
una distribución aproximadamente normal. En el Q-Q Plot, hay algunos datos en los extremos
que se desvían ligeramente de la línea, lo que podría señalar la presencia de algunos valores
atípicos.

44
 A pesar de las limitaciones del modelo en cuanto al diagnóstico de residuos y su
sensibilidad, y considerando que no se observan diferencias sustanciales en el valor de AIC, se
ha optado por seleccionar el modelo GLM. A continuación, se presenta el análisis del envelope
para este modelo seleccionado:
Figura 31.
Envelope Glm.

La mayoría de los puntos se sitúan dentro de las líneas del envelope, sin embargo, se observa que
algunos puntos, especialmente en los extremos inferiores de la distribución, se desvían
ligeramente fuera de las líneas del envelope. Estas desviaciones pueden indicar la presencia de
valores atípicos o efectos no capturados por el modelo, lo que podría ser un indicativo de que el
modelo no se ajusta completamente a las condiciones propuestas. Adicionalmente, se presenta
los puntos influyentes mediante la distancia de Cook:

45
Figura 32.
Distancia de Cook

En la gráfica, se destacan varias observaciones que superan el umbral establecido de

influencia, marcadas específicamente con sus números de identificación. Por ejemplo, las
observaciones 4, 5, 171, 45, 64, 1 y 86 exceden el umbral, lo que indica que estas observaciones
influyen en el ajuste del modelo.
Modelo 3
Debido a la presencia de los puntos influyentes, se estima nuevamente el modelo, pero
descartando los cuatro puntos influyentes más representativos. El resultado con esta corrección
es:

46
 El coeficiente asociado al intercepto y a la mayoría de las variables sugieren que sus
efectos sobre la variable respuesta son limitados, dada la ausencia de significancia estadística.
Sin embargo, variables como lugar 2, 3, 4, y tmax muestran una asociación positiva y
estadísticamente significativa con la respuesta, lo que indica que la extracción de los puntos
influyentes tiene un efecto apreciable sobre la captura por unidad de esfuerzo. La variable
periodoseco y marea 2, aunque no alcanza un nivel de significancia muy alta, sugiere una
tendencia podría ser de interés en análisis posteriores. El AIC es de 1346.3, demostrando una
mejora respecto a los anteriores modelos. Cuyo Envelope es el siguiente:

47
Figura 33.
Envelope, modelo considerando puntos influyentes

Se observa que los puntos del estudio exhiben una mayor consistencia en comparación
con otros modelos analizados anteriormente. Esta consistencia sugiere que el modelo es
adecuado para capturar las relaciones subyacentes en los datos. (Suponiendo que estos puntos
influyentes de glm sean similares para Gamlss se realizó la exclusión de los puntos influyentes
en el modelo Gamlss pero su resultado es desfavorable, dado que su AIC es 1415.542)
Modelos alternativos
Adicionalmente se realizó la estimación de un modelo de regresión robusta gamma, es
decir, el modelo ha sido ajustado mediante el método 'Mqle', que se enfoca en proporcionar
estimaciones robustas que minimizan la influencia de observaciones atípicas o anómalas en la
estimación de los parámetros, el resultado se representa en el modelo 4.

48
 Modelo 4

El intercepto y la mayoría de las variables, como viento, marea, temperatura máxima,

temperatura mínima e insolación, muestran coeficientes con valores p no significativos, lo que
indica que, bajo este modelo robusto, estas variables no tienen un efecto estadísticamente
significativo sobre CPUE. Sin embargo, la variable lugar en las áreas estudiadas, exhibe un
coeficiente positivo y significativo, sugiriendo que la localización influye significativamente en
la captura por unidad de esfuerzo. Para este modelo no se encuentra su AIC, lo que no permite la
comparación con los anteriores modelos.
Ahora, mediante la función Step se estiman diversos modelos mediante la asignación o
extracción de alguna variable que no es considerable para el modelo, se utilizó el método Both, el
cual comienza con el modelo inicial (con todas las variables) y en cada paso considera tanto la
adición de una nueva variable que no está en el modelo como la eliminación de una variable que
ya está incluida, seleccionando la opción que mejora el criterio de ajuste del modelo.

49
 Modelo 5
El modelo seleccionado es el siguiente:

Las variables Lugar, Marea y Tmax fueron las variables seleccionadas para el mejor
modelo, con un AIC de 1410.1 El Envelope y la distancia de Cook de los puntos influyentes
fueron:

Figura 34.
Envelope y Distancia de Cook modelo con selección de variables

50
 Se visualiza que aún hay puntos por fuera de las bandas de confianza, y se identifica los
puntos influyentes por la distancia de Cook. Extrayendo los 4 puntos influyentes más
representativos, se obtiene:

Las variables más significativas son lugar 2, 4 y Tmax, aunque el resto de variables son
también, relativamente significantes, el modelo logra un AIC de 1328,2 representando el mejor
ajuste encontrado hasta el momento.
Igualmente se realizó un proceso de selección de variables, pero para Galmss mediante la
función stepGAIC, el resultado fue la selección de las variables, Lugar, Marea, Temperatura
máxima. Pero su AIC es 1409.071, y puesto que no se puede realizar el análisis de sensibilidad o
residuos, este tipo de modelos no son los más adecuado para este caso de estudio.
Adicionalmente, se proponen diversas interacciones con todas las variables, para conocer
si alguna de ellas puede mejorar el modelo estimado, en la siguiente tabla se presenta algunos
modelos alternativos, sus respectivas estimaciones, Envelope y distancia de Cook.

51
Tabla 2.
Modelo con interacciones cuadráticas
Modelo 6
𝐿𝑜𝑔 (𝐶𝑃𝑈𝐸) = 𝛽0 + 𝛽1 𝑃𝑒𝑟𝑖𝑜𝑑𝑜 + 𝛽2 𝐿𝑢𝑔𝑎𝑟 + 𝛽3 𝑀𝑎𝑟𝑒𝑎
+ (𝛽4 𝑉𝑖𝑒𝑛𝑡𝑜 + 𝛽5 𝑇𝑚𝑎𝑥 + 𝛽6 𝑇𝑚𝑖𝑛 + 𝛽7 𝐼𝑛𝑠)2 + 𝜇
AIC 1440.3
Envelope Distancia de Cook

AIC corregido (Sin puntos influyentes) 1381.2

Se efectuaron todas las posibles interacciones de multiplicación entre dos variables, sin
eliminar ninguna. Las interacciones que mejoraron el AIC en comparación con el modelo base
sin interacciones (1416.5), son periodo*lugar, lugar*Tmin, marea*Tmax y marea*tmin, sus
comportamientos se presentan en la siguiente tabla:
Tabla 3.
Modelos con interacciones multiplicativas
Modelo 7
𝐿𝑜𝑔 (𝐶𝑃𝑈𝐸) = 𝛽0 + 𝛽1 𝑃𝑒𝑟𝑖𝑜𝑑𝑜 + 𝛽2 𝐿𝑢𝑔𝑎𝑟 + 𝛽3 𝑀𝑎𝑟𝑒𝑎 + 𝛽4 𝑉𝑖𝑒𝑛𝑡𝑜 + 𝛽5 𝑇𝑚𝑎𝑥 + 𝛽6 𝑇𝑚𝑖𝑛 + 𝛽7 𝐼𝑛𝑠
+ 𝛽8 (𝑝𝑒𝑟𝑖𝑜𝑑𝑜 ∗ 𝑙𝑢𝑔𝑎𝑟) + 𝜇
AIC 1415.3
Envelope Distancia de Cook

52
AIC corregido (Sin puntos influyentes) 1344.2
Modelo 8
𝐿𝑜𝑔 (𝐶𝑃𝑈𝐸) = 𝛽0 + 𝛽1 𝑃𝑒𝑟𝑖𝑜𝑑𝑜 + 𝛽2 𝐿𝑢𝑔𝑎𝑟 + 𝛽3 𝑀𝑎𝑟𝑒𝑎 + 𝛽4 𝑉𝑖𝑒𝑛𝑡𝑜 + 𝛽5 𝑇𝑚𝑎𝑥 + 𝛽6 𝑇𝑚𝑖𝑛 + 𝛽7 𝐼𝑛𝑠
+ 𝛽8 (𝑙𝑢𝑔𝑎𝑟 𝑋 𝑇𝑚𝑖𝑛) + 𝜇
AIC 1415
Envelope Distancia de Cook

AIC corregido (Sin puntos influyentes) 1346.9

Modelo 9

53
𝐿𝑜𝑔 (𝐶𝑃𝑈𝐸) = 𝛽0 + 𝛽1 𝑃𝑒𝑟𝑖𝑜𝑑𝑜 + 𝛽2 𝐿𝑢𝑔𝑎𝑟 + 𝛽3 𝑀𝑎𝑟𝑒𝑎 + 𝛽4 𝑉𝑖𝑒𝑛𝑡𝑜 + 𝛽5 𝑇𝑚𝑎𝑥 + 𝛽6 𝑇𝑚𝑖𝑛 + 𝛽7 𝐼𝑛𝑠
+ 𝛽8 (𝑀𝑎𝑟𝑒𝑎 𝑋 𝑇𝑚𝑎𝑥) + 𝜇
AIC 1414.2
Envelope Distancia de Cook

AIC corregido (Sin puntos influyentes) 1344.5

Modelo 10
𝐿𝑜𝑔 (𝐶𝑃𝑈𝐸) = 𝛽0 + 𝛽1 𝑃𝑒𝑟𝑖𝑜𝑑𝑜 + 𝛽2 𝐿𝑢𝑔𝑎𝑟 + 𝛽3 𝑀𝑎𝑟𝑒𝑎 + 𝛽4 𝑉𝑖𝑒𝑛𝑡𝑜 + 𝛽5 𝑇𝑚𝑎𝑥 + 𝛽6 𝑇𝑚𝑖𝑛 + 𝛽7 𝐼𝑛𝑠
+ 𝛽8 (𝑀𝑎𝑟𝑒𝑎 𝑋 𝑇𝑚𝑖𝑛) + 𝜇
AIC 1414.1
Envelope Distancia de Cook

54
 AIC corregido (Sin puntos influyentes) 1342.8

Adicionalmente, se realizó diferentes transformaciones de las variables, a continuación,

se presenta los modelos más interesantes por su AIC:
Tabla 4.
Modelo con transformaciones logarítmicas
Modelo 11
𝐿𝑜𝑔 (𝐶𝑃𝑈𝐸) = 𝛽0 + 𝛽1 𝑃𝑒𝑟𝑖𝑜𝑑𝑜 + 𝛽2 𝐿𝑢𝑔𝑎𝑟 + 𝛽3 𝑀𝑎𝑟𝑒𝑎 + 𝛽4 𝐿𝑜𝑔 (𝑉𝑖𝑒𝑛𝑡𝑜) + 𝛽5 𝐿𝑜𝑔(𝑇𝑚𝑎𝑥)
+ 𝛽6 𝐿𝑜𝑔 (𝑇𝑚𝑖𝑛) + 𝛽7 𝐼𝑛𝑠 + 𝜇
AIC 1416.5
Envelope Distancia de Cook

AIC corregido (Sin puntos influyentes) 1349.5

Tabla 5.
Modelo con transformaciones cuadráticas
Modelo 12
𝐿𝑜𝑔 (𝐶𝑃𝑈𝐸) = 𝛽0 + 𝛽1 𝑃𝑒𝑟𝑖𝑜𝑑𝑜 + 𝛽2 𝐿𝑢𝑔𝑎𝑟 + 𝛽3 𝑀𝑎𝑟𝑒𝑎 + 𝛽4 (𝑉𝑖𝑒𝑛𝑡𝑜)2 + 𝛽5 (𝑇𝑚𝑎𝑥)2
+ 𝛽6 (𝑇𝑚𝑖𝑛)2 + 𝛽7 (𝐼𝑛𝑠 2 ) + 𝜇
AIC 1416.1
Envelope Distancia de Cook

55
 AIC corregido (Sin puntos influyentes) 1348.7

Tabla 6.
Modelo con transformaciones de raíz cuadrada
Modelo 13
𝐿𝑜𝑔 (𝐶𝑃𝑈𝐸) = 𝛽0 + 𝛽1 𝑃𝑒𝑟𝑖𝑜𝑑𝑜 + 𝛽2 𝐿𝑢𝑔𝑎𝑟 + 𝛽3 𝑀𝑎𝑟𝑒𝑎 + 𝛽4 √𝑉𝑖𝑒𝑛𝑡𝑜 + 𝛽5 √𝑇𝑚𝑎𝑥 + 𝛽6 √𝑇𝑚𝑖𝑛
+ 𝛽7 √𝐼𝑛𝑠 + 𝜇
AIC 1415.963
Envelope Distancia de Cook
No disponible No disponible
AIC corregido (Sin puntos influyentes) No disponible

Se aplicaron diversas transformaciones adicionales a las variables, incluyendo

transformaciones recíprocas, que resultaron en un AIC de 1415.422; transformaciones con
funciones trigonométricas, con un AIC de 1433.466; y transformaciones polinómicas de segundo
orden, con un AIC de 1416.518. Sin embargo, estos modelos presentaron limitaciones, no fue
posible visualizar el Envelope o la distancia de Cook, ni corregir puntos influyentes. Aunque
estos modelos pueden ser relativamente mejores que el modelo base, enfrentan deficiencias en la
corrección de los puntos influyentes y en la interpretación de sus coeficientes.

56
 A manera de conclusión, se presenta la siguiente figura de comparación de AIC de los
modelos con mejor eficiencia.

Figura 35.
Comparación de AIC para los mejores modelos

En general, los modelos que integran correcciones (sin puntos influyentes) y seleccionan
variables o con alguna interacción muestran mejoras en el criterio de información de Akaike, su
diferencia es relativamente baja. El modelo que exhibe el AIC más bajo es el correspondiente al
GLM con selección de variables y corregido (modelo 5), con un valor de 1328.2

El modelo 5 constituye la estructura seleccionada para el análisis debido al cumplimiento

de los supuestos en el diagnóstico, formalmente se presenta de la siguiente manera:
Distribución variable de respuesta
𝑌𝑖𝑗𝑘𝑙 ~𝐺𝑎𝑚𝑚𝑎(𝜇𝑖𝑗𝑘𝑙 ; Φ)
Componente sistemático

Log (𝜇𝑖𝑘𝑙 ) = 𝛼 + 𝛾𝑘 + θ𝑙 + 𝛿1 + 𝛿2 + 𝜓1

57
 Donde: 𝑌𝑖𝑘𝑙 = 𝐶𝑎𝑝𝑡𝑢𝑟𝑎 𝑝𝑜𝑟 𝑢𝑛𝑖𝑑𝑎𝑑 𝑑𝑒 𝑒𝑠𝑓𝑢𝑒𝑟𝑧𝑜 𝑜𝑏𝑠𝑒𝑟𝑣𝑎𝑑𝑜 𝑒𝑛 𝐾𝑔 𝑝𝑎𝑟𝑎 𝑙𝑎 𝑖 −
é𝑠𝑖𝑚𝑎 𝑒𝑚𝑏𝑎𝑟𝑐𝑎𝑐𝑖ó𝑛 − é𝑠𝑖𝑚𝑜 𝑝𝑒𝑟𝑖𝑜𝑑𝑜 𝑑𝑒 𝑙𝑎 𝑘 − é𝑠𝑖𝑚𝑎 𝑚𝑎𝑟𝑒𝑎, 𝑒𝑛 𝑒𝑙 𝑙 − é𝑠𝑖𝑚𝑜 𝑙𝑢𝑔𝑎𝑟
𝜇𝑖𝑘𝑙 = 𝐶𝑎𝑝𝑡𝑢𝑟𝑎 𝑚𝑒𝑑𝑖𝑎 𝑑𝑒𝑙 𝑒𝑠𝑓𝑢𝑒𝑟𝑧𝑜 𝑜𝑏𝑠𝑒𝑟𝑣𝑎𝑑𝑜 𝑒𝑛 𝐾𝑔 𝑝𝑎𝑟𝑎 𝑙𝑎 𝑖
− é𝑠𝑖𝑚𝑎 𝑒𝑚𝑏𝑎𝑟𝑐𝑎𝑐𝑖ó𝑛 𝑑𝑒 𝑙𝑎 𝑘 − é𝑠𝑖𝑚𝑎 𝑚𝑎𝑟𝑒𝑎, 𝑒𝑛 𝑒𝑙 𝑙 − é𝑠𝑖𝑚𝑜 𝑙𝑢𝑔𝑎𝑟

𝑥𝑘 = 𝑀𝑎𝑟𝑒𝑎
𝑥𝑘 = 𝑀𝑎𝑟𝑒𝑎; k=1,2
𝑘 = 1; 𝐶𝑢𝑎𝑑𝑟𝑎𝑑𝑎
𝑘 = 2; 𝑆𝑖𝑠𝑖𝑔𝑖𝑎

𝑥𝑙 = 𝐿𝑢𝑔𝑎𝑟
𝑥𝑙 = 𝐿𝑢𝑔𝑎𝑟; l=1,2,3,4
𝑙 = 1; Á𝑟𝑒𝑎 1
𝑙 = 2; Á𝑟𝑒𝑎 2
𝑙 = 3; Á𝑟𝑒𝑎 3
𝑙 = 4; Á𝑟𝑒𝑎 4

𝑥1𝑖 = 𝑉𝑖𝑒𝑛𝑡𝑜 𝑑𝑒 𝑙𝑎 𝑖 − é𝑠𝑖𝑚𝑎 𝑒𝑚𝑏𝑎𝑟𝑐𝑎𝑐𝑖ó𝑛

𝑥2𝑖 = 𝑇𝑒𝑚𝑝𝑒𝑟𝑎𝑡𝑢𝑟𝑎 𝑚á𝑥𝑖𝑚𝑎 𝑑𝑒 𝑙𝑎 𝑖 − é𝑠𝑖𝑚𝑎 𝑒𝑚𝑏𝑎𝑟𝑐𝑎𝑐𝑖ó𝑛
𝑥3𝑖 = 𝑇𝑒𝑚𝑝𝑒𝑟𝑎𝑡𝑢𝑟𝑎 𝑚í𝑛𝑖𝑚𝑎 𝑑𝑒 𝑙𝑎 𝑖 − é𝑠𝑖𝑚𝑎 𝑒𝑚𝑏𝑎𝑟𝑐𝑎𝑐𝑖ó𝑛
𝑥4𝑖 = 𝐼𝑛𝑠𝑜𝑙𝑎𝑐𝑖ó𝑛 de la i-ésima embarcación
𝜓1 = 𝑥𝑙 ∗ 𝑥3𝑖 = Interacción entre lugar y temperatura mínima
𝛼 = 𝑖𝑛𝑡𝑒𝑟𝑐𝑒𝑝𝑡𝑜
𝛾𝑘 = Mide el cambio que la marea ejerce sobre la captura por unidad de esfuerzo manteniendo todo lo demás constante

θ𝑙 = Mide el cambio que el lugar ejerce sobre la captura por unidad de esfuerzo manteniendo todo lo demás constante

𝛿1 = Mide el cambio que el viento ejerce sobre la captura por unidad de esfuerzo manteniendo todo lo demás constante
𝛿2
= Mide el cambio que la temperatura máxima ejerce sobre la captura por unidad de esfuerzo manteniendo todo lo demás constante

𝜓1
= Mide el cambio que la interacción entre lugar y temperatura mínima ejerce sobre la captura por unidad de esfuerzo

manteniendo todo lo demás constante

Interpretación
Para este modelo se tuvieron en cuenta los valores y el antilogaritmo y se procede a la
realización del análisis de cada uno de la siguiente manera

58
 El intercepto del modelo, con un coeficiente de 0.36587, no es significativo, lo que
sugiere que, en ausencia de los efectos de las demás variables, la tasa esperada de cpue no difiere
significativamente de uno. El coeficiente correspondiente a la variable lugar2, con un valor de
0.65963, se asocia con una tasa de cpue aproximadamente 1.9345 veces mayor en comparación
con el lugar de referencia. Esto indica un aumento significativo en la tasa de cpue al cambiar al
lugar2.
Para la variable lugar3, el coeficiente de 0.34485 implica un efecto marginal que
incrementa la tasa de cpue en aproximadamente 1.4117 veces (𝑒 0.34486 = 1.4117). Aunque
marginalmente significativo, este resultado sugiere una tendencia al aumento en la tasa de cpue.
En el caso de lugar4, el coeficiente de 0.64140 muestra un aumento en la tasa de cpue a
aproximadamente 1.8992 veces, reflejando una influencia considerable (𝑒 0.64140 = 1.8992).
La variable marea2, con un coeficiente de 0.21248, presenta un efecto marginal que
incrementa la tasa de cpue en aproximadamente 1.2368 veces (𝑒 0.21248 = 1.2368). Aunque este
efecto es marginalmente significativo, indica una relación positiva entre marea2 y la tasa de
cpue. Por último, el coeficiente para tmax, que es 0.06374, se traduce en un efecto marginal que
incrementa la tasa de cpue en aproximadamente 1.0658 veces por cada unidad adicional de
temperatura máxima(𝑒 0.06374 = 1.0658).
Finalmente, se muestra el Envelope del modelo con el AIC más bajo. No obstante, los
Envelopes de los modelos que han sido ajustados y que incluyen interacciones o
transformaciones exhiben comportamientos similares. Esto indica que, aunque se logra una
mejora en el ajuste del modelo con el AIC más bajo, las características generales de los residuos
no varían significativamente entre los diferentes modelos ajustados.

59
Figura 36.
Envelope, Modelo con mejor AIC.

60
 4. Apéndice. Códigos Utilizados.
###### PUNTO_1
# Cargamos la base de datos
wine <- [Link]("D:/3. Academia/Estadistica aplicada/3. Tercer Semestre/Modelos Lineales
generalizados/Evaluacion 1/[Link]", sep="")
View(wine)
attach(wine)
# variable región como factor

wine$regiao<- [Link](wine$regiao)
str(wine)

######################3
###### variable explicativa - qualidade
names(wine)
require(ggplot2)
library(plotly)

# Análisis Cullen & Frey

[Link]("fitdistrplus")
library(fitdistrplus)

# Test de distribución de Cullen & Frey

descdist(qualidade, boot = 1000)

png("Cullen_Frey_Plot.png")
descdist(qualidade, boot = 1000)

g0<-ggplot(wine,aes(x=qualidade))+geom_density()+theme_classic()+theme([Link] = element_text(hjust =
0.5))+ggtitle("Distribución de la variable qualidade") #3 está normal, pues, qué brutos, no, pero si
ggplotly(g0)
# Qualidade vs claridade

g1<-ggplot(wine,aes(x=claridade,y=qualidade))+
geom_point(size=2.8,colour="darkred")+
scale_shape_manual(values=c(18))+ylab("Calidad del vino")+xlab("Claridad del vino")+theme_classic()+
theme([Link] = element_text(hjust = 0.5))+
ggtitle("Claridad y calidad del vino")+stat_smooth(method = loess,se=FALSE);g1
ggplotly(g1)

# Qualidade vs aroma

g2<-ggplot(wine,aes(x=aroma,y=qualidade))+
geom_point(size=2.8,colour="darkred")+
scale_shape_manual(values=c(18))+ylab("Calidad del vino")+xlab("Aroma del vino")+theme_classic()+
theme([Link] = element_text(hjust = 0.5))+
ggtitle("Aroma y calidad del vino")+stat_smooth(method = loess,se=FALSE);g2

ggplotly(g2)
# Qualidade vs corpo

g3<-ggplot(wine,aes(x=corpo,y=qualidade))+
geom_point(size=2.8,colour="darkred")+
scale_shape_manual(values=c(18))+ylab("Calidad del vino")+xlab("cuerpo del vino")+theme_classic()+

61
 theme([Link] = element_text(hjust = 0.5))+
ggtitle("Calidad y cuerpo del vino")+stat_smooth(method = loess,se=FALSE);g3

ggplotly(g3)

# Qualidade vs sabor

g4<-ggplot(wine,aes(x=sabor,y=qualidade))+
geom_point(size=2.8,colour="darkred")+
scale_shape_manual(values=c(18))+ylab("Calidad del vino")+xlab("Sabor del vino")+theme_classic()+
theme([Link] = element_text(hjust = 0.5))+
ggtitle("Calidad y sabor del vino")+stat_smooth(method = loess,se=FALSE);g4

ggplotly(g4)

# Qualidade vs aromac

g5<-ggplot(wine,aes(x=claridade,y=qualidade))+
geom_point(size=2.8,colour="darkred")+
scale_shape_manual(values=c(18))+ylab("Calidad del vino")+xlab("AromaC del vino")+theme_classic()+
theme([Link] = element_text(hjust = 0.5))+
ggtitle("Calidad y AromaC del vino")+stat_smooth(method = loess,se=FALSE);g5

ggplotly(g5)

# Qualidade vs regiao

g6 <- ggplot(wine, aes(y = qualidade, x = regiao, fill = regiao)) +

geom_boxplot([Link] = "black") +
ylab("Calidad del vino") +
xlab("Región de procedencia") +
scale_fill_manual(values = c("darkred", "pink", "lightgoldenrodyellow")) +
theme_classic() +
theme([Link] = element_text(hjust = 0.5, size = 14)) + # Ajuste del título
ggtitle("Región y Calidad del vino");g6

ggplotly(g6)

### Acoplar

library(gridExtra)
library(grid)
[Link](g1, g2, g3, g4, g5, g6, ncol = 3, top = textGrob("Relación gráfica de variables y calidad del vino", gp
= gpar(fontsize = 20, font = 3)))

####################################
###################################
## correlación
###################################

[Link]("corrr")
library(corrr)
require(reshape2)
require(car)
require(emmeans)

62
## relación entre variables numércias
[Link](qualidade) # normal
[Link](claridade) # no normal
[Link](aroma) # normal
[Link](aromac) # normal
[Link](corpo) # normal
[Link](sabor) # normal

# Dado que no todas las variables son normales, se procede con correlaciones de pearson y spearman
(cor_matrix1 <- cor(wine[, sapply(wine, [Link])], method = "spearman"))## solo para claridade
(cor_matrix2 <- cor(wine[, sapply(wine, [Link])], method = "pearson"))## para demás variables

### Se identifican correlaciones de magnitud considerable entre aroma y calidad (r=0.71),

### cuerpo y calidad (r=0.55), sabor y calidad (0.79). Se procede con un gráfico de calor:

library(reshape2)

## gráfico para pearson

melted_cor_matrix <- melt(cor_matrix2)

gcor<-ggplot(data = melted_cor_matrix, aes(x = Var1, y = Var2, fill = value)) +

geom_tile(color = "white") +
scale_fill_gradient2(low = "blue", high = "red", mid = "white",
midpoint = 0, limit = c(-1,1), space = "Lab",
name="Correlación") +
theme_minimal() +
theme([Link].x = element_text(angle = 45, vjust = 1,
size = 12, hjust = 1)) +
coord_fixed() +
theme([Link] = element_text(hjust = 0.5, size = 14))+
geom_text(aes(Var1, Var2, label = round(value, 2)), color = "black", size = 4) +
labs(title = "Matriz de Correlación de Pearson",
x = "Variables",
y = "Variables");gcor
# gráfico interactivo

ggplotly(gcor)

## gráfico para spearman

melted_cor_matrix <- melt(cor_matrix1)

gcor<-ggplot(data = melted_cor_matrix, aes(x = Var1, y = Var2, fill = value)) +

geom_tile(color = "white") +
scale_fill_gradient2(low = "blue", high = "red", mid = "white",
midpoint = 0, limit = c(-1,1), space = "Lab",
name="Correlación") +
theme_minimal() +
theme([Link].x = element_text(angle = 45, vjust = 1,
size = 12, hjust = 1)) +
coord_fixed() +
theme([Link] = element_text(hjust = 0.5, size = 14))+
geom_text(aes(Var1, Var2, label = round(value, 2)), color = "black", size = 4) +
labs(title = "Matriz de Correlación de Spearman",
x = "Variables",
y = "Variables");gcor

63
# gráfico interactivo

ggplotly(gcor)

## relación var numérica vs categórica

a<-aov(qualidade~regiao, data = wine)

summary(a) # con un valor P menor a 0.05, existe evidencia para rechazar
# H0, considerando que existen diferencias
# en términos de la calidad del producto dada la región de procedencia.

# comprobación de supuestos:

[Link](residuals(a))# superior a 0.05, normal

leveneTest(a) #superior a 0.05, homogeneidad de varianzas

emmeans(a, pairwise ~ regiao, adjust = "tukey") ## las diferencias en las comparaciones dos a dos
# se identifican entre la región 1 y 2 (p=0.02), y entre región 1 y 3 (p<0.0001), así como
# entre las regiones 2 y tres (p<0.0001). Es decir, cada región imprime un sello diferencia en la calidad
# del vino.

## gráfico

library(DescTools)
library(randtests)
library(car)
library(easyanova)
library(agricolae)
library(lmtest)
library(performance)
library(see)
library(patchwork)
##

test_lsd1=[Link](a,"regiao",group=T)
test_lsd1
plot(test_lsd1)## el gráfico deja ver que el vino de mejor calidad se produce en la región 3,

# seguido de la región 1, y finalmente el vino de menor calidad proviene de la región 2.

test_tukey=TukeyHSD(a,"regiao")
test_tukey
plot(test_tukey)# la misma información en otra presentación.

######### DESARROLLO DEL PUNTO

# 1. hacer un modelo normal linear

modelo1=lm(qualidade~., data=wine)
summary(modelo1)

## el modelo lineal normal aporta en la explicación del fenómeno en un 79%, es significativo

## (p-value: 3.295e-10), pero solo las variables región (en su valor 2) y sabor son estadísticamente
## significativos. Se procede a identificar si el quitar variables influye en la adecuación
## del modelo.

64
## 2. Seleccionar un modelo con stepAIC

require(MASS)

stepAIC(modelo1)

## el mejor modelo parece ser el que solo incluye región, aromac y [Link] ajusta tal modelo.

modelo2<-lm(qualidade~aromac+sabor+regiao, data=wine)
summary(modelo2)
AIC(modelo2)

## El modelo mejora tanto en el valor de R2, como en la significacia de los predictores.

## La variable aromaC no es significativa. Se procede con el diagnóstico.

# 3. diagnóstico
# 3.1. Envelope y Dx normal.

[Link]<-modelo2 ## el envelope del modelo se comporta de manera adecuada.

## el Dx - hay 2 datos posiblemente influyentes, 25, 12 y 20.

[Link]()

# 3.2 DX con gamlss

require(gamlss)

modelo2gamlss=gamlss(qualidade~aromac+sabor+regiao,
family=NO,data = wine)
newpar<-par(mfrow=c(2,2),mar=par("mar")+c(0,1,0,0),[Link]="black",
col="blue4", [Link]="black",[Link]="black",pch="+",
cex=.6, [Link]=1.2, [Link]=1, [Link]=1.2,bg="white")
plot(modelo2gamlss,par=newpar)

[Link]()

## hay un patrón de cono en los residuos ajustados. Así mismo,

# parece existir un leve patrón en residuos versus index. Dados tales hallazgos, se propone un modelo gamma.

#### 4. Modelo GAMMA

# se inicia nuevamente con todas las variables

### Envelope GAMA

[Link] <- function([Link]) {
par(mfrow=c(1,1))

X <- [Link]([Link])
n <- nrow(X)
p <- ncol(X)
w <- [Link]$weights
W <- diag(w)
H <- solve(t(X) %*% W %*% X)
H <- sqrt(W) %*% X %*% H %*% t(X) %*% sqrt(W)
h <- diag(H)
ro <- resid([Link], type = "response")

65
 fi <- (n - p) / sum((ro / (fitted([Link])))^2)
td <- resid([Link], type = "deviance") * sqrt(fi / (1 - h))

e <- matrix(0, n, 100)

for (i in 1:100) {
resp <- rgamma(n, fi)
resp <- (fitted([Link]) / fi) * resp
fit <- glm(resp ~ X - 1, family = Gamma(link = "log"))
w <- fit$weights
W <- diag(w)
H <- solve(t(X) %*% W %*% X)
H <- sqrt(W) %*% X %*% H %*% t(X) %*% sqrt(W)
h <- diag(H)
ro <- resid(fit, type = "response")
phi <- (n - p) / sum((ro / (fitted(fit)))^2)
e[,i] <- sort(resid(fit, type = "deviance") * sqrt(phi / (1 - h)))
}

e1 <- numeric(n)
e2 <- numeric(n)

for (i in 1:n) {
eo <- sort(e[i, ])
e1[i] <- (eo[2] + eo[3]) / 2
e2[i] <- (eo[97] + eo[98]) / 2
}

med <- apply(e, 1, mean)

faixa <- range(td, e1, e2)
par(pty = "s")
qqnorm(td, xlab = "Percentil da N(0,1)", ylab = "Componente do Desvio", ylim = faixa, pch = 16, main = "",
cex = 0.5)
par(new = TRUE)
qqnorm(e1, axes = F, xlab = "", ylab = "", type = "l", ylim = faixa, lty = 1, main = "")
par(new = TRUE)
qqnorm(e2, axes = F, xlab = "", ylab = "", type = "l", ylim = faixa, lty = 1, main = "")
par(new = TRUE)
qqnorm(med, axes = F, xlab = "", ylab = "", type = "l", ylim = faixa, lty = 2, main = "")
}

#------------------------------------------------------------#
########## COOK Distancia ##########

[Link] <- function([Link]) {

cdi <- [Link]([Link])
plot(cdi, main = "Distancia de Cook", ylab = "Cooks Distance", xlab = "Observaciones")

threshold <- 4 / length(cdi)

selected <- which(cdi > threshold)

text(selected, cdi[selected], labels = selected, pos = 4, col = "black")

}

modelogamma=glm(qualidade~.,family=Gamma(link = log),
data = wine)

66
stepAIC(modelogamma) # el modelo solo debe tener sabor y región

modelogamma1<-glm(qualidade~sabor+regiao,family=Gamma(link = log),
data = wine)
summary(modelogamma1)

### 4.1 DX y envelope

[Link](modelogamma1)
[Link](modelogamma1)

###
[Link](modelogamma1) # valor estimado de 172.22560, SE de 39.47308
0.22085*172.22560 # = 38.03602

pchisq(38.04,34,[Link] = F) # No rechazo H0, es un indicio de que el modelo se ajusta a los datos.

[Link]<-modelogamma1 # el envelope tiene un mejor ajuste.

# Las observaciones distantes son la 12 y 20

modelogamma1<-glm(qualidade~sabor+regiao,family=Gamma(link = log),
data = wine, subset=c(-20,-12))
summary(modelogamma1)

## sin embargo, las observaciones no generan cambios inferenciales en los datos

# 4.2 Dx A través de gamlss

modelogammagamlss <- gamlss(qualidade ~ sabor

+ regiao, family = GA([Link] = "log"),
data = wine)
summary(modelogammagamlss)

newpar<-par(mfrow=c(2,2),mar=par("mar")+c(0,1,0,0),[Link]="black",
col="blue4", [Link]="black",[Link]="black",pch="+",
cex=.6, [Link]=1.2, [Link]=1, [Link]=1.2,bg="white")
plot(modelogammagamlss,par=newpar)
## el patrón de cono se sigue presentando. El modelo sigue con problemas en los residuos

#### Por tal razón, se modela sigma

modelogammagamlss2 <- gamlss(qualidade ~ sabor

+ regiao, [Link] = qualidade ~ sabor
+ regiao,
family = GA([Link] = "log",[Link] = "log"),
data = wine)

newpar<-par(mfrow=c(2,2),mar=par("mar")+c(0,1,0,0),[Link]="blue4",
col="blue4", [Link]="blue4",[Link]="blue4",pch="+",
cex=.45, [Link]=1.2, [Link]=1, [Link]=1.2,bg="white")
plot(modelogammagamlss2,par=newpar)

gamlss::wp(modelogammagamlss2)

67
summary(modelogammagamlss2) ## modelo para interpretación

###### PUNTO_2

geriatra <- [Link]("D:/3. Academia/Estadistica aplicada/3. Tercer Semestre/Modelos Lineales

generalizados/Evaluacion 1/[Link]", quote="\"", [Link]="")
View(geriatra)
str(geriatra)

names(geriatra)=c("caidas","intervencion","genero","equilibrio","fuerza")
head(geriatra)

attach(geriatra)

geriatra$intervencion<-factor(geriatra$intervencion)
geriatra$genero<-factor(geriatra$genero)

require(ggplot2)
### Análisis descriptivo
ggplot(geriatra, aes(x=caidas))+
geom_bar(binwidth=0.5, aes(fill=..count..), col='black') +
labs(title = "Distribución de la variable caidas",
x = "Número de caídas",
y = "Frecuencia") +
theme_bw()

p1<- ggplot(data = geriatra, aes(x = intervencion, y = caidas)) +

# Cambiar el color del gráfico de caja según la variable 'intervencion'
geom_boxplot(aes(fill = intervencion), color = "black", width = 0.7) +
# Cambiar los colores de relleno utilizando una paleta de RColorBrewer
scale_fill_brewer(palette = "Set2") +
labs(title = "Boxplot por Intervención",
x = "Intervención",
y = "Número de caídas") +
theme_bw()

p2<- ggplot(data = geriatra, aes(x = genero, y = caidas)) +

# Cambiar el color del gráfico de caja según la variable 'intervencion'
geom_boxplot(aes(fill = genero), color = "black", width = 0.7) +
# Cambiar los colores de relleno utilizando una paleta de RColorBrewer
scale_fill_brewer(palette = "Set3") +
labs(title = "Boxplot por Genero",
x = "Genero",
y = "Número de caídas") +
theme_bw()

p3 <- ggplot(geriatra, aes(x = equilibrio, y = caidas, color = equilibrio, size = caidas)) +

geom_point(alpha = 0.6) +
labs(title = "Número de caídas vs Equilibrio") +
theme_bw()

p4 <- ggplot(geriatra, aes(x = fuerza, y = caidas, color = fuerza, size = caidas)) +

geom_point(alpha = 0.6) +
labs(title = "Número de caídas vs Fuerza") +
theme_bw()

68
require(patchwork)
# Combinar los gráficos en una sola ventana
combined_plot <- (p1 | p2) / (p3 | p4)
combined_plot + plot_annotation(title = "Caidas vs las otras variables")

?fitDist
[Link]("fitdistrplus")
library(fitdistrplus)

### modelo poisson

modelogeriatraglm<-glm(caidas~intervencion+genero+equilibrio+fuerza,
family=poisson(link=log),
data=geriatra)
summary(modelogeriatraglm)

#### modelo con interacciones

modelogeriatraint1<-glm(formula = caidas ~ intervencion + equilibrio + fuerza +genero+
intervencion*equilibrio+
intervencion*fuerza + intervencion*genero+ equilibrio*fuerza+
equilibrio*genero+fuerza*genero ,
family = poisson(link = log),
data = geriatra)

### modelos con interacciones uno por uno

m1<-glm(formula = caidas ~ intervencion + equilibrio + fuerza +genero+ intervencion*equilibrio ,
family = poisson(link = log),
data = geriatra)
m2<-glm(formula = caidas ~ intervencion + equilibrio + fuerza +genero+ intervencion*fuerza ,
family = poisson(link = log),
data = geriatra)
m3<-glm(formula = caidas ~ intervencion + equilibrio + fuerza +genero+ intervencion*genero ,
family = poisson(link = log),
data = geriatra)
m4<-glm(formula = caidas ~ intervencion + equilibrio + fuerza +genero+ equilibrio*fuerza ,
family = poisson(link = log),
data = geriatra)
m5<-glm(formula = caidas ~ intervencion + equilibrio + fuerza +genero+ equilibrio*genero ,
family = poisson(link = log),
data = geriatra)
m6<-glm(formula = caidas ~ intervencion + equilibrio + fuerza +genero+ fuerza*genero ,
family = poisson(link = log),
data = geriatra)
mp<-glm(formula = caidas ~ intervencion + equilibrio + fuerza +genero ,
family = poisson(link = log),
data = geriatra)
modsel<-glm(formula = caidas ~ intervencion + equilibrio + fuerza,
family = poisson(link = log),
data = geriatra)
AIC(modsel,mp,m1,m2,m3,m4,m5,m6)
summary(modelogeriatraint1)
require(MASS)
stepAIC(modelogeriatraint1)
[Link]<-modsel

### modelo poisson gamlss

69
modelogeriatragamlss<-gamlss(caidas~intervencion+equilibrio+fuerza,
family=PO,
data = geriatra)
summary(modelogeriatragamlss)
newpar<-par(mfrow=c(2,2),mar=par("mar")+c(0,1,0,0),[Link]="blue4",
col="blue4", [Link]="blue4",[Link]="blue4",pch="+",
cex=.45, [Link]=1.2, [Link]=1, [Link]=1.2,bg="white")

plot(modelogeriatragamlss,par=newpar)

###warm plot
par(mfrow = c(2, 2), mar = c(4, 3.8, 1, 2), bg = "white")
[Link](modelogeriatragamlss, type = "wp", howmany = 4)

?[Link]
[Link](modelogeriatragamlss, what="mu", pages = 1, ask = FALSE, rug = TRUE)

####### Puntos influyentes

modsinpi<-glm(formula = caidas ~ intervencion + equilibrio + fuerza,
family = poisson(link = log),
data = geriatra,
subset = c(-37,-42,-52))
summary(modsel)
summary(modsinpi)

####### modelo nbin

[Link]<-[Link](caidas~intervencion+genero+equilibrio+fuerza,
data=geriatra)

summary([Link])
stepAIC([Link])

[Link]<-[Link](formula = caidas ~ intervencion + equilibrio + fuerza +genero+

intervencion*equilibrio+
intervencion*fuerza + intervencion*genero+ equilibrio*fuerza+
equilibrio*genero+fuerza*genero ,
data = geriatra)
[Link]<-stepAIC([Link])
[Link]<-[Link]
#########################################################################################

modselgamlssnb<-gamlss(caidas~intervencion+genero+equilibrio+fuerza,
family= NBI,
data = geriatra)
summary(modselgamlssnb)

###### PUNTO_3

rm(list=ls())

library(ggplot2)
library(dplyr)
library(gamlss)
library(leaps)
require(jtools)
require(interactions)

70
library(reshape2)
library(lme4)

# Base Raia
# 186 descargas de pesca en Bahía de Todos los Santos (costa noreste de Brasil), en el período de enero de 2013,
referentes a la captura de la raya blanca por el método artesanal.
#(i) período (período de pesca, seco o lluvioso),
# (ii) lugar (lugar de pesca, área1, área2, área3 y área4),
# (iii) marea (marea, cuadrado o sicigia)
# (iv) viento (velocidad del viento en m /s)
# (v)tmax (temperatura máxima, en C)
# (vi) tmin (temperatura mínima en C)
# (vii) ins (insolación en horas)
# (viii) cpue (captura por unidad de esfuerzo, en kg)

# las variables (iii) a (vii) se observaron en el sitio de pesca.

# El objetivo principal del estudio es relacionar la cpue (viii) media con las demás variables explicativas.

library(readxl)
baseraia <- read_excel("C:/Users/jrias/Downloads/[Link]")
View(baseraia)

baseraia$periodo <- [Link](baseraia$periodo)

baseraia$lugar <- [Link](baseraia$lugar)
baseraia$marea <- [Link](baseraia$marea)

head (baseraia)
attach(baseraia)

######## I. analisis exploratorio ############

# Test de distribución de Cullen & Frey

descdist(cpue, boot = 1000)

# Guardar el gráfico como una imagen

png("Cullen_Frey_Plot.png")
descdist(cpue, boot = 1000)

# periodo
baseraiap <- baseraia %>%
group_by(periodo) %>%
summarise(count = n()) %>%
mutate(percentage = count / sum(count) * 100)

g1 <- ggplot(baseraiap, aes(x = periodo, y = count, fill = periodo)) +

geom_bar(stat = "identity", [Link] = FALSE) +
geom_text(aes(label = sprintf("%.1f%%", percentage)),
vjust = -0.5,
color = "black") +
scale_y_continuous(expand = expansion(mult = c(0, 0.1)),
limits = c(0, max(baseraiap$count) * 1.1)) +
theme_minimal() +
labs(title = "Distribución de periodos",
x = "Periodo",

71
 y = "Frecuencia")
g1

# Lugar
baseraia %>%
count(lugar) %>%
mutate(pct = n / sum(n) * 100) -> baseraia_counts

g2<-ggplot(baseraia_counts, aes(x = lugar, y = n, label = sprintf("%.1f%%", pct))) +

geom_bar(stat = "identity", fill = "lightblue") +
geom_text(position = position_stack(vjust = 0.5), size = 3.5) +
theme_minimal() +
labs(title = "Distribución de Observaciones por Lugar", x = "Lugar", y = "Frecuencia")
g2

#Viento

g3<-ggplot(baseraia, aes(x = viento)) +

geom_histogram(bins = 10, fill = "darkgreen") +
theme_minimal() +
labs(title = "Histograma de Velocidad del Viento", x = "Velocidad del Viento", y = "Frecuencia")

g3
# Marea
baseraia %>%
count(marea) %>%
mutate(pct = n / sum(n) * 100) -> baseraia_countsm

g4<-ggplot(baseraia_countsm, aes(x = marea, y = n, label = sprintf("%.1f%%", pct))) +

geom_bar(stat = "identity", fill = "orange") +
geom_text(position = position_stack(vjust = 0.5), size = 3.5) +
theme_minimal() +
labs(title = "Distribución de Observaciones por marea", x = "marea", y = "Frecuencia")
g4

# Temperatura Maxima

g5 <- ggplot(baseraia, aes(x = tmax)) +

geom_histogram(aes(y = after_stat(density)), bins = 20, fill = "coral", alpha = 0.5) +
geom_density(color = "blue", size = 1, linetype = "dashed") +
theme_minimal() +
labs(title = "Histograma de Temperatura Máxima con Densidad",
x = "Temperatura Máxima",
y = "Densidad")

g5

## Temperatura Minima

g6 <- ggplot(baseraia, aes(x = tmin)) +

geom_histogram(aes(y = after_stat(density)), bins = 20, fill = "coral", alpha = 0.5) +
geom_density(color = "blue", size = 1, linetype = "dashed") +
theme_minimal() +
labs(title = "Histograma de Temperatura Minima con Densidad",
x = "Temperatura Minima",

72
 y = "Densidad")

g6

## Horas de insolacion

g7<-ggplot(baseraia, aes(x = ins)) +

geom_histogram(bins = 10, fill = "skyblue", color = "black") +
theme_minimal() +
labs(title = "Histograma de Intensidad Solar", x = "Intensidad Solar", y = "Frecuencia")
g7

## cpue (captura por unidad de esfuerzo, en kg)

g8<-ggplot(baseraia, aes(x = cpue)) +

geom_histogram(aes(y = ..density..), bins = 10, fill = "#458B74", color = "black") +
geom_density(alpha = .2, fill = "#FF6666") +
theme_minimal() +
labs(title = "Histograma de CPUE con Densidad", x = "CPUE", y = "Densidad")
g8

## relaciones entre variables

### Periodo vs CPUE

ggplot(baseraia, aes(x = cpue, fill = periodo)) +

geom_histogram(aes(y = ..density..), position = "identity", alpha = 0.6, bins = 10, color = "black") +
geom_density(alpha = .2, aes(color = periodo)) +
theme_minimal() +
labs(title = "Histograma de CPUE por Periodo", x = "CPUE", y = "Densidad") +
scale_fill_brewer(palette = "Set1") +
scale_color_brewer(palette = "Set1")

## Lugar Vs CPUE
ggplot(baseraia, aes(x = cpue)) +
geom_histogram(aes(y = ..density..), bins = 10, fill = "lightblue", color = "black") +
facet_wrap(~lugar) +
theme_minimal() +
labs(title = "Histograma de CPUE por Lugar", x = "CPUE", y = "Densidad")

## viento, CPUE, Marea

ggplot(baseraia, aes(x = viento, y = cpue)) +
geom_point(aes(color = factor(marea)), size = 2) +
theme_minimal() +
labs(title = "Relación entre Velocidad del Viento y CPUE",
x = "Velocidad del Viento",
y = "CPUE",
color = "Marea")
## temperatura maxima, insolacion, periodo y CPUE

ggplot(baseraia, aes(x = tmax, y = cpue)) +

geom_point(aes(color = periodo, size = ins)) +
theme_minimal() +
labs(title = "Relación entre Temperatura Máxima y CPUE",
x = "Temperatura Máxima (°C)",

73
 y = "CPUE",
color = "Periodo",
size = "Intensidad Solar")
## temperatura maxima, insolacion, periodo, lugar y CPUE
ggplot(baseraia, aes(x = tmax, y = cpue)) +
geom_point(aes(color = periodo, size = ins)) +
facet_wrap(~ lugar) +
theme_minimal() +
labs(title = "Temperatura Máxima vs. CPUE por Lugar",
x = "Temperatura Máxima (°C)",
y = "CPUE",
color = "Periodo",
size = "Intensidad Solar")

## marea, periodo y CPUE

ggplot(baseraia, aes(x = factor(marea), y = cpue, fill = periodo)) +

geom_boxplot() +
theme_minimal() +
labs(title = "Distribución de CPUE por Marea y Periodo",
x = "Marea",
y = "CPUE",
fill = "Periodo")

## Anlisis descriptivo numerico

resumen_categorico <- baseraia %>%

group_by(periodo, lugar, marea) %>%
summarise(
tmax_media = mean(tmax, [Link] = TRUE),
tmax_min = min(tmax, [Link] = TRUE),
tmax_max = max(tmax, [Link] = TRUE),
tmin_media = mean(tmin, [Link] = TRUE),
tmin_min = min(tmin, [Link] = TRUE),
tmin_max = max(tmin, [Link] = TRUE),
ins_media = mean(ins, [Link] = TRUE),
ins_min = min(ins, [Link] = TRUE),
ins_max = max(ins, [Link] = TRUE),
cpue_media = mean(cpue, [Link] = TRUE),
cpue_min = min(cpue, [Link] = TRUE),
cpue_max = max(cpue, [Link] = TRUE),
.groups = 'drop'
)

# Imprime la tabla de resumen

print(resumen_categorico)

####### II. Prponga un modelo gama con enlace logarítmic con todas las variables
# explicativas y por medio del criterio de AKAIKE seleccione un submodelo adecuado ##########

## modelo 1
m1gama=glm(cpue~., family = Gamma(link = log), data=baseraia)
summary(m1gama)
# AIC: 1416.5

# modelo 2

74
m1gamaga=gamlss(cpue~., family = GA([Link] ="log"), data=baseraia)
summary(m1gamaga)
# 1415.542
newpar<-par(mfrow = c(2, 2),
mar = par("mar") + c(0, 1, 0, 0),
[Link] = "blue4",
col = "tomato4",
[Link] = "blue4",
[Link] = "blue4",
pch = "+",
cex = 0.8,
[Link] = 1.2,
[Link] = 1,
[Link] = 1.2,
bg="white")
plot(m1gamaga,par=newpar)
plot(m1gamaga)
[Link](m1gama)
[Link](m1gama)

############# MODELO SIN PUNTOS INFLUYENTES #############

modelo2cook=glm(cpue ~.,family = Gamma(link = log), data=baseraia,

subset=-c(4, 5, 171, 64), )
summary(modelo2cook)
# AIC: 1346.3

[Link](modelo2cook)

modelo2cook2=gamlss(cpue~., family = GA([Link] ="log"), data=baseraia,

subset=-c(4, 5, 171, 64), )
summary(modelo2cook2)
# AIC: 1415.542

######### otros modelos ####

# Modelo total

m1gama <- glm(cpue ~ periodo + lugar + viento + marea + tmax + tmin + ins, family = Gamma(link = "log"),
data = baseraia)
summary(m1gama)

m1gamaga <- gamlss(cpue ~ periodo + lugar + viento + marea + tmax + tmin + ins, family = GA([Link] =
"log"), data = baseraia)
summary(m1gamaga)

### otras opciones

## [Link] ="log"
m2gamaga <- gamlss(cpue ~ periodo + lugar + viento + marea + tmax + tmin + ins, family = GA([Link]
="log"), data = baseraia)
summary(m2gamaga)

## [Link] = "log", [Link] ="log"

75
m3gamaga <- gamlss(cpue ~ periodo + lugar + viento + marea + tmax + tmin + ins, family = GA([Link] =
"log", [Link] ="log"), data = baseraia)
summary(m3gamaga)
# 1415.542

## Regresión Robusta Gamma

library(robustbase)
robust_glm <- glmrob(cpue ~ periodo + lugar + viento + marea + tmax + tmin + ins,
family = Gamma(link = "log"), data = baseraia)
summary(robust_glm)

# quitando variables
modelo_seleccionado_glm<- step(m1gama, direction = "both", trace = FALSE)

summary(modelo_seleccionado_glm)
# 1410.1

[Link](modelo_seleccionado_glm)
[Link](modelo_seleccionado_glm)

modelo_seleccionado_glm2=glm(cpue ~lugar+ marea + tmax, family = Gamma(link = log), data=baseraia,

subset=-c(4, 171, 18, 5, 64), )
summary(modelo_seleccionado_glm2)
# 1328.2

modelo_seleccionado_gamlss<- stepGAIC(m1gamaga, method = "backward")

summary(modelo_seleccionado_gamlss)
# 1409.071

############ III. Interacciones del modelo ###############

# 1. ^2
m1 <- glm(cpue ~ (periodo + lugar + viento + marea + tmax + tmin + ins)^2, family = Gamma(link = "log"),
data = baseraia)
summary(m1)
# AIC: 1440.3
[Link](m1)
[Link](m1)

m1 <- glm(cpue ~ (periodo + lugar + viento + marea + tmax + tmin + ins)^2, family = Gamma(link = "log"),
data = baseraia, subset=-c(127, 171, 1, 45),)
summary(m1)
# AIC: 1381.2

# 2. multiplicacion entre dos Variables

# Variables para interacciones
variables <- c("periodo", "lugar", "viento", "marea", "tmax", "tmin", "ins")

# Modelo inicial sin interacciones

modelo_base <- glm(cpue ~ (periodo + lugar + viento + marea + tmax + tmin + ins)^2, family = Gamma(link =
"log"), data = baseraia)
aic_base <- AIC(modelo_base)
modelos_mejores <- list(modelo_base = list(formula = modelo_base$formula, aic = aic_base, summary =
summary(modelo_base)))

# Probar interacciones de dos en dos

76
for (i in 1:(length(variables)-1)) {
for (j in (i+1):length(variables)) {
formula_interaccion <- [Link](paste("cpue ~", paste("periodo + lugar + viento + marea + tmax + tmin +
ins +", variables[i], "*", variables[j], sep="")))
modelo_interaccion <- gamlss(formula_interaccion, family = GA([Link] = "log"), data = baseraia)
aic_interaccion <- AIC(modelo_interaccion)

if (aic_interaccion < aic_base) {

# Almacenar modelo y AIC si es mejor
modelos_mejores[[paste(variables[i], variables[j], sep = "*")]] <- list(formula = formula_interaccion, aic =
aic_interaccion, summary = summary(modelo_interaccion))
aic_base <- aic_interaccion # Actualizar el AIC base si encuentra uno mejor
}
}
}

# mejores modelos
for (modelo in names(modelos_mejores)) {
cat("\n\nResumen del modelo con interacción", modelo, ":\n")
print(modelos_mejores[[modelo]]$summary)
cat("AIC del modelo:", modelos_mejores[[modelo]]$aic, "\n")
}

m2 <- glm(cpue ~ periodo + lugar + viento + marea + tmax + tmin + ins + periodo*lugar, family = Gamma(link
= "log"), data = baseraia)
summary(m2)

[Link](m2)
[Link](m2)

m2 <- glm(cpue ~ periodo + lugar + viento + marea + tmax + tmin + ins + periodo*lugar, family = Gamma(link
= "log"),
data = baseraia, subset=-c(4, 64, 45, 5), )
summary(m2)
# AIC: 1344.2

m3 <- glm(cpue ~ periodo + lugar + viento + marea + tmax + tmin + ins + lugar*tmin, family = Gamma(link =
"log"),
data = baseraia )
summary(m3)
[Link](m3)
[Link](m3)

m3 <- glm(cpue ~ periodo + lugar + viento + marea + tmax + tmin + ins + lugar*tmin, family = Gamma(link =
"log"),
data = baseraia, subset=-c(4, 5, 64, 118), )
summary(m3)
# AIC: 1346.9

m4 <- glm(cpue ~ periodo + lugar + viento + marea + tmax + tmin + ins + marea*tmax, family = Gamma(link =
"log"), data = baseraia)
summary(m4)

77
[Link](m4)
[Link](m4)

m4 <- glm(cpue ~ periodo + lugar + viento + marea + tmax + tmin + ins + marea*tmax, family = Gamma(link =
"log"),
data = baseraia, subset=-c(4, 5, 171, 64), )
summary(m4)
# AIC: 1344.5

m5 <- glm(cpue ~ periodo + lugar + viento + marea + tmax + tmin + ins + marea*tmin, family = Gamma(link =
"log"), data = baseraia)
summary(m5)
#AIC: 1414.1

[Link](m5)
[Link](m5)

m5 <- glm(cpue ~ periodo + lugar + viento + marea + tmax + tmin + ins + marea*tmin, family = Gamma(link =
"log"),
data = baseraia, subset=-c(4, 5, 171, 64), )
summary(m5)
# 1342.8

### Transformaciones Logarítmicas

m5 <- glm(cpue ~ periodo + lugar + log(viento) + marea + log(tmax) + log(tmin) + ins, family = Gamma(link =
"log"), data = baseraia)
summary(m5)
[Link](m5)
[Link](m5)
m5 <- glm(cpue ~ periodo + lugar + log(viento) + marea + log(tmax) + log(tmin) + ins, family = Gamma(link =
"log"),
data = baseraia, subset=-c(4, 5, 171, 45), )
summary(m5)
# AIC: 1349.3
# Modelo con transformaciones cuadráticas

m6 <- glm(cpue ~ periodo + lugar + I(viento^2) + marea + I(tmax^2) + I(tmin^2) + I(ins^2), family =
Gamma(link = "log"), data = baseraia)
summary(m6)

[Link](m6)
[Link](m6)
m6 <- glm(cpue ~ periodo + lugar + I(viento^2) + marea + I(tmax^2) + I(tmin^2) + I(ins^2), family =
Gamma(link = "log"),
data = baseraia, subset=-c(4, 5, 171, 45), )
summary(m6)
# AIC: 1348.1
# Modelo con transformaciones de raíz cuadrada
m7 <- gamlss(cpue ~ periodo + lugar + sqrt(viento) + marea + sqrt(tmax) + sqrt(tmin) + sqrt(ins),
family = GA([Link] = "log"), data = baseraia)
summary(m7)
# 1415.963
# Transformaciones Recíprocas
m8 <- gamlss(cpue ~ periodo + lugar + I(1/viento) + marea + I(1/tmax) + I(1/tmin) + ins,
family = GA([Link] = "log"), data = baseraia)
summary(m8)

78
# 1415.422
# Modelo # AIC
# GLM base 1416.5
# Galmss base 1415.542
# Glm corregido 1346.3
# Glm con selección de variables y corregido 1328,2
# con interacciones cuadráticas corregido 1381.2
# con interacciones multiplicativas (periodo*lugar) corregido 1344.2
# con interacciones multiplicativas (lugar*tmin ) corregido 1346.9
# con interacciones multiplicativas (marea*tmax ) corregido 1344.5
# con interacciones multiplicativas (marea*tmin ) corregido 1342.8
# con transformaciones logarítmicas corregido 1349.5
# con transformaciones cuadráticas 1348.7
data <- [Link](
Modelo = c("GLM base", "Galmss base", "Glm corregido", "Glm con selección de variables y corregido",
"con interacciones cuadráticas corregido", "con interacciones multiplicativas (periodo*lugar) corregido",
"con interacciones multiplicativas (lugar*tmin) corregido", "con interacciones multiplicativas
(marea*tmax) corregido",
"con interacciones multiplicativas (marea*tmin) corregido", "con transformaciones logarítmicas
corregido",
"con transformaciones cuadráticas corregido"),
AIC = c(1416.5, 1415.542, 1346.3, 1328.2, 1381.2, 1344.2, 1346.9, 1344.5, 1342.8, 1349.5, 1348.7)
)
data <- data[order(data$AIC), ]
library(ggplot2)
ggplot(data, aes(x = reorder(Modelo, AIC), y = AIC, fill = AIC)) +
geom_bar(stat = "identity", width = 0.8, color = "black", size = 0.3) +
geom_text(aes(label = sprintf("%.1f", AIC)), vjust = 1, color = "black", size = 3.5) +
coord_flip() +
labs(title = "Comparación de AIC entre modelos",
x = "Modelo",
y = "AIC") +
theme_minimal() +
theme([Link] = "none",
[Link].x = element_text(angle = 45, hjust = 1),
[Link] = element_text(hjust = 0.5)) +
scale_fill_gradient(low = "#2F4F4F", high = "#A2CD5A")
####### modelo definitivo
mx=glm(cpue ~lugar+ marea + tmax, family = Gamma(link = log), data=baseraia,
subset=-c(4, 171, 18, 5, 64), )
summary(mx)
# 1328.2
[Link](mx)

79