Paquetes utilizados para realizar el examen:
library(stargazer)
##
## Please cite as:
## Hlavac, Marek (2018). stargazer: Well-Formatted Regression and Summary Statistics Tables.
## R package version 5.2.2. https://CRAN.R-project.org/package=stargazer
library(lmtest)
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
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:base':
##
## format.pval, units
library(readxl)
library(sandwich)
library(AER)
## Loading required package: car
## Loading required package: carData
Datos para realizar el examen:
<- read.csv("C:\\Users\\Yudith Lopez Gaitan\\Downloads\\P1.csv")
P1 <- read_excel("C:/Users/Yudith Lopez Gaitan/Downloads/P2.xlsx")
P2 <- read.csv("C:\\Users\\Yudith Lopez Gaitan\\Downloads\\P3.csv", header = TRUE, sep = ";") P3
<- lm(log(JV) ~ log(U), P1)
modelo_uno stargazer(modelo_uno, type = "text")
##
## ===============================================
## Dependent variable:
## ---------------------------
## log(JV)
## -----------------------------------------------
## log(U) -1.612***
## (0.156)
##
## Constant 3.503***
## (0.283)
##
## -----------------------------------------------
## Observations 24
## R2 0.830
## Adjusted R2 0.822
## Residual Std. Error 0.285 (df = 22)
## F Statistic 107.361*** (df = 1; 22)
## ===============================================
## Note: *p<0.1; **p<0.05; ***p<0.01
Un intervalo de confianza está definido por:
\[Prob\left[ \beta_i-t_{\alpha/2}\cdot se(\hat{\beta}_i)\leq\hat{\beta}_i \leq\beta_i+t_{\alpha/2}\cdot se(\hat{\beta}_i)\right]=1-\alpha\]
Teniendo en cuenta que para este caso \(H_o:\beta_i=0\) y \(H_a:\beta_i\neq 0\). Tambien, se pueden obtener los errores estándar de la regresión anterior:
#El intervalo de B_0 será:
qt(0.995,22)*0.283
## [1] 0.797708
#El intervalo de B_1 será:
qt(0.995,22)*0.156
## [1] 0.4397259
Por ende,
Intervalo \(\beta_0\):
\[Prob\left[-0.7977 \leq \hat{\beta}_0 \leq0.7977\right]=0.99\]
Intervalo \(\beta_1\)
\[Prob\left[-0.4397 \leq \hat{\beta}_1 \leq0.4397\right]=0.99\]
Si no se desea rechazar la hipótesis nula, el coeficiente se debe ubicar dentro del intervalo de confianza en una zona de no rechazo con un nivel de confianza del 99%. Si no se encuentra dentro de esa área se rechaza y se puede concluir que el coeficiente es significativo.
Nosotros sabemos que \(\beta_0=3.503\) por lo que cae fuera de la zona de no rechazo y podemos decir que se rechaza la hipótesis nula de que \(H_0:\beta_0=0\) por lo que es estadísticamente diferente de 0. Asimismo, sabemos que \(\beta_1=−1.612\) por lo que también cae en la zona de rechazo, es decir se encuentra fuera de la zona de no rechazo, así que con un nivel de confianza del 99% rechazamos la hipótesis nula.
Este se calcula de la siguiente manera:
\[d=\frac{\sum\hat{u}_t^2+\sum\hat{u}_{t-1}^2-2\sum\hat{u}_t\hat{u}_{t-1}}{\sum\hat{u}_{t}^2}\]
dwtest(modelo_uno)
##
## Durbin-Watson test
##
## data: modelo_uno
## DW = 1.0896, p-value = 0.008821
## alternative hypothesis: true autocorrelation is greater than 0
En este caso, solo hay un variable (explicativa excluyendo el intercepto). Por lo tanto: \(d_L=1.29\) y \(d_U=1.45\). Por lo tanto se rechaza \(H_o\), porque existe evidencia de autocorrelación positiva de primer orden.
Aplicando el test de Breuch-Goodfrey
bgtest(modelo_uno, order = 1)
##
## Breusch-Godfrey test for serial correlation of order up to 1
##
## data: modelo_uno
## LM test = 4.6503, df = 1, p-value = 0.03105
bgtest(modelo_uno, order = 2)
##
## Breusch-Godfrey test for serial correlation of order up to 2
##
## data: modelo_uno
## LM test = 4.8224, df = 2, p-value = 0.08971
Esto nos dice que existe evidencia de autocorrelación de orden 1 y 2, ya que en ambos casos, se rechaza \(H_0\) bajo una significancia del 10%.
En los incisos c y d concluimos que existe autocorrelación, por lo que esto generaría un problema en las varianzas de los estimadores. La varianza de los estimadores será sobreestimada y eso nos va a llevar a que la pruebas t realizadas en el inciso b no tengan validez.
Teniendo en cuenta que el modelo AR(1) se define como:
\[e_t=ρ_et−1+v_t\]
Podremos saber que \(\hat{\rho}=0.447\) gracias a lo siguiente:
#El primer rezago de la serie de residuos y la serie de residuos será:
$resid <- resid(modelo_uno)
P1$resid_lag1 <- Lag(P1$resid, shift = 1)
P1$resid_lag1[1]<-0
P1
#El modelo con los residuos:
<-lm(resid ~ resid_lag1 - 1, P1)
modelo_zzstargazer(modelo_zz, type = "text")
##
## ===============================================
## Dependent variable:
## ---------------------------
## resid
## -----------------------------------------------
## resid_lag1 0.447**
## (0.191)
##
## -----------------------------------------------
## Observations 24
## R2 0.193
## Adjusted R2 0.158
## Residual Std. Error 0.250 (df = 23)
## F Statistic 5.505** (df = 1; 23)
## ===============================================
## Note: *p<0.1; **p<0.05; ***p<0.01
#p será:
<-modelo_zz$coefficients rho
Ahora, gracias a el valor obtenido se podrá transformar el modelo para aplicar MCG:
\[ρy_{t−1}=ρβ_0+ρβ_1x_{t−1}+ρu_{t−1}\]
Luego, restando el modelo original obtendremos:
\[y_t−ρy_{t−1}=β_0(1−ρ)+β_1(x_t−ρx_{t−1})+u_t−ρu_{t−1}\]
\[y^∗_t=\beta^∗_0+\beta_1x^∗_t+e_t\]
Sabiendo todo esto, podremos realizar el modelo de la siguiente manera:
$logU <- log(P1$U)
P1$logJV <- log(P1$JV)
P1
$constante_mccartney<-NA
P1$logU_mccartney<-NA
P1$logPV_mccartney<-NA
P1
$constant_mccartney[1]<-((1-rho^{2})^{0.5})
P1$logU_mccartney[1]<-((1-rho^{2})^{0.5})*P1$logU[1]
P1$logPV_mccartney[1]<-((1-rho^{2})^{0.5})*P1$logJV[1]
P1
$constant_mccartney[2:nrow(P1)] <- 1 - rho
P1
for (i in 2:nrow(P1)) {
$logU_mccartney[i] <- P1$logU[i] - rho*P1$logU[i-1]
P1
$logPV_mccartney[i] <- P1$logJV[i] - rho*P1$logJV[i-1]
P1
}
<- lm(logPV_mccartney ~ constant_mccartney + logU_mccartney - 1, P1)
modelo_autocorre_P1
stargazer(modelo_autocorre_P1, type = "text")
##
## ===============================================
## Dependent variable:
## ---------------------------
## logPV_mccartney
## -----------------------------------------------
## constant_mccartney 3.514***
## (0.238)
##
## logU_mccartney -1.616***
## (0.123)
##
## -----------------------------------------------
## Observations 24
## R2 0.909
## Adjusted R2 0.901
## Residual Std. Error 0.256 (df = 22)
## F Statistic 110.247*** (df = 2; 22)
## ===============================================
## Note: *p<0.1; **p<0.05; ***p<0.01
Luego que tenemos esto, nos falta verificar si este modelo tiene aucorrelación:
dwtest(modelo_autocorre_P1)
##
## Durbin-Watson test
##
## data: modelo_autocorre_P1
## DW = 2.0682, p-value = 0.604
## alternative hypothesis: true autocorrelation is greater than 0
bgtest(modelo_autocorre_P1, order = 1)
##
## Breusch-Godfrey test for serial correlation of order up to 1
##
## data: modelo_autocorre_P1
## LM test = 0.053385, df = 1, p-value = 0.8173
bgtest(modelo_autocorre_P1, order = 2)
##
## Breusch-Godfrey test for serial correlation of order up to 2
##
## data: modelo_autocorre_P1
## LM test = 0.086767, df = 2, p-value = 0.9575
Podemos observar que:
Para la prueba del Durbin Watson, el nuevo estadístico nos indica no autocorrelación. La prueba de Breusch-Godfrey tenemos un un p-value mayor a 0.5564, esto nos indica que vamos a rechazar que existe autocorrelación con un nivel de 95%. Por ende, la elasticidad de la vacantes laborales con respecto a la tasa de desempleo no cambio de manera significativa.
Podemos plantear este ejercicio por medio de la siguiente regresión:
\[ \text{log(avgprc)}= \beta_0 + \beta_1\text{mon} + \beta_2\text{tues} + \beta_3\text{wed} + \beta_4\text{thrus} + \beta_5\text{t} + \mu_t \]
Realizando esto, obtenemos la siguiente regresión:
<- lm(log(avgprc) ~ mon + tues + wed + thurs + t, P2)
modelo_dos stargazer(modelo_dos, type = "text")
##
## ===============================================
## Dependent variable:
## ---------------------------
## log(avgprc)
## -----------------------------------------------
## mon -0.010
## (0.129)
##
## tues -0.009
## (0.127)
##
## wed 0.038
## (0.126)
##
## thurs 0.091
## (0.126)
##
## t -0.004***
## (0.001)
##
## Constant -0.073
## (0.115)
##
## -----------------------------------------------
## Observations 97
## R2 0.085
## Adjusted R2 0.035
## Residual Std. Error 0.397 (df = 91)
## F Statistic 1.701 (df = 5; 91)
## ===============================================
## Note: *p<0.1; **p<0.05; ***p<0.01
Podemos observar que la serie de tendencia (t) es significativa a un p-value < 0.01, las otras variables no son significativas (según las pruebas t) para el modelo. Con esto llegamos a la conclusión de que los días no afectan la variabilidad del precio promedio que hay durante la semana.
Si añadimos las variables “wave2” y “wave3”, obtendremos:
\[ \text{log(avgprc)}= \beta_0 + \beta_1\text{mon} + \beta_2\text{tues} + \beta_3\text{wed} + \beta_4\text{thrus} + \beta_5\text{t} + + \beta_6 \text{wave2} + \beta_7 \text{wave3} + \mu_t \]
Entonces, obtenemos:
<- lm(log(avgprc) ~ mon + tues + wed + thurs + t + wave2 + wave3, P2)
modelo_dosb stargazer(modelo_dosb, type = "text")
##
## ===============================================
## Dependent variable:
## ---------------------------
## log(avgprc)
## -----------------------------------------------
## mon -0.018
## (0.114)
##
## tues -0.009
## (0.112)
##
## wed 0.050
## (0.112)
##
## thurs 0.123
## (0.111)
##
## t -0.001
## (0.001)
##
## wave2 0.091***
## (0.022)
##
## wave3 0.047**
## (0.021)
##
## Constant -0.920***
## (0.190)
##
## -----------------------------------------------
## Observations 97
## R2 0.309
## Adjusted R2 0.255
## Residual Std. Error 0.349 (df = 89)
## F Statistic 5.699*** (df = 7; 89)
## ===============================================
## Note: *p<0.1; **p<0.05; ***p<0.01
Con esto llegamos a la conclusión de que las variables integradas para este inciso son significativas a un 5%. Como esto sucede, sabemos que estas explican la variación de los precios en el mercado. Además, hay que tomar en cuenta que lo que podría estar ocasionando los aumentos de precios según el modelo puede ser las altas olas porque estas provocan un pésimo clima para los pescadores.
Podemos observar que:
\[wave2=\lambda_0+\lambda_1t_i\]
Además,
\[wave3=\delta_0+\delta_1t_i\] Entonces,
<- lm(wave2 ~ t, P2)
modelo_t2 <- lm(wave3 ~ t, P2)
modelo_t3 stargazer(modelo_t2, modelo_t3, type = "text")
##
## ==========================================================
## Dependent variable:
## ----------------------------
## wave2 wave3
## (1) (2)
## ----------------------------------------------------------
## t -0.019*** -0.024***
## (0.006) (0.007)
##
## Constant 6.008*** 6.209***
## (0.351) (0.368)
##
## ----------------------------------------------------------
## Observations 97 97
## R2 0.086 0.127
## Adjusted R2 0.077 0.118
## Residual Std. Error (df = 95) 1.717 1.798
## F Statistic (df = 1; 95) 8.991*** 13.842***
## ==========================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
Podemos observar que ambas variables son negativas y significativas al 5%, esto quiere decir que explican la variación de los precios en el mercado. Esto se puede deber a que conforme pasa el tiempo el oleaje puede ser más bajo mientras que la variable t continuará creciendo porque depende de la temporada.
<- resid(modelo_dosb)
residuo_modelox <- cor(residuo_modelox,P2$mon)
mon_corr <- cor(residuo_modelox,P2$tues)
tues_corr <- cor(residuo_modelox,P2$wed)
wed_corr <- cor(residuo_modelox,P2$thurs)
thurs_corr <- cor(residuo_modelox,P2$wave2)
wave2_corr <-cor(residuo_modelox,P2$wave3)
wave3_corr
<- c("mon", "tues", "wed", "thurs", "wave2", "wave3")
variables <- c(mon_corr, tues_corr, wed_corr, thurs_corr, wave2_corr, wave3_corr)
corr_var_res <- data.frame(variables, corr_var_res)
Corr_dos names(Corr_dos) <- c("Variable", "Correlación")
stargazer(Corr_dos, summary = F, type = "text")
##
## ======================
## Variable Correlación
## ----------------------
## 1 mon -0
## 2 tues -0
## 3 wed 0
## 4 thurs 0
## 5 wave2 0
## 6 wave3 0
## ----------------------
El oleaje del mar no puede ser explicado de manera sencilla, ya que este tiende a depender variables geográficas. Esto nos lleva a la conclusión de que todas la variables son exógenas.
dwtest(modelo_dosb)
##
## Durbin-Watson test
##
## data: modelo_dosb
## DW = 0.74523, p-value = 2.654e-12
## alternative hypothesis: true autocorrelation is greater than 0
Podemos ver que el p-value es menor a 5% (muy bajo), por lo que podemos rechazar la hipótesis nula, es decir que que sí hay autocerrelación. Ahora, se sabe que la correción de la matriz de Newey-West genera los mismos estimadores de MCO, pero los errores estándar serán diferentes:
<- sandwich::NeweyWest(modelo_dosb)
nwcov ::coeftest(modelo_dosb, vcov. = nwcov) lmtest
##
## t test of coefficients:
##
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.9202543 0.3908063 -2.3548 0.02073 *
## mon -0.0181580 0.0808454 -0.2246 0.82280
## tues -0.0085330 0.0835063 -0.1022 0.91884
## wed 0.0500407 0.0754772 0.6630 0.50905
## thurs 0.1225463 0.0555319 2.2068 0.02990 *
## t -0.0011575 0.0032460 -0.3566 0.72225
## wave2 0.0908906 0.0349925 2.5974 0.01099 *
## wave3 0.0473748 0.0213163 2.2225 0.02879 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
<- c(0.0217, 0.0208)
se_mco <- c(0.0349, 0.0213)
se_nw <- c("wave2", "wave3")
se_nombres
<- data.frame(se_nombres, se_mco, se_nw)
se_comp names(se_comp) <- c("Variable", "MCO", "Newey-West")
stargazer(se_comp, summary = F, type = "text", align = T, no.space = T, title = "Errores estándar", digits = 5)
##
## Errores estándar
## =============================
## Variable MCO Newey-West
## -----------------------------
## 1 wave2 0.02170 0.03490
## 2 wave3 0.02080 0.02130
## -----------------------------
Se puede observar que las varianzas del Newey-West son mayores, lo cual no tiende a ser lo más normal ante un problema de autocorrelación, ante este problema podríamos estar incurriendo en inferencias incorrectas pues hay valores dentro de la zona de no rechazo que se están rechazando por lo que aunque la varianza sea la mínima no es a la que se desea llegar para concluir acerca de los vectores.
En conclusión, sabemos que el problema de autocorrelación hace que los estimadores sean consistentes, lineales e insesgados, pero no ineficientes. Es por esto que \(\hat{\beta}_i\) no será MELI. Además, hay una probabilidad alta de cometer el error de tipo 1 ya que la varianza es la menor.
Nosotros sabemos que el estadístico F plantea la hípotesis nula de que todas las variables explicativas no expliquen estadísticamente el modelo presentado, por lo tanto:
\[H_o: \beta_1=\beta_2=\beta_3=0\]
Entonces,
<- lm(y ~ xper + cap + lab, P3)
modelo_tres stargazer(modelo_tres, type = "text")
##
## ===============================================
## Dependent variable:
## ---------------------------
## y
## -----------------------------------------------
## xper 1.451**
## (0.633)
##
## cap 0.440***
## (0.118)
##
## lab 0.237**
## (0.100)
##
## Constant 17.948*
## (10.503)
##
## -----------------------------------------------
## Observations 75
## R2 0.564
## Adjusted R2 0.546
## Residual Std. Error 27.490 (df = 71)
## F Statistic 30.652*** (df = 3; 71)
## ===============================================
## Note: *p<0.1; **p<0.05; ***p<0.01
Gracias a el p-value sabemos que es singinicativa con un nivel del 0.01, por ello la regresión conjunta explica la variabilidad de \(y_i\).
Creando tres intervalos de confianza:
<- data.frame("xper"=c(10,20,30),"cap"= mean(P3$cap),"lab"= mean(P3$lab))
predib <-predict.lm(modelo_tres, predib, level = 0.95, interval = c("confidence"))
predib_result stargazer(predib_result, type = "text")
##
## ========================
## fit lwr upr
## ------------------------
## 1 90.703 82.699 98.707
## 2 105.218 95.229 115.206
## 3 119.732 98.416 141.048
## ------------------------
Podemos observar que entre más años de experiencia tengan, la producción entre los intervalos será mayor. Además, es importante recalcar que la producción tiene una correlación positiva con los años de experiencia.
Los años de experiencia de un individuos se relaciona directamente con la edad (entre más años tenga es más posible que tenga más años de experiencia), ante el problema de la correlación de \(xper\) con \(u_i\), se debe a eso en concreto.
Recordando lo visto en clase, para que un instrumento sea válido debe ser relevante y debe tener correlación con la variable a la que será instrumento.
\[Corr(z_i,x_i)\neq0\]
Además, debe ser exógeno, no debe correlacionarse con los residuos.
\[Corr(z_i,u_i)=0\]
Entonces, tenemos:
<- resid(modelo_tres)
resid_tres <- cor(P3$age,P3$xper)
VI_relevan <- cor(P3$age,resid_tres)
VI_exogena stargazer(VI_relevan, type = "text", title = "Correlación contra variable")
##
## Correlación contra variable
## =====
## 0.322
## -----
stargazer(VI_exogena, type = "text", title = "Correlación contra residuos")
##
## Correlación contra residuos
## =====
## 0.235
## -----
Con esto concluimos que la correlación entre los residuos es muy baja significativamente, pero aún así comoa mayor edad se espera que las personas tengan más experiencia el intrumento se tomará como válido.
Para este caso, vamos a hacerlo de dos maneras:
<- ivreg(formula = y ~ xper + cap + lab | cap + lab + age, data=P3)
modelo_tres_IV
# Ahora, realizando esto mismo pero de manual
<- lm(xper ~ cap + lab + age, P3)
modelo_tres_endo <- fitted(modelo_tres_endo)
xper_hat <- lm(y ~ xper_hat + cap +lab, P3)
modelo_tres_VIM
stargazer(modelo_tres, modelo_tres_IV, modelo_tres_VIM, type = "text", title = "Modelo original, modelo corregido con ivreg Y modelo corregido manualmente")
##
## Modelo original, modelo corregido con ivreg Y modelo corregido manualmente
## ==============================================================
## Dependent variable:
## --------------------------------
## y
## OLS instrumental OLS
## variable
## (1) (2) (3)
## --------------------------------------------------------------
## xper 1.451** 5.123**
## (0.633) (2.201)
##
## xper_hat 5.123***
## (1.778)
##
## cap 0.440*** 0.330** 0.330**
## (0.118) (0.156) (0.126)
##
## lab 0.237** 0.241* 0.241**
## (0.100) (0.121) (0.098)
##
## Constant 17.948* -24.748 -24.748
## (10.503) (27.161) (21.942)
##
## --------------------------------------------------------------
## Observations 75 75 75
## R2 0.564 0.358 0.581
## Adjusted R2 0.546 0.331 0.563
## Residual Std. Error (df = 71) 27.490 33.369 26.957
## F Statistic (df = 3; 71) 30.652*** 32.823***
## ==============================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
Podemos observar cómo los coeficientes de las variables son distintos, en el modelo de variable instrumental xper_hat toma un valor significativo comparado con los betas estimados por MCO. En cambio, para cap y lab sus errores estándar casi no varían. Lo único que cambia es el peso que van obteniendo para explicar la variabilidad de \(y_i\) donde cap termina aumentando y lab disminuyendo. También se ve como realizando el método de regresión utilizando la variable instrumental aumenta el \(R^2\) por lo que este modelo explica mejor la variabilidad de de \(y_i\), Por último, la constante se vuelve negativa y aumenta en términos de valor absoluto, lo que nos indica es que si los años de experiencia promedio estimados sujeto a edad promedio de los trabajadores el índice indicador de la cantidad de capital utilizada además del índice indicador de la cantidad de trabajo utilizado resultan en 0, es por ello la eficiencia relativa de producción termina siendo negativa.
<- predict(modelo_tres_IV, newdata = data.frame(lab=mean(P3$lab), cap=mean(P3$cap), age=mean(P3$age), xper=10), interval = c("confidence"))
edad1 <- predict(modelo_tres_IV, newdata = data.frame(lab=mean(P3$lab), cap=mean(P3$cap), age=mean(P3$age), xper=20), interval = c("confidence"))
edad2 <- predict(modelo_tres_IV, newdata = data.frame(lab=mean(P3$lab), cap=mean(P3$cap), age=mean(P3$age), xper=30), interval = c("confidence")) edad3
Notemos que la variable de la edad es relevante en la explicación de la variabilidad de yi, esto porque al contar con el mismo promedio de edad se marca significativamente el papel que juega la variable xper, ya cuán mayor sea esta, mayor será la predicción obtenida de yi. Como conclusión obtenemos que, si se tiene la misma edad, mismo cap y mismo lab entre mayor sea xper mayor va a ser la eficiencia en la producción de vino de las fincas. Esto sucede porque las variables mencionadas están correlacionadas positivamente.
head(edad1)
## 1
## 76.45899
head(edad2)
## 1
## 127.6852
head(edad3)
## 1
## 178.9115
Realizando la regresión de la primera etapa cuando una predicción de cuanto va a ser los años de experiencia promedio de los administradores cuando se tiene que la edad promedio es cero:
<- 0
edadd
<- predict(modelo_tres_endo, newdata = data.frame(lab=(P3$lab), cap=(P3$cap), age=edadd, xper=10),
pre_10 se.fit = TRUE, level = 0.95, interval = c("confidence"))
<- predict(modelo_tres_endo, newdata = data.frame(lab=(P3$lab), cap=(P3$cap), age=edadd, xper=20),
pre_20 se.fit = TRUE, level = 0.95, interval = c("confidence"))
<- predict(modelo_tres_endo, newdata = data.frame(lab=(P3$lab), cap=(P3$cap), age=edadd, xper=30),
pre_30 se.fit = TRUE, level = 0.95, interval = c("confidence"))
head(pre_10)
## $fit
## fit lwr upr
## 1 4.240558 -1.00395843 9.485074
## 2 6.476685 1.85274521 11.100624
## 3 7.404305 2.08102121 12.727588
## 4 4.705266 -0.41055230 9.821085
## 5 6.522222 1.80809725 11.236348
## 6 7.429309 2.80259335 12.056024
## 7 6.044339 1.14786111 10.940816
## 8 6.530608 1.91462891 11.146587
## 9 6.615871 1.88163644 11.350105
## 10 5.751238 1.00199130 10.500485
## 11 6.538799 1.90487045 11.172728
## 12 4.757334 -0.30497238 9.819641
## 13 7.515841 2.70659944 12.325083
## 14 9.333924 4.20921540 14.458632
## 15 4.948272 -0.62123537 10.517779
## 16 6.445936 1.77311386 11.118758
## 17 7.791769 3.18005298 12.403485
## 18 6.592612 1.89188362 11.293340
## 19 5.576835 0.68996369 10.463707
## 20 7.123267 2.54572880 11.700806
## 21 6.025266 1.28824704 10.762286
## 22 7.725536 2.64323517 12.807838
## 23 6.697423 2.02654518 11.368301
## 24 8.897131 4.12051277 13.673748
## 25 9.179348 3.99869357 14.360003
## 26 7.618566 3.03262904 12.204503
## 27 5.743695 0.66039868 10.826991
## 28 7.767499 3.11287489 12.422124
## 29 8.683890 3.79461619 13.573165
## 30 7.230243 2.64611810 11.814368
## 31 7.126689 2.54352244 11.709855
## 32 8.078451 3.38799867 12.768904
## 33 4.505999 -0.70508445 9.717082
## 34 5.364404 0.52293909 10.205869
## 35 8.113140 3.44157279 12.784707
## 36 6.644249 1.99192690 11.296571
## 37 9.164653 4.30286074 14.026445
## 38 4.267665 -1.27353234 9.808863
## 39 6.761631 2.16712030 11.356142
## 40 7.347467 2.77148237 11.923451
## 41 8.349144 3.66196633 13.036322
## 42 6.576705 1.74690650 11.406503
## 43 6.804966 1.93875566 11.671176
## 44 7.100902 2.48370547 11.718098
## 45 7.602713 2.89938615 12.306040
## 46 4.947963 -0.03328176 9.929207
## 47 6.694208 1.99459141 11.393824
## 48 5.224958 0.22874650 10.221170
## 49 6.494056 1.84852804 11.139585
## 50 7.069247 2.18471219 11.953782
## 51 5.569420 -0.84543445 11.984275
## 52 4.768318 -0.34449649 9.881132
## 53 9.362059 4.44391306 14.280204
## 54 6.433093 1.77191508 11.094272
## 55 5.793980 0.44219862 11.145762
## 56 5.989367 0.36317627 11.615557
## 57 8.248342 3.40073909 13.095946
## 58 7.257039 2.42634499 12.087733
## 59 7.675477 3.02415511 12.326799
## 60 8.029888 2.73332290 13.326454
## 61 6.214982 1.55733669 10.872627
## 62 6.103263 1.42731992 10.779205
## 63 7.190814 2.55790710 11.823721
## 64 6.360102 1.68988227 11.030321
## 65 5.357807 -0.23373361 10.949348
## 66 4.208170 -1.04889215 9.465233
## 67 7.017108 2.18023182 11.853984
## 68 5.131619 -0.33602338 10.599262
## 69 5.940752 0.74460213 11.136902
## 70 7.687221 3.09517046 12.279271
## 71 5.718821 0.95817921 10.479464
## 72 9.034312 4.19940623 13.869218
## 73 8.185622 3.01138673 13.359858
## 74 9.605880 3.76919434 15.442565
## 75 4.738593 -0.37299420 9.850180
##
## $se.fit
## 1 2 3 4 5 6 7 8
## 2.630223 2.318992 2.669727 2.565679 2.364222 2.320385 2.455675 2.315000
## 9 10 11 12 13 14 15 16
## 2.374307 2.381836 2.324002 2.538842 2.411925 2.570137 2.793212 2.343508
## 17 18 19 20 21 22 23 24
## 2.312862 2.357503 2.450858 2.295721 2.375704 2.548869 2.342533 2.395563
## 25 26 27 28 29 30 31 32
## 2.598195 2.299933 2.549368 2.334381 2.452063 2.299025 2.298544 2.352350
## 33 34 35 36 37 38 39 40
## 2.613456 2.428086 2.342879 2.333227 2.438280 2.779015 2.304233 2.294942
## 41 42 43 44 45 46 47 48
## 2.350708 2.422234 2.440496 2.315610 2.358807 2.498188 2.356946 2.505694
## 49 50 51 52 53 54 55 56
## 2.329820 2.449686 3.217170 2.564172 2.466542 2.337668 2.684019 2.821640
## 57 58 59 60 61 62 63 64
## 2.431164 2.422684 2.332725 2.656327 2.335896 2.345073 2.323490 2.342203
## 65 66 67 68 69 70 71 72
## 2.804263 2.636515 2.425784 2.742125 2.605967 2.302999 2.387551 2.424796
## 73 74 75
## 2.594976 2.927207 2.563557
##
## $df
## [1] 71
##
## $residual.scale
## [1] 4.827578
head(pre_20)
## $fit
## fit lwr upr
## 1 4.240558 -1.00395843 9.485074
## 2 6.476685 1.85274521 11.100624
## 3 7.404305 2.08102121 12.727588
## 4 4.705266 -0.41055230 9.821085
## 5 6.522222 1.80809725 11.236348
## 6 7.429309 2.80259335 12.056024
## 7 6.044339 1.14786111 10.940816
## 8 6.530608 1.91462891 11.146587
## 9 6.615871 1.88163644 11.350105
## 10 5.751238 1.00199130 10.500485
## 11 6.538799 1.90487045 11.172728
## 12 4.757334 -0.30497238 9.819641
## 13 7.515841 2.70659944 12.325083
## 14 9.333924 4.20921540 14.458632
## 15 4.948272 -0.62123537 10.517779
## 16 6.445936 1.77311386 11.118758
## 17 7.791769 3.18005298 12.403485
## 18 6.592612 1.89188362 11.293340
## 19 5.576835 0.68996369 10.463707
## 20 7.123267 2.54572880 11.700806
## 21 6.025266 1.28824704 10.762286
## 22 7.725536 2.64323517 12.807838
## 23 6.697423 2.02654518 11.368301
## 24 8.897131 4.12051277 13.673748
## 25 9.179348 3.99869357 14.360003
## 26 7.618566 3.03262904 12.204503
## 27 5.743695 0.66039868 10.826991
## 28 7.767499 3.11287489 12.422124
## 29 8.683890 3.79461619 13.573165
## 30 7.230243 2.64611810 11.814368
## 31 7.126689 2.54352244 11.709855
## 32 8.078451 3.38799867 12.768904
## 33 4.505999 -0.70508445 9.717082
## 34 5.364404 0.52293909 10.205869
## 35 8.113140 3.44157279 12.784707
## 36 6.644249 1.99192690 11.296571
## 37 9.164653 4.30286074 14.026445
## 38 4.267665 -1.27353234 9.808863
## 39 6.761631 2.16712030 11.356142
## 40 7.347467 2.77148237 11.923451
## 41 8.349144 3.66196633 13.036322
## 42 6.576705 1.74690650 11.406503
## 43 6.804966 1.93875566 11.671176
## 44 7.100902 2.48370547 11.718098
## 45 7.602713 2.89938615 12.306040
## 46 4.947963 -0.03328176 9.929207
## 47 6.694208 1.99459141 11.393824
## 48 5.224958 0.22874650 10.221170
## 49 6.494056 1.84852804 11.139585
## 50 7.069247 2.18471219 11.953782
## 51 5.569420 -0.84543445 11.984275
## 52 4.768318 -0.34449649 9.881132
## 53 9.362059 4.44391306 14.280204
## 54 6.433093 1.77191508 11.094272
## 55 5.793980 0.44219862 11.145762
## 56 5.989367 0.36317627 11.615557
## 57 8.248342 3.40073909 13.095946
## 58 7.257039 2.42634499 12.087733
## 59 7.675477 3.02415511 12.326799
## 60 8.029888 2.73332290 13.326454
## 61 6.214982 1.55733669 10.872627
## 62 6.103263 1.42731992 10.779205
## 63 7.190814 2.55790710 11.823721
## 64 6.360102 1.68988227 11.030321
## 65 5.357807 -0.23373361 10.949348
## 66 4.208170 -1.04889215 9.465233
## 67 7.017108 2.18023182 11.853984
## 68 5.131619 -0.33602338 10.599262
## 69 5.940752 0.74460213 11.136902
## 70 7.687221 3.09517046 12.279271
## 71 5.718821 0.95817921 10.479464
## 72 9.034312 4.19940623 13.869218
## 73 8.185622 3.01138673 13.359858
## 74 9.605880 3.76919434 15.442565
## 75 4.738593 -0.37299420 9.850180
##
## $se.fit
## 1 2 3 4 5 6 7 8
## 2.630223 2.318992 2.669727 2.565679 2.364222 2.320385 2.455675 2.315000
## 9 10 11 12 13 14 15 16
## 2.374307 2.381836 2.324002 2.538842 2.411925 2.570137 2.793212 2.343508
## 17 18 19 20 21 22 23 24
## 2.312862 2.357503 2.450858 2.295721 2.375704 2.548869 2.342533 2.395563
## 25 26 27 28 29 30 31 32
## 2.598195 2.299933 2.549368 2.334381 2.452063 2.299025 2.298544 2.352350
## 33 34 35 36 37 38 39 40
## 2.613456 2.428086 2.342879 2.333227 2.438280 2.779015 2.304233 2.294942
## 41 42 43 44 45 46 47 48
## 2.350708 2.422234 2.440496 2.315610 2.358807 2.498188 2.356946 2.505694
## 49 50 51 52 53 54 55 56
## 2.329820 2.449686 3.217170 2.564172 2.466542 2.337668 2.684019 2.821640
## 57 58 59 60 61 62 63 64
## 2.431164 2.422684 2.332725 2.656327 2.335896 2.345073 2.323490 2.342203
## 65 66 67 68 69 70 71 72
## 2.804263 2.636515 2.425784 2.742125 2.605967 2.302999 2.387551 2.424796
## 73 74 75
## 2.594976 2.927207 2.563557
##
## $df
## [1] 71
##
## $residual.scale
## [1] 4.827578
head(pre_30)
## $fit
## fit lwr upr
## 1 4.240558 -1.00395843 9.485074
## 2 6.476685 1.85274521 11.100624
## 3 7.404305 2.08102121 12.727588
## 4 4.705266 -0.41055230 9.821085
## 5 6.522222 1.80809725 11.236348
## 6 7.429309 2.80259335 12.056024
## 7 6.044339 1.14786111 10.940816
## 8 6.530608 1.91462891 11.146587
## 9 6.615871 1.88163644 11.350105
## 10 5.751238 1.00199130 10.500485
## 11 6.538799 1.90487045 11.172728
## 12 4.757334 -0.30497238 9.819641
## 13 7.515841 2.70659944 12.325083
## 14 9.333924 4.20921540 14.458632
## 15 4.948272 -0.62123537 10.517779
## 16 6.445936 1.77311386 11.118758
## 17 7.791769 3.18005298 12.403485
## 18 6.592612 1.89188362 11.293340
## 19 5.576835 0.68996369 10.463707
## 20 7.123267 2.54572880 11.700806
## 21 6.025266 1.28824704 10.762286
## 22 7.725536 2.64323517 12.807838
## 23 6.697423 2.02654518 11.368301
## 24 8.897131 4.12051277 13.673748
## 25 9.179348 3.99869357 14.360003
## 26 7.618566 3.03262904 12.204503
## 27 5.743695 0.66039868 10.826991
## 28 7.767499 3.11287489 12.422124
## 29 8.683890 3.79461619 13.573165
## 30 7.230243 2.64611810 11.814368
## 31 7.126689 2.54352244 11.709855
## 32 8.078451 3.38799867 12.768904
## 33 4.505999 -0.70508445 9.717082
## 34 5.364404 0.52293909 10.205869
## 35 8.113140 3.44157279 12.784707
## 36 6.644249 1.99192690 11.296571
## 37 9.164653 4.30286074 14.026445
## 38 4.267665 -1.27353234 9.808863
## 39 6.761631 2.16712030 11.356142
## 40 7.347467 2.77148237 11.923451
## 41 8.349144 3.66196633 13.036322
## 42 6.576705 1.74690650 11.406503
## 43 6.804966 1.93875566 11.671176
## 44 7.100902 2.48370547 11.718098
## 45 7.602713 2.89938615 12.306040
## 46 4.947963 -0.03328176 9.929207
## 47 6.694208 1.99459141 11.393824
## 48 5.224958 0.22874650 10.221170
## 49 6.494056 1.84852804 11.139585
## 50 7.069247 2.18471219 11.953782
## 51 5.569420 -0.84543445 11.984275
## 52 4.768318 -0.34449649 9.881132
## 53 9.362059 4.44391306 14.280204
## 54 6.433093 1.77191508 11.094272
## 55 5.793980 0.44219862 11.145762
## 56 5.989367 0.36317627 11.615557
## 57 8.248342 3.40073909 13.095946
## 58 7.257039 2.42634499 12.087733
## 59 7.675477 3.02415511 12.326799
## 60 8.029888 2.73332290 13.326454
## 61 6.214982 1.55733669 10.872627
## 62 6.103263 1.42731992 10.779205
## 63 7.190814 2.55790710 11.823721
## 64 6.360102 1.68988227 11.030321
## 65 5.357807 -0.23373361 10.949348
## 66 4.208170 -1.04889215 9.465233
## 67 7.017108 2.18023182 11.853984
## 68 5.131619 -0.33602338 10.599262
## 69 5.940752 0.74460213 11.136902
## 70 7.687221 3.09517046 12.279271
## 71 5.718821 0.95817921 10.479464
## 72 9.034312 4.19940623 13.869218
## 73 8.185622 3.01138673 13.359858
## 74 9.605880 3.76919434 15.442565
## 75 4.738593 -0.37299420 9.850180
##
## $se.fit
## 1 2 3 4 5 6 7 8
## 2.630223 2.318992 2.669727 2.565679 2.364222 2.320385 2.455675 2.315000
## 9 10 11 12 13 14 15 16
## 2.374307 2.381836 2.324002 2.538842 2.411925 2.570137 2.793212 2.343508
## 17 18 19 20 21 22 23 24
## 2.312862 2.357503 2.450858 2.295721 2.375704 2.548869 2.342533 2.395563
## 25 26 27 28 29 30 31 32
## 2.598195 2.299933 2.549368 2.334381 2.452063 2.299025 2.298544 2.352350
## 33 34 35 36 37 38 39 40
## 2.613456 2.428086 2.342879 2.333227 2.438280 2.779015 2.304233 2.294942
## 41 42 43 44 45 46 47 48
## 2.350708 2.422234 2.440496 2.315610 2.358807 2.498188 2.356946 2.505694
## 49 50 51 52 53 54 55 56
## 2.329820 2.449686 3.217170 2.564172 2.466542 2.337668 2.684019 2.821640
## 57 58 59 60 61 62 63 64
## 2.431164 2.422684 2.332725 2.656327 2.335896 2.345073 2.323490 2.342203
## 65 66 67 68 69 70 71 72
## 2.804263 2.636515 2.425784 2.742125 2.605967 2.302999 2.387551 2.424796
## 73 74 75
## 2.594976 2.927207 2.563557
##
## $df
## [1] 71
##
## $residual.scale
## [1] 4.827578
Como podemos observar, los resultados son los mismos que en el inciso anterior ya que las proyecciones no cambian realmente. Esto se debe a que la información que se incluye por parte de la variable age ya está incorporada por la variable xper, por lo que age no podría inferir.