# Paquetes
library(lmtest)
## Warning: package 'lmtest' was built under R version 4.2.3
## Loading required package: zoo
## Warning: package 'zoo' was built under R version 4.2.3
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
library(car)
## Warning: package 'car' was built under R version 4.2.3
## Loading required package: carData
## Warning: package 'carData' was built under R version 4.2.3
library(plm)
# Creación de vectores de datos
cre_chi <- c(3.9071,9.2941,14.2162,13.8676,13.0522,10.9492,9.9284,9.2308,7.8376,7.6675,8.4915,8.3399,9.1306,10.0356,10.1112,11.3958,12.7195,14.2314,9.6543,9.3998,10.6361,9.5509,7.8596,7.7686,7.2995,6.9053,6.7367)
scf_chi <- c(23.99,25.70,30.34,37.09,34.44,32.34,31.64,31.00,32.88,32.53,32.57,33.45,35.05,38.25,39.52,39.42,38.72,37.89,39.06,43.81,43.92,43.86,44.24,44.51,43.85,42.09,41.55)
cae_chi <- c(24.27,25.94,30.14,36.05,35.76,34.27,33.81,34.53,32.42,33.52,39.41,38.52,42.74,51.80,59.50,62.20,64.47,62.19,57.61,45.18,50.71,50.74,48.26,46.74,44.90,39.46,36.89)
fch_chi <- c(2.01,2.32,2.7,2.5,2.5,2.4,2.3,2.32,2.56,2.5,2.6,2.8,3.01,2.91,2.85,2.8,3,3.25,3.35,3.5,3.6,3.9,4.4,4.15,4.1,4.25,4.2)
# Creación de series de tiempo
ts_cre_chi <- ts(cre_chi, start=1990)
ts_scf_chi <- ts(scf_chi, start=1990)
ts_cae_chi <- ts(cae_chi, start=1990)
ts_fch_chi <- ts(fch_chi, start=1990)
# Logaritmos de las series
lts_cre_chi <- log(ts_cre_chi)
lts_scf_chi <- log(ts_scf_chi)
lts_cae_chi <- log(ts_cae_chi)
lts_fch_chi <- log(ts_fch_chi)
# Creación de DataFrame para modelo
mod2_chi90_16 <- data.frame(cre_chi, scf_chi, cae_chi, fch_chi)
# Gráfica de las variables
plot(mod2_chi90_16, main="Gráfica de las Variables")
# Creación del modelo GLM
glm_mod3_chi <- glm(cre_chi ~ scf_chi + cae_chi + fch_chi, data=mod2_chi90_16)
summary(glm_mod3_chi)
##
## Call:
## glm(formula = cre_chi ~ scf_chi + cae_chi + fch_chi, data = mod2_chi90_16)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -4.4767 -1.0434 -0.2309 0.6162 5.9577
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 6.91203 3.02432 2.285 0.0318 *
## scf_chi 0.24124 0.20314 1.188 0.2471
## cae_chi 0.07139 0.05556 1.285 0.2116
## fch_chi -3.00905 1.36631 -2.202 0.0379 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for gaussian family taken to be 4.640085)
##
## Null deviance: 158.96 on 26 degrees of freedom
## Residual deviance: 106.72 on 23 degrees of freedom
## AIC: 123.73
##
## Number of Fisher Scoring iterations: 2
# Prueba de significancia
Psignificacion <- function(alpha, modelo) {
TabChi <- qchisq(1-alpha, modelo$df.residual)
Pv <- 1 - pchisq(modelo$deviance, modelo$df.residual)
R2 <- (1 - modelo$deviance/modelo$null.deviance) * 100
print(c("Chi-cuadrado", "p-value", "R^2"))
print(c(TabChi, Pv, R2))
}
Psignificacion(0.05, glm_mod3_chi)
## [1] "Chi-cuadrado" "p-value" "R^2"
## [1] 3.517246e+01 9.522383e-13 3.286116e+01
# Creación del modelo MCO
mco_mod2_chi <- lm(lts_cre_chi ~ lts_scf_chi + lts_cae_chi + lts_fch_chi)
summary(mco_mod2_chi)
##
## Call:
## lm(formula = lts_cre_chi ~ lts_scf_chi + lts_cae_chi + lts_fch_chi)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.63352 -0.09852 -0.04321 0.06387 0.63504
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.0117 1.8643 -1.079 0.2917
## lts_scf_chi 1.1255 0.8165 1.379 0.1813
## lts_cae_chi 0.3744 0.2916 1.284 0.2119
## lts_fch_chi -1.0924 0.4698 -2.325 0.0292 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2383 on 23 degrees of freedom
## Multiple R-squared: 0.3438, Adjusted R-squared: 0.2583
## F-statistic: 4.018 on 3 and 23 DF, p-value: 0.01952
# Prueba de autocorrelación Durbin-Watson
dwtest(mco_mod2_chi)
##
## Durbin-Watson test
##
## data: mco_mod2_chi
## DW = 1.1288, p-value = 0.001456
## alternative hypothesis: true autocorrelation is greater than 0
# Prueba de autocorrelación Breusch-Godfrey
bgtest(mco_mod2_chi)
##
## Breusch-Godfrey test for serial correlation of order up to 1
##
## data: mco_mod2_chi
## LM test = 2.679, df = 1, p-value = 0.1017
mco_mod2_chi <- lm(lts_cre_chi ~ lts_scf_chi + lts_cae_chi + lts_fch_chi, data=mod2_chi90_16)
# Resumen del modelo
summary(mco_mod2_chi)
##
## Call:
## lm(formula = lts_cre_chi ~ lts_scf_chi + lts_cae_chi + lts_fch_chi,
## data = mod2_chi90_16)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.63352 -0.09852 -0.04321 0.06387 0.63504
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.0117 1.8643 -1.079 0.2917
## lts_scf_chi 1.1255 0.8165 1.379 0.1813
## lts_cae_chi 0.3744 0.2916 1.284 0.2119
## lts_fch_chi -1.0924 0.4698 -2.325 0.0292 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2383 on 23 degrees of freedom
## Multiple R-squared: 0.3438, Adjusted R-squared: 0.2583
## F-statistic: 4.018 on 3 and 23 DF, p-value: 0.01952
# Errores estándar robustos (White)
library(sandwich)
coeftest(mco_mod2_chi, vcov = vcovHC(mco_mod2_chi, type="HC0"))
##
## t test of coefficients:
##
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.01171 2.81101 -0.7157 0.48141
## lts_scf_chi 1.12555 1.01991 1.1036 0.28119
## lts_cae_chi 0.37440 0.23800 1.5732 0.12934
## lts_fch_chi -1.09239 0.47039 -2.3223 0.02942 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Prueba de heterocedasticidad con Breusch-Pagan
bptest(mco_mod2_chi)
##
## studentized Breusch-Pagan test
##
## data: mco_mod2_chi
## BP = 10.533, df = 3, p-value = 0.01454
# Corrección de heterocedasticidad con errores estándar robustos de White
robust_se <- vcovHC(mco_mod2_chi, type="HC0")
coeftest(mco_mod2_chi, vcov=robust_se)
##
## t test of coefficients:
##
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.01171 2.81101 -0.7157 0.48141
## lts_scf_chi 1.12555 1.01991 1.1036 0.28119
## lts_cae_chi 0.37440 0.23800 1.5732 0.12934
## lts_fch_chi -1.09239 0.47039 -2.3223 0.02942 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Prueba de normalidad de residuos
shapiro.test(residuals(mco_mod2_chi))
##
## Shapiro-Wilk normality test
##
## data: residuals(mco_mod2_chi)
## W = 0.89404, p-value = 0.009798
# Q-Q plot de residuos
qqnorm(residuals(mco_mod2_chi))
qqline(residuals(mco_mod2_chi), col="red")
#Modelo con errores autorregresivos (AR1)
library(nlme)
modelo_gls <- gls(lts_cre_chi ~ lts_scf_chi + lts_cae_chi + lts_fch_chi,
data = mod2_chi90_16,
correlation = corAR1(form = ~1))
summary(modelo_gls)
## Generalized least squares fit by REML
## Model: lts_cre_chi ~ lts_scf_chi + lts_cae_chi + lts_fch_chi
## Data: mod2_chi90_16
## AIC BIC logLik
## -5.331186 1.48178 8.665593
##
## Correlation Structure: AR(1)
## Formula: ~1
## Parameter estimate(s):
## Phi
## 1
##
## Coefficients:
## Value Std.Error t-value p-value
## (Intercept) -3.630014 13857.702 -0.0002619 0.9998
## lts_scf_chi 0.217329 0.631 0.3442790 0.7338
## lts_cae_chi 0.907740 0.396 2.2938871 0.0313
## lts_fch_chi 1.340587 0.577 2.3226194 0.0294
##
## Correlation:
## (Intr) lts_s_ lts_c_
## lts_scf_chi 0.000
## lts_cae_chi 0.000 -0.396
## lts_fch_chi 0.000 -0.311 0.131
##
## Standardized residuals:
## Min Q1 Med Q3 Max
## -4.596817e-05 -1.112684e-05 1.200131e-05 4.599837e-05 8.079256e-05
##
## Residual standard error: 13857.7
## Degrees of freedom: 27 total; 23 residual
#(A) Corrección de Escala
mod2_chi90_16_scaled <- as.data.frame(scale(mod2_chi90_16))
modelo_gls_scaled <- gls(lts_cre_chi ~ lts_scf_chi + lts_cae_chi + lts_fch_chi,
data = mod2_chi90_16_scaled,
correlation = corAR1(form = ~1))
summary(modelo_gls_scaled)
## Generalized least squares fit by REML
## Model: lts_cre_chi ~ lts_scf_chi + lts_cae_chi + lts_fch_chi
## Data: mod2_chi90_16_scaled
## AIC BIC logLik
## -5.331186 1.48178 8.665593
##
## Correlation Structure: AR(1)
## Formula: ~1
## Parameter estimate(s):
## Phi
## 1
##
## Coefficients:
## Value Std.Error t-value p-value
## (Intercept) -3.630014 13857.702 -0.0002619 0.9998
## lts_scf_chi 0.217329 0.631 0.3442790 0.7338
## lts_cae_chi 0.907740 0.396 2.2938871 0.0313
## lts_fch_chi 1.340587 0.577 2.3226194 0.0294
##
## Correlation:
## (Intr) lts_s_ lts_c_
## lts_scf_chi 0.000
## lts_cae_chi 0.000 -0.396
## lts_fch_chi 0.000 -0.311 0.131
##
## Standardized residuals:
## Min Q1 Med Q3 Max
## -4.596817e-05 -1.112684e-05 1.200131e-05 4.599837e-05 8.079256e-05
##
## Residual standard error: 13857.7
## Degrees of freedom: 27 total; 23 residual
# Diferencias temporales
mod2_chi90_16_diff <- data.frame(
diff(lts_cre_chi),
diff(lts_scf_chi),
diff(lts_cae_chi),
diff(lts_fch_chi)
)
# Modelo MCO con diferencias
mco_diff <- lm(diff(lts_cre_chi) ~ diff(lts_scf_chi) + diff(lts_cae_chi) + diff(lts_fch_chi),
data = mod2_chi90_16_diff)
summary(mco_diff)
##
## Call:
## lm(formula = diff(lts_cre_chi) ~ diff(lts_scf_chi) + diff(lts_cae_chi) +
## diff(lts_fch_chi), data = mod2_chi90_16_diff)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.32719 -0.07477 0.01658 0.06492 0.59737
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.04689 0.04183 -1.121 0.2744
## diff(lts_scf_chi) 0.34100 0.63738 0.535 0.5980
## diff(lts_cae_chi) 0.95157 0.39547 2.406 0.0250 *
## diff(lts_fch_chi) 1.59867 0.61844 2.585 0.0169 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1875 on 22 degrees of freedom
## Multiple R-squared: 0.3994, Adjusted R-squared: 0.3175
## F-statistic: 4.877 on 3 and 22 DF, p-value: 0.009501
# Transformación Box-Cox
# Cargar librerías necesarias
library(MASS) # Para Box-Cox
library(car) # Para pruebas de normalidad
library(ggplot2) # Para visualización
## Warning: package 'ggplot2' was built under R version 4.2.3
# Ajustar el modelo de regresión lineal
modelo <- lm(lts_cre_chi ~ lts_scf_chi + lts_cae_chi + lts_fch_chi, data = mod2_chi90_16)
# Extraer los residuos del modelo
residuos <- residuals(modelo)
# Prueba de Shapiro-Wilk
shapiro.test(residuos)
##
## Shapiro-Wilk normality test
##
## data: residuos
## W = 0.89404, p-value = 0.009798
# Gráfico Q-Q
qqnorm(residuos)
qqline(residuos, col = "red")
# Histograma de residuos
ggplot(data = data.frame(residuos), aes(x = residuos)) +
geom_histogram(bins = 15, fill = "blue", alpha = 0.5) +
labs(title = "Histograma de Residuos", x = "Residuos", y = "Frecuencia")
#Transformación Box-Cox
library(MASS)
# Aplicar la transformación Box-Cox
boxcox(modelo, lambda = seq(-2, 2, by = 0.1))
# Determinar el mejor valor de lambda
lambda_optimo <- boxcox(modelo, lambda = seq(-2, 2, by = 0.1))$x[which.max(boxcox(modelo, lambda = seq(-2, 2, by = 0.1))$y)]
print(paste("El valor óptimo de lambda es:", lambda_optimo))
## [1] "El valor óptimo de lambda es: 2"
#Transformación Box-Cox
# Crear una nueva variable transformada
mod2_chi90_16$cre_chi_bc <- mod2_chi90_16$cre_chi^2
# Ajustar el modelo con la variable transformada
modelo_bc <- lm(cre_chi_bc ~ scf_chi + cae_chi + fch_chi, data = mod2_chi90_16)
# Resumen del nuevo modelo
summary(modelo_bc)
##
## Call:
## lm(formula = cre_chi_bc ~ scf_chi + cae_chi + fch_chi, data = mod2_chi90_16)
##
## Residuals:
## Min 1Q Median 3Q Max
## -66.665 -22.162 -11.681 9.298 125.645
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 67.893 62.591 1.085 0.2893
## scf_chi 3.570 4.204 0.849 0.4046
## cae_chi 1.448 1.150 1.259 0.2206
## fch_chi -53.103 28.277 -1.878 0.0731 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 44.58 on 23 degrees of freedom
## Multiple R-squared: 0.2813, Adjusted R-squared: 0.1875
## F-statistic: 3 on 3 and 23 DF, p-value: 0.05138
# Comprobar nuevamente la normalidad de los residuos
shapiro.test(residuals(modelo_bc))
##
## Shapiro-Wilk normality test
##
## data: residuals(modelo_bc)
## W = 0.86527, p-value = 0.002346
# Graficar Q-Q plot para verificar normalidad
qqnorm(residuals(modelo_bc))
qqline(residuals(modelo_bc), col = "red")
#Mejorar la normalidad de los residuos
mod2_chi90_16$cre_chi_sqrt <- sqrt(mod2_chi90_16$cre_chi) # Raíz cuadrada
mod2_chi90_16$cre_chi_log <- log(mod2_chi90_16$cre_chi) # Logaritmo
# Modelos con las nuevas transformaciones
modelo_sqrt <- lm(cre_chi_sqrt ~ scf_chi + cae_chi + fch_chi, data = mod2_chi90_16)
modelo_log <- lm(cre_chi_log ~ scf_chi + cae_chi + fch_chi, data = mod2_chi90_16)
# Prueba de normalidad en los residuos de ambos modelos
shapiro.test(residuals(modelo_sqrt))
##
## Shapiro-Wilk normality test
##
## data: residuals(modelo_sqrt)
## W = 0.91385, p-value = 0.02814
shapiro.test(residuals(modelo_log))
##
## Shapiro-Wilk normality test
##
## data: residuals(modelo_log)
## W = 0.89554, p-value = 0.01059
# Gráficos Q-Q para comparar normalidad
par(mfrow=c(1,2))
qqnorm(residuals(modelo_sqrt)); qqline(residuals(modelo_sqrt), col="red")
qqnorm(residuals(modelo_log)); qqline(residuals(modelo_log), col="blue")
par(mfrow=c(1,1))
#Mejorar la normalidad:
mod2_chi90_16$cre_chi_quart <- mod2_chi90_16$cre_chi^(1/4)
modelo_quart <- lm(cre_chi_quart ~ scf_chi + cae_chi + fch_chi, data = mod2_chi90_16)
# Prueba de normalidad
shapiro.test(residuals(modelo_quart))
##
## Shapiro-Wilk normality test
##
## data: residuals(modelo_quart)
## W = 0.90806, p-value = 0.02055
# Gráfico Q-Q
qqnorm(residuals(modelo_quart)); qqline(residuals(modelo_quart), col="blue")
#Modelo Robusto
library(MASS)
modelo_robusto <- rlm(cre_chi ~ scf_chi + cae_chi + fch_chi, data = mod2_chi90_16)
summary(modelo_robusto)
##
## Call: rlm(formula = cre_chi ~ scf_chi + cae_chi + fch_chi, data = mod2_chi90_16)
## Residuals:
## Min 1Q Median 3Q Max
## -4.41078 -0.72620 -0.04288 0.96010 6.08331
##
## Coefficients:
## Value Std. Error t value
## (Intercept) 7.2695 2.1815 3.3323
## scf_chi 0.1888 0.1465 1.2884
## cae_chi 0.0768 0.0401 1.9154
## fch_chi -2.6585 0.9856 -2.6975
##
## Residual standard error: 1.123 on 23 degrees of freedom
library(nortest)
ad.test(residuals(modelo_quart))
##
## Anderson-Darling normality test
##
## data: residuals(modelo_quart)
## A = 1.0063, p-value = 0.009967
mod2_chi90_16$cre_chi_fifth <- mod2_chi90_16$cre_chi^(1/5)
modelo_fifth <- lm(cre_chi_fifth ~ scf_chi + cae_chi + fch_chi, data = mod2_chi90_16)
shapiro.test(residuals(modelo_fifth)) # Revisa si mejora
##
## Shapiro-Wilk normality test
##
## data: residuals(modelo_fifth)
## W = 0.90601, p-value = 0.01841