Informe de aplicación de regresión lineal Multiple

Punto 1

Analisis exploratorio de datos caso BD OLX venta de autos mazda2

En la tabla se observa que el promedio de precios de los autos mazda2 en la pagina olx Colombia es de $ 52.000.000 con una desviación estandar de $ 22.700.000, los caules tienen un promedio de modelo de 2020 y un kilometraje de 72.600 Km con una variación muy alta de 48.9000 KM, se observa en la figura 1 que la mayoria de los autos se agrupan entre precios de 22 a 95 Millones, en cuanto al modelo la mayor concentracion estan entre 2008 y 2022 siento este el modelo más predominante, en la figura 3 se observa que la mayor cantidad de autos son 0 Km y los demás no superan los 280.000 km, en cuanto a la relación del precio y el tipo de transmisión figura 7 los autos con transmisión automatica tienden ha ser más costosos que los mecanicos. Los municipios donde hay más ofertas de autos son Cali, medellin y manizales.

Finalmente la Figura 5 nos muestra que existe una relación lineal exponencial entre precio y modelo entre más reciente es el modelo tiene myor precio y en la figura 6 entre precio y kilometraje tambien existe una relación exponnecial inversa a medida que el kilometraje aunmenta disminuye el precio de venta.

# Cargar los datos

bdolxcar =read_excel("C:/Mesa/Informe regresión múltiple/olxcar.xlsx", 
    col_types = c("text", "numeric", "numeric", 
        "numeric", "text", "text", "text", 
        "text"))

bdolxcar=bdolxcar[,2:8]

bdolxcar$tranmision=factor(bdolxcar$tranmision,levels = c("Automática","Automática Secuencial","Mecánica"))
head(bdolxcar)
## # A tibble: 6 × 7
##      precio modelo kilometraje tranmision            ciudad      Munic…¹ Depar…²
##       <dbl>  <dbl>       <dbl> <fct>                 <chr>       <chr>   <chr>  
## 1 136900000   2010       88000 Automática            Astorga     Medell… Antioq…
## 2  38000000   2013       66000 Automática            San Carlos  Tunjue… Bogotá 
## 3  64500000   2018       74000 Automática Secuencial La Base     Cali    Valle …
## 4  88650000   2022           0 Automática            Ciudad Cam… Cali    Valle …
## 5  88650000   2022           0 Automática            Ciudad Cam… Cali    Valle …
## 6  88650000   2022           0 Automática            Ciudad Cam… Cali    Valle …
## # … with abbreviated variable names ¹​Municipio, ²​Departamento
## Análisis univariado y bivariado
attach(bdolxcar)
summary(bdolxcar)
##      precio              modelo      kilometraje    
##  Min.   :  6800000   Min.   :1995   Min.   :     0  
##  1st Qu.: 35500000   1st Qu.:2012   1st Qu.: 34650  
##  Median : 42300000   Median :2015   Median : 79000  
##  Mean   : 52012161   Mean   :2015   Mean   : 72564  
##  3rd Qu.: 66000000   3rd Qu.:2019   3rd Qu.:103000  
##  Max.   :169800000   Max.   :2022   Max.   :280000  
##                  tranmision     ciudad           Municipio        
##  Automática           :130   Length:273         Length:273        
##  Automática Secuencial: 20   Class :character   Class :character  
##  Mecánica             :121   Mode  :character   Mode  :character  
##  NA's                 :  2                                        
##                                                                   
##                                                                   
##  Departamento      
##  Length:273        
##  Class :character  
##  Mode  :character  
##                    
##                    
## 
table1(~precio+modelo+kilometraje,data=bdolxcar)
Overall
(N=273)
precio
Mean (SD) 52000000 (22700000)
Median [Min, Max] 42300000 [6800000, 170000000]
modelo
Mean (SD) 2020 (4.55)
Median [Min, Max] 2020 [2000, 2020]
kilometraje
Mean (SD) 72600 (48900)
Median [Min, Max] 79000 [0, 280000]
g1=ggplot(bdolxcar,aes(x=precio))+geom_histogram()+theme_bw()
g2=ggplot(bdolxcar,aes(x=modelo))+geom_histogram()+theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5))
g3=ggplot(bdolxcar,aes(x=kilometraje))+geom_histogram()+theme_bw()
g4 =ggplot(bdolxcar, aes(Municipio))+geom_bar()+ggtitle("Figura 4. Hist de municipios")+theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5))

fig1 =subplot(ggplotly(g1), ggplotly(g2), ggplotly(g3), titleX = TRUE,titleY = TRUE) %>% layout(title = "  Fig 1. Hist de precios | Fig 2. Hist de modelo | Fig 3. Hist de Kilometraje")
fig1
ggplotly(g4)
g5=ggplot(data = bdolxcar,mapping = aes(x=modelo,y=precio))+geom_point()+geom_smooth()+ggtitle("Figura 5. Relación precio VS modelo")+theme_bw()

g6=ggplot(data = bdolxcar,mapping = aes(x= kilometraje,y=precio))+geom_point()+geom_smooth()+ggtitle("Figura 6. Relación precio VS kilometraje")+theme_bw()

g5

g6

g7=ggplot(bdolxcar,aes(x=tranmision,y=precio,fill=tranmision))+geom_boxplot()+ggtitle("Figura 7. Relación precio VS transmision")+theme_bw()
ggplotly(g7)

Modelo de regresión multiple

Se realiza el modelo con todas las variables y se evidencia que los p-value del modelo, transmisión y kilometraje son significativos siendo este ultimo el menos de los tres, de igual forma el municipio arrojo que no aporta mucho al modelo.

El \(R^2\) ajustado es medio alto (0,61) y significa que es capaz de explicar el modelo en un 61%

##bdolxcar2=subset(bdolxcar, tranmision !='Automática Secuencial')
mod_autos=lm(precio~modelo+kilometraje+tranmision+Municipio,data=bdolxcar)
summary(mod_autos)
## 
## Call:
## lm(formula = precio ~ modelo + kilometraje + tranmision + Municipio, 
##     data = bdolxcar)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -67474245  -4122186    -63399   2330645  97884343 
## 
## Coefficients:
##                                   Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                     -6.101e+09  7.248e+08  -8.417 5.13e-15 ***
## modelo                           3.057e+06  3.586e+05   8.523 2.59e-15 ***
## kilometraje                     -6.426e+01  3.185e+01  -2.018 0.044862 *  
## tranmisionAutomática Secuencial  2.863e+06  3.678e+06   0.778 0.437291    
## tranmisionMecánica              -8.021e+06  2.203e+06  -3.641 0.000339 ***
## MunicipioArmenia                 3.638e+06  6.884e+06   0.529 0.597654    
## MunicipioBarranquilla           -1.245e+06  6.990e+06  -0.178 0.858857    
## MunicipioBarrios Unidos         -8.863e+06  9.639e+06  -0.920 0.358847    
## MunicipioBello                   8.749e+05  7.322e+06   0.119 0.905000    
## MunicipioBosa                   -8.111e+05  1.511e+07  -0.054 0.957241    
## MunicipioBucaramanga            -8.688e+05  7.348e+06  -0.118 0.905991    
## MunicipioCaldas                 -6.560e+06  1.515e+07  -0.433 0.665409    
## MunicipioCali                    3.844e+06  5.377e+06   0.715 0.475381    
## MunicipioCesar                   3.601e+06  1.513e+07   0.238 0.812141    
## MunicipioChapinero              -5.473e+06  9.733e+06  -0.562 0.574495    
## MunicipioChía                    1.955e+07  9.754e+06   2.004 0.046265 *  
## MunicipioCórdoba                -9.539e+06  1.513e+07  -0.630 0.529134    
## MunicipioCúcuta                 -2.260e+06  9.832e+06  -0.230 0.818425    
## MunicipioCundinamarca            7.203e+06  9.648e+06   0.747 0.456099    
## MunicipioEngativá                1.664e+05  7.012e+06   0.024 0.981090    
## MunicipioEnvigado                4.928e+06  8.218e+06   0.600 0.549363    
## MunicipioFloridablanca           1.600e+06  7.684e+06   0.208 0.835280    
## MunicipioFontibón                6.582e+06  1.506e+07   0.437 0.662466    
## MunicipioFunza                   2.750e+06  1.514e+07   0.182 0.856013    
## MunicipioGirón                  -7.651e+05  1.508e+07  -0.051 0.959582    
## MunicipioHuila                  -4.964e+06  1.536e+07  -0.323 0.746797    
## MunicipioIbagué                  2.213e+06  9.693e+06   0.228 0.819651    
## MunicipioItagüí                  7.658e+05  8.849e+06   0.087 0.931112    
## MunicipioKennedy                -5.747e+06  8.430e+06  -0.682 0.496140    
## MunicipioLos Mártires           -8.764e+06  1.133e+07  -0.773 0.440118    
## MunicipioManizales               7.828e+05  6.212e+06   0.126 0.899836    
## MunicipioMedellín                3.949e+06  5.535e+06   0.713 0.476333    
## MunicipioMeta                   -1.026e+07  1.522e+07  -0.675 0.500681    
## MunicipioNeiva                   1.306e+06  7.679e+06   0.170 0.865136    
## MunicipioPereira                 3.635e+05  7.104e+06   0.051 0.959238    
## MunicipioPiedecuesta             1.938e+05  1.506e+07   0.013 0.989742    
## MunicipioPopayán                 2.683e+06  1.512e+07   0.177 0.859311    
## MunicipioPuente Aranda           1.706e+05  1.508e+07   0.011 0.990984    
## MunicipioRafael Uribe Uribe     -1.148e+05  1.507e+07  -0.008 0.993928    
## MunicipioSabaneta               -3.675e+06  1.517e+07  -0.242 0.808791    
## MunicipioSan Cristóbal          -5.835e+06  1.124e+07  -0.519 0.604204    
## MunicipioSan Juan de Pasto       1.508e+06  1.129e+07   0.134 0.893838    
## MunicipioSanta Fé                1.338e+06  1.516e+07   0.088 0.929783    
## MunicipioSanta Marta             3.566e+06  1.530e+07   0.233 0.815929    
## MunicipioSantander               3.151e+06  1.124e+07   0.280 0.779521    
## MunicipioSuba                    3.823e+06  7.677e+06   0.498 0.619048    
## MunicipioSucre                  -1.193e+07  1.555e+07  -0.767 0.443873    
## MunicipioTeusaquillo             1.369e+06  1.507e+07   0.091 0.927686    
## MunicipioTunjuelito             -9.806e+06  1.516e+07  -0.647 0.518495    
## MunicipioUsaquén                 1.964e+04  1.130e+07   0.002 0.998616    
## MunicipioValle del Cauca        -1.896e+06  8.843e+06  -0.214 0.830473    
## MunicipioValledupar              3.396e+05  1.526e+07   0.022 0.982259    
## MunicipioVillavicencio          -1.317e+06  8.182e+06  -0.161 0.872241    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 14170000 on 218 degrees of freedom
##   (2 observations deleted due to missingness)
## Multiple R-squared:  0.6862, Adjusted R-squared:  0.6114 
## F-statistic:  9.17 on 52 and 218 DF,  p-value: < 2.2e-16
par(mfrow=c(2,2))
plot(mod_autos)

Selección de los mejores predictores

En este caso se van a emplear la estrategia de stepwise mixto y el valor matemático empleado para determinar la calidad del modelo va a ser Akaike(AIC), con esto se reafirma la no significancia del municipio y finalmente el modelo queda \(precio_auto\)= \(-6.607e+09\)+ \(3.308e+06*modelo-4.856e+01*kilometraje+2.279e+06*tranmisionAutomática Secuencial-7.448e+06*tranmisionMecánica\), donde los valores \(\beta\) reafirman la evidencia exploratoria entre mnayor sea el modelo y el kilometraje menor el precio, los autos de transmision automatica tienen mayor precio que los de trasmision mecanica.

El modelo con las dos variables predictorias tiene un R2 alto (0,8253), lo que significa que es capaz de explicarlo en un 82,53% de la variablidad observada del peso, el p-value de las variables diametro y altura son significativos.

Validación del modelo

Media cero: Se cumple por defecto.

Varianza Constante: Se observa en la grafica 1 de residuales vs ajustados que el comportamiento es aleatorio aunque los valores se agrupan en la recta

Normalidad: Se observa en la grafica 2 que los datos se ajustan bien a la linea de normalidad en el qqplot

Independencia: Dado que estos registros no corresponden a datos en el tiempo no se tiene un orden temporal para realizar la validación de este supuesto.

step(object = mod_autos, direction = "both", trace = 1)
## Start:  AIC=8971.82
## precio ~ modelo + kilometraje + tranmision + Municipio
## 
##               Df  Sum of Sq        RSS    AIC
## - Municipio   48 3.7719e+15 4.7524e+16 8898.2
## <none>                      4.3752e+16 8971.8
## - kilometraje  1 8.1696e+14 4.4569e+16 8974.8
## - tranmision   2 3.0537e+15 4.6806e+16 8986.1
## - modelo       1 1.4579e+16 5.8332e+16 9047.8
## 
## Step:  AIC=8898.23
## precio ~ modelo + kilometraje + tranmision
## 
##               Df  Sum of Sq        RSS    AIC
## <none>                      4.7524e+16 8898.2
## - kilometraje  1 6.2559e+14 4.8150e+16 8899.8
## - tranmision   2 3.8093e+15 5.1334e+16 8915.1
## + Municipio   48 3.7719e+15 4.3752e+16 8971.8
## - modelo       1 2.3988e+16 7.1512e+16 9007.0
## 
## Call:
## lm(formula = precio ~ modelo + kilometraje + tranmision, data = bdolxcar)
## 
## Coefficients:
##                     (Intercept)                           modelo  
##                      -6.607e+09                        3.308e+06  
##                     kilometraje  tranmisionAutomática Secuencial  
##                      -4.856e+01                        2.279e+06  
##              tranmisionMecánica  
##                      -7.448e+06
mod_autos_final=lm(formula = precio ~ modelo + kilometraje + tranmision, data = bdolxcar)
summary(mod_autos_final)
## 
## Call:
## lm(formula = precio ~ modelo + kilometraje + tranmision, data = bdolxcar)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -66583782  -4676144  -1293276   2609282  99806216 
## 
## Coefficients:
##                                   Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                     -6.607e+09  5.768e+08 -11.455  < 2e-16 ***
## modelo                           3.308e+06  2.855e+05  11.587  < 2e-16 ***
## kilometraje                     -4.856e+01  2.595e+01  -1.871   0.0624 .  
## tranmisionAutomática Secuencial  2.279e+06  3.263e+06   0.698   0.4856    
## tranmisionMecánica              -7.448e+06  1.725e+06  -4.319 2.22e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 13370000 on 266 degrees of freedom
##   (2 observations deleted due to missingness)
## Multiple R-squared:  0.6592, Adjusted R-squared:  0.6541 
## F-statistic: 128.6 on 4 and 266 DF,  p-value: < 2.2e-16
par(mfrow=c(2,2))
plot(mod_autos_final)

validar_modelo=function(datos,n){
  var_predict=array(NA,n)
  for (i in 1:n) {
    datos_training=datos[-i,]
    datos_test=datos[i,]
    mod=lm(formula = precio ~ modelo + kilometraje + tranmision, data = datos_training)
    var_predict[i]=predict(mod,list((modelo=datos_test$modelo),(kilometraje=datos_test$kilometraje),(tranmision=datos_test$tranmision)))
  
  }
  return(var_predict)
}

Validación cruzada del modelo final de autos

El modelo ajustado haciendo la validación cruzada con la técnica MAE se obtiene un 14,74% de probabilidad de error de predicción, lo que constituye un regular modelo predictivo y el RMSE 29,317% la proporción del error versus el promedio del precio observado

precio_predict=validar_modelo(bdolxcar,dim(bdolxcar)[1])
res=data.frame(bdolxcar,precio_predict)
##head(res)
MAE=mean(abs(res$precio-na.omit(res$precio_predict)))
print(paste("MAE= ",MAE/mean(res$precio)*100))
## [1] "MAE=  14.7418563122248"
RMSE=sqrt(mean((res$precio-na.omit(res$precio_predict))^2))
print(paste("RMSE",RMSE/mean(res$precio)*100))
## [1] "RMSE 29.3174676370685"

Modelo como herramienta practica

Discutir potenciales usos de modelo como herramienta practica (como monetizar los resultados de este modelo)

  • Este tipo de desarrollos permiten agregar valor a las organizaciones hacia la toma de decisiones basada en datos, para el caso puntual inicialmente permite describir que pasa en el sector de ventas de vehículos de la referencia Mazda2 y poder determinar estadísticas de promedio nacional, por ciudad. Kilometraje y tipo de vehículo. Sumado a esto el modelo proporciona una referencia predictiva importante que le permite a la empresa comercializadora de compra de vehículos mejorar su poder de negociación frente a la compra y venta de autos y la oportunidad de estimar el precio partir de variables predictoras.
  • El modelo acompañado con desarrollo de una APP para analítica de dato sistematizando las capas de extracción, transformación y visualización, representa una herramienta importante que les permita procesar grandes volúmenes de datos, trabajarlos en forma dinámica y rápida, y transformarlos en acciones concretas que generen valor en la compra y venta de autos.

Punto 2

Dos tipos de moluscos A y B fueron sometidos a tres concentraciones distintas de agua de mar (100%, 75%, 50%) y se observo el consumo de oxigeno midiendo la proporcion de O2 por unidad peso seco del molusco.

Realice un analisis exploratorio que permita conocer como es el consumo de oxigeno en las distintas concentraciones de agua de mar. y si estas conclusiones son las mismas para cada tipo de molusco

require(ggplot2)
require(dplyr)


load("C:/Mesa/Informe regresión múltiple/moluscos.RData")
BaseMoluscos = data.frame(BD_moluscos)
BaseMoluscos$c_agua <- as.factor(BaseMoluscos$c_agua)
head(BaseMoluscos)
##   c_agua molusco cons_o
## 1    100       A   7.16
## 2    100       A   8.26
## 3    100       A   6.78
## 4    100       A  14.00
## 5    100       A  13.60
## 6    100       A  11.10
moluscoA = filter(BaseMoluscos, BaseMoluscos$molusco == "A")

moluscoB = filter(BaseMoluscos, BaseMoluscos$molusco == "B")

Análisis Exploratorio

summary(moluscoA$cons_o)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    5.20    7.18    9.74   10.00   11.28   18.80
summary(moluscoB$cons_o)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   1.800   5.685   8.060   8.609  10.700  17.700

Se aplica summary a la variable consumo de oxígeno (cons_o), para el caso de los moluscos tipo A tenemos una media de 9.74 mientras que su mínimo es de 5.20 y su máximo de 18.80. Para el caso de los moluscos tipo B tenemos una media de 8.060, un mínimo de 1.8000 y un máximo de 17.700

attach(BaseMoluscos)
summary(BaseMoluscos)
##  c_agua     molusco              cons_o      
##  50 :16   Length:48          Min.   : 1.800  
##  75 :16   Class :character   1st Qu.: 6.312  
##  100:16   Mode  :character   Median : 9.700  
##                              Mean   : 9.305  
##                              3rd Qu.:11.232  
##                              Max.   :18.800
ggplot(BaseMoluscos,aes(x=cons_o))+geom_histogram()+labs(x="Consumo de Oxigeno",y="Count")+ggtitle("Consumo de Oxigeno por Tipo de Molusco")+facet_grid(~molusco)

En cuanto al consumo de oxígeno independientemente de la cantidad de concentración de agua, podemos observar que es muy similar para ambos tipos de moluscos, sin embargo, para los moluscos de Tipo A podemos ver que su consumo mínimo parte desde 5 mientras que para los moluscos de tipo B aproximadamente de 2

g3=ggplot(data = BaseMoluscos,mapping = aes(x=c_agua,y=cons_o))+geom_point()+geom_smooth()+ggtitle("Relacion Concretacion de agua VS Consumo de oxigeno")+theme_bw()

g3

g4 <- ggplot(BaseMoluscos, aes(x = c_agua, y = cons_o)) +
  geom_point() + facet_grid(~molusco) +theme_bw() +ggtitle("Relacion Tipo de molusco VS Consumo de Oxgeno")+ 
  labs(x = "Concentracion agua", y = "Consumo de Oxigeno")
g4

g1=ggplot(BaseMoluscos,aes(x= c_agua, y=cons_o, fill=molusco))+
  geom_boxplot()+
  theme_bw()+
  geom_point(position=position_jitterdodge(),alpha=0.4)+
  ggtitle("Consumo de Oxigeno por tipo de Molusco en diferentes consetracciones de Agua")+
labs(x="Consentracion de Agua",y="Consumo de Oxigeno")
g1

Podemos observar que el consumo de oxígeno es más alto cuando los moluscos se encuentran sometidos a una concentración de 50% de agua de mar y el consumo disminuye cuando la concentración de agua es del 75%.

interaction.plot(
  x.factor = BaseMoluscos$c_agua,
  trace.factor = BaseMoluscos$molusco,
  response = BaseMoluscos$cons_o,
  fun = median,
  ylab = "Consumo de Oxigeno",
  xlab = "Concetracion de Agua",
  trace.label = "Tipo de Molusco",
  lwd = 1
)

Con base a la gráfica podemos identificar el consumo de oxígeno por cada tipo de molusco dependiendo de la concentración de agua de mar a la que esten expuestos.

BaseMoluscos=na.omit(BaseMoluscos)
BaseMoluscos$molusco=as.character(BaseMoluscos$molusco)
g2=ggplot(BaseMoluscos,aes(x=c_agua,y=cons_o,color=molusco))+geom_point()+theme_bw()+geom_smooth()+ ggtitle("Consumo de Oxigeno por tipo de Molusco en diferentes consetracciones de Agua")+
labs(x="Consentracion de Agua",y="Consumo de Oxigeno")
g2

Con base a esta última gráfica podemos observar que independientemente de la concentración de agua a la que esten expuestos los tipos de moluscos el molusco de tipo A fue quien obtuvo mayor consumo de oxígeno en los 3 escenarios a los que fueron sometidos.

  1. Estime el modelo de diseño de experimentos el cual permita evaluar el efecto de la concentración de agua de mar y los tipos de molusco sobre el consumo de oxigeno. interprete los coeficientes del modelo, el valor p y realice un post anova de considerarlo necesario para los factores
ModMoluscos = lm(formula = cons_o ~ c_agua + molusco, data=BaseMoluscos)
summary(ModMoluscos)
## 
## Call:
## lm(formula = cons_o ~ c_agua + molusco, data = BaseMoluscos)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.1750 -1.9877 -0.7019  2.1244  6.1450 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  12.9463     0.8521  15.193  < 2e-16 ***
## c_agua75     -5.2581     1.0436  -5.038 8.49e-06 ***
## c_agua100    -3.5794     1.0436  -3.430  0.00132 ** 
## moluscoB     -1.3913     0.8521  -1.633  0.10966    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.952 on 44 degrees of freedom
## Multiple R-squared:  0.3986, Adjusted R-squared:  0.3575 
## F-statistic: 9.719 on 3 and 44 DF,  p-value: 4.866e-05

Concentracion de agua de 75% (c_agua75) = 8.49e-06 Concentracion de agua de 100% (c_agua100) = 0.00132 moluscoB = 0.10966

Análisis: Como evidenciamos en los resultados del summary el valor que se obtuvo para p-value es 0.05, lo que sugeriría que existe una relación lineal entre las variables, pero solo para los casos donde las concentraciones de agua son 75% y 100%, es decir que el modelo se ve afectado solo cuando hay un cambio en las concentraciones de agua, mientras que si se cambia el tipo de molusco este no tendría un cambio

anova(ModMoluscos)
## Analysis of Variance Table
## 
## Response: cons_o
##           Df Sum Sq Mean Sq F value   Pr(>F)    
## c_agua     2 230.82 115.408 13.2457 3.14e-05 ***
## molusco    1  23.23  23.227  2.6658   0.1097    
## Residuals 44 383.37   8.713                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Al analizar los resultados obtenidos aplicando post anova, podemos darnos cuenta que la variable moluscos, como se indicó previamente no es significativa dentro del modelo, ya que su p-value es de 3.14e-05

por otra parte, podemos evidenciar que la variable c_agua es la más representativa, ya que si la concentración de agua varia, los resultados del modelo también se verán afectados. p-value es de 0.1097

require(agricolae)
Compara= LSD.test(ModMoluscos,c("c_agua","molusco"))
Compara
## $statistics
##    MSerror Df     Mean       CV  t.value      LSD
##   8.712897 44 9.304792 31.72303 2.015368 2.974442
## 
## $parameters
##         test p.ajusted         name.t ntr alpha
##   Fisher-LSD      none c_agua:molusco   6  0.05
## 
## $means
##         cons_o      std r       LCL       UCL  Min   Max     Q25    Q50     Q75
## 100:A  9.93625 2.747976 8  7.833002 12.039498 6.78 14.00  7.9850  9.295 11.7250
## 100:B  7.40625 2.844076 8  5.303002  9.509498 3.68 11.60  5.7225  6.140 10.1000
## 50:A  12.17500 3.090178 8 10.071752 14.278248 9.74 18.80 10.3100 11.110 12.5000
## 50:B  12.32625 3.517909 8 10.223002 14.429498 6.38 17.70 10.0575 12.850 14.5000
## 75:A   7.89000 2.739578 8  5.786752  9.993248 5.20 13.20  6.0775  7.180  8.8925
## 75:B   6.09500 2.739108 8  3.991752  8.198248 1.80  9.96  4.8300  5.595  7.3425
## 
## $comparison
## NULL
## 
## $groups
##         cons_o groups
## 50:B  12.32625      a
## 50:A  12.17500      a
## 100:A  9.93625     ab
## 75:A   7.89000     bc
## 100:B  7.40625     bc
## 75:B   6.09500      c
## 
## attr(,"class")
## [1] "group"

Usamos la libreria agricolae la cual ordena los valores de p de acuerdo con las medias del grupo,

ModMoluscos
## 
## Call:
## lm(formula = cons_o ~ c_agua + molusco, data = BaseMoluscos)
## 
## Coefficients:
## (Intercept)     c_agua75    c_agua100     moluscoB  
##      12.946       -5.258       -3.579       -1.391
par(mfrow=c(2,2))
plot(ModMoluscos)

Análisis Supuestos

- Linealidad: Tiene correlacion Lineal.

- Normalidad: No se cumple el supuesto de normalidad nuestro p-value = 0.04342 por lo que es menor a 0.05

shapiro.test(ModMoluscos$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  ModMoluscos$residuals
## W = 0.9509, p-value = 0.04342

- Homocedatidad (Varianza Constante): No se evidencia Homocedaticidad por lo que se cosidera un modelo heterocedatico el cual no tiene una varianza constante y cuya P corresponde a 0.41514

library(car)
ncvTest(ModMoluscos)
## Non-constant Variance Score Test 
## Variance formula: ~ fitted.values 
## Chisquare = 0.6640401, Df = 1, p = 0.41514

- Independencia: los errores son independientes p-value = 0.2163

lmtest::dwtest(ModMoluscos)
## 
##  Durbin-Watson test
## 
## data:  ModMoluscos
## DW = 1.8991, p-value = 0.2163
## alternative hypothesis: true autocorrelation is greater than 0

Punto 3

Estudiar la relación entre ciertas características del suelo y la producción de biomasa (gr) de una planta forrajera natural se obtuvieron 45 muestras en diferentes ambientes, y en cada muestra se estimó la biomasa (respuesta Y) y se registraron las características (covariables X) del suelo en el que crecía (pH, Salinidad, Zinc y Potasio).

Análisis de correlaciones

Del análisis preliminar se pueden extraer las siguientes conclusiones:

Las variables que tienen una mayor relación lineal con la Biomasa son: PH (r= 0.928), Zinc (r= -0.781)

Tambien se oberva que PH y Zinc están fuerte mente correlacionados (r= -0,720), al igual de sanialidaad con ZINC por lo que se sugiere dejar solo Zinc para el analisis

load("C:/Mesa/Informe regresión múltiple/Salinidad.RData")
head(Salinidad)
##    Biomasa   pH Salinidad    Zinc Potasio
## 1  765.280 5.00        33 16.4524 1441.67
## 2  954.017 4.70        35 13.9852 1299.19
## 3  827.686 4.20        32 15.3276 1154.27
## 4  755.072 4.40        30 17.3128 1045.15
## 5  896.176 5.55        33 22.3312  521.62
## 6 1422.836 5.50        33 12.2778 1273.02
library(corrplot)
round(cor(x = Salinidad, method = "pearson"), 3)
##           Biomasa     pH Salinidad   Zinc Potasio
## Biomasa     1.000  0.928    -0.067 -0.781  -0.073
## pH          0.928  1.000    -0.045 -0.720   0.032
## Salinidad  -0.067 -0.045     1.000 -0.427  -0.020
## Zinc       -0.781 -0.720    -0.427  1.000   0.079
## Potasio    -0.073  0.032    -0.020  0.079   1.000
multi.hist(x = Salinidad, dcol = c("blue", "red"), dlty = c("dotted", "solid"),
           main = "")

ggpairs(Salinidad, lower = list(continuous = "smooth"),
        diag = list(continuous = "barDiag"), axisLabels = "none")

Estimación del modelo de regresión lineal múltiple

Se realiza el metodo mixto con todas las variables, como predictoras y realizando la selección con los mejores predictores con el calculo del R2 ajustado para ver el porcentaje de explicación del mismo y lo valores de p-value si son significativos tando del modelo, como las variables predictoras y finalmente se relaiza la selección de los mejores predictores con la medición Akaike(AIC).

El modelo con las cuatro variables predictorias tiene un R2 ajustado alto (0,9154), lo que significa que es capaz de explicarlo en un 91,54% de la variablidad observada de la biomasa, el p-value de las variables PH,Sanilidad y Zinc son significativos, la variable potación su p-value no arrojo significancia sim embargo en la prueba de selección AIC se mantiene como predictora del modelo

Validación del modelo

Media cero: Se cumple por defecto.

Varianza Constante: Se observa en la grafica 1 de residuales vs ajustados que el comportamiento es aleatorio.

Normalidad: Se observa en la grafica 2 que los datos se ajustan bien a la linea de normalidad en el qqplot

Independencia: Dado que estos registros no corresponden a datos en el tiempo no se tiene un orden temporal para realizar la validación de este supuesto.

Modelo: \(Biomasa= 1492.808+262.883 *pH -33.500*Salinidad-28.973*Zinc-0.115*Potasio\)

mod_sanilidad=lm(Biomasa~ .,data = Salinidad)
summary(mod_sanilidad) 
## 
## Call:
## lm(formula = Biomasa ~ ., data = Salinidad)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -293.98  -88.83   -9.48   88.20  387.27 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 1492.8076   453.6013   3.291 0.002091 ** 
## pH           262.8829    33.7304   7.794 1.51e-09 ***
## Salinidad    -33.4997     8.6525  -3.872 0.000391 ***
## Zinc         -28.9727     5.6643  -5.115 8.20e-06 ***
## Potasio       -0.1150     0.0819  -1.404 0.167979    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 158.9 on 40 degrees of freedom
## Multiple R-squared:  0.9231, Adjusted R-squared:  0.9154 
## F-statistic:   120 on 4 and 40 DF,  p-value: < 2.2e-16
step(object = mod_sanilidad, direction = "both", trace = 1)
## Start:  AIC=460.84
## Biomasa ~ pH + Salinidad + Zinc + Potasio
## 
##             Df Sum of Sq     RSS    AIC
## <none>                   1009974 460.84
## - Potasio    1     49785 1059759 461.01
## - Salinidad  1    378486 1388460 473.17
## - Zinc       1    660588 1670562 481.49
## - pH         1   1533665 2543639 500.41
## 
## Call:
## lm(formula = Biomasa ~ pH + Salinidad + Zinc + Potasio, data = Salinidad)
## 
## Coefficients:
## (Intercept)           pH    Salinidad         Zinc      Potasio  
##    1492.808      262.883      -33.500      -28.973       -0.115
par(mfrow = c(2, 2))
plot(mod_sanilidad)