Tarea 4
ARODNAO
1/6/2020
Solución Tarea 4 Modelo Econométrico
Autor: ANGEL RODRIGUEZ NAOLA
Cargando Librerías de R:
library(xtable)
library(knitr)
library(car)## Loading required package: carDatalibrary(plyr)
library(Hmisc)## Loading required package: lattice## Loading required package: survival## Loading required package: Formula## Loading required package: ggplot2##
## Attaching package: 'Hmisc'## The following objects are masked from 'package:plyr':
##
## is.discrete, summarize## The following objects are masked from 'package:xtable':
##
## label, label<-## The following objects are masked from 'package:base':
##
## format.pval, unitslibrary(broom)
library(normtest)
library(RColorBrewer)
library(ggsci)
library(lmtest)## Loading required package: zoo##
## Attaching package: 'zoo'## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numericlibrary(corrplot)## corrplot 0.84 loadedlibrary(readxl)
library(dynlm)Seleccionando la Base de Datos Prestige
data<-Prestige
head(data)## education income women prestige census type
## gov.administrators 13.11 12351 11.16 68.8 1113 prof
## general.managers 12.26 25879 4.02 69.1 1130 prof
## accountants 12.77 9271 15.70 63.4 1171 prof
## purchasing.officers 11.42 8865 9.11 56.8 1175 prof
## chemists 14.62 8403 11.68 73.5 2111 prof
## physicists 15.64 11030 5.13 77.6 2113 prof#dataAlmacenamos el objeto data la base da datos Prestige, para luego poder trabajar con ella.
Veremos si hay valores perdidos
summary(data)## education income women prestige
## Min. : 6.380 Min. : 611 Min. : 0.000 Min. :14.80
## 1st Qu.: 8.445 1st Qu.: 4106 1st Qu.: 3.592 1st Qu.:35.23
## Median :10.540 Median : 5930 Median :13.600 Median :43.60
## Mean :10.738 Mean : 6798 Mean :28.979 Mean :46.83
## 3rd Qu.:12.648 3rd Qu.: 8187 3rd Qu.:52.203 3rd Qu.:59.27
## Max. :15.970 Max. :25879 Max. :97.510 Max. :87.20
## census type
## Min. :1113 bc :44
## 1st Qu.:3120 prof:31
## Median :5135 wc :23
## Mean :5402 NA's: 4
## 3rd Qu.:8312
## Max. :9517Existen 4 valores perdidos en la variable Type, lo omitiremos solos para el análisis descriptivo.
Cambiando de nombre a las variables
data<-rename(data,
replace = c(education="educ",
income="ing",
women="mujeres",
prestige="prestigio",
census="cod",
type="tipo"))
head(data)## educ ing mujeres prestigio cod tipo
## gov.administrators 13.11 12351 11.16 68.8 1113 prof
## general.managers 12.26 25879 4.02 69.1 1130 prof
## accountants 12.77 9271 15.70 63.4 1171 prof
## purchasing.officers 11.42 8865 9.11 56.8 1175 prof
## chemists 14.62 8403 11.68 73.5 2111 prof
## physicists 15.64 11030 5.13 77.6 2113 profEtiquetando a las variables
label(data$educ)="Promedio de años de estudio del trabajador nombrado"
label(data$ing)="Promedio de ingreso de dólares del nombrado"
label(data$mujeres)="Porcentaje de mujeres nombradas"
label(data$prestigio)="Valoración de prestigio de Pineo-Porter"
label(data$cod)="Código de ocupación en el censo"
label(data$tipo)="Tipo de ocupación en 3 niveles"
label(data)## educ
## "Promedio de años de estudio del trabajador nombrado"
## ing
## "Promedio de ingreso de dólares del nombrado"
## mujeres
## "Porcentaje de mujeres nombradas"
## prestigio
## "Valoración de prestigio de Pineo-Porter"
## cod
## "Código de ocupación en el censo"
## tipo
## "Tipo de ocupación en 3 niveles"Recodificando la variable Tipo de Educación
data$tipo=as.character(data$tipo)
data$tipo[data$tipo=="bc"]=1
data$tipo[data$tipo=="prof"]=2
data$tipo[data$tipo=="wc"]=3
head(data)## educ ing mujeres prestigio cod tipo
## gov.administrators 13.11 12351 11.16 68.8 1113 2
## general.managers 12.26 25879 4.02 69.1 1130 2
## accountants 12.77 9271 15.70 63.4 1171 2
## purchasing.officers 11.42 8865 9.11 56.8 1175 2
## chemists 14.62 8403 11.68 73.5 2111 2
## physicists 15.64 11030 5.13 77.6 2113 2Recodificando en factor la variable Tipo de Educación y etiquetando
data$tipo=as.factor(data$tipo)
a<-levels(data$tipo)
data$tipo=factor(data$tipo,
ordered = F,
levels = a,
labels = c("cuello azul",
"profesional",
"cuello blanco"))
head(data)## educ ing mujeres prestigio cod tipo
## gov.administrators 13.11 12351 11.16 68.8 1113 profesional
## general.managers 12.26 25879 4.02 69.1 1130 profesional
## accountants 12.77 9271 15.70 63.4 1171 profesional
## purchasing.officers 11.42 8865 9.11 56.8 1175 profesional
## chemists 14.62 8403 11.68 73.5 2111 profesional
## physicists 15.64 11030 5.13 77.6 2113 profesionalRealizando el análisis descriptivo
layout(matrix(c(1,3,2,4), nrow = 2,ncol= 2))
plot(data$educ,data$prestigio,col="blue",main="Prestigio vs Promedio años de estudio",
xlab="Estudio",ylab="Prestigio")
plot(data$ing,data$prestigio,col="red",main="Prestigio vs Promedio de ingreso",
xlab="Ingreso",ylab="Prestigio")
plot(data$mujeres,data$prestigio,,col="darkgreen",main="Prestigio vs % de mujeres nombradas",
xlab="Mujeres Nombradas",ylab="Prestigio")
plot(data$tipo,data$prestigio,col = pal_jco()(10),main="Prestigio vs Tipo de Ocupación",
xlab="Tipo de Ocupación",ylab="Prestigio")
Se observa que existe relación lineal entre la variable dependiente Nivel de Prestigio y las variables independientes Promedio de años de estudio y Promedio de ingreso, las demás variables no muestran relación alguna con el nivel de Prestigio
PREGUNTA 1 Crea el siguiente modelo:
\[\begin{equation} prestigio = \beta_0 + \beta_1 educ + \beta_2 \sqrt{ing} + \epsilon. \end{equation}\]
Creando una nueva base y el modelo de regresión lineal múltiple respectivo
data2 <- (Prestige)
data2<-rename(data2,
replace = c(education="educ",
income="ing",
women="mujeres",
prestige="prestigio",
census="cod",
type="tipo"))
modelo<-lm(prestigio~educ+sqrt(ing),data = data2)PREGUNTA 2 Realice una interpretación de los coeficientes del modelo.
summary(modelo)##
## Call:
## lm(formula = prestigio ~ educ + sqrt(ing), data = data2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -16.2402 -5.1954 -0.3229 4.2815 17.9898
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -18.17457 3.13510 -5.797 8.06e-08 ***
## educ 3.91380 0.33209 11.785 < 2e-16 ***
## sqrt(ing) 0.29018 0.03933 7.379 5.01e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 7.349 on 99 degrees of freedom
## Multiple R-squared: 0.8212, Adjusted R-squared: 0.8175
## F-statistic: 227.3 on 2 and 99 DF, p-value: < 2.2e-16Modelo de Regresión Lineal
\[\begin{equation} prestigio = -181.7457 + 3.91380*educ + 0.29018*\sqrt{ing} + \epsilon. \end{equation}\]
Los predictores en su media nos permite decir que el Prestigio promedio estimado es de -18.17457, cuando no consideramos el número promedio de años de educación y el Promedio de ingreso de dólares del nombrado del conjunto de datos.Es decir se espera en promedio generar un valor de -18.17 de Prestigio cuando no influye ningún factor.
Del resultado del modelo y el diagrama de dispersión podemos hacer algunas observaciones interesantes:
Para cualquier nivel de ingreso en una profesión, al aumentar un año del nivel de educación en una profesión determinada, el Prestigio promedio aumentará en 3.91380.
Del mismo modo, para cualquier nivel de educación, ver una mejora en la raíz cuadrada del ingreso en un dólar en una profesión determinada generará un Prestigio promedio adicional de 0.29018.
PREGUNTA 3 Evalúe e interprete la bondad de ajuste del modelo (R2).
m1<-glance(modelo)
l<-labels(m1)
l## [[1]]
## [1] "1"
##
## [[2]]
## [1] "r.squared" "adj.r.squared" "sigma" "statistic"
## [5] "p.value" "df" "logLik" "AIC"
## [9] "BIC" "deviance" "df.residual"comp<-matrix(data = c(l[[2]],m1), nrow = 11)
comp_2<-as.data.frame(comp)
names(comp_2)## [1] "V1" "V2"comp_2<-rename(comp_2,c(V1="CARACTERÍSTICAS",V2="M1"))
comp_2## CARACTERÍSTICAS M1
## 1 r.squared 0.8211581
## 2 adj.r.squared 0.8175451
## 3 sigma 7.348851
## 4 statistic 227.2807
## 5 p.value 9.936502e-38
## 6 df 3
## 7 logLik -346.6527
## 8 AIC 701.3054
## 9 BIC 711.8053
## 10 deviance 5346.556
## 11 df.residual 99Observando el R cuadrado ajustado como una métrica de variabilidad más apropiada ya que el R cuadrado ajustado aumenta solo si el nuevo término agregado termina mejorando el modelo más de lo que se esperaría. El R2 de 0.8212 o 82.12% explican el evento de predecir el Valor del Prestigio en función al nivel de educación (añños de estudio) y nivel de ingresos (en dólares). El R2 tiene un valor de 0.8212 por lo que las variables Ingreso y Educación explican solamente el 82.12% de la variablidad total de Prestigio.
PREGUNTA 4 Evalúe e interprete la prueba de significancia individual del modelo.
summary(modelo)##
## Call:
## lm(formula = prestigio ~ educ + sqrt(ing), data = data2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -16.2402 -5.1954 -0.3229 4.2815 17.9898
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -18.17457 3.13510 -5.797 8.06e-08 ***
## educ 3.91380 0.33209 11.785 < 2e-16 ***
## sqrt(ing) 0.29018 0.03933 7.379 5.01e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 7.349 on 99 degrees of freedom
## Multiple R-squared: 0.8212, Adjusted R-squared: 0.8175
## F-statistic: 227.3 on 2 and 99 DF, p-value: < 2.2e-16Ho: Bi=0 ( la variable xi no aporta al modelo)
H1: Bi<>0 ( la variable xi aporta al modelo)
alfa 0.05
Rechazamos la Ho cuando el p-value es menor 0.05,
T: 11.875, p=0.0000000000000002 < 0.05 (No se rechaza)/Años Promedio En Educación
T: 7.379, p=0.0000000000501 <0.05 (No se rechaza)/ Raíz cuadrada del ingreso en dólares
Podemos concluir con un nivel de confianza del 95%, que resultaron significativas las siguientes variables: Educación y Raíz Cuadrada del Ingreso.
PREGUNTA 5 Evalúe e interprete la prueba de significancia global del modelo.
summary(modelo)##
## Call:
## lm(formula = prestigio ~ educ + sqrt(ing), data = data2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -16.2402 -5.1954 -0.3229 4.2815 17.9898
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -18.17457 3.13510 -5.797 8.06e-08 ***
## educ 3.91380 0.33209 11.785 < 2e-16 ***
## sqrt(ing) 0.29018 0.03933 7.379 5.01e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 7.349 on 99 degrees of freedom
## Multiple R-squared: 0.8212, Adjusted R-squared: 0.8175
## F-statistic: 227.3 on 2 and 99 DF, p-value: < 2.2e-16Prueba de modelo (Anova)
Ho: B=0 ( modelo no es adecuado)
H1: B<>0( modelo es adecuado)
alfa 0.05
Rechazamos la Ho cuando el p-value es menor 0.05, F:227.3, p=0.00000000000000022 < 0.05.Podemos concluir con un nivel de confianza del 95%, que el modelo es adecuado.
PREGUNTA 6 Verifique la existencia de error de especificación en el modelo. De haber error corríjalo (puede hacer uso de variables dummies si usted cree que es necesario).Realice su interpretación, en cualquiera de los casos.
Una herramienta potencialmente útil para construir un modelo apropiado es la prueba de especificación de Ramsey, RESET. Este método agrega automáticamente términos polinomiales de orden superior a su modelo y prueba la hipótesis conjunta de que sus coeficientes son todos cero. Por lo tanto, la hipótesis nula de la prueba es
H0: “No son necesarios términos polinomiales de orden superior”;
Si rechazamos la hipótesis nula, debemos considerar incluir tales términos.
Para confirmar las observaciones visuales, se suele utilizar un test estadístico para diagnosticar una mala especificación funcional del modelo: RESET Test de Ramsey. La idea es justamente evaluar si es que existe un error de especificación de la ecuación de regresión. Este test lo que hace es volver a estimar el modelo pero incorporando los valores predichos del modelo original con alguna transformación no lineal de las variables. Luego, a partir de un Test-F se evalúa si el modelo con la especificación no lineal tiene un mejor ajuste que el modelo original sin la transformación no lineal.
La hipótesis nula postula que las nuevas variables (en este caso ing^(1/2) no aportan significativamente para explicar la variación de la variable dependiente; es decir, que su coeficiente es igual a cero ( β = 0).
Ho: Fórmula Funcional Correcta
Ha: Fórmula Funcional Incorrecta
La función R que realiza una prueba RESET es resettest, que requiere los siguientes argumentos: fórmula, la fórmula del modelo a probar o el nombre de un objeto lm ya calculado; potencia, un conjunto de enteros que indican las potencias de los términos polinomiales que se incluirán; tipo, que podría ser uno de “fitted”, “regressor” o “princomp”, que indica si los términos adicionales deben ser potencias de los regresores, valores ajustados o el primer componente principal de la matriz del regresor; y, finalmente, datos, que especifican el conjunto de datos que se utilizará si se ha proporcionado una fórmula y no un objeto modelo.
modelo <- lm(prestigio ~ educ + sqrt(ing), data = data2)
summary(modelo)##
## Call:
## lm(formula = prestigio ~ educ + sqrt(ing), data = data2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -16.2402 -5.1954 -0.3229 4.2815 17.9898
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -18.17457 3.13510 -5.797 8.06e-08 ***
## educ 3.91380 0.33209 11.785 < 2e-16 ***
## sqrt(ing) 0.29018 0.03933 7.379 5.01e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 7.349 on 99 degrees of freedom
## Multiple R-squared: 0.8212, Adjusted R-squared: 0.8175
## F-statistic: 227.3 on 2 and 99 DF, p-value: < 2.2e-16modelo1 <- lm(prestigio ~ educ + ing, data = data2)
summary(modelo1)##
## Call:
## lm(formula = prestigio ~ educ + ing, data = data2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -19.4040 -5.3308 0.0154 4.9803 17.6889
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -6.8477787 3.2189771 -2.127 0.0359 *
## educ 4.1374444 0.3489120 11.858 < 2e-16 ***
## ing 0.0013612 0.0002242 6.071 2.36e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 7.81 on 99 degrees of freedom
## Multiple R-squared: 0.798, Adjusted R-squared: 0.7939
## F-statistic: 195.6 on 2 and 99 DF, p-value: < 2.2e-16modelo2 <- lm(prestigio ~ educ + ing^2, data = data2)
summary(modelo2)##
## Call:
## lm(formula = prestigio ~ educ + ing^2, data = data2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -19.4040 -5.3308 0.0154 4.9803 17.6889
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -6.8477787 3.2189771 -2.127 0.0359 *
## educ 4.1374444 0.3489120 11.858 < 2e-16 ***
## ing 0.0013612 0.0002242 6.071 2.36e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 7.81 on 99 degrees of freedom
## Multiple R-squared: 0.798, Adjusted R-squared: 0.7939
## F-statistic: 195.6 on 2 and 99 DF, p-value: < 2.2e-16resettest(modelo)##
## RESET test
##
## data: modelo
## RESET = 1.8168, df1 = 2, df2 = 97, p-value = 0.168resettest(modelo1)##
## RESET test
##
## data: modelo1
## RESET = 3.4709, df1 = 2, df2 = 97, p-value = 0.035resettest(modelo2)##
## RESET test
##
## data: modelo2
## RESET = 3.4709, df1 = 2, df2 = 97, p-value = 0.035resettest(modelo, power=2, type=c("fitted"), data=data2)##
## RESET test
##
## data: modelo
## RESET = 1.1888, df1 = 1, df2 = 98, p-value = 0.2783resettest(modelo2, power=2, type=c("fitted"), data=data2)##
## RESET test
##
## data: modelo2
## RESET = 2.3298, df1 = 1, df2 = 98, p-value = 0.1301De acuerdo al resultado del Test F, confirmamos lo siguiente: el incorporar un término raíz cuadrada del ingreso mejora el ajuste de nuestra estimación. A esta conclusión llegamos observando el valor-p del test RESET de Ramsey mucho mayor que el nivel de significancia estadística del 5%, no rechazamos la hipótesis nula de que la incorporación del término raíz cuadrada e inclusive la incorporación del término cuadrático del ingreso no es necesario en la mejora del ajuste del modelo.
En nuestro caso, ambos valores p del modelo y modelo2 son superiores a 0,05, lo que indica que el modelo no falla marginalmente la prueba de especificación y que pueden no ser necesarios algunos términos de orden superior. Es decir nos indica que no hay error de especificación ya que no rechazamos la hipótesis nula que formula una forma funcional correcta.
PREGUNTA 7 Evalúe e interprete la multicolinealidad del modelo usando la matriz de correlaciones
head(data2)## educ ing mujeres prestigio cod tipo
## gov.administrators 13.11 12351 11.16 68.8 1113 prof
## general.managers 12.26 25879 4.02 69.1 1130 prof
## accountants 12.77 9271 15.70 63.4 1171 prof
## purchasing.officers 11.42 8865 9.11 56.8 1175 prof
## chemists 14.62 8403 11.68 73.5 2111 prof
## physicists 15.64 11030 5.13 77.6 2113 profdatoscorr<-data.frame(prestigio=data2$prestigio,ingreso=data2$ing,educacion=sqrt(data2$educ))
head(datoscorr)## prestigio ingreso educacion
## 1 68.8 12351 3.620773
## 2 69.1 25879 3.501428
## 3 63.4 9271 3.573514
## 4 56.8 8865 3.379349
## 5 73.5 8403 3.823611
## 6 77.6 11030 3.954744modelocorr=cor(datoscorr)
modelocorr## prestigio ingreso educacion
## prestigio 1.0000000 0.7149057 0.8402604
## ingreso 0.7149057 1.0000000 0.5641525
## educacion 0.8402604 0.5641525 1.0000000corrplot(modelocorr, method = "number")
cor.mtest <- function(mat, ...) {
mat <- as.matrix(mat)
n <- ncol(mat)
p.mat<- matrix(NA, n, n)
diag(p.mat) <- 0
for (i in 1:(n - 1)) {
for (j in (i + 1):n) {
tmp <- cor.test(mat[, i], mat[, j], ...)
p.mat[i, j] <- p.mat[j, i] <- tmp$p.value
}
}
colnames(p.mat) <- rownames(p.mat) <- colnames(mat)
p.mat
}
p.mat <- cor.mtest(modelocorr)
head(p.mat[, 1:3])## prestigio ingreso educacion
## prestigio 0.0000000 0.59880824 0.61355761
## ingreso 0.5988082 0.00000000 0.01474937
## educacion 0.6135576 0.01474937 0.00000000col <- colorRampPalette(c("#BB4444", "#EE9988", "#FFFFFF", "#77AADD", "#4477AA"))
corrplot(modelocorr, method="color", col=col(200),
type="upper", order="hclust",
addCoef.col = "black", # Add coefficient of correlation
tl.col="black", tl.srt=45, #Text label color and rotation
# Combine with significance
p.mat = p.mat, sig.level = 0.01, insig = "blank",
# hide correlation coefficient on the principal diagonal
diag=FALSE
)
matriz of the p-value de la correlación
La matriz de correlación que se muestra resalta la situación que encontramos con el resultado del modelo. Tenga en cuenta que la correlación entre la variable independiente educación y la variable dependiente prestigio es muy alta 0.84.Además también se observa que la correlación entre la variable independiente ingreso y la variable dependiente prestigio es muy alta en 0.71.
Esto revela que el nivel de educación de cada profesión y el nivel de ingreso está fuertemente alineado con el nivel de prestigio de cada profesión.También se observa que la correlación entre la variable educación e ingreso es baja (0.56) por lo que se descartaría la existencia de Multicolinealidad.
Dado que observamos de que las dos variables predictores (ingreso y educación) están asociados con el Prestigio, podemos considerar no eliminar ninguna variable predictora en el modelo.
PREGUNTA 8 Evalúe e interprete la multicolinealidad del modelo usando la matriz de dispersiones.
plot(datoscorr, pch = 16, col = "cornflowerblue")
El diagrama de matriz de dispersiones nos permite visualizar la relación entre todas las variables en una sola imagen.
Observamos que el nivel de prestigio aumenta cuando se incrementa el promedio de años de educación. También se observa que a medida que aumenta el ingreso promedio el nivel de prestigio aumenta. Concluimos de la gráfica mostrada que no existe multicolinealidad.
PREGUNTA 9 Evalúe e interprete la existencia de multicolinealidad del modelo usando el Factor de Inflación de Varianzas (VIF).
Para detectar si la multicolinealidad es problemática es necesario realizar un test de vif (variance inflation factors) porque ver correlación de a pares no nos ayuda a dilucidar si más de dos variables tienen una correlación lineal. Lo que nos dice este test vif es qué tanto se “agrandan” los errores de cada coeficiente en presencia de las demás variables (qué tanto se incrementa la varianza del error).
vif(modelo)## educ sqrt(ing)
## 1.535454 1.535454Luego, realizamos una consulta sobre si la raíz cuadrada de vif para cada variable es menor que 2 (la raíz cuadrada porque lo que me interesa es el error estándar y no la varianza). Vif: debería ser menor a 2, si es mayor a dos quiere decir que la varianza es demasiado alta y por tanto hay problema de multicolinealidad.
sqrt(vif(modelo)) ## educ sqrt(ing)
## 1.239134 1.239134De acuerdo a la consulta, parece ser que no tenemos un problema de multicolinealidad. Al ejecutar la prueba VIF(Indice de Variancia Inflada) nos arroja valores muy bajos menores a 2 y podemos considerar que hay ausencia de multicolinealidad.
PREGUNTA 10 Evalúe e interprete la existencia de heterocedasticidad del modelo usando cualquiera de los tests vistos en clase.
Para comprobar trabajaremos con los test de White y el test de Brush-Pagan.
El test de White se utiliza para ver si hay heteroestacidad lo mismo que el test de Brush-Pagan,los dos se encuentran en el paquete lmtest de R.
Ho: Si p-valor > 0 Existencia de Homocedasticidad Ha: Si p-valor < 0 Existencia de Heterocedasticidad
test de white
gqtest(modelo) ##
## Goldfeld-Quandt test
##
## data: modelo
## GQ = 0.68379, df1 = 48, df2 = 48, p-value = 0.9042
## alternative hypothesis: variance increases from segment 1 to 2Como podemos observar en esta prueba estadística nos dice que si el p-valor es mayor al nivel de significancia del 5% entonces no rechazamos la hipótesis nula de homocedasticidad dado que p-valor =0.9042 > 0.05, es decir la varianza del error es constante, el error (residuos) varía según los valores de las variables educación e ingreso.
test de Brusch-Pagan
bptest(modelo,studentize=T)##
## studentized Breusch-Pagan test
##
## data: modelo
## BP = 5.8367, df = 2, p-value = 0.05402Como podemos observar en esta prueba estadística nos dice que si el p-valor es mayor al nivel de significancia del 5% entonces no rechazamos la hipótesis nula de homocedasticidad dado que p-valor =0.05402 > 0.05 (ligeramente mayor), es decir la varianza del error es constante, el error (residuos) varía según los valores de las variables educación e ingreso.
Evaluación del modelo
Antes de aceptar el resultado de una regresión lineal, es importante evaluar su idoneidad para explicar los datos. Una de las muchas maneras de hacer esto es examinar visualmente los residuos. Si el modelo es apropiado, los errores residuales deben ser aleatorios y normalmente distribuidos. Además, eliminar un caso no debería afectar significativamente la idoneidad del modelo. R proporciona cuatro enfoques gráficos para evaluar un modelo utilizando el comando plot ().
layout(matrix(1:4, ncol= 2))
plot(modelo)
La gráfica en la esquina superior izquierda muestra los errores residuales graficados versus sus valores ajustados. Los residuos deben distribuirse aleatoriamente alrededor de la línea horizontal que representa un error residual de cero; es decir, no debería haber una tendencia distinta en la distribución de puntos.
El gráfico en la parte inferior izquierda es un gráfico QQ estándar, que debería sugerir que los errores residuales se distribuyen normalmente.
El gráfico de ubicación de escala en la esquina superior derecha muestra la raíz cuadrada de los residuos estandarizados (una especie de raíz cuadrada de error relativo) en función de los valores ajustados. Nuevamente, no debería haber una tendencia obvia en esta trama.
Finalmente, la gráfica en la esquina inferior derecha muestra el apalancamiento de cada punto, que es una medida de su importancia para determinar el resultado de la regresión. Superpuestos en el gráfico hay líneas de contorno para la distancia de Cook, que es otra medida de la importancia de cada observación para la regresión. Las distancias más pequeñas significan que eliminar la observación tiene poco efecto en los resultados de la regresión. Las distancias mayores de 1 son sospechosas y sugieren la presencia de un posible valor atípico o deficiente.
summary(modelo)##
## Call:
## lm(formula = prestigio ~ educ + sqrt(ing), data = data2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -16.2402 -5.1954 -0.3229 4.2815 17.9898
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -18.17457 3.13510 -5.797 8.06e-08 ***
## educ 3.91380 0.33209 11.785 < 2e-16 ***
## sqrt(ing) 0.29018 0.03933 7.379 5.01e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 7.349 on 99 degrees of freedom
## Multiple R-squared: 0.8212, Adjusted R-squared: 0.8175
## F-statistic: 227.3 on 2 and 99 DF, p-value: < 2.2e-16confint(modelo)## 2.5 % 97.5 %
## (Intercept) -24.395291 -11.953844
## educ 3.254854 4.572751
## sqrt(ing) 0.212147 0.368209residuales<-modelo$residuals
par(mfrow=c(1,2))
qqnorm(residuales)
residuales_r<-as.data.frame(residuales)
is.atomic(residuales_r)## [1] FALSEqqline(residuales_r)
hist(residuales)
par(mfrow=c(1,1))Los residuos tienen distribución normal
Aplicaremos el Test de Normalidad en residuos para comprobar lo anterior
El test de Jarque-Bera no requiere estimaciones de los parámetros que caracterizan la normal. El estadístico de Jarque-Bera cuantifica que tanto se desvían los coeficientes de asimetría y curtosis de los esperados en una distribución normal. Puede calcularse mediante la función jarque.bera.test() del paquete tseries.
jb.norm.test(residuales)##
## Jarque-Bera test for normality
##
## data: residuales
## JB = 0.35297, p-value = 0.8315Hipótesis
H0: La muestra proviene de una distribución normal.
H1: La muestra no proviene de una distribución normal.
Nivel de Significancia
El nivel de significancia que se trabajará es de 0.05. Alfa=0.05
Criterio de Decisión
Si P < Alfa Se rechaza Ho
Si p >= Alfa No se rechaza Ho
Con la prueba de Jarque Bera, observamos que el p-valor = 0.818 > 0.05 No se rechaza Ho, entonces los Residuos tienen distribución normal, concluyendo que el modelo es adecuado para su implementación.
PREGUNTA 11 Para la base de datos adjuntada a la tarea que representa la evolución de la inflación, desde 1950 hasta 2015, estime un modelo donde la inflación sea explicada por ella misma en un periodo rezagado, tal y como se hizo en los ejemplos. Luego de esto pruebe con los tres métodos propuestos en clase la existencia de autocorrelación. De existir autocorrelación, corríjala.
PREGUNTA_11 <- read_excel("C:/Users/ANGEL/Desktop/MODELO ECONOMETRICO EN RSTUDIO/SESION 6/PREGUNTA 11.xlsx")
head(PREGUNTA_11)## # A tibble: 6 x 2
## tiempo inf
## <dbl> <dbl>
## 1 1950 10.3
## 2 1951 9.95
## 3 1952 6.92
## 4 1953 9.09
## 5 1954 5.40
## 6 1955 4.70Creando la base Inflación
Inflacion=PREGUNTA_11
head(Inflacion)## # A tibble: 6 x 2
## tiempo inf
## <dbl> <dbl>
## 1 1950 10.3
## 2 1951 9.95
## 3 1952 6.92
## 4 1953 9.09
## 5 1954 5.40
## 6 1955 4.70attach(Inflacion) Aplicando el modelo de regresión lineal con rezago de la variable inflación
ModelReg=dynlm(inf~L(inf),data = Inflacion)
summary(ModelReg)## Warning in summary.lm(ModelReg): essentially perfect fit: summary may be
## unreliable##
## Time series regression with "numeric" data:
## Start = 1, End = 66
##
## Call:
## dynlm(formula = inf ~ L(inf), data = Inflacion)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.062e-12 5.380e-14 5.560e-14 6.030e-14 5.217e-13
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.399e-13 6.477e-14 -2.16e+00 0.0345 *
## L(inf) 1.000e+00 6.370e-17 1.57e+16 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 5.156e-13 on 64 degrees of freedom
## Multiple R-squared: 1, Adjusted R-squared: 1
## F-statistic: 2.464e+32 on 1 and 64 DF, p-value: < 2.2e-16Aplicando el test de Durbin Watson para probar existencia de autocorrelación
dwtest(ModelReg)##
## Durbin-Watson test
##
## data: ModelReg
## DW = 2.2369, p-value = 0.8194
## alternative hypothesis: true autocorrelation is greater than 0Ho: No Existe Autocorrelación
Ha: Existe Autocorrelación
El modelo no presenta autocorrelación ya que el p-valor del modelo es mayor del nivel de significancia del 5% (0.8194 > 0.05) ,no se rechaza la hipótesis nula de que no presenta autocorrelación, en el modelo de rezago.
Aplicando el test de Ljung-Box para probar la ausencia de autocorrelación en serie, hasta un retraso especificado k.
Box.test(residuals(ModelReg),type="Ljung-Box")##
## Box-Ljung test
##
## data: residuals(ModelReg)
## X-squared = 1.1054, df = 1, p-value = 0.2931Ho: No Existe Autocorrelación
Ha: Existe Autocorrelación
El modelo no presenta autocorrelación ya que el p-valor del modelo es mayor del nivel de significancia del 5% (0.2931 > 0.05) ,no se rechaza la hipótesis nula de que no presenta autocorrelación, en el modelo de rezago.
Aplicando el test de Breusch-Godfrey para la correlación serial de orden superior
bgtest(ModelReg)##
## Breusch-Godfrey test for serial correlation of order up to 1
##
## data: ModelReg
## LM test = 1.0575, df = 1, p-value = 0.3038Ho: No Existe Autocorrelación de n-ésimo orden
\[\begin{equation} \rho_1=\rho_2=\rho_3=\rho_4=......\rho_n=0 \end{equation}\]
Ha: Existe Autocorrelación de n-ésimo orden
El modelo no presenta autocorrelación de n-ésimo orden ya que el p-valor del modelo es mayor del nivel designificancia del 5% (0.3038 > 0.05) ,no se rechaza la hipótesis nula de que no presenta autocorrelación, en el modelo de rezago.