Análisis de Regresión: Parcial 3

PUNTO 1

Descripción de la Base de Datos

La base de datos “data.2” contiene las siguientes variables relacionadas con el documento “xslx” enviado, las cuales para el uso propuesto serán renombradas tal que:

pedidos: pdids.
Porcentaje Descuento Catálogo: PorDscCtalg.
PUP:pup.
Venta Neta Final: VentNetFin.
Demanda: Dem.

data.2 <- read_excel("C:/Users/Usuario/Downloads/Parcial 3 - II 2023.xlsx", sheet = "Hoja1")
x1<-as.numeric(data.2$PorDscCtalg)
x2<-as.numeric(data.2$pdids)
x3<-as.numeric(data.2$pup)
x4<-as.numeric(data.2$VentNetFin)

a)

Hacer la prueba de Park sobre cada variable independiente (usar un nivel de significancia del 10%). ¿Qué le dice la gráfica de los errores frente a la variable “Porcentaje Descuento Catalogo”?

Dado el modelo de “demanda” en términos de las variables “pedidos”, “PUP”, “Venta Neta Final” y “Porcentaje Descuento Catálogo” se tiene que:

\[ \text{Dem} = \beta_0 + \beta_1 \cdot \text{pdids} + \beta_2 \cdot \text{pup} + \beta_3 \cdot \text{VentNetFin} + \beta_4 \cdot \text{PorDscCtalg} + \epsilon \] En nuestro caso, el modelo está definido en R con la función lm(Dem ~ x1 + x2 + x3, data = data.2). Donde x1, x2, x3 son las “pdids”, “pup”, “VentNetFin”, “PorDscCtalg” respectivamente.

La prueba de Park nos indica la siguiente información sobre los residuos del modelo de regresión y además nos ayudará para evaluar la presencia de heteroscedasticidad en nuestro modelo de regresión; al tener un estadístico BP de 6.8031 con 3 grados de libertad, y además, un p-valor de 0.07845 que es menor al valor establecido para el nivel de significancia. Este resultado indica que no hay evidencia suficiente para rechazar la hipótesis nula de heteroscedasticidad , lo que sugiere que la varianza de los errores es constante a través de las predicciones del modelo.

modelo1<-lm(Dem ~ x1 + x2 + x3, data = data.2)
resultado_park <- bptest(modelo1)
print(resultado_park)

## 
##  studentized Breusch-Pagan test
## 
## data:  modelo1
## BP = 6.8031, df = 3, p-value = 0.07845

nivel_significancia <- 0.10
if (resultado_park$p.value <= nivel_significancia) {
        cat("Se rechaza la hipótesis nula de homocedasticidad.\n")
} else {
        cat("No hay suficiente evidencia para rechazar la hipótesis
            nula de homocedasticidad.\n")
}

## Se rechaza la hipótesis nula de homocedasticidad.

Dado lo anterior nos implica que hay evidencia estadística de que la varianza de los errores no es constante en todo el rango de los valores predichos, además, la rejeción de la hipótesis nula sugiere que existe heteroscedasticidad en el modelo, por otro lado, la presencia de heteroscedasticidad puede violar supuestos y afectar la validez de las inferencias estadísticas asociadas al modelo trabajado.

# Obtiendo los residuos del modelo
residuos <- residuals(modelo1)
plot(x1, residuos, main = "Grafico de Residuos vs. x1",
     xlab = "x1", ylab = "Residuos")
abline (0,0)
abline(h = c(-3000, 3000), col = "red", lty = 2)

El gráfico nos sugiere visualmente la ausencia en que la varianza de los errores para la variable “Porcentaje Descuento Catálogo” sea constante en todo el rango de nuestra variable. En consecuencia, podría afectar la eficiencia de los estimadores de los coeficientes de regresión, sobre esta variable y su interpretación.

b)

Hacer la regresión de mínimos cuadrados ponderados, usando como criterio de ponderación la variable: “Porcentaje Descuento Catalogo”. Hacer nuevamente la prueba de Park, pero sólo sobre la variable: Porcentaje Descuento Catalogo

# Ajustando la regresión de mínimos cuadrados ponderados
modelo_ponderado <- lm(Dem ~ x1 + x2 + x3, data = data.2, weights = x3)
summary(modelo_ponderado)

## 
## Call:
## lm(formula = Dem ~ x1 + x2 + x3, data = data.2, weights = x3)
## 
## Weighted Residuals:
##     Min      1Q  Median      3Q     Max 
## -709.31 -216.27  -21.79  156.98  711.67 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -1.845e+04  2.497e+03  -7.390 2.38e-09 ***
## x1           1.666e+04  6.227e+03   2.675   0.0103 *  
## x2           3.978e+00  2.226e-01  17.866  < 2e-16 ***
## x3           2.297e+05  1.938e+04  11.850 1.41e-15 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 325.8 on 46 degrees of freedom
## Multiple R-squared:  0.8965, Adjusted R-squared:  0.8897 
## F-statistic: 132.8 on 3 and 46 DF,  p-value: < 2.2e-16

En este escenario el modelo parece tener un buen ajuste y las variables independientes x1, x2, y x3 están relacionadas significativamente con la variable dependiente “Dem” en una comparación con el modelo anterior. Veamos la prueba de Park, pero únicamente sobre la variable: Porcentaje Descuento Catálogo:

# Realizando la prueba de Park únicamente sobre la variable 
#'Porcentaje Descuento Catalogo'
resultado_park_variable <- bptest(modelo_ponderado, ~ x3, studentize = FALSE)
print(resultado_park_variable)

## 
##  Breusch-Pagan test
## 
## data:  modelo_ponderado
## BP = 722.45, df = 1, p-value < 2.2e-16

Se observa que en esta prueba de Park arroga la información sobre los residuos del modelo de regresión para la variable ‘Porcentaje Descuento Catálogo’, luego al tener un estadístico BP de 722.45 con 1 grados de libertad, y además, un p-valor menor a 2.2e-16. Con este resultado se indica que al tener un BP tan grande, mayor es la evidencia en contra de la homoscedasticidad. Por otro lado, al tener un valor p muy pequeño (en este caso, prácticamente cero) sugiere que se tiene evidencia significativa en contra de la homoscedasticidad.

Como conclusión, se puede decir que los errores en el modelo con esta variable tratada no tienen una varianza constante, lo que podría afectar la validez de las inferencias realizadas a partir del modelo. Con la presencia de heteroscedasticidad, se podría considerar técnicas robustas o ajustes para abordar este problema.

PUNTO 2

datosP2 <- data.frame(pobl = c(500, 1200, 100, 400, 500, 300),
                    Reclmns = c(42, 37, 1, 101, 73, 14),
                    tamnoCar = factor(c("Pequeno", "Mediano", "Grande", "Pequeno",
                                        "Mediano", "Grande")),
                    grpEdad = factor(c(1, 1, 1, 2, 2, 2))
                    )

a)

Ajustar el modelo Poisson e interpretar los resultados (con y sin la variable offset).

# Con offset
datosP2 <- data.frame(pobl = c(500, 1200, 100, 400, 500, 300),
                    Reclmns = c(42, 37, 1, 101, 73, 14),
                    tamnoCar = factor(
                            c("Pequeno", "Mediano", "Grande",
                              "Pequeno", "Mediano", "Grande")),
                    grpEdad = factor(c(1, 1, 1, 2, 2, 2))
                    )
mdlo_poisson<- glm( Reclmns ~  tamnoCar + grpEdad, family = poisson, 
                    data = datosP2)
summary(mdlo_poisson)

## 
## Call:
## glm(formula = Reclmns ~ tamnoCar + grpEdad, family = poisson, 
##     data = datosP2)
## 
## Coefficients:
##                 Estimate Std. Error z value Pr(>|z|)    
## (Intercept)       1.4991     0.2747   5.458 4.81e-08 ***
## tamnoCarMediano   1.9924     0.2752   7.239 4.52e-13 ***
## tamnoCarPequeno   2.2548     0.2714   8.308  < 2e-16 ***
## grpEdad2          0.8544     0.1335   6.401 1.55e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 177.341  on 5  degrees of freedom
## Residual deviance:   5.749  on 2  degrees of freedom
## AIC: 43.856
## 
## Number of Fisher Scoring iterations: 5

Se observa que para la variable “TamañoCoche” cuando este es pequeño, según el modelo, indica que el riesgo relativo de tener un Reclamo por parte del cliente que tenga este tamaño de coche es 9 veces mayor al riesgo relativo de los reclamos por clientes con “TamañoCoche” distinto al tamaño pequeño.

Por parte de la variable “GrupoEdad”; el modelo nos indica que el riesgo relativo de que un cliente que pertenezca al grupo de edad “2” haga reclamos, es dos veces mayor al riesgo de realizar reclamos por parte del grupo “1”. En resumen, el modelo Poisson sin incluir la variable “offset”; sugiere que las variables predictoras: “TamañoCoche” cuando este es mediano y pequeño, “GrupoEdad” cuando es el grupo 2, son significativas para predecir la variable de respuesta, “Reclamaciones”, además, la “Null deviance” y “AIC” indican que el modelo ajustado tiene un buen ajuste a los datos trabajados.

# Sin offset
mdlo_poisson_offset <- glm(Reclmns ~ tamnoCar + grpEdad + offset(log(pobl)),
                           family = poisson, data = datosP2)
summary(mdlo_poisson_offset)

## 
## Call:
## glm(formula = Reclmns ~ tamnoCar + grpEdad + offset(log(pobl)), 
##     family = poisson, data = datosP2)
## 
## Coefficients:
##                 Estimate Std. Error z value Pr(>|z|)    
## (Intercept)      -4.4010     0.2868 -15.347  < 2e-16 ***
## tamnoCarMediano   1.0715     0.2784   3.848 0.000119 ***
## tamnoCarPequeno   1.7643     0.2724   6.478 9.32e-11 ***
## grpEdad2          1.3199     0.1359   9.713  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 175.1536  on 5  degrees of freedom
## Residual deviance:   2.8207  on 2  degrees of freedom
## AIC: 40.928
## 
## Number of Fisher Scoring iterations: 4

Para el modelo de regresión Poisson donde se tiene en cuenta la variable “offset”, se toma la variable “Población” y así se incorpora esta variable explicativa que tiene un efecto conocido y constante en la tasa de ocurrencia del evento, pero que no se ajusta en el modelo por medio de “offset(log(pobl))”. Así, con el modelo “mdlo_poisson_offset” se tiene que la probabilidad de que alguna persona que pertenezca al “Grupo de Edad 2”, realice algún reclamo es casi cuatro veces la probabilidad a que el reclamo lo realice una persona del “Grupo de Edad 1”. Respecto al “TamañoCoche” al considerar la variable offset y en símil con el modelo sin la variable offset, se evidencia que la variable “TamañoCoche” cuando este es pequeño; indica que el riesgo relativo de tener un Reclamo por parte del cliente, que tenga este tamaño de coche, es un poco menos de 6 veces mayor al riesgo relativo de los reclamos por clientes con “TamañoCoche” distinto al tamaño pequeño respecto a la población. Del mismo modelo se observa que los Pr(>|z|) para cada coeficiente, están indicando que los coeficientes son significativamente diferentes de cero.

Por último, y en comparación con el modelo “mdlo_poisson” el cual tenía un AIC de 43.856 existe una diferencia frente al AIC del modelo “mdlo_poisson_offset”, la cual indica que este último modelo con la variable “offset(log(pobl))”; tiene mejores ajustes respecto, al tener un AIC de 40.93. Y en general, este modelo ajustado con la variable offset parece tener un buen ajuste a los datos, con todos los coeficientes siendo significativos. La deviance residual y el AIC sugieren que este modelo es adecuado para describir los datos observados.

b)

Calcular la probabilidad que se den 50 reclamaciones para una población de 1000, en tamaños de coche pequeño en el grupo de edad 1 (con y sin la variable offset)

Pobla.ejmpl <- data.frame(pobl = 1000,
                          tamnoCar = factor("Pequeno"),
                          grpEdad = factor(1))
Modl.sin.offset <- predict(mdlo_poisson, 
                           newdata = Pobla.ejmpl,
                           type = "response")
Modl.con.offset <- predict(mdlo_poisson_offset,
                           newdata = Pobla.ejmpl,
                           type = "response")
cat("Probabilidad de 50 reclamaciones sin offset:", 
    Modl.sin.offset, "\n")

## Probabilidad de 50 reclamaciones sin offset: 42.68657

cat("Probabilidad de 50 reclamaciones con offset:",
    Modl.con.offset, "\n")

## Probabilidad de 50 reclamaciones con offset: 71.5978

Se tiene que la probabilidad de que se den 50 reclamaciones en una población de 1000 habitantes, cuyas reclamaciones son realizadas por personas que tienen coche pequeño en el grupo de edad 1 es de 42.7%, sin la variable offset. Incluyendo la variable offset, se tiene que la probabilidad en el mismo escenario será de 71.6%.