Punto 1 (30pts) Considerando la información suministrada en la base de datos de Vivienda. Realice lo siguiente:

Comenzamos leyendo la base de datos vivienda

## # A tibble: 8,322 × 12
##    Zona       piso  Estrato precio_millon Area_contruida parqueaderos Banos
##    <chr>      <chr>   <dbl>         <dbl>          <dbl> <chr>        <dbl>
##  1 Zona Sur   2           6           880          237   2                5
##  2 Zona Oeste 2           4          1200          800   3                6
##  3 Zona Sur   3           5           250           86   NA               2
##  4 Zona Sur   NA          6          1280          346   4                6
##  5 Zona Sur   2           6          1300          600   4                7
##  6 Zona Sur   3           6           513          160   2                4
##  7 Zona Sur   2           6           870          490   3                6
##  8 Zona Sur   5           5           310           82.5 1                2
##  9 Zona Sur   9           4           240           80   1                2
## 10 Zona Sur   6           6           690          150   2                5
## # … with 8,312 more rows, and 5 more variables: Habitaciones <dbl>, Tipo <chr>,
## #   Barrio <chr>, cordenada_longitud <dbl>, Cordenada_latitud <dbl>

(5pt) Filtre por barrio las variables precio por millón y área construída. El barrio de interés es ’el ingenio’

Se procede a realizar la seleccion de las variables:

## # A tibble: 202 × 3
##    precio_millon Area_contruida Barrio    
##            <dbl>          <dbl> <chr>     
##  1          1450           1200 el ingenio
##  2           290            100 el ingenio
##  3           360             99 el ingenio
##  4           550            197 el ingenio
##  5           410            136 el ingenio
##  6           390            198 el ingenio
##  7           300            147 el ingenio
##  8           250             97 el ingenio
##  9           370            104 el ingenio
## 10          1200            616 el ingenio
## # … with 192 more rows

(10pt) Repita en análisis de estas dos variables, usando RLS.

se procede a realizar analisis descriptivo:

Observando el histograma de la variable precio millon, se puede decir que la mayoria de las viviendas ubicadas en el barrio el ingenio tienen un costo entre 0 y 500 millones de pesos; sin embargo tambien se observa que algunas viviendas pueden llegar a costar hasta mas de 1000 millones de pesos en el mismo barrio.

Observando el histograma de la variable area construida, se dice que la mayoria de las viviendas ubicadas en el barrio el ingenio, presentan un area entre 0 y los 200 metros cuadrados; tambien se pueden observar algunas viviendas que presentan hasta el doble o el triple del area en el mismo barrio.

Observando el grafico de dispersión, se dice que la mayoría de los datos de las viviendas en el barrio el ingenio, se encuentran concentrados en un rango de 0 a 200 millones de pesos junto con un area de 0 a 400 metros cuadrados.

Ahora se procede a realizar un analisis por medio de RLS:

## 
## Call:
## lm(formula = precio_millon ~ Area_contruida, data = viviendanew)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -282.09  -61.15   -8.32   45.89  477.33 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    195.71496   11.63854   16.82   <2e-16 ***
## Area_contruida   1.21159    0.04632   26.16   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 105.2 on 200 degrees of freedom
## Multiple R-squared:  0.7738, Adjusted R-squared:  0.7727 
## F-statistic: 684.2 on 1 and 200 DF,  p-value: < 2.2e-16

Observando el modelo, se tiene que tanto el intercepto como la variable de area construida son significativos. Con respecto al coeficiente \(\beta_{1}\), se dice que por cada unidad de area construida va a aumentar en 1.21159 para la variable respuesta (precio millon). Ahora, observando el \(R^{2}\) del modelo indica que el modelo calculado explica el 77.38% de la variabilidad presente en la variable respuesta y el \(R^{2}\) ajustado fue de 77.27%.

Se procede a realizar la validación de los supuestos:

## Loading required package: zoo
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
## 
##  studentized Breusch-Pagan test
## 
## data:  modelo
## BP = 46.496, df = 1, p-value = 9.183e-12

Realizando el test de Breusch-pagan obtuvimos un valor p menor a 0.05 de tal manera que rechazamos la hipotesis nula y de esta manera no se cumple el supuesto de homocedasticidad.

###RESIDUOS
res.stud.base = studres(modelo)
modelo.fit.base= modelo$fitted.values
par(mfrow=c(1,2))
plot(modelo.fit.base,res.stud.base, ylab='residuos estudentizados',
     xlab='valores ajustados',main='(a)', col = "BLUE")
abline(h=0,lty=2)
lines(lowess(res.stud.base~modelo.fit.base), col = 2)

plot(modelo.fit.base,abs(res.stud.base), 
     ylab='valor absoluto de los residuos estudentizados',
     xlab='valores ajustados',main='(b)', col = "BLUE")
lines(lowess(abs(res.stud.base)~modelo.fit.base), col = 2)

Comprobando el supuesto de homocedasticidad y observando lo que nos arrojan las graficas podemos concluir que no se cumple este supuesto; ya que en los graficos a y b los datos no se destribuyen de manera equivalente, los datos tienden a estar en el lado izquierdo de las dos graficas.

car::qqPlot(modelo,xlab='cuantiles teóricos',ylab='residuos estudentizados', distribution = 'norm')

## [1]  75 125
ks.test(modelo$residuals,pnorm,mean(modelo$residuals),sd(modelo$residuals))
## Warning in ks.test.default(modelo$residuals, pnorm, mean(modelo$residuals), :
## ties should not be present for the Kolmogorov-Smirnov test
## 
##  Asymptotic one-sample Kolmogorov-Smirnov test
## 
## data:  modelo$residuals
## D = 0.10619, p-value = 0.02101
## alternative hypothesis: two-sided

Con la grafica de normalidad, se puede observar que algunos de los datos no siguen el patron lineal; por lo cual se puede decir que no se cumple el supuesto de normalidad. Para confirmar lo observado en la grafica, se procede a realizar el test de Kolmogorov-Smirnov, en el cual se obtiene un valor p de 0.02101, el cual es menor al nivel de significancia y por ende se concluye que el modelo no cumple el supuesto de normalidad.

como el modelo no cumple ninguno de los supuestos, se procede a aplicar una tranformación al modelo:

boxcox = MASS::boxcox(modelo,lambda=seq(-3,3,length.out = 1000),
                               ylab='log-verosimilitud')

boxcox$x[boxcox$y ==max(boxcox$y)] # valor que maximiza la log-verosimilitud
## [1] 0.2012012

Una vez obtenido el \(\lambda = 0.2012012\), se procede a realizar la tranformación del modelo

modelonew<-lm(I(precio_millon^0.2012012)~Area_contruida,data = viviendanew)
summary(modelonew)
## 
## Call:
## lm(formula = I(precio_millon^0.2012012) ~ Area_contruida, data = viviendanew)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.57302 -0.10417 -0.00491  0.09291  0.42248 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    3.032e+00  1.759e-02  172.31   <2e-16 ***
## Area_contruida 1.556e-03  7.002e-05   22.22   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1591 on 200 degrees of freedom
## Multiple R-squared:  0.7118, Adjusted R-squared:  0.7103 
## F-statistic: 493.9 on 1 and 200 DF,  p-value: < 2.2e-16

Cabe observar que con la transformación de box-cox realizada, el nuevo modelo sigue teniendo tanto el intercepto como la variable area construida como significativas; sin embargo el \(R^{2}\) ajustado disminuyo a un 71.03%.

Dicho lo anterior, se procede a realizar la validación de los supuestos para el nuevo modelo:

bptest(modelonew)
## 
##  studentized Breusch-Pagan test
## 
## data:  modelonew
## BP = 35.776, df = 1, p-value = 2.214e-09

Realizando el test de Breusch-pagan nuevamente, obtuvimos un valor p menor a 0.05 de tal manera que rechazamos la hipotesis nula y de esta manera no se cumple el supuesto de homocedasticidad.

###RESIDUOS
res.stud.base = studres(modelonew)
modelonew.fit.base= modelonew$fitted.values
par(mfrow=c(1,2))
plot(modelonew.fit.base,res.stud.base, ylab='residuos estudentizados',
     xlab='valores ajustados',main='(a)', col = "BLUE")
abline(h=0,lty=2)
lines(lowess(res.stud.base~modelonew.fit.base), col = 2)

plot(modelonew.fit.base,abs(res.stud.base), 
     ylab='valor absoluto de los residuos estudentizados',
     xlab='valores ajustados',main='(b)', col = "BLUE")
lines(lowess(abs(res.stud.base)~modelonew.fit.base), col = 2)  

Comprobando el supuesto de homocedasticidad por medio de las graficas, se observa nuevamente una tendencia para el costado izquierdo evitando ser equivalente tanto en la grafica a como en la b. Con estó y el test anterior se confirma el que el supuesto de homocedasticidad no se logra cumplir en el nuevo modelo.

car::qqPlot(modelonew,xlab='cuantiles teóricos',ylab='residuos estudentizados', distribution = 'norm')

## [1]  1 20
ks.test(modelonew$residuals,pnorm,mean(modelonew$residuals),sd(modelonew$residuals))
## Warning in ks.test.default(modelonew$residuals, pnorm,
## mean(modelonew$residuals), : ties should not be present for the Kolmogorov-
## Smirnov test
## 
##  Asymptotic one-sample Kolmogorov-Smirnov test
## 
## data:  modelonew$residuals
## D = 0.042428, p-value = 0.8603
## alternative hypothesis: two-sided

Observando la grafica de normalidad del nuevo modelo, se dice que los datos si presentan un patron lineal y por ende se concluye que se cumple el supuesto de normalidad. Por otro lado se confirma el supuesto al realizar el test de Kolmogorov-Smirnov, ya que se obtiene un valor p superior al nivel de significancia.

(10pt) Involucre la variable tipo de vivienda, y analice de nuevo usando RC.

Se procede a realizar nuevamente la selección de variables para asi agregar la variable tipo de vivienda y lograr realizar una regresión categórica:

viviendanew2 <- vivienda %>%
  filter(Barrio == "el ingenio")%>%
  select(precio_millon,Area_contruida,Tipo,Barrio)

modelo2<-viviendanew2
modelo2<- dummy_cols(modelo2,
                     select_columns = "Tipo")
modelo2
## # A tibble: 202 × 6
##    precio_millon Area_contruida Tipo        Barrio    Tipo_Apartamento Tipo_Casa
##            <dbl>          <dbl> <chr>       <chr>                <int>     <int>
##  1          1450           1200 Casa        el ingen…                0         1
##  2           290            100 Apartamento el ingen…                1         0
##  3           360             99 Apartamento el ingen…                1         0
##  4           550            197 Apartamento el ingen…                1         0
##  5           410            136 Apartamento el ingen…                1         0
##  6           390            198 Apartamento el ingen…                1         0
##  7           300            147 Apartamento el ingen…                1         0
##  8           250             97 Apartamento el ingen…                1         0
##  9           370            104 Apartamento el ingen…                1         0
## 10          1200            616 Casa        el ingen…                0         1
## # … with 192 more rows

Obtenida la nueva base de datos con la variable categorica tipo de vivienda y filtrada con el barrio “el ingenio”, Primero que todo se procede a realizar una grafica descriptiva de la nueva base:

Con este grafico de cajas se puede concluir que las casas situadas en el barrio el ingenio, son mas costosas con respecto a los apartamentos situados en el mismo barrio.

Ahora, se procede a realizar el nuevo modelo, observando el comportamiento cuando el tipo de vivienda es “casa”:

modelocasa<-lm(precio_millon~Area_contruida+Tipo_Casa, data = modelo2)
summary(modelocasa)
## 
## Call:
## lm(formula = precio_millon ~ Area_contruida + Tipo_Casa, data = modelo2)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -277.45  -48.20   -9.33   53.92  491.45 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    196.11964   10.91567  17.967  < 2e-16 ***
## Area_contruida   1.02197    0.05616  18.196  < 2e-16 ***
## Tipo_Casa       99.24986   18.63125   5.327 2.69e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 98.68 on 199 degrees of freedom
## Multiple R-squared:  0.802,  Adjusted R-squared:    0.8 
## F-statistic: 403.1 on 2 and 199 DF,  p-value: < 2.2e-16

Observando el modelo con el tipo de vivienda casa, se tiene que tanto el intercepto como la variable de area construida y el tipo casa son significativos para el modelo. Con respecto al coeficiente \(\beta_{1}\), se dice que por cada unidad de area construida va aumentar en 1.02197 para la variable respuesta (precio millon), y con respecto al coeficiente \(\beta_{2}\), se dice que si el tipo de vivienda es casa entonces tendra un incremento de 99.24986 en el precio de la vivienda. Ahora, observando el \(R^{2}\) del modelo, indica que el modelo calculado explica el 80.2% de la variabilidad presente en la variable respuesta, lo cual es algo muy bueno para el modelo y ademas, el \(R^{2}\) ajustado fue de 80%.

Se procede a realizar la validación de los supuestos:

## 
##  studentized Breusch-Pagan test
## 
## data:  modelocasa
## BP = 38.517, df = 2, p-value = 4.326e-09

Realizando el test de Breusch-pagan obtuvimos un valor p menor a 0.05 de tal manera que rechazamos la hipotesis nula y de esta manera no se cumple el supuesto de homocedasticidad.

###RESIDUOS
res.stud.base = studres(modelocasa)
modelocasa.fit.base= modelocasa$fitted.values
par(mfrow=c(1,2))
plot(modelocasa.fit.base,res.stud.base, ylab='residuos estudentizados',
     xlab='valores ajustados',main='(a)', col = "BLUE")
abline(h=0,lty=2)
lines(lowess(res.stud.base~modelocasa.fit.base), col = 2)

plot(modelocasa.fit.base,abs(res.stud.base), 
     ylab='valor absoluto de los residuos estudentizados',
     xlab='valores ajustados',main='(b)', col = "BLUE")
lines(lowess(abs(res.stud.base)~modelocasa.fit.base), col = 2)

Comprobando el supuesto de homocedasticidad y observando lo que nos arrojan las graficas podemos concluir que el modelo con tipo de vivienda casa en el barrio el ingenio no cumple este supuesto; ya que en los graficos a y b los datos no se distribuyen de manera equivalente, los datos tienden a estar en el lado izquierdo de las dos graficas.

car::qqPlot(modelocasa,xlab='cuantiles teóricos',ylab='residuos estudentizados', distribution = 'norm')

## [1]  75 125
ks.test(modelocasa$residuals,pnorm,mean(modelocasa$residuals),sd(modelocasa$residuals))
## Warning in ks.test.default(modelocasa$residuals, pnorm,
## mean(modelocasa$residuals), : ties should not be present for the Kolmogorov-
## Smirnov test
## 
##  Asymptotic one-sample Kolmogorov-Smirnov test
## 
## data:  modelocasa$residuals
## D = 0.085829, p-value = 0.102
## alternative hypothesis: two-sided

Con la grafica de normalidad, se puede observar que la mayoria de los datos siguen el patron lineal; por lo cual se puede decir que se cumple el supuesto de normalidad. Para confirmar lo observado en la grafica, se procede a realizar el test de Kolmogorov-Smirnov, en el cual se obtiene un valor p de 0.102, el cual es mayor al nivel de significancia y por ende se concluye que el modelo con el tipo de vivienda casa cumple el supuesto de normalidad.

Ahora se procede a realizar el modelo para el tipo de vivienda apartamento y asi realizar la comparación:

modeloapart<-lm(precio_millon~Area_contruida+Tipo_Apartamento, data = modelo2)
summary(modeloapart)
## 
## Call:
## lm(formula = precio_millon ~ Area_contruida + Tipo_Apartamento, 
##     data = modelo2)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -277.45  -48.20   -9.33   53.92  491.45 
## 
## Coefficients:
##                   Estimate Std. Error t value Pr(>|t|)    
## (Intercept)      295.36951   21.65885  13.637  < 2e-16 ***
## Area_contruida     1.02197    0.05616  18.196  < 2e-16 ***
## Tipo_Apartamento -99.24986   18.63125  -5.327 2.69e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 98.68 on 199 degrees of freedom
## Multiple R-squared:  0.802,  Adjusted R-squared:    0.8 
## F-statistic: 403.1 on 2 and 199 DF,  p-value: < 2.2e-16

Con respecto al modelo con el tipo de vivienda apartamento, se tiene que tanto el intercepto como la variable de area construida y el tipo apartamento tambien son significativos para el modelo. Con respecto al coeficiente \(\beta_{1}\), se dice que por cada unidad de area construida va aumentar en 1.02197 para la variable respuesta(precio millon) al igual que con el modelo anterior, y con respecto al coeficiente \(\beta_{2}\), se dice que si el tipo de vivienda es aparmento entonces tendra un decrecimiento de 99.24986 en el precio de la vivienda (es la diferencia si el tipo de vivienda es casa o apartamento). Ahora, observando el \(R^{2}\) del modelo, indica que el modelo calculado al igual que el anterior, explica el 80.2% de la variabilidad presente en la variable respuesta, lo cual es algo muy bueno para el modelo y su \(R^{2}\) ajustado fue de 80%.

Se procede a realizar la validación de los supuestos:

## 
##  studentized Breusch-Pagan test
## 
## data:  modeloapart
## BP = 38.517, df = 2, p-value = 4.326e-09

Realizando el test de Breusch-pagan obtuvimos un valor p menor a 0.05 de tal manera que rechazamos la hipotesis nula y de esta manera no se cumple el supuesto de homocedasticidad.

###RESIDUOS
res.stud.base = studres(modeloapart)
modeloapart.fit.base= modeloapart$fitted.values
par(mfrow=c(1,2))
plot(modeloapart.fit.base,res.stud.base, ylab='residuos estudentizados',
     xlab='valores ajustados',main='(a)', col = "BLUE")
abline(h=0,lty=2)
lines(lowess(res.stud.base~modeloapart.fit.base), col = 2)

plot(modeloapart.fit.base,abs(res.stud.base), 
     ylab='valor absoluto de los residuos estudentizados',
     xlab='valores ajustados',main='(b)', col = "BLUE")
lines(lowess(abs(res.stud.base)~modeloapart.fit.base), col = 2)

Comprobando el supuesto de homocedasticidad y observando lo que nos arrojan las graficas podemos concluir que el modelo con tipo de vivienda apartamento en el barrio el ingenio no cumple este supuesto; ya que en los graficos a y b los datos no se distribuyen de manera equivalente, los datos tienden a estar en el lado izquierdo de las dos graficas al igual que el modelo de casa.

car::qqPlot(modeloapart,xlab='cuantiles teóricos',ylab='residuos estudentizados', distribution = 'norm')

## [1]  75 125
ks.test(modeloapart$residuals,pnorm,mean(modeloapart$residuals),sd(modeloapart$residuals))
## Warning in ks.test.default(modeloapart$residuals, pnorm,
## mean(modeloapart$residuals), : ties should not be present for the Kolmogorov-
## Smirnov test
## 
##  Asymptotic one-sample Kolmogorov-Smirnov test
## 
## data:  modeloapart$residuals
## D = 0.085829, p-value = 0.102
## alternative hypothesis: two-sided

Con la grafica de normalidad, se puede observar que la mayoria de los datos siguen el patron lineal; por lo cual se puede decir que se cumple el supuesto de normalidad al igual que en modelo casa. Para confirmar lo observado en la grafica, se procede a realizar el test de Kolmogorov-Smirnov, en el cual se obtiene un valor p igual al modelo casa de 0.102, el cual es mayor al nivel de significancia y por ende se concluye que el modelo con el tipo de vivienda apartamento cumple el supuesto de normalidad.

Finalmente se concluye que los dos modelos, tanto el de casa como el de apartamento cumplen el supuesto de normalidad pero no cumplen el supuesto de homocedasticidad.

(5pt) Realice un ANOVA y compare con lo obtenido en el ítem anterior.

anova(modelocasa,modeloapart)
## Analysis of Variance Table
## 
## Model 1: precio_millon ~ Area_contruida + Tipo_Casa
## Model 2: precio_millon ~ Area_contruida + Tipo_Apartamento
##   Res.Df     RSS Df   Sum of Sq F Pr(>F)
## 1    199 1937949                        
## 2    199 1937949  0 -2.3283e-10

Como se puede observar se obtiene un valor p inferior al nivel de significancia, por lo cual se puede decir que si se encuentra dierencias entre el precio en el tipo de vivienda casa y apartamento.

Punto 2 (25pts) Usando la base de datos adjunta, datosME.txt. Donde se presenta la información de la masa y la edad de un grupo de personas, recopilada en un centro médico. Realice lo siguiente:

Comenzamos leyendo la base datosMe y vemos las primeras filas de la base:

##   Masa Edad
## 1  106   43
## 2  106   41
## 3   97   47
## 4  113   46
## 5   96   45
## 6  119   41

(10pt) Defina \(x_{i} = X_{i} − \bar{X}\) y genere un modelo de regresión cuadrático. Así mismo, adjunte una gráfica.

Primero se ajusta el modelo cuadratico con la variable sin centrar

mod=lm(Masa ~ Edad+I(Edad^2), data = base)
summary(mod)
## 
## Call:
## lm(formula = Masa ~ Edad + I(Edad^2), data = base)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -15.086  -6.154  -1.088   6.220  20.578 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 207.349608  29.225118   7.095 2.21e-09 ***
## Edad         -2.964323   1.003031  -2.955  0.00453 ** 
## I(Edad^2)     0.014840   0.008357   1.776  0.08109 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 8.026 on 57 degrees of freedom
## Multiple R-squared:  0.7632, Adjusted R-squared:  0.7549 
## F-statistic: 91.84 on 2 and 57 DF,  p-value: < 2.2e-16

Vemos en los resultados del modelo que el intercepto y la variable Edad tienen un aporte significativo pero el termino de la Edad al cuadrado, no. Además de que el \(R^{2}\) ajustado fue de 75.49%.

Se realiza la comprobacion de los supuestos:

## 
##  studentized Breusch-Pagan test
## 
## data:  mod
## BP = 7.379, df = 2, p-value = 0.02498

Con el test de Breusch-Pagan se obtiene un valor p menor a 0.05, por lo que hay problemas de heterocedasticidad.

###RESIDUOS
res.stud.base = studres(mod)
mod.fit.base= mod$fitted.values
par(mfrow=c(1,2))
plot(mod.fit.base,res.stud.base, ylab='residuos estudentizados',
     xlab='valores ajustados',main='(a)', col = "BLUE")
abline(h=0,lty=2)
lines(lowess(res.stud.base~mod.fit.base), col = 2)

car::qqPlot(mod,xlab='cuantiles teóricos',ylab='residuos estudentizados', distribution = 'norm',main='(b)')

## [1] 53 57
ks.test(mod$residuals,pnorm,mean(mod$residuals),sd(mod$residuals))
## Warning in ks.test.default(mod$residuals, pnorm, mean(mod$residuals),
## sd(mod$residuals)): ties should not be present for the Kolmogorov-Smirnov test
## 
##  Asymptotic one-sample Kolmogorov-Smirnov test
## 
## data:  mod$residuals
## D = 0.085408, p-value = 0.7739
## alternative hypothesis: two-sided

En la figura anterior, en la grafica (a) de los residuos vemos que estos se distribuyen alrededor de cero y parece que unos valores ajustados pequeños, tienen residuos un poco altos. En la grafica (b) vemos el qqplot donde los puntos estan cerca de la linea, por lo que puede haber normalidad en los residuos. Tambien se realiza el test de Kolmogorov-Smirnov, en el cual se obtiene un valor p de 0.7739 y se concluye que se cumple el supuesto de normalidad.

Como el modelo no cumple el supuesto de normalidad, se realiza una transformacion por el metodo Box Cox:

boxcox=MASS::boxcox(mod,lambda=seq(-3,3,length.out = 1000),ylab='log-verosimilitud')

boxcox$x[boxcox$y ==max(boxcox$y)] # valor que maximiza la log-verosimilitud
## [1] 1.750751

Entonces se ajusta el modelo con la transformacion:

modT=lm(I(Masa^1.750751) ~ Edad+I(Edad^2), data = base)
summary(modT)
## 
## Call:
## lm(formula = I(Masa^1.750751) ~ Edad + I(Edad^2), data = base)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -658.24 -304.50  -57.89  305.06  902.62 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 9562.8811  1393.0831   6.865 5.35e-09 ***
## Edad        -183.6967    47.8118  -3.842 0.000309 ***
## I(Edad^2)      1.0436     0.3983   2.620 0.011250 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 382.6 on 57 degrees of freedom
## Multiple R-squared:  0.7796, Adjusted R-squared:  0.7719 
## F-statistic: 100.8 on 2 and 57 DF,  p-value: < 2.2e-16

Vemos que el \(R^{2}\) ajustado fue de 77.19% y que todas las variables tienen un aporte significativo.

Se realiza la comprobacion de los supuestos:

## 
##  studentized Breusch-Pagan test
## 
## data:  modT
## BP = 2.425, df = 2, p-value = 0.2974

Con el test de Breusch-Pagan se obtiene un valor p mayor a 0.05, por lo que no hay problemas de heterocedasticidad.

###RESIDUOS
res.stud.base = studres(modT)
modT.fit.base= modT$fitted.values
par(mfrow=c(1,2))
plot(modT.fit.base,res.stud.base, ylab='residuos estudentizados',
     xlab='valores ajustados',main='(a)', col = "BLUE")
abline(h=0,lty=2)
lines(lowess(res.stud.base~modT.fit.base), col = 2)

car::qqPlot(modT,xlab='cuantiles teóricos',ylab='residuos estudentizados', distribution = 'norm',main='(b)')

## [1] 44 53
ks.test(modT$residuals,pnorm,mean(modT$residuals),sd(modT$residuals))
## Warning in ks.test.default(modT$residuals, pnorm, mean(modT$residuals), : ties
## should not be present for the Kolmogorov-Smirnov test
## 
##  Asymptotic one-sample Kolmogorov-Smirnov test
## 
## data:  modT$residuals
## D = 0.10061, p-value = 0.5781
## alternative hypothesis: two-sided

En la figura anterior, en la grafica (a) de los residuos vemos que estos se distribuyen alrededor de cero y no se nota un patron en los residuos. En la grafica (b) vemos el qqplot donde los puntos estan cerca de la linea, por lo que puede haber normalidad en los residuos. Tambien se realiza el test de Kolmogorov-Smirnov, en el cual se obtiene un valor p de 0.5781 y se concluye que se cumple el supuesto de normalidad.

vif(modT)
##      Edad I(Edad^2) 
##   128.256   128.256

Verificamos los VIF (Factor de Inflación de Varianza) y vemos que son muy altos, por lo que hay problemas de multicolinealidad.

Teniendo el modelo listo, ahora se centra la variable Edad restandole la media:

Edad.<-base$Edad - mean(base$Edad)

Se ajusta el modelo con la variable centrada

modC=lm(base$Masa ~ Edad.+I(Edad.^2))
summary(modC)
## 
## Call:
## lm(formula = base$Masa ~ Edad. + I(Edad.^2))
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -15.086  -6.154  -1.088   6.220  20.578 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 82.935749   1.543146  53.745   <2e-16 ***
## Edad.       -1.183958   0.088633 -13.358   <2e-16 ***
## I(Edad.^2)   0.014840   0.008357   1.776   0.0811 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 8.026 on 57 degrees of freedom
## Multiple R-squared:  0.7632, Adjusted R-squared:  0.7549 
## F-statistic: 91.84 on 2 and 57 DF,  p-value: < 2.2e-16

Vemos que se consigue el mismo \(R^{2}\) y \(R^{2}\) ajustado teniendo en cuenta el modelo sin variables centradas y sin transformacion. Al comprobar los supuestos se obtuvieron los mismos resultados y la misma transformacion, por lo que se ajusta el modelo con la variable centrada y la transformacion, y este modelo cumpliria los supuestos.

modCT=lm(I(base$Masa^1.750751) ~ Edad.+I(Edad.^2))
summary(modCT)
## 
## Call:
## lm(formula = I(base$Masa^1.750751) ~ Edad. + I(Edad.^2))
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -658.24 -304.50  -57.89  305.06  902.62 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 2299.1525    73.5576   31.26   <2e-16 ***
## Edad.        -58.4949     4.2249  -13.85   <2e-16 ***
## I(Edad.^2)     1.0436     0.3983    2.62   0.0112 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 382.6 on 57 degrees of freedom
## Multiple R-squared:  0.7796, Adjusted R-squared:  0.7719 
## F-statistic: 100.8 on 2 and 57 DF,  p-value: < 2.2e-16

Ahora se realiza una grafica de este modelo ajustado, donde vemos que la linea se ajusta a los puntos.

(5pt) Analice los estimadores de los parámetros \(β_{1}\) y \(β_{11}\).

Analizando los coeficientes del modelo, tenemos que todos fueron significativos, ademas si tambien vemos la grafica anterior, la edad tiene un efecto negativo sobre la masa, pues a mayor edad, menor masa tendra la persona.

(5pt) Ahora plantee la hipótesis de si el término cuadrático puede ser eliminado o no del modelo y justifique su respuesta.

Se define un modelo reducido, el cual no cuenta con el termino cuadratico

modCTR=lm(I(base$Masa^1.750751) ~ Edad.)

Para probar si el termino cuadratico puede ser eliminado se plantean las siguientes hipotesis en terminos del coeficiente asociado al termino cuadratico: \(H_{o}: \beta_{2} = 0 \quad vs \quad H_{a}: \beta_{2} \neq 0\)

Esto lo probamos usando la funcion anova, donde ingresamos un modelo sin el término cuadrático, seguido del modelo que contiene el término cuadrático

anova(modCTR,modCT)
## Analysis of Variance Table
## 
## Model 1: I(base$Masa^1.750751) ~ Edad.
## Model 2: I(base$Masa^1.750751) ~ Edad. + I(Edad.^2)
##   Res.Df     RSS Df Sum of Sq      F  Pr(>F)  
## 1     58 9346446                              
## 2     57 8341859  1   1004587 6.8644 0.01125 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Como el valor p da menor a 0.05, se cuenta con evidencia para rechazar la hipotesis nula, por lo que el coeficiente \(\beta_{2}\) es diferente de cero. Por lo tanto, el término cuadrático no debería ser eliminado de este modelo.

(5pt) ¿Justificó la transformación de la variable en el punto inicial (centrar la variable) para eliminar la multicolinealidad entre la variable y su forma cuadrática? Argumente usando las correlaciones.

Inicialmente vimos que los VIF del modelo sin la variable centrada eran mayores a 10 y esto denotaba problemas de multicolinealidad. Pero como se muestra a continuacion, al centrar la variable se consiguio eliminar el problema de multicolinealidad, pues los VIF son menores a 10.

vif(modCT)
##      Edad. I(Edad.^2) 
##   1.001473   1.001473

Tambien podemos ver las correlaciones, pues con la variable sin escalar se obtuvo una correlacion de 0.9960939.

cor(base$Edad,base$Edad^2)
## [1] 0.9960939

Pero cuando miramos la correlacion al centrar la variable, vemos que se reduce notoriamente, obteniendo -0.03835694.

cor(Edad.,Edad.^2)
## [1] -0.03835694

Finalmente, vimos que si se justifica el hecho de centrar la variable para eliminar la multicolinealidad entre dicha variable y su forma cuadrática.