#install.packages("devtools") # solo la primera vez
#install.packages("rlang")
#devtools::install_github("dgonxalex80/paqueteMOD", force =TRUE)
#install.packages("learnr")          # solo una vez
#learnr::run_tutorial("Tutorial101", "paqueteMOD")  # carga Tutorial101
library(paqueteMOD)
data("arboles")
head(arboles)
##   id  peso diametro altura
## 1  1 13.73      4.7    5.0
## 2  2 14.58      5.3    5.6
## 3  3 15.88      4.8    5.8
## 4  4  8.99      3.2    4.3
## 5  5  6.99      2.2    3.3
## 6  6 19.34      6.3    7.9

Problema 1 Para el grupo de datos contenidos en la base arboles (data(arboles)) determine que variable (altura o diametro) tiene más relación con la variable peso.

Hipótesis:A mayor altura del arbol mayor sera el peso del mismo Analisis univariado de los datos:

require(ggplot2)
## Loading required package: ggplot2
require(ggpubr)
## Loading required package: ggpubr
g1=ggplot(arboles,aes(x=altura))+geom_bar()+theme_bw()
ggarrange(g1, labels = c("A"),ncol = 1, nrow = 1)

Se puede visualizar que la altura de los arboles en el conjunto de datos analizados puede estar entre 4 y 8 metros de altura acorde al grafico.

Hipótesis:A mayor diametro del arbol mayor sera el peso del mismo Analisis univariado de los datos:

require(ggplot2)
require(ggpubr)
g1=ggplot(arboles,aes(x=diametro))+geom_bar()+theme_bw()
ggarrange(g1, labels = c("A"),ncol = 1, nrow = 1)

Acorde al grafico el dimaetro de los arboles analizados se encuentra entre 4 y 6 metros.

require(ggplot2)
require(ggpubr)
g1=ggplot(arboles,aes(x=peso))+geom_histogram(bins=30)+theme_bw()
g3=ggplot(arboles, aes(x=peso))+geom_boxplot(width=0.5)+theme_bw()
g2=ggplot(arboles,aes(x=diametro))+geom_histogram(bins=30)+theme_bw()
g4=ggplot(arboles, aes(x=diametro))+geom_boxplot(width=0.5)+theme_bw()
g5=ggplot(arboles,aes(x=altura))+geom_histogram(bins=30)+theme_bw()
g6=ggplot(arboles, aes(x=altura))+geom_boxplot(width=0.5)+theme_bw()
ggarrange(g1, g3, g2, g4, g5, g6, labels = c("A", "B", "C", "D", "E", "F"),ncol = 2, nrow = 3)

mean(arboles$peso) #promedio
## [1] 18.76611
median(arboles$peso) #mediana
## [1] 17.485
sd(arboles$peso) #desv
## [1] 8.157309
max(arboles$peso) #maximo ----
## [1] 47.87
min(arboles$peso) #minimo -----
## [1] 5.98

Para el caso de los graficos A y B se identifica que el peso de los árboles esta concentrado entre 5 a 25 toneladas. Mostrando asi que esta distribución tiene una asimetría positiva. A su vez el promedio del peso de los árboles es de 18.76 toneladas y la mediana que corresponde a 17.48 toneladas, por otra parte, la desviación estándar es de 8.157 toneladas. También podemos decir que el arbol mas pesado es de 47.87 toneladas, mientras que el menos pesado es de 5.98 toneladas.

mean(arboles$diametro) #promedio
## [1] 5.445556
median(arboles$diametro) #mediana
## [1] 5.4
sd(arboles$diametro) #desv
## [1] 1.451784
max(arboles$diametro) #maximo ----
## [1] 8.8
min(arboles$diametro) #minimo -----
## [1] 2.2

En cuanto a los diámetros, se puede observar en los gráficos C y D que esta variable tiene una distribución normal. los árboles tienen un diamtro entre 4.5 a 6.5 metros, con un promedio de diámetro de 5.44 metros y una desviación estándar de 1.451 metros. También podemos decir que el arbol mas ancho es de 8.8 metros, mientras que el menos ancho es de 2.2 metros.

mean(arboles$altura) #promedio
## [1] 6.634444
median(arboles$altura) #mediana
## [1] 6.45
sd(arboles$altura) #desv
## [1] 1.799386
max(arboles$altura) #maximo ----
## [1] 11.3
min(arboles$altura) #minimo -----
## [1] 3.3

Con respecto a la altura, se identifica en los grafico E y F que la altura de los arboles esta concentrada entre 5 y 7.5 metros, esta distribución también tiene una asimetría positiva. El promedio de altura de los arboles es de 6.63 metros y desviación estándar de 1.79 metros. También podemos decir que el arbol mas alto es de 11.3 metros, mientras que el mas bajo es de 3.3 metros.

Problema 2 Realice las pruebas de hipótesis con las que pueda verificar si las correlaciones entre las variables : peso, diametro y altura son diferentes de cero

cor(arboles$peso,arboles$diametro)
## [1] 0.908123
cor(arboles$peso,arboles$altura)
## [1] 0.8582009
cor.test(arboles$peso,arboles$diametro)
## 
##  Pearson's product-moment correlation
## 
## data:  arboles$peso and arboles$diametro
## t = 20.346, df = 88, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.8634081 0.9386817
## sample estimates:
##      cor 
## 0.908123
cor.test(arboles$peso, arboles$altura)
## 
##  Pearson's product-moment correlation
## 
## data:  arboles$peso and arboles$altura
## t = 15.684, df = 88, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.7918402 0.9045332
## sample estimates:
##       cor 
## 0.8582009

La correlaciòn entre el peso y el diamtro de un arbol es de 0.90 lo que quiere decir que los cambios en el peso estan directamente relacionados al dimatro que tenga el arbol. Podemos decir que hay una correlación positiva y fuerte. Y en el caso de la altura, también tiene una relación positiva y directa con el peso en un 0.85. Sin embargo podriamos decir que el diametro tiene mayor mayor relaciòn con el peso en los arboles que con la altura. Sin embargo, la altura también tiene una gran insidencia en la variación del peso ya que tanto para el diametro como para la altura la correlación de Pearson es cercana a 1.

require(ggplot2)
require(ggpubr)
g7=ggplot(arboles,aes(y=peso,x=diametro))+geom_point()+theme_bw()+geom_smooth()
g8=ggplot(arboles,aes(y=peso,x=altura))+geom_point()+theme_bw()+geom_smooth()
ggarrange(g7, g8, labels = c("A", "B"),ncol = 2, nrow = 1)
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

#learnr::run_tutorial("Tutorial102", "paqueteMOD")  # carga Tutorial102

Problema 3 La data biomasa contenida en paqueteMOD contiene información sobre la biomasa total y la altura de 90 árboles. Estime el modelo biomasa=𝛽0+𝛽1altura+𝜀 mediate el método de mínimos cuadrados ordinarios (MCO)

library(paqueteMOD)
data("biomasa")
modelo=lm(biomasa$bio_total ~ biomasa$altura)
summary(modelo)
## 
## Call:
## lm(formula = biomasa$bio_total ~ biomasa$altura)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -8.228 -1.969  0.572  2.377 15.106 
## 
## Coefficients:
##                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     -7.0456     1.7046  -4.133 8.14e-05 ***
## biomasa$altura   3.8906     0.2481  15.684  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4.211 on 88 degrees of freedom
## Multiple R-squared:  0.7365, Adjusted R-squared:  0.7335 
## F-statistic:   246 on 1 and 88 DF,  p-value: < 2.2e-16

Biomasa = -7.04 + 3.89*(Altura) Donde β0 toma el valor de -7.04 y β1 toma el valor de 3.89. Por lo que a partir de estos coeficientes podemos decir que cuando la altura es 0.00 la biomasa sera igual a β0 correspondiente a -7.04 y en el caso que la altura tome un valor diferente de 0.00, el β1 indica que por cada metro de altura el valor de la vivienda se incrementa 3.89.

El valor p indica si el coeficiente es significativo, para este caso presenta *** mostrando que la variable altura es significativa para el modelo y el r cuadrado muestra que el modelo logra explicar un 73% de la biomasa de los arboles. Sin embargo, no podriamos indicar que este sea el mejor modelo que explique la biomasa de los arboles. Por esto, se calcula un segundo modelo incorporando otra variable que puede afectar la biomasa de los arboles.

Problema 4 A partir de la data biomasa contenida en paqueteMOD contiene información sobre la biomasa total y la altura de 90 árboles. se puede estimar el modelo biomasa=𝛽0+𝛽1altura+𝜀 . A partir de los resultados obtenidos en la estimación valide los supuestos sobre los errores:

library(lmtest)
## Loading required package: zoo
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
require(zoo)
modelo=lm(biomasa$bio_total ~ biomasa$altura)
shapiro.test(modelo$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  modelo$residuals
## W = 0.95439, p-value = 0.003157
gqtest(modelo)
## 
##  Goldfeld-Quandt test
## 
## data:  modelo
## GQ = 1.1017, df1 = 43, df2 = 43, p-value = 0.3761
## alternative hypothesis: variance increases from segment 1 to 2
dwtest(modelo)
## 
##  Durbin-Watson test
## 
## data:  modelo
## DW = 0.90216, p-value = 4.234e-09
## alternative hypothesis: true autocorrelation is greater than 0
t.test(modelo$residuals)
## 
##  One Sample t-test
## 
## data:  modelo$residuals
## t = -2.9075e-16, df = 89, p-value = 1
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -0.8770047  0.8770047
## sample estimates:
##    mean of x 
## -1.28331e-16

Para el caso de la normalidad al realizar el test de Shapiro donde el p valor (0,0031) es menor al alfa de (0,05). Por lo cual, se rechaza la hipótesis nula (H0),y por lo tanto, la variable no presenta un comportamiento normal o paramétrico. Pero podemos realizar un ajuste al modelo para obtener un resultado mas consistente en relación a los residuos.Por lo tanto para ajustar el modelo se puede hacer una transformación utilizando un logaritmo.

para el Durbin-Watson test dado que este valor p es menor que 0.05, podemos rechazar la hipótesis nula y concluir que los residuos en este modelo de regresión están autocorrelacionados por lo cual podriamos analizar otras variables que afecten la biomasa de los arboles.

Problema 5

Para el modelo biomasa=𝛽0+𝛽1altura+𝜀 , obtenido a partir de la data biomasa contenida en paqueteMOD , verifique gráficamente el cumplimiento de los supuestos sobre los errores.

par(mfrow = c(2, 2))
plot(modelo)

- Media cero: En la primera gráfica (residuales vs ajustados) se observa que los errores no tienen un comportamiento aleatorio con media cero.

En base a la validación de supuestos descrita anteriormente, se decide hacer una transformación a las variables, con el fin de ajustar el modelo de regresión lineal y así tener un mejor rendimiento en la validación de los supuestos de los errores.

Problema 6

La data biomasa contenida en paqueteMOD contiene información sobre la biomasa total y la altura de 90 árboles. Verifique el valor de 𝜆 para la constante de la transformación Box-Cox para establecer la transformación apropiada de la variable dependiente en el modelo biomasa=𝛽0+𝛽1altura+𝜀

A partir del método de box-cox se calculo el λ optimo, para así identificar cual era la transformación que se debía realizar sobre la variable de respuesta Y. Se obtuvo un valor de λ de 0.0292 y en base a esto debemos transformar la variable respuesta con la función logaritmo.

library(MASS)
bc=boxcox(lm(biomasa$bio_total ~ biomasa$altura), lambda = -1:1)

(lambda <- bc$x[which.max(bc$y)])
## [1] 0.2929293

Modelo 2: se toman las dos variables tanto el diametro como la altura

modelo2=lm(biomasa$bio_total ~ diametro+altura, biomasa)
summary(modelo2)
## 
## Call:
## lm(formula = biomasa$bio_total ~ diametro + altura, data = biomasa)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -6.3083 -2.5121  0.1608  2.0088 11.7446 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  -9.1205     1.4305  -6.376 8.44e-09 ***
## diametro      4.7395     0.7128   6.649 2.49e-09 ***
## altura        0.3132     0.5751   0.544    0.587    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.449 on 87 degrees of freedom
## Multiple R-squared:  0.8253, Adjusted R-squared:  0.8213 
## F-statistic: 205.5 on 2 and 87 DF,  p-value: < 2.2e-16

Este modelo nos aroja el siguiente resultado

Biomasa = -9.12 + 4.73 * (diamtero) + 0.31 * (altura) Donde β0 toma el valor de -9.12, β1 toma el valor de 4.73 y β2 toma el valor de 0.31. Por lo que a partir de estos coeficientes podemos decir que cuando el diametro es 0.00 y la altura es 0.00 la biomasa sera igual a β0 correspondiente a -9.12 y en el caso que el diamtro y la altura tomen un valor diferente de 0.00, el β1 indica que por cada metro adicional en el diametro la biomasa se incrementa 4.73 toneladas, mientras que por cada metro adicional de altura la biomasa se incremenrta 0.31.

De acuerdo con el coeficiente de determinación ajustado, se puede concluir que el modelo explica en 82.27% la variabilidad del peso de los arboles. Este coeficiente indica que el modelo se ajusta bien a los datos, y que es un modelo confiable y tendría un buen desempeño para realizar predicciones.

El valor p indica si el coeficiente es significativo, para este caso presenta *** mostrando que la variable diametro es significativa para el modelo y el r cuadrado muestra que el modelo logra explicar un 82% de la biomasa de los arboles. Sin embargo, no podriamos indicar que este sea el mejor modelo que la explique ya que la variable altura no refleja significancia. Por lo cual ejecutamos un modelo en el cual se excluye la variable altura.

Modelo 3: Se excluye la variable altura que no se identifico como significativa para el modelo dos.

modelo3=lm(biomasa$bio_total ~ diametro, biomasa)
summary(modelo3)
## 
## Call:
## lm(formula = biomasa$bio_total ~ diametro, data = biomasa)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -6.3775 -2.6594  0.0237  1.8758 11.9876 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  -9.0203     1.4129  -6.384 7.86e-09 ***
## diametro      5.1026     0.2508  20.346  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.435 on 88 degrees of freedom
## Multiple R-squared:  0.8247, Adjusted R-squared:  0.8227 
## F-statistic:   414 on 1 and 88 DF,  p-value: < 2.2e-16

Este modelo nos aroja el siguiente resultado

Biomasa = -9.02 + 5.10 * (diamtero) Donde β0 toma el valor de -9.02, β1 toma el valor de 5.10. Por lo que a partir de estos coeficientes podemos decir que cuando el diametro es 0.00 la biomasa sera igual a β0 correspondiente a -9.02 y en el caso que el diametro tome un valor diferente de 0.00, el β1 indica que por cada metro adicional en el diametro la biomasa se incrementa 5.10 toneladas.

De acuerdo con el coeficiente de determinación ajustado, se puede concluir que el modelo explica en 82.27% la variabilidad del peso de los arboles. Este coeficiente indica que el modelo se ajusta bien a los datos, y que es un modelo confiable y tendría un buen desempeño para realizar predicciones.

El valor p indica si el coeficiente es significativo, para este caso presenta *** mostrando que la variable diametro es significativa para el modelo y el r cuadrado muestra que el modelo logra explicar un 82% la biomasa de los arboles.

par(mfrow = c(2, 2))
plot(modelo3)

require(zoo)
modelo3=lm(biomasa$bio_total ~ biomasa$altura)
shapiro.test(modelo3$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  modelo3$residuals
## W = 0.95439, p-value = 0.003157
gqtest(modelo3)
## 
##  Goldfeld-Quandt test
## 
## data:  modelo3
## GQ = 1.1017, df1 = 43, df2 = 43, p-value = 0.3761
## alternative hypothesis: variance increases from segment 1 to 2
dwtest(modelo3)
## 
##  Durbin-Watson test
## 
## data:  modelo3
## DW = 0.90216, p-value = 4.234e-09
## alternative hypothesis: true autocorrelation is greater than 0
t.test(modelo3$residuals)
## 
##  One Sample t-test
## 
## data:  modelo3$residuals
## t = -2.9075e-16, df = 89, p-value = 1
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -0.8770047  0.8770047
## sample estimates:
##    mean of x 
## -1.28331e-16

Validación de los supuestos de los errores

Debido a las falencias en la validación de supuestos se procede conla transformación a las variables para ajustar el modelo y así tener un mejor rendimiento en la validación de los supuestos de los errores.

Transformación Box-Cox Con el método de box-cox se calculo el λ optimo, para así identificar cual es la transformación que se debía realizar sobre la variable de respuesta Y (biomasa del arbol). Se obtuvo un valor de λ de 0.0292 y con esta información transfofrmamos la variable Y con la función logaritmo.

bc1=boxcox(lm(biomasa$bio_total ~ biomasa$diametro), lambda = -1:1)

(lambda <- bc$x[which.max(bc$y)])
## [1] 0.2929293

Transformación de la variable Y

Una vez realizada la transformación de la variable Y (logaritmo), se realiza nuevamente la prueba de shapiro-wilk sobre esta variable, y para este caso no se rechaza la hipótesis nula y se por lo tanto podemos concluir que los datos observados una vez han sido transformados siguen una distribución normal y para los gráficos de dispersión, se identifica una mejora en relación con la linea de las variables peso y diámetro.

biomasa2=biomasa
biomasa2$biomasa_transformado=log(arboles$peso)
shapiro.test(biomasa2$biomasa_transformado)
## 
##  Shapiro-Wilk normality test
## 
## data:  biomasa2$biomasa_transformado
## W = 0.99033, p-value = 0.7545

Para el caso de la normalidad se mantiene en mayor medida el mismo comportamiento centrado a la línea. Al realizar el test de Shapiro donde el p valor (0,75) es mayor a alfa (0,05) por lo cual no se rechaza la hipótesis nula (H0), y se indica que la variable presenta un comportamiento normal o paramétrico.

require(ggplot2)
require(ggpubr)
g9=ggplot(biomasa2,aes(y=biomasa_transformado,x=diametro))+geom_point()+theme_bw()+geom_smooth()
g10=ggplot(biomasa2,aes(y=biomasa_transformado,x=altura))+geom_point()+theme_bw()+geom_smooth()
ggarrange(g9, g10, labels = c("A", "B"),ncol = 2, nrow = 1)
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Modelo de regresión lineal a partir de la variable Y transformada

modelo_transformado=lm(biomasa2$biomasa_transformado ~ biomasa2$diametro+ biomasa2$altura)
summary(modelo_transformado)
## 
## Call:
## lm(formula = biomasa2$biomasa_transformado ~ biomasa2$diametro + 
##     biomasa2$altura)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -0.2722 -0.1168 -0.0100  0.1018  0.3198 
## 
## Coefficients:
##                   Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        1.34378    0.05917  22.710   <2e-16 ***
## biomasa2$diametro  0.33539    0.02949  11.374   <2e-16 ***
## biomasa2$altura   -0.04934    0.02379  -2.074   0.0411 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1427 on 87 degrees of freedom
## Multiple R-squared:  0.8919, Adjusted R-squared:  0.8894 
## F-statistic: 358.8 on 2 and 87 DF,  p-value: < 2.2e-16

Ya que las variables diámetro y altura son significativas accorde a los P-valores, donde se concluye que no se rechazan las hipótesis nulas con un 95% de confianza, el peso de los árboles depende de las variables: Diametro y Altura.

Se tiene que 𝛽0= 1.34. Este resultado indica que: Con un diámetro y altura de 0, el peso esperado de los árboles seria de 3.83 toneladas.

Se tiene que 𝛽1= 0.33. Este resultado indica que: Por cada metro adicional en el diámetro del árbol se espera un aumento promedio de 1.39 toneladas en el peso del árbol.

Se tiene que 𝛽2= -0.04. Este resultado indica que: Por cada metro adicional en la altura del árbol se espera un aumento promedio de 0.95 toneladas en el peso del árbol.

El coeficiente de determinación ajustado,indica que el modelo explica en 88.94% la variabilidad del peso de los arboles y por lo tanto el modelo se ajusta bien a los datos, y que es confiable por lo que tendría un buen desempeño para realizar predicciones.

par(mfrow = c(2, 2))
plot(modelo_transformado)

library(lmtest)
require(zoo)
shapiro.test(modelo_transformado$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  modelo_transformado$residuals
## W = 0.97487, p-value = 0.07832
gqtest(modelo_transformado)
## 
##  Goldfeld-Quandt test
## 
## data:  modelo_transformado
## GQ = 1.0399, df1 = 42, df2 = 42, p-value = 0.4499
## alternative hypothesis: variance increases from segment 1 to 2
dwtest(modelo_transformado)
## 
##  Durbin-Watson test
## 
## data:  modelo_transformado
## DW = 0.83921, p-value = 2.663e-10
## alternative hypothesis: true autocorrelation is greater than 0
t.test(modelo_transformado$residuals)
## 
##  One Sample t-test
## 
## data:  modelo_transformado$residuals
## t = 2.1938e-15, df = 89, p-value = 1
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -0.02954031  0.02954031
## sample estimates:
##    mean of x 
## 3.261461e-17

Validación de supuestos de los errores